frequency respons of halfband filter?

Started by alex65111 April 15, 2009
Let the halfband filter is defined

H=[0.033402 0 -0.015082	0 0.018968 0 -0.024143	0 0.031336 0 -0.042322	0 
0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 -0.042322
0 0.031336 0 -0.024143	0 0.018968 0 -0.015082	0 0.033402];

According to property of the halfband filter it has suppression on 1/4
sampling frequency 6dB.

Fx=fft(H,1024);
plot(20*log10(abs(Fx)))

Let's calculate the frequency response

s1=[zeros(1,128) 1 zeros(1,128)]; 
s_f=filter(H,1,s1); 
s_f=s_f(1:2:end); 
figure, plot(20*log10(abs(fft(s_f))))

On figure it is clearly visible, that on half of new sampling of frequency
instead of 3dB suppression actually it is more than  15dB suppression.


Why at calculation of the frequency response instead of expected 3dB
actually it appears much more?

If to generate a signal on frequency of 1/4 sampling frequencies, to pass
it through the halfband filter and then to decimate

t=0:1023;
s=sin(2*pi*0.24999*t)+0.1*randn(1,length(t));

s_f1=filter(H,1,s); 
s_f1=s_f1(1:2:end); 
figure, plot(20*log10(abs(fft(s_f1))))

that will be visible, that suppression of this frequency will be not such
what should be according to the frequency response. 
Why?


Does the sum of your coefficients equal 1?

-jim

alex65111 wrote:
> > Let the halfband filter is defined > > H=[0.033402 0 -0.015082 0 0.018968 0 -0.024143 0 0.031336 0 -0.042322 0 > 0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 -0.042322 > 0 0.031336 0 -0.024143 0 0.018968 0 -0.015082 0 0.033402]; > > According to property of the halfband filter it has suppression on 1/4 > sampling frequency 6dB. > > Fx=fft(H,1024); > plot(20*log10(abs(Fx))) > > Let's calculate the frequency response > > s1=[zeros(1,128) 1 zeros(1,128)]; > s_f=filter(H,1,s1); > s_f=s_f(1:2:end); > figure, plot(20*log10(abs(fft(s_f)))) > > On figure it is clearly visible, that on half of new sampling of frequency > instead of 3dB suppression actually it is more than 15dB suppression. > > Why at calculation of the frequency response instead of expected 3dB > actually it appears much more? > > If to generate a signal on frequency of 1/4 sampling frequencies, to pass > it through the halfband filter and then to decimate > > t=0:1023; > s=sin(2*pi*0.24999*t)+0.1*randn(1,length(t)); > > s_f1=filter(H,1,s); > s_f1=s_f1(1:2:end); > figure, plot(20*log10(abs(fft(s_f1)))) > > that will be visible, that suppression of this frequency will be not such > what should be according to the frequency response. > Why?
>Does the sum of your coefficients equal 1? > >-jim >
It unless concerns that that at calculation of the frequency response on frequency on 1/2 new sampling frequency, suppression appears distinct from 3dB? P.S. Coefficients have been calculated in fdatool (matlab).

alex65111 wrote:
> > >Does the sum of your coefficients equal 1? > > > >-jim > > > > It unless concerns that that at calculation of the frequency response on > frequency on 1/2 new sampling frequency, suppression appears distinct from > 3dB?
Sorry I don't speak that language. The question I asked could be answered with yes or no. -jim
> > P.S. Coefficients have been calculated in fdatool (matlab).
On Apr 15, 5:31&#2013266080;pm, "alex65111" <alex65...@list.ru>
wrote:
> Let the halfband filter is defined > > H=[0.033402 0 -0.015082 0 0.018968 0 -0.024143 &#2013266080;0 0.031336 0
-0.042322 &#2013266080;0
> 0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 -0.042322 > 0 0.031336 0 -0.024143 &#2013266080;0 0.018968 0 -0.015082
&#2013266080;0 0.033402];
> > According to property of the halfband filter it has suppression on 1/4 > sampling frequency 6dB. > > Fx=fft(H,1024); > plot(20*log10(abs(Fx))) > > Let's calculate the frequency response > > s1=[zeros(1,128) 1 zeros(1,128)]; > s_f=filter(H,1,s1); > s_f=s_f(1:2:end); > figure, plot(20*log10(abs(fft(s_f)))) > > On figure it is clearly visible, that on half of new sampling of frequency > instead of 3dB suppression actually it is more than &#2013266080;15dB
suppression. Where did you get this idea? If you decimate the filter response to kill the DC offset you get zero.
> > Why at calculation of the frequency response instead of expected 3dB > actually it appears much more? > > If to generate a signal on frequency of 1/4 sampling frequencies, to pass > it through the halfband filter and then to decimate > > t=0:1023; > s=sin(2*pi*0.24999*t)+0.1*randn(1,length(t)); > > s_f1=filter(H,1,s); > s_f1=s_f1(1:2:end); > figure, plot(20*log10(abs(fft(s_f1)))) > > that will be visible, that suppression of this frequency will be not such > what should be according to the frequency response. > Why?
================ TRY THIS ... AND THINK!!!! ============== clear all; close all; eps=10^-7; H1=[0.033402 0 -0.015082 0 0.018968 0 -0.024143 0 0.031336 0 -0.042322 0 ... 0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 -0.042322 ... 0 0.031336 0 -0.024143 0 0.018968 0 -0.015082 0 0.033402]; H2=H1; H2((length(H2)+1)/2)=0; % reset 0.5 to 0 %According to property of the halfband filter it has suppression on 1/4 %sampling frequency 6dB. Fx1=fft(H1,1024); f=(0:1023)/1024; figure(1); plot(f,20*log10(abs(Fx1)+eps)); grid on; title('Plot 1'); Fx2=fft(H2,1024); f=(0:1023)/1024; figure(2); plot(f,20*log10(abs(Fx2)+eps)); grid on; title('Plot 2 - makes sense'); %Let's calculate the frequency response of decimated filter response which % eliminates zeroes and the 0.5 offset. Value at new f/4 should be sum of % all coefficients except 0.5 . s1=[zeros(1,128) 1 zeros(1,128)]; s_f1=filter(H1,1,s1); s_f1=s_f1(1:2:end); fd=(0:length(s_f1)-1)/length(s_f1); figure(3); plot(fd, 20*log10(abs(fft(s_f1))+eps)); grid on; title('Plot 3 - makes sense'); ====================================== SplatSpam
Sorry, but I not understand your thought. Why decimated filter response 
eliminates  the 0.5 offset?

If I do so

Fx3_1=abs(Fx1(1:512));
Fx3_2=abs(Fx1(513:1024)); 
plot(20*log10(Fx3_1+eps))
hold on
plot(20*log10(Fx3_2+eps),'r')
hold off

that I receive that in a point 0.25 (0.5 after decimate) suppression
corresponds 6dB.

But if I try to construct the characteristic of system by giving on an
input of system of an impulse ([0 0 1 0 0]->filter->decimate->impulse
response->fft), that I receive that in a point 0.25 (0.5) suppression
corresponds >>6dB (approximately 20dB).

And if on a system input to give a signal on frequency 0.25 on an output
suppression there will be the same 6dB, but not 20dB.

Thus it turns out that at giving [0 1 0] on frequency 0.25(0.5 after
decimate)  suppression appears 20dB, and at giving sin(0.25) suppression
appears 6dB.

alex65111 wrote:
> Let the halfband filter is defined > > H=[0.033402 0 -0.015082 0 0.018968 0 -0.024143 0 0.031336 0 -0.042322 > 0 > 0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 > -0.042322 0 0.031336 0 -0.024143 0 0.018968 0 -0.015082 0 0.033402]; > > According to property of the halfband filter it has suppression on 1/4 > sampling frequency 6dB. > > Fx=fft(H,1024); > plot(20*log10(abs(Fx))) > > Let's calculate the frequency response > > s1=[zeros(1,128) 1 zeros(1,128)]; > s_f=filter(H,1,s1); > s_f=s_f(1:2:end); > figure, plot(20*log10(abs(fft(s_f)))) > > On figure it is clearly visible, that on half of new sampling of > frequency instead of 3dB suppression actually it is more than 15dB > suppression. > > > Why at calculation of the frequency response instead of expected 3dB > actually it appears much more? > > If to generate a signal on frequency of 1/4 sampling frequencies, to > pass it through the halfband filter and then to decimate > > t=0:1023; > s=sin(2*pi*0.24999*t)+0.1*randn(1,length(t)); > > s_f1=filter(H,1,s); > s_f1=s_f1(1:2:end); > figure, plot(20*log10(abs(fft(s_f1)))) > > that will be visible, that suppression of this frequency will be not > such what should be according to the frequency response. > Why?
Check your code. An FFT of the filter yields the correct response it appears. Fred
On Apr 16, 4:32&#2013266080;am, "alex65111" <alex65...@list.ru>
wrote:
> Sorry, but I not understand your thought. Why decimated filter response > eliminates &#2013266080;the 0.5 offset?
Because it eliminates it from the impulse response, you throw it away. The second fequency response you are trying to compute is not real meaningful because you are not considering that the combined frequency response and signal spectrum above the original fs/2 BW aliases into the new fs/2 BW (old fs/4 BW). You need to draw some pictures. Dirk
> > If I do so > > Fx3_1=abs(Fx1(1:512)); > Fx3_2=abs(Fx1(513:1024)); > plot(20*log10(Fx3_1+eps)) > hold on > plot(20*log10(Fx3_2+eps),'r') > hold off > > that I receive that in a point 0.25 (0.5 after decimate) suppression > corresponds 6dB. > > But if I try to construct the characteristic of system by giving on an > input of system of an impulse ([0 0 1 0 0]->filter->decimate->impulse > response->fft), that I receive that in a point 0.25 (0.5) suppression > corresponds >>6dB (approximately 20dB). > > And if on a system input to give a signal on frequency 0.25 on an output > suppression there will be the same 6dB, but not 20dB. > > Thus it turns out that at giving [0 1 0] on frequency 0.25(0.5 after > decimate) &#2013266080;suppression appears 20dB, and at giving sin(0.25)
suppression
> appears 6dB.
All thanks for attempts to help me but I have not understood where I am
mistaken.

Let's try to ask a question easier and particularly.

Let's assume at us there is a system consisting of the filter and
following for it decimation.

The impulse response of filter

H=[0.033402 0 -0.015082	0 0.018968 0 -0.024143	0 0.031336 0 -0.042322	0 
0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 
-0.042322 0 0.031336 0 -0.024143	0 0.018968 0 -0.015082	0 0.033402];

Decimation factor is 2.

How correctly to calculate the combined characteristic of system as a
whole? 

The calculated characteristic should show what suppression on frequency
0.5 (after decimation)?
On Apr 17, 12:58&#2013266080;am, "alex65111" <alex65...@list.ru>
wrote:
> All thanks for attempts to help me but I have not understood where I am > mistaken. > > Let's try to ask a question easier and particularly. > > Let's assume at us there is a system consisting of the filter and > following for it decimation. > > The impulse response of filter > > H=[0.033402 0 -0.015082 0 0.018968 0 -0.024143 &#2013266080;0 0.031336 0
-0.042322 &#2013266080;0
> 0.061361 0 -0.10468 0 0.31784 0.5 0.31784 0 -0.10468 0 0.061361 0 > -0.042322 0 0.031336 0 -0.024143 &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0 0.018968 0 -0.015082 &#2013266080;0 0.033402];
> > Decimation factor is 2. > > How correctly to calculate the combined characteristic of system as a > whole? > > The calculated characteristic should show what suppression on frequency > 0.5 (after decimation)?
Alex, You are not going to be able to get a frequency response of this system because of the aliasing below 0.5. In the aliased region the frequency response would not be defined in any meaningful way. You need to understand what the filter does. Draw pictures. Dirk
On Apr 30, 3:43&#2013266080;pm, "Fred Marshall"
<fmarshallx@remove_the_x.acm.org>
wrote:
> bellda2...@cox.net wrote: > > Hi Fred, > > > I meant the a half-band with f(0)=1.0, if the code to generate on is > > readily available. > > > Thanks, > > > Dirk > > Here is a quickie but the transition bands have bumps in them which is > unecessary. &#2013266080;This means it was an early version of my code.
&#2013266080;I don't know if
> the 1.0 at f=0 constraint is built into a later version. I can keep looking. > Anyway, you can look at the passband to see how it looks. &#2013266080;You
can also look
> in the stopband to see that there's a double zero at fs/2. > > This works out because: > > A halfband filter response is made up of a constant = 0.5 and an > antisymmetric component that is 0.5 in the passband and -0.5 in the > stopband. &#2013266080;When added together, the components yield 1.0 in the
passband and
> zero in the stopband (+/- the ripples of course). &#2013266080;If the
frequency response
> has a constant term, in this case 0.5, then the temporal response must be > 0.5 at "DC" or, really, at time zero. > So, we can design the antisymmetric filter and then just add the center > coefficient of 0.5. &#2013266080;That's a way to look at it. > > The approximation only needs to be done in the passband because the > structure guarantees that fs/4=0.5 (or zero) and guarantees the response > from fs/4 to fs/2 is antisymmetric. > > Thinking of it this way, if the frequency response is 1.0 at f=0 and the > components are as described above, then the frequency response at fs/2 must > be zero. &#2013266080;Indeed, you see a double zero there. > > Anyway, here is that imperfect (if not interesting) filter set of coeff's: > &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.318178627640009E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.624465197324753E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.656183436512947E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.963471177965403E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.100250272080302E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.103756394237280E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.913409423083067E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.710289599373937E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.520208990201354E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.379590759985149E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.529637048020959E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.101068196818233E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.219537038356066E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.439153611660004E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.934268832206726E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.313605815172195E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.500000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.313605815172195E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.934268832206726E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.439153611660004E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.219537038356066E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.101068196818233E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.529637048020959E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.379590759985149E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.520208990201354E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.710289599373937E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.913409423083067E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.103756394237280E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.100250272080302E-01
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.963471177965403E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.656183436512947E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.624465197324753E-02
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; &#2013266080;0.000000000000000E+00
> &#2013266080; &#2013266080; &#2013266080; &#2013266080;
&#2013266080; -0.318178627640009E-02
> > I'll see if I can't find better code than what created this example. > > Fred
Thanks Fred, It does just what you advertised (not that I doubted you, I just wanted to see the response). Very interesting, I will have to reread your method. Thanks again, Dirk
bellda2005@cox.net wrote:
> Hi Fred, > > I meant the a half-band with f(0)=1.0, if the code to generate on is > readily available. > > Thanks, > > Dirk
Here is a quickie but the transition bands have bumps in them which is unecessary. This means it was an early version of my code. I don't know if the 1.0 at f=0 constraint is built into a later version. I can keep looking. Anyway, you can look at the passband to see how it looks. You can also look in the stopband to see that there's a double zero at fs/2. This works out because: A halfband filter response is made up of a constant = 0.5 and an antisymmetric component that is 0.5 in the passband and -0.5 in the stopband. When added together, the components yield 1.0 in the passband and zero in the stopband (+/- the ripples of course). If the frequency response has a constant term, in this case 0.5, then the temporal response must be 0.5 at "DC" or, really, at time zero. So, we can design the antisymmetric filter and then just add the center coefficient of 0.5. That's a way to look at it. The approximation only needs to be done in the passband because the structure guarantees that fs/4=0.5 (or zero) and guarantees the response from fs/4 to fs/2 is antisymmetric. Thinking of it this way, if the frequency response is 1.0 at f=0 and the components are as described above, then the frequency response at fs/2 must be zero. Indeed, you see a double zero there. Anyway, here is that imperfect (if not interesting) filter set of coeff's: -0.318178627640009E-02 0.000000000000000E+00 0.624465197324753E-02 0.000000000000000E+00 -0.656183436512947E-02 0.000000000000000E+00 0.963471177965403E-02 0.000000000000000E+00 -0.100250272080302E-01 0.000000000000000E+00 0.103756394237280E-01 0.000000000000000E+00 -0.913409423083067E-02 0.000000000000000E+00 0.710289599373937E-02 0.000000000000000E+00 -0.520208990201354E-02 0.000000000000000E+00 0.379590759985149E-02 0.000000000000000E+00 -0.529637048020959E-02 0.000000000000000E+00 0.101068196818233E-01 0.000000000000000E+00 -0.219537038356066E-01 0.000000000000000E+00 0.439153611660004E-01 0.000000000000000E+00 -0.934268832206726E-01 0.000000000000000E+00 0.313605815172195E+00 0.500000000000000E+00 0.313605815172195E+00 0.000000000000000E+00 -0.934268832206726E-01 0.000000000000000E+00 0.439153611660004E-01 0.000000000000000E+00 -0.219537038356066E-01 0.000000000000000E+00 0.101068196818233E-01 0.000000000000000E+00 -0.529637048020959E-02 0.000000000000000E+00 0.379590759985149E-02 0.000000000000000E+00 -0.520208990201354E-02 0.000000000000000E+00 0.710289599373937E-02 0.000000000000000E+00 -0.913409423083067E-02 0.000000000000000E+00 0.103756394237280E-01 0.000000000000000E+00 -0.100250272080302E-01 0.000000000000000E+00 0.963471177965403E-02 0.000000000000000E+00 -0.656183436512947E-02 0.000000000000000E+00 0.624465197324753E-02 0.000000000000000E+00 -0.318178627640009E-02 I'll see if I can't find better code than what created this example. Fred
On Apr 30, 3:03&#2013266080;am, "Fred Marshall"
<fmarshallx@remove_the_x.acm.org>
wrote:
> bellda2...@cox.net wrote: > > > Interesting, Fred. &#2013266080;Could you post the coefficients of an
example
> > filter? > > > Thanks, > > > Dirk > > Dirk, > > Actually, it really does depend on which type so if you're a bit more > specific then I could see if I still have the code that will generate it. > The Hermann-Schussler work was done quite long ago. > > Do you mean a half-band with f(0)=1.0 ??? or something else? > > Fred
Hi Fred, I meant the a half-band with f(0)=1.0, if the code to generate on is readily available. Thanks, Dirk
Fred Marshall wrote:
> bellda2005@cox.net wrote: >> >> Interesting, Fred. Could you post the coefficients of an example >> filter? >> >> Thanks, >> >> Dirk > > Dirk, > > Actually, it really does depend on which type so if you're a bit more > specific then I could see if I still have the code that will generate > it. The Hermann-Schussler work was done quite long ago. > > Do you mean a half-band with f(0)=1.0 ??? or something else? > > Fred
In the mean time, here is some insight into equality constraints: In the standard Remez algorithm (and this would be true of any approach reaching minimax criteria), the ripple peaks alternate in sign and there are generally N+1 such peaks where N is the number of variables (e.g. filter coefficients). On rare occasions there are N+2. If you impose an equality constraint (such as 1.0 at f=0), then you might visualize this as "pushing" the function toward that value. And, you might think of this point as having infinite weight re: the error .. thus forcing the error to zero. So, if the weight is infinite at that point it's as if there's a peak in the error (even though it's very small in magnitude due to the weighting - or even zero .. which is fine). The function "pushes back" in some direction. As it turns out then, each such equality point takes up the place of an error peak and, although you can't readily "see" it, the sign of this error is as any of the others - that is, it's in the set of peaks whose signs alternate. Accordingly, since there's no perceptible error at this point, it appears that the adjacent error peaks "skip" an alternation - the adjacent error peaks have the same sign and the point of equality "uses up" the "missing" alternating peak. All this leads to a modified Remez algorithm .. with proofs, etc. which takes into consideration the "modified sign alternation property". Since the equality constraint must be represented in the set of equations that are being solved at each step, and because this equation has no error value to be equalized along with the others (it has zero error value), it uses up one degree of freedom. And, because of this, the equalized error (the ripple peaks) go up relative to not having imposed the constraint. From ripple magnitude perspective, it's like having a FIR filter with one less coefficient. Fred
bellda2005@cox.net wrote:
> > Interesting, Fred. Could you post the coefficients of an example > filter? > > Thanks, > > Dirk
Dirk, Actually, it really does depend on which type so if you're a bit more specific then I could see if I still have the code that will generate it. The Hermann-Schussler work was done quite long ago. Do you mean a half-band with f(0)=1.0 ??? or something else? Fred
On Apr 23, 11:44&#2013266080;pm, "Fred Marshall"
<fmarshallx@remove_the_x.acm.org>
wrote:
> bellda2...@cox.net wrote: > > On Apr 22, 10:53 pm, dbd <d...@ieee.org> wrote: > >> On Apr 22, 3:20 pm, bellda2...@cox.net wrote: > > >>> ... > > >>> Just curious, why do want to constrain the gain at DC to be
exactly
> >>> 1.0? > > .................................... > > > > > I was actually interested in Fred's response, since he brought it up; > > I thought he might have a different reason. > > > Thanks for your input. > > > Dirk > > Oh, I didn't mean to suggest there was (or wasn't) a good reason for doing > it, just that it could be done. &#2013266080;If you do a minimax design and
constrain
> the DC value to be 1.0 while the mid-point between the ripple peaks is also > 1.0 then you use up one degree of freedom and the ripple magnitude is larger > than it needs to be otherwise. > > In a typical halfband filter this also forces the response at fs/2 to be > zero. > > One of my "things" is equality-constrained minimax design so I find
things
> like this interesting. &#2013266080;And, I had done exactly this for a
while in my
> halfband design program. > > In general, the DC value (well all the band edges) will be a ripple peak. > And, you can force that value to be 1.0 by scaling. &#2013266080;Here's a
trick I've
> used in minimax design to get all-positive stopband ripple leading up to > Hermann-Schussler equiripple minimum phase FIR designs: > > Start with the stopband desired gain at zero. > At the next iteration of the Remez algorithm, increase the desired stopband > gain to be the same as the equalized error (which is not yet the peak > ripple). &#2013266080;The equalized error grows at each iteration so you
continue to do
> this until you reach convergence. &#2013266080;(The error peaks decrease
with each
> iteration at the same time). &#2013266080;When it's converged, the stopband
value will
> be equal to the ripple peak magnitude and the negative going peaks will be > zero-valued. > > [I've done this successfully and worked on a convergence proof but didn't > quite finish it. &#2013266080;It seems a useful extension to the Remez
algorithm and
> could stand having that proof. &#2013266080;Good grad student exercise I
should think.]
> > You could do the same in the passband but it's not quite as cool: > Start with the passband desired value at 1.0. > Then, in the next iteration, make the passband desired value equal to 1.0 > minus (or plus) the equalized error. > If the DC value is going to be at a positive peak then it will end up at > 1.0. > If the DC value is going to be at a negative peak then it will end up at 1.0 > minus 2x the peak ripple. > You could always run it twice: one with 1.0 minus... and one with 1.0 > plus... > Of course, all this does is scale the result but there's no loss of degrees > of freedom. &#2013266080;You might just as well divide all the coefficients
by (1.0 +
> peak) after the filter is designed to get the DC value to be 1.0 AND not > lose degrees of freedom. > But then, the ripple is all one-sided and deviates 2x if a passband of 1.0 > is what you wanted. &#2013266080;You may well be better off losing one
degree of freedom
> and having slightly larger peaks but not 2x peaks!! > > Fred
Interesting, Fred. Could you post the coefficients of an example filter? Thanks, Dirk
bellda2005@cox.net wrote:
> On Apr 22, 10:53 pm, dbd <d...@ieee.org> wrote: >> On Apr 22, 3:20 pm, bellda2...@cox.net wrote: >> >>> ... >> >>> Just curious, why do want to constrain the gain at DC to be exactly >>> 1.0?
....................................
> > I was actually interested in Fred's response, since he brought it up; > I thought he might have a different reason. > > Thanks for your input. > > Dirk
Oh, I didn't mean to suggest there was (or wasn't) a good reason for doing it, just that it could be done. If you do a minimax design and constrain the DC value to be 1.0 while the mid-point between the ripple peaks is also 1.0 then you use up one degree of freedom and the ripple magnitude is larger than it needs to be otherwise. In a typical halfband filter this also forces the response at fs/2 to be zero. One of my "things" is equality-constrained minimax design so I find things like this interesting. And, I had done exactly this for a while in my halfband design program. In general, the DC value (well all the band edges) will be a ripple peak. And, you can force that value to be 1.0 by scaling. Here's a trick I've used in minimax design to get all-positive stopband ripple leading up to Hermann-Schussler equiripple minimum phase FIR designs: Start with the stopband desired gain at zero. At the next iteration of the Remez algorithm, increase the desired stopband gain to be the same as the equalized error (which is not yet the peak ripple). The equalized error grows at each iteration so you continue to do this until you reach convergence. (The error peaks decrease with each iteration at the same time). When it's converged, the stopband value will be equal to the ripple peak magnitude and the negative going peaks will be zero-valued. [I've done this successfully and worked on a convergence proof but didn't quite finish it. It seems a useful extension to the Remez algorithm and could stand having that proof. Good grad student exercise I should think.] You could do the same in the passband but it's not quite as cool: Start with the passband desired value at 1.0. Then, in the next iteration, make the passband desired value equal to 1.0 minus (or plus) the equalized error. If the DC value is going to be at a positive peak then it will end up at 1.0. If the DC value is going to be at a negative peak then it will end up at 1.0 minus 2x the peak ripple. You could always run it twice: one with 1.0 minus... and one with 1.0 plus... Of course, all this does is scale the result but there's no loss of degrees of freedom. You might just as well divide all the coefficients by (1.0 + peak) after the filter is designed to get the DC value to be 1.0 AND not lose degrees of freedom. But then, the ripple is all one-sided and deviates 2x if a passband of 1.0 is what you wanted. You may well be better off losing one degree of freedom and having slightly larger peaks but not 2x peaks!! Fred
On Apr 22, 10:53&#2013266080;pm, dbd <d...@ieee.org> wrote:
> On Apr 22, 3:20 pm, bellda2...@cox.net wrote: > > > ... > > > Just curious, why do want to constrain the gain at DC to be exactly > > 1.0? > > > Dirk > > ... > > In anti-aliasing filters that are going to be used in multiple passes > of filter/desample it is desirable to use a reduced or zero ripple at > DC design so that error doesn't accumulate with each application of > the filter. That's why I've done it. > > Dale B. Dalrymple
Hi Dale, I have used halfband filters for longer than I care to say; almost every 'permanent'/contract/consulting position I have ever had has had applications for them for rate changing or multi-rate applications. Making DC gain EXACTLY 1.0 is not a problem I have ever encountered in my designs, I am usually worried about the entire passband ripple and the corresponding stopband ripple. If you use the middle coefficient as 0.5 or 1.0 as is commonly mentioned in texts (it used to be important to save the multiply) you don't get exactly 1.0 at DC. For interpolation with a series of halfband filters if you want a gain of ~1.0 from each interpolating filter, the filters won't have a nominal passband gain of 1.0 (of course you could use a passband gain of ~1.0 and adjust the scaling elsewhere to get the desired end result). I was actually interested in Fred's response, since he brought it up; I thought he might have a different reason. Thanks for your input. Dirk
On Apr 22, 3:20 pm, bellda2...@cox.net wrote:
> ... > > Just curious, why do want to constrain the gain at DC to be exactly > 1.0? > > Dirk > ...
In anti-aliasing filters that are going to be used in multiple passes of filter/desample it is desirable to use a reduced or zero ripple at DC design so that error doesn't accumulate with each application of the filter. That's why I've done it. Dale B. Dalrymple
On Apr 22, 11:05&#2013266080;am, "Fred Marshall"
<fmarshallx@remove_the_x.acm.org>
wrote:
> bellda2...@cox.net wrote: > > If you shift the halfband filter to remove the delay and look at the > > frequency response (real and signed, linear scale) you find that the > > response is symmetric about the point at fs/4. &#2013266080;So the
ripple at fs/2
> > will be the same magnitude but opposite sign of the ripple at DC. > > > Dirk > > I agree but don't quite get the point.... > And, this is only correct if you don't constrain the gain at DC to be > exactly 1.0 - which can be done.
Fred, The ripples around the center of the ripples will still be the same for real and signed frequency response on a linear scale. The symmetry displayed this way around the point at fs/4 is not dependent on the filter scaling or what you set the center tap to. However, if you scale the filter and set the center tap to something other than half the distance between the two ripples' centers (like other than 0.5 for ripple center of gains 1 and 0), then the MAGNITUDE response will not be symmetric about the fs/4 point, because you will have more ripple in the stopband than is necessary (i.e. could be reduced at no cost, so why wouldn't you?). If the design was an equiripple design, setting the center tap in this way would make alternating stopband peaks have different magnitudes, the larger stopband peaks larger than they need to be. Just curious, why do want to constrain the gain at DC to be exactly 1.0? Dirk
> > If you shift the halfband filter unit sample response to be centered at zero > delay it doesn't change the frequency magnitude response - except the > function becomes purely real. &#2013266080;The magnitude response in any
case is
> symmetrical around zero and fs/2. > > I have seen cases where the purely real response flips sign at fs/2 relative > to zero but don't recall the details. &#2013266080;Probably an x4 case
because of what I
> just said above. > > Fred