Pemodelan suara not gitar menggunakan algoritma Karplus-Strong dengan python

Memenuhi catatan referensi A dari oktaf pertama (440 Hz):





Kedengarannya menyakitkan, bukan? Apa lagi yang harus dikatakan tentang fakta bahwa not yang sama berbunyi berbeda pada alat musik yang berbeda. Kenapa gitu? Ini semua tentang adanya harmonik tambahan yang menciptakan timbre unik untuk setiap instrumen.



Tetapi kami tertarik dengan pertanyaan lain: bagaimana cara mensimulasikan timbre unik ini di komputer?

Algoritme Karplus-Kuat standar

Ilustrasi diambil dari situs ini .



Inti dari algoritma ini adalah sebagai berikut:



1) Buat array berukuran N dari bilangan acak (N berhubungan langsung dengan frekuensi suara dasar).



2) Tambahkan ke akhir larik ini nilai yang dihitung dengan rumus berikut:

$$

$$y(n)=\frac{y(n-N)+y(n-N-1)}{2},$$

Dimana $$$y$Apakah array kami.



3) Kami melaksanakan poin 2 beberapa kali.



Mari mulai menulis kode:



1) Impor perpustakaan yang diperlukan.

import numpy as np
import scipy.io.wavfile as wave

2) Kami menginisialisasi variabel.

frequency = 82.41     #     
duration = 1          #    
sample_rate = 44100   #  

3) Ciptakan kebisingan.

#  ,  frequency, ,        frequency .
#      sample_rate/length .
#  length = sample_rate/frequency.
noise = np.random.uniform(-1, 1, int(sample_rate/frequency))   

4) Buat array untuk menyimpan nilai dan menambahkan noise di awal.

samples = np.zeros(int(sample_rate*duration))
for i in range(len(noise)):
    samples[i] = noise[i]

5) Kami menggunakan rumus.

for i in range(len(noise), len(samples)):
    #   i   ,      .
    #  ,  i   ,       .
    samples[i] = (samples[i-len(noise)]+samples[i-len(noise)-1])/2

6) Kami menormalkan dan menerjemahkan ke dalam tipe data yang diinginkan.

samples = samples / np.max(np.abs(samples))  
samples = np.int16(samples * 32767)     

7) Simpan ke file.

wave.write("SoundGuitarString.wav", 44100, samples)

8) Mari mendesain semuanya sebagai fungsi. Sebenarnya, itu semua kodenya.

import numpy as np
import scipy.io.wavfile as wave
 
def GuitarString(frequency, duration=1., sample_rate=44100, toType=False):
    #  ,  frequency, ,        frequency .
    #      sample_rate/length .
    #  length = sample_rate/frequency.
    noise = np.random.uniform(-1, 1, int(sample_rate/frequency))      #  
 
    samples = np.zeros(int(sample_rate*duration))
    for i in range(len(noise)):
        samples[i] = noise[i]
    for i in range(len(noise), len(samples)):
        #   i   ,      .
        #  ,  i   ,       .
        samples[i] = (samples[i-len(noise)]+samples[i-len(noise)-1])/2
 
    if toType:
        samples = samples / np.max(np.abs(samples))  #   -1  1
        return np.int16(samples * 32767)             #     int16
    else:
        return samples
 
 
frequency = 82.41
sound = GuitarString(frequency, duration=4, toType=True)
wave.write("SoundGuitarString.wav", 44100, sound)

9) Ayo lari dan dapatkan:





Untuk membuat senar terdengar lebih baik, mari kita sedikit meningkatkan rumus:

$$

$$y(n)=0.996\frac{y(n-N)+y(n-N-1)}{2}$$

Senar keenam terbuka (82,41 Hz) berbunyi seperti ini:





String pertama yang terbuka (329,63 Hz) berbunyi seperti ini:





Kedengarannya bagus, bukan?



Anda dapat tanpa henti memilih koefisien ini dan menemukan rata-rata antara suara dan durasi yang indah, tetapi lebih baik langsung ke algoritma Advanced Karplus-Strong.

Sedikit tentang Z-transform

Biarlah $$$x$ Merupakan larik nilai input, dan $$$y$- array nilai keluaran. Setiap elemen dalam y diekspresikan dengan rumus berikut:

$$

$$y(n)=x(n)+x(n-1).$$

Jika indeks berada di luar larik, maka nilainya adalah 0. Yaitu $$$x(0-1)=0$... (Lihat kode sebelumnya, di sana secara implisit digunakan).



Rumus ini dapat ditulis dalam Z-transform yang sesuai:

$$

$$H(z)=1+z^{-1}.$$

Jika rumusnya seperti ini:

$$

$$y(n)=x(n)+x(n-1)-y(n-1).$$

Artinya, setiap elemen dari array input bergantung pada elemen sebelumnya dari array yang sama (kecuali elemen nol, tentu saja). Kemudian Z-transform yang sesuai terlihat seperti ini:

$$

$$H(z)=\frac{1+z^{-1}}{1+z^{-1}}.$$

Proses sebaliknya: dapatkan rumus untuk setiap elemen dari Z-transform. Misalnya,

$$

$$H(z)=\frac{1+z^{-1}}{1-z^{-1}}.$$

$$

$$H(z)=\frac{Y(z)}{X(z)} =\frac{1+z^{-1}}{1-z^{-1}}.$$

$$

$$Y(z)*(1-z^{-1} )=X(z)*(1+z^{-1} ).$$

$$

$$Y(z)*1-Y(z)*z^{-1}= X(z)*1+X(z)*z^{-1}.$$

$$

$$y(n)-y(n-1)=x(n)+x(n-1).$$

$$

$$y(n)=x(n)+x(n-1)+y(n-1).$$

Jika seseorang tidak mengerti, maka rumusnya adalah: $$$Y(z)*α*z^{-k}=α*y(n-k)$dimana $$$α$- bilangan riil apapun.



Jika Anda perlu mengalikan dua transformasi Z satu sama lain, maka$$$z^{-a}*z^{-b}=z^{-a-b}.$

Algoritma Karplus-Kuat yang diperluas

Ilustrasi diambil dari situs ini .



Berikut penjelasan singkat dari setiap fitur.

Bagian I. Fungsi yang mengubah kebisingan awal

1) Filter lowpass arah pilih (low pass filter)$$$H_p (z)$...

$$

$$H_p (z)=\frac{1-p}{1-pz^{-1}},p ∈ [0,1).$$

Rumus yang sesuai:

$$

$$y(n)=(1-p)x(n)+py(n-1).$$

Kode:

buffer = np.zeros_like(noise)
buffer[0] = (1 - p) * noise[0]
for i in range(1, N):
    buffer[i] = (1-p)*noise[i] + p*buffer[i-1]
noise = buffer

Anda harus selalu membuat larik lain untuk menghindari kesalahan. Mungkin tidak dapat digunakan di sini, tetapi di filter berikutnya Anda tidak dapat melakukannya tanpanya.



2) Filter sisir pilih posisi (filter sisir)$$$H_β (z)$...

$$

$$H_β (z)=1-z^{-int(β*N+1/2)},β∈(0,1).$$

Rumus yang sesuai:

$$

$$y(n)=x(n)-x(n-int(β*N+1/2)).$$

Kode:

pick = int(beta*N+1/2)
if pick == 0:
    pick = N   #      
buffer = np.zeros_like(noise)
for i in range(N):
    if i-pick < 0:
        buffer[i] = noise[i]
    else:
        buffer[i] = noise[i]-noise[i-pick]
noise = buffer

Pada paragraf pertama di halaman 13 dokumen ini , tertulis berikut ini (tidak secara harfiah, tetapi dengan menjaga artinya): koefisien β meniru posisi string yang dipetik. Jika$$$β=1/2$, maka ini berarti pemetikan dilakukan di tengah-tengah benang. Jika$$$β=1/10$ - pemetikan dilakukan pada sepersepuluh dari tali jembatan.

Bagian II. Fungsi terkait dengan bagian utama algoritme

Ada jebakan di sini yang harus kita atasi. Misalnya, filter peredam tali$$$H_d (z)$ ditulis seperti ini: $$$H_d (z)=(1-S)+Sz^{-1}$... Tetapi gambar tersebut menunjukkan bahwa dia mengambil makna dari mana dia memberikannya. Artinya, ternyata sinyal masukan dan keluaran untuk filter ini adalah satu dan sama. Artinya, setiap filter tidak dapat diterapkan secara terpisah, seperti di bagian sebelumnya, semua filter harus diterapkan secara bersamaan. Ini dapat dilakukan, misalnya dengan mencari produk dari setiap filter. Tetapi pendekatan ini tidak rasional: saat menambah atau mengubah filter, Anda harus menggandakan semuanya lagi. Itu mungkin untuk dilakukan, tetapi itu tidak masuk akal. Saya ingin mengubah filter dalam satu klik, dan tidak menggandakan semuanya lagi dan lagi.

Karena sinyal keluaran dari filter dianggap sebagai masukan untuk filter lain, saya mengusulkan untuk menulis setiap filter sebagai fungsi terpisah yang memanggil fungsi filter sebelumnya di dalam dirinya sendiri.

Saya pikir kode contoh akan menjelaskan apa yang saya maksud.

1) Filter Garis Tunda $$$z^{-N}.$

$$

$$H(z)=z^{-N}.$$

Rumus yang sesuai:

$$

$$y(n)=x(n-N).$$

Kode:

#    ,     samples  0.
#    n-N<0   0,    .
def DelayLine(n):
    return samples[n-N]

2) Filter peredam tali $$$H_d (z)$...

$$

$$H_d (z)=(1-S)+Sz^{-1},S∈[0,1].$$

Dalam algoritma aslinya $$$S=0.5.$

Rumus yang sesuai:

$$

$$y(n)=(1-S)x(n)+Sx(n-1).$$

Kode:

# String-dampling filter.
# H(z) = 0.996*((1-S)+S*z^(-1)).    S = 0.5. S ∈ [0, 1]
# y(n)=0.996*((1-S)*x(n)+S*x(n-1))
def StringDampling_filter(n):
    return 0.996*((1-S)*DelayLine(n)+S*DelayLine(n-1))

Dalam hal ini, filter ini adalah filter peredam String Satu Nol. Ada opsi lain, Anda dapat membacanya di sini .



3) Filter allpass kekakuan tali$$$H_s (z)$...

Tidak peduli seberapa banyak saya melihat, sayangnya, saya tidak dapat menemukan sesuatu yang spesifik. Di sini filter ditulis secara umum. Tetapi itu tidak berhasil karena bagian tersulit adalah menemukan peluang yang tepat. Ada hal lain dalam dokumen ini di halaman 14, tetapi saya tidak memiliki dasar matematis yang cukup untuk memahami apa yang terjadi di sana dan bagaimana menggunakannya. Jika ada yang bisa, beri tahu saya.



4) Filter allpass penyetelan string urutan pertama$$$H_ρ (z)$...

Halaman 6, kiri bawah dokumen ini :

$$

$$H_ρ (z)=\frac{C+z^{-1}}{1+Cz^{-1}},C∈(-1,1).$$

Rumus yang sesuai:

$$

$$y(n)=Cx(n)+x(n-1)-Cy(n-1).$$

Kode:

# First-order string-tuning allpass filter
# H(z) = (C+z^(-1))/(1+C*z^(-1)). C ∈ (-1, 1)
# y(n) = C*x(n)+x(n-1)-C*y(n-1)
def FirstOrder_stringTuning_allpass_filter(n):
    #        ,    ,  
    #    ,          samples.
    return C*(StringDampling_filter(n)-samples[n-1])+StringDampling_filter(n-1)

Harus diingat bahwa jika Anda menambahkan lebih banyak filter setelah filter ini, maka Anda harus menyimpan nilai lampau, karena nilai itu tidak akan lagi disimpan dalam larik sampel.

Karena panjang derau awal adalah bilangan bulat, kita membuang bagian pecahan saat menghitung. Ini menyebabkan kesalahan dan ketidakakuratan. Sebagai contoh, jika sample rate 44100 dan panjang noise 133 dan 134, maka frekuensi sinyal yang sesuai adalah 331.57 Hz dan 329.10 Hz. Dan frekuensi not E oktaf pertama (string terbuka pertama) adalah 329,63 Hz. Di sini perbedaannya ada pada sepersepuluh, tetapi, misalnya, untuk fret ke-15, perbedaannya mungkin sudah beberapa Hz. Filter ini ada untuk mengurangi kesalahan ini. Ini tidak dapat digunakan jika frekuensi sampling tinggi (sangat tinggi: beberapa ratus ribu Hz, atau bahkan lebih) atau frekuensi fundamental rendah, seperti, misalnya, untuk string bass.

Ada variasi lain, Anda bisa membaca semuanya di sana .



5) Kami menggunakan fungsi kami.

def Modeling(n):
    return FirstOrder_stringTuning_allpass_filter(n)
 
for i in range(N, len(samples)):
    samples[i] = Modeling(i)

Bagian III. Filter Lowpass Tingkat Dinamis $$$H_L (z).$

$$$\check{ω}=ω\frac{T}{2}=2πf\frac{T}{2}=π\frac{f}{F_s}$dimana $$$f$ - frekuensi dasar, $$$F_s$- frekuensi sampling.

Pertama kita temukan arraynya$$$y$ dengan rumus sebagai berikut:

$$

$$H(z)=\frac{\check{ω}}{1+\check{ω}}\frac{1+z^{-1}}{1-\frac{1-\check{ω}}{1+\check{ω}} z^{-1}}$$

Rumus yang sesuai:

$$

$$y(n)=\frac{\check{ω}}{1+\check{ω}}(x(n)+x(n-1))+\frac{1-\check{ω}}{1+\check{ω}}y(n-1)$$

Kemudian kami menerapkan rumus berikut:

$$

$$x(n)=L^{\frac{4}{3}}x(n)+(1-L)y(n),L∈(0,1)$$

Kode:

# Dynamic-level lowpass filter. L ∈ (0, 1/3)
w_tilde = np.pi*frequency/sample_rate
buffer = np.zeros_like(samples)
buffer[0] = w_tilde/(1+w_tilde)*samples[0]
for i in range(1, len(samples)):
    buffer[i] = w_tilde/(1+w_tilde)*(samples[i]+samples[i-1])+(1-w_tilde)/(1+w_tilde)*buffer[i-1]
samples = (L**(4/3)*samples)+(1.0-L)*buffer

Parameter L mempengaruhi nilai penurunan volume. Dengan nilainya sama dengan 0,001, 0,01, 0,1, 0,32, volume sinyal berkurang masing-masing sebesar 60, 40, 20 dan 10 dB.



Mari kita rancang semuanya sebagai sebuah fungsi. Sebenarnya, itu semua kodenya.

import numpy as np
import scipy.io.wavfile as wave
 
 
def GuitarString(frequency, duration=1., sample_rate=44100, p=0.9, beta=0.1, S=0.5, C=0.1, L=0.1, toType=False):
    N = int(sample_rate/frequency)            #    
 
    noise = np.random.uniform(-1, 1, N)   #  
 
    # Pick-direction lowpass filter (  ).
    # H(z) = (1-p)/(1-p*z^(-1)). p ∈ [0, 1)
    # y(n) = (1-p)*x(n)+p*y(n-1)
    buffer = np.zeros_like(noise)
    buffer[0] = (1 - p) * noise[0]
    for i in range(1, N):
        buffer[i] = (1-p)*noise[i] + p*buffer[i-1]
    noise = buffer
 
    # Pick-position comb filter ( ).
    # H(z) = 1-z^(-int(beta*N+1/2)). beta ∈ (0, 1)
    # y(n) = x(n)-x(n-int(beta*N+1/2))
    pick = int(beta*N+1/2)
    if pick == 0:
        pick = N   #      
    buffer = np.zeros_like(noise)
    for i in range(N):
        if i-pick < 0:
            buffer[i] = noise[i]
        else:
            buffer[i] = noise[i]-noise[i-pick]
    noise = buffer
 
    #    .
    samples = np.zeros(int(sample_rate*duration))
    for i in range(N):
        samples[i] = noise[i]
 
    #    ,     samples  0.
    #    n-N<0   0,    .
    def DelayLine(n):
        return samples[n-N]
 
    # String-dampling filter.
    # H(z) = 0.996*((1-S)+S*z^(-1)).    S = 0.5. S ∈ [0, 1]
    # y(n)=0.996*((1-S)*x(n)+S*x(n-1))
    def StringDampling_filter(n):
        return 0.996*((1-S)*DelayLine(n)+S*DelayLine(n-1))
 
    # First-order string-tuning allpass filter
    # H(z) = (C+z^(-1))/(1+C*z^(-1)). C ∈ (-1, 1)
    # y(n) = C*x(n)+x(n-1)-C*y(n-1)
    def FirstOrder_stringTuning_allpass_filter(n):
        #        ,    ,  
        #    ,          samples.
        return C*(StringDampling_filter(n)-samples[n-1])+StringDampling_filter(n-1)
 
    def Modeling(n):
        return FirstOrder_stringTuning_allpass_filter(n)
 
    for i in range(N, len(samples)):
        samples[i] = Modeling(i)
 
    # Dynamic-level lowpass filter. L ∈ (0, 1/3)
    w_tilde = np.pi*frequency/sample_rate
    buffer = np.zeros_like(samples)
    buffer[0] = w_tilde/(1+w_tilde)*samples[0]
    for i in range(1, len(samples)):
        buffer[i] = w_tilde/(1+w_tilde)*(samples[i]+samples[i-1])+(1-w_tilde)/(1+w_tilde)*buffer[i-1]
    samples = (L**(4/3)*samples)+(1.0-L)*buffer
 
    if toType:
        samples = samples/np.max(np.abs(samples))   #   -1  1
        return np.int16(samples*32767)              #     int16
    else:
        return samples
 
 
frequency = 82.51
sound = GuitarString(frequency, duration=4, toType=True)
wave.write("SoundGuitarString.wav", 44100, sound)

Senar keenam terbuka (82,41 Hz) berbunyi seperti ini:





Dan string pertama yang terbuka (329,63 Hz) berbunyi seperti ini:





Senar pertama tidak terdengar sangat bagus, secara halus. Lebih seperti lonceng daripada tali. Untuk waktu yang sangat lama saya mencoba mencari tahu apa yang salah dalam algoritme. Pikir itu filter yang tidak digunakan. Setelah beberapa hari bereksperimen, saya menyadari bahwa saya perlu meningkatkan sample rate menjadi setidaknya 100.000:





Kedengarannya lebih baik, bukan?



Pengaya seperti memainkan glissando atau simulasi string simpatik dapat dibaca dalam dokumen ini (hlm. 11-12).



Inilah pertarungan untuk Anda:





Perkembangan akor: CG # Am F. Strike: Enam. Penundaan antara dua pemetikan string yang berurutan adalah 0,015 detik; penundaan antara dua serangan berturut-turut dalam pertempuran adalah 0,205 detik; penundaan itu sendiri dalam pertempuran adalah 0,41 detik. Algoritme telah mengubah nilai L menjadi 0,2.



Terima kasih telah membaca artikelnya. Semoga berhasil!