Matlab では、Welch の方法 ( pwelch) を使用してパワー スペクトルを頻繁に計算し、それを両対数プロットに表示します。によって推定された周波数pwelchは等間隔ですが、対数対数プロットには対数間隔の点の方が適しています。特に、プロットを PDF ファイルに保存する場合、高頻度でポイントが過剰になるため、ファイル サイズが非常に大きくなります。
線形間隔の周波数から対数間隔の周波数まで、スペクトルをリサンプリング (再ビン化) する効果的なスキームは何ですか? または、過度に大きなファイル サイズを生成せずに PDF ファイルに高解像度スペクトルを含める方法は何ですか?
行うべき明白なことは、単に使用することinterp1です:
rate = 16384; %# sample rate (samples/sec)
nfft = 16384; %# number of points in the fft
[Pxx, f] = pwelch(detrend(data), hanning(nfft), nfft/2, nfft, rate);
f2 = logspace(log10(f(2)), log10(f(end)), 300);
Pxx2 = interp1(f, Pxx, f2);
loglog(f2, sqrt(Pxx2));
ただし、これはスペクトル内の電力を節約しないため、望ましくありません。たとえば、新しい周波数ビンの 2 つの間に大きなスペクトル線がある場合、結果のログ サンプリングされたスペクトルから単純に除外されます。
これを修正するには、代わりにパワー スペクトルの積分を補間します。
df = f(2) - f(1);
intPxx = cumsum(Pxx) * df; % integrate
intPxx2 = interp1(f, intPxx, f2); % interpolate
Pxx2 = diff([0 intPxx2]) ./ diff([0 F]); % difference
これはかわいらしく、ほとんどの場合機能しますが、ビンの中心が正しくなく、周波数グリッドがより細かくサンプリングされる可能性のある低周波数領域をインテリジェントに処理しません。
その他のアイデア:
- 新しい周波数ビニングを決定
accumarrayし、リビニングを行うために使用する関数を作成します。 - 補間を行う前に、スペクトルに平滑化フィルターを適用します。問題: スムージング カーネル サイズは、目的の対数スムージングに適応する必要があります。
- この
pwelch関数は周波数ベクトル引数 を受け入れますf。この場合、Goetzel アルゴリズムを使用して目的の周波数で PSD を計算します。pwelchそもそも対数間隔の周波数ベクトルで呼び出すだけで十分かもしれません。（これは多かれ少なかれ効率的ですか？） - PDF のファイル サイズの問題については、スペクトルのビットマップ イメージを含めます (見栄えが悪い - きれいなベクター グラフィックスが必要です!)。
- または、スペクトルを示す単純な線分ではなく、領域(多角形/信頼区間) を表示することもできます。