- ベストアンサー
FFTの周波数分解能を上げる方法とは?
- 通常DFTの周波数分解能(δf)を小さくするには分析対象データの時間長を延ばす必要があります。
- しかし、分析対象データの時間長を延ばさずにΔfを小さくする方法は数学的に可能ですか?
- それらがFFTのように直交性や重複演算の排除を利用して計算量を減らすことができるのかという質問です。
- みんなの回答 (6)
- 専門家の回答
質問者が選んだベストアンサー
- ベストアンサー
#2です。δfについて了解しました。以下、参考意見です。 時間窓長を伸ばさないという条件がありましたが、実際に実行するのではなく、時間窓長を伸ばしたら何が出来るのかを見てみるのも悪くないと思いました。 時間窓長をLとして正規化された三角関数系、 L^(1/2),(L/2)^(1/2)×cos(2πt/L),(L/2)^(1/2)×cos(2πt×2/L),・・・ (L/2)^(1/2)×sin(2πt/L),(L/2)^(1/2)×sin(2πt×2/L),・・・ で考えます。対象時系列をf(t)とし、f(t)をm回繰り返した時系列をF(t)とします。波数jのcos成分の振幅a(j)は、 L/2×a(j)=∫f(t)cos(2πt×j/L)dt ※0≦t≦L. などとなりますが、これをF(t)で考えた場合、f(t)がm回繰り返されただけなので、 m×L/2×a(j)=∫F(t)cos(2πt×j/L)dt=m×∫f(t)cos(2πt×j/L)dt ※F(t)で0≦t≦m×L,f(t)で0≦t≦L. となり、cos(2πt×j/L)のa(j)はf(t)でもF(t)でも同じになります。一方、時間窓長mLでは、cos(2πt×j/L/m)なども正規化された三角関数系に含まれるので(当然直交する)、同じように振幅a(j/m)を計算できす。 しかしF(t)は、0≦t≦Lでf(t)に一致しなければならないので、フーリエ分解の一意性を考えると、a(j/m)=0になると思えます。もっとも三角関数系は完全系でもあり、関数空間C[0,L]の基底ですので、[0,L]では直交しないcos(2πt×j/L/m)などをcos(2πt×j/L)の線形和で表す事は、原理的には可能でしょう。cos(2πt×j/L/m)などを{cos(2πt×j/L),sin(2πt×j/L)}でフーリエ分解して(基底変換して)グチャグチャやれば良い訳です。 ただ、そこまでして得た{cos(2πt×j/L/m),sin(2πt×j/L/m)}のスペクトル表現に、いったいどんな意味があるのか?、という疑問はあります。 例えばf(t)=sin(2πt/L)だった場合、それを何回繰り返したところで1/Lに輝線スペクトルが1本でるだけで、それ以外の事が起こったら、かなりびっくりします。「1/Lに輝線スペクトル一つだけ」というスペクトル表現以外のものが欲しいのだろうか?、という疑問です。
その他の回答 (5)
- noname_deadbeef
- ベストアンサー率50% (10/20)
>δfの刻みが小さくなったことによって設けられたところに大きな成分があれば見つかることを期待しています。すなわち見えなかったピークが見えるようになる事を期待しています。 少なくともエネルギー的には見落とされることはありません。パーセバルの定理より。 データ全長に1.5次相当の正弦波が書かれていたとしてDFTすると直流成分+1次+2次+…と分解されてしまいますが、そのデータ全長しか情報がないのであればこれでも一つの正しい答えです。データの後ろに0を追加して1.5次相当も基底に含まれるようにすれば1.5次のところにピークが出るでしょうが…鋭いピークにはなりません。
お礼
ご回答ありがとうございます。DFTというのはアナログフーリエ変換結果の0.5~1.5次の平均振幅(平均パワー)を1次、1.5~2.5次の平均を2次の振幅とする計算なのでおっしゃる通りエネルギー的に見逃されることはないと思います。 その上で波形の特徴として分解能を上げる事で卓越周期成分が検出できないかと思っています。大人しくWaveletかけたほうが早い、というのは抜きにして。
#4です。実際に2個つないでやったんですか?(^^)。 >・・・Matlabとかそのたぐいのフリーソフト(Scilabとか)を調べれば一発だったりしますか? 自分でプログラムを組んじゃう傾向があるので、Matlabとかは使わないのでわかりません。でもそういう需要があれば(需要がないとは言い切れない)、Matlabとかに該当する機能はあるような気はします。 以下も参考意見です。 リーケージは確かにあるんですよね。それにエイリアジングはデジタル処理では防げないし、元データが周期関数でない事から発生する不連続ギャップの影響もあります。これらにより、同じ時系列を繰り返しても微妙に違ったスペクトル分布になるのは、その通りだと思います。ただそんなに違ったものにはならないんじゃないか?、というのが自分の経験でした。またそういう事を気にし出すと、ラニングスペクトルなんか作れませんし・・・(^^;)。 後続に0を足す方向ですが、FFTではデータ数を2^mにするために良くやります。その時思うんですが、厳密には元データと実用的には問題ないくらいに「近い」別データを扱っている、と考えた方が妥当なんじゃないのかな?、と感じています。 それでFFTではデータ数を合わせる不便さもあり、また最近のPCが高速な事もあって、個人的には時間窓長に制限のないDFT(演算時間n^2の方)をけっこう使います。もとデータに忠実であった方が良いと感じるからです。 真実がどうかは誰か教えて下さい状態ですが(^^;)、0を足す方向を理想化すれば、前にも後にも0を無限個足して、DFTでなく連続的なフーリエ変換に持ち込むのが一番妥当かも知れませんね。もちろん途中で計算を打ち切らざるえないですが、こっちの方が疑問の余地がない気がします。
お礼
ご回答ありがとうございます。 > #4です。実際に2個つないでやったんですか?(^^)。 はい。正確にはそうなる事はわかってて、やってそうなることを確認しましたと言ったところです。 実は2つつないで1サンプルでも足すとリーケージの具合が変わってくるのでA(2j-1)≠0にすることができ、スペクトルの形もあまり変わりません。 後続のゼロ埋めで何に使うのかにもよりますが、建築・土木の振動解析なんかではさほど影響がないなどと教科書には書かれています。画像変換や通信では何か顕著な問題が起こるかもかも。 そもそも後続のゼロがどうこう言うレベルであればリーケージ対策の窓関数なんて論外なのと、もしも短いデータで(窓内に収まらない長いサイン波を使ってDFTできたとして)結果は窓関数をかけたのと同じことをしているような気もしますけど。
すいません、不注意で間違いました。#3です。 まず「正規化しない」三角関数系、 1,cos(2πt/L),cos(2πt×2/L),・・・ sin(2πt/L),sin(2πt×2/L),・・・ で考えてます。それに正規化するんだったら、L^(1/2)と(L/2)^(1/2)を掛けるんじゃなくて、それらで割らなきゃ駄目ですよね(^^;)。
お礼
ちなみにこういうのはMatlabとかそのたぐいのフリーソフト(Scilabとか)を調べれば一発だったりしますか?
質問の趣旨をわかってないのかも知れませんが、周波数分解能がナイキスト周波数の事だとすれば、時系列の幅が長かろうが短かろうが、周波数分解能はサンプリング周波数で決まると思います。 #1さんの仰るように、分解能を上げたければ同一時系列幅のサンプリング点数を増やす(サンプリング周波数を上げる)しかない気がします。だってサンプリング点数以上の情報は、元データから採取しない訳ですから。
お礼
ご回答ありがとうございます。周波数分解能はナイキスト周波数(サンプル周波数の半分)ではなく、スペクトルの横軸の周波数刻みΔfの事です。
- bran111
- ベストアンサー率49% (512/1037)
お礼
ご回答ありがとうございます。 >周波数分解能は、時間窓長で決まり、Δf = 1/T となります。 この時間窓長を延ばさないという条件でお願いします。
お礼
ご回答ありがとうございます。 実はf(t)を複数回つないでF(t)にして時間窓長を延ばすという発想はDFTをする人には普通にあります。というのもFFTは2の累乗個のサンプル数で行うのが一般的で、サンプル数が足りない場合に後ろをゼロうめすることはよくあるわけです。 これはスペクトルの分解能を上げるためではなくただのフーリエ変換を高速フーリエ変換にするためです。 >これをF(t)で考えた場合、f(t)がm回繰り返されただけなので、 > m×L/2×a(j)=∫F(t)cos(2πt×j/L)dt=m×∫f(t)cos(2πt×j/L)dt > ※F(t)で0≦t≦m×L,f(t)で0≦t≦L. >となり、cos(2πt×j/L)のa(j)はf(t)でもF(t)でも同じになります。 で、f(t)を複数回つないでF(t)にして時間窓長を延ばすという発想ですが、cos(2πt×j/L)のa(j)はf(t)とF(t)で同じになるのはf(t)が狭帯域である場合です。広帯域波である場合にはそうはなりません。 すなわちF(t)がf(t)を2つつないだものであるとして、a(j)とA(2j)は狭帯域波であれば同じでA(2j-1)=0だと思いますが、広帯域波の場合にはそもそもリーケージ(漏れ誤差)があるのでa(j)=(A(2j-1)+A(2j))/2、A(2j-1)≠A(2j)です。 > ただ、そこまでして得た{cos(2πt×j/L/m),sin(2πt×j/L/m)}の > スペクトル表現に、いったいどんな意味があるのか?、という疑問はあります。 δfの刻みが小さくなったことによって設けられたところに大きな成分があれば見つかることを期待しています。すなわち見えなかったピークが見えるようになる事を期待しています。
補足
すみません。下記間違いです。 >広帯域波の場合にはそもそもリーケージ(漏れ誤差)があるのでa(j)=(A(2j-1)+A(2j))/2、A(2j-1)≠A(2j)です。 元の信号が何だろうが窓で切り取って離散化したものを窓内で周期的なサイン波に分解したのだから、それを2個つなぎ合わせてもやはり窓内で周期的(すなわち(A(2j-1)=0)になりますね。 やはり伸ばすなら後続のゼロを足してやる方法でしょう。