JPH10247838A - Adaptive control method - Google Patents

Abstract

PROBLEM TO BE SOLVED: To reduce an arithmetic operation amount without decreasing the convergent speed. SOLUTION: An autocorrelation vector RP+ N2-2 (k1 ), a smoothing coefficient vector SP+n2 (k1 ), an approximation filter coefficient Z(k), and an if error signal vector EP(k1 ) of an input x(k1 ) are initialized (S0 ), and if k1 mod N1 =1 is true (C1 ), then inner product is calculated as y'(k1 +i)=XL(k2 +i).Z(k2 ) (S2 ), and then RP+ N2-2 (k1 ), y(k1 ), e(k1 ), EP(k1 ), and a pre-filter coefficient vector G(k1 ), SP+n2 (k1 ) are calculated (S3 ), k1 =k1 +1 and n2 =n2 +1 are calculated respectively (S4 ). If k1 mod N2 =1 is true (L2 ), Z(k2 +N2 ) is calculated, k2 =k2 +N2 and n2 =0 are calculated and the step returns to the step C1 (S6 ). The y'(k1 +i) is calculated for each N1 time, Z(k2 ) is calculated for each time N2 (N2 <N1 ) and the others are calculated each time.

Description

DETAILED DESCRIPTION OF THE INVENTION

[0001]

BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention relates to an adaptive control method for minimizing the power of an estimation error of an adaptive filter in an acoustic echo canceller, active noise control or the like.

[0002]

2. Description of the Related Art First, in this section, an adaptive filter and an adaptive algorithm used for adaptive control will be described. Next, the projection method in the adaptive algorithm for correcting the adaptive filter will be described. Furthermore, as a method of reducing the amount of calculation of the projection method, the high-speed projection method with a modified operation procedure and the adaptive filter correction are not performed every time. The projection method based on the block processing will be described. Finally, problems of the projection method based on block processing will be described. 3.2.1 Adaptive Filter and Adaptive Algorithm FIG. 6A shows a block diagram of an adaptive control system using the adaptive filter. In the figure, an input signal is represented by x (k), an output signal is represented by y (k), a desired signal is represented by d (k), and an error signal is represented by e (k). However, in this specification, the expression of time represents time by an integer k as discrete time. The adaptive filter is an FIR type filter, and has a time-varying coefficient ｛h
₁ (k), h ₂ (k),..., H _L (k)}. Here, L is the number of filter coefficients. The filter coefficients are expressed as vectors as H (k) = [h ₁ (k), h ₂ (k),..., H _L (k)] ^T (1). Here, [] ^T represents transposition of a vector or a matrix. The adaptive filter adjusts the filter coefficient H (k) based on the error signal e (k) and the input signal x (k) to minimize the power of the error signal e (k). The procedure of coefficient correction at this time is known as an adaptive algorithm.

Next, an adaptive algorithm will be described. The adaptive algorithm means that the input signal x (k) is supplied to the filtering unit 11 as shown in FIG.
Output signal y (k) of filtering section 11 and desired signal d
An error signal e (k) with respect to (k) is obtained by a difference circuit 13, the error signal e (k) and the input signal x (k) are input to a filter coefficient update unit 12, and a filter coefficient H (K + 1) is updated. An adaptive algorithm is an operation procedure inside such an adaptive filter,
At time k, the filter coefficient H (k) is corrected by the filter coefficient update unit 12 to obtain a coefficient H (k + 1) that minimizes the power of the error signal e (k). The correction is expressed by the following equation.

H (k + 1) = H (k) + μΔH (k) (2.1) Here, ΔH (k) is a correction vector, and μ is an amount for controlling the magnitude of correction and is called a step size. Currently, the most widely used adaptive algorithm is LMS
Methods and learning identification methods are known. In the LMS method, the correction vector is an error signal e (k) and an input vector x
Δ (k) = e (k) X _L (k) (2.2) as a product of (k), and the correction formula is H (k + 1) = H (k) + μe (k) X _L (k (2.3). However, X _L (k) = [x (k), x (k−1),..., X (k−L + 1)] ^T (3.1), and this is called an input (signal) vector or simply an input. Also, the filter output y (k + 1) at the next time k + 1
Are synthesized by the inner product operation by the following equation in the filtering unit 11 of FIG. 6A.

Y (k + 1) = X _L (k + 1) ^T H (k + 1) (3.2) The filter coefficient is corrected based on equation (2.3), and the next time is corrected based on equation (3.2). The product-sum operation amount required for the operation of the adaptive filter for synthesizing outputs is 2L (L is the number of filter coefficients). This computation amount of 2L is the minimum computation amount among the currently known adaptive algorithms for correcting the hourly filter coefficient. 3.2.2 Projection method 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, it has a drawback that the convergence speed is slow (it takes time until an appropriate filter coefficient is obtained) for a colored signal, periodic noise, and the like.

A projection method was proposed in 1984 as a method for overcoming the disadvantages of the LMS method and the learning identification method. However, the projection method has a disadvantage that the amount of calculation is large. However, in recent years, a high-speed projection method has been proposed in which the disadvantage of the projection method is improved and the method is executed with a small amount of calculation. Here, the outline will be described below. In projection method, in order to reduce the error, the correction of the filter coefficients, the next number p of formula _{d (k) = X L (} k) T H (k + 1) d (k + 1) = X L (k + 1) T H (k + 1) (4)... D (k−p + 1) = X _L (k−p + 1) ^T H (k + 1) Here, p (p <L) is called a 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 _L (k−i +
1) convolution with ^T result _{X L (k-i + 1} ) T H (k +
1) matches the past desired signal d (ki + 1),
Is established for the past p inputs,
H (k + 1) is determined. As described above, if the filter coefficient is determined so that the same value as the desired signal is output with respect to the past input, the input X _{L at the} future time can be obtained.
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, the future error signal e (k + 1) = d (k + 1) -y
(K + 1) is also considered to be a small value. Now,
Solving for formula (2.1) are substituted in equation (4) ΔH (k), ΔH (k) = X L, p (k) [X L, p (k) T X L, p (k )] ^-1 _Ep (k) (5). Where _{XL, p} (k) = [ _XL (k), _XL (k-1), ..., _XL (k-p + 1)] (6) _Ep (k) = [d (k ), d (k-1) , ..., a d (k-p + 1) ] T -X L, p (k) T H (k) (7.1). The subscripts L and p of X indicate that X is L rows and p columns. This E _p (k) is expressed by the following equation:
It can also be expressed recursively.

[0007] Here, _Ep-1 (k-1) = [e (k-1), (1-.mu.) E (k-2), ..., (1-.mu.) ^P- 2e (k-p + 1)] ^T (7.3).

From equations (5) and (2.1), the input X
_The filter coefficient is corrected using _L (k) and the error signal e (k) to obtain a coefficient H (k + 1).
The procedure for synthesizing the output y (k + 1) of +1 can be expressed as the following equation. _{G (k) = [X L} , p (k) T X L, p (k)] -1 E p (k) (8) H (k + 1) = H (k) + μX Lp (k) G (k) (9) y (k + 1 ) = X L (k + 1) T H (k + 1) (10) However, G (k) represents the second half of the equation (5), called pre-filter coefficients. Equations (8), (9), and (10) are equations representing the basic procedure of the conventional projection method, and FIG. 6B shows a functional block diagram thereof. In FIG. 6B, the pre-filter coefficient calculation unit 21 performs an operation of Expression (8) from the input signal x (k) and the error signal e (k), and obtains the obtained pre-filter coefficient G (k) and input signal x (k). ), The filter coefficient updating unit 22 calculates the equation (9) to correct the filter coefficient, and the filtering unit 11 executes the calculation of the equation (10) using the corrected filter coefficient and the input signal.

Here, the amount of calculation of this system will be considered. In most DSPs, addition and multiply-accumulate operations are executed in one clock, respectively, and addition and multiply-accumulate operations occupy the largest part of the calculation. I do. First, the calculation of the equation (8) performed by the pre-filter coefficient calculation unit 21 is performed by a matrix X of p rows and p columns.
_{L, P} (k) ^T XL Since the calculation of the inverse matrix of _{L, P} (k) is included, the product-sum operation amount (hereinafter referred to as O
(Denoted as (p ³ )). This means that as p increases, the amount of calculation becomes very large. Next, the calculation of the expression (9) performed by the filter coefficient update unit 22 is such that X _{L, p} (k) is L
Since the matrix of rows and p columns and G (k) is a p-order vector, pL times of product-sum operations are required. Furthermore, since the operation of Expression (10) performed by the filter execution unit 11 is an inner product of an L-th order vector, L product-sum operations are required. Therefore,
The sum of (p + 1) L + O (p ³ ) is required. For example, if p = 20 and L = 1000, 2
More than 1,000 product-sum operations are required.

Next, a description will be given of a high-speed projection method with a devised operation and a projection method based on block processing in which adaptive filter updating is not performed every time but at certain block intervals, as a method of reducing the amount of calculation of the projection method. 3.2.3 High-speed projection method On the other hand, in the recently proposed high-speed projection method, the amount of product-sum operation can be reduced to 2L + 20p by devising the operation. Hereinafter, the high-speed projection method will be briefly described.

First, in the high-speed projection method, the O (p ³ ) of the O (p ³ ) is calculated by calculating the pre-filter coefficient G (k) using the linear prediction method without directly performing the inverse matrix operation of the equation (8). The product-sum operation amount is reduced to a product-sum operation amount of 16p. Since the reduction in the amount of calculation is not directly related to the present invention, the details are omitted. Next, in the high-speed projection method, the smoothing coefficient S
_{By using i} (k) and the approximate filter coefficient Z (k), the product-sum operation of (p + 1) L performed by the filter coefficient update unit 22 and the filtering unit 11 in FIG.
L + 4p product-sum operation amount is reduced. This method will be described with reference to FIG. 7A.

Referring to FIG. 7A, a pre-filter coefficient smoothing unit 31 obtains a smoothed coefficient vector by smoothing the pre-filter coefficient, and an approximate filter coefficient updating unit 32 uses an element of the smoothed coefficient vector as a weighting factor for the input signal vector. The approximate filter coefficient vector is corrected in the direction, the inner product of the approximate filter coefficient vector and the input signal vector is performed by the filtering unit 36, and the autocorrelation vector of the input signal vector is obtained by the correlation calculation unit 33, The inner product calculation with the vector is performed in the inner product calculator 35, and this output and the output of the filtering unit 11 are output from the output signal synthesizer 3
The sum is added at 4 to obtain an output signal. The following operation is performed in each of these components. That is, the pre-filter coefficient smoothing unit 31 performs the operation of Expression (11), the approximate filter coefficient update unit 32 performs the operation of Expression (12), the correlation calculation unit 33 performs the operation of Expression (13), and the inner product operation calculation unit 35 calculates Equation (14)
The inner product operation of the second term on the right side of
The addition of the first and second terms on the right side of 4) is performed. FIG. 8 shows a calculation procedure of the high-speed projection method corresponding to FIG. 7A.

[0013] Z (k + 1) = Z (k) + X L (k-p + 1) s p (k) (12) R p-1 (k + 1) = R p-1 (k) + x (k + 1) X p-1 (k) −x (k−L + 1) X _p−1 (k−L) (13) y (k + 1) = X _L (k + 1) ^T Z (k + 1) + R _p−1 (k + 1) ^T S _p−1 (k) ( 14) However, S _P (k) and S _p-1 (k-1) are each a p-order and p-1 order smoothing coefficient _{vector, S p (k) = [} s 1 (k) _{, s 2 (k), ...} , s p (k)] T (15) S p-1 (k) = [s 1 (k), s 2 (k), ..., s p-1 (k)] ^It is represented as ^T (16). As can be seen from equation (11), S _p (k)
Is recursively calculated using the pre-filter coefficient vector G (k). Further, Z (k) is an approximate filter coefficients vector, recursively calculated using the _{X L (k-p + 1} ) and s _p (k). R _p-1 (k + 1) is a correlation vector represented by the following equation.

R _p−1 (k + 1) = [ _XL (k + 1) ^T _XL (k), _XL (k + 1) ^T _XL (k−1),..., _XL (k _{^{+1) T X L (kp-}} 2)] T (17) in this case, a mathematical proof is omitted, the configuration of FIG. 7A,
Alternatively, the output signal y is calculated using the input signal x (k) and the error signal e (k) by the calculation of the equations (11) to (14).
(K + 1) can be synthesized. Expression (13)
Although the time parameter of (14) is k + 1, it is k in FIG. 7A. This means that at time k, equation (1)
3), (14), (11), and (12) are calculated in this order.

Now, the following two points can be pointed out. First, regarding the calculation amount of the above method, p is used in equation (11), L is used in equation (12), 2p is used in equation (13), and equation (1) is used.
In 4), it can be executed with L + p, that is, a total of 2L + 4p product-sum operations. This is compared with the case where the operations of Expressions (9) and (10) performed by the filter coefficient updating unit 22 and the filtering unit 11 of the conventional projection method shown in FIG. 6B require the product-sum operation of (p + 1) L. Has been greatly reduced. For example, p
= 20 and L = 1000, the projection method is 21000, whereas the high-speed projection method is 2080, which is about 1/10
Is the calculation amount.

The second point is that, in the high-speed projection method shown in FIG. 7A, the output signal y (k) is synthesized without explicitly obtaining the filter coefficient H (k + 1), as shown in FIG. 6B. This is a major difference from the conventional projection method. More specifically,
Instead of the filter coefficient H (k + 1), an approximate filter coefficient Z (k + 1) is used, and H (k + 1) and Z (k + 1) are used.
In order to correct the difference from (k + 1), the correlation vector R
_p-1 (k + 1) is used. This has the effect of reducing the amount of computation. 3.2.4 Conventional block projection method In the high-speed projection method described so far, the filter coefficient H
The calculation amount can be reduced by updating the approximate filter coefficient Z (k) every time instead of (k). Separately, the amount of calculation can be reduced by updating the filter coefficient H (k) every N (> 1) time (block processing).

The block processing method is classified into two types in terms of whether or not the output y (k) is equal to the projection method for correcting the hourly filter coefficient. When the output y (k) is different from the projection method for correcting the filter coefficient every time, the amount of calculation is reduced, but the convergence speed is reduced. On the other hand, in the case of the block processing having the same output y (k) as the projection method for correcting the hourly filter coefficient, the degree of reduction in the amount of operation is different from the output y (k) different from the projection method for correcting the hourly filter coefficient. However, there is an advantage that the convergence speed is not deteriorated, although it is smaller than that of the block processing method having.

<A> A block projection method in which the output y is different from the projection method in which the filter coefficient is corrected every time. In the conventional block projection method in which block processing is applied to the projection method, Ep (k) and G (k) are calculated. And H (k +
N) is updated in the same manner as in equations (7.1), (8), and (9).

[0019] _{E p (k + N-1} ) = [d (k + N-1), d (k + N-2), ..., d (k + N-p)] T -X L, p (k + N-1) T H (k ) (18) G (k + N-1) = [X L, p (k + N-1) T X L, p (k + N-1)] -1 E p (k + N-1) (19) H (k + N) = H (K) + μX _{L, p} (k + N−1) G (k + N−1) (20) The calculation of the output y is the output y (k +
1) Instead of the calculation, N (y + k),
, Y (k + 2N-1) are calculated at once.

[0020] _{Y N (k + 2N-1} ) -X L, N (k + 2N-1) T H (k + N) (21) Some or all of the elements of the second term of Formula (18) Formula (2
In 1), it is calculated by an equation in which k + 2N is replaced by k + N. Equation (1
8), (19), (20), and (21) are conventional output y
(K) is a calculation procedure of the block projection method different from the projection method of correcting the filter coefficient for each time, and its functional block diagram is shown in FIG. 7B. In FIG. 7B, a pre-filter coefficient calculation unit 21, a filter coefficient update unit 22, a block filtering unit 41, and a block error calculation unit 42 execute the operations of Expressions (19), (20), (21), and (18), respectively. I do. FIG. 9 shows a calculation procedure of the conventional block projection method corresponding to FIG. 7B.

There are two aspects to the reduction in the amount of calculation by block processing. One is a filter coefficient H (k) of the operation amount.
Is corrected every time in the projection method, but is reduced to once at N (> 1) time in the block projection method. The other is X _{L, P} (k + N) on the right side of equation (20).
-1) G (k + N-1) and XL _{, N} (k
^{+ 2N-1) T H (} k + N) is a reduction of the calculation amount by computable at high speed. XL _{, N} (k + 2N) of the latter
-1) calculation of ^T H (k + N) will be described.

[0022] X _L, to calculate N a ^{(k + 2N-1) T} H (k + N) efficiently, as in the following _{equation, X L, N (k +} 2
N-1) is decomposed into a partial matrix of N rows and N columns, and H (k + N) is decomposed into a partial vector of N elements. However, H is added with an element having a value of 0 so that L = NM (M is an integer). X _{L, N} (k + 2N−1) ^T = [X _{N, N} (k + 2N−1) ^T , X _{N, N} (k + N−1),..., X _{N, N} (k + 2N−1− (M−1) N ^{) T] T (22) H} (k + N) T = [H (k + N) 0 T, H (k + N) 1 T, ..., H (k + N) M-1 T] T is decomposed in this manner, X _{L, N} (K + 2N-1) ^T H
(K + N) can be calculated as the sum of the product between the submatrix of X and the subvector of H.

[0023] Each in _{X L, N (k + 2N} -1) T H (k + N) = Σ i = 0 M-1 X N, N (k + 2N-1-iN) T H (k + N) i (23) Σ Square matrix X _{N, N} (k + 2N−1−iN)
Is the product of the input signal vector X _2N (k + 2N-1-iN) and the filter coefficient H (k + N) _i.
N) can be considered as a convolution of _i . As a method of performing high-speed convolution, a method using FFT, Winog
rad based fast convolution method (for details, see Reference B
lahut, FAST ALGORITHMS FOR
DIGITAL SIGNAL PROCESSIN
G, Addison-Westley, 1987, Chapters 3, 9). These fast convolution methods are expressed by the following equation (2).
By applying the method to 3), the calculation amount of Y _N (k + N) can be reduced.

The principle of fast convolution of the method using FFT will be briefly described. As is well known, convolution in the time domain is a product in the frequency domain. F [x * h] = F [x] ★ F [h] (24) where F [] is Fourier transform, * is convolution,
★ represents the product for each element. Using this relationship, the convolution of x and h can be calculated by the following equation.

X * h = F ⁻¹ [F [x * h]] = F ⁻¹ [F [x] ★ F [h]] (25) where F ⁻¹ [] represents an inverse Fourier transform. . By using the fast Fourier transform (FFT) for the Fourier transform, it is possible to perform high-speed convolution. Equation (20),
When applied to (21), the following equation is obtained.

H (k + N)_m = H (k)_m+ ΜF^-1[F [X_2p(K + N−1−mN)] F [G (k + N−1)]] m = 0, 1,..., L / p−1 (26) Y_N(K + 2N-1) = Σ_{m = 0} ^M-1F^-1[F [X_2N(K + 2N-1-mN)] F [H (k + N)_m]] = F^-1[Σ_{m = 0} ^M-1F [X_2N(K + 2N-1-mN)] F [H (k + N)_m]] (27) Taking as an example the equation (27) corresponding to the equation (21),
The 2N FFT is X _2N(K + 2N-1) once (X_2N
(K + 2N-1-mN), the FFT of m> 1 is the time k + 2
N−1−mN, so do not count), H (k
+ N)_mM times, 2N inverse FFT once, and then
Thus, 2NM complex product sum operations are required. This operation
When the quantity is replaced by the sum of the addition of real numbers and the product-sum operation,
Data length is 2N FFT, inverse FFT is about 8Nlog_Two(2
N) The product-sum operation of one complex number is 4. These
Thus, the amount of calculation required for equation (27) is about (M + 2) 8Nlo
g_TwoN + 8 NM, almost MNlog_TwoProportional to N
You. On the other hand, equation (27) is directly calculated as in equation (21).
The amount of calculation when NL (= N^TwoM), so FF
It can be seen that the use of T reduces the amount of calculation.

<B> Block processing projection method having the same output as the projection method for correcting the filter coefficient every time In the block projection method described above in 3.2.4 <A>, only one block length N is used. Since the adaptive filter is not modified, there is a disadvantage that the convergence of the adaptive filter becomes slow as the block length N is increased.

For this disadvantage, in the high-speed projection method,
Correlation vector R _p-1 (k) and smoothing coefficient vector S
_The solution is to extend the calculation range of _p (k) and to reduce the amount of computation of the projection method by block processing without lowering the convergence speed of the adaptive filter. This is a block process having the same output as the projection method for correcting the time filter coefficient.

Hereinafter, this method will be described. In the high-speed projection method for modifying the above-mentioned adaptive filter at each time, Z (k + i + 1) and (i> = 0) are obtained from Expression (12).
(K). Z (k + 1 + i) = Z (k + i) + s p (k + i) X L (k + 1 + i-p) = Z (k + i-1) + s p (k + i-1) X L (k + i-p) + s p (k + i) X L ( k + 1 + i-p) ... = Z (k) + Σ j = 0 i s p (k + i-j) X L (k + 1 + i-p-j) = Z (k) + X L, i + 1 (k + 1 + i-p) S p + _i (k + i) _p to _{p + i} (28), and especially when i = N−1, the update formula is an update from the approximate filter coefficient Z (k) to Z (k + N).

Z (k + N) = Z (k) + X_{L, N}(K + N-p) S_{p + N-1}(K + N-1)_p~_{p + N-1} (29) S in equation (28)_{p + i}(K + i) is the value of equation (15)
S of degree p_p(K) with order p + i
This is called an extended smoothing coefficient vector, and S _{p + i}(K
+ I)_p~_{p + i}Is the extended smoothing vector S_{p + i}(K + i)
From the p-th element to the p + i-th element
Means S_{p + i}(K + i)_p~_{p + i}Is recursive by
Is calculated.

[0031] Further, the extended smoothing coefficient vector S _{p + i} (k + i) is represented by S _p−1 (k + i) and S _{p + i} (k + i) _p in Expression (16).
It is calculated by the following equation by combining _{p + i} .

[0032] By substituting the equation (30) and the regression equation (11) of the smoothing coefficient vector into the equation (31) in order, the extended smoothed vector S
The regression equation of _{p + i} (k + i) is obtained.

[0033] By the way, in the equation (14), by setting k to k + i−1, the following equation y (k + 1 + i) is obtained.

Y (k + i) = X_L(K + i)^TZ (k + i) + R_p-1(K + i)^TS_p-1(K + i-1) (33) From the right side of the equation (28), Z (k + i) is converted to Z (k)
The following expression is obtained. Z (k + i) = Z (k) + X_{L, i}(K + ip) S_{p + i-1}(K + i-1)_p~_{p + i-1} (34) By substituting equation (34) into equation (33), the following equation is obtained. y (k + i) = X_L(K + i)^T[Z (k) + X_{L, i}(K + ip) S_{p + i-1}(K + i-1)_p~_{p + i-1}] + R_p-1(K + i)^TS_p-1(K + i-1) = X_L(K + i)^TZ (k) + X_L(K + i)^TX_{L, i}(K + ip) S_{p + i-1}(K + i-1)_p~_{p + i-1} + R_p-1(K + i)^TS_p-1(K + i-1) = X_L(K + i)^TZ (k) + R_{p + i-1}(K + i)^T _p~_{p + i-1} S_{p + i-1}(K + i-1)_p~_{p + i-1} + R_p-1(K + i)^TS_p-1(K + i-1) (35) where R_{p + i-1}(K + i) is called an extended correlation vector.
And R_{p + i-1}(K + i) _p~_{p + i-1}Is the extended correlation vector
A vector consisting of the p-th element to the (p + i-1) -th element
Means torr. R_{p + i-1}(K + i)_p~_{p + i-1}Is R
_{N + p-2}(K + i)_p~_{N + p-2}Some or all elements of
And R_{N + p-2}(K + i)_p~_{N + p-2}Is recursive by
Is calculated.

R _{N + p−2} (k + i) ^T _p to _{N + p−2} = _XL (k + i) ^T _{XL, N−1} (k + ip) = R _{N + p−2} (k + i−1) ^T _{p NN + p−2} + x (k + i) X _N−1 (k + ip) −x (k + i−L) X _N−1 (k + i−p−L) 17) R _p-1 (k +
It is calculated by combining i) and R _{p + i-1} (k + i) _p to _{p + i-1} .

[0036] The time update of the extended correlation vector is performed by using the equations (13) and (3).
6) is substituted into equation (37) and multiplied as follows.

RN _{+ p-2} (k + i) = _{RN + p-2} (k + i-1) + x (k + i) _{XN + p-2} (k + i-1) -x (k + i-L) _{XN + p −2} (k + i−L−1) (38) When the extended smoothing coefficient vector and the extended correlation vector are used, Expression (35) is multiplied as follows. y (k + i) = y '(k + i) + R _{i + p-1} (k + i) ^T S _{i + p-1} (k + i-1) (39) y' (k + i) = _XL (k + i) ^T Z (k) (40) Since y (k), _Ep (k), G (k) and _Sp (k) can be calculated by the procedure of the high-speed projection method, the equation (39) is set as i = 0.
Can be used to calculate y (k). Using equation (7.2), e (k) and E _p (k) can be calculated from y (k).
_{From p} (k), G (k) can be calculated by equation (8). G
From (k), S _{p + 1} (k) can be calculated using equations (11), (30), and (31).
(K + 1) can be calculated. Thus, y (k +
i), (i = 0, 1,..., N−1) can be calculated, and at time k + N−1, Z (k +
N) can be calculated.

In summary, as shown in FIG. 11, first, as initialization, R _{p + N−1} (0) = _XL (0) ^T [ _XL (−
1), X _L (−2),..., X _L (−p−N + 1)] ^T, and E _p−1 (0) = S _p−1 (0) = X _L (0) = [0,
0,..., 0] ^T (S ₁₀₀ ), and then k = 1 (S
_101), then the inner product computation y the approximate vector and the input vector '(k + i) = X L (k + i) T Z (k), i = 0,
,..., N−1 are performed ( _S102 ). Then, for i = 0, 1,..., N−1, equations (38), (39), (7.2), (8),
The calculation of (32) is performed for each time, and when the calculation of i = N−1 is completed (S ₁₀₃ ), the expression (29) is calculated and Z
(K) is updated to Z (k + N) ( _S104 ). next,
The process returns to step _S102 with k set to k + N ( _S105 ).

A functional block diagram corresponding to this procedure is shown in FIG.
0, and the same reference numerals as those in FIGS. 7A and 7B denote the same parts. Expression (13) is calculated from the input signal vector by the correlation calculation unit 33, and Expression (36) is calculated from the input signal vector by the correlation calculation expansion unit 52. The correlation vector from the correlation calculation unit 33 and the correlation calculation expansion The extended vector of the correlation vector from the section 52 is connected to the vector of the expanded section by the equation (37) in the correlation vector connecting section 39 to obtain an expanded correlation vector. In short, the correlation calculation unit 33, the correlation calculation expansion unit 52, and the correlation vector connection unit 53 calculate Expression (38), and constitute an extended correlation vector calculation unit 56.

The pre-filter coefficient calculator 21 calculates equation (8) from the input signal vector and the error signal vector to calculate a pre-filter coefficient. The pre-filter coefficient is calculated by the pre-filter coefficient smoother 31 as equation (30). ), And the vector of the extended portion and the smoothing coefficient vector from the smoothing unit 31 are connected by the operation of the equation (31) by the smoothing coefficient vector connecting unit 55 to obtain an extended smoothing coefficient vector. In short, the pre-filter coefficient smoothing unit 31, the smoothing coefficient expanding unit 54, and the smoothing coefficient vector connecting unit 55 are expressed by the formula (3)
An extended smoothing coefficient vector calculator 57 for calculating 2) is configured.

The inner product operation of the extended correlation vector from the extended correlation vector operation unit 56 and the extended smooth coefficient vector from the extended smooth coefficient vector calculation unit 57, that is, equation (39)
Is performed on the second term on the right side of. The extended partial vector from the smoothing coefficient extending unit 54 is input to the approximate filter coefficient block updating unit 51, and the equation (30) is calculated. An inner product operation between the approximate filter coefficient of the operation result and the input vector is performed by the block filtering unit 41, and the operation result and the operation result from the inner product operation calculation unit 35 are added by the output synthesis unit 34, and the output signal and Is done. Block filtering section 41, approximate filter coefficient block updating section 5
1 is calculated every N times, the other parts are performed every time, and the extended smoothing coefficient vector part 57 is calculated by the block filtering part 41, and the number of characters is changed from p to primary at each time until the next calculation. This is repeated, and this is repeated. In the extended correlation vector calculation unit 56, a correlation vector of the same order as the maximum order p + N-1 of the extended smoothing coefficient vector is always calculated.
The number of characters of the same degree as the number of characters of the extended smoothing coefficient vector at that time from the extended smoothing coefficient vector calculator 57 is input.

As described above, the filter coefficient for each time is corrected.
The complexity of the block projection method with the same output as the projection method
Ask. First, S in FIG.
₁₀₂As described in 3.2.4 <A>,
The high-speed convolution used can be used, and the amount of computation is about 8 (M
+2) Nlog_Two2N + 8NM, per hour
Is 8 (M + 2) log_Two2N + 8M. S₁₀₃so
Is (2 + N−2) times in equation (38) per time,
Equation (39) is given by p + N / 2 times, Equation (7.2) is given p times,
(8) requires 15p computations, and equation (32) requires p computations.
It is important, and the total is 2.5N + 20p. S₁₀₄
Is high using FFT as described in 3.2.4 <A>.
It can be calculated quickly. That is, equation (29) is the same as equation (20).
, So that Z (k + N)_m= Z (k)_m+ ΜF^-1[F [X_2N(K + N-p-mN)] ★ F [S (k + N-1)_{p, N + p-1}]] M = 0, 1,..., L / N-1 (41) X having a data length of 2N_2NTo FFT (k + N-p)
Once, a vector obtained by adding N zeros to S (k + N-1)
One FFT of 2N length data and one inverse FFT
L / N (= M) times are required. In addition, the product-sum operation of complex numbers
Is required 2MN times. These FFT, inverse FFT,
The operation amount required for the product-sum operation of prime numbers is about 8 (M + 2) Nlo
g_Two2N + 8NM, and 8 (M +
2) log _Two2N + 8M. From the above, S₁₀₂, S
₁₀₃And S₁₀₄The amount of computation per time required for
Thus, T = 16 (M + 2) log_Two(2N) + 2.5N + 16M + 20p (42)

[0043]

The output y is the same as that of the conventional projection method for correcting the hourly filter coefficient (hereinafter, the conventional method).
In the block projection method having (k), there is a block length that minimizes the calculation amount T, and a longer block length has a problem that the effect of reducing the calculation amount is reduced or eliminated. The reason for this is that the amount of computation 2N required for updating the extended correlation vector R _{N + p−2} calculated by the equation (38), the extended correlation vector calculated by the second term on the right side of the equation (39), and the extended smoothing coefficient vector The amount of computation required for the product of
Since 5N is proportional to the block length N, the block length N
Is to increase with lengthening. Using the equation, it can be seen that the differential coefficient is positive with respect to the block length N of T in equation (42).

∂T / ∂N = 16 ((−L / N^Two) Log_Two(2N) + ((L / N) +2) / (Nln2)) + 2.5-16 (L / N^Two) = N^-2(16 (−Llog_Two(2N) + (L + 2N) / ln2) + 2.5N^Two-16L) (43) N^-2Is always positive and 16 (-Llog_Two(2N)-+
(L + 2N) / ln2) + 2.5N^Two-16L is L >>
N takes a negative value, but when N increases, N, N ^TwoSection
Log compared to_Two(2N) term can be ignored,
Since it can be regarded as a quadratic expression of N, the value is positive. That is,
The complexity decreases at N << L, but increases at a certain N or more.
I do.

When the filter length L is given, N at which the amount of operation T is the minimum is the block length N at which the right side of equation (43) becomes zero. An equation in which the right side of equation (43) is set to 0 can be obtained by recursive calculation regarding N as shown below. First, an equation in which the right side of equation (43) is set to 0 is regarded as a quadratic equation of N. N = [16 / ln2 + 4√ (1 / ln2) 2 + 2.5L (log 2 (2N) -1 / ln2 + 1)] /2.5 = 9 + 1.6 √2.5Llog 2 (2N) = 9 + 2.5 √ Llog ₂ (2N) (44) The recursive calculation is performed by substituting the block length N on the left side of equation (44) for the right side. If the block length N obtained by the i-th recursive calculation is expressed as N (i), equation (44) becomes the following equation.

N (i + 1) = 9 + 2.5√Llog ₂ (2N (i)) (45)

[0047]

In order to solve the above-mentioned problems, the present invention solves the above problems by calculating the output calculated by the equation (40) and updating the approximate filter coefficient Z calculated by the equation (29) in the conventional method. Is to solve the problem by using different block lengths, thereby making it possible to reduce the amount of computation and the amount of delay that cannot be obtained by the conventional method.

According to the present invention, the input signal is subjected to the variable filter processing, the error signal between the output signal and the desired signal is obtained, and the characteristic of the variable filter is adaptively adjusted so as to minimize the power of the error signal. A first step of performing an operation corresponding to a product of an inverse matrix of the autocorrelation matrix of the input signal vector and the error signal vector to obtain a prefilter coefficient vector; A second step of smoothing to obtain a smoothed coefficient vector, a third step of correcting an approximate filter coefficient vector in the direction of the input signal vector using the elements of the smoothed coefficient vector as weighting coefficients, A fourth step of performing an inner product operation with an approximate filter coefficient vector, and calculating an autocorrelation vector of the input signal vector A fifth step of calculating the inner product of the smoothing coefficient vector and the autocorrelation vector, a result of the inner product operation obtained in the fourth step, and an inner product operation obtained in the sixth step And a seventh step of adding the result to obtain the output signal. The correction of the third step and the calculation of the fourth step are performed at N ₁ and N ₂ times (N ₁ and N _{2 respectively).} 2
For each of the above integers, N ₁ ≠ N ₂ ), the third and fourth
Each step other than the step performs an operation at each time,
Every time the calculation of the third step is performed, the order of the smoothing coefficient vector obtained in the second step is increased by one at each time until the calculation of the third step next, and the next third The operation at the operation time of the step returns to the original expression. In the fifth step, the autocorrelation vector having the same order as the maximum order of the smoothing coefficient vector is calculated. In the sixth step, the order of the smoothing coefficient vector is calculated. Is used from the autocorrelation vector obtained in the fifth step by an order equal to

[0049]

Embodiments of the present invention will be described below. Generally, in an adaptive algorithm that performs block processing, there are two calculation steps that can reduce the amount of calculation by block processing. In block projection method having the same output as the projection algorithm to correct each time the adaptive filter mentioned above, wherein the output calculation and 11 in S ₁₀₄ as calculated by the equation (40) in S ₁₀₂ in FIG. 11 ( This is the update of the approximate filter coefficient Z calculated in 29). By the way, the reason why the operation reduction effect becomes smaller with respect to the long block length N, which has been a problem in the above-described conventional method, is that the extended correlation vector R _{N + p−2} calculated by Expression (38) is updated over time. Since the required operation amount 2N and the operation amount 0.5N required for multiplying the extended correlation vector calculated by the second term on the right side of the equation (39) and the extended smoothing coefficient vector are proportional to the block length N, the block length N is increased. It was to increase with time. The calculation that increases in proportion to the block length occurs because the update of the approximate filter is suspended only for the block length. Therefore, the above problem can be solved by using different block lengths in the above two block processes.

The block length in the two block processing is
Being independent, the process of the present invention is different from the conventional process.
May differ in order or timing of processing
You. The block length in the output calculation is N₁, Approximate filter
Block length in coefficient update is N _TwoEquation (40),
(39), (7.2), (8), (32) and (2)
Corresponding to the processing by the conventional method given in 9), the method of the present invention
Is as follows.

Y ′ (k₁+ I) = X_L(K_Two+ I)^TZ (k_Two), I = 0, ..., N₁-1 (46) R_{N2 + p-2}(K₁) = R_{N2 + p-2}(K-1) + x (k₁) X_{N2 + p-2}(K₁-1) -x (k₁-L) X_{N2 + p-2}(K₁−L-1) (47) y (k₁) = Y '(k₁) + R_{n2 + p-1}(K₁)^TS_{n2 + p-1}(K₁-1) (48) e (k₁) = D (k₁) -Y (k₁) (49) G (k₁) = [X_Lp(K₁)^TX_Lp(K₁)]^-1E_p(K₁) (51) Z (k_Two+ N_Two) = Z (k_Two) + X_{L, N2}(K_Two+ N_Two-P) S_{N2 + p-1}(K + N_Two-1)_{p-N2 + p-1} (53) where k_TwoIs related to the modification of the approximate filter coefficient Z
Time index, k₁Is the time indicator for the other variables.
You. Also, n_Two= K₁mod N_TwoVariable with a value of -1
Is a number. FIG. 1 shows a calculation procedure according to the present invention.
First, as initialization, R_{p + N2-1}(0) = X_L(0)^T[X_L(-1), X_L(-2),
…, X_L(-p-N_Two+1)]^T S_p-1(0) = [0,0, ..., 0]^T Z_L-1(1) = [0,0, ..., 0]^T E_p-1(0) = [y (0), y (-1), ..., y (-p
+2)]^T And (S₀), And then set k1 = k2 = 1 and n2 = 0 (S
₁), Then, if k₁mod N₁Judgment whether = 1
Set (C₁), Judgment C₁Is true, the input vector
Product operation y '(k₁+ I) = X_L(K_Two+ I)^TZ (k_Two), I = 0, ..., N₁-1 (S_Two)I do. Thereafter, equations (47), (48), and (4)
9), (50), (51), and (52) are sequentially calculated.
No (S_Three), Then k₁And n_TwoIs increased by 1 (S _Four),
Then, if k₁mod N_Two= 1
(C_Two), Judgment C _TwoIs true, approximate filter coefficients
Modification of Z (S_Five) Z (k_Two+ N_Two) = Z (k_Two) + X_{L, N2}(K_Two+ N_Two−
p) S_{N2 + p-1}(K + N_Two-1)_{p-N2 + p-1} And k_TwoTo N_TwoIncrease by n_TwoTo 0
(S₆). Next, the process returns to the determination C1.

A functional block diagram corresponding to this procedure is shown in FIG.
7A and 7B, and are denoted by the same reference numerals. The unit in FIG. 2 has the same function as the unit in FIG. 11 which is a functional block diagram of the conventional method, but the order of operation is different. Estimate the operation amount of the present invention. First, S ₂ of Figure 1 showing a procedure of the method, as described in the conventional method, FFT fast convolution are available with, the calculation amount is about _{T 0 = 8 (M 1 +2} ) N ₁ log ₂ 2N ₁ +
8N ₁ M ₁ , and per time T ₀ / N ₁ = 8
A _{(M 1 +2) log 2 2N} 1 + 8M 1. here,
L = N ₁ M ₁ . In S ₃ , the expression (4)
7) 2 (N ₂ + p−2) times, and equation (48) shows p + N ₂ /
The calculation needs to be performed twice, p times in equation (50), 15p times in equation (51), and p times in equation (52), and the total is 2.5N ₂ + 20p (T _c = N ₂ (2.5N ₂ +
20p)). Calculation amount required for the S _5, as described in the conventional method, if L = N ₂ M _2, about T _u =
8 (M ₂ +2) a _{N 2 log 2 2N 2 + 8N} 2 M 2, 1 in per _{time, T u / N 2 = 8} (M 2 +2) l
og ₂ 2N ₂ + 8M ₂ . As described above, when the amount of operation required to execute the present invention is described as T, T is the sum of the amount of operation per time required for S ₂ , S ₃ , and S ₅ , and T = T ₀ / N ₁ + T _C / N ₂ + T _u / N ₂ = 8 (M ₁ +2) log ₂ (2N ₁ ) +8 (M ₂ +2) log ₂ (2N ₂ ) + 2.5N ₂ + 8M ₁ + 8M ₂ + 20p (54) FIG. 3 shows an example of the relationship between the calculation amount T calculated by the equation (42) and the block lengths N ₁ and N ₂ .
The conditions are: filter length L = 4096, block length N ₁ ,
N ₂ is a 128,256,512,1024,2048. In the figure, Δ is the amount of calculation possible by the conventional method of N ₁ = N ₂ , the block length is 512 and the amount of calculation is the minimum,
The calculation amount increases whether the block length is shorter or longer. Meanwhile, □ is a computation amount by N ₁ ≠ N ₂ becomes the present invention method, with a longer block length N _1, the amount of computation is reduced, the block length N ₁ = 1024, 2048, conventional methods of equal block length It is possible to make it smaller than the minimum calculation amount in the case of. When N ₂ is changed, N ₂ = 25
6 is the minimum for any N ₁ . From FIG. 3, it can be seen that the operation amount of the present invention is smaller than that of the conventional method when the block length N ₁ is longer than the block length N which gives the minimum operation amount of the conventional method calculated by the equation (45). . L = 40
In this example of 96, the recursive calculation of equation (45) is given by N (0) =
Starting with an initial value of 512, we get N = 487.

Since the amount of calculation is usually defined and a filter length as long as possible within the range is implemented, it may be more convenient to express the filter length as a function of the amount of calculation.
By substituting L = N ₁ M ₁ and L = N ₂ M ₂ into equation (42) and solving for the filter length L, the following equation is obtained. L = [T-16log ₂ (2N ₁ ) -16log ₂ (2N ₂ ) -2.5N ₂ -20p] / [8 (log ₂ (2N ₁ ) +1) / N ₁ +8 (log ₂ (2N ₂ ) +1) / N ₂ ] (55) FIG. 4 shows an example of the relationship between the filter length L and the block lengths N ₁ and N ₂ calculated by the equation (55). The condition is that the calculation amount is T = 8000 and the block lengths N ₁ and N ₂ are 128,2.
56, 512, 1024, 2048, and 4096.
In the figure, Δ is a filter length that can be obtained by the conventional method of N ₁ = N _{2. The} block length N ₁ is 1024, the filter length is the longest, and the calculation amount increases even if the block length is shorter or longer. I do. In particular, the block length N ₁ does not L> 0 and becomes solution at 4096 or higher. On the other hand, □ is a filter length according to the present invention, which satisfies N ₁ ≠ N _2. As the block length N ₁ increases, the filter length increases, and the block length N ₁
In the case of = 2048 or more, it is possible to implement a filter length longer than the longest filter length in the case of the conventional method having an equal block length.

Next, consider the time when the filter output y (k) is delayed by performing the block processing. If this delay time is represented by D, it can be seen from FIG. D = N ₁ + T ₀ / T, N ₁ > N ₂ = (2T ₀ + T _c + T _u ) / T, N ₁ < N ₂ (56) The output delay time D and the block calculated by the equation (56) in FIG. an example of a relationship between the length N _2. The conditions are filter length L = 4096, and block lengths N ₁ and N ₂ are 128 and 25.
6, 512, 1024, and 2048. In the figure, Δ is the output delay time according to the conventional method of N ₁ = N ₂ , □ is the output delay time according to the present invention of N ₁ ≠ N ₂ ,
The block length N ₁ is 128 from the bottom for each N ₂ ,
256, 512, 1024, and 2048. N ₂ is 2
Above 56, it can be seen that the method of the present invention can achieve a shorter output delay time than the conventional method.

The feature of the present invention described above is that N ₁ Ｎ
By setting N ₂ , it is possible to execute with a smaller amount of computation than the block projection method having the same output as the conventional projection method for modifying the time-of-day filter. In other words, a longer filter length can be used with the same computation amount. , The power of the error signal is further reduced. Also, the output delay time can be made shorter than that of the conventional method. In this case, since N ₂ > N ₁ and the update cycle of the filter coefficient becomes long, in the case of an echo canceller, abnormalities such as double talk may occur. Upon detection, it is also possible to perform a process of canceling the update at that time.

Embodiment 1 The first example is an acoustic echo canceller. FIG. 12A shows a functional configuration of an acoustic echo canceller to which the method of the present invention is applied. As the input signal x (k), the audio signal is converted by the filtering unit 11, the filter coefficient updating unit 12, the speaker 6
1 and the output of the microphone 62 is
(K) is supplied to the difference circuit 13. A signal reaching the microphone 62 from the speaker 61 through an acoustic signal path, that is, a so-called echo transmission path 63 is an acoustic echo, and an output signal y of the adaptive filter (filtering unit) 11.
(K) corresponds to the estimated echo, and the estimation error e (k) corresponds to the residual echo.

The acoustic echo canceller is a device for preventing howling and echo in a voice communication system such as a video conference system. In the case of an acoustic echo canceller, the adaptive filter 11 transmits a signal from the speaker 61 to the microphone 62.
Is simulated. The output of the adaptive filter 11 is an estimated echo, and the echo is removed by subtracting it from the echo received by the microphone 62. When the transmission path 63 between the speaker 61 and the microphone 62 is simulated by the adaptive filter 11, the length of the filter required is several hundreds to several thousands, and a large amount of calculation is required.
In addition, since the transmission path 63 between the speaker 61 and the microphone 62 changes due to movement of a person, opening and closing of a door, and the like, a high-speed adaptive algorithm that can follow the change is required.

The second application example is noise control. FIG.
Shows the principle of noise control. The purpose of noise control is to control the noise source 8
The noise observed through the noise transmission path 82 from No. 1 is converted from the speaker 86 into a negative sound (when the sound pressure is expressed as y (k),
The sound pressure represented by −y (k) is called a negative sound with respect to y (k). In this example, the adaptive filter (filtering unit) 1
1 is a noise obtained by observing the noise of the noise source 81 with the monitor microphone 84, and the desired signal d
(K) is the noise at the observation point (microphone 83), the output y (k) of the adaptive filter 11 is the estimated noise, and the error signal e (k) is the estimated noise radiated from the speaker 86 at the observation point. And the sound pressure of the noise arriving from the noise source are synthesized and obtained as the output of the microphone 83. The adaptive filter 11 is connected to the observation microphone 83 from the noise source 81.
If the transmission path 82 to the transmission path 82 is well simulated, the sound pressure -y (k) radiated from the speaker 86 causes the noise d.
(K) is deleted.

The two application examples have been described above. In either case, in order to achieve the purpose of canceling the echo and the noise, it is sufficient to obtain the estimated echo and the estimated noise, and it is not necessary to obtain the transmission path of the echo and the noise. Therefore, the present invention method using an approximate filter coefficient can be applied to both application examples, and the amount of calculation can be reduced as compared with the conventional apparatus. Although the calculation amount is further increased, the output y
(K) can be obtained quickly. Also,
The present invention has a basic function of minimizing the power of the error signal e (k) in the configuration as shown in FIG. 6A. Therefore, the present invention can be applied to all applications in which the problem to be solved can be modeled as the minimization of the power of the error signal shown in FIG. 6A.

[0060]

As described above, according to the present invention, two different block lengths are introduced in a block projection method having the same output as a projection method for correcting a time filter, which is one of adaptive algorithms, and an operation is performed. By changing the procedure, the amount of calculation can be reduced without lowering the convergence speed. As a result, the present invention does not reduce the convergence speed with respect to the acoustic echo canceller device in which a method such as the conventional block projection method having a slow convergence is used.
The amount of calculation can be reduced. In other words, it is possible to use a higher-order projection method that requires a little computational effort but to use a higher-order projection method, and to increase the filter length with the same hardware scale, thereby enabling faster echo reduction and steady-state echo reduction. Is achieved. According to the present invention, the amount of calculation cannot be reduced, but the delay time for obtaining the output y (k) can be shortened. In this case, the update time of the filter becomes longer, and during this time, abnormalities such as double talk may occur. When the occurrence is detected, a process of stopping the next update of the filter if detected can be performed.

[Brief description of the drawings]

FIG. 1 is a flowchart showing an embodiment of the present invention.

FIG. 2 is a diagram showing a functional configuration of an apparatus to which an embodiment of the present invention is applied.

FIG. 3 shows a block length N ₁ when a filter length is fixed,
Diagram showing the relationship N ₂ and calculation amount T.

FIG. 4 shows block lengths N ₁ and N ₂ when the calculation amount is fixed.
FIG. 6 is a diagram showing a relationship between the filter length and the filter length L.

FIGS. 5A and 5B are diagrams showing the operation execution time and the output delay time, respectively, and C is a diagram showing the relationship between the block lengths N ₁ and N ₂ and the output delay time D;

6A is a diagram illustrating an adaptive control system using an adaptive filter, and FIG. 6B is a diagram illustrating a functional configuration of control using a conventional projection method.

7A is a diagram showing a functional configuration of control using a conventional high-speed projection method, and FIG. 7B is a diagram showing a functional configuration of control using a conventional block projection method.

FIG. 8 is a flowchart showing a calculation procedure of control using a conventional high-speed projection method.

FIG. 9 is a flowchart showing a calculation procedure of control using a conventional block method.

FIG. 10 is a diagram showing a functional configuration of control using a block projection method having the same output as a conventional projection method for correcting a time filter.

FIG. 11 is a flowchart showing a calculation procedure of control using a block projection method having the same output as a conventional projection method for correcting a time filter.

12A is a diagram illustrating a functional configuration in which the present invention is applied to an acoustic echo canceller, and FIG. 12B is a diagram illustrating a functional configuration in which the present invention is applied to noise control.

Claims (1)

[Claims] 1. A method of performing variable filter processing on an input signal, obtaining an error signal between the output signal and a desired signal, and adaptively controlling characteristics of the variable filter so as to minimize the power of the error signal. A first step of obtaining a prefilter coefficient vector by performing an operation corresponding to a product of an inverse matrix of an autocorrelation matrix of the input signal vector and the error signal vector; and smoothing the prefilter coefficient vector by smoothing. A second step of obtaining a coefficient vector, a third step of correcting an approximate filter coefficient vector in the direction of the input signal vector using the elements of the smooth coefficient vector as weighting coefficients, and the input signal vector and the approximate filter coefficient A fourth step of calculating an inner product of the input signal vector and a fifth step of obtaining an autocorrelation vector of the input signal vector
A step of performing an inner product operation of the smoothing coefficient vector and the autocorrelation vector; a result of the inner product operation obtained in the fourth step;
A seventh step of adding the inner product operation result obtained in the step to obtain the output signal, wherein the modification of the third step and the operation of the fourth step are performed at N ₁ and N ₂ times, respectively. (N ₁ and N ₂ are integers of 2 or more,
N ₁ ≠ N ₂ ), the steps other than the third and fourth steps are performed at each time, and each time the third step is performed, the third step is performed. At each time until the calculation, the order of the smoothing coefficient vector obtained in the second step is increased by one, and the calculation at the calculation time of the next third step is returned to the original expression. The autocorrelation vector of the order equal to the maximum order of the smoothing coefficient vector is calculated, and the sixth step uses an order equal to the order of the smoothing coefficient vector from the autocorrelation vector obtained in the fifth step. Adaptive control method.