c - Caculating the PSD using FFTW -
c - Caculating the PSD using FFTW -
i'm sound file i've record using `alsa' lib using next setups :
fs = 96000; // sample frequency channelnumber = 1 ; format =int16 ; length = 5sec;
meaning 480000 16bit value. want calculate psd of set of :
what i'm trying tto save result bunch of double value in info can plot them of evaluating them ( i'm not sure if that's correct) :
#include <stdio.h> #include <stdlib.h> #include <math.h> #include <fftw3.h> int main(){ char filename[] = "sound.raw"; char magnfile[] = "data.txt"; file* inp = null; file* oup = null; float* info = null; fftwf_complex* out; int index = 0; fftwf_plan plan; double var =0; short wert = 0; float r,i,magn; int n = 512; info =(float*)fftwf_malloc(sizeof(float)*n); out = (fftwf_complex*) fftwf_malloc(sizeof(fftwf_complex)*n); //allocating memory input info plan = fftwf_plan_dft_r2c_1d(n,data,out, fftw_measure); // opening file reading inp = fopen(filename,"r"); oup = fopen(magnfile,"w+"); if(inp== null){ printf(" couldn't open file \n "); homecoming -1; } if(oup==null){ printf(" couldn't open output file \n"); } while(!feof(inp)){ if(index < n){ fread(&wert,sizeof(short),1,inp); //printf(" wert %d \n",wert); data[index] = (float)wert; //printf(" wert %lf \n",data[index]); index = index +1; } else{ index = 0; fftwf_execute(plan); //printf("new plan \n"); //printf(" real \t imag \t magn \t \n"); for(index = 0 ; index<n; index++){ r=out[index][0]; =out[index][1]; magn = sqrt((r*r)+(i*i)); printf("%.10lf \t %.10lf \t %.10lf \t \n",r,i,magn); //fwrite(&magn,sizeof(float),1,oup); //fwrite("\n",sizeof(char),1,oup); fprintf(oup,"%.10lf\n ", magn); } index = 0 ; fseek(inp,n,seek_cur); } } fftwf_destroy_plan(plan); fftwf_free(data); fftwf_free(out); fclose(inp); fclose(oup); homecoming 0 ; }
the problem have how can implement winding function in code ? , don't think result accurate, since i'get lot of 0 in magnitude values ? ? if has illustration i'll thankful .
here simple illustration of applying "hanning" window info prior fft:
for (int = 0; < n; ++i) { data[i] *= 0.5 * (1.0 + cos(2.0 * m_pi * (double)i / (double)(n - 1))); }
c signal-processing fftw
Comments
Post a Comment