FFTのピーク周波数​以外の周波数のレベル​を下げることはできま​すか？

100kHz正弦波を周波数1MHzでサンプリングした信号を2048点でFFTした場合に、振幅1、10^-5、10^-10の正弦波いずれも、ピーク高-30～40dBあたりに下辺のレベルが形成されて同じような形状の波形が形成されます。 これを、共通した最小のレベルが形成されて、そこから振幅に比例した高さのピークが立つような波形にすることは可能でしょうか。 使用したコードは以下です。

Fs = 1e+6;                      % サンプリング周波数
Ts = 1/Fs;                      % サンプリング周期
Ls = 2048;                      % サンプル数
t = (0:Ls-1)*Ts;                % 時刻ベクトル

Fi  = Fs * 0.1;                 % 入力正弦波周波数
x1 = (1e-0)*sin(2*pi*Fi*t);     %           振幅1
x2 = (1e-5)*sin(2*pi*Fi*t);     %           振幅10^-5
x3 = (1e-10)*sin(2*pi*Fi*t);    %           振幅10^-10

if 0
for i=1:Ls
    %%ハン窓
    x1(i) = x1(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    x2(i) = x2(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    x3(i) = x3(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    %%ブラックマン・ハリス窓
    x1(i) = x1(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
    x2(i) = x2(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
    x3(i) = x3(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
end
end
y1 = fft(x1);
P21 = 10*log10(abs(y1/Ls));
P11 = P21(1:Ls/2+1);
y2 = fft(x2);
P22 = 10*log10(abs(y2/Ls));
P12 = P22(1:Ls/2+1);
y3 = fft(x3);
P23 = 10*log10(abs(y3/Ls));
P13 = P23(1:Ls/2+1);

f = Fs*(0:(Ls/2))/Ls;
hold off;
plot(f,P11);
hold on;
plot(f,P12);
hold on;
plot(f,P13);
xlabel('f (Hz)')
ylabel('|P| (dB)')

途中の窓関数を掛けても全体の症状は変わりません。宜しくお願い致します。

Answer 1

[Takuji Fukumoto]

まず、下記についてですが、

>ピーク高-30～40dBあたりに下辺のレベルが形成されて

データの長さLsが周期の整数倍でないことで別の成分が見えてしまっています。

Fi  = Fs * 0.1;                 % 入力正弦波周波数

なのでLsが10の倍数であればノイズフロアは量子化誤差レベルまでさがります。

その上で3つの入力信号に対して同じホワイトノイズをかけることで同じノイズフロアを指定することができます。 Lsとx1,x2,x3の定義のみを書き換えたものです

Fs = 1e+9;                      % サンプリング周波数
Ts = 1/Fs;                      % サンプリング周期
Ls = 2000;                      % サンプル数
t = (0:Ls-1)*Ts;                % 時刻ベクトル

Fi  = Fs * 0.1;                 % 入力正弦波周波数
Ti = 1/Fi;
x1 = (1e-0)*sin(2*pi*Fi*t)+1e-12*rand([1 Ls]);     %           振幅1
x2 = (1e-5)*sin(2*pi*Fi*t)+1e-12*rand([1 Ls]);     %           振幅10^-5
x3 = (1e-10)*sin(2*pi*Fi*t)+1e-12*rand([1 Ls]);    %           振幅10^-10

if 0
for i=1:Ls
    %%ハン窓
    x1(i) = x1(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    x2(i) = x2(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    x3(i) = x3(i) * (0.5-0.5*cos(2*pi*(i-1)/(Ls-1)));
    %%ブラックマン・ハリス窓
    x1(i) = x1(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
    x2(i) = x2(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
    x3(i) = x3(i) * (0.35875-0.48829*cos(2*pi*(i-1)/(Ls-1))+0.14128*cos(4*pi*(i-1)/(Ls-1))-0.01168*cos(6*pi*(i-1)/(Ls-1)));
end
end
y1 = fft(x1);
P21 = 10*log10(abs(y1/Ls));
P11 = P21(1:Ls/2+1);
y2 = fft(x2);
P22 = 10*log10(abs(y2/Ls));
P12 = P22(1:Ls/2+1);
y3 = fft(x3);
P23 = 10*log10(abs(y3/Ls));
P13 = P23(1:Ls/2+1);

f = Fs*(0:(Ls/2))/Ls;
hold off;
plot(f,P11);
hold on;
plot(f,P12);
hold on;
plot(f,P13);
xlabel('f (Hz)')
ylabel('|P| (dB)')

[Shuichi Watanabe]

ご教示ありがとうございます。

ノイズフロアを量子化誤差レベルまで下げればホワイトノイズでノイズフロア一定にし所望の形にできる、反対に量子化誤差レベルまで下がってもそのままでは一定にはならない(元の振幅に依存)ということですね。

ノイズフロアを量子化誤差レベルまで下げるにはデータ長を周期(サンプリング周期/入力信号周期)の整数倍にするということですが、入力信号周波数が未知、変動、あるいは複数の重ね合せなどによってその条件が確保できない場合に、それでも今の「ピーク-30～40dB」の幅を広げるための何か決まった方法などあるのでしょうか。窓関数がそれにあたりますか。

宜しくお願い致します。

[Takuji Fukumoto]

はい、窓関数の種類によっても変わりますし、 同じ窓関数ノイズフロアを下げる方法としては、

・ノイズを分散させる→高い周波数成分まで利用する→サンプリング周波数Fsを大きくする。

・周波数分解能を上げる→基本周波数を低くする→長時間のデータをとる→ｔが大きくなるようにLsを大きくする

ことでダイナミックレンジを広くとることができます。 トレードオフでデータ量が大きくなるので、計算に時間がかかるようになります。