filtfilt

ゼロ位相デジタル フィルター処理

構文

y = filtfilt(b,a,x)

y = filtfilt(sos,g,x)

y = filtfilt(d,x)

y = filtfilt(B,A,x,"ctf")

y = filtfilt({B,A,g},x)

y = filtfilt(___,Name=Value)

説明

y = filtfilt(b,a,x) は、入力データ x を順方向と逆方向の両方で処理することにより、ゼロ位相デジタル フィルター処理を実行します。データを順方向でフィルター処理した後、関数は初期条件を一致させて起動時と終了時の過渡現象を最小限に抑え、フィルター処理したシーケンスを反転し、その反転したシーケンスを逆方向にフィルター処理します。結果には、以下の特性があります。

• ゼロ位相の歪み

• 元のフィルター伝達関数の振幅の 2 乗と等しい、フィルター伝達関数

• b と a により指定されるフィルターの次数の倍のフィルター次数

微分器フィルターやヒルベルト FIR フィルターの演算は位相応答に大きく依存しているため、関数 filtfilt をこれらのフィルターと共に使うことは避けてください。

y = filtfilt(sos,g,x) ゼロ位相は、行列 sos およびスケール値 g で表される 2 次セクション (双二次) フィルターを使用して入力データ x をフィルター処理します。

y = filtfilt(d,x) ゼロ位相は、デジタル フィルター d を使用して入力データ x をフィルター処理します。d を周波数応答仕様に基づいて生成するには、関数 designfilt を使用します。

y = filtfilt(B,A,x,"ctf") は、分子係数 B と分母係数 A によって定義されるCascaded Transfer Functions (CTF) を使用し、入力データ x をゼロ位相でフィルター処理します。 (R2024b 以降)

メモ

A をスカラーまたはベクトルとして指定する場合、B が 2 次セクション行列入力 sos から取得した 6 列の CTF 分子行列であることを明確にするため、"ctf" オプションを指定します。

y = filtfilt({B,A,g},x) は、入力データ x のフィルター処理にフィルターのスカラー値 g を含めます。 (R2024b 以降)

y = filtfilt(___,Name=Value) は、名前と値の引数を使用して追加オプションを指定します。 (R2026a 以降)

例

心電図波形のゼロ位相フィルター処理

ゼロ位相フィルター処理は、フィルター処理された時間波形の特徴を、フィルター処理されていない信号で現れる場所にそのまま保持する上で役立ちます。

合成心電図 (ECG) 波形をゼロ位相フィルター処理します。波形を生成する関数は、例の最後にあります。QRS 群は心電図波形の重要な特徴です。この例では時点 160 あたりから始まっています。

wform = ecg(500);

plot(wform)
axis([0 500 -1.25 1.25])
text(155,-0.4,"Q")
text(180,1.1,"R")
text(205,-1,"S")

心電図を加法性ノイズで乱します。再現可能な結果が必要な場合は、乱数発生器をリセットします。ローパス FIR 等リップル フィルターを作成し、ゼロ位相のフィルターと通常のフィルターの両方を使用してノイズの多い波形をフィルター処理します。

rng("default")

x = wform' + 0.25*randn(500,1);
d = designfilt("lowpassfir", ...
    PassbandFrequency=0.15,StopbandFrequency=0.2, ...
    PassbandRipple=1,StopbandAttenuation=60, ...
    DesignMethod="equiripple");
y = filtfilt(d,x);
y1 = filter(d,x);

tiledlayout("flow")
nexttile
plot([y y1])
title("Filtered Waveforms")
legend(["Zero-phase Filtering" "Conventional Filtering"])

nexttile
plot(wform)
title("Original Waveform")

ゼロ位相フィルター処理では信号のノイズが低減し、QRS 群は元々の発生と同じ時点に保持されます。通常のフィルター処理では信号のノイズは低減しますが、QRS 群に遅延が生じます。

バタワース 2 次セクション フィルターを使用してフィルター処理を繰り返します。

d1 = designfilt("lowpassiir",FilterOrder=12, ...
    HalfPowerFrequency=0.15,DesignMethod="butter");
y = filtfilt(d1,x);

figure
plot(x)
hold on
plot(y,LineWidth=3)
hold off
legend(["Noisy ECG" "Zero-Phase Filtering"])

この関数は ECG 波形を生成します。

function x = ecg(L)
%ECG Electrocardiogram (ECG) signal generator.
%   ECG(L) generates a piecewise linear ECG signal of length L.
%
%   EXAMPLE:
%   x = ecg(500).';
%   y = sgolayfilt(x,0,3); % Typical values are: d=0 and F=3,5,9, etc. 
%   y5 = sgolayfilt(x,0,5); 
%   y15 = sgolayfilt(x,0,15); 
%   plot(1:length(x),[x y y5 y15]);

%   Copyright 1988-2002 The MathWorks, Inc.

a0 = [0,1,40,1,0,-34,118,-99,0,2,21,2,0,0,0]; % Template
d0 = [0,27,59,91,131,141,163,185,195,275,307,339,357,390,440];
a = a0/max(a0);
d = round(d0*L/d0(15)); % Scale them to fit in length L
d(15)=L;

for i=1:14
       m = d(i):d(i+1)-1;
       slope = (a(i+1)-a(i))/(d(i+1)-d(i));
       x(m+1) = a(i)+slope*(m-d(i));
end

end

線形チャープのゼロ位相カスケード型フィルター処理

R2024b 以降

カスケード型伝達関数を使用して、ゼロ位相フィルター処理を実行します。

通過帯域リップルと阻止帯域の減衰量がそれぞれ 0.1 dB と 40 dB である楕円フィルターを設計します。サンプル レート 2000 Hz を指定します。フィルターの周波数応答をプロットします。

Fs = 2000;
[B,A] = ellip(20,0.1,40,[0.3 0.7],"ctf");
freqz(B,A,2048,Fs,"ctf")

1 秒でナイキスト周波数に達する線形スイープ チャープ信号をフィルター処理します。入力信号と出力信号のスペクトルを比較します。

t = 0:1/Fs:1;
x = chirp(t,0,t(end),Fs/2)';
y = filtfilt(B,A,x,"ctf");
pspectrum([x y],Fs,Leakage=1,FrequencyResolution=1)

振動信号からの正弦波アーティファクトのフィルター処理

R2024b 以降

伝達関数によるゼロ位相フィルター処理を使用して、ノイズの多い正弦波アーティファクトを再作成します。CTF とスケール値を使用して、振動信号をフィルター処理します。

正規分布ノイズと 3 つの正弦波形で構成される信号 u を作成します。サンプル レートは 1 kHz です。

rng("default")
Fs = 1e3;
ts = (0:1/Fs:2)';

a0 = [3 2 1];
f0 = [0.1 0.5 0.9]*Fs/2;
p0 = [0 pi/4 pi/2];

u = 0.1*randn(size(ts)) + 0.1*sin(2*pi*f0.*ts+p0)*a0';

3 次のバタワース バンドストップ デジタル フィルターで n0 をフィルター処理して信号 v を作成し、ノイズの多い正弦波アーティファクトを再作成します。

[b,a] = butter(3,[0.15 0.85],"stop");
v = filtfilt(b,a,u);

u と v を比較します。両方の信号が同位相であることを確認します。

tiledlayout("flow")
nexttile
strips([u(ts<0.1) v(ts<0.1)],0.1,Fs)
legend(["u" "v"],Location="southeast")
xlabel("Time (seconds)")
nexttile
pspectrum([u v],Fs)
legend(["u" "v"],Location="southeast")

電圧制御された振動信号 x を作成します。信号 v で表されるノイズの多い正弦波アーティファクトを追加します。

vo = exp(-2*abs(ts-1)).*sin(8*pi*ts);
x = vco(vo,[0.25 0.75]*Fs/2,Fs) + v;

信号 x を 24 次のチェビシェフ II 型フィルターでフィルター処理します。CTF 形式とスケール値 (B,A,g) を使用します。

[B,A,g] = cheby2(24,50,[0.2 0.8],"ctf");
y = filtfilt({B,A,g},x);

短時間フーリエ変換の振幅の 2 乗を比較します。阻止帯域で振幅が急激に減少するのを観察します。

tiledlayout("flow")
nexttile
stft(x,Fs,Window=bohmanwin(128),OverlapLength=120, ...
    FFTLength=512,FrequencyRange="onesided")
title("Input x")
nexttile
stft(y,Fs,Window=bohmanwin(128),OverlapLength=120, ...
    FFTLength=512,FrequencyRange="onesided")
title("Output y")

複数の高調波を含む正弦波トーンのゼロ位相フィルター処理

R2026a 以降

600 Hz まで 60 Hz 間隔で複数の高調波を含む複素数値の 30 Hz 正弦波トーンに対して、ゼロ位相フィルター処理を適用します。

サンプル レート 600 Hz で 2 秒間の複素数値信号を生成します。この信号は、10 V、30 Hz の正弦波トーンで構成され、60 Hz から 600 Hz まで等間隔に分布する 10 個の高調波を含み、各高調波の振幅は 1 V です。

Fs = 600;
t = (0:1/Fs:2)';
x = 10*exp(1i*2*pi*30*t) + sum(exp(1i*2*pi*60*t*(0:9)),2);

トーンは 30 Hz で振動し、不要成分は 60 Hz 間隔でその 10 倍まで周波数のピークをもつため、対象信号であるトーンを復元するには、10 次の IIR 櫛形ノッチ フィルターが適しています。

b と a を、品質係数 35 の 10 次 IIR 櫛形ノッチ フィルターの分子係数および分母係数を表すベクトルとして定義します。特定の次数および品質係数をもつ IIR 櫛形フィルターの設計の詳細については、iircomb (DSP System Toolbox)を参照してください。

b = [0.957 zeros(1,9) -0.957];
a = [1 zeros(1,9) -0.914];

信号の実数部を 0.9 秒から 1.1 秒までの範囲でプロットします。フィルター応答を 0 Hz から 600 Hz までの範囲でプロットします。

tiledlayout("flow")
nexttile
plot(t,real(x))
xlim([0.9 1.1])
xlabel("Time (s)")
title("real(x)")
nexttile
[h,f] = freqz(b,a,8192,"whole",Fs);
plot(f,mag2db(abs(h)))
xlabel("Frequency (Hz)")
title("Magnitude Response (b,a)")

位相を保持したまま入力信号をフィルター処理します。Gustafsson の手法を使用して、フィルター状態の初期条件を推定します。過渡長を信号長に等しいと仮定し、入力信号はそのサンプルをミラーリングしてパディングします。

y = filtfilt(b,a,x, ...
    InitialStatesMethod="gustafsson", ...
    TransientLength=length(x), ...
    PaddingPattern="reflect");

入力信号とフィルター処理後の信号を、時間領域と周波数領域で比較します。両方の信号の実数部を 0.9 秒から 1.1 秒までの範囲で、ウェルチのパワー スペクトルを 0 Hz から 600 Hz までの範囲でプロットします。フィルター処理後の信号は高調波が除去された 30 Hz のトーンを示し、フィルター処理後の信号 y は入力信号 x と同相です。

figure;
tiledlayout("vertical")
nexttile
plot(t,real(x),t,real(y),"--")
legend("real(" + ["x" "y"] +")",Location="southeast")
xlabel("Time (s)")
xlim([0.9 1.1])
nexttile
[p,f] = pwelch([x y],[],[],8192,Fs);
pdb = pow2db(abs(p));
plot(f,pdb(:,1),"-",f,pdb(:,2),"--")
xlabel("Frequency (Hz)")
legend("PSD of " + ["x" "y"],Location="southeast")

入力引数

b, a — 伝達関数の係数

ベクトル

伝達関数の係数。ベクトルで指定します。全極フィルターを使用する場合は、b を 1 に設定します。全零 (FIR) フィルターを使用する場合は、a を 1 に設定します。

例: b = [1 3 3 1]/6 および a = [3 0 1 0]/3 は、正規化された 3 dB の周波数 0.5π ラジアン/サンプルをもつ 3 次のバタワース フィルターを指定します。

データ型: double | single

x — 入力信号

ベクトル | 行列 | N 次元配列

入力信号。実数値または複素数値のベクトル、行列、または N 次元配列として指定します。x は有限の値でなければなりません。x の長さは、max(length(B)-1, length(A)-1) として定義されるフィルター次数の 3 倍より大きくなければなりません。関数は、x が行ベクトルでなければ、x の最初の配列次元に沿って動作します。x が行ベクトルである場合は、2 番目の次元に沿って動作します。

例: cos(pi/4*(0:159))+randn(1,160) は単一チャネルの行ベクトル信号です。

例: cos(pi./[4;2]*(0:159))'+randn(160,2) は 2 チャネル信号です。

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

sos — 2 次セクションの係数

行列

2 次セクションの係数。行列として指定します。sos は L 行 6 列の行列で、セクション数 L が 2 以上でなければなりません。セクション数が 2 未満の場合、関数は入力を分子ベクトルとして扱います。sos の各行は 2 次 (双二次) フィルターの係数に対応しています。sos の i 行目は [bi(1) bi(2) bi(3) ai(1) ai(2) ai(3)] に対応しています。

例: s = [2 4 2 6 0 2;3 3 0 6 0 0] は、正規化された 3 dB の周波数 0.5π ラジアン/サンプルをもつ 3 次のバタワース フィルターを指定します。

データ型: double | single

g — スケール値

スカラー | ベクトル

スケール値。実数値のスカラー、または L + 1 個の要素をもつ実数値のベクトルとして指定します。ここで、L は CTF セクションの数です。スケール値は、カスケード フィルター表現のセクション全体にわたるフィルター ゲインの分布を表します。

filtfilt 関数は、g の指定方法に応じて、scaleFilterSections 関数を使用してフィルター セクションにゲインを適用します。

• スカラー — この関数は、すべてのフィルター セクションにわたってゲインを均一に配分します。

• ベクトル — この関数は、最初の L 個のゲイン値を対応するフィルター セクションに適用し、最後のゲイン値をすべてのフィルター セクションに均一に配分します。

データ型: double | single

d — デジタル フィルター

digitalFilter オブジェクト

デジタル フィルター。digitalFilter オブジェクトで指定します。デジタル フィルターを周波数応答仕様に基づいて生成するには、関数 designfilt を使用します。

例: d = designfilt("lowpassiir",FilterOrder=3,HalfPowerFrequency=0.5) は、正規化された 3 dB の周波数 0.5π ラジアン/サンプルをもつ 3 次のバタワース フィルターを指定します。

B,A — カスケード伝達関数 (CTF) 係数

スカラー | ベクトル | 行列

カスケード伝達関数 (CTF) 係数。スカラー、ベクトル、または行列として指定します。B と A は、それぞれカスケード伝達関数の分子係数と分母係数をリストします。

B のサイズは L 行 (m + 1) 列、A のサイズは L 行 (n + 1) でなければなりません。ここで、

• L はフィルター セクションの数を表します。

• m はフィルターの分子の次数を表します。

• n はフィルターの分母の次数を表します。

カスケード伝達関数の形式と係数行列の詳細については、CTF 形式によるデジタル フィルターの指定を参照してください。

メモ

A(:,1) のいずれかの要素が 1 と等しくない場合、filtfilt はフィルター係数を A(:,1) で正規化します。この場合、A(:,1) はゼロ以外でなければなりません。

例: B = [2 4 2;3 3 0] と A = [6 0 2;6 0 0] は、正規化された 3 dB の周波数 0.5π ラジアン/サンプルをもつ 3 次のバタワース フィルターを指定します。

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

名前と値の引数

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

例: y = filtfilt(d,x,InitialStatesMethod="gustafsson") は、デジタル オブジェクト d を使用して入力信号 x をゼロ位相フィルター処理し、最適な初期状態を計算するために Gustafsson の順方向-逆方向最小二乗手法を適用します。

InitialStatesMethod — フィルター状態の初期条件の推定方法

"const" (既定値) | "gustafsson" | "zero"

R2026a 以降

フィルター状態の初期条件 ("初期状態") の推定方法。次のいずれかとして指定します。

• "const" — filtfilt は、フィルター係数から導出される線形システムを解くことで初期状態を推定します。

  • この方法は、入力信号が一定であり、システムが定常状態で動作していることを前提としています。したがって、フィルター処理後の信号は、定数入力信号に対応する定常値から開始します。

  • x が定数でない場合、関数はその最初のサンプルを使用し入力信号が一定であると仮定して、この方法を使用します。

• "gustafsson" — filtfilt は、Gustafsson のフォワード-バックワード最小二乗手法を使用して初期状態を推定します[1]。この方法は、フォワード-バックワードおよびバックワード-フォワードにフィルター処理した出力間の差分のエネルギーを最小化し、対称性を達成するとともに両端での境界アーティファクトを低減します。

• "zero" — filtfilt はすべての初期状態をゼロに設定します。この方法は計算量が最も少なく、入力信号やフィルター応答に依存しません。

TransientLength — 過渡応答の長さ

"filtord" (既定値) | "stepzlen" | 非負の整数

R2026a 以降

過渡応答の長さ ("過渡長")。"filtord"、"stepzlen"、または length(x) 以下の非負の整数として指定します。

過渡長は、filtfilt がフィルター処理後の信号の開始時および終了時における過渡応答のモデル化および抑制に使用するサンプル数を表します。

• "filtord" — filtfilt は、フィルター次数の 3 倍の過渡長を使用します。

• "stepzlen" — filtfilt は、フィルターのステップ応答の整定時間に基づく過渡長を使用します。

  ステップ応答が s[k]、定常値が s_final であるフィルターに対して、関数は |s[k] – s_final| < 0.02×max(|s[k] – s_final|) となる最小のインデックス k を過渡長として決定します。

  k が信号長よりも大きい場合、関数は信号長 length(x) を過渡長として使用します。

• length(x) 以下の非負の整数 — filtfilt は、指定された過渡長を使用します。

PaddingPatern を "none" として指定した場合、filtfilt では、TransientLength に指定した値が無視され、フィルター処理後の信号にモデル化または抑制を行うべき過渡は存在しないと仮定されます。

PaddingPattern — フィルター処理前に入力信号をパディングするパターン

"linear" (既定値) | "reflect" | "none"

R2026a 以降

フィルター処理前に入力信号をパディングするパターン。"linear"、"reflect"、または "none" として指定します。

関数は、TransientLength に指定した値を使用して、操作次元に沿って入力信号 x をパディングするサンプル数を決定します。関数は、パディングしたサンプルを使用して過渡応答をモデル化し、フィルター処理後の信号の両端におけるモデル化された過渡応答を抑制してから、y を返します。

以下の表に、パターンの名称と説明、入力データ x = [a b c d ... j k l m] を各パターンがどのようにパディングするかの例を示します。

パディング パターン	説明	パディングされたデータ
"linear"	関数は、境界における傾きの連続性を維持するように入力信号の両端をパディングします。 このパディング パターンは、信号が遷移領域に滑らかに続くようにし、境界での過渡を抑制するのに役立ちます。
"reflect"	関数は、端点を重複させずに、各端でサンプルをミラーリングすることで入力信号をパディングします。
"none"	関数は入力信号をパディングしません。 この場合、filtfilt は TransientLength に指定した値も無視します。

出力引数

y — フィルター処理された信号

ベクトル | 行列 | N 次元配列

フィルター処理された信号。ベクトル、行列、または N 次元配列として返されます。

• filtfilt 関数は、x と同じサイズの y を返します。

• 入力引数を single 型として指定した場合、filtfilt は単精度演算を使用してフィルター処理演算を実行し、single 型として y を返します。それ以外の場合、filtfilt は double 型として y を返します。

詳細

カスケード伝達関数

IIR デジタル フィルターをカスケード セクションに分割すると、数値安定性が向上し、係数の量子化誤差の影響を受けにくくなります。L 個の伝達関数 H₁(z)、H₂(z)、…、H_L(z) による伝達関数 H(z) のカスケード形式は次のとおりです。

$$H(z)={\displaystyle \prod _{l=1}^L{H}_l(z)}=H_1(z)\times H_2(z)\times \cdots \times H_L(z).$$

CTF 形式によるデジタル フィルターの指定

解析、可視化、信号フィルター処理で使用するため、CTF 形式でデジタル フィルターを指定できます。係数 B と A をリストしてフィルターを指定します。スカラーまたはベクトルの g を指定することで、セクション全体にフィルターのスケーリング ゲインを含めることもできます。

フィルター係数

係数を L 行の行列

$$B=\left[\begin{array}{cccc}{b}_{11}& b_{12}& \cdots & b_{1,m+1}\\ b_{21}& b_{22}& \cdots & b_{2,m+1}\\ ⋮& ⋮& \ddots & ⋮\\ b_{L1}& b_{L2}& \cdots & b_{L,m+1}\end{array}\right],A=\left[\begin{array}{cccc}{a}_{11}& a_{12}& \cdots & a_{1,n+1}\\ a_{21}& a_{22}& \cdots & a_{2,n+1}\\ ⋮& ⋮& \ddots & ⋮\\ a_{L1}& a_{L2}& \cdots & a_{L,n+1}\end{array}\right],$$

として指定した場合、L 個のカスケード伝達関数のシーケンスとしてフィルターを指定したものとみなされます。このとき、フィルターの完全な伝達関数は次のようになります。

$$H\left(z\right)=\frac{b_{11}+b_{12}{z}^{-1}+\cdots +b_{1,m+1}{z}^{-m}}{a_{11}+a_{12}{z}^{-1}+\cdots +a_{1,n+1}{z}^{-n}}\times \frac{b_{21}+b_{22}{z}^{-1}+\cdots +b_{2,m+1}{z}^{-m}}{a_{21}+a_{22}{z}^{-1}+\cdots +a_{2,n+1}{z}^{-n}}\times \cdots \times \frac{b_{L1}+b_{L2}{z}^{-1}+\cdots +b_{L,m+1}{z}^{-m}}{a_{L1}+a_{L2}{z}^{-1}+\cdots +a_{L,n+1}{z}^{-n}},$$

ここで、m ≥ 0 はフィルターの "分子次数"、n ≥ 0 は "分母次数" です。

• B と A の両方をベクトルとして指定した場合、基になるシステムは 1 つのセクションから成る IIR フィルター (L = 1) であるとみなされます。このとき、B は伝達関数の分子を表し、A はその分母を表します。

• B がスカラーの場合、フィルターは全極 IIR フィルターのカスケードであるとみなされます。このとき、各セクションの全体的なシステム ゲインは B に等しくなります。

• A がスカラーの場合、フィルターは FIR フィルターのカスケードであるとみなされます。このとき、各セクションの全体的なシステム ゲインは 1/A に等しくなります。

メモ

• 2 次セクションの行列をカスケード伝達関数に変換するには、sos2ctf 関数を使用します。

• 零点-極-ゲイン フィルター表現をカスケード伝達関数に変換するには、zp2ctf 関数を使用します。

係数とゲイン

全体的なスケーリング ゲイン、または係数値から共通因数としてくくりだされた複数のスケーリング ゲインがある場合、係数とゲインを {B,A,g} 形式の cell 配列として指定できます。フィルター セクションのスケーリングは、固定小数点演算を使用して各フィルター セクションの出力の振幅レベルを確実に同等にする場合に特に重要であり、数値精度の制限によるフィルター応答の不正確さを回避するのに役立ちます。

ゲインは、全体的なゲインを表すスカラー、またはセクションのゲインを表すベクトルとして指定できます。

• ゲインがスカラーの場合、その値はすべてのカスケード フィルター セクションに均一に適用されます。

• ゲインがベクトルの場合、カスケード内のフィルター セクションの数 L よりも 1 つ多い要素をもたなければなりません。最初の L 個のスケール値はそれぞれ対応するフィルター セクションに適用され、最後の値はすべてのカスケード フィルター セクションに均一に適用されます。

係数行列とゲイン ベクトルを

$$B=\left[\begin{array}{cccc}{b}_{11}& b_{12}& \cdots & b_{1,m+1}\\ b_{21}& b_{22}& \cdots & b_{2,m+1}\\ ⋮& ⋮& \ddots & ⋮\\ b_{L1}& b_{L2}& \cdots & b_{L,m+1}\end{array}\right],A=\left[\begin{array}{cccc}{a}_{11}& a_{12}& \cdots & a_{1,n+1}\\ a_{21}& a_{22}& \cdots & a_{2,n+1}\\ ⋮& ⋮& \ddots & ⋮\\ a_{L1}& a_{L2}& \cdots & a_{L,n+1}\end{array}\right],g=\left[\begin{array}{ccccc}{g}_1& g_2& \cdots & g_L& g_{\text{S}}\end{array}\right],$$

として指定した場合、フィルター システムの伝達関数は次であるとみなされます。

$$H\left(z\right)=g_{\text{S}}\left(g_1\frac{b_{11}+b_{12}{z}^{-1}+\cdots +b_{1,m+1}{z}^{-m}}{a_{11}+a_{12}{z}^{-1}+\cdots +a_{1,n+1}{z}^{-n}}\times g_2\frac{b_{21}+b_{22}{z}^{-1}+\cdots +b_{2,m+1}{z}^{-m}}{a_{21}+a_{22}{z}^{-1}+\cdots +a_{2,n+1}{z}^{-n}}\times \cdots \times g_L\frac{b_{L1}+b_{L2}{z}^{-1}+\cdots +b_{L,m+1}{z}^{-m}}{a_{L1}+a_{L2}{z}^{-1}+\cdots +a_{L,n+1}{z}^{-n}}\right).$$

ヒント

スケーリング ゲインを含むフィルターを CTF 形式で取得できます。butter、cheby1、cheby2、ellip などのデジタル IIR フィルター設計関数の出力を使用します。これらの関数で、"ctf" フィルタータイプ引数を指定し、スケール値を取得するために B、A、g を返すように指定します。 (R2024b 以降)

参照

[1] Gustafsson, F. “Determining the initial states in forward-backward filtering.” IEEE^® Transactions on Signal Processing. Vol. 44, April 1996, pp. 988–992. https://doi.org/10.1109/78.492552.

[2] Lyons, Richard G. Understanding Digital Signal Processing. Upper Saddle River, NJ: Prentice Hall, 2004.

[3] Mitra, Sanjit K. Digital Signal Processing. 2nd Ed. New York: McGraw-Hill, 2001.

[4] Oppenheim, Alan V., and Ronald W. Schafer, with John R. Buck. Discrete-Time Signal Processing. 2nd Ed. Upper Saddle River, NJ: Prentice Hall, 1999.

拡張機能

C/C++ コード生成

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

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

digitalFilter オブジェクトは、コード生成でサポートされません。

GPU コード生成

GPU Coder™ を使用して NVIDIA® GPU のための CUDA® コードを生成します。

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

digitalFilter オブジェクトは、コード生成でサポートされません。

GPU 配列

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

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

バージョン履歴

R2006a より前に導入

R2026a: ゼロ位相デジタル フィルター処理を実行する新しいオプション

filtfilt 関数は、追加のオプションを指定して信号にゼロ位相フィルター処理を適用できる名前と値の引数をサポートしています。

• InitialStatesMethod — フィルター状態の初期条件を推定するアルゴリズムを選択します。

• TransientLength — フィルター処理後の信号の過渡応答を抑制するためにサンプル数を指定するか、推定方法を選択します。

• PaddingPattern — フィルター処理前に入力信号の両端をパディングするパターンを選択します。

R2024b: カスケード伝達関数を使用した信号のフィルター処理

filtfilt 関数は、カスケード伝達関数 (CTF) 形式の入力をサポートします。

参考

ctffilt | designfilt | digitalFilter | fftfilt | filter | filter2 | paddata

トピック

• 信号の 60 Hz ハムの削除
• 実践に即したデジタル フィルター処理の紹介
• 非因果的ゼロ位相フィルターの実装