Digital all-pass filter design method based on L-M algorithm
阅读说明：本技术 一种基于l-m算法的数字全通滤波器设计方法 (Digital all-pass filter design method based on L-M algorithm ) 是由 赵禹 孟婕 叶芃 杨扩军 张沁川 黄武煌 于 2021-06-22 设计创作，主要内容包括：本发明公开了一种基于L-M算法的数字全通滤波器设计方法,先根据设计要求设置数字全通滤波器的目标群时延,再通过计算目标群时延和GD-(sum),进而设置数字全通滤波器的阶数L；然后计算增加固定偏置后的数字全通滤波器的时延设计目标通过将划分成L个频率区间FB-(l),确定各个FB-(l)的边界值接着,通过每个FB-(l)区间对应数字全通滤波器二阶节的极点模值及相角构造向量,进而计算出数字全通滤波器的群时延最后以最小均方误差E(U)准则,利用Levenberg-Marquardt算法对A-(l)和θ-(l)的参数进行优化,完成数字全通滤波器设计。(The invention discloses a digital all-pass filter design method based on an L-M algorithm, which comprises the steps of firstly setting a target group delay of a digital all-pass filter according to design requirements, and then calculating the target group delay and GD sum Further setting the order L of the digital all-pass filter; then calculating the time delay design target of the digital all-pass filter after adding the fixed offset By mixing Divided into L frequency intervals FB l Determining each FB l Boundary value of Then, through each FB l The interval corresponds to the pole module value and the phase angle construction vector of the second order section of the digital all-pass filter, and then the group delay of the digital all-pass filter is calculated Finally, using Levenberg-Marquardt algorithm to pair A according to the minimum mean square error E (U) criterion l And theta l The parameters are optimized to complete the design of the digital all-pass filter.)
1. A digital all-pass filter design method based on an improved graphical method is characterized by comprising the following steps:
(1) and a target group delay GD with a digital all-pass filtergoal(omega) or GDgoal(ωk) Wherein ω represents a uniform distribution in [0, π]Inner continuous digital angular frequency, omegakRepresents a uniform distribution in [0, π]K discrete digital angular frequencies within;
(2) calculating target group delay and GDsum；
Or
Wherein Δ ω is the frequency interval of the discrete digital angular frequency;
(3) setting the order L of the digital all-pass filter, wherein the selection of the order L satisfies the following conditions: GD (GD) devicesum≤2πL；
(4) Calculating group delay and GD according to order L of digital all-pass filtersumFixed offset value GD ofextra；
(5) Calculating the time delay design target of the digital all-pass filter after adding the fixed offset
(6) Using single second order digital all-pass filter group time delay at [0, pi]The characteristic that the internal integral is constant to 2 pi isDivided into L frequency intervals FBlThe integration or summation of the group delays in each frequency interval is 2 pi and the individual FBs are determinedlBoundary value of
(7) Calculate each FBlInterval corresponds to pole modulus of second order section of digital all-pass filterAnd phase angleWhereinAnd
(8) define vector U ═ a1,…,Al,…,AL,θ1,…,θl,…,θL]Using F (A)l) Replacing pole modulus Ml，AlAn argument, F (A), representing the pole modulus of the mapped ith second order sectionl)∈(0,1)，-∞<Al＜∞；
The pole modulus M of each second order section of the digital all-pass filterlBy using the lower partCarrying out mapping conversion on the formula;
wherein F (-) represents a mapping function;
(9) after the mapping of the above formula, the stable interval of the second order section of the digital all-pass filter is formed by Ml<1 conversion to- ∞<Al<Infinity to calculate the group delay of the digital all-pass filter
(10) Constructing a nonlinear minimum mean square error equation E (U);
(11) using Levenberg-Marquardt algorithm to A according to the minimum mean square error E (U) criterionlAnd thetalOptimizing the parameters;
(11.1) setting the maximum iteration number P, initializing the current iteration number P to be 1, and setting the initial iteration starting point of the Levenberg-Marquardt algorithm to be the variable U solved in the step (7), namely, setting U to be the variable U1＝U；
(11.2) calculating the update vector in the p iterationWherein, the matrixλpIteration control parameters of a Levenberg-Marquardt algorithm; j. the design is a squarepIs E (U)p) The jacobian matrix of (a) is,W＝diag{W(ω1),W(ω2),…,W(ωK) Is a weighting matrix, Dp＝diag{HpIs composed of HpA diagonal matrix of diagonal elements;
(11.3) updating the vector U after the p iterationp+1＝Up+Δp；
(11.4) calculating the minimum mean square error E (U) after the p iterationp+1) And with E (U)p) Comparing the sizes, if E (U)p+1)≤E(Up) Then accept the vector Up+1And update λp+1＝λp10; otherwise, the vector U is heldpUnchanged and as the vector U after the p-th iterationp+1Update λp+1＝λp×10；
(11.5) judging whether the current iteration number P reaches the maximum iteration number P, if not, making P equal to P +1, and then returning to the step (11.2); otherwise, entering the step (11.6);
(11.6) outputting the vector UPReading the vector UPIn (A) corresponds tolAnd thetalThen, using the formula (4) to convert AlConversion to MlThereby constructing a form H (z) of digital all-pass filter second-order cascade;
wherein the content of the first and second substances,is the pole of the first second order section, al1、al2The filter coefficient of the first second order section of the all-pass filter is obtained;
so far, the design of the digital all-pass filter is completed.
Technical Field
The invention belongs to the technical field of digital filters, and particularly relates to a digital all-pass filter design method based on an L-M (Levenberg-Marquardt, Levenberg Marquardt) algorithm.
Background
In an analog-digital hybrid electronic system, when an electric signal is transmitted in the system, certain nonlinear phase distortion is introduced due to the non-ideal characteristics of a channel or an analog circuit, different phase shifts are added to signals with different frequencies, so that the system output is greatly influenced, and the signals are distorted in severe cases. Therefore, for systems requiring strict linearity in phase, such as phased array radar, digital oscilloscope, etc., non-linearity compensation of phase is crucial.
The digital all-pass filter can change the phase characteristics of the digital signals, and further effectively solves the problem of nonlinear distortion of the phase of an electronic system. The gain of the digital all-pass filter is 0dB for all frequency components of the input signal, so that the phase characteristics of the input signal are changed without attenuating any frequency signal. The digital all-pass filter plays an important role in application occasions such as phase compensation, group delay equalization and the like to meet the requirement of a system on linear phase-frequency response.
Disclosure of Invention
The invention aims to overcome the defects of the prior art and provides a digital all-pass filter design method based on an L-M algorithm, which realizes the coefficient design of an all-pass filter with specific group delay characteristics by adjusting the pole distribution of each cascade second-order node.
In order to achieve the above object, the present invention provides a method for designing a digital all-pass filter based on an L-M algorithm, comprising the steps of:
(1) and a target group delay GD with a digital all-pass filtergoal(omega) or GDgoal(ωk) Wherein ω represents a uniform distribution in [0, π]Inner continuous digital angular frequency, omegakRepresents a uniform distribution in [0, π]K discrete digital angular frequencies within;
(2) calculating target group delay and GDsum；
Or
Wherein Δ ω is the frequency interval of the discrete digital angular frequency;
(3) setting the order L of the digital all-pass filter, wherein the selection of the order L satisfies the following conditions: GD (GD) devicesum≤2πL；
(4) Calculating group delay and GD according to order L of digital all-pass filtersumFixed offset value GD ofextra；
(5) Calculating the time delay design target of the digital all-pass filter after adding the fixed offset
(6) Using single second order digital all-pass filter group time delay at [0, pi]The characteristic that the internal integral is constant to 2 pi isDivided into L frequency intervals FBlThe integration or summation of the group delays in each frequency interval is 2 pi and the individual FBs are determinedlBoundary value of
(7) Calculate each FBlInterval corresponds to pole modulus of second order section of digital all-pass filterAnd phase angleWhereinAndbeta is a set threshold;
(8) define vector U ═ a1,…,Al,…,AL,θ1,…,θl,…,θL]Using F (A)l) Replacing pole modulus Ml，AlAn argument, F (A), representing the pole modulus of the mapped ith second order sectionl)∈(0,1)，-∞<Al＜∞；
The pole modulus M of each second order section of the digital all-pass filterlCarrying out mapping conversion by using the following formula;
wherein F (-) represents a mapping function;
(9) after the mapping of the above formula, the stable interval of the second order section of the digital all-pass filter is formed by Ml<1 conversion to- ∞<Al<Infinity to calculate the group delay of the digital all-pass filter
(10) Constructing a nonlinear minimum mean square error equation E (U);
(11) using Levenberg-Marquardt algorithm to A according to the minimum mean square error E (U) criterionlAnd thetalOptimizing the parameters;
(11.1) setting the maximum iteration number P, initializing the current iteration number P to be 1, and setting the initial iteration starting point of the Levenberg-Marquardt algorithm to be the variable U solved in the step (7), namely, setting U to be the variable U1＝U；
(11.2) calculating an update vector delta at the p iterationp＝(Hp+λpDp)-1JpWRpWherein, the matrixλpIteration control parameters of a Levenberg-Marquardt algorithm; j. the design is a squarepIs E (U)p) Jacobian moment, W ═ diag { W (ω)1),W(ω2),…,W(ωK) Is a weighting matrix, Dp＝diag{HpThe method is a diagonal matrix;
(11.3) updating the vector U after the p iterationp+1＝Up+Δp；
(11.4) calculating the minimum mean square error E (U) after the p iterationp+1) And with E (U)p) Comparing the sizes, if E: (Up+1)≤E(Up) Then accept the vector Up+1And update λp+1＝λp10; otherwise, the vector U is heldpUnchanged and as the vector U after the p-th iterationp+1Update λp+1＝λp×10；
(11.5) judging whether the current iteration number P reaches the maximum iteration number P, if not, making P equal to P +1, and then returning to the step (11.2); otherwise, entering the step (11.6);
(11.6) outputting the vector UPReading the vector UPIn (A) corresponds tolAnd thetalThen, using the formula (4) to convert AlConversion to MlThereby constructing a form H (z) of digital all-pass filter second-order cascade;
wherein, Z represents a Z domain,is the pole of the first second order section, al1、al2The filter coefficient of the first second order section of the all-pass filter is obtained;
so far, the design of the digital all-pass filter is completed.
The invention aims to realize the following steps:
the invention relates to a digital all-pass filter design method based on an L-M algorithm, which comprises the steps of firstly setting a target group delay of a digital all-pass filter according to design requirements, and then calculating the target group delay and GDsumFurther setting the order L of the digital all-pass filter; then calculating the time delay design target of the digital all-pass filter after adding the fixed offsetBy mixingDivided into L frequency intervals FBlDetermining each FBlBoundary value ofThen, through each FBlThe interval corresponds to the pole module value and the phase angle construction vector of the second order section of the digital all-pass filter, and then the group delay of the digital all-pass filter is calculatedFinally, using Levenberg-Marquardt algorithm to pair A according to the minimum mean square error E (U) criterionlAnd thetalThe parameters are optimized to complete the design of the digital all-pass filter.
Meanwhile, the digital all-pass filter design method based on the L-M algorithm further has the following beneficial effects:
(1) the characteristic that the integral of the second-order nodal group time delay of the all-pass filter is constant to 2 pi in [0, pi ] is utilized, the target group time delay of the all-pass filter is cut into a plurality of frequency bands, and the module value and the phase angle value of the pole of the all-pass filter are obtained through simple calculation on the basis of the frequency bands, so that the design difficulty of all-pass filtering is simplified.
(2) The introduced mapping mechanism transforms the pole modulus MlIs mapped as AlSo as to change the stable interval of the second order section of the digital all-pass filter from Ml<1 conversion to- ∞<Al<Infinity, then the all-pass filter can be designed using an unconstrained optimization algorithm with assurance that the all-pass filter is stable.
(3) The introduction of unconstrained optimization algorithms such as LM greatly improves the design precision of the all-pass filter, so that the digital all-pass filter can obtain higher compensation precision under the same order; the introduction of the weighting function enables the algorithm to perform an accurate optimization process for the concerned frequency range, and the flexibility of the algorithm is increased.
Drawings
FIG. 1 is a flow chart of the digital all-pass filter design method based on the L-M algorithm of the present invention;
FIG. 2 is a graph of target group delayA frequency slicing diagram of (a);
FIG. 3 is a graph illustrating a second order group delay bound value;
FIG. 4 is a schematic diagram of target group delay frequency slicing after adding a fixed group delay;
FIG. 5 is a diagram of the pole modulus M of an all-pass filterlA mapping process schematic diagram;
fig. 6 is a schematic diagram of a second order cascade structure of an all-pass filter.
Detailed Description
The following description of the embodiments of the present invention is provided in order to better understand the present invention for those skilled in the art with reference to the accompanying drawings. It is to be expressly noted that in the following description, a detailed description of known functions and designs will be omitted when it may obscure the subject matter of the present invention.
Examples
FIG. 1 is a flow chart of the digital all-pass filter design method based on the L-M algorithm.
In this embodiment, as shown in fig. 1, the method for designing a digital all-pass filter based on an L-M algorithm of the present invention includes the following steps:
s1, firstly, according to the actually designed digital all-pass filter, the target group delay GD of the digital all-pass filter needs to be setgoal(omega) or GDgoal(ωk) Wherein ω represents a uniform distribution in [0, π]Inner continuous digital angular frequency, omegakRepresents a uniform distribution in [0, π]K discrete digital angular frequencies within;
s2, calculating the target group delay and GDsum；
Or
Wherein Δ ω is the frequency interval of the discrete digital angular frequency;
s3, setting the order L of the digital all-pass filter, wherein the selection of L satisfies the following conditions: GD (GD) devicesum≤2πL；
S4, in order to obtain higher design precision of the all-pass filter group delay, the group delay and GD can be calculated according to the order L of the digital all-pass filtersumFixed offset value GD ofextra；
In this embodiment, GD due to increased fixed group delay offset, as shown in fig. 4extraIs constant, i.e. introduces an additional linear phase difference, which does not change the nonlinear phase characteristics of the system.
S5, calculating the time delay design target of the digital all-pass filter with the fixed offset added
S6, as shown in figure 2, using a single second order digital all-pass filter group delay at [0, π]The characteristic that the internal integral is constant to 2 pi isDivided into L frequency intervals FBlReferring to equation (1), the integration or summation of the group delays in each frequency interval is 2 π and the individual FBs are determinedlAs shown in FIG. 2, the l-th frequency interval FBlIn the range of
S7, making each frequency interval correspond to a second-order section of the all-pass filter, each FB can be calculatedlInterval corresponds to pole modulus of second order section of digital all-pass filterAnd phase angleWhereinAndbeta is a shape parameter, and the value of beta determines the delay peak point tau of the ith all-pass second order nodal groupl(θl) Group delay corresponding to frequency interval boundaryAndthe proportional relationship of (a) is shown in FIG. 3, and according to experience, the value range of beta is [0.75,0.95 ]]；
S8, definition vector U ═ a1,L,Al,L,AL,θ1,L,θl,L,θL]Using F (A)l) Replacing pole modulus Ml，AlAn argument, F (A), representing the pole modulus of the mapped ith second order sectionl)∈(0,1)，-∞<Al＜∞；
As shown in fig. 5, the pole modulus M of each second order section of the digital all-pass filter is measuredlCarrying out mapping conversion by using the following formula;
wherein F (-) represents a mapping function;
s9, after the mapping of the above formula is carried out, the stable interval of the second order section of the digital all-pass filter is formed by Ml<1 conversion to- ∞<Al<Infinity to calculate the group delay of the digital all-pass filter
At this time, no matter AlAnd thetalThe stability of the digital all-pass filter can be ensured by taking which value.
S10, constructing a nonlinear minimum mean square error equation E (U);
s11, due to the non-linear function such as cos in the formula (5), the non-linear function is represented by a non-linear equation system with the argument of U, so that U can be further optimized and solved by using a Levenberg-Marquardt (L-M) algorithm by means of the minimum mean square error criterion, namely AlAnd thetalThe parameters are optimized, and the specific optimization process is as follows:
s11.1, setting the maximum iteration number P, initializing the current iteration number P to be 1, and enabling the initial iteration starting point of the L-M algorithm to be the variable U solved in the step (7), namely enabling U to be obtained1＝U；
S11.2, calculating an update vector delta in the p iterationp＝(Hp+λpDp)-1JpWRpWherein, the matrixλpIs an iterative control parameter, lambda, of the L-M algorithmpWhen the overall search is larger, the L-M algorithm is close to the gradient descent algorithm, and has higher overall search capability, namely lambdapWhen the local search algorithm is small, the L-M algorithm is close to a Gauss-Newton iteration method, and the local search capability is high; j. the design is a squarepIs E (U)p) The jacobian moment of; rpAs a residual vector, specifically satisfying: rp＝[Rp(ω1),Rp(ω2),…,Rp(ωK)]T，Dp＝diag{HpThe method is a diagonal matrix; w ═ diag { W (ω)1),W(ω2),…,W(ωK) Is a weighting matrix;
s11.3, updating the vector U after the p iterationp+1＝Up+Δp；
S11.4, calculating the minimum mean square error E (U) after the p iterationp+1) And with E (U)p) Comparing the sizes, if E (U)p+1)≤E(Up) Then accept the vector Up+1And update λp+1＝λp10; otherwise, the vector U is heldpUnchanged and as the vector U after the p-th iterationp+1Update λp+1＝λp×10；
S11.5, determining whether the current iteration number P reaches the maximum iteration number P, if not, making P equal to P +1, and returning to step S11.2; otherwise, go to step S11.6;
s11.6, outputting the vector UPReading the vector UPIn (A) corresponds tolAnd thetalThen, using the formula S4 to convert AlConversion to MlThereby constructing a form H (z) of digital all-pass filter second-order cascade;
wherein, Z represents a Z domain,is the pole of the first second order section, al1、al2The filter coefficient of the first second order section of the all-pass filter is obtained; the Z-domain expression of the digital all-pass filter is expanded into the form of a second order cascade as shown in fig. 6.
So far, the design of the digital all-pass filter is completed.
Although illustrative embodiments of the present invention have been described above to facilitate the understanding of the present invention by those skilled in the art, it should be understood that the present invention is not limited to the scope of the embodiments, and various changes may be made apparent to those skilled in the art as long as they are within the spirit and scope of the present invention as defined and defined by the appended claims, and all matters of the invention which utilize the inventive concepts are protected.