体系 › 数値計算
·無料·論文·20分で読めます

ゼロから学ぶフーリエ変換 — 信号を純音の重ね合わせとして読む

時間軸に描かれた信号は、本当に知りたいことを隠している——そこにどの周波数が、どれだけ入っているのか。この記事はフーリエ変換をたった1つの考えから組み立てる——内積で信号を純音と照らし合わせる。そこから連続の積分を記号ごとに分解し、DFTにも同じことをし、2つの正弦波を手で分離し、最後に現場で必ず刺される3つの罠——エイリアシング、有限の分解能、スペクトル漏れ——で締める。

要点

  • 正弦波は振幅A・周波数f・位相φの3つの数で完全に記述され、複素指数 e^{2πift} として扱える。
  • 長さTの窓の上では、周波数の差がk/Tの音どうしは直交し、他の周波数は掛けて平均すると打ち消される。
  • フーリエ変換 X(f)=∫x(t)e^{-2πift}dt は信号と周波数fの探り針との内積で、|X(f)|が振幅、偏角が位相になる。
  • DFTではビン番号kが周波数f_k=k f_s/Nに対応し、ビン間隔Δf=f_s/N=1/Tは録音の長さだけで決まる。
  • f_sがf_maxの2倍以下だとエイリアシングが起き、Δfより近い音は分離できず、ビン中心を外れた音は漏れる。
対象textタスクnumerical

Théorie analytique de la chaleur

一次資料 — この記事の根拠

この解説の公開 2026-09-23

Théorie analytique de la chaleurFourier 1822論文ページ
Discrete-Time Signal Processing"Oppenheim & Schafer 2010
https://www.pearson.com/en-us/subject-catalog/p/discrete-time-signal-processing/P200000003226"3rd ed.

手元にあるグラフは、見たいグラフではない

ピアノの和音を1秒録音してプロットしてみてください。得られるのはゆらぐ1本の線、つまり時間に対する気圧です。そのゆらぎのどこかに、押された3つの音、部屋の暗騒音、マイクのヒスノイズが入っている——けれど、そのどれも目には見えません。このプロットが答えているのは「各瞬間の圧力はいくらだったか?」であって、あなたが実際に持っている問いは「ここにどの音が、どれだけの大きさで入っているのか?」です。

フーリエ変換は、1番目の答えを2番目の答えに変換する機械です。時間の関数を受け取り、周波数の関数を返す。そして本格的な信号処理はすべて——オーディオコーデック、JPEG、MRIの再構成、音声モデルが朝飯前に食べるスペクトログラム——まずこれを走らせることから始まります。

この記事は、この考えが実際に組み上がる順序どおりに進みます。2つの前提(純音と、2つの音が「直交する」とはどういうことか)、すべての仕事をこなすたった1つの測定器(内積)、連続の式を記号ごとに分解、実際に計算できる離散版、実例、そして最後にうまくいかない3つのこと。

前提1: 純音

この分野全体の原子が正弦波——ただ1つの純音です。

x(t)=Asin⁡(2πft+ϕ)x(t) = A \sin(2\pi f t + \phi)

3つの数がこれを完全に記述します。AA は振幅、どれだけ大きいか。ff は周波数、1秒あたりの周期数(ヘルツ)で、音の高さ——2πft2\pi f t は周期をラジアンに換算するので、tt が 1/f1/f 秒進むと角度はちょうど1回転進みます。ϕ\phi は位相、t=0t = 0 のときに波が周期のどこから始まるか。

代数を扱いやすくする約束が2つあります。1つめは、オイラーの公式を使って、音をサインではなく複素指数で書くこと。

eiθ=cos⁡θ+isin⁡θe^{i\theta} = \cos\theta + i\sin\theta

絵として読むと、eiθe^{i\theta} は単位円上の角度 θ\theta の点で、コサインがその横座標、サインが縦座標です。つまり e2πifte^{2\pi i f t} は、tt が進むにつれて1秒あたり ff 回、単位円を回る点。この1つの対象がコサインとサインを同時に持ち運ぶので、両者を置き換えられます——実数の音2つの代わりに、複素の音1つ。

2つめは、振幅と位相を1つの複素数に同居させること。複素数 c=∣c∣eiϕc = |c| e^{i\phi} は大きさと角度を持ちます——音が必要とするものがちょうど2つ揃っています。ですから「周波数 ff がどれだけ入っているか」は、大きさが振幅、角度が位相である1つの複素数として出てきます。これは小細工ではなく、変換の出力が複素数値である理由そのものです。

前提2: 直交性 — なぜ音どうしは測定を汚染しないのか

この手法全体を可能にしている事実がこれです。異なる周波数は、互いの測定を邪魔しない。 正確に述べると、長さ TT の窓の上で、

∫0Te2πif1t e2πif2t‾  dt=∫0Te2πi(f1−f2)t dt={Tf1=f20f1−f2=k/T, k≠0\int_0^T e^{2\pi i f_1 t}\,\overline{e^{2\pi i f_2 t}}\; dt = \int_0^T e^{2\pi i (f_1 - f_2) t}\, dt = \begin{cases} T & f_1 = f_2 \\ 0 & f_1 - f_2 = k/T,\ k \neq 0 \end{cases}
(1)

項ごとに読むと、上線は複素共役(虚部の符号を反転すること。eiθe^{i\theta} なら e−iθe^{-i\theta} になります)。e2πif1te^{2\pi i f_1 t} に e−2πif2te^{-2\pi i f_2 t} を掛けると指数が足され、差の周波数 f1−f2f_1 - f_2 の音が1つだけ残る。そして、ある音を自分自身の周期の整数回ぶん積分すると0になります——上がったぶんは必ず同じだけ下がるからです。差の周波数が0のときだけ、被積分関数が定数 11 になり、その積分は窓の長さ TT になります。

言葉にすると、信号にある音を掛けて平均すると、その周波数以外はすべて自分で打ち消えて消える。窓の中に整数回の周期が収まる周波数 fk=k/Tf_k = k/T は直交する族をなし、これは垂直な座標軸の連続版です。フーリエ変換とは、その軸に沿った座標を読み取るだけのもの。そして式の「f1−f2=k/Tf_1 - f_2 = k/T」という条件は、のちに漏れ(leakage)として噛みついてくる細字の注意書きです。

直感: 信号は重なった層で、1層ずつ重さを量る

さて中心の主張、これがフーリエの1822年の洞察です。まともな信号はどれも、純音の和として書ける。 近似できる、ではなく、書ける。矩形波、ドラムの一撃、話し声の母音——それぞれが正弦波の、あるレシピそのものです。信号は透明な層の重なりで、層1枚が周波数1つ。オシロスコープに見えているのは、その積み重ねを真上から見た姿です。

それが正しいなら、信号を「解析する」とは1つのこと——レシピを見つけることです。そして直交性のおかげで、全部を同時に解かずに材料を1つずつ見つけられます。「この信号に440 Hzはどれだけ入っているか?」と問うには、信号に440 Hzの音を掛けて平均する。440 Hzの層は自分自身と掛かるのでいたるところ非負になり、大きな値に平均されます。それ以外の層は、440 Hzの音と掛かると0でない差の周波数の音になり、0に平均されます。1つの数が落ちてきて、それが440 Hzの量です——ほかに何が入っていようと関係なく。

同じことの幾何的な言い方: 信号を、軸が周波数である空間のベクトルだと思ってください。掛けて平均する操作は、1つの軸との内積、つまり射影です。v⋅e^xv \cdot \hat{e}_x が yy や zz を知らなくても3次元ベクトルから xx 座標を取り出すのと同じように、e2πifte^{2\pi i f t} との内積は信号から ff 座標を取り出します。フーリエ変換は、その座標の全リストです。

FIG 1内積は「片方がもう片方の向きにどれだけ入っているか」を測る。直角(直交)にすると0になる——違う周波数の音どうしが互いの測定を汚さないのは、これと同じ理由

連続フーリエ変換を記号ごとに

内積のレシピを取り、窓を時間軸全体まで伸ばして、書き下します。

X(f)=∫−∞∞x(t) e−2πift dtX(f) = \int_{-\infty}^{\infty} x(t)\, e^{-2\pi i f t}\, dt
(2)

左から右へ読むと、

そして出力の中身: ∣X(f)∣|X(f)| は周波数 ff がどれだけ入っているか(振幅)、arg⁡X(f)\arg X(f) はその成分が周期のどこにいるか(位相)。∣X(f)∣|X(f)| を ff に対してプロットすれば、冒頭で欲しかった絵——和音の音が別々の棒として立ち上がる図——が得られます。

この変換は可逆で、それが「信号は層の重なりそのものだ」の本当の意味です。

x(t)=∫−∞∞X(f) e+2πift dfx(t) = \int_{-\infty}^{\infty} X(f)\, e^{+2\pi i f t}\, df
(3)

言葉にすると、信号を組み立て直すには、すべての周波数の音 e2πifte^{2\pi i f t} を取り、その周波数の係数 X(f)X(f) で大きさを合わせ、全部足し上げる。符号が ++ に変わっていることに注意してください——解析は共役へ射影し、合成は音そのものを足す。どちらの向きにも何も失われません。時間の絵と周波数の絵は、1つの対象の2つの記述です。

すぐ読める帰結が1つ: f=0f = 0 では探り針が e0=1e^0 = 1 なので、X(0)=∫x(t) dtX(0) = \int x(t)\,dt ——直流成分(DC)、信号全体のゼロからのずれです。いちばん低い「周波数」とは、平均値のことです。

離散版: コンピュータが実際に走らせるもの

コンピュータには無限の時間軸も、連続無限のサンプルもありません。あるのは NN 個の数 x0,x1,…,xN−1x_0, x_1, \dots, x_{N-1} で、1秒あたり fsf_s サンプルのサンプリング周波数で取られたもの。ですからサンプル nn は時刻 t=n/fst = n/f_s にあり、記録全体は T=N/fsT = N/f_s 秒続きます。積分を和で置き換え、連続無限の周波数を直交族 fk=k/Tf_k = k/T で置き換えると、式は離散フーリエ変換(DFT)になります。

Xk=∑n=0N−1xn e−2πikn/N,k=0,1,…,N−1X_k = \sum_{n=0}^{N-1} x_n\, e^{-2\pi i k n / N}, \qquad k = 0, 1, \dots, N-1
(4)

どの部分も連続版の何かに対応しています。

この最後の行が、応用フーリエ解析でもっとも役に立つ事実です。ビン番号 kk はヘルツに fk=kfs/Nf_k = k f_s / N で換算され、隣り合うビンの間隔は

Δf=fsN=1T\Delta f = \frac{f_s}{N} = \frac{1}{T}

間隔は、どれだけ長く録ったかだけで決まり、どれだけ速くサンプリングしたかには依りません。440 Hzと441 Hzを区別したい? ならば Δf≤1\Delta f \le 1 Hz が必要で、つまり少なくとも1秒ぶんの音が必要です。どれだけオーバーサンプリングしても、その代わりにはなりません。

逆変換は、符号を反転して先頭に 1/N1/N を付けた同じ和です。

xn=1N∑k=0N−1Xk e+2πikn/Nx_n = \frac{1}{N} \sum_{k=0}^{N-1} X_k\, e^{+2\pi i k n / N}

1/N1/N があるのは、離散の直交関係が k=k′k = k' のとき ∑n=0N−1e2πi(k−k′)n/N=N\sum_{n=0}^{N-1} e^{2\pi i (k - k')n/N} = N になるからです(それ以外は 00——公比が1の NN 乗根である等比級数だからです)。長さの2乗が NN の基底ベクトルへ射影すると NN 倍に数え過ぎるので、逆変換がそれを割り戻します。この定数をどちらに置くかは純粋に慣習で、だからライブラリによって規格化が違い、NN 倍のバグの温床になります。

書いたとおりに計算すると、DFTは NN 個のビンそれぞれに NN 回の掛け算で O(N2)O(N^2)。FFT は同じ数を O(Nlog⁡N)O(N \log N) で得る分割統治の組み替えです。別の変換ではなく、この式への速い道筋にすぎません。そのアルゴリズムはゼロから学ぶFFTの主題です。

実例: 2つの音を元に戻す

大きくて低い音と、小さくて高い音、2つの正弦波を混ぜます。

x(t)=sin⁡(2π⋅5t)+0.5sin⁡(2π⋅12t)x(t) = \sin(2\pi \cdot 5t) + 0.5 \sin(2\pi \cdot 12t)

fs=128f_s = 128 Hz でちょうど1秒サンプリングすると、N=128N = 128、Δf=128/128=1\Delta f = 128/128 = 1 Hz。各ビンの幅は1ヘルツで、どちらの音もビンの中心にぴったり乗ります——5 Hzはビン5、12 Hzはビン12に着地します。

import numpy as np

fs, N = 128, 128
t = np.arange(N) / fs
x = np.sin(2*np.pi*5*t) + 0.5*np.sin(2*np.pi*12*t)

X = np.fft.fft(x)
mag = np.abs(X)

for k in np.argsort(mag)[::-1][:4]:
    print(f"bin {k:3d}  =  {k*fs/N:5.1f} Hz   |X| = {mag[k]:.1f}")

ちょうど4つのピークが出てきて、そのすべてが予測できます。

bin   5  =    5.0 Hz   |X| = 64.0
bin 123  =  123.0 Hz   |X| = 64.0
bin  12  =   12.0 Hz   |X| = 32.0
bin 116  =  116.0 Hz   |X| = 32.0

読み取るべきことが3つ。1つめ、大きさ。 振幅 AA の実正弦波は高さ NA/2NA/2 のピークを作ります。sin⁡θ=(eiθ−e−iθ)/2i\sin\theta = (e^{i\theta} - e^{-i\theta})/2i が音を大きさ A/2A/2 の複素指数2つに分け、NN サンプルにわたる和がそれぞれを NN 倍するからです。ですから A=1A = 1 は 128×1/2=64128 \times 1/2 = 64、A=0.5A = 0.5 は 3232。ピークの比が振幅の比です——変換がレシピを復元しました。

2つめ、鏡像。 ビン123と116は e−iθe^{-i\theta} 側の半分です。123=128−5123 = 128 - 5、116=128−12116 = 128 - 12。実数値の入力に対しては常に XN−k=Xk‾X_{N-k} = \overline{X_k} なので、スペクトルの上半分は下半分の反射で、新しい情報を持ちません。プロットが fs/2f_s/2 で止まるのも、np.fft.rfft が存在するのもこのためです。

3つめ、位相。 生の値を印字すると X5=−64iX_5 = -64i、X12=−32iX_{12} = -32i——純虚数で、実部がありません。これはコサインではなくサインである署名です。コサインなら純実数で出てきます。位相を捨てて ∣X∣|X| だけを残すとこれが失われる。だから、振幅スペクトログラムだけから位相復元の工夫なしに音声を再構成できないのです。

3つの罠を、刺される順に

エイリアシング: 遅すぎるサンプリングは周波数を偽造する

標本化定理(ナイキスト–シャノン)は、fmax⁡f_{\max} より上の周波数を含まない信号は fs>2fmax⁡f_s > 2 f_{\max} でサンプリングすれば完全に決まる、と言います。サンプリング周波数の半分 fs/2f_s/2 がナイキスト周波数——そのサンプル列が表現できる最高の周波数です。

これを破ると、エネルギーは消えません。姿を変えて戻ってきます。 ナイキストより上の周波数 ff の音は、サンプルの瞬間においては、範囲内に収まる整数 mm に対する ∣f−mfs∣|f - m f_s| の音と区別がつきません。上の例に fs=128f_s = 128 のまま np.sin(2*np.pi*100*t) を足してみると、100 Hzはナイキスト(64 Hz)を越えており、ピークはビン28に現れます——128−100=28128 - 100 = 28 だからです。完全にきれいな100 Hzの音が、完全にきれいな28 Hzの音として報告され、それが偽物だと示すものはデータの中に何もありません。対処は予防しかありません——サンプラーの前に置くローパスのアンチエイリアシングフィルタで、だからどのオーディオインターフェースもそれをハードウェアに持っています。

分解能: 短い窓では近い周波数を分けられない

Δf=1/T\Delta f = 1/T は硬い床です。1ビンより近い2つの音は同じビンに落ちて1つのピークに融合し、後処理では分離できません。定番の失敗は、記録をゼロ埋めして、滑らかになったスペクトルを見て喜ぶこと。ゼロ埋めは同じ情報を細かいグリッドに補間するので、ピークの位置は読みやすくなりますが、すでに融合した音を分解しはしません。そして TT はスペクトルが平均される区間でもあるので、周波数の細かさと時間の細かさを両方欲しがると正面から衝突します——スペクトログラムに窓長を1つ選んでそれと付き合うことを強いる、時間–周波数のトレードオフです。

漏れ: ビン中心でない音はスペクトル全体に滲む

直交性の条件に戻ってください。それは、周波数の差が「窓あたり整数回の周期」であるときだけ成り立ちました。窓の中で端数の周期を回る音は何とも直交せず、そのエネルギーはすべてのビンにこぼれます。同じ N=128N = 128 で sin⁡(2π⋅5.5t)\sin(2\pi \cdot 5.5 t) を取ってみましょう——窓の中で5.5周期、ビンのちょうど間です。

bin:     3      4      5      6      7      8
|X|:  10.5   15.7   42.6   39.0   12.0    6.7

64の棒1本の代わりに、約41の肩が2つと、長く伸びた裾があります。物理的には、DFTは記録が永遠に繰り返されると仮定するので、波を周期の途中で切ると継ぎ目に不連続ができ、その人工的な段差が広帯域の成分を持つのです。

標準の処置が窓関数です。両端で滑らかに0へ落ちる曲線をサンプルに掛け、継ぎ目を取り除きます。同じ信号にハン窓を掛けると、

bin:     3      4      5      6      7      8
|X|:  0.81   5.60  27.01  27.03   5.58   0.79

遠くの裾が1桁以上落ちました。その代償に注意してください。メインローブは依然2ビンに広がり(しかも前より少し広く)、窓が信号のエネルギーを削ったのでピークの高さは縮んでいます。これが窓関数の永久の取引です——サイドローブを下げる代わりに、メインローブが広がり、自分で振幅補正をかけることになる。

現場ではこう出てくる

DFTの和を手で書くことはまずありません。呼ぶのは np.fft.rfft、torch.stft、librosa.melspectrogram です。実際に時間を使うのは、それらの関数が露出させている4つのパラメータで、そのどれもが上の罠がキーワード引数の服を着たものです。

n_fft は周波数分解能です。 これが Δf=fs/nfft\Delta f = f_s / n_{\mathrm{fft}} を決めます。16 kHzの音声で n_fft=512 ならビン間隔は31 Hz——フォルマントには十分ですが、低い男声の倍音を分けるには粗すぎます。上げれば周波数の細かさを買い、時間の細かさで払います。

hop_length と window は漏れの管理です。 スペクトログラムは音声を重なり合うフレームに切ります。各フレームは有限の窓なので、各フレームが漏れます。ハン窓が既定なのには理由があり、慣習的な50〜75 %のオーバーラップは、あるフレームでテーパーで削られた端を次のフレームの中央が覆うために存在します。

規格化の慣習は違い、症状はきれいな NN 倍として出ます。 NumPyは 1/N1/N を逆変換側に置きますが、両側に 1/N1/\sqrt{N} を置くライブラリもあります。rfft は N/2+1N/2 + 1 個のビンを返すので、折り返された半分に2の係数を掛けずにパワーを合計すると過小になります。測った値が怪しくきりのいい係数だけずれているときは、物理より先に慣習を確認してください。

ナイキストより下に報告された音も、エイリアスでありうる。 説明できない音がスペクトルに見えたら、信号源の中を探す前にアンチエイリアシングの経路を確かめる。ソフトウェアでのリサンプリングにも同じ要件があります——先にローパスせずにダウンサンプリングすると、捨てられるはずの最上オクターブが可聴帯域へ折り返されてきます。

まとめ

O(Nlog⁡N)O(N \log N) のアルゴリズムと、そこから従う畳み込み定理はゼロから学ぶFFTへ。「直交基底へ射影する」という線形代数としての読み方はAIのための線形代数を参照してください。

確認問題

  1. fs=8000f_s = 8000 Hz で0.25秒の音を録音しました。ビン間隔 Δf\Delta f はいくらで、DFTはいくつのビンを返し、1000 Hzはどのビン番号に対応しますか。
  2. ある信号を fs=1000f_s = 1000 Hz でサンプリングし、アンチエイリアシングフィルタなしで700 Hzの純音が含まれていました。ピークはどの周波数に現れますか。また、なぜサンプル値だけからそれが偽物だと検出できないのでしょうか。
  3. 同僚は3 Hz離れた2つのエンジン倍音を区別する必要があります。彼は非常に高いサンプリング周波数で0.2秒録音し、65536点までゼロ埋めしたうえで「それでもピークが融合する」と報告しました。彼のやり方の何がまずいのか、何を変える必要がありますか。
  4. fs=64f_s = 64、N=64N = 64 で cos⁡(2π⋅10t)\cos(2\pi \cdot 10 t) のDFTを取ると X10=32X_{10} = 32(純実数)になり、同じ実験を sin⁡\sin で行うと X10=−32iX_{10} = -32i になります。大きさ32と、位相の違いの両方を説明してください。
解答
  1. Δf=1/T=1/0.25=4\Delta f = 1/T = 1/0.25 = 4 Hz。記録には N=fsT=8000×0.25=2000N = f_s T = 8000 \times 0.25 = 2000 サンプル入るので、DFTは2000ビンを返します(実信号なら独立なのはそのうち1001個)。1000 Hzは k=f/Δf=1000/4=250k = f/\Delta f = 1000/4 = 250 にあります。
  2. ナイキストは fs/2=500f_s/2 = 500 Hzで、700 Hzはその上なので、音は ∣700−1000∣=300|700 - 1000| = 300 Hzへ折り返されます。700 Hzの音と300 Hzの音のサンプル値は、これらの瞬間においては数値として完全に同一なので、サンプル値をどう解析しても区別できません——情報は処理で失われたのではなく、サンプリングの時点で破壊されています。防げるのはサンプラーの前のフィルタだけです。
  3. 分解能は記録長で決まります。Δf=1/0.2=5\Delta f = 1/0.2 = 5 Hzで、これは3 Hzの間隔より粗いので、2つの音は1つのビンに落ちます。サンプリング周波数は表現できる最高周波数を上げるだけで分解能は上げません。ゼロ埋めもスペクトルを細かいグリッドに補間するだけで、窓がそもそも捉えなかった区別を取り戻すことはできません。彼は長く録る必要があります——少なくとも T>1/3≈0.34T > 1/3 \approx 0.34 秒、実用上はきれいに分けるためにその数倍。
  4. 大きさ: 振幅 AA の実正弦波は大きさ A/2A/2 の複素指数2つに分かれ、それぞれが自分のビンで NA/2N A/2 に合計されます——ここでは 64×1/2=3264 \times 1/2 = 32。位相: cos⁡θ=(eiθ+e−iθ)/2\cos\theta = (e^{i\theta} + e^{-i\theta})/2 は正の半分が2つなので係数は実数になりますが、sin⁡θ=(eiθ−e−iθ)/2i\sin\theta = (e^{i\theta} - e^{-i\theta})/2i は 1/i=−i1/i = -i の因子を持ち、係数を負の虚軸へ回します。振幅が同じなので大きさは同じ。偏角が90°違うのは、それがまさにサインとコサインの間の位相差だからです。
🎬 この内容の解説動画(Manim・字幕つき)

参考文献

  1. Théorie analytique de la chaleur. Fourier 1822論文ページ
  2. Discrete-Time Signal Processing. "Oppenheim & Schafer 2010
  3. https://www.pearson.com/en-us/subject-catalog/p/discrete-time-signal-processing/P200000003226". 3rd ed.

本記事は上記論文の本文にもとづいて執筆しています。数値・主張は原典を優先してください。

運営: TechAthletes · 編集方針 · 出典の付け方 · 訂正の記録

コメント

コメントにはログインが必要です