• ベストアンサー

FFTの基底の選び方について教えてください。

FFTの基底の選び方について教えてください。 NETLIB-NAPACKから任意個数のFFTをダウンロードし、個数が素数の計算をしたところ計算時間がものすごくかかります。例えば、N=16411で複素数の場合、そのままだと約15秒、後にゼロをつけて2のべき乗(N=32768)にすると0.02秒、N=17496(=2^3*3^7、試行で求めた)にすると0.03秒となります。そこで後につけるゼロの個数をできるだけ減らし、かつ元個数に最も近いような個数（例えば今回の17496)を、2,3,5,7のべき乗で表す個数（今回の17496)を求めるアルゴリズム、もしくはプログラムは無いものでしょうか。

質問者が選んだベストアンサー

  • ベストアンサー
  • Tacosan
  • ベストアンサー率23% (3656/15482)
回答No.1

2, 3, 5, 7 のべきだけで与えられた N 以上の最も小さい数を作ればいいわけだから, O((log N)^4) のアルゴリズムは一瞬で考えられる. もちろん, 2^p 3^q 5^r 7^s の p, q, r, s を全部調べるだけ.

pipiruru11
質問者

お礼

ありがとうございます。 言われてみれば目からうろこが落ちたような感じです。 今、実際にFFTの前段にプログラミング（Fortrannです）してますが、 てこずっています。

その他の回答 (3)

  • Tacosan
  • ベストアンサー率23% (3656/15482)
回答No.4

4乗時間でよければ, 単純に 4重ループを組むだけです. 3乗時間にするなら, そのうち内側の 2重ループを工夫すればいい.

pipiruru11
質問者

お礼

何度もご親切な回答、ありがとうございます。ご教示の方法でうまくいきました。もう一つ、素因数分解で2,3,5,7以外の数値が出る場合は元の個数を1ずつ増やし繰り返す方法も試しましたが、こちらもOKでした。

  • Tacosan
  • ベストアンサー率23% (3656/15482)
回答No.3

O((log N)^3) まではちょっとがんばれば落ちる. もっと落ちるかなぁ?

  • f272
  • ベストアンサー率46% (7998/17100)
回答No.2

#1の言うように作ってみると，やっぱり一瞬で 2^4 * 3^1 * 5^0 * 7^3 = 16464 が見つかる。 2^4 * 3^1 * 5^0 * 7^3 = 16464 2^5 * 3^1 * 5^2 * 7^1 = 16800 2^0 * 3^0 * 5^0 * 7^5 = 16807 2^1 * 3^5 * 5^1 * 7^1 = 17010 2^1 * 3^0 * 5^2 * 7^3 = 17150 2^7 * 3^3 * 5^1 * 7^0 = 17280 2^3 * 3^7 * 5^0 * 7^0 = 17496

pipiruru11
質問者

補足

ありがとうございます。 今、実際にFFTの前段にプログラミング（FortranとVBAです）してますが、 てこずっています。最速時間となるNでなくても、ご教示のN=16464でも良いので、総当たり戦であたっているのですが、うまくいかず、四苦八苦してます。

  1. カテゴリ
  2. 学問・教育
  3. 数学・算数

関連するQ&A

  • MATLABでのFFTについて

    MATLAB という汎用数値解析プログラムを使って、人でとった心電図のRR間隔の周波数解析を試みております。R波の検出し、512ポイントのRR間隔(RRI)のデータを得ました{R波をとるときの時間をtとし、RRI(n) = t(n+1)-t(n)}。ここからデータを補完してリサンプリングし、データポイントを2のべき乗にそろえ、FFTに持っていきたいと思います。 そこで質問は、＃1．RRIを得る時点で、512ポイントにする意味はありますか？ ＃2．最終的にFFTに行く前に2のべき乗のポイント数になれば問題ないでしょうか？

  • MDCT→FFTへの変換？（オーディオ信号）

    いまオーディオ圧縮アルゴリズムの勉強をしています。 IMDCTの式が 　　 y[n]= (2/N)*Σ｛x[k]*cos((2π/N)*(n+a)*(k+0.5))｝ 　　　　　　　　k ただし0≦n<N,　aは定数, Σはk=0,1,…N/2-1 とあるのですが、これをFFTで計算するには どのような処理をすればいいのでしょうか？ また、cos(x)の計算は何回すればいいのでしょうか？ ＃cos(x)の値をテーブルで持っておくには、 ＃何種類の値をもっていればいいのでしょうか？ 本やネットで調べたら、DCTの式はあるのですが、 MDCTについてはあまり触れられていないので… よろしくお願いします。

  • 高速フーリエ変換でデータ数が２のべき乗でない時

    こんにちは。現在、フーリエ変換について勉強しているのですが、ちょっとわからないことがあったので質問させていただきました。 質問内容は高速フーリエ変換についてで、cooley＆tukeyのアルゴリズムを利用すると、データが２の冪乗個のときは計算量をО(NlogN)に減らせる事ができるというものでした。 しかしデータが2の冪乗個でないとき。例えばN=5000くらいのときはデータを切り取って無理やりN=4096（=2^12)みたいな感じにすれば良いんですよね？ やっぱりその時って、N=5000で通常の離散フーリエ変換したときと周波数値に誤差が出ると思うのですが、それはどうやったら計算できるのでしょうか。。。 どなたかご教授していただければ幸いです。

  • 画像をFFTした際のスペクトル分布について

    現在、森北出版から出版されている「画像処理とパターン認識入門 基礎からVC#/VC++.NETによるプロジェクト作成まで」という本の8章を参考にテクスチャーマッチングに用いるスペクトル分布を求めようとしています。 この本のソースコードをダウンロードし、スペクトル分布を計算してグラフ化すると本に掲載されている通りになります。 私は現在FFTを行うためにMIST(Media Integration Standard Toolkit)というライブラリを利用させていただいています。 このライブラリを用いてFFTを行い、パワースペクトルを画像化すると、本のソースを用いて作成したパワースペクトルの画像と見た目は一致します。 しかし、スペクトル分布はまったく違ったものになります。 本に書かれているしているスペクトル分布は2種類あり、本から引用すると「動径方向分布」と「角度方向分布」です。以降「」内は本からの引用です。 「動径方向分布」は、「中心からr=√k^2+l^2の距離に存在する環状領域内のパワースペクトルの和」のことです。 本ではこのグラフがでこぼこしているのに対して、私が作成したプログラムでは右肩上がりのグラフや右下下がりのグラフになってしまいます。 「角度方向分布」は、「水平軸から角度θの線形領域内のパワースペクトルの和」のことです。 本のグラフではピークが出ているの対して、私が作成したプログラムで出力した値を用いたグラフではピークがでません。 グラフ化にはエクセルを使用しています。 スペクトル分布の計算部分は基本的に本のソースを使用し、ライブラリによって違いが出る点を修正して使用しております。 FFTを行うアルゴリズムは異なると思うのですが、スペクトル分布がアルゴリズムによって異なることはあるのでしょうか？ また、これらの違いの原因にはどのようなことが考えられるでしょうか？ 何か気づかれた点などありましたら、よろしくお願い致します。

  • 2次元データの複素フーリエ変換するコードの作成

    数値計算等の2次元や3次元の空間データ（実数）をFFTによって複素フーリエ変換する実際のプログラム化についてお尋ねします。プログラムの実装ということなので実際的な質問で長文になっています。すみません。 まず、手持ちに1次元のFFTプログラムがあるということを前提とします（逆フーリエ変換すると、元の実数の系列が出ることは確認済のコード）。そして2次元配列の実数のデータがあるとします。この2次元のデータを2次元の複素フーリエ成分に変換することが目的です。（私の分野では波数空間への展開ということになり、複素数ですから位相情報も含まれることになります。） 例えば、x,y方向に16x16のデータあるとすると、 do j=1,16 ここでjを固定してi:1～16の実数データについて1次元のFFTをかける。 このとき、FFTにかける16個の実数データを複素数の実部に入れて、虚部はゼロとする。 FFTの出力も複素数となっている。 ここで出てくる複素フーリエ変換の結果は実部・虚部で前半（0～7）であり、後半（8～15）はその対称とか点対称（符号が逆）とかになっている（虚部をゼロとしているから）。それを複素数の2次元配列として保存する。 enddo 次いで、 do i=1,16 iを固定してj方向にFFTをかける。このとき、FFTに放り込むデータは上記の複素フーリエ変換の出力結果である2次元データを使う。具体的には複素数の2次元データをj方向の1次元の複素数配列にコピーしてFFTをかけて、その出力結果を新たな2次元配列の複素数に保存する。 enddo この結果、得られた2次元の複素数のデータが、私の所望のデータである、ということです。 式が指し示すとおりのことをすればいいのだ、ということに尽きるのだろうと思いますが、アルゴリズム的にアンバランスのように見えてこれでいいのかなと思えてしまいます。最初に虚部をゼロにするというようなこととかです。そのため確信が持てません。また、結果を見てもわかりにくい面があります。 このような考え方で実装するということいいのでしょうか。全く間違っているでしょうか。もしその場合、考え方の間違いを指摘して頂けると助かりますが（根本的な間違いだったら指摘しようがないということにもなりますが。） また、例えば、始めから実数の2次元配列をすぐに2次元複素数の実部に入れて、虚部をゼロとしてそこからコード方がすっきりするのかなと思いますが。 この辺が確定すると、3次元は同じことということになります。 サンプルコードがネットに出ているという面もありますが、自分でやる方が組み込みやすいのでお尋ねしました。 長文で申し訳ありませんが、よろしくお願いします。

  • FFT（高速フーリエ変換）の考え方

    ///////////質問(1)//////////////////// 離散フーリエ変換が以下の式で表せることは理解できます。 　 　　 N-1　　 　 　　nk 　 　　 　 　 -2πi / N Ａ ＝　 Σ　　ａ 　 Ｗ ，　　 Ｗ ＝ ｅ 　 　n 　 k=0　　　k 　　 　　 　 　 　 　 　 n：基本周波数のn倍 N：データ数 ここで、実際のプログラム上では、 (1)akに サンプリングしたデータとサンプリング周期を掛けたもの を入れ、 (2)それをWと掛けあわせて行くと思います。 では、Wには何を代入すればいいのでしょうか？ Wを実部・虚部にわける、もしくは、 絶対値・位相にわけて計算する方法で あっているのでしょうか？ Anは周波数n成分の面積を表していると思うのですが、面積に虚数概念がで てくるのも変な話なような気がして、質問させていただきました。 ///////////質問(2)//////////////////// 　　　N/2-1　　　　2nk　　　n　N/2-1　　　　　2nk Ａ　＝　Σ　ａ　　Ｗ　＋　Ｗ　　Σ　　ａ　　Ｗ 　n　　k=0　　2k　　　　　　　　k=0　　2k+1 Cooley-Tukey 型 FFTは質問(1)の式を上記のように２項に分解することで、 計算を減らすことができるとのことです。 しかし、私には理解できません。 左側の項でN/2回の掛け算を行い、右側の項でもN/2回の掛け算を行うのでは、 結局式を分解しただけで、何も変わっていないように思えます。 どのように考えればいいのでしょうか？ アドバイスよろしくお願いします。

  • ハミング窓関数とFFT(高速フーリエ変換)を用いて音(WAVE)ファイルを周波数解析したいのですが

    FFTを用いてドレミの周波数解析のプログラムを作ろうと勉強しています。 そこで理論的にわかっていないため質問させていただきました。 ハミング窓の使い方についてです。 ドレミファソラシドと1secごとに変化していくWAVEファイルがあった場合。(0-1秒＝ド　1-2秒＝レ) １：0～50msハミング窓を与えてやる ２：その区間をFFTする。(データ数は512で考える) ３：求まった基本周波数からドレミの判定を行う ４：50msずらす ５：１に戻る ※汚いですが手書きで自分の↑のイメージを書いて見ました。http://www.jpdo.com/link/1/img/49864.jpg このようにして9秒まで分割して計算すればドレミの解析が出来るはずだと思っているのですが、どうやってハミング窓を0～50msなどという風に使えばいいのかがわかりません。 このあたりの概念がさっぱりなので、 どうすれば0～50msの区間をハミング窓で得られるのか。 ハミング窓で得た区間をFFTするにはどう考えるのか。 を是非教えてください。

  • アルゴリズム

    アルゴリズムの勉強をしていて、時間計算量に関する問題があり、解いたのですが、解答が載ってなく困ってます。正誤の判断と、もし間違っているなら、何が間違っているのかを教えて頂けると助かります。 ある問題において、大きさ(データ量)nに対して、アルゴリズムA、B、Cの時間計算量が、それぞれn、n^2(nの2乗)、2^n(2のn乗)であるとする。 (1)アルゴリズムAを用いて10秒間にn=100の問題が解けた。20秒かけるとき扱える問題の大きさnの値を求めよ。 解) n=100*2 =200 (2)ある計算機を用いてアルゴリズムBで10秒間にn=100の問題が解けた。100倍早い計算機を用いたとき、10秒間に扱える問題の大きさを求めよ。 解) 求める問題の大きさをxとおくと n=(100^2)*100 =10000*100 =1000000 (3)アルゴリズムCを用いて1時間にn=20の問題が解けた。n=40の問題を解くのに何時間かかるか。 解) 2^40=(2^20)*(2^20) =1*(2^20) =2^20[時間]

  • 反復試行の確率と最大値の問題について

    式が読みにくくなるため、べき乗は「^」, 割り算は「/」で表します。 数学IAの反復試行の確率と最大値の問題について途中計算について質問があります。 Pn＝n-1C2（1/5）^2(4/5)^n-3×1/3 ちなみに、n-1C2は確率に出てくる“組み合わせの数”のことなので、“シーのエヌ引く１・２”のことです。 　 (n-1)!　　　　　 　4^n-3 ＝-------　×　----------------- 2! (n-3)! 5^2×5^n-3×5 　 (n-1)(n-2)4^n-3 解答の途中式＝--------------- 2・5^n 上記の解答の途中式までの途中計算を省かずに教えて頂けると嬉しいです。 よろしくお願いします。

  • 数学的帰納法

    問題 任意の自然数nに対して5・2^n+(-4)^n-1をある素数pで割った時の余りが常に1になるとする時のpの値を求めよ。 解説は添付の資料の通りです。 n=1,2を代入してp=5であるところまでは出来ます。その後帰納法を使った証明で、 途中の解答に、 突然、なぜ漸化式が出てきたのかがわかりません。漸化式を使う必要性はなんですか？ 計算が簡単だから？ 通常の帰納法のように解答するという方法はダメなのでしょうか？今回の証明は特別に漸化式を使わないと解けない問題だということでしょうか。

注目のQ&A

カテゴリ

専門家に質問してみよう

職業から探して質問する

あなたにピッタリな商品が見つかる！ OKWAVEセレクト

質問する