フィルター処理によって生じる遅延の補正

デジタルフィルターを使用すると信号に遅延が生じます。フィルター特性によって，遅延が全周波数にわたって一定の场合と周波数に応じて変化する场合があります。遅延を补正するために行う处理は，遅延のタイプによって异なります。关数grpdelayを使用すると，フィルター遅延を周波数关数としてとらえることができます。この关数出力を见ることによって，フィルター遅延が一定であるか，周波数に応じて変化する（つまり，遅延が周波数に依存する）かを特定することができます。

フィルター遅延が全周波数にわたって一定の場合,信号を時間シフトすることによって容易に補正できます。冷杉フィルターでは通常,遅延が一定です。一方,遅延が周波数に応じて変化する場合,位相歪みが生じて信号波形が著しく変わることがあります。周波数に依存する遅延の補正は,遅延が一定の場合ほど簡単ではありません。IIRフィルターでは,周波数に依存する遅延が生じます。

定遅延フィルターの补正

前述のように，フィルターの群遅延を测定して周波数关数が一定であることを确认できます。关数grpdelayを使用してフィルター遅延Dを測定し,入力信号にD個の0を追加して出力信号をDサンプルだけ時間シフトすると,この遅延を補正できます。

ノイズを含んだ心电図信号について，75赫兹を上回る高周波ノイズを除去するためにフィルター处理する例を考えてみます.FIRローパスフィルターをかけてから，フィルター遅延を补正します。これにより，ノイズを含んだ信号とフィルター处理后の信号が正しく整列され，比较できるようにプロットが互いに重ね合わせられます。

FS = 500;%采样频率，单位为HzN = 500;%信号采样个数RNG.默认；x =心电图(N) + 0.25 * randn (N, 1);％噪声波形T =（0：N-1）/ FS;%的时间向量

75赫兹のカットオフ周波数をもつ70次ローパスFIRフィルターを设计します。

Fnorm = 75 / (Fs / 2);%归一化频率DF = designfilt（'lowpassfir'，'FilterOrder', 70,“CutoffFrequency”, Fnorm);

フィルターの群遅延をプロットして,全周波数で群遅延が一定であり,このフィルターが線形位相であることを確認します。群遅延を使用してフィルター遅延を測定します。

grpdelay（DF，2048，FS）图群延迟

D =意味着(grpdelay (df))采样滤波器延迟%

d = 35

フィルター处理の前に，入力データベクトルXの最后にd个の0を追加します。これにより，有用なサンプルはすべてフィルターから除外され，また，入力信号と遅延补正された出力信号が同じ长さになることが确実になります。データのフィルター处理を行い，出力信号をdサンプルだけ时间シフトして遅延を补正します。この最后の手顺で，フィルターの过渡特性が效果的に除去されます。

y =过滤器(df, [x;0 (D, 1)]);%在输入数据中添加D个零Y = Y（d + 1：结束）;％移位数据以补偿延迟情节(t t, x,, y,“r”，“线宽”, 1.5)标题(过滤后的波形的）Xlabel（“时间(s)”） 传奇（“原始噪声信号”，“滤波信号”） 网格在轴紧

周波数依存の遅延の补正

周波数に依存する遅延は信号の位相歪みを発生させます。このタイプの遅延を補正することは,遅延が一定である場合ほど簡単ではありません。オフライン処理が可能なアプリケーションの場合,関数filtfiltを使用したゼロ位相フィルター处理によって周波数依存の遅延を除去できます。filtfiltは入力データの处理を顺方向と逆方向の両方向に行うことにより，ゼロ位相フィルター处理を行います。主な效果はゼロ位相の歪みを得ること，つまり，データは0サンプルの定遅延をもつ同等のフィルターによってフィルター处理されます。その他の效果は，元のフィルターの伝达关数の振幅の2乘に等しいフィルター伝达关数と，元のフィルターの次数の2倍のフィルター次数が得られることです。

前節で定義した心電図信号を考えてみましょう。この信号を,遅延の補正がある場合とない場合に分けてフィルター処理します。75赫兹のカットオフ周波数をもつ7次ローパスIIR楕円フィルターを設計します。

Fnorm = 75 / (Fs / 2);%归一化频率DF = designfilt（“lowpassiir”，......'PassbandFrequency'Fnorm,......'FilterOrder'7......“PassbandRipple”，1，......'StopbandAttenuation'、60);

フィルターの群遅延をプロットします。群遅延が周波数に応じて変化し，フィルター遅延が周波数に依存しています。

grpdelay (df, 2048,'一半'Fs)

データをフィルター处理して，时间信号に対するそれぞれのフィルターの実装の效果を见てみます。ゼロ位相フィルター处理がフィルター遅延を效果的に除去しています。

Y1 =滤波器（DF，X）;％非线性相位滤波器 - 没有延迟补偿Y2 = filtfilt（DF，X）;％零阶段实施 - 延迟补偿图（T，X）保持在图（T，Y，“r”，“线宽”1.5)情节(t, y2,“线宽”, 1.5)标题(过滤后的波形的）Xlabel（“时间(s)”） 传奇（原始信号的，“非线性相位IIR输出”，......“零相IIR输出”xlim([0.25 0.55])网格在

因果性のない顺方向/逆方向のフィルター处理が可能で，フィルター応答を元の応答の2乘に変更できるアプリケーションの场合，ゼロ位相フィルター处理は优れたツールです。

定遅延が発生するフィルターは线形位相フィルターです。周波数依存の遅延が発生するフィルターは非线形位相フィルターです。

信号の望ましくないスペクトル成分の除去

フィルターは,信号からの望ましくないスペクトル成分の除去によく使用されます。これを行うため,さまざまな種類のフィルターからの選択が可能です。高周波数成分を除去する場合はローパスフィルターを選び,低周波数成分を除去する場合はハイパスフィルターを選びます。また,中間周波数帯域をそのまま残して低周波数成分も高周波数成分も除去する場合は,バンドパスフィルターを選びます。与えられた帯域の周波数成分を除去する場合は,バンドストップフィルターを選びます。

電源供給ラインからのハムおよびホワイトノイズを含むオーディオ信号について考えてみます。電源供給ラインからのハムは60 Hzトーンが原因となって引き起こされます。ホワイトノイズはオーディオ帯域全体にわたって存在する信号です。

オーディオ信号を読み込みます。サンプルレートを44.1千赫に指定します。

Fs = 44100;y = audioread (“noisymusic.wav”）;

信号のパワースペクトルをプロットします。赤い三角形のマーカーは,オーディオ信号に干渉する強い60 Hzトーンを示します。

[P F] = pwelch (y)的(8192 1),8192/2,8192 Fs,'力量'）;helperFilterIntroductionPlot1 (F P [60 60], [-9.365 - -9.365],......｛“原始信号功率谱”，“60 Hz的语气”}）

最初に,できるだけ多くのホワイトノイズスペクトル成分をローパスフィルターで除去します。フィルターの通過帯域の値は,ノイズ除去と高周波数成分の損失によるオーディオ品質の劣化との適切なトレードオフが得られる設定にすべきです。60 Hzハムを除去する前にローパスフィルターをかけると,帯域制限された信号のダウンサンプリングが可能になり,非常に便利です。信号のレートを低くすると,よりシャープで狭い60 Hzバンドストップフィルターを低いフィルター次数で設計できます。

通过帯域周波数1千赫，阻止帯域周波数1.4千赫のローパスフィルターを设计します。最小次数设计を选択します。

FP = 1E3;％通带频率以HzFST = 1.4e3;停止频带频率，单位为Hz美联社= 1;％通带波纹，以dBAST = 95;阻带衰减百分比，单位为dBDF = designfilt（'lowpassfir'，'PassbandFrequency'《外交政策》,......“StopbandFrequency”置,“PassbandRipple”，鸭，......'StopbandAttenuation'，AST，'采样率'Fs);

フィルター応答を解析します。

fvtool（DF，'fs'Fs,'FrequencyScale'，“日志”，......“FrequencyRange”，“指定freq.向量”，“FrequencyVector”，F）

データをフィルター处理し，遅延を补正します。

d =平均值（grpdelay（DF））;%过滤器延迟ylp =过滤器(df, [y;0 (D, 1)]);ylp = ylp (D + 1:结束);

ローパスフィルター处理を行った信号のスペクトルを见てみます0.1400赫兹を超える周波数成分が除去されています。

Flp (Plp) = pwelch (ylp (8192 1), 8192/2, 8192 Fs,'力量'）;Flp helperFilterIntroductionPlot1 (F P、Plp......｛“原始信号”，“低通滤波的信号”}）

上のパワースペクトルのプロットより,ローパスフィルター処理を行った信号の無視できない周波数成分の最大値は1400 Hzであることがわかります。サンプリング理論により,サンプルレート $2\times 1400＝2800$は信号を正確に表すのに十分です。しかし,必要以上のサンプルを処理する必要があって無駄な44100 Hzのサンプルレートを使用しています。信号をダウンサンプリングしてサンプルレートを下げ,処理サンプル数を減らすことで計算負荷を軽減できます。サンプルレートを低くすると,60赫兹ノイズ除去に必要なよりシャープで狭いバンドストップフィルターを低いフィルター次数で設計することができます。

サンプルレートFs / 10 = 4.41 kHzを得るために,ローパスフィルター処理を行った信号を係数10でダウンサンプリングします。ダウンサンプリング前後の信号のスペクトルをプロットします。

Fs = f / 10;码= downsample (ylp 10);(Pds, Fds) = pwelch(码,(8192 1),8192/2,8192 Fs,'力量'）;helperFilterIntroductionPlot1（F，P，FDS，PDS，......｛“在44100赫兹的信号采样”，'下采样信号，Fs = 4410hz '}）

次に,IIRバンドストップフィルターで60 Hzトーンを除去します。阻止帯域が中心60 Hz,幅4赫兹になるようにします。シャープな周波数ノッチおよび小さい通過帯域リップルを比較的低い次数で実現できるIIRフィルターを選択します。

DF = designfilt（“bandstopiir”，'PassbandFrequency1'55，......“StopbandFrequency1”58，“StopbandFrequency2”, 62,......'PassbandFrequency2'65，'PassbandRipple1'，1，......'StopbandAttenuation'现年60岁的'PassbandRipple2'，1，......'采样率'Fs,'DesignMethod'，“ellip”）;

振幅応答を解析します。

fvtool（DF，'fs'Fs,'FrequencyScale'，“日志”，......“FrequencyRange”，“指定freq.向量”，“FrequencyVector”，FDS（FDS> F（2）））

位相歪みを防ぐために,ゼロ位相フィルター処理を実行します。関数filtfiltを使用してデータを処理します。

YBS = filtfilt（DF，码）;

最後に,信号をアップサンプリングして,オーディオサウンドカードに準拠した元のオーディオサンプルレート44.1 kHzに戻します。

yf =插值函数(yb, 10);Fs = f * 10;

最後に,元の信号のスペクトルとフィルター処理後のスペクトルを見てみます。フィルターによって,高周波ノイズフロアと60 Hzトーンが減衰しています。

[Pfinal, Ffinal] = pwelch (yf, (8192 1), 8192/2, 8192 Fs,'力量'）;helperFilterIntroductionPlot1（F，P，Ffinal，Pfinal，......｛“原始信号”，“最终过滤信号”}）

コンピューターにサウンドカードが搭載されているなら,処理前後の信号を再生できます。上述したように,最終的にオーディオファイルの60 Hzハムと高周波ノイズを効果的に減衰させることができました。

％为了发挥原有的信号，取消接下来的两行％hplayer = audioplayer（Y，FS）;% (hplayer)玩要播放降噪信号，取消注释下面两行% hplayer = audioplayer(yf, Fs);％播放（hplayer）;

信号の微分

MATLAB関数diffによる信号の微分には,出力時ノイズレベル増大の可能性という欠点があります。このため,より良い選択として微分器フィルターを使用します。微分器フィルターは,目的の帯域で微分器として,他の全周波数では減衰器として動作し,高周波ノイズを効果的に除去します。

例として，地震発生时の建物の床の変位速度を解析します。変位またはドリフトの测定値は，地震の条件下で3阶建ての试験构造物の1阶で记录されてquakedrift.matファイルに保存されました。データベクトルの长さは10E3，サンプルレートは1千赫，测定単位は厘米です。

変位データを微分して,地震発生時の建物の床の速度と加速度の推定値を求めます。diffおよび冷杉微分器フィルターを使用した結果を比較します。

负载quakedrift.matFS = 1000;％ 采样率dt = 1 / f;％时间差t =(0:长度(漂移)1)* dt;%的时间向量

ほとんどの信号エネルギーが存在する帯域である100 Hzの通過帯域周波数の50次微分器フィルターを設計します。フィルターの阻止帯域周波数を120 Hzに設定します。

DF = designfilt（“differentiatorfir”，'FilterOrder', 50岁,......'PassbandFrequency', 100,“StopbandFrequency”, 120,......'采样率'Fs);

関数diffは,応答 $$H（Z）＝1-Z^{-1}$$をもつ1次FIRフィルターと见なすことができます0.50次微分器FIRフィルターの振幅応答と关数diffの応答を比較するために,FVToolを使用します。明らかに,0 ~ 100 Hzの通過帯域領域で両方の応答が等価です。しかし,阻止帯域領域では,50次微分器フィルターが成分を減衰させているのに対し,diffの応答は成分を増幅させています。これが実質的に,高周波ノイズのレベルを上げています。

hfvt = fvtool（DF，[1 -1]，1，'MagnitudeDisplay'，“零相”，'fs'Fs);传奇(hfvt“50阶FIR微分”，“DIFF函数的响应”）;

関数diffを使用して微分します为.diffの演算によって欠损しているサンプルを补正するために，0を追加します。

v1 = diff / dt(漂移);a1 = diff / dt (v1);v1 = [0;v1);a1 = [0;0;a1];

50次微分器冷杉フィルターを使用して微分し,遅延を補正します。

d =平均值（grpdelay（DF））;%过滤器延迟v2 =过滤器(df,漂移;0 (D, 1)]);v2 = v2 (D + 1:结束);a2 =过滤器(df, v2;0 (D, 1)]);a2 = a2 (D + 1:结束);v2 = v2 / dt;a2 = a2 / dt ^ 2;

床の変位データ点をいくつかプロットします。diffおよび50 次微分器 FIR フィルターによって算出した速度と加速度のデータ点もいくつかプロットします。diffを使用した場合,速度の推定値ではノイズがわずかに増幅し,加速度の推定値ではノイズが大幅に増幅していることがわかります。

helperFilterIntroductionPlot2 (t,漂移,v1、v2, a1, a2)

信号の积分

漏洩積分器フィルターは,伝達関数 $$H（Z）＝1/［1-c{Z}^{-1}］$$を使用する全極フィルターです。 $$c$$は定数で,フィルターの安定性を確保するためには1未満でなければなりません。当然ですが, $$c$$が1に近づくにつれて，漏泄积分器がdiffの伝达关数の逆关数に近づきます。前节で求めた加速度と速度の推定値を漏泄积分器で积分して，速度とドリフトをそれぞれ取り戻します。关数diffを使用して求めた推定値のほうがノイズが多いので,こちらを使用します。

漏洩積分器を $$\mathrm{一个}＝0．999$$で使用します。漏泄积分器フィルターの振幅応答をプロットします。フィルターは，效果的に高周波ノイズを除去するローパスフィルターとして动作します。

fvtool（1，[1 -.999]'fs'Fs)

漏泄积分器で，速度と加速度をフィルター处理します。时间差で乘算します。

v_original = v1;a_original = a1;D_leakyint = filter(1，[1 -0.999]，v_original);V_leakyint = filter(1，[1 -0.999]，a_original);D_leakyint = D_leakyint * dt;V_leakyint = V_leakyint * dt;

変位と速度の推定値をプロットし,元の信号と比較します。

helperFilterIntroductionPlot3 (t,漂移,d_leakyint、v_original v_leakyint)

関数浓汤およびcumtrapzを使用しても信号を積分できます。結果は,漏洩積分器で求めた結果と同様になります。

まとめ

この例では,線形位相フィルターおよび非線形位相フィルターについて,また,それぞれのフィルタータイプによって生じる位相遅延の補正方法について説明しました。さらに,フィルターをかけて望ましくない周波数成分を信号から取り除く方法やローパスフィルターで帯域を制限した後で信号をダウンサンプリングする方法について説明しました。最後に,デジタルフィルター設計を使用した,信号の微分および積分方法を説明しました。この例の全体を通して,解析ツールを使用してフィルターの応答と群遅延を確認する方法についても説明しました。

フィルター适用の详细については，信号处理工具箱™ドキュメンテーションを参照してください。デジタルフィルターの设计方法の详细については，実践に即したデジタルフィルター設計の紹介の例を参照してください。

参考文献

Proakis, j.g.，和d.g. Manolakis。数字信号处理:原理，算法和应用.Englewood峭壁，NJ：普伦蒂斯霍尔，1996年。

Orfanidis, s . J。信号处理概论.Englewood峭壁，NJ：普伦蒂斯霍尔，1996年。

付录

この例では,次の補助関数が使用されています。