JP3435675B2 - Adaptive control method - Google Patents

Description

【発明の詳細な説明】Detailed Description of the Invention

【０００１】[0001]

【発明の属する技術分野】この発明は、音場制御や能動
騒音制御、振動制御などに適用され、制御点における誤
差信号のパワーを最小にするように、可変フィルタの特
性を適応的に制御して可変フィルタの出力信号を制御
し、特に可変フィルタの出力端と制御点までの伝達特性
を可変フィルタの入力信号に付与した信号に基づいて制
御を行う適応的制御方法に関する。BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention is applied to sound field control, active noise control, vibration control, etc., and adaptively controls the characteristics of a variable filter so as to minimize the power of an error signal at a control point. The present invention relates to an adaptive control method for controlling an output signal of a variable filter, and in particular, performing control based on a signal in which a transfer characteristic from an output end of the variable filter to a control point is added to an input signal of the variable filter.

【０００２】[0002]

【従来の技術】[Prior art]

適応フィルタと filtered-x 構成 図５Ａは適応フィルタを用いた適応制御系のブロック図
を示す。また、適応フィルタ１１はＦＩＲ形フィルタで
あって、時変係数｛ｈ₁ （ｋ），ｈ₂ （ｋ），... ，ｈ
_L （ｋ）｝を持つものとする。ただし、Ｌはフィルタ係
数の数である。また、フィルタ係数はベクトルとして、 （１）ｈ(k) ＝[h₁(k)，ｈ₂(k)，... ，ｈ_L (k)]^T と表す。ただし、[ ] ^T はベクトルまたは行列の転置を
表す。またこの明細書では、時間の表現は離散時間とし
て整数ｋで時間を表すものとする。Adaptive Filter and Filtered-x Configuration FIG. 5A shows a block diagram of an adaptive control system using an adaptive filter. Further, the adaptive filter 11 is an FIR type filter and has time-varying coefficients {h ₁ (k), h ₂ (k), ..., H.
_L (k)}. However, L is the number of filter coefficients. The filter coefficient is expressed as a vector (1) h (k) = [h ₁ (k), h ₂ (k), ..., H _L (k)] ^T. However, [] ^T represents the transpose of a vector or matrix. Further, in this specification, time is represented by an integer k as discrete time.

【０００３】適応フィルタ１１は、誤差信号ｅ（ｋ）と
入力信号ｘ（ｋ）に基づいて、フィルタ係数ｈ（ｋ）
を調整して、誤差信号ｅ（ｋ）のパワーを最小にする。
この時の係数修正の手順は適応アルゴリズムとして知ら
れている。図５Ｂに図５Ａの適応制御系を騒音制御に直
接適用した図を示す。適応フィルタ１１の出力は位相反
転器１２を通じてスピーカ１３へ供給され、スピーカ１
３よりの音響信号ｙ’（ｋ）は制御点１４に到達させ、
制御点１４における騒音ｎ（ｋ）を抑圧するようにす
る。制御点１４にはマイクロホン１５が設置され、誤差
信号ｅ（ｋ）をモニタする。この適応制御系の目的は、
スピーカ１３から出た音ｙ’（ｋ）で騒音ｎ（ｋ）を消
去することである。The adaptive filter 11 has a filter coefficient h (k) based on the error signal e (k) and the input signal x (k).
Is adjusted to minimize the power of the error signal e (k).
The procedure of coefficient correction at this time is known as an adaptive algorithm. FIG. 5B shows a diagram in which the adaptive control system of FIG. 5A is directly applied to noise control. The output of the adaptive filter 11 is supplied to the speaker 13 through the phase inverter 12, and the speaker 1
The acoustic signal y ′ (k) from 3 reaches the control point 14,
The noise n (k) at the control point 14 is suppressed. A microphone 15 is installed at the control point 14 and monitors the error signal e (k). The purpose of this adaptive control system is
The noise n (k) is eliminated by the sound y '(k) emitted from the speaker 13.

【０００４】図５Ｃに図５Ｂの等価ブロック図を示し、
図５Ｂと共通要素は同一の記号を付けてある。スピーカ
１３の入力端から制御点１４までの伝達関数Ｃ（スピー
カ１３の電気音響変換特性とスピーカ１３から制御点１
４までの音響特性）を持った伝達系１６が適応フィルタ
１１と制御点１４との間に介在されていると表せる。ま
た、図５Ｂ中の騒音ｎ（ｋ）は所望信号ｄ（ｋ）と表
し、位相反転器１２による演算は制御点１４における差
演算として表した。図５Ｃを図５Ａと比較した場合、図
５Ａにおいてはフィルタ出力ｙ（ｋ）が直接、所望信号
ｄ（ｋ）と差演算されているのに対して、図５Ｃではフ
ィルタ出力ｙ（ｋ）は伝達系１６を経た後、ｄ（ｋ）と
差演算されている点が相違点である。通常の適応アルゴ
リズムでは、フィルタ出力ｙ（ｋ）を直接演算して誤差
信号ｅ（ｋ）を計算することが仮定されているので、図
５Ｃの系では適応アルゴリズムは良好に動作しない。FIG. 5C shows an equivalent block diagram of FIG. 5B.
Elements common to FIG. 5B are labeled with the same symbols. Transfer function C from the input end of the speaker 13 to the control point 14 (electroacoustic conversion characteristics of the speaker 13 and the control point 1 from the speaker 13)
It can be said that a transmission system 16 having acoustic characteristics up to 4) is interposed between the adaptive filter 11 and the control point 14. Further, the noise n (k) in FIG. 5B is expressed as the desired signal d (k), and the calculation by the phase inverter 12 is expressed as the difference calculation at the control point 14. When FIG. 5C is compared with FIG. 5A, the difference between the filter output y (k) and the desired signal d (k) is directly calculated in FIG. 5A, whereas the filter output y (k) is calculated in FIG. 5C. The difference is that after passing through the transmission system 16, a difference calculation with d (k) is performed. Since it is assumed that the normal adaptive algorithm directly calculates the filter output y (k) to calculate the error signal e (k), the adaptive algorithm does not work well in the system of FIG. 5C.

【０００５】そこで、図６Ａに示す系を考える。この系
においては図５Ｃ中の適応フィルタ１１と伝達系１６と
の順序が入れ替わっている。適応フィルタ１１の特性の
変化がゆるやかであると仮定すれば、このような順序の
変更を行っても、図５Ｃおよび図６Ａ中の信号ｘ（ｋ）
とｙ’（ｋ）の関係は維持される。従って、図５Ｃ中の
ｙ’（ｋ）と図６Ａ中のｙ’（ｋ）は同じ信号となる。
図５Ｃと図６Ａとで異なっている点は、適応フィルタ１
１に対する入力が図５Ｃではｘ（ｋ）であるのに対し、
図６Ａでは、ｘ（ｋ）をＣの特性を持った伝達系１６に
通した信号ｘ’（ｋ）となっている点である。また、図
６Ａでは、誤差信号ｅ（ｋ）を得るのに適応フィルタ１
１の出力ｙ’（ｋ）を直接利用しており、このようにす
ることで従来の適応アルゴリズムを利用することが可能
となる。Therefore, consider the system shown in FIG. 6A. In this system, the order of the adaptive filter 11 and the transfer system 16 in FIG. 5C is exchanged. Assuming that the characteristic change of the adaptive filter 11 is gentle, even if the order is changed, the signal x (k) in FIGS. 5C and 6A is changed.
And y '(k) is maintained. Therefore, y ′ (k) in FIG. 5C and y ′ (k) in FIG. 6A are the same signal.
The difference between FIG. 5C and FIG. 6A is that the adaptive filter 1
While the input for 1 is x (k) in FIG. 5C,
In FIG. 6A, the point is that the signal x ′ (k) is obtained by passing x (k) through the transmission system 16 having the characteristic of C. Also, in FIG. 6A, the adaptive filter 1 is used to obtain the error signal e (k).
The output y ′ (k) of 1 is directly used, and by doing so, the conventional adaptive algorithm can be used.

【０００６】しかし、図５Ｂに示した騒音制御系などに
おいて、スピーカ１３と適応フィルタ１１との順序を入
れ替えることは物理的に不可能である。そこで、図６Ｂ
の系を考える。騒音制御系におけるスピーカ１３の入力
端から制御点１４までの伝達系１６の特性Ｃはあらかじ
め測定または推定しておき、これをフィルタ係数計算部
１７に含まれるフィルタ１８の特性として与える。この
系ではまず、フィルタ実行部１９において、入力信号ｘ
（ｋ）を可変フィルタ１１に通すことで出力ｙ（ｋ）を
合成し、その合成出力を伝達系１６を経て制御点１４に
発生させた音と雑音ｄ（ｋ）とから誤差信号ｅ（ｋ）を
得る。以上の動作は図５Ｃの系と同様である。However, in the noise control system and the like shown in FIG. 5B, it is physically impossible to switch the order of the speaker 13 and the adaptive filter 11. Therefore, FIG. 6B
Consider the system. The characteristic C of the transmission system 16 from the input end of the speaker 13 to the control point 14 in the noise control system is measured or estimated in advance, and this is given as the characteristic of the filter 18 included in the filter coefficient calculation unit 17. In this system, first, in the filter execution unit 19, the input signal x
The output y (k) is combined by passing (k) through the variable filter 11, and the combined output is generated from the sound and noise d (k) generated at the control point 14 via the transmission system 16 and the error signal e (k). ) Get. The above operation is similar to that of the system of FIG. 5C.

【０００７】一方、フィルタ係数計算部１７では、入力
信号ｘ（ｋ）を特性Ｃを持ったフィルタ１８を通した
ｘ’（ｋ）と誤差信号ｅ（ｋ）を用いて適応部２１にお
いてフィルタ係数を計算する。この動作は図６Ａに示し
たものと同様である。図５Ｃと図６Ａとの比較として先
に説明したのと同様に、図６Ｂにおける伝達系１６の出
力ｙ’（ｋ）と、適応部２１の出力ｙ’（ｋ）は同じと
なる。従って、誤差信号ｅ（ｋ）は適応部２１の出力を
使って計算したものとみなすことができ、図６Ａの系と
の等価性が説明される。従って、適応部２１では従来の
適応アルゴリズムが使用可能である。適応部２１で計算
されたフィルタ係数ｈ（ｋ）は、毎時刻、可変フィル
タ１１に転送される。On the other hand, in the filter coefficient calculation section 17, the input signal x (k) is passed through the filter 18 having the characteristic C and x '(k) and the error signal e (k) are used in the adaptation section 21 to obtain the filter coefficient. To calculate. This operation is similar to that shown in FIG. 6A. As described above as a comparison between FIG. 5C and FIG. 6A, the output y ′ (k) of the transmission system 16 and the output y ′ (k) of the adaptation unit 21 in FIG. 6B are the same. Therefore, the error signal e (k) can be regarded as calculated using the output of the adaptation unit 21, and the equivalence with the system of FIG. 6A will be explained. Therefore, the adaptation unit 21 can use the conventional adaptation algorithm. The filter coefficient h (k) calculated by the adaptation unit 21 is transferred to the variable filter 11 every time.

【０００８】この図６Ｂに示したように、入力信号をフ
ィルタリングした後、適応アルゴリズムを使用する方法
は、 filtered-x 構成と呼ばれる。 filtered-x 構成
は、騒音制御や振動制御などのように、適応フィルタ１
１の出力点と制御点１４との間に伝達系１６を含む制御
系においては代表的な系構成である。 適応アルゴリズム 次に、適応アルゴリズムについて説明する。The method of using the adaptive algorithm after filtering the input signal as shown in FIG. 6B is called a filtered-x configuration. The filtered-x configuration is an adaptive filter 1 such as noise control and vibration control.
This is a typical system configuration in the control system including the transmission system 16 between the output point 1 of 1 and the control point 14. Adaptive Algorithm Next, the adaptive algorithm will be described.

【０００９】適応アルゴリズムとは、例えば図５Ａのよ
うに結線された適応フィルタ１１の内部の演算手順であ
って、時刻ｋにおいて、フィルタ係数ｈ（ｋ）を修正
して、誤差信号ｅ（ｋ）がなるべく小さくなるような係
数ｈ（ｋ＋１）を得る方法である。その修正は、次式
で表される。 （2.1) ｈ(k+1) ＝ｈ(k) ＋μ・Δｈ(k) ここで、Δｈ(k) は修正ベクトル、μは修正の大きさ
を制御する量でステップサイズと呼ばれる。The adaptive algorithm is a calculation procedure inside the adaptive filter 11 connected, for example, as shown in FIG. 5A. At the time k, the filter coefficient h (k) is corrected and the error signal e (k) is corrected. This is a method of obtaining a coefficient h (k + 1) such that is as small as possible. The correction is expressed by the following equation. (2.1) h (k + 1) = h (k) + μ · Δh (k) where Δh (k) is the correction vector and μ is the amount that controls the size of the correction and is called the step size.

【００１０】現在、最も多く利用されている適応アルゴ
リズムとしては、ＬＭＳ法や学習同定法などが知られて
いる。ＬＭＳ法においては、修正ベルトルは誤差信号ｅ
（ｋ）と入力ベクトルｘ（ｋ）との積として、 （2.2) Δｈ(k) ＝ｅ(k) ｘ(k) と表され、修正式は、 （2.3) ｈ(k+1) ＝ｈ(k) ＋μ・ｅ(k) ｘ(k) となる。ただし、 （3.1) ｘ(k) ＝[x(k),ｘ(k-1), ... ,ｘ(k-L+1)]^T で、これを入力ベクトルまたは単に入力と呼ぶ。また、
次の時刻ｋ＋１のフィルタ出力ｙ（ｋ＋１）は、次式に
よる内積演算により合成される。At present, the LMS method and the learning identification method are known as the most widely used adaptive algorithms. In the LMS method, the corrected belt is the error signal e
The product of (k) and the input vector x (k) is expressed as (2.2) Δh (k) = e (k) x (k), and the correction formula is (2.3) h (k + 1) = h (k) + μ · e (k) x (k). However, (3.1) x (k) = [x (k), x (k-1), ..., x (k-L + 1)] ^T , which is called an input vector or simply an input. Also,
The filter output y (k + 1) at the next time k + 1 is synthesized by the inner product calculation according to the following equation.

【００１１】 （3.2) ｙ（ｋ＋１）＝ｘ^T (k+1) ｈ(k+1) 式（２．３）に基づいてフィルタ係数を修正して、式
（３．２）に基づいて次の時刻の出力を合成するという
適応フィルタ１１の動作に必要な積和演算量は２Ｌ（Ｌ
はフィルタ係数の数）である。この２Ｌという演算量は
現在知られている適応アルゴリズムのなかでは最小の演
算量となっている。(3.2) y (k + 1) = x ^T (k + 1) h (k + 1) The filter coefficient is corrected based on the equation (2.3), and then based on the equation (3.2), The product-sum calculation amount required for the operation of the adaptive filter 11 to synthesize the outputs at the time
Is the number of filter coefficients). This calculation amount of 2L is the minimum calculation amount among currently known adaptive algorithms.

【００１２】このように、ＬＭＳ法や学習同定法は演算
量が少ないという長所を持っている。しかし、その反
面、有色信号や周期性雑音などに対しては、収束速度が
遅い（適切なフィルタ係数を得るまでの時間がかかる）
という欠点を持っている。ＬＭＳ法や学習同定法の欠点
を克服する方法として、１９８４年に射影法が提案され
たが、射影法では演算量が多くなるという欠点を持って
いた。しかし近年、この射影法の欠点を改善して、低演
算量で実行する高速射影法が提案された。その詳細は特
開平７−９２９８０号公報に記されているが、ここでは
その概要を以下に説明する。ただし、以下では、説明の
簡単のためにμ＝１を仮定する。As described above, the LMS method and the learning identification method have an advantage that the amount of calculation is small. However, on the other hand, the convergence speed is slow for colored signals and periodic noise (it takes time to obtain appropriate filter coefficients).
Has the drawback. The projection method was proposed in 1984 as a method for overcoming the drawbacks of the LMS method and the learning identification method, but the projection method had the drawback of increasing the amount of calculation. However, in recent years, a high-speed projection method has been proposed which improves the shortcomings of this projection method and executes with a small amount of calculation. Details thereof are described in Japanese Patent Application Laid-Open No. 7-92980, and the outline thereof will be described below. However, in the following, μ = 1 is assumed for simplicity of description.

【００１３】射影法において、誤差を小さくするために
は、フィルタ係数の修正は、次のｐ個の式 （４） d(k)＝ｘ^T (k) ｈ(k+1) d(k-1)＝ｘ^T (k-1) ｈ(k+1) ・・・ d(k-p+1)＝ｘ^T (k-p+1) ｈ(k+1) を満たすように行われる。ただし、ｐ（＜Ｌ）は射影の
次数と呼ばれる。この式（４）の意味は、修正を行った
後の係数ｈ（ｋ＋ｌ）を過去の入力ｘ^T （ｋ−ｉ＋
１）と畳み込んだ結果ｘ^T （ｋ−ｉ＋１）ｈ（ｋ＋
ｌ）が、過去の所望信号ｄ（ｋ−ｉ＋１）一致する、と
いう関係が過去ｐ個の入力に対して成立するように、
ｈ（ｋ＋ｌ）を決定するというものである。このよう
に、過去の入力に対して所望信号と同一の値を出力する
ようにフィルタ係数を定めれば、未来の時刻の入力ｘ
（ｋ＋ｌ）に対する出力ｙ（ｋ＋ｌ）も所望信号ｄ（ｋ
＋ｌ）に近い値が出力できるものと考えられる。その結
果、未来の誤差信号ｅ（ｋ＋ｌ）＝ｄ（ｋ＋ｌ）−ｙ
（ｋ＋ｌ）も小さな値となるものと考えられる。In the projection method, in order to reduce the error, the filter coefficient is corrected by the following p equations (4) d (k) = x ^T (k) h (k + 1) d (k- ^{1) = x T (k-} 1) h (k + 1) ··· d (k-p + 1) = x T (k-p + 1) is performed so as to satisfy the h (k + 1). However, p (<L) is called the projection order. The meaning of this equation (4) is that the coefficient h (k + 1) after the correction is applied to the past input x ^T (k−i +).
1) The result of convolution with x ^T (k-i + 1) h (k +
l) is the same as the desired signal d (k−i + 1) in the past, so that the relationship holds for the past p inputs.
It is to determine h (k + 1). In this way, if the filter coefficient is determined so as to output the same value as the desired signal with respect to the past input, the input x at the future time x
The output y (k + 1) for (k + 1) is also the desired signal d (k
It is considered that a value close to +1) can be output. As a result, future error signal e (k + 1) = d (k + 1) -y
It is considered that (k + 1) also has a small value.

【００１４】さて、式（２．１）を式（４）に代入して
解けば、 （５） Δｈ(k) ＝Ｘ(k)(Ｘ^T (k) Ｘ(k))^-1ｅ_p (k) と表される。ただし、 （６）Ｘ(k) ＝[ ｘ(k),ｘ(k-1), ... ,ｘ(k-p+1)] （7.1)ｅ_p (k) ＝[d(k), d(k-1), ... ,d(k-p+1)] ^T −Ｘ^T (k）ｈ(k) である。By substituting the equation (2.1) into the equation (4) and solving it, (5) Δh (k) = X (k) (X ^T (k) X (k)) ^-1 e _p It is expressed as (k). However, (6) X (k) = [x (k), x (k-1), ..., x (k-p + 1)] (7.1) e _p (k) = [d (k) , d (k-1), ..., a d (k-p + 1) ] T -X T (k) h (k).

【００１５】式（５）の右辺中の行列（Ｘ^T (k) Ｘ
(k) ）はＸ(k) を構成する列ベクトルの自己相関行列
と呼ばれる。また、このｅ_p (k) は、次式のように、
再帰的にも表すことができる。 ただし、 （7.3)ｅ_p-1(k-1)＝[e(k-1),(1−μ)e(k-2),...,(1−μ）^p-2 e(k-p+1)] ^T である。The matrix (X ^T (k) X in the right side of equation (5)
(k)) is called the autocorrelation matrix of the column vectors that make up X (k). Also, this e _p (k) is expressed by the following equation:
It can be expressed recursively. However, (7.3) e _p-1 (k-1) = [e (k-1), (1−μ) e (k-2), ..., (1−μ) ^p-2 e (k -p + 1)] ^T.

【００１６】式（５），（２．１）より、入力ｘ(k)
と誤差信号ｅ（ｋ）を用いてフィルタ係数を修正して、
係数ｈ（ｋ＋１）を求め、これを用いて時刻ｋ＋１の
出力ｙ（ｋ＋１）を合成する手順は、次式のように表す
ことができる。 （８）ｇ(k) ＝（Ｘ^T (k) Ｘ(k))^-1ｅ_p (k) （９）ｈ(k+1) ＝ｈ(k) ＋μ・Ｘ(k) ｇ(k) （10）ｙ(k+1) ＝ｘ^T (k+1) ｈ(k+1) ただし、ｇ（ｋ）は式（５）の後半を表し、プレフィ
ルタ係数と呼ばれる。この式（８）（９）（１０）が従
来の射影法の原理的手順を表す式であり、これを機能ブ
ロック図で表したものを図７Ａに示す。即ち、プレフィ
ルタ係数計算部２４で入力ベクトルＸ（ｋ）の自己相
関行列の逆行列と、誤差信号ベクトル（ｋ）との積、つ
まり式（８）を演算し、得られたプレフィルタ係数ｇ
（ｋ）と入力ベクトルＸ（ｋ）の積にステップサイズ
μを乗算し、これによりそれまでのフィルタ係数ｈ
（ｋ）を修正し、つまり式（９）を演算し、このフィル
タ係数を適応フィルタ１１の係数に設定して、その適応
フィルタ１１に入力ｘ（ｋ）を通し、つまり式（１０）
を演算して、出力ｙ（ｋ＋１）を得る。From equations (5) and (2.1), the input x (k)
And the error signal e (k) is used to modify the filter coefficient,
The procedure for obtaining the coefficient h (k + 1) and synthesizing the output y (k + 1) at time k + 1 using this can be expressed as in the following equation. (8) g (k) = (X T (k) X (k)) -1 e p (k) (9) h (k + 1) = h (k) + μ · X (k) g (k) (10) y (k + 1 ) = x T (k + 1) h (k + 1) However, g (k) represents the second half of the equation (5), called pre-filter coefficients. Equations (8), (9), and (10) are equations that represent the principle procedure of the conventional projection method, and a functional block diagram showing these is shown in FIG. 7A. That is, the prefilter coefficient calculator 24 calculates the product of the inverse matrix of the autocorrelation matrix of the input vector X (k) and the error signal vector (k), that is, the expression (8), and obtains the prefilter coefficient g
The product of (k) and the input vector X (k) is multiplied by the step size μ, which gives the filter coefficient h so far.
(K) is modified, that is, equation (9) is calculated, this filter coefficient is set to the coefficient of the adaptive filter 11, and the input x (k) is passed through the adaptive filter 11, that is, equation (10).
To obtain the output y (k + 1).

【００１７】ここで、この系の計算量を考えてみる。ま
ず、プレフィルタ係数計算部２４で行う式（８）の演算
は、ｐ行ｐ列の行列Ｘ^T （ｋ）Ｘ（ｋ）の逆行列の
計算を含んでいるので、ｐの３乗に比例する積和演算量
（以降、Ｏ（ｐ³)と表す）が必要となる。これは、ｐが
大きくなると、たいへん大きな演算量となる。次に、フ
ィルタ係数更新部２５で行う式（９）の演算はＸ
（ｋ）がＬ行ｐ列の行列、ｇ（ｋ）がｐ次のベクトル
であることより、ｐＬ回の積和演算が必要である。ま
た、フィルタ実行部１１で行う式（１０）の演算はＬ次
のベクトルの内積であるので、Ｌ回の積和演算が必要で
ある。従って、合計（ｐ＋１）Ｌ＋Ｏ（ｐ³)の積和演算
量が必要となる。例えば、ｐ＝２０，Ｌ＝１０００の場
合には、２1,０００回以上もの積和演算が必要となる。Now, let us consider the calculation amount of this system. First, since the calculation of the equation (8) performed by the pre-filter coefficient calculation unit 24 includes the calculation of the inverse matrix of the matrix ^XT (k) X (k) of p rows and p columns, it is proportional to the cube of p. The product-sum calculation amount (hereinafter referred to as O (p ³ )) is required. This is a very large calculation amount as p increases. Next, the calculation of the equation (9) performed by the filter coefficient updating unit 25 is X
Since (k) is a matrix of L rows and p columns and g (k) is a p-th order vector, pL times of product-sum operations are required. Further, since the calculation of the expression (10) performed by the filter execution unit 11 is an inner product of Lth-order vectors, L-time product-sum calculation is required. Therefore, a total of (p + 1) L + O (p ³ ) product-sum operation amount is required. For example, when p = 20 and L = 1000, 21,000 or more product-sum operations are required.

【００１８】これに対して近年提案された高速射影法で
は、積和演算量を２Ｌ＋２０ｐにまで低減することが可
能である。以下、特開平７−９２９８０号公報に従っ
て、高速射影法を簡単に説明する。まず、高速射影法で
は、式（８）の逆行列演算を直接行わずに線形予測法を
利用してプレフィルタ係数ベクトルｇ（ｋ）を計算す
ることで、Ｏ（ｐ³)の積和演算量を１６ｐの積和演算量
に低減している。この演算量の低減に関してはこの発明
とは直接関連しないので詳細は省略する。更に高速射影
法では、平滑化プレフィルタ係数ベクトルｓ_i （ｋ）と
近似的フィルタ係数ベクトルｚ（ｋ）を利用すること
で、図７Ａのフィルタ係数更新部２５およびフィルタ実
行部１１で行っている（ｐ＋１）Ｌの積和演算を２Ｌ＋
４ｐの積和演算量にまで低減している。この方法につい
て図７Ｂに基づいて説明を行う。即ち、プレフィルタ係
数計算部２４よりのプレフィルタ係数ベクトルｇ
（ｋ) はプレフィルタ係数平滑部２６で平滑化され、そ
の平滑化されたプレフィルタ係数ベクトルを加重係数と
して、入力ベクトルＸ（ｋ）の方向に近似的フィルタ
係数ベクトルを修正することが近似的フィルタ係数更新
部２７で行われる。入力ｘ（ｋ）の自己相関が相関計算
部２８で計算され、その出力はフィルタ演算部２９で入
力と近似的フィルタ係数ベクトルとの内積演算と、自己
相関ベクトルと平滑プレフィルタ係数ベクトルとの内積
演算と、これら内積演算の加算演算とが行われる。即
ち、プレフィルタ係数平滑部２６では式（１１）の演
算、近似的フィルタ係数更新部２７では式（１２）の演
算、相関計算部２８では式（１３）の演算、フィルタ演
算部２９では式（１４）の演算がそれぞれ行われる。な
お、フィルタ演算部２９内にはタイミングを合わせるた
めのバッファを備えている。On the other hand, in the recently proposed high-speed projection method, the product-sum calculation amount can be reduced to 2L + 20p. The high-speed projection method will be briefly described below according to Japanese Patent Laid-Open No. 7-92980. First, in the high-speed projection method, the pre-filter coefficient vector g (k) is calculated by using the linear prediction method without directly performing the inverse matrix operation of Expression (8), and thus the product sum operation of O (p ³ ) is performed. The amount is reduced to the product sum calculation amount of 16p. Since the reduction of the calculation amount is not directly related to the present invention, the details will be omitted. Furthermore, in the high-speed projection method, the smoothed pre-filter coefficient vector s _i (k) and the approximate filter coefficient vector z (k) are used to perform the filter coefficient update unit 25 and the filter execution unit 11 in FIG. 7A. The product sum operation of (p + 1) L is 2L +
It has been reduced to a product sum calculation amount of 4p. This method will be described with reference to FIG. 7B. That is, the pre-filter coefficient vector g from the pre-filter coefficient calculation unit 24
(K) is smoothed by the pre-filter coefficient smoothing unit 26, and it is approximate to correct the approximate filter coefficient vector in the direction of the input vector X (k) using the smoothed pre-filter coefficient vector as a weighting coefficient. This is performed by the filter coefficient updating unit 27. The autocorrelation of the input x (k) is calculated by the correlation calculation unit 28, and the output is calculated by the filter calculation unit 29 as the inner product of the input and the approximate filter coefficient vector and the inner product of the autocorrelation vector and the smoothed prefilter coefficient vector. The operation and the addition operation of these inner product operations are performed. That is, the pre-filter coefficient smoothing unit 26 calculates the formula (11), the approximate filter coefficient updating unit 27 calculates the formula (12), the correlation calculating unit 28 calculates the formula (13), and the filter calculating unit 29 calculates the formula (12). The calculation of 14) is performed. A buffer for adjusting the timing is provided in the filter calculation unit 29.

【００１９】 (11) ｓ_p (k) ＝「 0 ｜ ＋μ・ｇ（ｋ） ｜ｓ_p-1(k-1)」 (12) ｚ(k+1) ＝ｚ(k) ＋ｘ(k-p+1) ｓ_p (k) (13) ｒ_p-1(k+1)＝ｒ_p-1(k)＋ｘ(k+1) ｘ_p-1(k) −ｘ(k-L+1) ｘ_p-1(k-L) (14) ｙ(k+1) ＝ｘ^T (k+1) ｚ(k+1) ＋ｒ_p-1 ^T (k+1) ｓ_p-1(k) ただし、ｓ_p (k) およびｓ_p-1(ｋ−１）は、それぞ
れｐ次およびｐ−１次の平滑化係数ベクトルであって、 (15) ｓ_p (k) ＝[ ｓ₁(k)，ｓ₂(k), ... , ｓ_p (k)] (16) ｓ_p-1(k)＝[ ｓ₁(k)，ｓ₂(k), ... , ｓ_p-1(k)] と表される。式（１１）からわかるように、ｓ_p (k)
はプレフィルタ係数ベクトルｇ(k) を用いて再帰的に
計算される。また、ｚ(k) は、近似的フィルタ係数ベ
クトルであり、ｘ（ｋ−ｐ＋１）とｓ_p (k) を用いて
再帰的に計算される。また、ｒ_p-1(ｋ＋１）は、次式
で表される相関ベクトルである。[0019] _{(11) s p (k)} = "0 | + μ · g (k ) | s p-1 (k-1) " (12) z (k + 1 ) = z (k) + x (k- _{p + 1) s p (k} ) (13) r p-1 (k + 1) = r p-1 (k) + x (k + 1) x p-1 (k) -x (k-L + 1 ) x _p-1 (kL) (14) y (k + 1) = x ^T (k + 1) z (k + 1) + r _p-1 ^T (k + 1) sp _-1 (k) where s _p (k) and s _p-1 (k-1) are smoothing coefficient vectors of p-th order and p-1st order, respectively, and are (15) s _p (k) = [s ₁ (k), s ₂ (k), ..., s _p (k)] (16) s _p-1 (k) = [s ₁ (k), s ₂ (k), ..., s _p-1 (k )] Is represented. As can be seen from equation (11), s _p (k)
Is recursively calculated using the prefilter coefficient vector g (k). Furthermore, z (k) is an approximate filter coefficients vector, recursively calculated using the x (k-p + 1) and s _p (k). Further, r _p-1 (k + 1) is a correlation vector represented by the following equation.

【００２０】 (17) ｒ_p-1(k+1)＝[ ｘ^T (k+1) ｘ(k),ｘ^T (k+1) ｘ(k-1）， ... , ｘ(k+1) ^T ｘ(k-p-2)]^T このとき、数学的な証明は特開平７−９２９８０号公報
に記述されているので省略するが、図７Ｂの構成、また
は式（１１）−（１４）の演算により、入力信号ｘ
（ｋ）および誤差信号ｅ（ｋ）を用いて出力信号ｙ（ｋ
＋１）を合成することができる。なお、式（１３）（１
４）の時間パラメータがｋ＋１となっているのに図７Ｂ
ではｋとなっている。これは時刻ｋにおいて式（１３）
（１４）（１１）（１２）の順に計算することを表した
ものである。(17) r _p-1 (k + 1) = [x ^T (k + 1) x (k), x ^T (k + 1) x (k-1), ..., x (k +1) ^T x (kp-2)] ^T At this time, the mathematical proof is described in Japanese Patent Laid-Open No. 7-92980, and therefore omitted, but the configuration of FIG. 7B or equations (11)-(14) ), The input signal x
(K) and the error signal e (k), the output signal y (k
+1) can be synthesized. Note that equation (13) (1
7B even though the time parameter of 4) is k + 1.
Is k. This is formula (13) at time k.
The calculation is performed in the order of (14), (11), and (12).

【００２１】さてここで、以下の２点が指摘できる。ま
ず第１は、以上の方法の演算量に関して、式（１１）で
はｐ，式（１２）ではＬ，式（１３）では２ｐ，式（１
４）ではＬ＋ｐ，合計２Ｌ＋４ｐの積和演算量で実行で
きる。このことは、図７Ａに示した従来の射影法のフィ
ルタ係数更新部２５およびフィルタ実行部１１で行う式
（９）（１０）の演算が（ｐ＋１）Ｌの積和演算を必要
とするのに比べて大幅な低減となっている。例えばｐ＝
２０，Ｌ＝１０００の場合、従来法が２１０００である
のに対して、高速射影法では２０８０と約１／１０の演
算量となっている。The following two points can be pointed out here. First, regarding the calculation amount of the above method, p is used in Expression (11), L is used in Expression (12), 2p is used in Expression (13), and Expression (1) is used.
In 4), L + p can be executed with a total of 2L + 4p sum of products. This means that the calculation of the equations (9) and (10) performed by the filter coefficient updating unit 25 and the filter executing unit 11 of the conventional projection method shown in FIG. 7A requires the product-sum calculation of (p + 1) L. It is a significant reduction. For example p =
In the case of 20, L = 1000, the conventional method has 21000, whereas the high-speed projection method has 2080, which is a calculation amount of about 1/10.

【００２２】第２の点は、図７Ｂに示した高速射影法で
は、フィルタ係数ｈ（ｋ＋１）を陽に求めることなく
出力信号ｙ（ｋ）を合成している点が、図７Ａに示した
従来の射影法との大きな相違点である。より具体的に
は、フィルタ係数ｈ（ｋ＋１）の代わりに、近似的フ
ィルタ係数ベクトルｚ（ｋ＋１）を利用し、またｈ
（ｋ＋１）とｚ（ｋ＋１）との差を補正するために、
相関ベクトルｒ_p-1(ｋ＋１）を利用している。そし
て、このことが演算量削減に効果をなしている。The second point is that in the high speed projection method shown in FIG. 7B, the output signal y (k) is synthesized without explicitly obtaining the filter coefficient h (k + 1), which is shown in FIG. 7A. This is a big difference from the conventional projection method. More specifically, the approximate filter coefficient vector z (k + 1) is used instead of the filter coefficient h (k + 1), and h
In order to correct the difference between (k + 1) and z (k + 1),
The correlation vector r _p-1 (k + 1) is used. This is effective in reducing the calculation amount.

【００２３】[0023]

【発明が解決しようとする課題】さてここで、この高速
射影法を図６Ｂに示した filtered-x 構成に適用するこ
とを考える。filtered-x 構成では、図６Ｂに示したよ
うに、適応部２１で計算したフィルタ係数ｈ（ｋ）を
可変フィルタ１１に転送して使用する必要がある。しか
し、高速射影法では、フィルタ係数ｈ（ｋ）を陽に計
算することなく出力信号を計算することが特徴であっ
た。ただし、図７Ｂ中のｘ（ｋ），ｙ（ｋ）は、それぞ
れ図６Ｂ中のｘ’（ｋ），ｙ’（ｋ）に対応しており、
図７Ｂに示した系のみでは、図６Ｂ中のｙ（ｋ）は合成
できない点に留意が必要である。Now, consider application of this fast projection method to the filtered-x configuration shown in FIG. 6B. In the filtered-x configuration, as shown in FIG. 6B, it is necessary to transfer the filter coefficient h (k) calculated by the adaptation unit 21 to the variable filter 11 for use. However, the fast projection method is characterized in that the output signal is calculated without explicitly calculating the filter coefficient h (k). However, x (k) and y (k) in FIG. 7B correspond to x ′ (k) and y ′ (k) in FIG. 6B, respectively,
It should be noted that y (k) in FIG. 6B cannot be synthesized only with the system shown in FIG. 7B.

【００２４】従って、高速射影法を図６Ｂに示した fil
tered-x 構成に適用するためには、フィルタ係数を陽に
計算する手順を追加する必要がある。具体的には、フィ
ルタ係数ｈ（ｋ＋１）は近似的フィルタ係数ベクトル
ｚ（ｋ＋１）を修正することで次式に基づいて導出す
ることができる。 (18) ｈ(k+1) ＝ｚ(k＋１) ＋[ ｘ(k),ｘ(k-1), ... , ｘ(k-p＋２)]ｓ_p-1(k) 以上の手順を組み込んだ filtered-x 構成の高速射影法
の計算手順を図８に示す。図８において手順（ａ）〜
（ｄ）までは初期化の手順である。（ｆ）〜（ｌ）の手
順は毎時刻に繰り返し演算される。各手順での実行する
式の先頭にかっこで囲んだ数字は対応する式番号を表
す。手順（ｇ）は図７Ａ中のフィルタ実行部１１での式
（１０）の演算であり、次の手順（ｇ）は、通常、例え
ば図６Ｂ中の制御点において物理的に演算がなされる。
手順（ｉ）における式（８）中の行列Ｘ（ｋ）は式
（６）における入力ベクトルｘ（ｋ），ｘ（ｋ−
１），... の代わりに、伝達特性Ｃのフィルタを通した
信号ベクトルｘ’（ｋ），ｘ’（ｋ−１），... を構
成要素とするものであり、手順（ｋ）の式（１２）中の
入力ベクトルｘ（ｋ），ｘ（ｋ−１），... もｘ’
（ｋ），ｘ（ｋ−１），... が用いられる。式（８）を
演算する手順（ｉ）は、この式を直接計算する以外に
も、逆行列の補助定理を用いたり、また特開平７−９２
９８０号公報に示されるように線形予測係数を利用する
ことにより低演算量で計算ができる。つまり、式（８）
に相当する演算をしてもよい。Therefore, the fast projection method is shown in FIG.
In order to apply it to the tered-x configuration, it is necessary to add a procedure to explicitly calculate the filter coefficient. Specifically, the filter coefficient h (k + 1) can be derived based on the following equation by modifying the approximate filter coefficient vector z (k + 1). (18) h (k + 1) = z (k + 1) + [x (k), x (k-1), ..., x (k-p + 2)] s _p-1 (k) Figure 8 shows the calculation procedure of the built-in filtered-x high-speed projection method. In FIG. 8, procedure (a)-
The procedure up to (d) is the initialization procedure. The procedure of (f) to (l) is repeatedly calculated at every time. The number in parentheses at the beginning of the expression executed in each step represents the corresponding expression number. The procedure (g) is the calculation of the equation (10) in the filter execution unit 11 in FIG. 7A, and the next procedure (g) is usually physically calculated, for example, at the control point in FIG. 6B.
The matrix X (k) in equation (8) in step (i) is the input vector x (k), x (k- in equation (6).
Instead of 1), ..., Signal vectors x '(k), x' (k-1), ... The input vectors x (k), x (k-1), ... Also in equation (12) of
(K), x (k-1), ... Are used. The procedure (i) for calculating the equation (8) uses an inverse matrix lemma in addition to the direct calculation of the equation, and is disclosed in Japanese Patent Laid-Open No. 7-92.
As shown in Japanese Patent Laid-Open No. 980, calculation can be performed with a low calculation amount by using a linear prediction coefficient. That is, equation (8)
May be calculated.

【００２５】図９に図８に示した処理手順を機能ブロッ
クとして示す。即ち、フィルタ処理部１１で入力ｘ
（ｋ）とフィルタ係数計算部３１よりのフィルタ係数と
に手順（ｆ），つまり式（１０）が演算され、この演算
結果出力ｙ（ｋ）が制御点１４に与えられ、手順（ｇ）
が自動的になされ、その誤差ｅ（ｋ）がプレフィルタ係
数計算部２４に入力される。一方、入力ｘ（ｋ）はフィ
ルタ実行部１１の出力から制御点１４までの伝達特性Ｃ
がフィルタ１８で付与されてｘ’（ｋ）としてプレフィ
ルタ計数計算部２４に供給され、プレフィルタ計数計算
部２４で手順（ｉ），（ｊ）の式（7.２）と式（８）と
の計算がなされ、これにより得られたプレフィルタ係数
ベクトルｇ（ｋ）がプレフィルタ係数平滑部２６で手
順（ｊ）の平滑化処理を式（１１）により行い、近似的
フィルタ係数更新部２７で、その平滑化プレフィルタ係
数ベクトルと、伝達特性Ｃが付与された入力ｘ’（ｋ）
とにより手順（ｋ）の式（１２）を実行して近似的フィ
ルタ係数ベクトルを順次更新し、得られた近似的フィル
タ係数ベクトルと、ｘ’（ｋ）と、平滑化プレフィルタ
係数ベクトルとを用いてフィルタ係数計算部３１で手順
（ｌ）の式（１８）を演算してフィルタ係数を求める。FIG. 9 shows the processing procedure shown in FIG. 8 as a functional block. That is, the input x in the filter processing unit 11
The procedure (f), that is, the equation (10) is calculated for (k) and the filter coefficient from the filter coefficient calculation unit 31, this calculation result output y (k) is given to the control point 14, and the procedure (g)
Is automatically performed, and the error e (k) is input to the pre-filter coefficient calculation unit 24. On the other hand, the input x (k) is the transfer characteristic C from the output of the filter execution unit 11 to the control point 14.
Is added by the filter 18 and supplied to the pre-filter count calculation unit 24 as x ′ (k), and in the pre-filter count calculation unit 24, the formulas (7.2) and (8) of the steps (i) and (j) are calculated. And the pre-filter coefficient vector g (k) thus obtained is subjected to the smoothing process of the procedure (j) in the pre-filter coefficient smoothing unit 26 by the equation (11), and the approximate filter coefficient updating unit 27 Then, the smoothed pre-filter coefficient vector and the input x ′ (k) to which the transfer characteristic C is added
And by successively updating the formula approximate filter coefficient vector by running (12) of the steps (k), and the resulting approximate filter coefficient vector, x 'and (k), and a smoothing pre-filter coefficient vector Using the filter coefficient calculation unit 31, the expression (18) of the procedure (1) is calculated to obtain the filter coefficient.

【００２６】しかし、この手順では、従来の高速射影法
に加えて手順（１）式（１８）の演算が必要となるた
め、新たに（ｐ−１）Ｌ回の積和演算が必要となり、高
速射影法の特徴である低演算量の利点は失われる。この
ように、従来技術において演算量の大幅な増加を伴うこ
となく高速射影法を filtered-x 構成に適用することは
できない、という問題点があった。However, in this procedure, in addition to the conventional high-speed projection method, the calculation of the procedure (1) expression (18) is required, and therefore (p-1) L product-sum calculations are newly required, The advantage of low computational complexity, which is a feature of the fast projection method, is lost. As described above, the conventional technique has a problem that the fast projection method cannot be applied to the filtered-x configuration without significantly increasing the amount of calculation.

【００２７】[0027]

【課題を解決するための手段】この発明では、上記の問
題点に対して相互相関ベクトルを導入することで解決を
図り、演算量の増加を伴うことなく、高速射影法を fil
tered-x 構成に適用することを可能とするものである。
即ち、この発明によれば、入力信号ベクトルと、これに
伝達特性Ｃを付与した信号ベクトルとの相互相関ベクト
ルを求め、この相互相関ベクトルと平滑化されたプレフ
ィルタ係数ベクトルとの内積演算と、入力信号ベクトル
と近似的フィルタ係数ベクトルとの内積演算とを行い、
これら内積演算結果を合成して出力信号とする。In the present invention, the above problem is solved by introducing a cross-correlation vector, and the high-speed projection method is applied to the high-speed projection method without increasing the calculation amount.
It is applicable to tered-x configuration.
That is, according to the present invention, the cross-correlation vector between the input signal vector and the signal vector having the transfer characteristic C added thereto is obtained, and the inner product operation of the cross-correlation vector and the smoothed pre-filter coefficient vector is performed. Do the inner product operation of the input signal vector and the approximate filter coefficient vector,
These inner product calculation results are combined and used as an output signal.

【００２８】まず最初にこの発明の原理を説明する。図
６Ｂに示した filtered-x 構成の系における適応部２１
に、図７Ｂで示した高速射影法を適用したと考える。こ
のとき、信号の表示としては、図７Ｂにおけるｘ（ｋ）
は図６Ｂにおいては、ｘ’（ｋ）となり、また、ｙ
（ｋ）はｙ’（ｋ）となる。図６Ｂにおける可変フィル
タ１１の出力ｙ（ｋ）は、 （19） ｙ(k) ＝ｘ^T (k) ｈ(k) と表される。一方、フィルタ係数ｈ（ｋ）は、式（１
８）の時刻ｋをｋ−１として、ｘをｘ’とおいて、 （20) ｈ(k) ＝ｚ(k) ＋[ ｘ'(ｋ−1), ｘ'(ｋ−２）， ... , ｘ'(ｋ−ｐ＋１)]ｓ_p-1(ｋ−１） と表される。式（２０）を式（１９）に代入すれば、 （21) ｙ(k) ＝ｘ^T (k) ｚ(k) ＋ｘ^T (k)[ｘ'(ｋ−1), ｘ'(ｋ−２）， ... , ｘ'(ｋ−ｐ＋１)]ｓ_p-1(ｋ−１） と表される。ここで、相互相関ベクトルｐ（ｋ）を次
式で定義する。First, the principle of the present invention will be described. The adaptation unit 21 in the system having the filtered-x configuration shown in FIG. 6B.
It is considered that the fast projection method shown in FIG. At this time, the signal is displayed as x (k) in FIG. 7B.
Becomes x ′ (k) in FIG. 6B, and y
(K) becomes y '(k). The output y of the variable filter 11 in FIG. 6B (k) is represented as (19) y (k) = x T (k) h (k). On the other hand, the filter coefficient h (k) is calculated by the equation (1
8) time k is k−1, x is x ′, and (20) h (k) = z (k) + [x ′ (k−1), x ′ (k-2) ,. ., x '(k-p + 1)] s _p-1 (k-1). Substituting equation (20) into equation (19), (21) y (k) = x T (k) z (k) + x T (k) [x '(k-1), x' (k- 2), ..., X '(k-p + 1)] s _p-1 (k-1). Here, the cross-correlation vector p (k) is defined by the following equation.

【００２９】 （22）ｐ(k) ＝（ｘ^T (k)[ｘ'(ｋ−1), ｘ'(ｋ−２）， ... , ｘ'(ｋ−ｐ＋１)]) ^T すると、式（２１）は、次式のように表される。 （23） ｙ(k) ＝ｘ^T (k) ｚ(k) ＋ｐ^T (k) ｓ_p-1(ｋ−１） ここで、ｐ（ｋ）はｐ（ｋ−１）より再帰的に計算
できることを示す。まず、式（３．１）に示した定義よ
り、 （24）ｘ'(k-1)＝[ ｘ'(k-1), ｘ'(k-2), ... , ｘ'(k-L)] ^T ｘ'(k-2)＝[ ｘ'(k-2), ｘ'(k-3), ... , ｘ'(k-L-1)] ^T ・・・ ｘ'(k-p+1)＝[ ｘ'(k-p+1)，ｘ'(k-p), ... , ｘ'(k-p-L+2)] ^T である。次に、ｐ−１個の入力信号を要素とした入力信
号ベクトルｘ’_p-1(ｋ−１）を次式で定義する。[0029] (22) p (k) = (x T (k) [x '(k-1), x' (k-2), ..., x '(k-p + 1)]) T Then, Expression (21) is expressed as the following expression. (23) y (k) = x T (k) z (k) + p T (k) s p-1 (k-1) where, p (k) is computed recursively from p (k-1) Show what you can do. First, according to the definition shown in Expression (3.1), (24) x '(k-1) = [x' (k-1), x '(k-2), ..., x' (kL )] ^T x '(k-2) = [x' (k-2), x '(k-3), ..., x' (kL-1)] ^T ... x '(k-p +1) = [x '(k-p + 1), x' (kp), ..., X '(kp-L + 2)] ^T. Next, an input signal vector _x'p-1 (k-1) having _p-1 input signals as elements is defined by the following equation.

【００３０】 （25）ｘ' _p-1(k-1)＝[ ｘ'(k-1), ｘ'(k-2), ... , ｘ'(k-p+1)] ^T これを用いれば、次式の関係が成立する。 （26) [ ｘ'(k-1), ｘ'(k-2), ... , ｘ'(k-p+1)] ＝[ ｘ' _p-1(k-1)，ｘ' _p-1(k-2), ... , ｘ' _p-1(k-L)] ^T これを式（２２）に代入すれば （27）ｐ(k) ＝（ｘ^T (k)[ｘ' _p-1(k-1)，ｘ' _p-1(k-2), ... ,ｘ' _p-1(k-L)] ^T ）^T ＝ｘ(k) ｘ' _p-1(k-1)＋ｘ(k-1) ｘ' _p-1(k-2) ...,＋x(k-L+2)ｘ' _p-1(k-L+1)＋x(k-L+1)ｘ' _p-1(k-L) の関係が得られる。式（２７）のｋにｋ−１を代入すれ
ば、 （28）ｐ(k-1) ＝ｘ(k-1) ｘ' _p-1(k-2)＋x(k-2)ｘ' _p-1(k-3) ...,＋x(k-L+1)ｘ' _p-1(k-L)＋ｘ(k-L) ｘ' _p-1(k-L-1) 式（２８）と式（２７）とを比較すれば、次の再帰式が
得られる。(25) _x'p-1 (k-1) = [x '(k-1), x' (k-2), ..., x '(k-p + 1)] ^T Using, the following relationship holds. (26) [x '(k-1), x' (k-2), ..., x '(k-p + 1)] = [x' _p-1 (k-1), _{x'p -1} (k-2), ..., x'p _-1 (kL)] ^T Substituting this into equation (22) gives (27) p (k) = (x ^T (k) [x ' _{p -1} (k-1), _x'p-1 (k-2), ..., x'p _-1 (kL)] ^T ) ^T = x (k) _x'p-1 (k-1) + X (k-1) _x'p-1 (k-2) ..., + x (k-L + 2) x ' _p-1 (k-L + 1) + x (k-L + 1) x' The relation of _p-1 (kL) is obtained. Substituting k-1 for k in equation (27), (28) p (k-1) = x (k-1) _x'p-1 (k-2) + x (k-2) x ' _{p -1 (k-3) ...,} + x (kL + 1) x 'p-1 (kL) + x (kL) x' p-1 (kL-1) equation (28) and (27) Comparing with, the following recursive expression is obtained.

【００３１】 （29）ｐ(k) ＝ｐ(k-1)+x(k)ｘ' _p-1(k-1)−ｘ(k-L) ｘ' _p-1(k-L-1) この再帰式（２９）に基づいて相互相関ベクトルｐ
（ｋ）を計算し、式（２３）を実行すれば、陽にフィル
タ係数ｈ（ｋ）を計算することなくfiltered-x構成の
出力信号ｙ（ｋ）を算出することができる。[0031] (29) p (k) = p (k-1) + x (k) x 'p-1 (k-1) -x (kL) x' p-1 (kL-1) this recursive formula The cross-correlation vector p based on (29)
By calculating (k) and executing Expression (23), the output signal y (k) of the filtered-x configuration can be calculated without explicitly calculating the filter coefficient h (k).

【００３２】[0032]

【発明の実施の形態】課題を解決するための手段項に示
したこの発明の原理から、この発明の filtered-x 構成
の高速射影法の手順は図１に示すようになる。図１にお
いて、初期化手順（ａ）〜（ｄ）において、相互相関ベ
クトルｐ（ｏ）＝（Ｘ^T （０）［Ｘ’（−１），
Ｘ’（−２），... ，Ｘ’（−ｐ＋１）］）^T と
し、平滑化フィルタ係数ベクトルｓ_p-1 （０），誤差
ベクトルｅ_p-1(０）をそれぞれ［０，０，... ，０］
^T ，近似的フィルタ係数ベクトルｚ₁(１）を［０，
０，... ，０］とする。DESCRIPTION OF THE PREFERRED EMBODIMENTS From the principle of the present invention described in the means for solving the problems, the procedure of the fast projection method of the filtered-x configuration of the present invention is as shown in FIG. In FIG. 1, in the initialization procedures (a) to (d), the cross-correlation vector p (o) = (X ^T (0) [X ′ (− 1),
X ′ (− 2), ..., X ′ (− p + 1)]) ^T , the smoothing filter coefficient vector s _p-1 (0) and the error vector e _p-1 (0) are respectively [0, 0 , ..., 0]
^T , the approximate filter coefficient vector z ₁ (1) is [0,
0, ..., 0].

【００３３】次に時刻をｋ＝１とし、手順（ｅ）に移
る。この手順（ｅ）から手順（ｋ）は毎時刻に繰り返し
演算され、つまり手順（ｅ）〜（ｋ）の演算が終了する
と、時刻ｋを＋１して手順（ｅ）〜（ｋ）の演算を再び
実行することを繰り返す。この繰り返し演算において、
手順（ｅ）では式（２９）により入力ｘ（ｋ）のベクト
ルと、伝達特性Ｃが付与された入力ｘ’（ｋ）のベクト
ルとの相互相関ベクトルｐ（ｋ）を計算し、手順
（ｆ′）で式（２３），つまり、入力ベクトルｘ
^T （ｋ）と近似的フィルタ係数ベクトルｚ（ｋ）との
内積演算と、相互相関ベクトルｐ（ｋ）^T と、平滑化
プレフィルタ係数ベクトルｓ_p-1(ｋ−１）との内積と
を加算して出力ｙ（ｋ）を得る。Next, the time is set to k = 1, and the procedure moves to step (e). The procedure (e) to the procedure (k) are repeatedly calculated every time, that is, when the calculations of the procedures (e) to (k) are completed, the time k is incremented by 1 and the calculations of the procedures (e) to (k) are performed. Repeat to execute again. In this iterative operation,
In step (e), the cross-correlation vector p (k) between the vector of the input x (k) and the vector of the input x ′ (k) to which the transfer characteristic C is added is calculated by the equation (29), and the step (f ′) In equation (23), that is, the input vector x
The inner product operation of ^T (k) and the approximate filter coefficient vector z (k) and the inner product of the cross-correlation vector p (k) ^T and the smoothed pre-filter coefficient vector sp _-1 (k-1) Addition is performed to obtain the output y (k).

【００３４】手順（ｇ）は図８の場合と同様に、出力ｙ
（ｋ）が制御点１４に与えられて、所望信号との誤差を
得、手順（ｈ）乃至（ｋ）は図８中のそれと同様であ
る。ただし、手順（ｉ）における行列Ｘ（ｋ）は式
（６）における入力ベクトルｘ（ｋ），ｘ（ｋ−
１），... の代わりに特性Ｃのフィルタを通した信号ベ
クトルｘ’（ｋ），ｘ’（ｋ−１），... を構成要素と
したものであり、また手順（ｋ）においてもｘ（ｋ）の
代わりにｘ’（ｋ）を用いる。つまり、この図１では図
８中のフィルタ係数ｈ（ｋ＋１）を陽に計算する手順
（ｌ）と、そのフィルタ係数ｈ（ｋ）を用いて出力
ｙ’（ｋ）を計算する手順（ｆ）の代わりに、相互相関
ベクトルｐ（ｋ）を計算する手順（ｅ），およびその
相互相関ベクトルを用いて出力ｙ（ｋ）を計算する手順
（ｆ′）が用いられている。なお、手順（ｉ）は図８の
場合と同様に、逆行列の補助定理を用いたり、線形予測
係数を利用することにより、低演算量で計算してもよ
い。The procedure (g) is the same as in the case of FIG.
(K) is given to the control point 14 to obtain an error from the desired signal, and steps (h) to (k) are the same as those in FIG. However, the matrix X (k) in step (i) is the input vector x (k), x (k- in equation (6).
1), ..., Instead of 1), ..., Signal vectors x '(k), x' (k-1), ... Also uses x '(k) instead of x (k). That is, in FIG. 1, the procedure (1) for explicitly calculating the filter coefficient h (k + 1) in FIG. 8 and the procedure (f) for calculating the output y ′ (k) using the filter coefficient h (k). Instead of, the procedure (e) for calculating the cross-correlation vector p (k) and the procedure (f ') for calculating the output y (k) using the cross-correlation vector are used. Note that the procedure (i) may be calculated with a low calculation amount by using the lemma of the inverse matrix or by using the linear prediction coefficient, as in the case of FIG.

【００３５】従来技術である図８における手順（ｆ）
（１）を実行するのに必要な演算量は約ｐＬである。こ
れに対して図１における手順（ｅ）（ｆ′）を実行する
のに必要な演算量は約Ｌ＋３ｐと大幅な低減となってい
る。例えばｐ＝２０，Ｌ＝１０００の場合、図８におけ
る手順（ｆ）（１）を実行するのに必要な演算量は約２
００００，これに対して図１における手順（ｅ）
（ｆ′）を実行するのに必要な演算量は約１０６０と、
約１／２０に低減している。Procedure (f) in FIG. 8 which is a conventional technique.
The amount of calculation required to execute (1) is about pL. On the other hand, the amount of calculation required to execute the steps (e) and (f ') in FIG. 1 is about L + 3p, which is a large reduction. For example, when p = 20 and L = 1000, the calculation amount required to execute the procedure (f) (1) in FIG. 8 is about 2
0000, in contrast, the procedure (e) in FIG.
The amount of calculation required to execute (f ') is about 1060,
It has been reduced to about 1/20.

【００３６】一方、以上のこの発明を従来の高速射影法
と比べた場合、高速射影法における式（１３）（１４）
の演算を式（２３）（２９）の演算（手順（ｅ）
（ｆ′））に入れ替えたものとなっている。これらの式
の演算量は等しいので、この発明は演算量を増加させる
ことなく、従来の高速射影法を filtered-x 構成に適用
可能とする方法と言うこともできる。On the other hand, when the present invention described above is compared with the conventional high-speed projection method, equations (13) and (14) in the high-speed projection method are used.
Of the equations (23) and (29) (procedure (e)
(F ')) has been replaced. Since the calculation amounts of these equations are equal, the present invention can be said to be a method that makes it possible to apply the conventional high-speed projection method to the filtered-x configuration without increasing the calculation amount.

【００３７】図１の処理手順を機能的に示すと図２に示
すようになる。図２に図９と対応する部分に同一符号を
付けてあり、図９と異なる点は、図９中からフィルタ実
行部１１，フィルタ係数計算部３１が除去され、代わっ
て相関計算部３３が設けられ、入力ベクトルｘ（ｋ）
と、フィルタ１８からの伝達特性が付与された信号ベク
トルｘ’（ｋ）とから手順（ｅ）の式（２９）により
相互相関ベクトルｐ（ｋ）が計算され、かつフィルタ
実行部３４でその相互相関ベクトルｐ（ｋ）と、平滑
化プレフィルタ係数ベクトルｓ_p-1(ｋ−１）との内積
演算と、入力ベクトルｘ（ｋ）と近似的フィルタ係数ベ
クトルｚ（ｋ）との内積演算と、これら内積演算の結
果を合成して出力信号ｙ（ｋ）を得ること、つまり手順
（ｆ′）の式（２３）の計算がなされる。The functional procedure of the processing procedure shown in FIG. 1 is as shown in FIG. 2 are denoted by the same reference numerals as those in FIG. 9, and the difference from FIG. 9 is that the filter execution unit 11 and the filter coefficient calculation unit 31 are removed from FIG. 9, and a correlation calculation unit 33 is provided instead. Input vector x (k)
And the signal vector x ′ (k) to which the transfer characteristic from the filter 18 is added, the cross-correlation vector p (k) is calculated by the equation (29) of the procedure (e), and the filter execution unit 34 calculates the cross-correlation vector a correlation vector p (k), and the inner product computation of the smoothed pre-filter coefficient vector _{s p-1 (k-1} ), approximate filter coefficients base and the input vector x (k)
The inner product operation with the cutler z (k) and the result of these inner product operations are combined to obtain the output signal y (k), that is, the calculation of the equation (23) of the procedure (f ′) is performed.

【００３８】この図２における動作は次のように行われ
る。まず、入力信号ｘ（ｋ）をフィルタ１８を通して信
号ｘ’（ｋ）を合成する。次にプレフィルタ係数計算部
２４において、信号ｘ’（ｋ）と誤差ｅ（ｋ）を用い
て、式（７．２）（８）に基づいてプレフィルタ係数ベ
クトルｇ（ｋ）を計算する。なお、ｇ（ｋ）は、式
（８）を直接計算する以外にも、逆行列の補助定理を用
いたり、また特開平７−９２９８０号公報に示されるよ
うに線形予測係数を利用することにより低演算量で計算
ができる。The operation in FIG. 2 is performed as follows. First, the input signal x (k) is passed through the filter 18 to synthesize the signal x '(k). Next, in pre-filter coefficient calculating unit 24, using the signal x '(k) and error e a (k), pre-filter coefficients base on the basis of the equation (7.2) (8)
Calculate the ktor g (k). Note that g (k) can be calculated by using the inverse matrix lemma other than directly calculating equation (8), or by using a linear prediction coefficient as disclosed in JP-A-7-92980. Can be calculated with a low amount of calculation.

【００３９】次にプレフィルタ係数ベクトルｇ（ｋ）
はプレフィルタ係数平滑部２６に転送され、式（１１）
に基づいて平滑化プレフィルタ係数ベクトルｓ
_p （ｋ）を合成する。次にｓ_p （ｋ）の第ｐ番目の要
素であるｓ_p （ｋ）を近似的フィルタ係数更新部２７に
転送して、同時に入力されるｘ’（ｋ）と共に式（１
２）を用いて近似的フィルタ係数ベクトルｚ（ｋ＋
１）を計算する。Next, the pre-filter coefficient vector g (k)
Is transferred to the pre-filter coefficient smoothing unit 26, and equation (11)
Smoothing prefilter coefficient vector s
Synthesize _p (k). Next, s _p (k), which is the p-th element of s _p (k), is transferred to the approximate filter coefficient updating unit 27, and x ′ (k) that is input at the same time is used in the equation (1
2) using the approximate filter coefficient vector z (k +
Calculate 1).

【００４０】ここで時間パラメータを一つ更新する。即
ち、ｋをｋ＋１とする。そして、時間ｋ＋１における入
力ｘ（ｋ＋１）とフィルタ出力ｘ’（ｋ＋１）を用い
て、相互相関計算部３３において、式（２９）に基づい
て、相互相関ベクトルｐ（ｋ＋１）を計算する。ｐ
（ｋ＋１）は、時刻ｋにおいて計算された近似的フィル
タ係数ベクトルｚ（ｋ＋１）および平滑化プレフィル
タ係数ベクトルｓ_p （ｋ），および入力ｘ（ｋ＋
１）と共にフィルタ実行部３４に入力される。フィルタ
実行部３４においては式（２３）に基づいて出力ｙ（ｋ
＋１）を合成する。なお、図１の各手順における時刻ｋ
と、図２における時刻ｋとを合わせるために、例えばフ
ィルタ実行部３４中に近似的フィルタ係数ベクトルｚ
（ｋ），平滑化プレフィルタ係数ベクトルｓ_p （ｋ）
を一時保持するバッファが設けられることになる。ただ
図２はあくまでも機能構成を示すものであり、これらの
各部は通常はそれぞれハードウェアとして独立に設けら
れるのではなく、例えば一つ乃至、複数のＤＳＰ（ディ
ジタルシグナルプロセッサ）により全体の処理がなされ
る。Here, one time parameter is updated. That is, k is set to k + 1. Then, using the input x (k + 1) and the filter output x ′ (k + 1) at time k + 1, the cross-correlation calculation unit 33 calculates the cross-correlation vector p (k + 1) based on Expression (29). p
(K + 1) is the approximate filter coefficient vector z (k + 1) and the smoothing prefill calculated at time k.
Data coefficient vector s _p (k), and the input x (k +
It is input to the filter execution unit 34 together with 1). In the filter execution unit 34, the output y (k
+1) is synthesized. The time k in each procedure in FIG.
And the time k in FIG. 2, the approximate filter coefficient vector z
(K), smoothed pre-filter coefficient vector s _p (k)
Will be provided with a buffer for temporarily holding. However, FIG. 2 merely shows the functional configuration, and each of these units is not normally provided independently as hardware, but the entire processing is performed by, for example, one or a plurality of DSPs (digital signal processors). It

【００４１】図２に示したこの発明はｘ（ｋ），ｅ
（ｋ）を入力としてｙ（ｋ）を出力する。一方、図５Ｂ
に示した騒音制御系における適応フィルタ１１も同様に
ｘ（ｋ),ｅ（ｋ）を入力としてｙ（ｋ）を出力するもの
である。従って、図５Ｂにおける適応フィルタ１１を直
接、図２に示した構成で置き換えることで、この発明を
騒音制御系に適用することが実現できる。The present invention shown in FIG. 2 has x (k), e
Y (k) is output with (k) as an input. On the other hand, FIG.
Similarly, the adaptive filter 11 in the noise control system shown in FIG. 6 also inputs x (k) and e (k) and outputs y (k). Therefore, the present invention can be applied to the noise control system by directly replacing the adaptive filter 11 in FIG. 5B with the configuration shown in FIG.

【００４２】図３にこの発明を振動制御系に適用した場
合を示す。この制御系は以下のように動作する。壁面４
１に固定された板４２は第１の加震機４３によって加震
される。この加震は望ましいものではないものとする。
例えば、この加震機４３は発電機であって、発電の際に
発生する振動が板４２に伝わって板４２を振動させるも
のとする。この制御系の目的は第２の加震機４４を用い
て、板４２の振動を減衰させることである。そのために
は、板４２に第１の振動検出器４５を取り付けて、この
出力を誤差信号（パワーを最小にすべき信号）ｅ（ｋ）
として取り出し、これを適応フィルタ４６に入力する。
また、第１の加震機４３には第２の振動検出器４７を取
り付けて振動の原信号ｘ（ｋ）を取り出し、適応フィル
タ４６に入力する。適応フィルタ４６は誤差信号ｅ
（ｋ），即ち板４２の振動が最小となるような出力信号
ｙ（ｋ）を合成して第２の加震機４４に供給する。FIG. 3 shows a case where the present invention is applied to a vibration control system. This control system operates as follows. Wall 4
The plate 42 fixed to 1 is shaken by the first shaker 43. This shaking is not desirable.
For example, it is assumed that the shaker 43 is a generator and that vibration generated during power generation is transmitted to the plate 42 and vibrates the plate 42. The purpose of this control system is to dampen the vibration of the plate 42 using the second shaker 44. For that purpose, the first vibration detector 45 is attached to the plate 42, and this output is used as an error signal (a signal whose power should be minimized) e (k).
And inputs it to the adaptive filter 46.
Further, a second vibration detector 47 is attached to the first shaker 43 to take out the original vibration signal x (k) and input it to the adaptive filter 46. The adaptive filter 46 outputs the error signal e
(K), that is, the output signal y (k) that minimizes the vibration of the plate 42 is synthesized and supplied to the second shaker 44.

【００４３】この制御系は第１の振動検出器４５の設置
された点が制御点である。そして、フィルタ４６の出力
端から制御点までの間には伝達系（第２の加震機４４の
電気−振動変換特性など）が含まれている。従って、適
応フィルタ４６は filtered-x 構成を持つ必要がある。
そこで、この制御系にはこの発明が適用できる。図２に
示したように、この発明は入力信号ｘ（ｋ）および誤差
信号ｅ（ｋ）を用いて出力信号ｙ（ｋ）を合成する fil
tered-x 構成の適応フィルタである。従って、図２の系
を図３の適応フィルタ４６としてそのまま利用すること
で、この発明を振動制御系に適用することができる。In this control system, the point where the first vibration detector 45 is installed is the control point. Further, a transmission system (electric-vibration conversion characteristic of the second shaker 44, etc.) is included between the output end of the filter 46 and the control point. Therefore, the adaptive filter 46 needs to have a filtered-x configuration.
Therefore, the present invention can be applied to this control system. As shown in FIG. 2, the present invention uses the input signal x (k) and the error signal e (k) to synthesize the output signal y (k) fil
It is an adaptive filter with a tered-x configuration. Therefore, the present invention can be applied to the vibration control system by directly using the system of FIG. 2 as the adaptive filter 46 of FIG.

【００４４】[0044]

【発明の効果】図４はこの発明の効果を示すために行っ
た騒音制御シミュレーションの結果を示す図である。実
験は残響時間が約４００msの室内データを用いて行っ
た。図４の縦軸は制御点における誤差信号レベルを表
し、横軸は時間をサンプル数で表した。ただしサンプリ
ング周波数は１kHz とした。適応制御法としては、この
発明と同程度の演算量を持つ学習同定法（従来技術）、
および射影次数ｐ＝２，４，８としたこの発明を用い
た。室内の伝達特性は時刻ｋ＝１０００で変化させた。
従来技術である学習同定法を用いた場合には、誤差の収
束が遅く、誤差レベルを−２０dBに減衰させるために８
０００サンプル程度（ｋ＝１０００から９０００まで）
必要としている。これに対してこの発明でｐ＝８の場合
には、５００サンプル程度（ｋ＝１５００まで）で−２
０dBの誤差レベルに到達している。このように、この発
明では従来技術に比べて約１５倍の収束速度の向上が得
られた。FIG. 4 is a diagram showing the result of the noise control simulation performed to show the effect of the present invention. The experiment was performed using room data with a reverberation time of about 400 ms. The vertical axis of FIG. 4 represents the error signal level at the control point, and the horizontal axis represents the time by the number of samples. However, the sampling frequency was 1 kHz. As an adaptive control method, a learning identification method (prior art) having a calculation amount equivalent to that of the present invention,
The invention was used with the projection orders p = 2, 4, and 8. The transfer characteristic in the room was changed at time k = 1000.
When the learning identification method, which is a conventional technique, is used, the convergence of the error is slow, and in order to reduce the error level to −20 dB, 8
About 000 samples (k = 1000 to 9000)
In need of. On the other hand, in the case of p = 8 in the present invention, about 500 samples (up to k = 1500) is −2.
The error level of 0 dB has been reached. As described above, in the present invention, the convergence speed is improved about 15 times as compared with the conventional technique.

【００４５】以上説明したように、この発明では、相互
相関ベクトルｐ（ｋ＋１）を導入し、これを利用する
ことにより、フィルタ係数ｈ（ｋ＋１）を陽に計算す
ることなく filtered-x 構成の出力信号ｙ（ｋ）を算出
することができる。そして、その結果、この発明は演算
量を増加させることなく、高速射影法を filtered-x構
成に適用可能とする。その結果、従来技術ではＬＭＳ法
などの収束の遅い手法が利用されていた騒音制御装置に
対して、演算量を増加させることなく、即ち、同一のハ
ードウェア規模で収束の速い射影法の利用が可能とな
り、高速な騒音低減効果や振動低減効果の実現が達成さ
れる。As described above, according to the present invention, the cross-correlation vector p (k + 1) is introduced, and by utilizing this, the output of the filtered-x configuration can be obtained without explicitly calculating the filter coefficient h (k + 1). The signal y (k) can be calculated. As a result, the present invention enables the fast projection method to be applied to the filtered-x configuration without increasing the amount of calculation. As a result, it is possible to use a projection method with a fast convergence without increasing the amount of calculation, that is, with respect to a noise control device in which a slow convergence method such as the LMS method is used in the related art. It becomes possible and high-speed noise reduction effect and vibration reduction effect are achieved.

【００４６】なお、この発明は filtered-x 構成と類似
の構成、即ち、適応部に対する入力信号とフィルタ部に
対する入力信号が異なるようなあらゆる構成の適応制御
装置にも直接適応が可能である。The present invention can be directly applied to an adaptive control device having a configuration similar to the filtered-x configuration, that is, any configuration in which the input signal to the adaptive unit and the input signal to the filter unit are different.

【図面の簡単な説明】[Brief description of drawings]

【図１】この発明の方法の実施例を示す流れ図。FIG. 1 is a flow chart showing an embodiment of the method of the present invention.

【図２】この発明の方法を機能的に示すブロック図。FIG. 2 is a block diagram functionally showing the method of the present invention.

【図３】この発明の方法を適応することができる振動制
御系の例を示す図。FIG. 3 is a diagram showing an example of a vibration control system to which the method of the present invention can be applied.

【図４】この発明の方法と従来方法とによる誤差の収束
状態を計算機シミュレーションの結果を示すグラフ。FIG. 4 is a graph showing the results of computer simulation of the convergence state of the error between the method of the present invention and the conventional method.

【図５】Ａは適応制御系を示すブロック図、Ｂは適応制
御系と騒音制御系に適用した例を示すブロック図、Ｃは
その等価構成図である。5A is a block diagram showing an adaptive control system, FIG. 5B is a block diagram showing an example applied to an adaptive control system and a noise control system, and C is an equivalent configuration diagram thereof.

【図６】Ａは図５Ｃに示した制御系の順序を入れ替えた
図、Ｂは図５Ｂ中の伝達系の特性を入力信号に付与し、
その信号と誤差信号とでフィルタ係数を制御する構成を
示すブロック図である。6A is a diagram in which the order of the control system shown in FIG. 5C is changed, and B is the characteristic of the transmission system in FIG.
It is a block diagram which shows the structure which controls a filter coefficient with the signal and an error signal.

【図７】Ａは従来の射影法を適用した制御方法を示す機
能ブロック図、Ｂは従来の高速射影法を適用した制御方
法を示す機能ブロック図である。7A is a functional block diagram showing a control method applying a conventional projection method, and FIG. 7B is a functional block diagram showing a control method applying a conventional high-speed projection method.

【図８】filtered-x 構成に高速射影法を適用した処理
手順を示す図。FIG. 8 is a diagram showing a processing procedure in which a fast projection method is applied to a filtered-x configuration.

【図９】図８の手順を機能的に示すブロック図。9 is a block diagram functionally showing the procedure of FIG.

───────────────────────────────────────────────────── フロントページの続き (56)参考文献 特開 平７−92980（ＪＰ，Ａ) 田中雅史，外２名著，Ｆｉｌｔｅｒｅ ｄ−ｘ 構成に対する高速射影アルゴリ ズムの適用について，日本音響学会平成 ７年度秋季研究発表会講演論文集−Ｉ −，社団法人日本音響学会，1995年 ９ 月27日，ｐ．545−546 (58)調査した分野(Int.Cl.⁷，ＤＢ名) G06F 17/00 - 17/18 G05B 1/00 - 7/04 G05B 11/00 - 13/04 G05B 17/00 - 17/02 G05B 21/00 - 21/02 G10K 11/00 - 13/00 H03H 15/00 - 15/02 H03H 17/00 - 21/00 ─────────────────────────────────────────────────── ─── Continuation of the front page (56) Reference Japanese Patent Laid-Open No. 7-92980 (JP, A) Masafumi Tanaka, 2 authors, Application of high-speed projection algorithm to Filtere dx construction, Acoustical Society of Japan, 1995 Proceedings of Autumn Research Conference-I-, The Acoustical Society of Japan, September 27, 1995, p. 545-546 (58) Fields investigated (Int.Cl. ⁷ , DB name) G06F 17/00-17/18 G05B 1/00-7/04 G05B 11/00-13/04 G05B 17/00-17 / 02 G05B 21/00-21/02 G10K 11/00-13/00 H03H 15/00-15/02 H03H 17/00-21/00

Claims (3)

(57)【特許請求の範囲】(57) [Claims] 【請求項１】 制御点における誤差信号のパワーを最小
とするように、 可変フィルタの特性を適応的に制御して上記可変フィル
タの出力信号を制御する方法であって、 上記可変フィルタの出力端と、上記制御点までの伝達特
性を上記可変フィルタの入力信号に付与した信号に基づ
いて制御を行う方法において、 上記伝達特性を入力信号に付与した信号ベクトルと、上
記入力信号ベクトルとの相互相関ベクトルを求める第１
ステップと、 上記伝達特性を入力信号に付与した信号ベクトルの自己
相関行列の逆行列と、上記誤差信号ベクトルとの積に相
当する演算を行ってプレフィルタ係数ベクトルを求める
第２ステップと、 上記プレフィルタ係数ベクトルを平滑化する第３ステッ
プと、 上記平滑化されたプレフィルタ係数ベクトルを加重係数
として、上記伝達特性を入力信号に付与した信号ベクト
ルの方向に近似的フィルタ係数ベクトルを修正する第４
ステップと、 上記相互相関ベクトルと上記平滑化されたプレフィルタ
係数ベクトルとの内積演算を行い、かつ上記入力信号ベ
クトルと、上記近似的フィルタ係数ベクトルとの内積演
算を行い、これらを合成して出力信号とする第５ステッ
プと、 を有することを特徴とする適応的制御方法。1. A method of adaptively controlling characteristics of a variable filter to control an output signal of the variable filter so as to minimize power of an error signal at a control point, the output terminal of the variable filter. And a method of controlling the transfer characteristic up to the control point based on a signal added to the input signal of the variable filter, wherein a cross-correlation between the signal vector added with the transfer characteristic to the input signal and the input signal vector First to find a vector
A second step of obtaining a pre-filter coefficient vector by performing an operation corresponding to the product of the inverse matrix of the autocorrelation matrix of the signal vector in which the transfer characteristic is given to the input signal and the error signal vector; A third step of smoothing the filter coefficient vector, and a fourth step of modifying the approximate filter coefficient vector in the direction of the signal vector in which the transfer characteristic is given to the input signal, using the smoothed pre-filter coefficient vector as a weighting coefficient.
Step, inner product operation of the cross-correlation vector and the smoothed pre-filter coefficient vector is performed, and inner product operation of the input signal vector and the approximate filter coefficient vector is performed, and these are combined and output. A fifth step of making a signal, and the adaptive control method. 【請求項２】 上記入力信号をｘ（ｋ）、上記伝達特性
を付与した入力信号ベクトルをｘ’_p-1(ｋ−１）とす
ると、上記第１ステップでの上記相互相関ベクトル ｐ
（ｋ）は ｐ(k) ＝ｐ(k-1)+ｘ(k) ｘ’_p-1(k-1)−ｘ(k-L) ｘ' _p-1(k-L-1) を演算して求めることを特徴とする請求項１記載の適応
的制御方法。Wherein said input signal x (k), when an input signal vector obtained by applying the transfer characteristic x 'and _{p-1 (k-1)} , the cross-correlation vector p in the first step
(K) is determined by calculating the p (k) = p (k -1) + x (k) x 'p-1 (k-1) -x (kL) x' p-1 (kL-1) The adaptive control method according to claim 1, wherein: 【請求項３】 上記近似的フィルタ係数ベクトルは、高
速射影法に基づいて計算されることを特徴とする請求項
１記載の適応的制御方法。3. The approximate filter coefficient vector is calculated based on a fast projection method.
1. The adaptive control method described in 1 .