プロが教える店舗&オフィスのセキュリティ対策術

高速フーリエ変換(FFT)のプログラム
高速フーリエ変換のプログラムが必要になり、visual c++で作成しています。
http://www.kiso.tsukuba.ac.jp/~watanabe/FFT.htm
こちらのサイトの一番下のプログラムを参考にして作成しています。
ここで質問なのですが信号データが8点のプログラムが上記のサイトにはあり、私は作成したプログラムが正しいかどうかチェックするため4点の信号データのFFTを求めたいのですが、以下のようにx0=2,x1=2,x3=8,x4=-4の値をフーリエ変換するプログラムにしたところ、この4点フーリエ変換の正規の値と、このプログラムの出力値が違ってしまいます。私の知識不足も否めませんがどこが間違っているのかわからずこまっています。
よろしければ教えてください。

void four1(float datar[],float datai[],unsigned long nn, int isign){
unsigned long n,mmax,m,j,istep,i;
double wtemp,wr,wpr,wpi,wi,theta;
float tempr,tempi;

n=nn;
j=0;
for (i=0;i<n ;i++) {
if (j > i) {
SWAP(datar[j], datar[i]);
SWAP(datai[j], datai[i]);
}
m=n >> 1;
while (m >= 1 && j > m) {
j -= m;
m >>= 1;
}
j += m;
}
mmax=1;
while (n > mmax) {
istep=mmax << 1;
theta=isign*(6.28318530717959/istep);
wtemp=sin(0.5*theta);
wpr = -2.0*wtemp*wtemp;
wpi=sin(theta);
wr=1.0;
wi=0.0;
for (m=0;m<mmax;m++) {
for (i=m;i<n;i+=istep) {
j=i+mmax;
tempr=wr*datar[j]-wi*datai[j];
tempi=wr*datai[j]+wi*datar[j];
datar[j]=datar[i]-tempr;
datai[j]=datai[i]-tempi;
datar[i] += tempr;
datai[i] += tempi;
}
wr=(wtemp=wr)*wpr-wi*wpi+wr;
wi=wi*wpr+wtemp*wpi+wi;
}
mmax=istep;
}

}
void realft()
{
float datar[4]={2.0,2.0,8.0,-4.0};
float datai[4]={0.0,0.0,0.0,0.0};
int i;

four1(datar,datai,4,-1);

for(i=-1;i<3;i++){


printf("i = %d:\tReal: %f, \tImag: %f\n", i, datar[i]/4,datai[i]/4);

}

}

A 回答 (1件)

while (m >= 1 && j > m) {


じゃなくて
while (m >= 1 && j >= m) {
だろうね。
それから
for(i=-1;i<3;i++){
これは
for(i=0;i<4;i++){
だろうね。

この回答への補足

早速の回答ありがとうございます。
値が合いました。ありがとうございます。
一つ質問なのですが、x0~x3までの信号データをフーリエ変換して求まる周波数スペクトル値はX-1~X2で与えられると教科書等に書いてあります。
iの値が0~3をとるとするとX-1の値がi=3のときに与えられますよね。
この場合X-1=1/4*{x0+jx1-x2-jx3}=X3 なのでしょうか。X-1がi=3の場合に対応するだけなのでしょうか。

 理解が浅く意味がわからない質問になってしまっているかもしれませんがよろしくお願いします。

補足日時:2010/11/08 17:58
    • good
    • 0

お探しのQ&Aが見つからない時は、教えて!gooで質問しましょう!