NumPy fft.fftshift()徹底解説:FFTスペクトルを「中央寄せ」する基本
なぜfftshift
が必要なのか?
FFTは、信号を時間領域から周波数領域に変換します。NumPyのfft.fft()
関数を使ってFFTを実行すると、結果の配列(スペクトル)は通常、以下のような構造になります。
- 後半の半分
負の周波数成分(高い負の周波数から0 Hzに近い負の周波数へ) - 中央付近
ナイキスト周波数(サンプリング周波数の半分) - 最初の半分
低周波数成分(0 Hzに近い成分)から正の周波数成分へ
つまり、DC成分(0 Hz成分)が配列の最初のインデックスに位置し、正の周波数が続き、その後負の周波数が並びます。これは数学的な処理としては正しいのですが、人間がスペクトルを視覚的に解釈する際には不便です。通常、私たちは0 Hzを中心として、左右に対称的に正負の周波数が広がっている形式を好みます。
fftshift
の機能
fft.fftshift(x)
は、入力配列x
のゼロ周波数成分(DC成分)を配列の中心に移動させ、それに伴い他の周波数成分も適切に並べ替えます。具体的には、配列の左半分と右半分を入れ替えるような操作を行います。
これにより、FFT結果の配列は以下のような構造になります。
- 右半分
正の周波数成分 - 中央
DC成分(0 Hz) - 左半分
負の周波数成分
この配置は、グラフでプロットする際に、0 を中心とする見慣れた周波数スペクトルになります。
import numpy as np
import matplotlib.pyplot as plt
# サンプリング周波数
fs = 100
# 信号の長さ
t = np.arange(0, 1, 1/fs) # 0秒から1秒まで、fsHzでサンプリング
# 信号 (10Hzの正弦波と20Hzの正弦波の組み合わせ)
f1 = 10
f2 = 20
signal = np.sin(2 * np.pi * f1 * t) + 0.5 * np.sin(2 * np.pi * f2 * t)
# FFTを実行
N = len(signal)
yf = np.fft.fft(signal)
# 周波数軸の生成
xf = np.fft.fftfreq(N, 1/fs)
# プロット1:fftshiftなし
plt.figure(figsize=(12, 5))
plt.subplot(1, 2, 1)
plt.plot(xf, np.abs(yf))
plt.title('FFT結果 (fftshiftなし)')
plt.xlabel('周波数 (Hz)')
plt.ylabel('振幅')
plt.grid(True)
# fftshiftを適用
yf_shifted = np.fft.fftshift(yf)
xf_shifted = np.fft.fftshift(xf) # 周波数軸もfftshiftに合わせて並べ替える
# プロット2:fftshiftあり
plt.subplot(1, 2, 2)
plt.plot(xf_shifted, np.abs(yf_shifted))
plt.title('FFT結果 (fftshiftあり)')
plt.xlabel('周波数 (Hz)')
plt.ylabel('振幅')
plt.grid(True)
plt.tight_layout()
plt.show()
この例では、fftshift
を適用しない場合と適用した場合のFFT結果の振幅スペクトルを比較しています。fftshift
を適用しない左側のグラフでは、0 Hzが左端にあり、10 Hzと20 Hzのピークが正の周波数側に現れ、負の周波数成分は右側に現れています。一方、fftshift
を適用した右側のグラフでは、0 Hzが中央に移動し、10 Hzと20 Hzのピークが正の周波数側に、対応する負の周波数成分が負の周波数側に、対称的に現れているのがわかります。
- 周波数軸も
fftfreq
で生成した後、fftshift
で並べ替えることで、スペクトルと周波数軸が一致するようになります。 - 通常、
fftshift
はFFT結果のプロットや分析を行う前に適用されます。 - これにより、FFTスペクトルを視覚的に解釈しやすくなります。
fft.fftshift()
は、NumPyのFFT結果のゼロ周波数成分(DC成分)を配列の中心に移動させます。
fftshiftの適用漏れ、または間違った理解
よくあるエラー
FFT結果をプロットした際に、ゼロ周波数成分(DC成分)がグラフの中心に表示されず、端に表示される。これはfftshift
を適用し忘れているか、その役割を理解していないために起こります。
トラブルシューティング
-
周波数軸もシフトする
スペクトルを正しくプロットするためには、周波数軸もfftshift
に合わせて調整する必要があります。np.fft.fftfreq()
で生成した周波数軸もnp.fft.fftshift()
でシフトしてください。import numpy as np import matplotlib.pyplot as plt # 信号の生成 (例) fs = 100 # サンプリング周波数 t = np.arange(0, 1, 1/fs) signal = np.sin(2 * np.pi * 5 * t) # 5 Hzの正弦波 yf = np.fft.fft(signal) xf = np.fft.fftfreq(len(signal), 1/fs) # 誤ったプロット (fftshiftなし) # plt.plot(xf, np.abs(yf)) # 0Hzが左端に来る # 正しいプロット (fftshiftあり) yf_shifted = np.fft.fftshift(yf) xf_shifted = np.fft.fftshift(xf) plt.plot(xf_shifted, np.abs(yf_shifted)) plt.title('FFT結果 (fftshift適用済み)') plt.xlabel('周波数 (Hz)') plt.ylabel('振幅') plt.grid(True) plt.show()
-
fftshiftを適用する
FFTによって得られた配列(例:yf = np.fft.fft(signal)
)をプロットする前に、必ずnp.fft.fftshift(yf)
を適用してください。
ifftshiftとの混同
よくあるエラー
FFTとIFFT(逆フーリエ変換)を組み合わせた処理で、fftshift
とifftshift
の使い分けを誤り、結果がシフトしたり、期待する値にならない。
トラブルシューティング
-
「シフトしてからFFTし、FFT結果をシフトしてからIFFTし、IFFT結果をシフトする」という慣用句
FFTの演算は、0番目の要素をDC成分として扱いますが、時間領域の信号をt=0を中心として配置したい場合(例えば、フィルタ設計などで中心対称なインパルス応答を考える場合など)には、FFTの入力にifftshift
を適用する必要があります。そして、FFT結果を可視化する際にはfftshift
を適用します。IFFTの場合も同様です。一般的な正しい使用例 (中心にt=0とf=0が来るようにしたい場合)
# 時間領域の信号 x が t=0 を中心に配置されている場合 # FFT: X = fftshift(fft(ifftshift(x))) # IFFT: x = fftshift(ifft(ifftshift(X)))
これは特に2次元の画像処理などで、空間フィルタリングを行う際によく見られます。フィルタの原点(中心)が画像の中心に来るようにするために、FFTの入力に
ifftshift
を適用し、FFT後のスペクトルを可視化するためにfftshift
を適用します。多くの場合は、単にFFT結果をプロットするために
fftshift
を使うだけです。fft(signal)
の結果はDC成分が先頭にあるため、それを中心に移動させるためにfftshift(fft_result)
を使います。入力信号をシフトする必要があるかどうかは、その信号が時間領域でどのように定義されているか、そしてどのような位相特性を期待するかによります。もし、FFTの入力が既にt=0が先頭にあるような一般的な時系列データであれば、
fftshift
はFFT結果をプロットする目的でのみ使われることがほとんどです。
多次元配列(画像など)での軸の指定ミス
よくあるエラー
2D以上の配列(例: 画像)に対してfftshift
を適用する際に、シフトしたい軸(axes
引数)を誤って指定する。
トラブルシューティング
-
axes引数を確認する
fftshift
はデフォルトでは全ての軸に対してシフトを行いますが、特定の軸のみをシフトしたい場合はaxes
引数を使用します。import numpy as np # 2D配列の例 img = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]]) # 全ての軸をシフト (デフォルト) shifted_all = np.fft.fftshift(img) print("全ての軸をシフト:\n", shifted_all) # 結果例: # [[9 7 8] # [3 1 2] # [6 4 5]] # 0番目の軸のみをシフト shifted_axis0 = np.fft.fftshift(img, axes=0) print("\n0番目の軸のみをシフト:\n", shifted_axis0) # 結果例: # [[7 8 9] # [1 2 3] # [4 5 6]] # 1番目の軸のみをシフト shifted_axis1 = np.fft.fftshift(img, axes=1) print("\n1番目の軸のみをシフト:\n", shifted_axis1) # 結果例: # [[3 1 2] # [6 4 5] # [9 7 8]]
画像処理では、通常は水平・垂直両方向でシフトしたいので、
axes=None
(デフォルト) を使うことが多いでしょう。
奇数長の配列と偶数長の配列の違い
よくあるエラー
配列の長さが奇数か偶数かによって、fftshift
の挙動が微妙に異なることを考慮しない。特にDC成分の位置が偶数長の場合はlen(x)//2
番目、奇数長の場合はlen(x)//2
番目(どちらもPythonの0ベースインデックスで)に移動しますが、その「中心」の定義が少し変わります。
トラブルシューティング
-
プロットの際にはfftfreqと組み合わせる
どちらの場合でも、np.fft.fftfreq
で生成した周波数軸とfftshift
を組み合わせることで、正しい周波数とスペクトルの対応が得られます。import numpy as np # 偶数長の例 a_even = np.array([0, 1, 2, 3, 4, 5]) shifted_even = np.fft.fftshift(a_even) print(f"偶数長 {a_even} -> {shifted_even}") # 結果: [3 4 5 0 1 2] (半分を入れ替える) # 奇数長の例 a_odd = np.array([0, 1, 2, 3, 4]) shifted_odd = np.fft.fftshift(a_odd) print(f"奇数長 {a_odd} -> {shifted_odd}") # 結果: [3 4 0 1 2] (中央の0が残る)
-
数学的な定義を理解する
fftshift
は、配列を半分に分割して入れ替える操作です。偶数長の配列の場合、ちょうど真ん中で分割できますが、奇数長の場合は中心の要素が残ります。このため、シフト後の中心の定義が少し異なります。しかし、fftshift
自体は内部でこの違いを適切に処理してくれるため、ユーザーが直接何かを調整する必要はほとんどありません。
よくあるエラー
周波数領域でフィルタリングを行う際に、fftshift
の適用順序を誤ると、フィルタが正しく適用されないか、結果がシフトしてしまう。
トラブルシューティング
-
フィルタはシフト後のスペクトルに合わせる
フィルタを設計する際、通常はDC成分が中心にあるような周波数空間を想定して設計します。したがって、FFT結果をfftshift
でシフトしてからフィルタを適用し、その後のIFFTのためにifftshift
で元に戻す、という手順を踏むのが一般的です。# 例: ローパスフィルタリング # 1. 信号のFFT yf = np.fft.fft(signal) # 2. FFT結果をシフト (DC成分を中央へ) yf_shifted = np.fft.fftshift(yf) # 3. フィルタを設計 (DC成分を中心に設計) # 例: 中心が0となるローパスフィルタ rows, cols = yf_shifted.shape if yf_shifted.ndim > 1 else (1, len(yf_shifted)) center_row, center_col = rows // 2, cols // 2 # ... (フィルタマスクの生成ロジック。例えば、ガウシアンフィルタや円形フィルタなど) # filter_mask = ... (center_row, center_colを中心とする) # 4. フィルタ適用 filtered_yf_shifted = yf_shifted * filter_mask # 5. IFFTのために元の順序に戻す filtered_yf = np.fft.ifftshift(filtered_yf_shifted) # 6. IFFT実行 restored_signal = np.fft.ifft(filtered_yf)
フィルタリングにおいては、
fftshift
とifftshift
のペアが非常に重要です。フィルタ設計の基準となる「中心」と、FFTの出力における「DC成分の位置」を一致させるために使われます。
fftshift
は、FFTの結果を人間が直感的に理解しやすい形式に変換するためのツールです。主なトラブルシューティングのポイントは以下の通りです。
- プロットする際は常に
fftshift
とfftfreq
の組み合わせを忘れない。 - FFT-IFFTの往復変換で位相がずれる場合は、
ifftshift(input)
->fft
->fftshift(output)
のパターンを検討する。 - 多次元配列では
axes
引数を正しく指定する。 - フィルタリングなど周波数領域での操作を行う際は、フィルタ設計の中心と
fftshift
の適用順序を一致させる。
基本的な使用例:1次元信号のFFTとプロット
最も一般的な使用例は、1次元の時系列信号のFFT結果を可視化する際です。
import numpy as np
import matplotlib.pyplot as plt
# 1. 信号の生成
# サンプリング周波数 (Hz)
fs = 100
# 信号の期間 (秒)
T = 1
# 時間軸
t = np.linspace(0, T, int(fs * T), endpoint=False)
# 信号 (10 Hzの正弦波 + 25 Hzの余弦波)
f1 = 10
f2 = 25
signal = 2 * np.sin(2 * np.pi * f1 * t) + 1 * np.cos(2 * np.pi * f2 * t)
# 2. FFTの実行
N = len(signal) # 信号の長さ
yf = np.fft.fft(signal) # FFT計算
# 周波数軸の生成 (fftshift適用前)
xf = np.fft.fftfreq(N, 1/fs)
# 3. FFT結果のプロット (fftshiftなし)
plt.figure(figsize=(14, 5))
plt.subplot(1, 2, 1)
plt.plot(xf, np.abs(yf)) # 振幅スペクトルをプロット
plt.title('FFT結果 (fftshiftなし)')
plt.xlabel('周波数 (Hz)')
plt.ylabel('振幅')
plt.grid(True)
# 4. fftshiftを適用
yf_shifted = np.fft.fftshift(yf) # FFT結果をシフト
xf_shifted = np.fft.fftshift(xf) # 周波数軸もシフト
# 5. FFT結果のプロット (fftshiftあり)
plt.subplot(1, 2, 2)
plt.plot(xf_shifted, np.abs(yf_shifted)) # 振幅スペクトルをプロット
plt.title('FFT結果 (fftshiftあり)')
plt.xlabel('周波数 (Hz)')
plt.ylabel('振幅')
plt.grid(True)
plt.tight_layout()
plt.show()
解説
- 右側のグラフ (
fftshift
あり) では、0 Hzがグラフの中心に移動し、正負の周波数成分が0 Hzを中心に左右対称に配置されています。これにより、周波数スペクトルが直感的に理解しやすくなります。10 Hzと25 Hzのピークが正の周波数側に、−10 Hzと−25 Hzのピークが負の周波数側に明確に現れています。 - 左側のグラフ (
fftshift
なし) では、0 Hzがグラフの左端に位置し、正の周波数成分(10 Hzと25 Hzのピーク)がその右に、負の周波数成分がグラフの右端付近に現れています。これはFFTの数学的な出力順序です。
2次元配列(画像)に対する使用例
画像処理において、2D-FFT(2次元高速フーリエ変換)を実行した際にもfftshift
は非常に役立ちます。周波数成分を画像の中心に配置することで、高周波成分(エッジなど)と低周波成分(滑らかな領域)の関係を視覚的に把握しやすくなります。
import numpy as np
import matplotlib.pyplot as plt
from scipy.ndimage import gaussian_filter # 例としてガウシアンフィルタを使用
# 1. 画像の生成 (例: 中央に正方形がある画像)
rows, cols = 64, 64
image = np.zeros((rows, cols), dtype=np.float32)
image[rows//4 : 3*rows//4, cols//4 : 3*cols//4] = 1 # 中央に白い正方形
# 2. 画像のFFT
image_fft = np.fft.fft2(image)
# 3. FFT結果の可視化 (対数スケールで振幅を表示)
# FFT結果は複素数なので、振幅(np.abs)を取り、可視化のために対数スケールに変換する
magnitude_spectrum = np.log(np.abs(image_fft) + 1e-9) # +1e-9 はlog(0)を避けるため
plt.figure(figsize=(15, 6))
# FFT結果 (fftshiftなし)
plt.subplot(1, 3, 1)
plt.imshow(magnitude_spectrum, cmap='gray')
plt.title('FFT結果 (fftshiftなし)')
plt.colorbar()
# 4. fftshiftを適用
image_fft_shifted = np.fft.fftshift(image_fft)
magnitude_spectrum_shifted = np.log(np.abs(image_fft_shifted) + 1e-9)
# FFT結果 (fftshiftあり)
plt.subplot(1, 3, 2)
plt.imshow(magnitude_spectrum_shifted, cmap='gray')
plt.title('FFT結果 (fftshiftあり)')
plt.colorbar()
# 5. フィルタリングの例 (低周波パスフィルタ)
# FFTシフトされたスペクトルをガウシアンフィルタで平滑化(低周波成分を残す)
filtered_fft_shifted = gaussian_filter(image_fft_shifted.real, sigma=5) + 1j * gaussian_filter(image_fft_shifted.imag, sigma=5)
# 6. IFFTのためにfftshiftを元に戻す (ifftshift)
filtered_fft = np.fft.ifftshift(filtered_fft_shifted)
# 7. IFFTを実行して画像に戻す
restored_image = np.fft.ifft2(filtered_fft).real # 虚数部が小さいので実数部のみ
# 8. フィルタリングされた画像の表示
plt.subplot(1, 3, 3)
plt.imshow(restored_image, cmap='gray')
plt.title('フィルタリング後の画像 (低周波パス)')
plt.colorbar()
plt.tight_layout()
plt.show()
解説
- **右の画像 (
フィルタリング後):**
fftshiftでDC成分を中央に移動させた後、中央(低周波成分)を残すようなフィルタを適用し、その後
ifftshiftで元に戻してから
ifft2`を実行しています。結果として、元の画像のシャープなエッジがぼやけ、滑らかな画像になっていることがわかります。これは高周波成分が除去されたためです。 - 中央の画像 (fftshiftあり)
fftshift
を適用することで、DC成分が画像の中央に移動し、画像の周波数構成がより直感的に理解できるようになります。中央が低周波成分、外側が高周波成分を示します。正方形のエッジに対応する十字型のパターンが見られます。 - 左の画像 (fftshiftなし)
FFTスペクトルのDC成分(低周波の中心)が左上隅にあり、そこから放射状に高周波成分が広がっています。
fft.fftshift()
は、axes
引数を使って特定の軸のみをシフトすることもできます。これは、多次元配列の一部の次元のみを操作したい場合に便利です。
import numpy as np
# 3D配列の例
arr_3d = np.arange(2*3*4).reshape((2, 3, 4))
print("元の3D配列:\n", arr_3d)
# デフォルト (全ての軸をシフト)
shifted_all_axes = np.fft.fftshift(arr_3d)
print("\n全ての軸をシフト:\n", shifted_all_axes)
# 0番目の軸 (深さ方向) のみシフト
shifted_axis0 = np.fft.fftshift(arr_3d, axes=0)
print("\n0番目の軸のみシフト:\n", shifted_axis0)
# 1番目の軸 (行方向) のみシフト
shifted_axis1 = np.fft.fftshift(arr_3d, axes=1)
print("\n1番目の軸のみシフト:\n", shifted_axis1)
# 2番目の軸 (列方向) のみシフト
shifted_axis2 = np.fft.fftshift(arr_3d, axes=2)
print("\n2番目の軸のみシフト:\n", shifted_axis2)
# 複数の軸を指定 (例: 0番目と2番目の軸をシフト)
shifted_axes02 = np.fft.fftshift(arr_3d, axes=(0, 2))
print("\n0番目と2番目の軸をシフト:\n", shifted_axes02)
解説
axes
引数にタプルで軸のインデックスを指定することで、必要な軸に対してのみシフト操作を行うことができます。これは、例えば、時空間データ(時間と空間の2次元FFTなど)で、時間軸だけをシフトしたい場合などに役立ちます。
Therefore, providing "alternative methods" would either involve explaining how fftshift
itself works (which is not an "alternative" to using it) or describing entirely different approaches to frequency analysis that do not rely on standard FFT output conventions (which would be outside the scope of "alternatives to fft.fftshift()
").