前回紹介した窓関数に関するウィキペディアですが、下のような図が載っています。
この図の意味するところが最初はピンと来なかったのですが、要は観測時間を限定して切り出すとその瞬間にこういうスペクトルが重畳されて現れると云うことのようです。
早い話が、直流を時間で区切って矩形波にするとこういうスペクトラムになると考えて良いでしょう。
ウィキペディアにいくつもの種類の窓関数が紹介されていますが、何がお勧めというわけではなくその時の都合で選べ、という感じです。
何を選ぶかの基準は、周波数分解能（周波数 0Hz が尖塔である）とダイナミックレンジ（ノイズスペクトラムが小さい）ということになります。
前の解説で FFT は切り取られた部分が延々に繰り返されるという前提で分析する、と書きましたが、切り取ったという時点でこのようなノイズが乗るということになります。その時に観測時間の中で整数倍で繰り返される信号しかなければ、上述の矩形波パワースペクトラムをみての通り、整数の周波数（BIN といいます）ではレベルが 0 ですから正しく観測できます。ところが端数の周波数に対してはノイズが加算されるというか、そもそも端数の周波数の存在自体を DFT = FFT は認めていないのですから、そういった信号が入ってくれば観測結果が狂ってきます。これは前々回確認しました。
で、今回は矩形波で切り取った場合とカイザー関数といわれる柔軟性があるのでよく使われる窓関数で切り取った場合の比較をしてみることにします。カイザー窓の場合はこのような図が載っています。
まずはそれぞれの窓関数の形と特性を以下のスクリプトで比較してみます。Wikipeda に載っている図の手元確認というところです。
今までと違うのは端数の周波数でも窓関数を掃引していることです。
xdel(winsid()) //ウィンドをクリア
clear // Data をクリア
N=2048
t=0:1/N:1-1/N
// x=sin(2*%pi*t*f1);
// y=cos(2*%pi*t*f1);
Signal=1;
Rect=window('re',N)
Kaiser=window('kr',N,%pi*3) // Kaiser window を選択。αは 3
scf(1);plot(1:N,Rect,1:N,Kaiser)
mtlb_axis([0,N,0,1.5])
legend('Rectungular','Kaiser')
count=0
fmax=10
fstep=0.1
for f1=0:fstep:fmax // 調査周波数
count=count+1
Rect_combo_sin(count,:)=sin(2*%pi*t*f1).*Rect;
Rect_combo_cos(count,:)=cos(2*%pi*t*f1).*Rect;
Rect_sin_level(count)=sum(Rect_combo_sin(count,:))./N*2
Rect_cos_level(count)=sum(Rect_combo_cos(count,:))./N*2
Rect_level(count)=((Rect_sin_level(count))^2+(Rect_cos_level(count))^2)^0.5
Kaiser_combo_sin(count,:)=sin(2*%pi*t*f1).*Kaiser;
Kaiser_combo_cos(count,:)=cos(2*%pi*t*f1).*Kaiser;
Kaiser_sin_level(count)=sum(Kaiser_combo_sin(count,:))./N*2
Kaiser_cos_level(count)=sum(Kaiser_combo_cos(count,:))./N*2
Kaiser_level(count)=((Kaiser_sin_level(count))^2+(Kaiser_cos_level(count))^2)^0.5
end
scf(1000);clf
plot(0:fstep:fmax,(Rect_level'),0:fstep:fmax,(Kaiser_level'))
legend('Rectungular','Kaiser')
scf(1001);clf
wt_axis={0,fmax,-100,20]
plot(0:fstep:fmax,20*log10(Rect_level'),0:fstep:fmax,20*log10(Kaiser_level'))
mtlb_axis(wt_axis)
legend('Rectungular','Kaiser')
窓の形はこうなっています。
青線が矩形波です。当然のことながら観測区間では DC に見えてしまいます。
緑線がカイザー窓です。時間軸の両端をゼロにすることで窓内の信号に周期性がなくても強制的に周期性があるように見せます。
周波数特性はこうなります。青線は矩形波窓、緑線はカイザー窓です。デシベル表示です。
スクリプトを見ていただければ分かるように、今までの手製のスペクトラムアナライザをそのまま使っています。で、ちょっと気になったのが矩形波窓の場合の 0Hz のレベルが 6dB なってしまっていることで、これは計算式を追っていけば自明なのですが今までの評価結果と整合性がとれてません。何か勘違いしているのかもともと別のものを評価しているのか今すぐには説明出来ません。
これらから読み取れることは、矩形波窓だと観測周期に整数倍の周波数はノイズレスで出力されますが、端数の周波数を含むと大きなノイズを含む、カイザー窓の場合は 0 ～ 2 Hz にはノイズが現れますが、それより高い周波数では端数の周波数が入ってきてもノイズとしては小さそうです。
次回は色々な周波数をこれらの窓を通して FFT を掛けて測定してみたいと思います。
この図の意味するところが最初はピンと来なかったのですが、要は観測時間を限定して切り出すとその瞬間にこういうスペクトルが重畳されて現れると云うことのようです。
早い話が、直流を時間で区切って矩形波にするとこういうスペクトラムになると考えて良いでしょう。
ウィキペディアにいくつもの種類の窓関数が紹介されていますが、何がお勧めというわけではなくその時の都合で選べ、という感じです。
何を選ぶかの基準は、周波数分解能（周波数 0Hz が尖塔である）とダイナミックレンジ（ノイズスペクトラムが小さい）ということになります。
前の解説で FFT は切り取られた部分が延々に繰り返されるという前提で分析する、と書きましたが、切り取ったという時点でこのようなノイズが乗るということになります。その時に観測時間の中で整数倍で繰り返される信号しかなければ、上述の矩形波パワースペクトラムをみての通り、整数の周波数（BIN といいます）ではレベルが 0 ですから正しく観測できます。ところが端数の周波数に対してはノイズが加算されるというか、そもそも端数の周波数の存在自体を DFT = FFT は認めていないのですから、そういった信号が入ってくれば観測結果が狂ってきます。これは前々回確認しました。
で、今回は矩形波で切り取った場合とカイザー関数といわれる柔軟性があるのでよく使われる窓関数で切り取った場合の比較をしてみることにします。カイザー窓の場合はこのような図が載っています。
まずはそれぞれの窓関数の形と特性を以下のスクリプトで比較してみます。Wikipeda に載っている図の手元確認というところです。
今までと違うのは端数の周波数でも窓関数を掃引していることです。
xdel(winsid()) //ウィンドをクリア
clear // Data をクリア
N=2048
t=0:1/N:1-1/N
// x=sin(2*%pi*t*f1);
// y=cos(2*%pi*t*f1);
Signal=1;
Rect=window('re',N)
Kaiser=window('kr',N,%pi*3) // Kaiser window を選択。αは 3
scf(1);plot(1:N,Rect,1:N,Kaiser)
mtlb_axis([0,N,0,1.5])
legend('Rectungular','Kaiser')
count=0
fmax=10
fstep=0.1
for f1=0:fstep:fmax // 調査周波数
count=count+1
Rect_combo_sin(count,:)=sin(2*%pi*t*f1).*Rect;
Rect_combo_cos(count,:)=cos(2*%pi*t*f1).*Rect;
Rect_sin_level(count)=sum(Rect_combo_sin(count,:))./N*2
Rect_cos_level(count)=sum(Rect_combo_cos(count,:))./N*2
Rect_level(count)=((Rect_sin_level(count))^2+(Rect_cos_level(count))^2)^0.5
Kaiser_combo_sin(count,:)=sin(2*%pi*t*f1).*Kaiser;
Kaiser_combo_cos(count,:)=cos(2*%pi*t*f1).*Kaiser;
Kaiser_sin_level(count)=sum(Kaiser_combo_sin(count,:))./N*2
Kaiser_cos_level(count)=sum(Kaiser_combo_cos(count,:))./N*2
Kaiser_level(count)=((Kaiser_sin_level(count))^2+(Kaiser_cos_level(count))^2)^0.5
end
scf(1000);clf
plot(0:fstep:fmax,(Rect_level'),0:fstep:fmax,(Kaiser_level'))
legend('Rectungular','Kaiser')
scf(1001);clf
wt_axis={0,fmax,-100,20]
plot(0:fstep:fmax,20*log10(Rect_level'),0:fstep:fmax,20*log10(Kaiser_level'))
mtlb_axis(wt_axis)
legend('Rectungular','Kaiser')
窓の形はこうなっています。
青線が矩形波です。当然のことながら観測区間では DC に見えてしまいます。
緑線がカイザー窓です。時間軸の両端をゼロにすることで窓内の信号に周期性がなくても強制的に周期性があるように見せます。
周波数特性はこうなります。青線は矩形波窓、緑線はカイザー窓です。デシベル表示です。
スクリプトを見ていただければ分かるように、今までの手製のスペクトラムアナライザをそのまま使っています。で、ちょっと気になったのが矩形波窓の場合の 0Hz のレベルが 6dB なってしまっていることで、これは計算式を追っていけば自明なのですが今までの評価結果と整合性がとれてません。何か勘違いしているのかもともと別のものを評価しているのか今すぐには説明出来ません。
これらから読み取れることは、矩形波窓だと観測周期に整数倍の周波数はノイズレスで出力されますが、端数の周波数を含むと大きなノイズを含む、カイザー窓の場合は 0 ～ 2 Hz にはノイズが現れますが、それより高い周波数では端数の周波数が入ってきてもノイズとしては小さそうです。
次回は色々な周波数をこれらの窓を通して FFT を掛けて測定してみたいと思います。