"""
FFT_func.py - 関数データの高速フーリエ変換（FFT）および逆高速フーリエ変換（IFFT）を実行するスクリプト。
概要:
このスクリプトは、ユーザー定義の関数 `func_str` に基づいてデータを生成し、
そのデータに対してFFTまたはIFFTを実行し、結果をExcelファイルとして保存し、matplotlibで可視化します。
コマンドライン引数でモード、関数文字列、x範囲、データ点数を指定できます。
詳細説明:
- `cparams` オブジェクトを通じて、計算パラメータ（モード、関数文字列、xの最小値、最大値、データ点数など）を管理します。
- `common_func` 関数で、指定された関数文字列とx範囲・点数に基づいて入力データ (x, y) を生成します。
- `FFT` 関数または `IFFT` 関数が選択されたモードに応じて実行され、それぞれ `numpy.fft.fft` または `numpy.fft.ifft` を使用して変換を行います。
- 変換結果はExcelファイルに保存され、元のデータと変換後のデータの両方がグラフとして表示されます。
- `tklib` ライブラリのユーティリティ関数 (`tkApplication`, `tkParams`, `getarg` など) を利用しています。
関連リンク:
:doc:`FFT_func_usage`
"""
import sys
import csv
import numpy as np
from numpy import sin, cos, tan, pi
from numpy import exp, log, sin, cos, tan, arcsin, arccos, arctan, sinh, cosh, tanh, sqrt, abs
from random import random as rand
import pandas as pd
import matplotlib.pyplot as plt
from tklib.tkutils import replace_path, getarg, getintarg, getfloatarg, pint, pfloat
from tklib.tksci.tksci import kB, e
from tklib.tkvariousdata import tkVariousData
from tklib.tkapplication import tkApplication
from tklib.tkparams import tkParams
from tklib.tksci.tksci import asin, acos, atan, degcos, degsin, degtan, degacos, degasin, degatan, arcsin, arccos, arctan
from tklib.tksci.tksci import factorial, log10, gamma, combination, eVTonm, nmToeV, Bn
from tklib.tksci.tksci import h, h_bar, hbar, e, kB, NA, c, pi, pi2, torad, todeg, basee
from tklib.tksci.tksci import me, mp, mn
from tklib.tksci.tksci import u0, e0, a0, R, F, g, G
from tklib.tksci.tksci import HartreeToeV, RyToeV, KToeV, eVToK, JToeV, eVToJ, Debye
app = tkApplication()
cparams = tkParams()
cparams.func_str = 'sin(10.0 * pi * x) + 0.5 * rand()'
cparams.mode = 'FFT'
cparams.xmin = -1.0
cparams.xmax = 1.0
cparams.nx = 1024
cparams.figsize = [12, 8]
cparams.fontsize = 12
cparams.legend_fontsize = 10
logfile = 'fft_func-out.txt' #app.replace_path(cparams.infile)
print(f"Open logfile [{logfile}]")
app.redirect(targets = ["stdout", logfile], mode = 'w')
[ドキュメント]
def usage(app):
"""
スクリプトの利用方法を標準出力に表示します。
概要:
スクリプトのコマンドライン引数の正しい使い方を示します。
:param app: tkApplicationのインスタンス。メッセージ表示に使用されます。
:type app: tklib.tkapplication.tkApplication
"""
print(f"Usage: python {sys.argv[0]} func_str x0 x1 nx")
[ドキュメント]
def update_vars(app):
"""
コマンドライン引数からグローバルパラメータ `cparams` を更新します。
概要:
スクリプト実行時に指定されたコマンドライン引数を解析し、
計算に利用する各種パラメータ (`mode`, `func_str`, `xmin`, `xmax`, `nx`) を更新します。
引数が指定されていない場合は、`cparams` のデフォルト値が使用されます。
:param app: tkApplicationのインスタンス。現在この関数では直接使用されていませんが、引数として渡されます。
:type app: tklib.tkapplication.tkApplication
"""
cparams.mode = getarg(1, cparams.mode)
cparams.func_str = getarg(2, cparams.func_str)
cparams.xmin = getfloatarg(3, cparams.xmin)
cparams.xmax = getfloatarg(4, cparams.xmax)
cparams.nx = getintarg (5, cparams.nx)
[ドキュメント]
def func(x):
"""
グローバルパラメータ `cparams.func_str` に定義された関数を評価します。
概要:
与えられた `x` 値に対して、`cparams.func_str` に文字列として定義された数学関数を計算します。
`eval` 関数を用いて動的に関数を評価します。
:param x: 関数を評価する入力値。
:type x: float or numpy.ndarray
:returns: 入力 `x` に対応する関数の計算結果。
:rtype: float or numpy.ndarray
"""
y = eval(cparams.func_str, globals(), {"x": x})
return y
[ドキュメント]
def common_func():
"""
FFT/IFFT処理の共通前処理として、入力データ (x, y) を生成します。
概要:
`cparams` に設定された範囲 (`xmin`, `xmax`) と点数 (`nx`) に基づいて、
等間隔の `x` データポイントと、それに対応する `func_str` から計算された `y` データポイントを生成します。
:returns: `x` データポイントのリストと `y` データポイントのリストのタプル。
:rtype: tuple[list[float], list[float]]
"""
print("")
print("========================================")
print(f"mode: {cparams.mode}")
print(f"func: {cparams.func_str}")
print(f"xrange : {cparams.xmin} - {cparams.xmax} for {cparams.nx} points")
print("")
print(f"Build function [{cparams.func_str}]")
cparams.xstep = (cparams.xmax - cparams.xmin) / (cparams.nx - 1)
x = [cparams.xmin + i * cparams.xstep for i in range(cparams.nx)]
y = [func(_x) for _x in x]
return x, y
[ドキュメント]
def IFFT():
"""
指定された関数データの逆高速フーリエ変換 (IFFT) を実行し、結果を保存・プロットします。
概要:
`common_func` で生成された入力データ `y` に対して `numpy.fft.ifft` を適用し、
逆フーリエ変換された結果を計算します。変換されたデータと元のデータをExcelファイルに保存し、
matplotlibを使用して両者をグラフで可視化します。
詳細説明:
1. `common_func` を呼び出し、元の `x`, `y` データを取得します。
2. `numpy.fft.ifft` を使用して、`y` のIFFTを計算します。
3. IFFT後の周波数軸 (`xf`) を計算します。
4. 変換結果の実部 (`Fr`) と虚部 (`Fi`) を抽出し、`pd.DataFrame` にまとめます。
5. データフレームを `fft_func-ifft.xlsx` というファイル名でExcelに保存します。
6. `matplotlib.pyplot` を使用して、元のデータ (`x` vs `y`) とIFFT結果 (`xf` vs `Fr`, `Fi`, `Fr`) を2つのサブプロットに描画します。
7. グラフのタイトル、ラベル、凡例を設定し、表示後にアプリケーションを終了します。
"""
x, y = common_func()
F = np.fft.ifft(y)
# calculate FFTed x axis
xfrange = 1.0 / cparams.xstep
dxf = xfrange / (cparams.nx - 1)
xf = [i * dxf for i in range(cparams.nx)]
print("")
print("x range : {} - {}".format(min(x), max(x)))
print("x(fft) range: {} - {}".format(min(xf), max(xf)))
Fr = [F[i].real for i in range(cparams.nx)]
Fi = [F[i].imag for i in range(cparams.nx)]
Fr = [F[i].real for i in range(cparams.nx)]
df = pd.DataFrame(np.array([x, y, xf, Fr, Fi, Fr]).T,
columns = ['x', 'y', 'x(fft)', 'y(fft).r', 'y(fft).i', '|y(fftr)|'])
outfile = 'fft_func-ifft.xlsx' #app.replace_path(cparams.infile, template = ["{dirname}", "{filebody}-ifft.xlsx"])
print("")
print(f"Save to [{outfile}]")
df.to_excel(outfile, index = False)
print("")
print("Plot")
fig, axes = plt.subplots(1, 2, figsize = cparams.figsize)
axes = axes.flatten()
# ax2 = axes[1].twinx()
axes[0].tick_params(labelsize = cparams.fontsize)
axes[1].tick_params(labelsize = cparams.fontsize)
axes[0].plot(x, y, label = 'input data', color = 'black', linewidth = 0.5)
axes[0].set_xlabel('x', fontsize = cparams.fontsize)
axes[0].set_ylabel(f'y = {cparams.func_str}', fontsize = cparams.fontsize)
axes[0].legend(fontsize = cparams.legend_fontsize)
ins1 = axes[1].plot(xf, Fr, label = 'real', color = 'black', linewidth = 0.5)
ins2 = axes[1].plot(xf, Fi, label = 'imag', color = 'red', linewidth = 0.5)
ins3 = axes[1].plot(xf, Fr, label = 'abs', color = 'blue', linewidth = 0.5)
axes[1].set_xlabel('k', fontsize = cparams.fontsize)
axes[1].set_ylabel('Inverse $FT^*$', fontsize = cparams.fontsize)
axes[1].legend(fontsize = cparams.legend_fontsize)
plt.tight_layout()
plt.pause(0.001)
app.terminate("\n", usage = usage, pause = True)
[ドキュメント]
def FFT():
"""
指定された関数データの高速フーリエ変換 (FFT) を実行し、結果を保存・プロットします。
概要:
`common_func` で生成された入力データ `y` に対して `numpy.fft.fft` を適用し、
フーリエ変換された結果を計算します。変換されたデータと元のデータをExcelファイルに保存し、
matplotlibを使用して両者をグラフで可視化します。
詳細説明:
1. `common_func` を呼び出し、元の `x`, `y` データを取得します。
2. `numpy.fft.fft` を使用して、`y` のFFTを計算します。
3. FFT後の周波数軸 (`xf`) を計算します。
4. 変換結果の実部 (`Fr`) と虚部 (`Fi`) を抽出し、`pd.DataFrame` にまとめます。
5. データフレームを `fft_func-fft.xlsx` というファイル名でExcelに保存します。
6. `matplotlib.pyplot` を使用して、元のデータ (`x` vs `y`) とFFT結果 (`xf` vs `Fr`, `Fi`, `Fr`) を2つのサブプロットに描画します。
7. グラフのタイトル、ラベル、凡例を設定し、表示後にアプリケーションを終了します。
"""
x, y = common_func()
F = np.fft.fft(y)
# calculate FFTed x axis
xfrange = 1.0 / cparams.xstep
dxf = xfrange / (cparams.nx - 1)
xf = [i * dxf for i in range(cparams.nx)]
print("")
print("x range : {} - {}".format(min(x), max(x)))
print("x(fft) range: {} - {}".format(min(xf), max(xf)))
Fr = [F[i].real for i in range(cparams.nx)]
Fi = [F[i].imag for i in range(cparams.nx)]
Fr = [F[i].real for i in range(cparams.nx)]
df = pd.DataFrame(np.array([x, y, xf, Fr, Fi, Fr]).T,
columns = ['x', 'y', 'x(fft)', 'y(fft).r', 'y(fft).i', '|y(fftr)|'])
outfile = 'fft_func-fft.xlsx' #app.replace_path(cparams.infile, template = ["{dirname}", "{filebody}-fft.xlsx"])
print("")
print(f"Save to [{outfile}]")
df.to_excel(outfile, index = False)
print("")
print("Plot")
fig, axes = plt.subplots(1, 2, figsize = cparams.figsize)
axes = axes.flatten()
# ax2 = axes[1].twinx()
axes[0].tick_params(labelsize = cparams.fontsize)
axes[1].tick_params(labelsize = cparams.fontsize)
axes[0].plot(x, y, label = 'input data', color = 'black', linewidth = 0.5)
axes[0].set_xlabel('x', fontsize = cparams.fontsize)
axes[0].set_ylabel(f'y = {cparams.func_str}', fontsize = cparams.fontsize)
axes[0].legend(fontsize = cparams.legend_fontsize)
ins1 = axes[1].plot(xf, Fr, label = 'real', color = 'black', linewidth = 0.5)
ins2 = axes[1].plot(xf, Fi, label = 'imag', color = 'red', linewidth = 0.5)
ins3 = axes[1].plot(xf, Fr, label = 'abs', color = 'blue', linewidth = 0.5)
axes[1].set_xlabel('k', fontsize = cparams.fontsize)
axes[1].set_ylabel('$FT^*$', fontsize = cparams.fontsize)
axes[1].legend(fontsize = cparams.legend_fontsize)
plt.tight_layout()
plt.pause(0.001)
app.terminate("\n", usage = usage, pause = True)
[ドキュメント]
def main():
"""
スクリプトのメイン処理を実行します。
概要:
スクリプトの開始点であり、初期設定、コマンドライン引数によるパラメータ更新、
および選択されたモード（FFTまたはIFFT）に応じた処理の呼び出しを行います。
詳細説明:
1. 起動メッセージを表示します。
2. `update_vars` を呼び出してコマンドライン引数からパラメータを更新します。
3. `cparams.print_parameters()` を実行して現在のパラメータ設定を表示します。
4. `cparams.mode` の値に応じて `FFT()` または `IFFT()` 関数を呼び出します。
5. 未実装のモードが指定された場合はエラーメッセージを表示し、スクリプトを終了します。
"""
print("")
print( "#######################################")
print(f"# {sys.argv[0]}")
print( "#######################################")
update_vars(app)
cparams.print_parameters()
if cparams.mode == 'FFT':
FFT()
elif cparams.mode == 'IFFT':
IFFT()
else:
app.terminate(f"\nmode [{cparams.mode}] is not implemented.\n", usage = usage, pause = True)
if __name__ == "__main__":
main()