Modellierung des Klangs von Gitarrennoten mit dem Karplus-Strong-Algorithmus in Python

Treffen Sie die Referenznote A der ersten Oktave (440 Hz):





Klingt schmerzhaft, nicht wahr? Was kann man noch dazu sagen, dass dieselbe Note auf verschiedenen Musikinstrumenten unterschiedlich klingt? Wieso ist es so? Es geht um das Vorhandensein zusätzlicher Harmonischer , die für jedes Instrument ein einzigartiges Timbre erzeugen.



Wir interessieren uns jedoch für eine andere Frage: Wie kann man dieses einzigartige Timbre auf einem Computer simulieren?

Standard-Karplus-Strong-Algorithmus

Illustration von dieser Seite genommen .



Das Wesentliche des Algorithmus ist wie folgt:



1) Erstellen Sie ein Array der Größe N aus Zufallszahlen (N steht in direktem Zusammenhang mit der Grundschallfrequenz).



2) Fügen Sie am Ende dieses Arrays den nach der folgenden Formel berechneten Wert hinzu:

$$

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

Wo $$$y$Ist unser Array.



3) Wir führen Punkt 2 so oft aus, wie es erforderlich ist.



Beginnen wir mit dem Schreiben des Codes:



1) Importieren Sie die erforderlichen Bibliotheken.

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

2) Wir initialisieren die Variablen.

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

3) Lärm erzeugen.

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

4) Erstellen Sie ein Array, um die Werte zu speichern und am Anfang Rauschen hinzuzufügen.

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

5) Wir verwenden die Formel.

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

6) Wir normalisieren und übersetzen in den gewünschten Datentyp.

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

7) In Datei speichern.

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

8) Lassen Sie uns alles als Funktion entwerfen. Eigentlich ist das der ganze Code.

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) Lass uns rennen und holen:





Um den Klang der Saite zu verbessern, verbessern wir die Formel leicht:

$$

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

Eine offene sechste Saite (82,41 Hz) klingt folgendermaßen:





Die offene erste Saite (329,63 Hz) klingt folgendermaßen:





Hört sich gut an, nicht wahr?



Sie können diesen Koeffizienten endlos auswählen und den Durchschnitt zwischen schönem Klang und Dauer ermitteln. Es ist jedoch besser, direkt zum Advanced Karplus-Strong-Algorithmus zu wechseln.

Ein wenig über Z-Transformation

Lassen $$$x$ Ist ein Array von Eingabewerten und $$$y$- ein Array von Ausgabewerten. Jedes Element in y wird durch die folgende Formel ausgedrückt:

$$

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

Wenn sich der Index außerhalb des Arrays befindet, ist der Wert 0. Das heißt $$$x(0-1)=0$... (Schauen Sie sich den vorherigen Code an, dort wurde er implizit verwendet).



Diese Formel kann in die entsprechende Z-Transformation geschrieben werden:

$$

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

Wenn die Formel so ist:

$$

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

Das heißt, jedes Element des Eingabearrays hängt vom vorherigen Element desselben Arrays ab (außer natürlich dem Nullelement). Dann sieht die entsprechende Z-Transformation folgendermaßen aus:

$$

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

Umgekehrter Prozess: Holen Sie sich die Formel für jedes Element aus der Z-Transformation. Zum Beispiel,

$$

$$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).$$

Wenn jemand nicht verstanden hat, lautet die Formel: $$$Y(z)*α*z^{-k}=α*y(n-k)$wo $$$α$- eine beliebige reelle Zahl.



Wenn Sie zwei Z-Transformationen miteinander multiplizieren müssen, dann$$$z^{-a}*z^{-b}=z^{-a-b}.$

Erweiterter Karplus-Strong-Algorithmus

Illustration von dieser Seite genommen.



Hier finden Sie eine kurze Beschreibung der einzelnen Funktionen.

Teil I. Funktionen, die das Anfangsrauschen transformieren

1) Tiefpassfilter in Auswahlrichtung (Tiefpassfilter)$$$H_p (z)$...

$$

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

Entsprechende Formel:

$$

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

Der Code:

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

Sie sollten immer ein anderes Array erstellen, um Fehler zu vermeiden. Vielleicht hätte es hier nicht verwendet werden können, aber im nächsten Filter können Sie nicht darauf verzichten.



2) Pick-Position-Kammfilter (Kammfilter)$$$H_β (z)$...

$$

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

Entsprechende Formel:

$$

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

Der Code:

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

Im ersten Absatz auf Seite 13 dieses Dokuments wird Folgendes geschrieben (nicht wörtlich, aber unter Beibehaltung der Bedeutung): Der Koeffizient β ahmt die Position der gezupften Zeichenfolge nach. Wenn$$$β=1/2$, dann bedeutet dies, dass der Zupf in der Mitte der Saite gemacht wurde. Wenn$$$β=1/10$ - Der Zupf wurde an einem Zehntel der Schnur von der Brücke gemacht.

Teil II. Funktionen, die sich auf den Hauptteil des Algorithmus beziehen

Hier gibt es eine Falle, die wir umgehen müssen. Zum Beispiel String-Dampling-Filter$$$H_d (z)$ so geschrieben: $$$H_d (z)=(1-S)+Sz^{-1}$... Aber die Abbildung zeigt, dass er den Sinn dort nimmt, wo er ihn gibt. Das heißt, es stellt sich heraus, dass die Eingangs- und Ausgangssignale für diesen Filter ein und dasselbe sind. Dies bedeutet, dass nicht jeder Filter separat angewendet werden kann, da im vorherigen Abschnitt alle Filter gleichzeitig angewendet werden müssen. Dies kann beispielsweise durch Auffinden des Produkts jedes Filters erfolgen. Dieser Ansatz ist jedoch nicht rational: Wenn Sie einen Filter hinzufügen oder ändern, müssen Sie alles erneut multiplizieren. Es ist möglich, dies zu tun, aber es macht keinen Sinn. Ich möchte den Filter mit einem Klick ändern und nicht alles immer wieder multiplizieren.

Da das Ausgangssignal des Filters als Eingang für ein anderes Filter betrachtet wird, schlage ich vor, jedes Filter als separate Funktion zu schreiben, die die Funktion des vorherigen Filters in sich selbst aufruft.

Ich denke, der Beispielcode wird deutlich machen, was ich meine.

1) Verzögerungsleitungsfilter $$$z^{-N}.$

$$

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

Entsprechende Formel:

$$

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

Der Code:

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

2) String-Dämpfungsfilter $$$H_d (z)$...

$$

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

Im ursprünglichen Algorithmus $$$S=0.5.$

Entsprechende Formel:

$$

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

Der Code:

# 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))

In diesem Fall ist dieser Filter der One Zero String-Dämpfungsfilter. Es gibt andere Optionen, Sie über sie lesen hier .



3) Allpassfilter mit Saitensteifigkeit$$$H_s (z)$...

Egal wie sehr ich aussah, leider konnte ich nichts Bestimmtes finden. Hier ist der Filter allgemein geschrieben. Aber das funktioniert nicht, da es am schwierigsten ist, die richtigen Chancen zu finden. In diesem Dokument auf Seite 14 befindet sich noch etwas anderes , aber ich habe nicht genügend mathematischen Hintergrund, um zu verstehen, was dort passiert und wie man es verwendet. Wenn jemand kann, lass es mich wissen.



4) String-Tuning-Allpassfilter erster Ordnung$$$H_ρ (z)$...

Seite 6, unten links in diesem Dokument:

$$

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

Entsprechende Formel:

$$

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

Der Code:

# 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)

Es ist zu beachten, dass Sie, wenn Sie nach diesem Filter weitere Filter hinzufügen, den vorherigen Wert speichern müssen, da er nicht mehr im Beispielarray gespeichert wird.

Da die Länge des Anfangsrauschens eine ganze Zahl ist, werfen wir den Bruchteil beim Zählen weg. Dies führt zu Fehlern und Ungenauigkeiten. Wenn beispielsweise die Abtastrate 44100 und die Rauschlänge 133 und 134 beträgt, betragen die entsprechenden Signalfrequenzen 331,57 Hz und 329,10 Hz. Die Frequenz der E-Noten der ersten Oktave (der ersten offenen Saite) beträgt 329,63 Hz. Hier liegt der Unterschied in Zehnteln, aber zum Beispiel für den 15. Bund kann der Unterschied bereits mehrere Hz betragen. Um diesen Fehler zu reduzieren, ist dieser Filter vorhanden. Sie kann weggelassen werden, wenn die Abtastfrequenz hoch ist (wirklich hoch: mehrere hunderttausend Hz oder sogar mehr) oder die Grundfrequenz niedrig ist, wie zum Beispiel für Basssaiten.

Es gibt noch andere Variationen, können Sie über sie alle lesen gibt .



5) Wir nutzen unsere Funktionen.

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

Teil III. Dynamischer Pegel-Tiefpassfilter $$$H_L (z).$

$$$\check{ω}=ω\frac{T}{2}=2πf\frac{T}{2}=π\frac{f}{F_s}$wo $$$f$ - fundamentale Frequenz, $$$F_s$- Abtastfrequenz.

Zuerst finden wir das Array$$$y$ mit folgender Formel:

$$

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

Entsprechende Formel:

$$

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

Dann wenden wir die folgende Formel an:

$$

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

Der Code:

# 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

Der Parameter L beeinflusst den Wert für die Lautstärkeverringerung. Mit Werten von 0,001, 0,01, 0,1, 0,32 nimmt das Signalvolumen um 60, 40, 20 bzw. 10 dB ab.



Lassen Sie uns alles als Funktion entwerfen. Eigentlich ist das der ganze Code.

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)

Eine offene sechste Saite (82,41 Hz) klingt folgendermaßen:





Und die offene erste Saite (329,63 Hz) klingt so:





Die erste Saite klingt, gelinde gesagt, nicht sehr gut. Eher wie eine Glocke als eine Schnur. Ich habe sehr lange versucht herauszufinden, was im Algorithmus falsch war. Ich dachte, es wäre ein unbenutzter Filter. Nach Tagen des Experimentierens wurde mir klar, dass ich die Abtastrate auf mindestens 100.000 erhöhen musste:





Klingt besser, nicht wahr?



Add-Ons wie das Spielen von Glissando oder das Simulieren einer sympathischen Zeichenfolge können in diesem Dokument gelesen werden (S. 11-12).



Hier ist ein Kampf für dich:





Akkordfolge: CG # Am F. Strike: Sechs. Die Verzögerung zwischen zwei aufeinanderfolgenden Zupfen der Saite beträgt 0,015 Sekunden. Die Verzögerung zwischen zwei aufeinanderfolgenden Treffern in einem Kampf beträgt 0,205 Sekunden. Die Verzögerung selbst im Kampf beträgt 0,41 Sekunden. Der Algorithmus hat den Wert von L auf 0,2 geändert.



Vielen Dank für das Lesen des Artikels. Viel Glück!