instfreq

瞬時周波数の推定

構文

ifq = instfreq(x,fs)

ifq = instfreq(x,t)

ifq = instfreq(xt)

ifq = instfreq(tfd,fd,td)

ifq = instfreq(___,Name=Value)

[ifq,t] = instfreq(___)

instfreq(___)

説明

ifq = instfreq(x,fs) は、fs のレートでサンプリングした信号の瞬時周波数、x を推定します。x が行列の場合、この関数は各列の瞬時周波数を個別に推定し、結果を対応する列 ifq で返します。

ifq = instfreq(x,t) は、t に格納されている時間値でサンプリングした x の瞬時周波数を推定します。

ifq = instfreq(xt) は、MATLAB^® timetable xt に格納されている信号の瞬時周波数を推定します。関数は、timetable 内のすべての変数および各変数内のすべての列を個別に処理します。

ifq = instfreq(tfd,fd,td) は、その時間-周波数分布、tfd が fd に格納されている周波数値および td に格納されている時間値でサンプリングされる信号の瞬時周波数を推定します。

ifq = instfreq(___,Name=Value) は、名前と値の引数を使用して、上記の任意の構文に追加オプションを指定します。瞬時周波数の推定に使用するアルゴリズム、または計算に使用する周波数範囲を指定できます。

また [ifq,t] = instfreq(___) は、ifq に対応するサンプル時間のベクトル t も返します。

出力引数なしで instfreq(___) を使用すると、推定された瞬時周波数がプロットされます。

例

非定常信号の瞬時周波数

5 kHz で 4 秒間サンプリングされた信号を生成します。信号は、振動する振幅の領域および増加トレンドで変動する周波数の領域によって分離された、持続時間が減少していく一連のパルスで構成されています。信号をプロットします。

fs = 5000;
t = 0:1/fs:4-1/fs;

s = besselj(0,1000*(sin(2*pi*t.^2/8).^4));

% To hear, type sound(s,fs)

plot(t,s)

信号の時間依存周波数をパワー スペクトログラムの最初のモーメントとして推定します。パワー スペクトログラムをプロットして瞬時周波数を重ね合わせます。

instfreq(s,fs)

複素数値の信号の瞬時周波数

正弦関数的に変化する周波数成分をもつチャープからなる複素数値の信号を生成します。信号は 3 kHz で 1 秒間サンプリングされ、ホワイト ガウス ノイズに組み込まれます。

fs = 3000;
t = 0:1/fs:1-1/fs;
x = exp(2j*pi*100*cos(2*pi*2*t))+randn(size(t))/100;

信号の時間依存周波数をパワー スペクトログラムの最初のモーメントとして推定します。これは、instfreq が複素数値の信号をサポートする唯一のメソッドです。パワー スペクトログラムをプロットして瞬時周波数を重ね合わせます。

instfreq(x,t)

マルチチャネル信号の瞬時周波数

2 つの電圧制御発振器で構成される、1 kHz で 2 秒間サンプリングされた 2 チャネルの信号を生成します。

• 1 つのチャネルでは、周期の 75% で最大になるノコギリ波として、瞬時周波数が時間と共に変化します。

• もう 1 つのチャネルでは、デューティ比 30% の方形波として、瞬時周波数が時間と共に変化します。

2 つのチャネルのスペクトログラムをプロットします。ノコギリ チャネルには 0.1 秒の時間分解能を指定し、方形チャネルには 10 Hz の周波数分解能を指定します。

fs = 1000;
t = (0:1/fs:2)';
x = vco(sawtooth(2*pi*t,0.75),[0.1 0.4]*fs,fs);
y = vco(square(2*pi*t,30),[0.1 0.3]*fs,fs);

subplot(1,2,1)
pspectrum(x,fs,'spectrogram','TimeResolution',0.1)
subplot(1,2,2)
pspectrum(y,fs,'spectrogram','FrequencyResolution',10)

信号を timetable に保存します。瞬時周波数を計算して表示します。

xt = timetable(seconds(t),x,y);

clf
instfreq(xt)

解析信号を使用して計算を繰り返します。

instfreq(xt,'Method','hilbert')

チャープの瞬時周波数

ガウス変調された 2 次チャープを生成します。2 kHz のサンプルレートと 4 秒の信号持続時間を指定します。

fs = 2000;
t = 0:1/fs:4-1/fs;

q = chirp(t-1,0,1/2,20,'quadratic',100,'convex').*exp(-1.7*(t-2).^2);
plot(t,q)

既定の設定で関数 pspectrum を使用して信号のパワー スペクトルを推定します。推定を使用して瞬時周波数を計算します。

[p,f,t] = pspectrum(q,fs,'spectrogram');

instfreq(p,f,t)

シンクロスクイーズド フーリエ変換を使用して計算を繰り返します。500 サンプルのハン ウィンドウを使用して信号をセグメントに分割し、ウィンドウを適用します。

[s,sf,st] = fsst(q,fs,hann(500));

instfreq(abs(s).^2,sf,st)

2 つの異なる方法を使用して検出された瞬時周波数を比較します。

[psf,pst] = instfreq(p,f,t);
[fsf,fst] = instfreq(abs(s).^2,sf,st);

plot(fst,fsf,pst,psf)

正弦波の瞬時周波数

1 kHz で 0.3 秒間サンプリングされ、分散 1/16 のホワイト ガウス ノイズに組み込まれた正弦波信号を生成します。200 Hz の正弦波周波数を指定します。信号の瞬時周波数を推定して表示します。

fs = 1000;
t = (0:1/fs:0.3-1/fs)';

x = sin(2*pi*200*t) + randn(size(t))/4;

instfreq(x,t)

もう一度信号の瞬時周波数を推定しますが、ここでは入力として 25 Hz のおおまかな周波数分解能の時間-周波数分布を使用します。

[p,fd,td] = pspectrum(x,t,'spectrogram','FrequencyResolution',25);

instfreq(p,fd,td)

条件スペクトル モーメントとしての瞬時周波数と帯域幅

300 Hz と 1,200 Hz の間で周波数が正弦関数的に変化するチャープで構成される信号を生成します。信号は 3 kHz で 2 秒間サンプリングされます。

fs = 3e3;
t = 0:1/fs:2;
y = chirp(t,100,1,200,"quadratic");
y = vco(cos(2*pi*t),[0.1 0.4]*fs,fs);

instfreq を使用して、信号の瞬時周波数と対応するサンプル時間を計算します。tfsmoment (Predictive Maintenance Toolbox)で計算された信号の時間-周波数分布の非集中 1 次条件スペクトル モーメントに出力が対応していることを確認します。

[z,tz] = instfreq(y,fs);
[a,ta] = tfsmoment(y,fs,1,Centralize=false);

plot(tz,z,ta,a,'.')
legend("instfreq","tfsmoment")

instbw を使用して、信号の瞬時帯域幅と対応するサンプル時間を計算します。スケール係数を 1 に指定します。信号の時間分布の集中 2 次条件スペクトル モーメントの平方根に出力が対応していることを確認します。つまり、instbw は標準偏差を生成し、tfsmoment は分散を生成します。

[w,tw] = instbw(y,fs,ScaleFactor=1);
[m,tm] = tfsmoment(y,fs,2);

plot(tw,w,tm,sqrt(m),'.')
legend("instfreq","tfsmoment")

入力引数

x — 入力信号

ベクトル | 行列

ベクトルまたは行列として指定される入力信号。x がベクトルの場合、instfreq では単一チャネルとして取り扱います。x が行列の場合、この関数は各列の瞬時周波数を個別に計算し、結果を対応する列 ifq で返します。

例: sin(2*pi*(0:127)/16)+randn(1,128)/100 は、ノイズを含む正弦波を指定します

例: [2 1].*sin(2*pi*(0:127)'./[16 64]) は、2 チャネルの正弦波を指定します。

データ型: single | double
複素数のサポート: あり

fs — サンプルレート

正のスカラー

サンプルレート。正のスカラーで指定します。サンプルレートは単位時間あたりのサンプル数です。時間の単位が秒の場合、サンプルレートの単位は Hz です。

データ型: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

t — サンプル時間

実数ベクトル | duration スカラー | duration 配列 | datetime 配列

サンプル時間。実数ベクトル、duration スカラー、duration 配列、または datetime 配列で指定します。

• duration スカラー — x の連続するサンプル間の時間間隔。

• 実数ベクトル、duration 配列、または datetime 配列 — x の各要素に対応する時点。

例: seconds(1) は信号の連続測定値間に 1 秒の中断を指定します。

例: seconds(0:8) は信号が 1 Hz で 8 秒間サンプリングされることを指定します。

データ型: single | double | duration | datetime

xt — 入力 timetable

timetable

入力 timetable。xt に増加する有限の行時間を含めなければなりません。

timetable が欠損している場合や時間点が重複している場合、欠損または重複する時間および非等間隔の時間をもつ timetable の整理のヒントを使用して修正できます。

例: timetable(seconds(0:4)',randn(5,1)) は 1 Hz で 4 秒間サンプリングされたランダム過程を指定します。

例: timetable(seconds(0:4)',randn(5,3),randn(5,4)) にはいずれも 1 Hz で 4 秒間サンプリングされた、3 チャネルのランダム過程と 4 チャネルのランダム過程が含まれます。

データ型: single | double
複素数のサポート: あり

tfd — 時間-周波数分布

行列

時間-周波数分布。fd に格納されている周波数および td に格納されている時間値でサンプリングされる行列として指定します。この入力引数は、Method を "tfmoment" に設定した場合にのみサポートされます。

例: [p,f,t] = pspectrum(sin(2*pi*(0:511)/4),4,"spectrogram") は、4 Hz で 128 秒間サンプリングされた 1 Hz の正弦波の時間-周波数分布と、それを計算する周波数および時間も指定します。

データ型: single | double

fd, td — 時間-周波数分布の周波数および時間値

ベクトル

時間-周波数分布の周波数および時間値。ベクトルで指定します。これらの入力引数は、Method を "tfmoment" に設定した場合にのみサポートされます。

例: [p,f,t] = pspectrum(sin(2*pi*(0:511)/4),4,"spectrogram") は、4 Hz で 128 秒間サンプリングされた 1 Hz の正弦波の時間-周波数分布と、それを計算する周波数および時間も指定します。

データ型: single | double

名前と値の引数

オプションの引数のペアを Name1=Value1,...,NameN=ValueN として指定します。ここで、Name は引数名で、Value は対応する値です。名前と値の引数は他の引数の後に指定しなければなりませんが、ペアの順序は重要ではありません。

R2021a より前は、名前と値をそれぞれコンマを使って区切り、Name を引用符で囲みました。

例: 'Method','tfmoment','FrequencyLimits',[25 50] は、時間-周波数分布の最初の条件スペクトル モーメントを検出することによって、25 Hz から 50 Hz の範囲で入力の瞬時周波数を計算します。

FrequencyLimits — 周波数範囲

[0 fs/2] (実数値の信号の既定値) | [-fs/2 fs/2] (複素数値の信号の既定値) | 2 要素ベクトル (Hz 単位)

周波数範囲。Hz 単位の 2 要素ベクトルとして指定します。FrequencyLimits を指定していない場合、この引数の既定値は実数値の信号では [0 fs/2]、複素数値の信号では [-fs/2 fs/2] になります。この引数は、Method を "tfmoment" に設定した場合にのみサポートされます。

データ型: single | double

Method — 計算法

"tfmoment" (既定値) | "hilbert"

計算法。"tfmoment" または "hilbert" として指定します。

• "tfmoment" — x の時間-周波数分布の最初の条件スペクトル モーメントとして瞬時周波数を計算します。x が不等間隔サンプルの場合、関数は信号を等間隔グリッドに内挿して、瞬時周波数を計算します。

• "hilbert" — ヒルベルト変換を使用して検出された x の解析信号の位相の微分として瞬時周波数を計算します。このメソッドは等間隔サンプルの実数値の信号のみを受け入れ、時間-周波数分布の入力はサポートしません。

出力引数

ifq — 瞬時周波数

ベクトル | 行列 | timetable

瞬時周波数。入力と同じ次元のベクトル、行列または timetable として返されます。

データ型: single | double

t — 周波数推定値の時間

実数ベクトル | duration 配列 | datetime 配列

周波数推定値の時間。実数ベクトル、duration 配列、または datetime 配列として返されます。

データ型: single | double

メモ

t、fd、または td が単精度で、入力ベクトルまたは行列が倍精度の場合、関数は出力を単精度で返します。出力データ型はサンプル レートのデータ型に依存しません。

詳細

瞬時周波数

非定常信号の "瞬時周波数" は、信号が徐々に変化する際、その中にある周波数の平均に関連する時変パラメーターです [1]、[2]。

• 'Method' が 'tfmoment' に設定されている場合、instfreq は入力信号の時間-周波数分布の最初の条件スペクトル モーメントとして瞬時周波数を推定します。関数は以下を実行します。

  1. 関数 pspectrum を使用して入力のスペクトログラム パワー スペクトル P(t,f) を計算し、このスペクトルを時間-周波数分布として使用します。

  2. 以下を使用して瞬時周波数を推定します。

     $$f_{\text{inst}}(t)=\frac{{\displaystyle {\int }_0^{\infty }fP(t,f)df}}{{\displaystyle {\int }_0^{\infty }P(t,f)df}}.$$

• 'Method' が 'hilbert' に設定されている場合、instfreq は入力の解析信号の位相の微分として瞬時周波数を推定します。関数は以下を実行します。

  1. 関数 hilbert を使用して入力の解析信号 x_A を計算します。

  2. 以下を使用して瞬時周波数を推定します。

     $$f_{\text{inst}}(t)=\frac{1}{2\pi }\frac{d\varphi }{dt},$$

     ここで ϕ は入力の解析信号の位相です。

参照

[1] Boashash, Boualem. “Estimating and Interpreting the Instantaneous Frequency of a Signal. I. Fundamentals.” Proceedings of the IEEE^® 80, no. 4 (April 1992): 520–538. https://doi.org/10.1109/5.135376.

[2] Boashash, Boualem. "Estimating and Interpreting The Instantaneous Frequency of a Signal. II. Algorithms and Applications." Proceedings of the IEEE 80, no. 4 (May 1992): 540–568. https://doi.org/10.1109/5.135378.

拡張機能

C/C++ コード生成

MATLAB® Coder™ を使用して C および C++ コードを生成します。

使用上の注意および制限:

• 名前と値の引数はコンパイル時の定数でなければなりません。

• コード生成では、timetable はサポートされていません。

GPU 配列

Parallel Computing Toolbox™ を使用してグラフィックス処理装置 (GPU) 上で実行することにより、コードを高速化します。

この関数は、GPU 配列を完全にサポートします。詳細については、GPU での MATLAB 関数の実行 (Parallel Computing Toolbox)を参照してください。

バージョン履歴

R2018a で導入

R2022b: gpuArray オブジェクトの使用

関数 instfreq は gpuArray オブジェクトをサポートしています。この機能を使用するには、Parallel Computing Toolbox™ が必要です。

参考

hilbert | instbw | pspectrum | tfmoment (Predictive Maintenance Toolbox) | tfsmoment (Predictive Maintenance Toolbox) | tftmoment (Predictive Maintenance Toolbox)

トピック

• ヒルベルト変換および瞬時周波数
• 複素チャープの瞬時周波数

×

例を開く

この例の変更されたバージョンがあります。編集された方の例を開きますか?

いいえ、変更されたバージョンを上書きします。 はい