Red de conocimiento informático - Conocimiento del nombre de dominio - Cómo programar en C para implementar cualquier número de FFT

Cómo programar en C para implementar cualquier número de FFT

Esto no se puede explicar claramente en sólo una o dos frases.

FFT requiere que el número de puntos de sincronización de entrada sea una potencia entera de 2, como 1024, 2048...

int get_pow_2(int NP)

{

int j = 1, k, m = 1;

mientras ( m lt; NP ) { m = m * 2; >

return m;

}

Debe asegurarse de que esté remuestreando, cambie el número de puntos de tiempo originales a la potencia entera de 2 o agréguelos más tarde. El método de suma es 0 completo, o se supone que la secuencia de tiempo se repite, es decir, la secuencia de puntos de la secuencia de tiempo original se suma hasta que se satisface la potencia entera de 2.

Para encontrar el promedio de la serie temporal, elimine el valor promedio de todos los puntos (es decir, desplace la línea central).

Filtrado de tiempos. Consideración de la ventana del filtro: rectangular, en forma de cos u otra.

Determinar el tiempo dt, df transformado, grados de libertad o banda de paso

Hacer FFT

#define SWAP(a, b) tempr=(a) ; (a)=(b);(b)=tempr

void jfour1(float ya[], unsigned long nn, int isign)

{

unsigned long n, mmax, m, j, istep, i;

doble wtemp, wr, wpr, wpi, wi, theta;

n=nn lt; lt; 1;

j=1;

para (i=1; ilt; n; i =2) {

if (j gt; i) {

SWAP(ya[j], ya[i]);

SWAP(ya[j 1], ya[i 1 ]

}

m=n gt;

mientras (m gt; = 2 amp; j gt; m) {

j -= m;

m gt; = 1

}; >

p>

};

mmax=2;

mientras (n gt; mmax) {

istep=mmax lt; ; 1;

theta=isign*(6.28318530717959/mmax);

wtemp=sin(0.5*theta);

wpr = -2.0*wtemp* wtemp;

p>

wpi=sin(theta);

wr=1.0;

wi=0.0;

for (m=1;mlt;mmax ; m =2) {

for (i=m; ilt; =n; i =istep) {

j=i mmax;

tempr = wr * ya[j]- wi * ya[j 1];

tempi = wr * ya[j 1] wi * ya[j];

ya[j] = ya[i] - tempr;

ya[j 1] = ya[i 1] - tempi;

ya[i] = tempr ;

ya [i 1] = tempi;

};

wr = (wtemp=wr) * wpr - wi * wpi wr;

wi = wi * wpr wtemp * wpi wi

};

mmax=istep

}

}

# undef SWAP

void jrealft(float ya[], unsigned long n, int isign)

{

void jfour1 (float ya[], unsigned long nn , int isign);

unsigned long i, i1, i2, i3, i4, np3, n05;

float c1=0.5, c2; , h1r, h1i, h2r, h2i

>doble wr, wi, wpr, wpi, wtemp, theta;

n05 = n gt;

theta=3.141592653589793/(doble) (n05);

p>

if (isign == 1) {

c2 = -0.5;

jfour1(ya,n05,1);

} más {

c2=0.5;

theta = -theta;

};

wtemp=sin(0.5 *theta);

wpr = -2.0*wtemp*wtemp;

wpi=sin(theta);

wr=1.0 wpr;

wi= wpi;

np3=n 3;

for (i=2;ilt;=(ngt;gt;2);i) {

i4= 1 (i3=np3-(i2=1 (i1=i i-1)));

h1r = c1 * (ya[i1] ya[i3]);

h1i = c1 * (ya[i2] - ya[i4]);

h2r = -c2* (ya[i2] ya[i4]);

h2i = c2 * (ya[i1] - ya[i3]);

ya[i1] = h1r wr * h2r - wi * h2i;

ya[i2] = h1i wr * h2i wi * h2r;

ya[i3] = h1r - wr * h2r wi * h2i;

ya[i4] = -h1i wr * h2i wi * h2r;

wr= (wtemp=wr) * wpr - wi * wpi wr;

wi=wi * wpr wtemp * wpi wi;

};

if (isign == 1) {

ya[1] = (h1r=ya[1]) ya[2];

ya[2] = h1r-ya[2 ];

} más {

ya[1] = c1 * ((h1r=ya[1]) ya[2]);

ya[ 2]=c1 * (h1r - ya[2]);

jfour1(ya,n05,-1);

}

}

Filtrado espectral (eliminando la fuga de espectro fuera de la frecuencia de corte).