Detekcja częstotliwoci podstawowej przy pomocy funkcji AMDF

Š 1999, 2000, Marek Szczerba

Inicjalizacja

Wczytanie bibliotek

Biblioteka obsługi dwięku

<< Miscellaneous`Audio`

Biblioteki statystyczne

<< Statistics`DescriptiveStatistics`

Akwizycja sygnału

SetDirectory["D:\\MIDI\\"];

snd = ReadSoundfile["mobs.wav", PrintHeader -> True];

Parametry definiowane przez użytkownika

Odstęp w próbkach pomiędzy poszczególnymi punktami analizy

Skok=200;

Częstotliwoć próbkowania

CzP=22050;

Granice pasma analizy

CzMin=100;
CzMax=2000;

Wartoć parametru AMDF

K=0.9547;

Wartoć parametru wzmocnienia niskich częstotliwoci

Tłumik=0.9547;

Parametry obliczane automatycznie

Konwersja pasma analizy

\(OkrMin = Ceiling[CzP\/CzMax]; 
OkrMax = Floor[CzP\/CzMin]; \)

Liczba punktów analizy

\(LRamek = 0; 

    For[i = 1, i <= Length[snd] - 2  OkrMax, i += Skok, \(LRamek++\)]; \)

Definicje funkcji przetwarzania sygnału

Definicja funkcji AMDF

\(\(AMDF = 
	
      Function[{dane, K}, 
		fL = Length[dane]\  - \ OkrMax; 
		
        Wynik = Table[0, \ {fi, \ OkrMax}]; 
		
        For[fi = OkrMin, \ fi\  <= \ OkrMax, \ \(fi++\), 
			
          ρ = 
            ∑\+\(fn = 1\)\%fL Abs[
                dane[[fn]] - 
                  dane[[fn + fi]]]; 

          			Wynik[[fi]] = ρ; ]; 
        
\ \ \ \ \ \ \ \ Wynik]; \)\)

Pętla funkcji AMDF

\(\(OblAMDF = 
	
      Function[snd, 
		
        Wyniki = Table[Table[0, \ {fx, \ OkrMax}], \ {fy, \ LRamek}]; 
		
        Numer = 0; 
		
        For[fk = 1, fk <= Length[snd] - 2\ OkrMax, fk += Skok, 
			
          \(Numer++\); 
			
          Ramka = Table[
              snd[[fj]], {fj, fk, 
                fk + 2\ OkrMax - 1}]; 
			AKRamka\  = \ AMDF[Ramka, K]; \ 
          
			For[fm\  = \ OkrMin, \ fm\  <= \ OkrMax, \ \(fm++\), 

            				\(Wyniki[[Numer, fm
                ]] = 
              AKRamka[[fm]]; \)]; ]; 

        		Wyniki]; \)\)

Poszukiwanie maksimum funkcji autokorelacji

\(\(MinAMDF = 
	
      Function[dane, 
		Wynik = {}; 
		
        For[fi = 1, fi <= Length[dane], \(fi++\), \ 
			
          Minim = 100000000; 
			Lok = 0; 
			
          For[fj = OkrMin + 1, 
            fj <= Length[dane[[fi]]], 
            \(fj++\), 
				
            \(If[dane[[fi, fj]] < Minim, 
              
					
              Minim = dane[[fi, fj]]; 

              					Lok = fj; ]; \)]; 
			
          AppendTo[Wynik, {Lok, Minim}]; ]; 
		Wynik]; \)\)

Funkcja tłumienia składowych o niskiej częstotliwoci

\(\(Tłumienie = 
      Function[{dane, k}, 
		Wynik = dane; 
		M = OkrMax - OkrMin; 

        		For[fn = 1, fn <= Length[dane], \(fn++\), 
			
          For[fm = 1, 
            fm <= Length[dane[[fn]]], 
            \(fm++\), 
				
            \(Wynik[[fn, fm]] = 
              dane[[fn, fm]]*
                \((1 - \(k \((M - \ fm)\)\)\/M)\); \)]]; 
		Wynik]; \)\)

Funkcja detekcji częstotliwoci podstawowej

\(\(DetF0 = 
      Function[dane, 
		Wy = Table[0, {fx, Length[dane]}]; 
		
        For[fi = 1, fi <= Length[dane], \(fi++\), 
			
          Wy[[fi]] = 
            CzP\/dane[[fi, 1]]]; 
		
        Wy]; \)\)

Obliczenie poprawnoci detekcji

\(\(BłšdProc = 
      Function[dane, 
		OK = 0; 
		fj = 1; 
		
        Błędzik = {}; 
		Liczba = Length[dane]; 
		
        For[fi = 1, fi <= Liczba, \(fi++\), 
			
          While[\(\((fi - 1)\) Skok + 1\)\/CzP > 
                Fraza[[fj, 2]] && 
              fj < Length[Fraza], \(fj++\)]; 
			
          GrDolna = 
            Fraza[[fj, 1]] 
              \@\(2\^\(-1\)\)\%12; 
			
          GrGórna = 
            Fraza[[fj, 1]] \@2\%12; 

          			AppendTo[Błędzik, 
            Fraza[[fj, 1]] - 
              dane[[fi]]]; 
			
          If[GrDolna < dane[[fi]] && 
              dane[[fi]] < 
                GrGórna, \(OK++\)]; 
			
          If[dane[[fi]] == 0, 
            \(Liczba--\)]]; 
		100  OK\/Liczba // N]; \)\)

Obliczenia

Obliczenie funkcji autokorelacji dla poszczególnych ramek sygnału

Timing[AM = OblAMDF[snd]; ]

Tłumienie składowych o niższej częstotliwoci

TłAMDF = Tłumienie[AM, Tłumik];

Detekcja minimów funkcji AMDF

MAMDF = MinAMDF[TłAMDF];

Detekcja częstotliwoci podstawowej

FS = DetF0[MAMDF];

InS = 
      Table[{(((i - 1) Skok + 1)/CzP), 
          FS[[i]]}, {i, Length[FS]}];

a = ListPlot[InS, PlotStyle -> RGBColor[0, 0, 1]];

Show[a, Graph,  Axes -> False,  Frame -> True];

Obliczenie błędu detekcji

Błšd = BłšdProc[FS]

Wykres błędu

\(\(ListPlot[Błędzik, PlotStyle -> RGBColor[0, 0, 1], 
      PlotJoined -> True, PlotRange -> All]; \)\)

Optymalizacja wartoci parametru tłumienia

\(\(P = Tłumik; 
R = 2; 
Niegotowe = True; 

    While[Niegotowe, 
	Parametry = {P - 1\/R, P, P + 1\/R}; 
	
      Dobre = {0, 0, 0}; 
	Szumy = {0, 0, 0}; 
	
      For[i = 1, i <= Length[Parametry], \(i++\), 
		
        TłAMDF = 
          Tłumienie[AM, 
            Parametry[[i]]]; 
		
        MAMDF = MinAMDF[TłAMDF]; 
		FS = DetF0[MAMDF]; 
		
        OT = BłšdProc[FS]; 
		
        Dobre[[i]] = OT; 
	]; 
	
      Print["P=", N[P], "   R=", N[R], "   Dobre=", Dobre]; 
	
      If[Max[Dobre] == Dobre[[2]], 
        R += \(1\/2\) R, 
		
        If[Dobre[[1]] > 
            Dobre[[3]], 
          P = Parametry[[1]], 
          P = Parametry[[3]]]]; 
	
      If[R > 200, Niegotowe = False]; ]; 			\)\)