HSPポータル
サイトマップ お問い合わせ


HSPTV!掲示板


未解決 解決 停止 削除要請

2014
0220
down周波数を取得したい9解決


down

リンク

2014/2/20(Thu) 18:53:15|NO.60151

題名のとうり周波数を取得したいです。
fftのプラグイン等みたのですが、どうすれば周波数を取得するのか全く分かりませんでした。
サンプルでも良いのでスクリプトとか教えていただきたいです。
(fftのプラグインを使っていただいていると嬉しいです。)
お願いします。



この記事に返信する


pippi

リンク

2014/2/21(Fri) 03:44:04|NO.60159

最近FFTの質問多いですね。
どうぞ。

beki=10
LIST_SIZE=1<<beki;(2^10=1024)
LIST_SIZE2=LIST_SIZE>>1
LIST_SIZE8=LIST_SIZE<<3
ddim x,LIST_SIZE;実数部
ddim y,LIST_SIZE;虚数部
周波数その1=50.0;Hz
周波数その2=100.0;Hz
周波数その3=150.0;Hz
周波数その4=200.0;Hz
周波数その5=250.0;Hz
周波数その6=300.0;Hz
周波数その7=350.0;Hz

//////////////;実数部にwaveデータ
repeat LIST_SIZE
x.cnt=1.0*sin(周波数その1*6.283185307179586*cnt/LIST_SIZE)+2.0*sin(周波数その2*6.283185307179586*cnt/LIST_SIZE)
x.cnt+=1.0*sin(周波数その3*6.283185307179586*cnt/LIST_SIZE)+2.0*sin(周波数その4*6.283185307179586*cnt/LIST_SIZE)
x.cnt+=1.0*sin(周波数その5*6.283185307179586*cnt/LIST_SIZE)+2.0*sin(周波数その6*6.283185307179586*cnt/LIST_SIZE)
x.cnt+=1.0*sin(周波数その7*6.283185307179586*cnt/LIST_SIZE)
loop
//////////////

////////////wave波形表示
screen 1,LIST_SIZE,128
pos 0,0
repeat LIST_SIZE
line cnt,64+int(x.cnt*8.0)
loop
////////////


fft x,y,1

///////////;√(実数部^2+虚数部^2) が特定周波数の振幅を表す
screen 0,LIST_SIZE,128
color 128,128,128:boxf LIST_SIZE/2,0,LIST_SIZE,128
color 0,0,0
pos 0,0
repeat LIST_SIZE
line cnt,64+int(-sqrt(y.cnt*y.cnt+x.cnt*x.cnt)/16.0)
loop
///////////////
stop







#deffunc fft array A,array B,int inv
ddim c,LIST_SIZE
ddim d,LIST_SIZE

repeat beki;fft
bekii=LIST_SIZE>>(cnt+1)
repeat LIST_SIZE2
t2 = cnt\bekii;
t0 = (cnt/bekii)*bekii*2 + t2;
t1 = t0+bekii;
r4=A(t1)
i4=B(t1);
rm0=A(t0)-r4
im0=B(t0)-i4;
pii=3.14159265358979323846264338328*t2/bekii*inv;
dsin=sin(pii)
dcos=cos(pii);
A(t0) += r4;
B(t0) += i4;
A(t1) = rm0*dcos-im0*dsin;
B(t1) = rm0*dsin+im0*dcos;
loop
loop


repeat LIST_SIZE;ビット逆順
ic=cnt/65536
ic = (ic&0x00005555)<<1 | (ic&0x0000AAAA)>>1;
ic = (ic&0x00003333)<<2 | (ic&0x0000CCCC)>>2;
ic = (ic&0x00000F0F)<<4 | (ic&0x0000F0F0)>>4;
ic = (ic&0x000000FF)<<8 | (ic&0x0000FF00)>>8;
iic=cnt\65536
iic = (iic&0x00005555)<<1 | (iic&0x0000AAAA)>>1;
iic = (iic&0x00003333)<<2 | (iic&0x0000CCCC)>>2;
iic = (iic&0x00000F0F)<<4 | (iic&0x0000F0F0)>>4;
iic = (iic&0x000000FF)<<8 | (iic&0x0000FF00)>>8;
ic=(iic*32768+ic/2)>>(31-beki);
C(cnt)=A(ic)
D(cnt)=B(ic)
loop

if inv=-1{
repeat LIST_SIZE
a.cnt=c.cnt/LIST_SIZE
b.cnt=d.cnt/LIST_SIZE
loop
}else{
memcpy a,c,LIST_SIZE8
memcpy b,d,LIST_SIZE8
}
return




サンプリング点数が1024ならば、1/2の512Hzまでが観測できる最大の周波数となります。



down

リンク

2014/2/21(Fri) 21:15:51|NO.60172

有難うございます!
それで何ですが、ddim x,LIST_SIZE ○○(waveファイル名)
とすればよいのですか?
ご回答宜しくお願いします!!



pippi

リンク

2014/2/22(Sat) 06:56:20|NO.60184

waveファイルの構造はご存じですか?
.wavは整数で波形を保存しているのに対しこのソースではdouble型で波形を扱っています。
.wavを読み込んだ後なんとか波形をdouble型に変換してください
fft解析はそれからです



down

リンク

2014/2/22(Sat) 16:11:35|NO.60193

ありがとうございました。



down

リンク

2014/2/22(Sat) 16:12:47|NO.60194

ありがとうございました。



なたで

リンク

2014/2/27(Thu) 20:36:48|NO.60302

少し補足です。

下記は恐らく誤りです。
>サンプリング点数が1024ならば、1/2の512Hzまでが観測できる最大の周波数
正しくは、次のようになっています。
サンプリング周波数が1024Hzならば、1/2の512Hzまでが観測できる最大の周波数

また、「double型にしたい」という質問でもあったのですが、
http://hsp.tv/play/pforum.php?mode=all&num=60214

ここでいう、double型に変換というのは、
単純に整数から実数へ変換というものではなくて、
一般的には、-1≦x<1の間にするという意味だと思います。

これは、8ビットのwavの場合は、-128〜127
16ビットのwavの場合は、-32768〜32767のため、
単純に、実数へ変換してしまうと、量子化ビット数によって
音量が変わってしまうからです。



down

リンク

2014/2/27(Thu) 21:30:06|NO.60307

"一般的には、-1≦x<1の間にするという意味だと思います。"
これは、例えば10とかだったら、
−1以上1未満の数に、何らかの計算によってするという事ですか?
ご回答お願いします。



kanahiron

リンク

2014/2/27(Thu) 22:40:19|NO.60310

(自分の勉強を兼ねて)実験がてら

screen 0,250,800 wds = 1024*1024*10 whs = 44 sr = 48000 sdim wavedata,wds+1 sdim waveheader,whs+1 lpoke waveheader,0,0x46464952 ;"RIFF" lpoke waveheader,4,wds+whs-8 ;データサイズ+ヘッダーサイズ-8 lpoke waveheader,8,0x45564157 ;"WAVE" lpoke waveheader,12,0x20746D66 ;"fmt " lpoke waveheader,16,16 ;fmt チャンクのバイト数 wpoke waveheader,20,1 ;フォーマットID (01) wpoke waveheader,22,2 ;モノラルは1,ステレオは2 lpoke waveheader,24,sr ;サンプリングレート lpoke waveheader,28,sr*4 ;データ速度(サンプリングレート*4) wpoke waveheader,32,4 ;ブロックサイズ (Byte/sample*チャンネル数) wpoke waveheader,34,16 ;サンプルあたりのビット数 lpoke waveheader,36,$61746164 ;"data" lpoke waveheader,40,wds ;波形データ rad = 0.0 rad2 = 0.0 repeat wds/4 rad += 1000.0/48000.0*2 rad2 += 1000.0/48000.0*2 a = int(sin(rad )*32768) b = int(sin(rad2)*32768) wpoke wavedata,cnt*4,a wpoke wavedata,cnt*4+2,b loop bsave "テスト.wav",waveheader,wds+whs,0 bsave "テスト.wav",wavedata,wds+whs,whs ;exec "テスト.wav" ,16 ;再生できる ;/* repeat wds/4 as = wpeek(wavedata,cnt*4) //処理 if as > 32768{ d = as - 65536 } else { d = as } // title strf("%05d,%05d",as,d) redraw 0 color 255,255,255 boxf color pos 10,380 mes "wpeekで読み込んだデータ" pos 130,410 mes "処理したデータ" line 0,400,250,400 pos 100,as/250+400 mes "●" pos 130,d/250+400 mes "●" redraw 1 await 16 loop ;*/ stop //もし実数にするなら ddim data,wds/4 repeat wds/4 as = wpeek(wavedata,cnt*4) if as > 32768{ d = as - 65536 } else { d = as } data(cnt) = 0.000030517578125*d;double(d)/32768でも同じ loop stop



down

リンク

2014/3/1(Sat) 16:07:20|NO.60354

すみません···
色々それについて考えていて、それから返信使用と思ったので···
ありがとうございます。(文章がぐちゃぐちゃでスミマセン)



なたで

リンク

2014/2/28(Fri) 22:43:37|NO.60331

wavフォーマットについて調べているかとおもいますが、
wav形式では、音声波形を8ビットの整数か、16ビットの整数などで書き込みます。

そのため、その10という数字が、
16ビットの整数で表されていたとしたなら 1/(2の15乗)
8ビットの整数で表されていたとしたなら 1/(2の7乗)となります。
音声を信号処理するときは、変換を行うことが多いです。

kanahiron様のスクリプトを参考にするといいです。

なお、整数を直接実数とみなしたままでもFFTはできます。

wavのフォーマットからFFTや、エフェクターとかに興味があるなら
青木さんの「C言語ではじめる音のプログラミング」がおすすめです。
http://floor13.sakura.ne.jp/book03/book03.html



ONION software Copyright 1997-2021(c) All rights reserved.