Transformata fourier și transformata Fourier Digitală
Definitie și formulă matematică
Transformata Fourier a F(w) a unui semnal f(t), în condiţiile formulei Fourier de inversare (transformata Fourier inversă), se mai numeşte spectru de frecventă al semnalului f(t) pentru orice frecvenţă w (măsurată în Hz) avem:
iar fkse numeşte eşantionul spectrului lui x pe frecvenţa k
Digital fourier transformation – DFT:
Modul de calcul al spectrului unui semnal in timp discret x[n], n=0,…,N-1, se bazează pe utilizarea funcţiilor sinsi cos din biblioteca matematică a limbajului C. După ce se aplică formulele prezentate mai jos rezultatul este o mărime complexă, z.
Ceea ce ne interesează în cazul analizei de frecvență a unui semnal digital este modulul acestui număr complex, având semnificaţia fizica de amplitudine a semnalului la diferite frecvenţe, practic ne interesează doar P[k].
Practic eu am decupat semnalul sonor folosind o fereastră de tip Hamming, și am calculat spectrul de frecvență pe această bucățică din semnal. Apoi am deplasat fereastra cu câteva sample-uri la dreapta și am repetat procesul:
Rezultatul analizei de frecvență pentru un singur cuvânt:
Implementare
Clasa DFT_Spectrum de mai jos contine pe lângă metoda de determinare a spectrului (dft) și toate metodele pentru decupare a semnalului acustic, MelScale , convolutie și filtrare :
<br />//------------------------clasa pentru DSP analiza spectrala si compresie-------<br />class DFT_Spectrum<br />{<br />public:<br /> double s[1024];<br /> Complex spec[1024];<br /> double mels[1024];<br /> double cep[1024];<br /> double filtreTriunghiulare[1024][16];<br /> double mas;<br /> double esantioane[80000];<br /> double Window[3000];<br /> double WindowS[1024];<br /> double Convolutia[3000];<br /> //-----------metodele--------------------------------------------<br /> void dft(int ns)<br /> {<br /> double suma_real=0,suma_imag=0,Magnitudine;<br /> int h=0,k=0,j=0;<br /> double re,im;<br /> double c1,c2;<br /> c1=2.0*M_PI/double(ns);<br /> c2=2.0/double(ns);<br /> Complex z;<br /> h=0;<br /> for (k=0; k<ns/2; k++)<br /> {<br /> h++;<br /> suma_real = 0;<br /> suma_imag = 0;<br /> for(j=0; j<ns; j++)<br /> {<br /> suma_real += double(esantioane[j])*cos(double(k)*double(j)*c1);<br /> suma_imag += double(esantioane[j])*sin(double(k)*double(j)*c1);<br /> }<br /> spec[h].real =c2*double(suma_real);<br /> spec[h].imaginar=c2*double(suma_imag);<br /> Magnitudine =z.modul(spec[h]);<br /> s[h] =20.0*log10(1+Magnitudine);<br /><br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void SetEsant_Double(double Esant[80000],int nr)<br /> {<br /> for(int i=0; i<nr; i++)<br /> {<br /> esantioane[i]=Esant[i];<br /> }<br /> }<br /> //--------------------------------------------------------------<br /> void DCT(int ns)<br /> {<br /> static double suma=0;<br /> static int h=0,k=0,j=0;<br /> static double c1;<br /> c1= M_PI/double(ns);<br /> h=0;<br /> for(k=0; k<1024; k++)<br /> {<br /> h++;<br /> suma = 0;<br /> for(j=0; j<ns; j++)<br /> {<br /> suma+= mels[j]*cos((double(k)-0.5)*double(j)*c1);<br /> }<br /> cep[h]=suma;<br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void Logaritmare(int n)<br /> {<br /> for(int i=0; i<n; i++)<br /> {<br /> mels[i]=20.0*log10(1 + s[i]);<br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void GenHamming(int wn)<br /> {<br /> for(int i=0; i<wn-1; i++)<br /> Window[i]=0.54-0.46*cos(2.0*double(i)*1.0/wn*M_PI);<br /> }<br /> //------------------------------------------------------------<br /> void GenRandomSpectrum()<br /> {<br /> randomize();<br /> for(int i=0; i<1024; i++)<br /> s[i]=random(1000)/100.0;<br /> }<br /> //-------------------------------------------------------------<br /> void GenBlackman(int wn)<br /> {<br /> for(int i=0; i<wn-1; i++)<br /> Window[i]=0.35875-0.48829*cos(2.0*double(i)*1.0/wn*M_PI)<br /> +0.14128*cos(4.0*double(i)*1.0/wn*M_PI)<br /> -0.01168*cos(6.0*double(i)*1.0/wn*M_PI);<br /> }<br /> //-------------------------------------------------------------<br /> void GenHann(int wn)<br /> {<br /> for(int i=0; i<wn-1; i++)<br /> Window[i]=0.5*(1-cos(2.0*double(i)*1.0/wn*M_PI));<br /> }<br /> //-------------------------------------------------------------<br /> void GenXHann(int wn)<br /> {<br /> for(int i=0; i<wn-1; i++)<br /> {<br /> if(i<=wn/4-1)<br /> Window[i]=0.5*(1-cos(2.0*double(i)*1.0/(wn/2)*M_PI));<br /> else<br /> {<br /> if((i>=wn/4)&&(i<=3*wn/4-1))<br /> Window[i]=1.0;<br /> if(i>=3*wn/4)<br /> Window[i]=0.5*(1-cos(2.0*double(i)*1.0/(wn/2)*M_PI));<br /> }<br /> }<br /> }<br /> //------------------------------------------------------------<br /> void GenRect(int wn)<br /> {<br /> for(int i=0; i<wn-1; i++)<br /> Window[i]=1.0;<br /> }<br /> //-------------------------------------------------------------<br /> void GenTriunghi(int wn)<br /> {<br /> int h=0;<br /> for(int i=0; i<wn-1; i++)<br /> if(i<wn/2)<br /> {<br /> Window[i]=double(i)/double(wn);<br /> h++;<br /> }<br /> else<br /> {<br /> Window[i]=double(h)/double(wn);<br /> h--;<br /> }<br /> }<br /> //------------------------------------------------------------------<br /> void GenFiltreTriungiulare(int wn)<br /> {<br /> int h=1,lungime=0,nrFer=0,a32=0,depl=wn;<br /> wn=1;<br /> a32=0;<br /> do<br /> {<br /> h=1;<br /> if(nrFer==0)<br /> {<br /> for(int i=1; i<=depl; i++)<br /> {<br /> filtreTriunghiulare[i][nrFer]=1.0;<br /> }<br /> a32=depl;<br /> }<br /> else<br /> {<br /> wn=wn*2;<br /> for(int i=1; i<=wn; i++)<br /> {<br /> if(i<wn/2)<br /> {<br /> filtreTriunghiulare[i+a32][nrFer]=double(i)/double(wn);<br /> h++;<br /> }<br /> else<br /> {<br /> filtreTriunghiulare[i+a32][nrFer]=double(h)/double(wn);<br /> h--;<br /> }<br /> }<br /> }<br /> if(nrFer==0)<br /> {<br /> a32=depl;<br /> nrFer++;<br /> }<br /> else<br /> {<br /> nrFer++;<br /> a32=wn+depl;<br /> }<br /> }<br /> while(nrFer<11);<br /> }<br /> //-------------------------------------------------------------<br /> void FiltruTriunghiular(int nrFiltru)<br /> {<br /> for(int i=0; i<512; i++)<br /> {<br /> WindowS[i] = s[i]* filtreTriunghiulare[i][nrFiltru];<br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void BubbleSort()<br /> {<br /> double aux=0.0;<br /> bool gata;<br /> do<br /> {<br /> gata=true;<br /> for(int i=1; i<1024; i++)<br /> if(WindowS[i+1]>WindowS[i])<br /> {<br /> aux=WindowS[i+1];<br /> WindowS[i+1]=WindowS[i];<br /> WindowS[i]=aux;<br /> gata=false;<br /> }<br /> }<br /> while(gata!=true);<br /> }<br /><br /> //-------------------------------------------------------------<br /> void Convolutie(int n,int indx)<br /> {<br /> static double s=0;<br /> static double tau=0;<br /> static double tx=0;<br /> for(int tau=0; tau<(n-1); tau++)<br /> {<br /> s=0;<br /> for(int i=0; i<(n-1); i++)<br /> {<br /> tx=tau-i;<br /> if(tx<0)<br /> {<br /> tx=tx+n;<br /> }<br /> s+=esantioane[i+indx]*Window[i];<br /><br /> }<br /> Convolutia[tau]=s;<br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void PseudoConvolutie(int n,int indx)<br /> {<br /> for(int i=0; i<n-1; i++)<br /> {<br /> Convolutia[i]=0;<br /> Convolutia[i]=esantioane[i+indx]*Window[i];<br /><br /> }<br /> }<br /> //-------------------------------------------------------------<br /> void SumConvolutie(int n,int indx)<br /> {<br /> for(int i=0; i<(n-1); i++)<br /> {<br /> Convolutia[i]=esantioane[i+indx]+Window[i];<br /><br /> }<br /> }<br /> //-------------------------------------------------------------<br /> int getMelscale(int numar)<br /> {<br /> //-------compresia MelScale-----<br /> double nr_frecv=2,suma=0,frecventa_prag=1;<br /> int j=0,index=1,indexp=0,h=0;<br /> //am 512 frecvente<br /> do<br /> {<br /> if(frecventa_prag>1000)<br /> nr_frecv=8;<br /> suma=0;<br /> for(int i=1; i<=nr_frecv; i=i+1)<br /> {<br /> suma+=s[i+j+int(nr_frecv)]; //suma a n frecvente <br /> index=(i+j)/nr_frecv;<br /> }<br /> mels[index]=suma/nr_frecv; //media aritmetica a n frecvente<br /> j=j+nr_frecv;<br /> frecventa_prag = frecventa_prag + 8000.0/(numar);<br /> }<br /> while(j<numar/2);<br /> return index;<br /> }<br /><br />}; Documentație
- https://ro.wikipedia.org/wiki/Transformata_Fourier
- https://en.wikipedia.org/wiki/Fast_Fourier_transform
An nou fericit !
Pentru întrebari și/sau consultanță tehnică vă stau la dispozitie pe blog sau pe email simedruflorin@automatic-house.ro. O seară/zi plăcută tuturor !