【名称】 mcfft - multi channnel FFT program 【用法】 mcfft [オプション] [入力ファイル] ( > 出力ファイル ) 【説明】 テキスト形式の時系列データに対して高速フーリエ変換を行う。 トリガ方式の指定と平均化、および周波数伝達関数の算出も可 能。 入力は、タブ・スペース・コンマを区切り文字とするテキスト ファイル。入力ファイルを指定しない場合は標準入力を読む。 入力の１行目はタイトル行で２行目以降に時間刻み一定の時刻 歴の数値データ。入力には少なくとも時間と信号の２列が必要。 列数と行数は(メモリの許す限り)無制限。 出力は後述の出力オプションで指定された項目をタブ区切りテ キストとして標準出力に書き出す。 【オプション】 [入力オプション] -t 列番号 時刻の列番号 デフォルト 1 -i 列番号 伝達関数の入力(基準)信号の列番号 デフォルト 2 -n 行数 最大読み込み可能データ行数 デフォルト 262144 (=2^18) [分析オプション] -b 指数 ブロックサイズ 2^n の n デフォルト 10 (2^10=1024) 但し、以下は例外 0: 冒頭から 2^n として取れる最大点数 1: 末尾にパディングして全体をカバーできる最小の 2^n 点 パディングするのは (0 ではなく) 先頭と末尾の平均値 -a 回数 最大平均回数 デフォルト 0 0 なら元データから取り込める最大回数 -s 文字列 トリガ方法を指定する文字列 デフォルト F,0 文字列内の各要素はスペースを挟まずコンマで区切る。 オーバーラップとプレトリガはブロックサイズに対する % で指定。 -sF,オーバーラップ量 (フリートリガでの連続取り込み) -sT,列番号,値,勾配,プレトリガ量 (トリガ信号による同期取り込み) 勾配は + (増加方向)または - (減少方向)の１文字で指定。 列番号の信号が指定値を指定勾配の方向に横切った時刻から プレトリガ量だけ前の時刻をブロックの開始点とする。 プレトリガ量が負の値の場合はポストトリガとして扱う。 -w 文字 窓関数 デフォルト r r: 矩形 h: 標準ハニング ( 0.5 * (1 + cos(2πt/T)) ) H: ２倍ハニング ( 1 + cos(2πt/T) ) [出力オプション] -d データチェック デフォルト 無し 取り込み開始時刻にはブロックサイズを、それ以外の時刻には 0 を記入した block_head という項目名の列を元データの最後 に加え、時系列データ全体を出力して終了する。 -c 分析条件確認情報出力 デフォルト 無し 周波数(-d 指定の場合は時刻)列の項目名に、ブロックサイズ と平均回数を追記する。(例: Freq[Hz]_N=1024_ave=16 ) -l 周波数 最低出力周波数(l は L の小文字) デフォルト 0.0 -h 周波数 最高出力周波数 デフォルト 分析最大周波数 -o 文字列 出力項目指定文字(列) デフォルト p 出力項目を示す f,s,p,e,x,S,P,E の文字の任意の並び (例 -o xpe) f: フーリエスペクトル (振幅と位相) s: パワー スペクトル PS ( 実効値^2 ) p: パワー スペクトル密度 PSD ( 実効値^2/Δf) e: エネルギスペクトル密度 ESD (T * 実効値^2/Δf) x: 伝達関数 (ゲイン、位相、コヒーレンス) S: パワー スペクトル PS ( 振幅^2 ) P: パワー スペクトル密度 PSD ( 振幅^2/Δf) E: エネルギスペクトル密度 ESD (T * 振幅^2/Δf) 出力の１列目は周波数 Freq[Hz] (-d 指定の場合以外) ２列目以降に(時刻列を除く)元データの全列について列毎に、 出力項目指定文字列の文字の並び順にその文字に対応する分析 結果を出力 【出力データ詳細】 説明のための例として、入力データの列の並びは Time chA chB の３列であるとする。また、周波数 f Hz の成分の時間軸波形 を x(t) = A * cos(2πft+θ) とする。更に、周波数分解能Δf、およびブロック長 T を Δf ( = サンプリング周波数 / ブロックサイズ ) T = 1 / Δf とする。 f: フーリエスペクトル (振幅と位相) 振幅 A, 位相 θ を、それぞれ項目名に _amp, _phase[deg] を付けた列に出力する。 平均はフーリエスペクトルの実部、虚部それぞれのアンサ ンブル平均である。 従って、トリガを用いて着目する現象に同期した平均をと っている場合以外では、一般に平均化した結果には意味が ない。（各回の位相が揃っていないため、振幅は平均回数 が増すにつれて０に近づく。） s: パワー スペクトル PS ( 実効値^2 ) p: パワー スペクトル密度 PSD ( 実効値^2/Δf) e: エネルギスペクトル密度 ESD (T * 実効値^2/Δf) パワスペクトルまたはパワスペクトル密度またはエネルギ スペクトル密度を、項目名に _PS または _PSD または _ESD を付けて出力する。 パワーは実効値の２乗、すなわち A^2/2 で定義しており、 PS : (A^2 / 2) PSD : (A^2 / 2) / Δf ESD : T * (A^2 / 2) / Δf である。 なお、ここで使っている定義の片側パワスペクトルは Parseval の等式を満たしている。 x: 伝達関数 (ゲイン、位相、コヒーレンス) 基準信号を chA とした場合、例えば chB 列に相当する分 析結果として周波数伝達関数のゲイン、位相、およびコヒー レンスの３列を以下のような項目名で出力する。 chB/chA_gain chaB/chaA_phase[deg] chB/chA_coh なお、出力項目の列の数を数えやすくするために、基準信 号 chA 自身に対する伝達関数（ゲイン=1, 位相=0, コヒー レンス=1)も同じ形式で３列で出力する。 内部的にはクロススペクトルのアンサンブル平均を入力の パワスペクトルで割って求めている。従って、出力に重畳 しているノイズの影響を平均化によって低減することが期 待できる。 S: パワー スペクトル PS ( 振幅^2 ) P: パワー スペクトル密度 PSD ( 振幅^2/Δf) E: エネルギスペクトル密度 ESD (T * 振幅^2/Δf) 項目名に _PS(amp^2) または _PSD(amp^2) または _ESD(amp^2) を付けて出力する。 パワーとして実効値の２乗(= A^2/2)ではなく単純に振幅の ２乗(= A^2)を使っていること以外は上述の s, p, e 指定 の場合と同じである。 言い換えると、数値としては s, p, e 指定の場合に比べて ２倍の値になる。(これを dB 表示すると 3 dBの違いになる。) なお、この定義のパワーは Parseval の等式を満たさない。 【用例】 (1) mcfft infile.csv > outfile.txt １列目が時刻データである infile.csv を全てデフォルト 条件で分析し、各列のパワスペクトル密度を outfile.txt に保存する。 (2) mcfft -sT,2,5,+,20 -b9 -a50 -i3 -oxp infile.csv > outfile.txt トリガ信号を 2 列目として、それが値 5 を増加方向に横 切る時刻を基準に、その時刻よりブロックサイズ長の 20 % 前の時刻からデータを取り込む。 ブロックサイズは 512 (=2^9) で、平均化は 50 回以内と する。 出力は 3 列目を基準信号とする周波数伝達関数とパワスペ クトル密度とする。 (3) mcfft -sT,2,5,+,20 -b9 -a50 -i3 -d -c infile.csv > outfile.txt 上記の (2) と同じ分析条件の指定で、時系列データのどの 部分が実際に取り込まれて平均化に使われるかを確認する ための block_head という項目名の列を追加した時系列デー タを出力する。 時刻列(この例ではデフォルトで 1 列目)の項目名の後ろに ブロックサイズと平均化回数を付け足す。 【作者】 尾川敏也 t_ogw[at]nifty.ne.jp : [at] を @ に書き替えて下さい。 参考:2020年6月以前のアドレスは