On Mar 19, 8:01�am, "Akshaya Mukund"
<akshaya.mukund@n_o_s_p_a_m.gmail.com> wrote:
> I have been trying to code a normal C program for a radix 3 fft problem. My
> program works for 3,6,9,27 points but conks off for the 81 point fft
> onwards. I would like to know if there is a mistake in my program or is it
> because VC can't support such intensive calculations ;P. Well I really
> don't know how to go ahead.I thought maybe I had calculated my twiddle
> factors wrong and I have even tried writing out the complete equation(81
> terms) for one particular output. Nothing wrong with my algo (but I could
> definitely be wrong as well). Can anyone please check it out. Here is my
> code.
>
> The possible sources of error could be my butterfly function (you can chk
> it out below), or maybe the fact that I have messed up some type casting.
>
> If anyone gets lucky with this please let me know. thank you for your
> time.
>
> Cheers
> Akshaya
>
> /////////////////////////////////////////////////////////////////
>
> #include <stdio.h>
> #include <stdlib.h>
> #include <math.h>
>
> struct dft
> {
> � � � � double r,i;
>
> };
>
> #define pi (double)2*asin((double)1)
>
> double* data_input(char*, double*); �//used to input data
> int trittoint(int *trit,int m); � � �//convert trit to int
> void tritrev(int *trit, int m); � � �//reverse the trit
> void tritconvert(int m, int i,int* trit); � �//convert to trit
> double twiddle_geni(int N,int k,int m); � //imag part of twiddle factor
> double twiddle_genr(int N,int k,int m); � //real part
> void butterfly(struct dft *a,int aidx,int stage,int bno,int sets,int N);
> � � � //butterfly structure (radix-3)and also twiddle factor mult
>
> int main()
> {//to do a 3,6,9,27 just change N to no of pts and m to no. of trits
> � � � � FILE *fid, *fid1;
> � � � � int i,j,k,l,N=81,m=4,trit[4],limit,limit1,aidx=0; � � � � � � � � � � � � � � � //3^m = N (no of
> trits)
> � � � � struct dft a[81],b[81],*temp;
> � � � � char *path = "data2.txt";
> � � � � double temp1,t2[81],*t1;
>
> � � � � //init data for fft
> � � � � for(i=0;i<N;i++)
> � � � � {
> � � � � � � � � b[i].r = b[i].i = 0;
> � � � � � � � � a[i].r = 0;
> � � � � � � � � a[i].i = 0;
> � � � � � � � � c[i].r=c[i].i=0;
> � � � � }
>
> � � � � //data entered into input struct var from file
> � � � � //t1 = data_input(path,t2);
> � � � � for(i=0;i<N;i++) � � � � � � � �// i'm entering arbitrary data here
> � � � � � � � � a[i].r �= i;
>
> � � � � //reordering loop
> � � � � //reordered into proper order trit conversion then reversal
> � � � � // this part is fine but you are welcome to chk it out
>
> � � � � for(i=1;i<2*N/3;i++)
> � � � � {
> � � � � � � � � tritconvert(i,m,trit); � � � � �//convert to trits
> � � � � � � � � tritrev(trit,m); � � � � � � � � � � � �//reverse the trit
> � � � � � � � � j = trittoint(trit,m); � � � � �//convert to integer
>
> � � � � � � � � if(j>i) � � � � � � � � � � � � � � � � � � � � � � �//swap
> � � � � � � � � {
> � � � � � � � � � � � � temp1 = a[j].r;
> � � � � � � � � � � � � a[j].r = a[i].r;
> � � � � � � � � � � � � a[i].r = temp1;
> � � � � � � � � }
> � � � � }
>
> � � � � //fft loop
>
> � � � � for(i=0;i<m;i++) � � //loop takes care of the stage
> � � � � {
> � � � � � � � � limit = (int)pow((double)3,(double)(m-1-i)); //sets of
> � � � � � � � � � � � � � � � � � � � � � � � � � � � //radix 3
> butterflies
> � � � � � � � � aidx = 0;
> � � � � � � � � for(j=0;j<limit;j++) //no of sets of butterflies
> � � � � � � � � {
> � � � � � � � � � � � � limit1 = (int)pow((double)3,(double)i);
> � � � � � � � � � � � � for(k=0;k<limit1;k++) � � � �//no of butterflies per set
> � � � � � � � � � � � � {
> � � � � � � � � � � � � � � � � butterfly(a,aidx,i,k,limit,N);
> � � � � � � � � � � � � � � � � aidx++;
> � � � � � � � � � � � � }
> � � � � � � � � � � � � aidx += (int)pow((double)3,(double)(i+1))-limit1;
> � � � � � � � � }
> � � � � }
>
> � � � � //writing the fft onto text file
> � � � � fid = fopen("datafftreal.txt","wb");
> � � � � fid1 = fopen("datafftimag.txt","wb");
> � � � � for(i=0;i<27;i++)
> � � � � {
> � � � � � � � � temp = &(a[i].r);
> � � � � � � � � fwrite(temp,sizeof(double),1,fid);
> � � � � � � � � temp = &(a[i].i);
> � � � � � � � � fwrite(temp,sizeof(double),1,fid1);
> � � � � }
>
> � � � � fclose(fid);
> � � � � fclose(fid1);
>
> � � � � return (0);
>
> }
>
> double* data_input(char* path, double* input)
> {
> � � � � FILE *fid;
> � � � � int i;
>
> � � � � fid = fopen(path,"rb");
> � � � � fread(input,8,27,fid);
>
> � � � � for(i=0;i<27;i++)
> � � � � {
> � � � � � � � � //svbarve
> � � � � � � � � input[i] = 00.0;
> � � � � � � � � printf("%lf\t",input[i]);
> � � � � }
>
> � � � � fclose(fid);
>
> � � � � return (input);
>
> }
>
> void butterfly(struct dft *a,int aidx,int stage,int bno,int sets,int N)
> {
> � � � � int l,k,i,loopinc,n=3;
> � � � � struct dft temp[3];
>
> � � � � loopinc = (int)pow((double)3,(double)stage);
>
> � � � � for(l=aidx,i=0;i<3;i++)
> � � � � {
> � � � � � � � � temp[i].r = a[l].r;
> � � � � � � � � temp[i].i = a[l].i;
> � � � � � � � � a[l].r=0;
> � � � � � � � � a[l].i=0;
> � � � � � � � � l += loopinc;
> � � � � }
>
> � � � � k = aidx;
> � � � � for(i=0;i<3;i++)
> � � � � {
> � � � � � � � � for(l=0;l<3;l++)
> � � � � � � � � {
> � � � � � � � � � � � � a[k].r += (temp[l].r*twiddle_genr(n,l,i) -
> temp[l].i*twiddle_geni(n,l,i))*twiddle_genr(N,l,bno*sets);//doubt
> � � � � � � � � � � � � a[k].r -= (temp[l].r*twiddle_geni(n,l,i) +
> temp[l].i*twiddle_genr(n,l,i))*twiddle_geni(N,l,bno*sets);
> � � � � � � � � � � � � a[k].i += (temp[l].r*twiddle_geni(n,l,i) +
> temp[l].i*twiddle_genr(n,l,i))*twiddle_genr(N,l,bno*sets);//doubt
> � � � � � � � � � � � � a[k].i += (temp[l].r*twiddle_genr(n,l,i) -
> temp[l].i*twiddle_geni(n,l,i))*twiddle_geni(N,l,bno*sets);
> � � � � � � � � }
> � � � � � � � � k += loopinc;
> � � � � }
>
> }
>
> void tritconvert(int i, int m, int* trit)
> //used inreordering loop (conv to trit)
> {
> � � � � int j=0;
> � � � � for(j=0;j<m;j++)
> � � � � � � � � trit[j]=0;
>
> � � � � j=0;
> � � � � while(i>0)
> � � � � {
> � � � � � � � � trit[j]=i%3;
> � � � � � � � � i = i/3;
> � � � � � � � � j++;
> � � � � }
>
> }
>
> void tritrev(int *trit,int m) � //used in reordering loop (rev the trits)
> {
> � � � � int temp,i,j;
> � � � � for(i=0,j=m-1;i<m/2;i++,j--)
> � � � � {
> � � � � � � � � temp = trit[i];
> � � � � � � � � trit[i] = trit[j];
> � � � � � � � � trit[j] = temp;
> � � � � }
>
> }
>
> int trittoint(int *trit,int m) �//used in reordering loop(convert to int)
> {
> � � � � int i,j=0,z=1;
> � � � � for(i=0;i<m;i++)
> � � � � {
> � � � � � � � � j = j + trit[i]*z;
> � � � � � � � � z*= 3;
> � � � � }
> � � � � return (j);
>
> }
>
> double twiddle_genr(int N,int k,int m)
> {
> � � � � double te;
> � � � � te = cos((-2*k*m*pi)/N);
> � � � � return(te);}
>
> double twiddle_geni(int N,int k,int m)
> {
> � � � � double te;
> � � � � te = sin((-2*k*m*pi)/N);
> � � � � return(te);
>
>
>
> } � � �- Hide quoted text -
>
> - Show quoted text -
I'm sure the problem is in your code and not a limitation of VC.
Clay