HSPTV!掲示板

あるサイトからダウンロードさせていただいた、
;SlowFFT (^^)
;
;2^13 = 8192ぐらいまでなら実用できそう。
;それ以上は要検討
;
;code by natade
;http://twitter.com/natadea
;修正BSDライセンス

#ifndef _getBitSize@SlowFFT
#module SlowFFT bitrv, samplelength, samplebit, fft_re2, fft_im2, fft_exp_re, fft_exp_im
#define ctype _getBitSize(%1) int(logf(%1)/logf(2)+0.5)

#modinit int _samplelength
ddim fft_re2, _samplelength
ddim fft_im2, _samplelength
ddim fft_exp_re, _getBitSize(_samplelength), _samplelength
ddim fft_exp_im, _getBitSize(_samplelength), _samplelength
_makeBitReverseTable bitrv, _samplelength
_makeFFTTable fft_exp_re,fft_exp_im
return

#modterm
bitrv = 0:fft_re2 = 0:fft_im2 = 0:fft_exp_re = 0:fft_exp_im = 0
return

#modfunc fft array fft_re1,array fft_im1,local len
_FFT fft_re1, fft_im1, fft_exp_re, fft_exp_im
len = length2(fft_exp_re)
repeat len
fft_re2(cnt) = fft_re1(bitrv(cnt))
fft_im2(cnt) = fft_im1(bitrv(cnt))
loop
memcpy fft_re1, fft_re2, len * 8
memcpy fft_im1, fft_im2, len * 8
return

#deffunc _FFT array f_re, array f_im, array exp_re, array exp_im,local stage,local blocks,local i,local j,local center,local blocklength,local pointlength,local re,local im
stage = length(exp_re) - 1
pointlength = length2(exp_re)
center = length2(exp_re) / 2
blocklength = 1
repeat length(exp_re)
repeat blocklength
i = cnt * pointlength
j = i + center
repeat center
re = f_re(i)
im = f_im(i)
f_re(i) += f_re(j)
f_im(i) += f_im(j)
re -= f_re(j)
im -= f_im(j)
f_re(j) = re * exp_re(stage, cnt) - im * exp_im(stage, cnt)
f_im(j) = im * exp_re(stage, cnt) + re * exp_im(stage, cnt)
i++
j++
loop
loop
stage--
blocklength *= 2
pointlength /= 2
center /= 2
loop
return

#modfunc ifft array fft_re1,array fft_im1,local len,local inv_samplelength
_IFFT fft_re1, fft_im1, fft_exp_re, fft_exp_im
len = length2(fft_exp_re)
repeat len
fft_re2(cnt) = fft_re1(bitrv(cnt))
fft_im2(cnt) = fft_im1(bitrv(cnt))
loop
inv_samplelength = 1.0 / length(fft_re2)
repeat len
fft_re1(cnt) = inv_samplelength * fft_re2(cnt)
fft_im1(cnt) = inv_samplelength * fft_im2(cnt)
loop
return

#deffunc _IFFT array f_re, array f_im, array exp_re, array exp_im,local stage,local blocks,local i,local j,local center,local blocklength,local pointlength,local re,local im
stage = length(exp_re) - 1
pointlength = length2(exp_re)
center = length2(exp_re) / 2
blocklength = 1
repeat length(exp_re)
repeat blocklength
i = cnt * pointlength
j = i + center
repeat center
re = f_re(i)
f_re(i) += f_re(j)
re = re - f_re(j)
im = f_im(i)
f_im(i) += f_im(j)
im = im - f_im(j)
f_re(j) = re * exp_re(stage, cnt) + im * exp_im(stage, cnt)
f_im(j) = im * exp_re(stage, cnt) - re * exp_im(stage, cnt)
i++
j++
loop
loop
stage--
blocklength *= 2
pointlength /= 2
center /= 2
loop
return

#deffunc _makeFFTTable array exp_re, array exp_im,local i,local s,local N
repeat length(exp_re)
i = cnt
s = - 6.28318531 / (2 << i)
repeat (2 << i)
exp_re(i,cnt) = cos(s*cnt)
exp_im(i,cnt) = sin(s*cnt)
loop
loop
return

#deffunc _makeBitReverseTable array out,int len,local i,local j,local k
dim out, len
repeat len
out(cnt) = cnt
loop
i = 0:j = 1
repeat len - 2
k = len >> 1
repeat
if k > i:break
i -= k
k >>= 1
loop
i += k
if(j < i) {
k = out(j)
out(j) = out(i)
out(i) = k
}
j++
loop
return

#global
#endif

font "ＭＳ ゴシック", 10

; フーリエ変換したい値の初期化
size = 1 << 4
mes "length -> "+size
newmod fftobj,SlowFFT,size
ddim f_re, size
ddim f_im, size
repeat size
f_re(cnt) = double(cnt)
loop
repeat size
mes ""+f_re(cnt) + " + " + f_im(cnt) + "j"
loop

; フーリエ変換
pos 0, 200
mes "start"
fft fftobj,f_re,f_im
mes "end"
repeat size
mes ""+f_re(cnt) + " + " + f_im(cnt) + "j"
loop

; 逆フーリエ変換
pos 300, 200
mes "start"
ifft fftobj,f_re,f_im
mes "end"
repeat size
mes ""+f_re(cnt) + " + " + f_im(cnt) + "j"
loop
というプログラムについてなんですが
これで声の音程を調べたいのです。
何を追加すればよいかなど教えていただければ幸いです。
お願いします。

何を追加すればいいかというと答えにくいのだが、
まずは離散フーリエ変換(DFT)そのものを理解するのが先ではないのか？
(DFTの高速化されたアルゴリズムがFFT)

別にプログラマーはそのアルゴリズムを覚える必要はなく、
ただどんな値を入れると、どんな値が返ってくるかということを知っていればいい。

checkさんが言うようにFFT（高速フーリエ変換）は、
計算時間が多いDFT(離散フーリエ変換)を高速に行うものです。

離散フーリエ変換は、結局のところコンピュータが扱うような
離散的な数値をフーリエ変換するというものです。
というわけで、まずはフーリエ変換が何なのかを知った方がいいと思います。
ちなみに、計算式は覚えなくていいです。

音の高さは周期性がある波です。
（波の幅が長いほど、低い音に聞こえる）
周期的な波の幅を調べるのが。フーリエ変換です。
それで、音声波形から、
音の高さへ変換するという流れです。

これらの説明は長くなるので、
本買って勉強するのが近道です。
音声の分析なら「ディジタル・サウンド処理入門」がおすすめです。
http://www.cqpub.co.jp/hanbai/books/30/30901.htm

そして、この本に書いてある、
「3.6　ケプストラム法」を用いるのが、
正確に声の音程を調べるのに適しています。

結構勉強大変です。

ご回答ありがとうございます。
それでなんですが、
最後の方のフーリエ変換の所を何かすればできるようになりませんかね？
そこに声のファイルを読み込むようにして・・・
だとできませんかね？
FFTは結構調べましたので詳しくはわかりませんが、
大体の意味は分かります。
回答宜しくお願いします！！

周波数を知りたいのです。
お願いです。サンプルを少し教えてください。