Metoda cepstralna

Š Marek Szczerba, 2000

Inicjalizacja

Wczytanie bibliotek

\(<< Miscellaneous`Audio`
<< CCRMA`SCMTheory`
<< CCRMA`Windows`

      \(<< Statistics`MultiDescriptiveStatistics`; \)\)

Akwizycja sygnału

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

    snd = ReadSoundfile["moboe2.wav"]; \)

Parametry okrelane przez użytkownika

Częstotliwoć próbkowania

CzP=22050;

Rozmiar okna czasowego

OknoN=512;

Zakładka

Zakładka=(1/2);

Współczynnik zeropaddingu

wZP=3;

Zakres częstotliwoci

CzMin=100;
CzMaks=1500;

Próg szumowy

PrSz=-∞;

Detekcja maksymalnej wartoci bezwzględnej próbki w sygnale

MaksPróbka = Max[Abs[snd]];

Zakres częstotliwoci podstawowej (w próbkach)

\(FMin = Floor[CzMin \((\((1 + wZP)\)\ OknoN)\)/CzP]; 

    FMax = Ceiling[CzMaks \((\((1 + wZP)\)\ OknoN)\)/CzP]; \)

Definicje funkcji przetwarzania sygnału

Definicja okna czasowego

\(\(FOkna = Function[dane, 
		dane*HammingW[Length[dane]]]; \)\)

Zeropadding

\(\(ZP = 
      Function[{dane, wZP}, 
		Wy = dane; 
		
        zera = Table[0, {x, wZP\ Length[dane]}]; 
		AppendTo[Wy, zera]; 

        		Flatten[Wy]]; \)\)

Obliczenie widm

\(\(OblWidma = 
      Function[dane, 
		TempWy = Abs[Fourier[dane]]; 
		
        Wy = Table[
            TempWy[[fx]], {fx, 
              Length[TempWy]\/2}]; 
		
        For[fi = 1, fi <= FMin, \(fi++\), 
          \(Wy[[fi]] = 0; \)]; 
		Wy

        	]; \)\)

Tłumienie cepstrum

\(\(Tłumienie = 
      Function[{dane, par}, 
		Wy = dane; 
		
        Tłumik = 
          Table[Log[fx]\^par, {fx, 1, 
              Length[dane[[1]]]}]; 
		
        For[fi = 1, fi <= Length[dane], \(fi++\), 
			
          For[fj = 1, 
            fj <= Length[dane[[fi]]], 
            \(fj++\), 
				
            Wy[[fi]] = 
              dane[[fi]] 
                Tłumik]]; 
		Wy]; \)\)

Definicja funkcji błędu

\(\(BłšdProc = 
      Function[dane, 
		OK = 0; 
		fj = 1; 
		
        Błędzik = {}; 
		Liczba = Length[dane]; 
		
        For[fi = 1, fi <= Liczba, \(fi++\), 
			
          While[\(\((fi - 1)\) OknoN \((1 - Zakładka)\) + 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 logarytmu widma

\(i = 0; 
Widma = {}; 

    While[i <= Length[snd] - OknoN, 
	
      Ramka = Table[
          snd[[x]], {x, i + 1, 
            i + OknoN}]; 
	
      If[20  Log[10, Max[Abs[Ramka]]\/MaksPróbka] >= PrSz, 
		
        Okno = FOkna[Ramka]; 
		ZeroP = ZP[Okno, wZP]; 
		
        AppendTo[Widma, N[Log[Abs[OblWidma[ZeroP]]]]]]; 

      i += Floor[OknoN \((1 - Zakładka)\)]; ]; \)

"Lustrane odbicie" logarytmu widma

\(PreCeps = {}; 

    Lram = Length[Widma[[1]]]; 

    For[i = 1, i <= Length[Widma], \(i++\), 
	
      For[j = 1, j <= Lram, \(j++\), 
        If[Widma[[i, j]] < 0, 
          Widma[[i, j]] = 0]]; 
	
      PreCFR = Table[0, {x, 2  Lram}]; 
	
      For[j = 1, j <= Lram, \(j++\), 
        PreCFR[[j]] = 
          Widma[[i, Lram - j + 1]]]; 
	
      For[j = Lram + 1, j <= 2  Lram, \(j++\), 
        PreCFR[[j]] = 
          Widma[[i, j - Lram]]]; 
	
      AppendTo[PreCeps, PreCFR]]; \)

Obliczenie cepstrum

\(Cepstrum = {}; 

    For[i = 1, i <= Length[PreCeps], \(i++\), 
	
      AppendTo[Cepstrum, 
        Abs[InverseFourier[
            PreCeps[[i]]]]]]; 

    For[i = 1, i <= Length[Cepstrum], \(i++\), 
	
      Cepstrum[[i]] = 
        Table[Cepstrum[[i, x]], {x, 
            Lram}]]; \)

Detekcja częstotliwoci podstawowej (z interpolacjš wielomianowš)

\(FS = {}; 

    For[i = 1, i <= Length[Cepstrum], \(i++\), 
	Maks = {0, 0}; 
	
      For[j = FMin, j <= FMax, \(j++\), \ 
        If[Cepstrum[[i, j]] > 
            Maks[[1]], 
          Maks = {Cepstrum[[i, j]], 
              j}]]; 
	
      If[Maks[[2]] != FMin && 
          Maks[[2]] != FMax && 
          Maks != {0, 0}, 
		
        dopol = {
            Cepstrum[[i, 
              Maks[[2]] - 1
                ]], 
            Cepstrum[[i, 
              Maks[[2
                ]]]], 
            Cepstrum[[i, 
              Maks[[2]] + 1
                ]]}; 
		
        poly = \(-Fit[dopol, {1, x, x\^2}, x]\); 
		
        korekta = 
          x /. \(FindMinimum[poly, {x, 1}]\)[[2
              ]]; 
		
        Maks[[2]] += korekta - 2]; 
	
      If[Maks == {0, 0}, 
        AppendTo[FS, FS[[Length[FS]]]], 
        AppendTo[FS, 
          CzP\/\(\ Maks[[2]]\)]]]; \)

Konwersja czasu

InS = 
      Table[{(((i - 1) OknoN (1 - Zakładka) + 1)/CzP), 
          FS[[i]]}, {i, Length[FS]}];

a = 
      ListPlot[InS, PlotRange -> {0, 1200}, PlotJoined -> True, 
        PlotStyle -> RGBColor[0, 0, 1]];

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

Obliczenie błędu detekcji

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

Wykres błędu

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

Wartoć rednia miary błędu

\(Mean[Abs[Błędzik]]\)

Tłumienie wartoci cepstrum dla wysokich częstotliwoci

\(Wyniki = {}; 

    For[k = 1, k <= 2, k += 1\/10, 
	
      Tłumik = 
        Table[Log[x + 1]\^k, {x, 1, \ 
            Length[Cepstrum[[1]]]}]; 
	
      FS = {}; 
	
      For[i = 1, i <= Length[Cepstrum], \(i++\), 
		Maks = {0, 0}; 
		
        For[j = FMin, j <= FMax, \(j++\), \ 
          If[Cepstrum[[i, j]] 
                Tłumik[[j]] > 
              Maks[[1]], 
            Maks = {Cepstrum[[i, j]] 
                  Tłumik[[j]], 
                j}]]; 
		
        If[Maks[[2]] != FMin && 
            Maks[[2]] != FMax && 
            Maks != {0, 0}, 
			
          dopol = {
              Cepstrum[[i, 
                  Maks[[2]] - 1
                    ]] 
                Tłumik
                  [[Maks[[2
                      ]] - 1]], 
              Cepstrum[[i, 
                  Maks[[2
                    ]]]] 
                Tłumik
                  [[Maks[[2
                    ]]]], 
              Cepstrum[[i, 
                  Maks[[2]] + 1
                    ]] 
                Tłumik
                  [[Maks[[2
                      ]] + 1]]}; 

          			poly = \(-Fit[dopol, {1, x, x\^2}, x]\); 
			
          korekta = 
            x /. \(FindMinimum[poly, {x, 1}]\)[[2
                ]]; 
			
          Maks[[2]] += korekta - 2]; 

        		If[Maks == {0, 0}, 
          AppendTo[FS, 
            FS[[Length[FS]]]], 
          AppendTo[FS, 
            CzP\/\(\ Maks[[2]]\)]]]; 
	
      AppendTo[Wyniki, {BłšdProc[FS], 
          Mean[Abs[Błędzik]]}]; 
	
      Print[N[k], "   ", BłšdProc[FS], "   ", 
        Mean[Abs[Błędzik]]]]; \)

Wyniki // N