Ron Mayer FFT convolution algorithm -

does know mayer fft implementation (without me having take quite bit of time study code)?

i attempting perform convolution, , ifft appears produce call "mirrored" output. in other words, kernel+signal length limited n/2 , whatever occupies n=0...n/2 mirrored n=n...n/2. sort of looks expect fft in terms of negative frequencies...except it's mirror in negative time.

here convolution code:

   void convolve(struct cxtype* data,  struct cxtype* kernel, int size)      {     int i,j;     int wrksz = size;     float gain = 1.0f/((float) wrksz);       mayer_fft(wrksz, data->re, data->im);     mayer_fft(wrksz, kernel->re, kernel->im);      for(i=0;i<wrksz;i++)     {     data->re[i]*=kernel->re[i]*gain;     data->im[i]*=kernel->im[i]*gain;     }      mayer_ifft(wrksz, data->re, data->im); } 

doing same thing using gnu octave (matlab syntax equivalent unfamiliar) produces expected results, including allowing me occupy m+n-1 in signal output:

fs=48000; ts = 1/fs; nn =  1024 sincsz = floor(nn/4); sigstart = floor(nn/16); sigend = floor(nn/2); dpi=2*pi;  %window func tau=(1:sincsz)/sncsz; window=0.5*(1.0 - cos(dpi*tau)); %plot(tau,window)  %sinc filter kernel fc=5050; wc=dpi*fc; ts=1/fs; segtime=ts*sincsz; t0=-segtime/2; t=ts*((1:sincsz) - 1) + t0 ;  s=zeros(1,sincsz); s=window.*sin(wc*t)./(wc*t); s(sincsz/2+1) = 1;  %plot(t,s) fund = 1650; tt = (1:nn)*ts; signal = sin(dpi*tt*fund) + sin(dpi*tt*2*fund) + sin(dpi*tt*3*fund) + sin(dpi*tt*4*fund) + sin(dpi*tt*5*fund); signal(1:sigstart) = signal(1:sigstart)*0; signal(sigend:nn) = signal(sigend:nn)*0; %plot(tt,signal)  h=zeros(1,nn); h(1:sincsz) = s(1:sincsz);  h=fft(h); x=fft(signal);  y=h.*x;  y=ifft(y);  plot(real(y)) 

the equivalent signal , fir kernel synthesis implemented in c (not shown). using gnuplot display results of c implementation of this, know filter kernel , signal implemented identically did octave.

the difference in being done far can tell fft implementation. know if these results fundamental problem understanding of fft algorithm in general, or fht-based implementation archaic code written ron mayer? can go archived web site code using: http://reocities.com/researchtriangle/8869/fft_summary.html

now, if perform fft on block of data , perform ifft, original data expect. if modify data spectrally in way, results don't expect.

i once tried use mayer fft algorithm replace used in s. bernsee's pitch shift algorithm, , did not work @ all. used fftw3, , code worked expected. curious try same basic algorithm fftw3 see happens.

what don't know if misunderstand fundamental causing me mis-apply rmayer's implementation, or if simple artifact have work around (that is, using fft size double of expect).

doh! 1 of things forgetting put semicolon @ end of line. convolution complex multiplication in frequency domain -- doing naive point-by-point multiplication. here corrected code, shows complex multiplication. of course, there c , c++ constructs , macros/routines doing this, here brute force method pedagogical purposes: assume struct cxtype defined as:

struct cxtype { float* re; float* im; };  //and such struct should have mem allocated before sending convolve()      void convolve(struct cxtype* data,  struct cxtype* kernel, int size)          {         int i,j;         int wrksz = size;         float gain = 1.0f/((float) wrksz);         float a,b,c,d;           mayer_fft(wrksz, data->re, data->im);         mayer_fft(wrksz, kernel->re, kernel->im);          for(i=0;i<wrksz;i++)         {         a=data->re[i];         b=data->im[i];         c=kernel->re[i];         d=kernel->im[i];         data->re[i]=(a*c - b*d)*gain;         data->im[i]=(a*d + b*c)*gain;         }          mayer_ifft(wrksz, data->re, data->im);     } 

and above code snippet works when called , not produce these strange mirroring effects talking about. matlab code, works because octave/matlab hide detail of complex multiplication behind handy syntax, h.*x.

i able reproduce problem in octave multiplying , combining real , imaginary parts separately mimic mistake in c, , fixed doing similar above. said, there nothing unique rmayer's implementation of fft...only implementation of convolution amiss.