CN105719255B - Bandwidth varying linear-phase filter method based on Laplace structure - Google Patents
Bandwidth varying linear-phase filter method based on Laplace structure Download PDFInfo
- Publication number
- CN105719255B CN105719255B CN201610044267.XA CN201610044267A CN105719255B CN 105719255 B CN105719255 B CN 105719255B CN 201610044267 A CN201610044267 A CN 201610044267A CN 105719255 B CN105719255 B CN 105719255B
- Authority
- CN
- China
- Prior art keywords
- filter
- fliplr
- work
- composite filter
- square estimation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 12
- 239000002131 composite material Substances 0.000 claims abstract description 84
- 238000005457 optimization Methods 0.000 claims abstract description 23
- 230000004044 response Effects 0.000 claims description 60
- 210000003462 vein Anatomy 0.000 claims description 6
- 238000001914 filtration Methods 0.000 description 6
- 238000001228 spectrum Methods 0.000 description 6
- 238000004088 simulation Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 210000001367 artery Anatomy 0.000 description 2
- 230000015572 biosynthetic process Effects 0.000 description 2
- 230000006835 compression Effects 0.000 description 2
- 238000007906 compression Methods 0.000 description 2
- 230000007812 deficiency Effects 0.000 description 2
- 239000011159 matrix material Substances 0.000 description 2
- 238000003786 synthesis reaction Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 238000000354 decomposition reaction Methods 0.000 description 1
- 238000012938 design process Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000003708 edge detection Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 235000015170 shellfish Nutrition 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- H—ELECTRICITY
- H03—ELECTRONIC CIRCUITRY
- H03H—IMPEDANCE NETWORKS, e.g. RESONANT CIRCUITS; RESONATORS
- H03H17/00—Networks using digital techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20172—Image enhancement details
- G06T2207/20182—Noise reduction or smoothing in the temporal domain; Spatio-temporal filtering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
-
- H—ELECTRICITY
- H03—ELECTRONIC CIRCUITRY
- H03H—IMPEDANCE NETWORKS, e.g. RESONANT CIRCUITS; RESONATORS
- H03H17/00—Networks using digital techniques
- H03H2017/0072—Theoretical filter design
- H03H2017/0081—Theoretical filter design of FIR filters
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Filters That Use Time-Delay Elements (AREA)
- Complex Calculations (AREA)
Abstract
The present invention settles sth. according to policy or law a kind of design method of the filter suitable for laplacian pyramid structure, mainly solves the filter of existing laplacian pyramid structure, in addition to Ha Er filters, other filters cannot have the problem of linear phase and orthogonal property simultaneously.Its technical solution is：First set filter length N, cut-off frequecy of passband ωpWith stopband initial frequency ωs, decimation factor M；Then according to these parameters, the firpm functions in MatLab is called to generate ptototype filter p；To the half coefficient of the ptototype filter, fmincon functions are called to optimize, and composite filter is obtained according to the result of optimization；Relationship is finally overturn according to time domain and can be obtained resolution filter.The present invention can not only make filter have orthogonal and linear phase characteristic, and adaptive-bandwidth, wider application conditions are provided for the filter of laplacian pyramid structure by providing a kind of rational condition of orthogonal constraints.
Description
Technical field
The invention belongs to signal processing technology field, more particularly to the linear-phase filter design side of a kind of adaptive-bandwidth
Method can be used for image co-registration, increase, compression, denoising and edge detection.
Background technology
Laplacian pyramid structure LP is proposed at first by Burt et al., and image coding and compression are used primarily for.Later,
Laplacian pyramid structure LP is extended to a kind of multiple dimensioned, multiresolution analysis tool by people, is applied to image procossing
Various fields, such as image co-registration, image denoising.
Fig. 1 gives laplacian pyramid structure chart.In Fig. 1, for given input signal x, output c is x
Close approximation, d are the residual errors between original signal and its prediction signal p.In this structure, resolution filter h and synthetic filtering
Device g is typically a pair of of low-pass filter, and analysis matrix H and synthetical matrix G then constitute one group of transformation pair.
In traditional laplacian pyramid structure LP designs, resolution filter h and composite filter g the two filtering
Device takes the low frequency channel of multichannel perfect reconstruction filter group, filter bank structure as shown in Figure 2.In each channel,
Filter and its upper and lower decimation factor constitute positive inverse transformation pair.Therefore, in filter group the filter of low frequency channel usually by
It directly applies in laplacian pyramid structure LP.It is carried although this way is the structure of laplacian pyramid structure LP
It has supplied to help, but has also limited the performance of laplacian pyramid structure LP median filters.
Laplacian pyramid structure LP is played an important role in the various fields of image procossing, but in its design side
There is also a critical issues in face, i.e.,：Filter in existing laplacian pyramid structure LP is taken from complete in multichannel
Low frequency channel in full reconfigurable filter group.This makes the filter in laplacian pyramid structure LP be filtered by Perfect Reconstruction
The constraint of wave device group itself, certain characteristics are restricted, including the very important linear phase in image procossing
Characteristic.In the design of perfect reconstruction filter group, although linear phase characteristic is easy to realize in two channels filter group,
But be often difficult to realize in multichannel perfect reconstruction filter group, cause be except bandwidthTwo channels filter group with
Outside, the problem of perfect reconstruction filter group of other bandwidth is difficult with linear phase characteristic.So in laplacian pyramid
In structure LP, it is taken from this indirect filter design method of perfect reconstruction filter group, it will cause Laplce golden
Filter in word tower structure LP is limited by perfect reconstruction filter group and is difficult while reaching bandwidth varying and linearly
The deficiency of phase.
Invention content
It is an object of the invention to the deficiencies for above-mentioned prior art, propose a kind of bandwidth varying based on Laplace structure
Linear-phase filter method, to break away from the design constraint of perfect reconstruction filter group so that the frequency band of filter is variable, and
With linear phase property while with orthogonal property.
The technical proposal of the invention is realized in this way：
One, technical principles
The decomposition part of laplacian pyramid structure is as shown in Figure 1, when H and G meet HG=I relationships, LP structures
In this channel constitute the subsystem of a Perfect Reconstruction, can realize the Perfect Reconstruction to sub- spacing wave.According to H
This relationship of G=I can derive that resolution filter h and composite filter g need the condition met, can indicate as follows,
Wherein, N is the length of filter, and M is decimation factor, and the bandwidth of filter isM is shift count, and value is
Arbitrary integer, δ (m) are dirac sequences,
Above-mentioned condition enables to the filter of bandwidth varying while having orthogonal and linear phase characteristic.
Two, technical solutions
According to above-mentioned principle, technical scheme is as follows：
Technical solution 1：A kind of low pass filter design method of the bandwidth varying linear phase based on Laplace structure, including
To low pass composite filter gLDesign and to low pass resolution filter hLDesign：
(1) setting low pass composite filter gLLength NL, cut-off frequecy of passbandStopband initial frequencyBandwidth isWherein M is decimation factor, M and NLIt is integer, M >=2, NLFor the integral multiple of M,
(2) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the lowpass prototype filter p of linear phaseL
Least square estimation pL(n), n=0,1,2 ... NL-1；
(3) lowpass prototype filter p is determinedLHalf coefficient qL(k)：
Work as NLFor odd number when,
Work as NLFor even number when,
(4) by qL(k) initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gL0(n),gL0(n-mM)>=δ (m)
Wherein, gL0(n) it is low-pass filter impulse response free variable, works as NLFor even number when, gL0(n)=[qL(k),
fliplr(qL(k))], work as NLFor odd number when, gL0(n)=[qL(k),fliplr(qL(k-1)], fliplr is sequence in MatLab
Left and right overturning function, GL(ejw) it is gL0(n) frequency response, α are weights, and value is 0 ＜ α ＜ 1, and m is shift count, and value is
Arbitrary integer, δ (m) are dirac sequences,
In the case where meeting constraints, weight α is adjusted, as object function φLWhen value reaches minimum, after obtaining optimization
Low-pass reconstruction filters gLHalf coefficient QL(k)；
(5) according to the half coefficient Q of the low pass composite filter after optimizationL(k), low pass composite filter g is obtainedLArteries and veins
Rush response sequence gL(n)：
Work as NLFor odd number when, gL(n)=[QL(k),fliplr(QL(k-1))]；
Work as NLFor even number when, gL(n)=[QL(k),fliplr(QL(k))]；
(6) according to the Least square estimation g of low pass composite filterL(n) with the impulse response sequence of low pass resolution filter
Arrange hL(n) meet time domain inverse relation：hL(n)=gL(NL- 1-n), by the Least square estimation g of low pass composite filterL(n) i.e.
The Least square estimation h of low pass resolution filter can be found outL(n)。
Technical solution 2：A kind of Design of Bandpass method of the bandwidth varying linear phase based on Laplace structure, including
To band logical composite filter gBDesign and to band logical resolution filter hBDesign：
1) setting band logical composite filter gBLength NB, cut-off frequecy of passband is respectivelyStopband initial frequency
RespectivelyDecimation factor is M, and bandwidth isWherein M and NBIt is integer, M >=2, NBIt is the integral multiple of M,
2) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the band logical ptototype filter p of linear phaseB
Least square estimation pB(n), n=0,1,2 ... NB-1；
3) band logical ptototype filter p is determinedBHalf coefficient qB(k)：
Work as NBFor odd number when,
Work as NBFor even number when,
4) by qB(k) initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gB0(n),gB0(n-mM)>=δ (m)
Wherein, gB0(n) it is bandpass filter impulse response free variable, works as NBFor even number when, gB0(n)=[qB(k),
fliplr(qB(k))], work as NBFor odd number when, gB0(n)=[qB(k),fliplr(qB(k-1)], fliplr is sequence in MatLab
Left and right overturning function, GB(ejw) it is gB0Frequency response, α is weight, and value is 0 ＜ α ＜ 1, and m is shift count, and value is to appoint
Meaning integer, δ (m) is dirac sequence,
In the case where meeting constraints, weight α is adjusted, as object function φBValue reaches minimum, after being optimized
The half coefficient Q of low pass composite filterB(k)；
5) according to the half coefficient Q of the band logical composite filter after optimizationB(k), band logical composite filter g is obtainedBPulse
Response sequence gB(n)：
Work as NBFor odd number when, gB(n)=[QB(k),fliplr(QB(k-1))]；
Work as NBFor even number when, gB(n)=[QB(k),fliplr(QB(k))]；
6) according to the Least square estimation g of band logical composite filterB(n) with the Least square estimation of band logical resolution filter
hB(n) meet time domain inverse relation：hB(n)=gB(NB- 1-n), by the Least square estimation g of band logical composite filterB(n)
Find out the Least square estimation h of band logical resolution filterB(n)。
Technical solution 3：A kind of high-pass filter design method of the bandwidth varying linear phase based on Laplace structure, including
To high pass composite filter gHDesign and to high pass resolution filter hHDesign：
(A) setting high pass composite filter gHLength NH, cut-off frequecy of passband isStopband initial frequency isBandwidth isWherein M is decimation factor, M and NHIt is integer, M >=2, NHFor the integral multiple of M,
(B) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the high pass ptototype filter p of linear phaseH
Least square estimation pH(n), n=0,1,2 ... NH-1；
(C) high pass ptototype filter p is determinedHHalf coefficient qH(k)：
Work as NHFor odd number when,
Work as NHFor even number when,
(D) by qHThe initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gH0(n),gH0(n-mM)>=δ (m)
Wherein, gH0(n) it is high-pass filter impulse response free variable：Work as NHFor even number when, gH0(n)=[qH(k),-
fliplr(qH], (k)) work as NHWhen being strange, gH0(n)=[qH(k),fliplr(qH(k-1)], fliplr is that sequence is left in MatLab
Right overturning function, GH(ejw) it is gH0(n) frequency response, α are weights, and value is 0 ＜ α ＜ 1, and m is shift count, and value is to appoint
Meaning integer, δ (m) is dirac sequence,
In the case where meeting constraints, weight α is adjusted, as object function φHWhen value reaches minimum, after obtaining optimization
High pass composite filter half coefficient QH(k)；
(E) according to the half coefficient Q of the high pass composite filter after optimizationH(k), high pass composite filter g is obtainedHArteries and veins
Rush response sequence gH(n)：
Work as NHFor odd number when, gH(n)=[QH(k),fliplr(QH(k-1))],
Work as NHFor even number when, gH(n)=[QH(k),-fliplr(QH(k))]；
(F) according to the Least square estimation g of high pass composite filterH(n) with the impulse response sequence of high pass resolution filter
Arrange hH(n) meet time domain inverse relation：hH(n)=gH(NH- 1-n), by the Least square estimation g of high pass composite filterH(n) i.e.
The Least square estimation h of high pass resolution filter can be found outH(n)。
Compared with the prior art, the present invention has the following advantages：
First, the filter that the present invention designs can have linear phase and orthogonal property simultaneously.
The existing western Daubechies filters of more shellfishes all only have orthogonal property in addition to Ha Er Haar filters, without
With linear phase characteristic；Though biorthogonal Biorthogonal filters have linear phase characteristic, do not have orthogonal property.Its
Its all filter designed by perfect reconstruction filter group all only has orthogonal property and line in addition to Haar filters
One of property phase characteristic, and cannot have the two characteristics simultaneously.The present invention breaks away from setting for perfect reconstruction filter group
Meter constraint directly designs the filter in laplacian pyramid structure LP.When resolution filter h and reconfigurable filter g meet
Condition：When, filter can have linear phase and orthogonal property simultaneously；
Second, adaptive-bandwidth of the invention.
Existing Ha Er Haar filter bandwidhts areThe present invention filter bandwidht can beI.e. spectral limit is
Third, design process of the present invention are simple.
Existing Ha Er Haar filters be composite filter g and resolution filter h are designed, and the present invention
Filter is only designed with to composite filter g, is then filtered according to the Least square estimation h (n) of resolution filter h and synthesis
Relationship h (n)=g (N-1-n) of the Least square estimation g (n) of wave device g, you can acquire resolution filter h, greatly reduce
The complexity of design.
Description of the drawings
Fig. 1 is that existing Laplce's structure LP decomposes part-structure figure；
Fig. 2 is existing multichannel perfect reconstruction filter group structure chart；
Fig. 3 is the design flow diagram of the present invention；
Fig. 4 is that the present invention is based on the bandwidth of pulling structureThe performance simulation figure for the low-pass filter that length is 24；
Fig. 5 is that the present invention is based on the bandwidth of pulling structureThe performance simulation figure for the bandpass filter that length is 60；
Fig. 6 is that the present invention is based on the bandwidth of pulling structureThe performance simulation figure for the high-pass filter that length is 60.
Specific implementation mode
Technical solutions and effects of the present invention is described in further detail below in conjunction with attached drawing.
The filter that the present invention designs is to be based on laplacian pyramid structure, as shown in Figure 1, its filter includes decomposing
Filter and composite filter.Input signal x passes through resolution filter h, upper and lower decimation factor M and composite filter g successively
Afterwards, the approximation signal p of input signal x is obtained, the residual error between input signal and its approximation signal p is d.
Example 1：Designing the bandwidth based on Laplace structure isThe low-pass filter that length is 24.
This example includes to low pass composite filter gLDesign and to low pass resolution filter hLDesign two parts.One,
Design low-pass reconstruction filters gL
Step 1, setting low pass composite filter gLParameter.
If length NL=24, decimation factor M=4, cut-off frequecy of passbandStopband initial frequencyWherein rLIt is low pass intermediate zone adjustment parameter, by changing rLSize, adjustment intermediate zone can be facilitated
Width, this example takes rL=0.4.
Step 2, lowpass prototype filter p is determinedL。
Lowpass prototype filter pLA substantially low-pass impulse response sequence pL(n), this example is by calling MatLab
In firpm functions obtain linear phase low pass ptototype filter pLLeast square estimation pL(n), n=0,1 ... NL-1。
The firpm functions are a functions for being designed with limit for length's Least square estimation FIR filter,《Digital signal
Processing》It is described in book, the filter that firpm functions are designed all has linear phase characteristic.
Step 3, lowpass prototype filter p is determinedLThe Least square estimation q of halfL(k)。
Due to lowpass prototype filter pLWith linear phase characteristic, i.e. its Least square estimation pL(n) it is central symmetry
, in order to reduce computation complexity, therefore only to pL(n) Least square estimation of half, which optimizes, can obtain low pass conjunction
At filter gLLeast square estimation gL(n)。
Take the Least square estimation p of lowpass prototype filterL(n) half Least square estimation qL(k), wherein working as NLIt is
When even number,Work as NLWhen being odd number,For example, setting NL=4, pL(n)=[pL(0)
pL(1) pL(2) pL(3)] p, is takenL(n) Least square estimation of half is：N is taken in this exampleL=24, therefore the arteries and veins of lowpass prototype filter
Rush response sequence pL(n) half Least square estimation
Step 4, the parameter of setting majorized function fmincon.
If the initial value of majorized function fmincon is qL(k), constraints is ＜ gL0(n),gL0(n-mM) ＞=δ (m), mesh
Scalar functions areWherein：gL0(n) it is that low-pass filter pulse is rung
It is shift count to answer free variable, m, and value is arbitrary integer, gL0(n-mM) it is gL0(n) m shift sequence, α are weights, 0
＜ α ＜ 1, GL(ejω) it is gL0(n) frequency response, δ (m) are dirac sequences,
The majorized function fmincon is common majorized function in a kind of MatLab；
gL0(n) value is according to NLParity depending on：
Work as NLWhen being even number, gL0(n)=[qL(k),fliplr(qL(k))],
Work as NLWhen being odd number, gL0(n)=[qL(k),fliplr(qL(k-1))],
Fliplr is that sequence or so overturns function in MatLab, such as sets qL(k)=[qL(0) qL(1) qL(2) qL
(3)], NL=4, by sequence qL(k) left and right overturning, is represented by fliplr (qL(k))=[qL(3) qL(2) qL(1) qL
(0)]；
Constraints ＜ gL0(n),gL0(n-mM) ＞=δ (m) shows gL0(n) and gL0(n-mM) sequence is orthogonal, i.e.,
Object function is passband and stopband to filter while optimizing, to the passband of filter when weight α determines to optimize
With the limited degree of stopband, α show more greatly to passband limit it is bigger, i.e., passband is more smooth, and stop band ripple is bigger.
Step 5, the half impulse response q of the impulse response of lowpass prototype filterL(k) low pass composite filter g is determinedL
Half coefficient QL(k)。
To the impulse response p of lowpass prototype filterL(n) half impulse response qL(k), pass through majorized function fmincon
Optimization adjusts the size of weight α value, this example takes α=0.8, in the case where meeting constraints, as object function φLMost
Hour, function return value is low pass composite filter gLHalf coefficient QL(k)。
Step 6, low pass composite filter g is determinedLLeast square estimation gL(n)。
Due to low pass composite filter gLIt is linear-phase filter, Least square estimation gL(n) be it is centrosymmetric,
It therefore can be by low pass composite filter gLHalf coefficient QL(k) low pass composite filter g is determinedLLeast square estimation gL
(n)：
Work as NLWhen being even number, gL(n)=[QL(k),fliplr(QL(k))],
Work as NLWhen being odd number, gL(n)=[QL(k),fliplr(QL(k-1))]。
Two, design low pass resolution filter hL
Step 7, the Least square estimation h of low pass resolution filter is determinedL(n)。
Set the impulse response g of low pass composite filterL(n) with the Least square estimation h of low pass resolution filterL(n) full
Sufficient time domain inverse relation：hL(n)=gL(NL- 1-n), the pulse that low pass resolution filter can be obtained according to the time domain inverse relation is rung
Answer sequences hL(n), i.e. low pass resolution filter hL。
Such as set gL(n)=[gL(0) gL(1) gL(2) gL(3)],NL=4, gL(n) and hL(n) meet time domain reversion to close
It is, then hL(n)=gL(NL- n-1)=[gL(3) gL(2) gL(1) gL(0)]。
By above step, can obtain be based on the bandwidth of Laplace structureThe low-pass filter that length is 24, frequency spectrum branch
Support is ranging fromThe performance of the low-pass filter is as shown in Figure 4.Wherein：
Fig. 4 (a) is low pass composite filter gLTime-domain pulse response,
Fig. 4 (b) is low pass composite filter gLFrequency domain response,
Fig. 4 (c) is low pass resolution filter hLTime-domain pulse response,
Fig. 4 (d) is low pass resolution filter hLFrequency domain response.
It can be seen that the low-pass filter belongs to linear-phase filter, demonstrates this example design from Fig. 4 (a) and 4 (c)
Low pass composite filter Least square estimation gL(n) meet ＜ gL(n),gL(n-mM) ＞=δ (m), therefore the low pass synthetic filtering
Utensil has orthogonality, and due to the impulse response g of low pass composite filterL(n) with the impulse response sequence of low pass resolution filter
Arrange hL(n) meet time domain inverse relation：hL(n)=gL(NL- 1-n), therefore the impulse response h of low pass resolution filterL(n) also full
Sufficient ＜ hL(n),hL(n-mM) ＞=δ (m), i.e. the low pass resolution filter have orthogonality.It can be seen that from Fig. 4 (b) and 4 (d), it should
Low-pass filter frequency spectrum support range beI.e. bandwidth is
Example 2 designs the bandwidth based on Laplace structureThe bandpass filter that length is 60.
This example includes to band logical composite filter gBDesign and to band logical resolution filter hBDesign two parts.1、
Design band logical composite filter gB
Step 1, setting band logical composite filter gBParameter.
If length NB=60, decimation factor M=3, cut-off frequecy of passband are respectively Stopband initial frequency is respectivelyWherein rBIt is
Band logical intermediate zone adjustment parameter, by changing rBSize, can facilitate adjustment intermediate zone width, this reality rB=0.3.
Step 2 determines band logical ptototype filter pB。
Band logical ptototype filter pBA substantially band logical Least square estimation pB(n), this example is by calling MatLab
In firpm functions obtain linear phase, band and lead to ptototype filter pBLeast square estimation pB(n), n=0,1 ... NB-1。
Step 3 determines band logical ptototype filter pBThe Least square estimation q of halfB(k)。
Due to band logical ptototype filter pBWith linear phase characteristic, i.e. its Least square estimation pB(n) center is symmetrical
, in order to reduce computation complexity, therefore only to pB(n) Least square estimation of half, which optimizes, can obtain band logical synthesis
Filter gBLeast square estimation gB(n)。
Take the Least square estimation p of band logical ptototype filterB(n) half Least square estimation qB(k), wherein working as NBIt is
When even number,Work as NBWhen being odd number,For example, setting NB=4, pB(n)=[pB(0)
pB(1) pB(2) pB(3)] p, is takenB(n) Least square estimation of half is：N is taken in this exampleB=60, therefore the arteries and veins of band logical ptototype filter
Rush response sequence pB(n) half Least square estimation
Step 4, the parameter of setting majorized function fmincon.
If the initial value of majorized function fmincon is qB(k), constraints is<gB0(n),gB0(n-mM) ＞=δ (m), mesh
Scalar functions areWherein：gB0
(n) it is bandpass filter impulse response free variable, m is shift count, and value is arbitrary integer, gB0(n-mM) it is gB0(n)
M shift sequence, α are weight, 0 ＜ α ＜ 1, GB(ejω) it is gB0(n) frequency response, δ (m) are dirac sequences,
gB0(n) value is according to NBParity depending on：
Work as NBWhen being even number, gB0(n)=[qB(k),fliplr(qB(k))],
Work as NBWhen being odd number, gB0(n)=[qB(k),fliplr(qB(k-1))],
Fliplr is that sequence or so overturns function, constraints in MatLab<gB0(n),gB0(n-mM)>=δ (m) shows
gB0(n) and gB0(n-mM) sequence is orthogonal, i.e.,
Object function is passband and stopband to filter while optimizing, to the passband of filter when weight α determines to optimize
With the limited degree of stopband, α show more greatly to passband limit it is bigger, i.e., passband is more smooth, and stop band ripple is bigger.
Step 5, the half impulse response q of the impulse response of band logical ptototype filterB(k) low pass composite filter g is determinedB
Half coefficient QL(k)。
Optimized by majorized function fmincon, to the impulse response p of band logical ptototype filterB(n) half impulse response
qB(k), the size of weight α value is adjusted, this example takes α=0.7, in the case where meeting constraints, as object function φBMost
Hour, function return value is band logical composite filter gBHalf coefficient QB(k)。
Step 6 determines band logical composite filter gBLeast square estimation gB(n)。
Due to band logical composite filter gBIt is linear-phase filter, Least square estimation gB(n) be it is centrosymmetric,
It therefore can be by band logical composite filter gBHalf coefficient QB(k) band logical composite filter g is determinedBLeast square estimation gB
(n)：
Work as NBWhen being even number, gB(n)=[QB(k),fliplr(QB(k))],
Work as NBWhen being odd number, gL(n)=[QB(k),fliplr(QB(k-1))]。
2, design band logical resolution filter hB
Step 7 determines the Least square estimation h of band logical resolution filterB(n)。
According to the impulse response g of setting band logical composite filterB(n) with the Least square estimation h of band logical resolution filterB
(n) meet time domain inverse relation：hB(n)=gB(NB- 1-n), it can obtain the Least square estimation h of resolution filterB(n), i.e. band
Logical resolution filter hB。
Such as set gB(n)=[gB(0) gB(1) gB(2) gB(3)],NB=4, gB(n) and hB(n) meet time domain reversion to close
It is, then hB(n)=gB(NB- n-1)=[gB(3) gB(2) gB(1) gB(0)]。
By above step, can obtain be based on the bandwidth of Laplace structureThe bandpass filter that length is 60, frequency spectrum branch
Support is ranging fromThe performance of the bandpass filter is as shown in Figure 5.Wherein：
Fig. 5 (a) is band logical composite filter gBTime-domain pulse response,
Fig. 5 (b) is band logical composite filter gBFrequency domain response,
Fig. 5 (c) is band logical resolution filter hBTime-domain pulse response,
Fig. 5 (d) is band logical resolution filter hBFrequency domain response.
It can be seen that the bandpass filter belongs to linear-phase filter, demonstrates this example design from Fig. 5 (a) and 5 (c)
Band logical composite filter Least square estimation gB(n) meet<gB(n),gB(n-mM)>=δ (m), therefore the band logical synthetic filtering
Utensil has orthogonality, and because of the impulse response g of band logical composite filterB(n) with the Least square estimation h of band logical resolution filterB
(n) meet time domain inverse relation：hB(n)=gB(NB- 1-n), therefore the impulse response h of band logical resolution filterB(n) also meet<
hB(n),hB(n-mM)>=δ (m), i.e. the band logical resolution filter have orthogonality.It can be seen that from Fig. 5 (b) and 5 (d), the band logical
Filter spectrum support range beI.e. bandwidth is
Example 3. designs the bandwidth based on Laplace structureThe high-pass filter that length is 60.
This example includes to high pass composite filter gHDesign and to high pass resolution filter hHDesign two parts.1)、
Design high pass composite filter gH
Step A, setting high pass composite filter gHParameter.
If length NH=60, decimation factor M=3, cut-off frequecy of passbandStopband initial frequencyWherein rHIt is high pass intermediate zone adjustment parameter, by changing rHSize, adjustment transition can be facilitated
The width of band, this example take rH=0.2.
Step B determines high pass ptototype filter pH。
High pass ptototype filter pHA substantially high pass Least square estimation pH(n), this example is by calling MatLab
In firpm functions obtain linear phase high pass ptototype filter pHLeast square estimation pH(n), n=0,1 ..., NH-1。
Step C determines high pass ptototype filter pHThe Least square estimation q of halfH(k)。
Due to high pass ptototype filter pHWith linear phase characteristic, i.e. its Least square estimation pH(n) it is central symmetry
, in order to reduce computation complexity, therefore only to pH(n) Least square estimation of half, which optimizes, can obtain high pass conjunction
At filter gHLeast square estimation gH(n)。
Take the Least square estimation p of high pass ptototype filterH(n) half Least square estimation qH(k), wherein working as NHIt is
When even number,Work as NHWhen being odd number,For example, setting NH=4, pH(n)=[pH(0)
pH(1) pH(2) pH(3)] p, is takenH(n) Least square estimation of half is：N is taken in this exampleH=24, therefore high pass ptototype filter
Least square estimation pH(n) half Least square estimation
Step D, the parameter of setting majorized function fmincon.
If the initial value of majorized function fmincon is qH(k), constraints is<gH0(n),gH0(n-mM)>=δ (m), mesh
Scalar functions areWherein：gH0(n) it is high-pass filter pulse
Free variable is responded, m is shift count, and value is arbitrary integer, gH0(n-mM) it is gH0(n) m shift sequence, α are power
Weight, 0 ＜ α ＜ 1, GH(ejω) it is gH0(n) frequency response, δ (m) are dirac sequences,
gH0(n) value is according to NHParity depending on：
Work as NHWhen being even number, gH0(n)=[qH(k),-fliplr(qH(k))],
Work as NHWhen being odd number, gH0(n)=[qH(k),fliplr(qH(k-1))],
Fliplr is that sequence or so overturns function in MatLab, such as sets qH(k)=[qH(0) qH(1) qH(2) qH
(3)], NH=4, by sequence qH(k) left and right overturning, is represented by fliplr (qH(k))=[qH(3) qH(2) qH(1) qH
(0)]；
Constraints<gH0(n),gH0(n-mM)>=δ (m) shows gH0(n) and gH0(n-mM) sequence is orthogonal, i.e.,
Object function is passband and stopband to filter while optimizing, to the passband of filter when weight α determines to optimize
With the limited degree of stopband, α show more greatly to passband limit it is bigger, i.e., passband is more smooth, and stop band ripple is bigger.
Step E, according to the half impulse response q of the impulse response of high pass ptototype filterH(k) high pass synthetic filtering is determined
Device gHHalf coefficient QH(k)。
To the impulse response p of high pass ptototype filterH(n) half impulse response qH(k), pass through majorized function fmincon
Optimization adjusts the size of weight α value, this example takes α=0.8, in the case where meeting constraints, as object function φHMost
Hour, function return value is high pass composite filter gHHalf coefficient QH(k)。
Step F determines high pass composite filter gHLeast square estimation gH(n)。
Due to high pass composite filter gHIt is linear-phase filter, Least square estimation gH(n) be it is centrosymmetric,
It therefore can be by high pass composite filter gHHalf coefficient QH(k) it can determine high pass composite filter gHLeast square estimation gH
(n)：
Work as NHWhen being even number, gH(n)=[QH(k),-fliplr(QH(k))],
Work as NHWhen being odd number, gH(n)=[QH(k),fliplr(QH(k-1))]。
Step G, according to the impulse response g of high pass composite filterH(n) high pass resolution filter h is obtainedH。
Due to setting the impulse response g of high pass composite filterH(n) with the Least square estimation h of high pass resolution filterH
(n) meet time domain inverse relation：hH(n)=gH(NH- 1-n), therefore can obtain high pass resolution filter by the time domain inverse relation
Least square estimation hH(n), i.e. high pass resolution filter hH。
Such as set gH(n)=[gH(0) gH(1) gH(2) gH(3)],NH=4, gH(n) and hH(n) meet time domain reversion to close
It is, then hH(n)=gH(NH- n-1)=[gH(3) gH(2) gH(1) gH(0)]。
By above step, can obtain be based on the bandwidth of Laplace structureThe high-pass filter that length is 60, frequency spectrum branch
Support is ranging fromThe performance of the high-pass filter is as shown in Figure 6.Wherein：
Fig. 6 (a) is high pass composite filter gHTime-domain pulse response,
Fig. 6 (b) is high pass composite filter gHFrequency domain response,
Fig. 6 (c) is high pass resolution filter hHTime-domain pulse response,
Fig. 6 (d) is high pass resolution filter hHFrequency domain response.
It can be seen that the high-pass filter belongs to linear-phase filter, demonstrates this example design from Fig. 6 (a) and 6 (c)
High pass composite filter Least square estimation gH(n) meet<gH(n),gH(n-mM)>=δ (m), therefore the high pass synthetic filtering
Utensil has orthogonality, and the impulse response g of high pass composite filterH(n) with the Least square estimation h of high pass resolution filterH
(n) meet time domain inverse relation：hH(n)=gH(NH- 1-n), therefore the impulse response h of high pass resolution filterH(n) also meet<
hH(n),hH(n-mM)>=δ (m), i.e. the high pass resolution filter have orthogonality.It can be seen that from Fig. 6 (b) and 6 (d), the high pass
Filter spectrum support range beI.e. bandwidth is
Claims (3)
1. a kind of low pass filter design method of the bandwidth varying linear phase based on Laplace structure, including low pass is synthesized and is filtered
Wave device gLDesign and to low pass resolution filter hLDesign：
(1) setting low pass composite filter gLLength NL, cut-off frequecy of passbandStopband initial frequencyBandwidth is
Wherein M is decimation factor, M and NLIt is integer, M >=2, NLFor the integral multiple of M,
(2) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the lowpass prototype filter p of linear phaseLArteries and veins
Rush response sequence pL(n), n=0,1,2 ... NL-1；
(3) lowpass prototype filter p is determinedLHalf coefficient qL(k)：
Work as NLFor odd number when,
Work as NLFor even number when,
(4) by qL(k) initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gL0(n),gL0(n-mM)>=δ (m)
Wherein, gL0(n) it is low-pass filter impulse response free variable, works as NLFor even number when, gL0(n)=[qL(k),fliplr
(qL(k))], work as NLFor odd number when, gL0(n)=[qL(k),fliplr(qL(k-1)], fliplr is that sequence or so is turned in MatLab
Turn function, GL(ejω) it is gL0(n) frequency response, α are weights, and value is 0 ＜ α ＜ 1, and m is shift count, and value is arbitrary
Integer, δ (m) are dirac sequences,
In the case where meeting constraints, weight α is adjusted, as object function φLIt is low after being optimized when value reaches minimum
Logical composite filter gLHalf coefficient QL(k)；
(5) according to the half coefficient Q of the low pass composite filter after optimizationL(k), low pass composite filter g is obtainedLPulse ring
Answer sequence gL(n)：
Work as NLFor odd number when, gL(n)=[QL(k),fliplr(QL(k-1))]；
Work as NLFor even number when, gL(n)=[QL(k),fliplr(QL(k))]；
(6) according to the Least square estimation g of low pass composite filterL(n) with the Least square estimation h of low pass resolution filterL
(n) meet time domain inverse relation：hL(n)=gL(NL- 1-n), by the Least square estimation g of low pass composite filterL(n)
Find out the Least square estimation h of low pass resolution filterL(n)。
2. a kind of Design of Bandpass method of the bandwidth varying linear phase based on Laplace structure, including band logical is synthesized and is filtered
Wave device gBDesign and to band logical resolution filter hBDesign：
1) setting band logical composite filter gBLength NB, cut-off frequecy of passband is respectivelyStopband initial frequency
RespectivelyDecimation factor is M, and bandwidth isWherein M and NBIt is integer, M >=2, NBIt is the integral multiple of M,
2) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the band logical ptototype filter p of linear phaseBPulse
Response sequence pB(n), n=0,1,2 ... NB-1；
3) band logical ptototype filter p is determinedBHalf coefficient qB(k)：
Work as NBFor odd number when,
Work as NBFor even number when,
4) by qB(k) initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gB0(n),gB0(n-mM)>=δ (m)
Wherein, gB0(n) it is bandpass filter impulse response free variable, works as NBFor even number when, gB0(n)=[qB(k),fliplr
(qB(k))], work as NBFor odd number when, gB0(n)=[qB(k),fliplr(qB(k-1)], fliplr is that sequence or so is turned in MatLab
Turn function, GB(ejω) it is gB0Frequency response, α is weight, and value is 0 ＜ α ＜ 1, and m is shift count, and value is arbitrary whole
Number, δ (m) is dirac sequence,
In the case where meeting constraints, weight α is adjusted, as object function φBValue reaches minimum, the low pass after being optimized
The half coefficient Q of composite filterB(k)；
5) according to the half coefficient Q of the band logical composite filter after optimizationB(k), band logical composite filter g is obtainedBImpulse response
Sequence gB(n)：
Work as NBFor odd number when, gB(n)=[QB(k),fliplr(QB(k-1))]；
Work as NBFor even number when, gB(n)=[QB(k),fliplr(QB(k))]；
6) according to the Least square estimation g of band logical composite filterB(n) with the Least square estimation h of band logical resolution filterB(n)
Meet time domain inverse relation：hB(n)=gB(NB- 1-n), by the Least square estimation g of band logical composite filterB(n) it can find out
The Least square estimation h of band logical resolution filterB(n)。
3. a kind of high-pass filter design method of the bandwidth varying linear phase based on Laplace structure, including high pass is synthesized and is filtered
Wave device gHDesign and to high pass resolution filter hHDesign：
(A) setting high pass composite filter gHLength NH, cut-off frequecy of passband isStopband initial frequency isBandwidth isWherein M is decimation factor, M and NHIt is integer, M >=2, NHFor the integral multiple of M,
(B) according to above-mentioned parameter, the firpm functions in MatLab is called to generate the high pass ptototype filter p of linear phaseHArteries and veins
Rush response sequence pH(n), n=0,1,2 ... NH-1；
(C) high pass ptototype filter p is determinedHHalf coefficient qH(k)：
Work as NHFor odd number when,
Work as NHFor even number when,
(D) by qHThe initial value of function fmincon as an optimization is optimized according to following optimization formula：
s.t.<gH0(n),gH0(n-mM)>=δ (m)
Wherein, gH0(n) it is high-pass filter impulse response free variable：Work as NHFor even number when, gH0(n)=[qH(k),-fliplr
(qH], (k)) work as NHFor odd number when, gH0(n)=[qH(k),fliplr(qH(k-1)], fliplr is that sequence or so is turned in MatLab
Turn function, GH(ejω) it is gH0(n) frequency response, α are weights, and value is 0 ＜ α ＜ 1, and m is shift count, and value is arbitrary
Integer, δ (m) are dirac sequences,
In the case where meeting constraints, weight α is adjusted, as object function φHWhen value reaches minimum, the height after being optimized
The half coefficient Q of logical composite filterH(k)；
(E) according to the half coefficient Q of the high pass composite filter after optimizationH(k), high pass composite filter g is obtainedHPulse ring
Answer sequence gH(n)：
Work as NHFor odd number when, gH(n)=[QH(k),fliplr(QH(k-1))],
Work as NHFor even number when, gH(n)=[QH(k),-fliplr(QH(k))]；
(F) according to the Least square estimation g of high pass composite filterH(n) with the Least square estimation h of high pass resolution filterH
(n) meet time domain inverse relation：hH(n)=gH(NH- 1-n), by the Least square estimation g of high pass composite filterH(n)
Find out the Least square estimation h of high pass resolution filterH(n)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610044267.XA CN105719255B (en) | 2016-01-22 | 2016-01-22 | Bandwidth varying linear-phase filter method based on Laplace structure |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610044267.XA CN105719255B (en) | 2016-01-22 | 2016-01-22 | Bandwidth varying linear-phase filter method based on Laplace structure |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105719255A CN105719255A (en) | 2016-06-29 |
| CN105719255B true CN105719255B (en) | 2018-08-17 |
Family
ID=56153847
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201610044267.XA Active CN105719255B (en) | 2016-01-22 | 2016-01-22 | Bandwidth varying linear-phase filter method based on Laplace structure |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105719255B (en) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP3236626B1 (en) * | 2016-04-21 | 2020-09-23 | Institut Mines Telecom / Telecom Bretagne | Filter for linear modulation based communication systems |
| CN108427032B (en) * | 2018-01-29 | 2020-12-11 | 中国电子科技网络信息安全有限公司 | Frequency spectrum decomposition method and frequency-time inversion method |
| CN111061150B (en) * | 2019-10-23 | 2020-11-27 | 南京大学 | Hardware implementation method of Laplace frequency response |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101546992A (en) * | 2008-03-29 | 2009-09-30 | 华为技术有限公司 | Filtering method and filter |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| EP2417704B1 (en) * | 2009-04-10 | 2017-08-02 | Hittite Microwave LLC | Fractional-n frequency synthesizer having reduced fractional switching noise |
-
2016
- 2016-01-22 CN CN201610044267.XA patent/CN105719255B/en active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101546992A (en) * | 2008-03-29 | 2009-09-30 | 华为技术有限公司 | Filtering method and filter |
Non-Patent Citations (1)
| Title |
|---|
| 两通道LPBFB滤波器对称性与长度研究;张池军等;《信号处理》;20090831;第25卷(第8期);1325-1327 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105719255A (en) | 2016-06-29 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Unser et al. | Cardinal exponential splines: Part I-Theory and filtering algorithms | |
| CN101635563B (en) | Method and device for processing time domain signal of real-value | |
| Bayram et al. | Frequency-domain design of overcomplete rational-dilation wavelet transforms | |
| CN105719255B (en) | Bandwidth varying linear-phase filter method based on Laplace structure | |
| Blu et al. | Wavelets, fractals, and radial basis functions | |
| Bayram et al. | Overcomplete discrete wavelet transforms with rational dilation factors | |
| DE102015116269B4 (en) | SCRATCH TRANSFORMER, ANALOG-DIGITAL TRANSFORMER WITH A SCREEN TRANSFORMER, AND METHOD FOR CONVERTING A DATA FLOW FROM A DATA RATE TO ANOTHER DATA RATE | |
| CN101977033B (en) | Digital filtering method for underground instrument signal transmission | |
| Patil et al. | On the design of FIR wavelet filter banks using factorization of a halfband polynomial | |
| CN107294512B (en) | Non-uniform filter bank filtering method based on tree structure | |
| CN105811920B (en) | A kind of FRM narrow transition band filter group structure | |
| Wu et al. | Design of discrete Fourier transform modulated filter bank with sharp transition band | |
| CN107239623B (en) | Optimal design method of M-channel oversampling image filter bank based on convex optimization | |
| Lin et al. | New perspectives and improvements on the symmetric extension filter bank for subband/wavelet image compression | |
| Salgado et al. | Power and area efficient comb-based decimator for sigma-delta ADCs with high decimation factors | |
| CN105099398B (en) | The construction method of non-homogeneous DFT modulated filters group based on phase-modulation | |
| CN103973254A (en) | Transimpedance type integrated band-pass filter design method | |
| CN114331926B (en) | Two-channel graph filter bank coefficient design optimization method based on element changing idea | |
| Tay et al. | Biorthogonal filter banks constructed from four halfband filters | |
| Lu et al. | A mapping-based design for nonsubsampled hourglass filter banks in arbitrary dimensions | |
| Sreelekha et al. | Development of variable 2D FIR filter structures using farrow approximation and row-wise polyphase decomposition | |
| CN111010146A (en) | Signal reconstruction structure based on fast filter bank and design method thereof | |
| Lai | Popular Wavelet Families and Filters and Their Use. | |
| Benkrid et al. | Handling finite length signals borders in two-channel multirate filter banks for perfect reconstruction | |
| CN2826831Y (en) | Parallel superposed frequency domain digital array filter |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |