Filter design with 'freqsamp' method?

  • Follow


I'm new to Matlab, so this is probably a really dumb question.
I'm trying to design an FIR where I can set the points.
I want a symmetric FIR.
I'm assuming that an N tap FIR gives me N/2 degrees of freedom
for an even number of taps, and (N + 1)/2 degrees of freedom
for an odd number of taps.
This should allow me to set the frequency response at a
number of points equal to the degrees of freedom, and
have the frequency response travel through these points
exactly. If I don't specify transitions that are too sharp
relative to the number of taps, then the frequency response
should be reasonably well controlled between the set points.

I had assumed that I could use the filter toolbox's freqsamp
method to do this, but it doesn't seem to work as I expected.
What am I doing wrong? I'd normally do this in Mathematica
or C, but I'm trying to get into the spirit of Matlab.

I would want to solve the point setting problem anyway, but
perhaps there's an easier way to do the immediate task that
I wanted it for. I want to implement a Butterworth amplitude
response in an FIR filter (to emulate some legacy equipment)
but I can't find an easy way of doing that in Matlab (hence the
point setting experiments).

Thanks

Pete 


0
Reply Pete 7/1/2010 1:21:26 AM


Pete Fraser wrote:

> I'm new to Matlab, so this is probably a really dumb question.
> I'm trying to design an FIR where I can set the points.
> I want a symmetric FIR.
> I'm assuming that an N tap FIR gives me N/2 degrees of freedom
> for an even number of taps, and (N + 1)/2 degrees of freedom
> for an odd number of taps.
> This should allow me to set the frequency response at a
> number of points equal to the degrees of freedom, and
> have the frequency response travel through these points
> exactly. If I don't specify transitions that are too sharp
> relative to the number of taps, then the frequency response
> should be reasonably well controlled between the set points.
> 
> I had assumed that I could use the filter toolbox's freqsamp
> method to do this, but it doesn't seem to work as I expected.
> What am I doing wrong? I'd normally do this in Mathematica
> or C, but I'm trying to get into the spirit of Matlab.
> 
> I would want to solve the point setting problem anyway, but
> perhaps there's an easier way to do the immediate task that
> I wanted it for. I want to implement a Butterworth amplitude
> response in an FIR filter (to emulate some legacy equipment)
> but I can't find an easy way of doing that in Matlab (hence the
> point setting experiments).

But why would anyone want to emulate Butterworth magnitude response in a 
linear phase FIR filter whereas it is very straightforward to do the 
(**exact) Butterworth as IIR filter?

** frequency warping



Vladimir Vassilevsky
DSP and Mixed Signal Design Consultant
http://www.abvolt.com


0
Reply Vladimir 7/1/2010 1:35:09 AM


"Vladimir Vassilevsky" <nospam@nowhere.com> wrote in message 
news:0Jidnc3knONVbbbRnZ2dnUVZ_oSdnZ2d@giganews.com...

> But why would anyone want to emulate Butterworth magnitude response in a 
> linear phase FIR filter whereas it is very straightforward to do the 
> (**exact) Butterworth as IIR filter?

Butterworth will be one of several possible responses
from the same hardware.

I'd rather not deal with Butterworth's phase issues.

Pete 


0
Reply Pete 7/1/2010 2:06:54 AM

On 6/30/2010 10:06 PM, Pete Fraser wrote:
> "Vladimir Vassilevsky"<nospam@nowhere.com>  wrote in message
> news:0Jidnc3knONVbbbRnZ2dnUVZ_oSdnZ2d@giganews.com...
>
>> But why would anyone want to emulate Butterworth magnitude response in a
>> linear phase FIR filter whereas it is very straightforward to do the
>> (**exact) Butterworth as IIR filter?
>
> Butterworth will be one of several possible responses
> from the same hardware.
>
> I'd rather not deal with Butterworth's phase issues.

What is so good about a Butterworth amplitude response that you want to 
imitate it?

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply Jerry 7/1/2010 2:39:29 AM

On 06/30/2010 06:21 PM, Pete Fraser wrote:
> I'm new to Matlab, so this is probably a really dumb question.
> I'm trying to design an FIR where I can set the points.
> I want a symmetric FIR.
> I'm assuming that an N tap FIR gives me N/2 degrees of freedom
> for an even number of taps, and (N + 1)/2 degrees of freedom
> for an odd number of taps.
> This should allow me to set the frequency response at a
> number of points equal to the degrees of freedom, and
> have the frequency response travel through these points
> exactly. If I don't specify transitions that are too sharp
> relative to the number of taps, then the frequency response
> should be reasonably well controlled between the set points.
>
> I had assumed that I could use the filter toolbox's freqsamp
> method to do this, but it doesn't seem to work as I expected.
> What am I doing wrong? I'd normally do this in Mathematica
> or C, but I'm trying to get into the spirit of Matlab.
>
> I would want to solve the point setting problem anyway, but
> perhaps there's an easier way to do the immediate task that
> I wanted it for. I want to implement a Butterworth amplitude
> response in an FIR filter (to emulate some legacy equipment)
> but I can't find an easy way of doing that in Matlab (hence the
> point setting experiments).

Make a vector that "draws" the frequency response that you want, then 
take its FFT.

It may not be the entirely kosher way to do it, but it certainly works.

-- 

Tim Wescott
Wescott Design Services
http://www.wescottdesign.com

Do you need to implement control loops in software?
"Applied Control Theory for Embedded Systems" was written for you.
See details at http://www.wescottdesign.com/actfes/actfes.html
0
Reply Tim 7/1/2010 2:39:33 AM

Tim Wescott  <tim@seemywebsite.com> wrote:

>Make a vector that "draws" the frequency response that you want, then 
>take its FFT.
>
>It may not be the entirely kosher way to do it, but it certainly works.

Doesn't this give you a complex filter?

Can you just take the real part afterwards?

Steve
0
Reply spope33 7/1/2010 3:10:52 AM

On Jun 30, 10:39=A0pm, Tim Wescott <t...@seemywebsite.com> wrote:
> On 06/30/2010 06:21 PM, Pete Fraser wrote:
>
>
>
> > I'm new to Matlab, so this is probably a really dumb question.
> > I'm trying to design an FIR where I can set the points.
> > I want a symmetric FIR.
> > I'm assuming that an N tap FIR gives me N/2 degrees of freedom
> > for an even number of taps, and (N + 1)/2 degrees of freedom
> > for an odd number of taps.
> > This should allow me to set the frequency response at a
> > number of points equal to the degrees of freedom, and
> > have the frequency response travel through these points
> > exactly. If I don't specify transitions that are too sharp
> > relative to the number of taps, then the frequency response
> > should be reasonably well controlled between the set points.
>
> > I had assumed that I could use the filter toolbox's freqsamp
> > method to do this, but it doesn't seem to work as I expected.
> > What am I doing wrong? I'd normally do this in Mathematica
> > or C, but I'm trying to get into the spirit of Matlab.
>
> > I would want to solve the point setting problem anyway, but
> > perhaps there's an easier way to do the immediate task that
> > I wanted it for. I want to implement a Butterworth amplitude
> > response in an FIR filter (to emulate some legacy equipment)
> > but I can't find an easy way of doing that in Matlab (hence the
> > point setting experiments).
>
> Make a vector that "draws" the frequency response that you want, then
> take its FFT.
>
> It may not be the entirely kosher way to do it, but it certainly works.

Pete might want to use the "window method".  with an iFFT *larger*
than N, draw in the frequency response, he'll get an impulse response
that is longer than N, then window it to length N.  (then maybe FFT
back to see how bad the window mangled things.)

r b-j

0
Reply robert 7/1/2010 3:11:40 AM

On 06/30/2010 08:10 PM, Steve Pope wrote:
> Tim Wescott<tim@seemywebsite.com>  wrote:
>
>> Make a vector that "draws" the frequency response that you want, then
>> take its FFT.
>>
>> It may not be the entirely kosher way to do it, but it certainly works.
>
> Doesn't this give you a complex filter?
>
> Can you just take the real part afterwards?

If the desired magnitude response is real and symmetric around the 
frequency = 0 bin (element 1 in Matlab's indexing scheme) then the FFT 
will be largely real, with only leftover numerical noise in the 
imaginary part.  If the imaginary part has any weight to it at all, then 
you've done something wrong.

Note that you want to be symmetric around frequency = 0 modulo the whole 
thing; i.e. for a 256 element vector element 2 should equal element 256 
with Matlab-style indexing, element 3 should equal element 255, etc.

-- 

Tim Wescott
Wescott Design Services
http://www.wescottdesign.com

Do you need to implement control loops in software?
"Applied Control Theory for Embedded Systems" was written for you.
See details at http://www.wescottdesign.com/actfes/actfes.html
0
Reply tim177 (4404) 7/1/2010 3:34:42 AM

Tim Wescott  <tim@seemywebsite.com> wrote:

>On 06/30/2010 08:10 PM, Steve Pope wrote:

>> Tim Wescott<tim@seemywebsite.com>  wrote:

>>> Make a vector that "draws" the frequency response that you want, then
>>> take its FFT.

>>> It may not be the entirely kosher way to do it, but it certainly works.

>> Doesn't this give you a complex filter?

>> Can you just take the real part afterwards?

>If the desired magnitude response is real and symmetric around the 
>frequency = 0 bin (element 1 in Matlab's indexing scheme) then the FFT 
>will be largely real, with only leftover numerical noise in the 
>imaginary part.  If the imaginary part has any weight to it at all, then 
>you've done something wrong.

That makes sense.  Thanks.

Steve
0
Reply spope33 7/1/2010 3:54:01 AM

"Jerry Avins" <jya@ieee.org> wrote in message 
news:BdTWn.3387$OU6.1391@newsfe20.iad...

> What is so good about a Butterworth amplitude response that you want to 
> imitate it?

Nothing.
It's just that, under some circumstances, equipment that is
later in the signal chain expects a Butterworth response. 


0
Reply Pete 7/1/2010 4:17:15 AM

"Tim Wescott" <tim@seemywebsite.com> wrote in message 
news:26Odnbyf_5N1YrbRnZ2dnUVZ_gadnZ2d@web-ster.com...

> Make a vector that "draws" the frequency response that you want, then take 
> its FFT.

That's typically how I have done it in the past (well, DFT anyway).
It's just that I'm new to Matlab, trying to learn it, and assumed that
was what freqsamp was supposed to do. 


0
Reply Pete 7/1/2010 4:19:28 AM

On 7/1/2010 12:17 AM, Pete Fraser wrote:
> "Jerry Avins"<jya@ieee.org>  wrote in message
> news:BdTWn.3387$OU6.1391@newsfe20.iad...
>
>> What is so good about a Butterworth amplitude response that you want to
>> imitate it?
>
> Nothing.
> It's just that, under some circumstances, equipment that is
> later in the signal chain expects a Butterworth response.

Usually, that includes the Butterworth phase response. If magnitude 
alone is enough, I suspect that it's just window dressing. YMMV

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply Jerry 7/1/2010 4:27:45 AM

Pete Fraser wrote:
> I'm new to Matlab, so this is probably a really dumb question.
> I'm trying to design an FIR where I can set the points.
> I want a symmetric FIR.
> I'm assuming that an N tap FIR gives me N/2 degrees of freedom
> for an even number of taps, and (N + 1)/2 degrees of freedom
> for an odd number of taps.
> This should allow me to set the frequency response at a
> number of points equal to the degrees of freedom, and
> have the frequency response travel through these points
> exactly. If I don't specify transitions that are too sharp
> relative to the number of taps, then the frequency response
> should be reasonably well controlled between the set points.
> 
> I had assumed that I could use the filter toolbox's freqsamp
> method to do this, but it doesn't seem to work as I expected.
> What am I doing wrong? I'd normally do this in Mathematica
> or C, but I'm trying to get into the spirit of Matlab.
> 
> I would want to solve the point setting problem anyway, but
> perhaps there's an easier way to do the immediate task that
> I wanted it for. I want to implement a Butterworth amplitude
> response in an FIR filter (to emulate some legacy equipment)
> but I can't find an easy way of doing that in Matlab (hence the
> point setting experiments).
> 
> Thanks
> 
> Pete 
> 
> 

Hi Pete!  How have you been?

Well, I looked into it a bit and found that I'm pretty far removed from 
the filter design toolbox / notation, etc.

What have you done so far?  any code samples?
What were the results?
Maybe that will get me going further.

Fred
0
Reply Fred 7/1/2010 4:51:38 AM

robert bristow-johnson wrote:
> On Jun 30, 10:39 pm, Tim Wescott <t...@seemywebsite.com> wrote:
>> On 06/30/2010 06:21 PM, Pete Fraser wrote:
>>
>>
>>
>>> I'm new to Matlab, so this is probably a really dumb question.
>>> I'm trying to design an FIR where I can set the points.
>>> I want a symmetric FIR.
>>> I'm assuming that an N tap FIR gives me N/2 degrees of freedom
>>> for an even number of taps, and (N + 1)/2 degrees of freedom
>>> for an odd number of taps.
>>> This should allow me to set the frequency response at a
>>> number of points equal to the degrees of freedom, and
>>> have the frequency response travel through these points
>>> exactly. If I don't specify transitions that are too sharp
>>> relative to the number of taps, then the frequency response
>>> should be reasonably well controlled between the set points.
>>> I had assumed that I could use the filter toolbox's freqsamp
>>> method to do this, but it doesn't seem to work as I expected.
>>> What am I doing wrong? I'd normally do this in Mathematica
>>> or C, but I'm trying to get into the spirit of Matlab.
>>> I would want to solve the point setting problem anyway, but
>>> perhaps there's an easier way to do the immediate task that
>>> I wanted it for. I want to implement a Butterworth amplitude
>>> response in an FIR filter (to emulate some legacy equipment)
>>> but I can't find an easy way of doing that in Matlab (hence the
>>> point setting experiments).
>> Make a vector that "draws" the frequency response that you want, then
>> take its FFT.
>>
>> It may not be the entirely kosher way to do it, but it certainly works.
> 
> Pete might want to use the "window method".  with an iFFT *larger*
> than N, draw in the frequency response, he'll get an impulse response
> that is longer than N, then window it to length N.  (then maybe FFT
> back to see how bad the window mangled things.)
> 
> r b-j
> 

r b-j,

Well, maybe I've had it wrong all these years but I'd say that the 
windowing method starts with N frequency samples where N is the length 
of the filter you want.  In fact, first order, it gives you a pretty 
good idea of how sharp the transitions can be, etc.  And, in reverse, 
that suggests the length of the filter you need.

If you start with some higher freq resolution then you will be *more* 
surprised if you window down to the desired N.  I'd not recommend doing 
it that way - but maybe there's some trick that I don't understand yet.

Anyway, having selected the frequency samples [the assumption is that 
they are on a regular grid from 0 to Fs-1) there are a couple of options 
of course:

1)  IFFT the samples / Window in time (i.e. not just the no-op trivial 
rectangular) to get rid of intersample "trillies" in frequency / FFT the 
result with lots of zero-padding to see what you've really got after 
windowing.

or

2)  Convolve the freq samples in frequency with the sinc-like FFT of the 
temporal window you have in mind to use.  If the "window" function is a 
lot longer than the filter sample sequence then you'll also get the 
interpolated values between samples to see what you've really got in 
frequency after windowing.

Fred
0
Reply Fred 7/1/2010 2:59:36 PM


Pete Fraser wrote:

> "Vladimir Vassilevsky" <nospam@nowhere.com> wrote in message 
> news:0Jidnc3knONVbbbRnZ2dnUVZ_oSdnZ2d@giganews.com...
> 
> 
>>But why would anyone want to emulate Butterworth magnitude response in a 
>>linear phase FIR filter whereas it is very straightforward to do the 
>>(**exact) Butterworth as IIR filter?
> 
> 
> Butterworth will be one of several possible responses
> from the same hardware.
> 
> I'd rather not deal with Butterworth's phase issues.

Check the FDLS method of Greg Berchin. That is one of the most elegant 
methods of designing to a prototype.


Vladimir Vassilevsky
DSP and Mixed Signal Design Consultant
http://www.abvolt.com
0
Reply Vladimir 7/1/2010 3:55:46 PM

"Fred Marshall" <fmarshallx@remove_the_xacm.org> wrote in message 
news:BIednY5EMdtRg7HRnZ2dnUVZ_qmdnZ2d@centurytel.net...

> Hi Pete!  How have you been?

I'm doing well thanks.
Working on a couple of interesting projects.

>
> Well, I looked into it a bit and found that I'm pretty far removed from 
> the filter design toolbox / notation, etc.
>
> What have you done so far?  any code samples?
> What were the results?
> Maybe that will get me going further.

I think I'll just do it the brute force way (using Matlab),
then post the code I used.

Hopefully then folks can tell me how I _should_ have done it.

Pete 


0
Reply Pete 7/1/2010 5:02:28 PM

"Vladimir Vassilevsky" <nospam@nowhere.com> wrote in message 
news:OP-dnf44bO7hJ7HRnZ2dnUVZ_t-dnZ2d@giganews.com...

> Check the FDLS method of Greg Berchin. That is one of the most elegant 
> methods of designing to a prototype.

Thanks.
That looks like a very useful tool.
I'll play with the code this morning.

I'm still curious about freqsamp though.
WhereGreg talks about the pseudoinverse, I had intended
to give Matlab the exact number of points to correctly
determine the inverse. Obviously I can do this with an explicit
matrix inversion, but had assumed that freqsamp was intended
to do that. Did I just assume wrong, or am I using it incorrectly?

Thanks

Pete 


0
Reply Pete 7/1/2010 5:23:44 PM

Pete Fraser <pfraser@covad.net> wrote:

>I'm still curious about freqsamp though.

Near as I can tell there is not sufficient information in the
matlab documentation to state exactly (or perhaps, even generally)
what they are doing here.

If you're lucky someone from Mathworks will see your question
and answer it for you.  I believe you can also ask this on the Mathworks
site, after registering.


Steve
0
Reply spope33 7/1/2010 5:42:51 PM

"Steve Pope" <spope33@speedymail.org> wrote in message 
news:i0ik2r$85g$2@blue.rahul.net...
> Near as I can tell there is not sufficient information in the
> matlab documentation to state exactly (or perhaps, even generally)
> what they are doing here.

I thought perhaps I was looking in the wrong area.
The Mathworks help is a model of clarity, with lots of
well-thought-out examples, but when I try to dig a
little deeper I come up short.

>
> If you're lucky someone from Mathworks will see your question
> and answer it for you.  I believe you can also ask this on the Mathworks
> site, after registering.

I'll try that. Thanks. 


0
Reply Pete 7/1/2010 6:05:57 PM

"Vladimir Vassilevsky" <nospam@nowhere.com> wrote in message 
news:OP-dnf44bO7hJ7HRnZ2dnUVZ_t-dnZ2d@giganews.com...

> Check the FDLS method of Greg Berchin. That is one of the most elegant 
> methods of designing to a prototype.

I think I must be doing something dumb.
I started off with a really simple half-band FIR, just
as a test.

The input file was:

11     100000
0              1     0         1
5000        1     -72      1
10000     1      -144    1
15000     1      -216    1
20000     1      -288    1
25000     .5     -360    1
30000     0      -432    1
35000     0      -504    1
40000     0      -576    1
45000     0      -648    1
50000     0      -720    1


N=8
D=0
delta=0

and I get B =


   -3.55068459956424e-018
       -0.0454056666613377
     2.42666919829597e-016
         0.292936161872951
                       0.5
         0.326331656363899
    -2.00804673259988e-016
       -0.0874864446451108
    -1.76118166519034e-017

so it's doing a lot of things right.

The center is 0.5, four are close to 0,
but 0.29 is not 0.32, and -0.04 is not -0.08.

What am I doing wrong?

Thanks

Pete 


0
Reply Pete 7/2/2010 1:48:22 AM

On Thu, 1 Jul 2010 18:48:22 -0700, "Pete Fraser" <pfraser@covad.net> wrote:

>I think I must be doing something dumb.
>I started off with a really simple half-band FIR, just
>as a test.
....
>The center is 0.5, four are close to 0,
>but 0.29 is not 0.32, and -0.04 is not -0.08.
>
>What am I doing wrong?

You are expecting a generic least-squares fit to produce a symmetric filter,
without telling the algorithm to expect symmetry.  If you know that your filter
is symmetric, you should constrain FDLS to guarantee it:

Standard FDLS formulation (set up for zero-order denominator):

                      -1            -n
               b + b z   + ... + b z
        y(z)    0   1             n
        ---- = -----------------------
        u(z)              1


Constrained formulation:

                      -1            -(n-1)      -n
               b + b z   + ... + b z       + b z
        y(z)    0   1             1           0
        ---- = -----------------------------------
        u(z)                  1

Greg
0
Reply Greg 7/2/2010 11:05:36 AM

In fact, if you know that your filter is halfband (every other coefficient is
zero), you should also constrain for that:

Standard FDLS formulation (set up for zero-order denominator):

                      -1            -n
               b + b z   + ... + b z
        y(z)    0   1             n
        ---- = -----------------------
        u(z)              1

Constrained formulation:

                      -2            -(n-2)      -n
               b + b z   + ... + b z       + b z
        y(z)    0   2             2           0
        ---- = -----------------------------------
        u(z)                  1

Always try to use as much information as as you have.

Greg
0
Reply Greg 7/2/2010 11:19:34 AM

"Greg Berchin" <gberchin@comicast.net.invalid> wrote in message 
news:uchr26djlcjnmrmvi4k94ue126i7sf9vov@4ax.com...

> You are expecting a generic least-squares fit to produce a symmetric 
> filter,
> without telling the algorithm to expect symmetry.

I thought I was doing that be specifying a linear phase response but,
of course, the FDLS is just using the phase response as a goal.
It can do a better match for my requested amplitude response by
using an asymmetric kernel, so it does.
Is that correct?

> If you know that your filter
> is symmetric, you should constrain FDLS to guarantee it:

I'm using version 2.0 (10/21/2006) and I don't see any way
of constraining the filter to be symmetrical (or will I need to
tweak the code)?

Thanks for a nice tool and tutorial.

Pete 


0
Reply Pete 7/2/2010 4:46:15 PM

On Fri, 2 Jul 2010 09:46:15 -0700, "Pete Fraser" <pfraser@covad.net> wrote:

>I thought I was doing that be specifying a linear phase response but,
>of course, the FDLS is just using the phase response as a goal.
>It can do a better match for my requested amplitude response by
>using an asymmetric kernel, so it does.
>Is that correct?

It's just a pseudoinverse and it's just data.  Truncation and rounding,
numerical ill-conditioning, etc., all affect the result.

Also, you are trying to identify nine coefficients with only eleven data points.
It should be sufficient, but give it some more data to work with anyway.

>I'm using version 2.0 (10/21/2006) and I don't see any way
>of constraining the filter to be symmetrical (or will I need to
>tweak the code)?

The Matlab code is only intended to be a starting point, an introductory
example.  The possibilities for variations on FDLS are limitless.  So modify the
code to fit your specific problem.  I won't mind!

>Thanks for a nice tool and tutorial.

You're welcome.

Greg
0
Reply Greg 7/2/2010 4:51:22 PM

>
>
>Pete Fraser wrote:
>
>> "Vladimir Vassilevsky" <nospam@nowhere.com> wrote in message 
>> news:0Jidnc3knONVbbbRnZ2dnUVZ_oSdnZ2d@giganews.com...
>> 
>> 
>>>But why would anyone want to emulate Butterworth magnitude response in a

>>>linear phase FIR filter whereas it is very straightforward to do the 
>>>(**exact) Butterworth as IIR filter?
>> 
>> 
>> Butterworth will be one of several possible responses
>> from the same hardware.
>> 
>> I'd rather not deal with Butterworth's phase issues.
>

Maybe this could help:
http://www.icsi.berkeley.edu/~storn/fiwiz.html

In the 'Z plane design' option, you pretty much can give it your desired
magnitude response and press 'search'...

At the very least, it's quite interesting to look at the genetic algorithm
converging!

Dave
0
Reply gretzteam 7/2/2010 7:33:12 PM

Pete Fraser wrote:

> I'd rather not deal with Butterworth's phase issues.

A good thing to remember is that any filter that is as
selective as a given Butterworth filter will have a 
similar RMS delay spread as the Butterworth filter.
So no free lunch on phase issues in many cases.


Steve
0
Reply spope33 7/2/2010 7:39:47 PM

"Steve Pope" <spope33@speedymail.org> wrote in message 
news:i0lfa3$rk9$1@blue.rahul.net...

> A good thing to remember is that any filter that is as
> selective as a given Butterworth filter will have a
> similar RMS delay spread as the Butterworth filter.
> So no free lunch on phase issues in many cases.

But I'm doing a symmetric FIR with the Butterworth
amplitude response. 


0
Reply Pete 7/2/2010 8:33:04 PM

Pete Fraser <pfraser@covad.net> wrote:

>"Steve Pope" <spope33@speedymail.org> wrote in message 

>> A good thing to remember is that any filter that is as
>> selective as a given Butterworth filter will have a
>> similar RMS delay spread as the Butterworth filter.
>> So no free lunch on phase issues in many cases.

>But I'm doing a symmetric FIR with the Butterworth
>amplitude response. 

Linear phase is overrated.  :-)


S.
0
Reply spope33 7/2/2010 8:42:58 PM

On Jul 2, 4:42=A0pm, spop...@speedymail.org (Steve Pope) wrote:
> Pete Fraser <pfra...@covad.net> wrote:
> >"Steve Pope" <spop...@speedymail.org> wrote in message
> >> A good thing to remember is that any filter that is as
> >> selective as a given Butterworth filter will have a
> >> similar RMS delay spread as the Butterworth filter.
> >> So no free lunch on phase issues in many cases.
> >But I'm doing a symmetric FIR with the Butterworth
> >amplitude response.

so then you'll get even *more* phase shift delay.  (Butterworths are
minimum phase and linear phase is not minimum phase unless it's a wire
or a simple gain.)

> Linear phase is overrated. =A0:-)

yeah, and what Steve said.

r b-j
0
Reply robert 7/2/2010 11:05:34 PM

On Jul 1, 10:59=A0am, Fred Marshall <fmarshallx@remove_the_xacm.org>
wrote:
>
> r b-j,
>
> Well, maybe I've had it wrong all these years but I'd say that the
> windowing method starts with N frequency samples where N is the length
> of the filter you want.

so then, since the DFT and iDFT are bijective (i love using fancy-
pants words), why window?  if h[n] has N samples and N degrees of
freedom, so does H[k].  you specify your N frequency samples and you
can hit it perfectly with no windowing.

so, that seems curious to me.

r b-j
0
Reply robert 7/2/2010 11:09:39 PM

Steve Pope wrote:
> A good thing to remember is that any filter that is as
> selective as a given Butterworth filter will have a
> similar RMS delay spread as the Butterworth filter.
> So no free lunch on phase issues in many cases.

Pete Fraser replied:
>But I'm doing a symmetric FIR with the Butterworth
>amplitude response.

Then robert bristow-johnson said:
> so then you'll get even *more* phase shift delay.  (Butterworths are
> minimum phase and linear phase is not minimum phase unless it's a wire
> or a simple gain.)

True, but Steve and I were talking about RMS delay spread.

Pete 


0
Reply Pete 7/2/2010 11:55:05 PM

On Jul 2, 7:55=A0pm, "Pete Fraser" <pfra...@covad.net> wrote:
>
> ... but Steve and I were talking about RMS delay spread.
>

"spread" from what?  zero?  or some mean delay?  how is the mean
defined?

"RMS delay spread" hasn't been a parameter i've been familiar with.
is it a measure of overall phase nonlinearity?

just curious.

r b-j


0
Reply robert 7/3/2010 3:37:14 AM

robert bristow-johnson  <rbj@audioimagination.com> wrote:

>On Jul 2, 7:55�pm, "Pete Fraser" <pfra...@covad.net> wrote:

>> ... but Steve and I were talking about RMS delay spread.

>"spread" from what?  zero?  or some mean delay?  how is the mean
>defined?

>"RMS delay spread" hasn't been a parameter i've been familiar with.

>just curious.

Interestingly some of the literature defintions of RMS delay spread
I disagree with, but it basically goes something like this: if the impulse
response is:

      (sum over i) (h(i) * (t - tau(i)))

then the average delay is

     ave = (sum over i)(h(i)*tau(i) / (sum over i)(h(i))

and the RMS delay spread is

    sqrt((sum over i)(h(i)*((tau(i) - ave)^2)).

>is it a measure of overall phase nonlinearity?

No, the phase can be linear but you can have a large RMS delay spread.

I'm interested if anyone doesn't like the expression I just gave above.
I've seen expressions that were not equivalent to this (I'm going to
say by Spencer).


Steve
0
Reply spope33 7/3/2010 3:52:56 AM

Oops, let me fix this.  I left out some absolute values.

the impulse response is:

      (sum over i) (h(i) * (t - tau(i)))

then the average delay is

     ave = (sum over i)(abs(h(i))*tau(i) / (sum over i)(abs(h(i)))

and the RMS delay spread is

    sqrt((sum over i)(abs(h(i))*((tau(i) - ave)^2)).



S.
0
Reply spope33 7/3/2010 3:58:44 AM

On Fri, 2 Jul 2010 20:37:14 -0700 (PDT), robert bristow-johnson
<rbj@audioimagination.com> wrote:

>> ... but Steve and I were talking about RMS delay spread.
>>
>
>"spread" from what?  zero?  or some mean delay?  how is the mean
>defined?
>
>"RMS delay spread" hasn't been a parameter i've been familiar with.
>is it a measure of overall phase nonlinearity?

Possibly useful:  see US Patent 5375067 for a discussion of Central Time and RMS
Delay.

Greg
0
Reply Greg 7/3/2010 10:48:57 AM

On Sat, 03 Jul 2010 06:48:57 -0400, Greg Berchin <gberchin@comicast.net.invalid>
wrote:

>Possibly useful:  see US Patent 5375067 for a discussion of Central Time and RMS
>Delay.

I gotta stop typing before my brain is awake.

Make that Central Time and RMS *Duration*.
0
Reply Greg 7/3/2010 11:10:26 AM

On Jul 3, 6:48=A0am, Greg Berchin <gberc...@comicast.net.invalid> wrote:
> On Fri, 2 Jul 2010 20:37:14 -0700 (PDT), robert bristow-johnson
>
> <r...@audioimagination.com> wrote:
> >> ... but Steve and I were talking about RMS delay spread.
>
> >"spread" from what? =A0zero? =A0or some mean delay? =A0how is the mean
> >defined?
>
> >"RMS delay spread" hasn't been a parameter i've been familiar with.
> >is it a measure of overall phase nonlinearity?
>
> Possibly useful: =A0see US Patent 5375067 for a discussion of Central Tim=
e and RMS
> Duration.
>


Greg, i didn't know about this patent of yours.  good show!

well, i'm still trying to pin these guys down.  if their r.m.s. delay
parameter ("spread" or whatever it's called) is a measure of deviation
from phase linearity, then nothing can beat a phase linear FIR.  but
if it's a measure of phase delay or group delay (from zero), then a
minimum-phase filter will beat linear-phase.

that's all i wanted to say about it.

r b-j


0
Reply robert 7/3/2010 6:34:44 PM

robert bristow-johnson  <rbj@audioimagination.com> wrote:
>well, i'm still trying to pin these guys down.  if their r.m.s. delay
>parameter ("spread" or whatever it's called) is a measure of deviation
>from phase linearity, then nothing can beat a phase linear FIR.  but
>if it's a measure of phase delay or group delay (from zero), then a
>minimum-phase filter will beat linear-phase.

It is neither.

Steve
0
Reply spope33 7/3/2010 6:35:28 PM

"robert bristow-johnson" <rbj@audioimagination.com> wrote in message 
news:c2f65e9b-74e2-41e3-85f0-61684db85f48@k39g2000yqd.googlegroups.com...

> well, i'm still trying to pin these guys down.  if their r.m.s. delay
> parameter ("spread" or whatever it's called) is a measure of deviation
> from phase linearity, then nothing can beat a phase linear FIR.

I don't know if there's a formal definition.
I just assumed it was a measure of the dispersiveness
(frequency dependent delay) of the filter.

> but if it's a measure of phase delay or group delay (from zero),
> then minimum-phase filter will beat linear-phase.

True, but I'm not sure why anybody would refer to that
as "spread".

Pete 


0
Reply Pete 7/3/2010 6:55:44 PM

Steve Pope <spope33@speedymail.org> wrote:

>robert bristow-johnson  <rbj@audioimagination.com> wrote:

>>well, i'm still trying to pin these guys down.  if their r.m.s. delay
>>parameter ("spread" or whatever it's called) is a measure of deviation
>>from phase linearity, then nothing can beat a phase linear FIR.  but
>>if it's a measure of phase delay or group delay (from zero), then a
>>minimum-phase filter will beat linear-phase.

>It is neither.

Here's a reference:

Ibnkahla, _Signal Processing for Mobile Communications Handbook_,
section 1.2.2.1.4.  It is viewable on Google Books.

He used the same expression for average delay I stated, but
his expression for RMS Delay Spread is a little different,
and I don't like it as much, but I think his expression is
as close to a standard definition as you will get.

RMS delay spread must be used carefully; there are responses with
large RMS delay spreads which do not have bad phase qualities for
a most signals, i.e.  something like

      h(t) = 1*(t) + 100*(t - 20) - 100*(t - 20.0001)



S.
0
Reply spope33 7/3/2010 7:04:53 PM

robert bristow-johnson wrote:
> On Jul 1, 10:59 am, Fred Marshall <fmarshallx@remove_the_xacm.org>
> wrote:
>> r b-j,
>>
>> Well, maybe I've had it wrong all these years but I'd say that the
>> windowing method starts with N frequency samples where N is the length
>> of the filter you want.
> 
> so then, since the DFT and iDFT are bijective (i love using fancy-
> pants words), why window?  if h[n] has N samples and N degrees of
> freedom, so does H[k].  you specify your N frequency samples and you
> can hit it perfectly with no windowing.
> 
> so, that seems curious to me.
> 
> r b-j

You window because the frequency response between those points can be 
nasty.  Either you don't window and convolve those samples with a 
matching sinc er.. Dirichlet or you window and convolve those samples 
with something else.

The trade is that the Dirichlet matches all the samples exactly.  Other 
windows don't.  Some match them all except the adjacent two as in the 
1/2  1  1/2 sum of adjacent Dirichlets.  It still has zeros in the sum 
at all the other sample points - so the convolution doesn't perturb 
their values.  And, the values between points are better behaved.

And making N larger to begin with doesn't affect any of this in my way 
of looking at it.

Fred
0
Reply Fred 7/3/2010 10:21:10 PM

Pete Fraser <pfraser@covad.net> wrote:

>"robert bristow-johnson" <rbj@audioimagination.com> wrote in message 

>> well, i'm still trying to pin these guys down.  if their r.m.s. delay
>> parameter ("spread" or whatever it's called) is a measure of deviation
>> from phase linearity, then nothing can beat a phase linear FIR.

>I don't know if there's a formal definition.
>I just assumed it was a measure of the dispersiveness
>(frequency dependent delay) of the filter.

So as to save you guys some googling I have uploaded Ibnkahla's
definition of RMS delay spread to the following link which
I will leave in place for a few days.

Hopefully this makes things less ambiguous.

http://www.rahul.net/spp/IbnkahlaDelaySpread.bmp

Steve
0
Reply spope33 7/4/2010 4:15:35 PM

On Sun, 4 Jul 2010 16:15:35 +0000 (UTC), spope33@speedymail.org (Steve Pope)
wrote:

>So as to save you guys some googling I have uploaded Ibnkahla's
>definition of RMS delay spread to the following link which
>I will leave in place for a few days.
>
>Hopefully this makes things less ambiguous.
>
>http://www.rahul.net/spp/IbnkahlaDelaySpread.bmp

I only gave it a cursory look, but comparing that with the discussion of moments
that I put into US Patent 5375067, Ibnkahla's definition of RMS delay spread
appears to be the same as RMS Duration.

Greg
0
Reply Greg 7/4/2010 5:35:10 PM

On Jul 4, 1:35=A0pm, Greg Berchin <gberc...@comicast.net.invalid> wrote:
> On Sun, 4 Jul 2010 16:15:35 +0000 (UTC), spop...@speedymail.org (Steve Po=
pe)
> wrote:
>
> >So as to save you guys some googling I have uploaded Ibnkahla's
> >definition of RMS delay spread to the following link which
> >I will leave in place for a few days.
>
> >Hopefully this makes things less ambiguous.

might be if i knew what P(t) is.

> >http://www.rahul.net/spp/IbnkahlaDelaySpread.bmp
>
> I only gave it a cursory look, but comparing that with the discussion of =
moments
> that I put into US Patent 5375067, Ibnkahla's definition of RMS delay spr=
ead
> appears to be the same as RMS Duration.

can you define that in terms of impulse response, h(t)?

otherwise i'm still confused.

r b-j
0
Reply robert 7/5/2010 12:30:48 AM

robert bristow-johnson  <rbj@audioimagination.com> wrote:

>> On Sun, 4 Jul 2010 16:15:35 +0000 (UTC), spop...@speedymail.org (Steve Pope)

>> >So as to save you guys some googling I have uploaded Ibnkahla's
>> >definition of RMS delay spread to the following link which
>> >I will leave in place for a few days.

>> >Hopefully this makes things less ambiguous.

>might be if i knew what P(t) is.

For the full details, enough of Ibnakahla's text (see my previous
post) is on Google books to extract those.

I think he wants P(tau) to be the power delay profile, but another
possibility is the just magnitude of the impulse response.

(So, yes, still ambiguous.)

S.
0
Reply spope33 7/5/2010 12:50:49 AM

On Sun, 4 Jul 2010 17:30:48 -0700 (PDT), robert bristow-johnson
<rbj@audioimagination.com> wrote:

>can you define that in terms of impulse response, h(t)?

In the case of the RMS Duration, it is simply the normalized second moment of a
generic time domain waveform, so h(t) will do nicely.

In the case of Ibnkahla's RMS delay spread, I interpreted P(tau) to be something
like h(t-tau), where tau is the temporal centroid.  If P(t) is something else,
then I don't know.

Greg
0
Reply Greg 7/5/2010 1:03:15 AM

Fred Marshall wrote:
> robert bristow-johnson wrote:
>> On Jul 1, 10:59 am, Fred Marshall <fmarshallx@remove_the_xacm.org>
>> wrote:
>>> r b-j,
>>>
>>> Well, maybe I've had it wrong all these years but I'd say that the
>>> windowing method starts with N frequency samples where N is the length
>>> of the filter you want.
>>
>> so then, since the DFT and iDFT are bijective (i love using fancy-
>> pants words), why window?  if h[n] has N samples and N degrees of
>> freedom, so does H[k].  you specify your N frequency samples and you
>> can hit it perfectly with no windowing.
>>
>> so, that seems curious to me.
>>
>> r b-j
> 
> You window because the frequency response between those points can be 
> nasty.  Either you don't window and convolve those samples with a 
> matching sinc er.. Dirichlet or you window and convolve those samples 
> with something else.
> 
> The trade is that the Dirichlet matches all the samples exactly.  Other 
> windows don't.  Some match them all except the adjacent two as in the 
> 1/2  1  1/2 sum of adjacent Dirichlets.  It still has zeros in the sum 
> at all the other sample points - so the convolution doesn't perturb 
> their values.  And, the values between points are better behaved.
> 
> And making N larger to begin with doesn't affect any of this in my way 
> of looking at it.
> 
> Fred

And, I'm sure you know ..

If you use the basis 1/2 1 1/2 sum of Dirichlets (raised cosine in time) 
which I will refer to as "RC" then you do have to convolve and can't 
solve a system of equations that matches the samples exactly at each 
point.  If you did that, then it would revert to using a single 
Dirichlet I do believe.

Since the stopbands are zero then all you get in them are the 
fast-decaying tails of RC components from the passbands by 1/f^3 [rather 
than 1/f for a single Dirichlet].

Fred
0
Reply Fred 7/5/2010 4:54:07 PM

On Jul 3, 6:21=A0pm, Fred Marshall <fmarshallx@remove_the_xacm.org>
wrote:
> robert bristow-johnson wrote:
> > On Jul 1, 10:59 am, Fred Marshall <fmarshallx@remove_the_xacm.org>
> > wrote:
> >> r b-j,
>
> >> Well, maybe I've had it wrong all these years but I'd say that the
> >> windowing method starts with N frequency samples where N is the length
> >> of the filter you want.
>
> > so then, since the DFT and iDFT are bijective (i love using fancy-
> > pants words), why window? =A0if h[n] has N samples and N degrees of
> > freedom, so does H[k]. =A0you specify your N frequency samples and you
> > can hit it perfectly with no windowing.
>
> > so, that seems curious to me.
>
> > r b-j
>
> You window because the frequency response between those points can be
> nasty. =A0Either you don't window and convolve those samples with a
> matching sinc er.. Dirichlet or you window and convolve those samples
> with something else.

you're always windowing it (unless your impulse response is an
infinitely-repeating periodic function.  you draw your N (or more)
discrete points in the frequency domain, you iDFT it and what you have
there is still a periodic sequence.  if you use just one cycle of that
(and zero the rest, as you would for an FIR), then you've applied the
rectangular window.

> The trade is that the Dirichlet matches all the samples exactly. =A0Other
> windows don't. =A0Some match them all except the adjacent two as in the
> 1/2 =A01 =A01/2 sum of adjacent Dirichlets. =A0It still has zeros in the =
sum
> at all the other sample points - so the convolution doesn't perturb
> their values. =A0And, the values between points are better behaved.

the values between points exhibit less "high frequency" behavior or
are less smooth.  that's normally thought of as "better".  and the
reason is that in the other domain (where the window is applied),
you've reduced the amplitude of the points further from h[0] (the
"high frequency" points, but they're really the high time-displacement
points that fuel the "high frequency" wiggling in the frequency
domain.

> And making N larger to begin with doesn't affect any of this in my way
> of looking at it.

it allows you to draw the frequency response more densely.  so you can
explicitly specify how "wiggly" you want it between the sparser points
you have with your smaller "N".  but to do so, the impulse response is
very large, too large.  so you have to shorten it and that is
windowing.  then, starting with the ideal (but too long) impulse
response, you can try shortening it in a variety of different ways
(using different windows) and see how well you do.  you can also try
shortening it to different lengths (using whatever window that looks
best) to make a tradeoff decision on FIR length and its performance.

it's very similar in philosophy to what we do with the windowed-sinc
LPF design.  start with the ideal (which is too long) and see what
happens when you settle for something shorter.

r b-j


0
Reply robert 7/5/2010 5:18:54 PM

robert bristow-johnson wrote:
> On Jul 3, 6:21 pm, Fred Marshall <fmarshallx@remove_the_xacm.org>
> wrote:
>> robert bristow-johnson wrote:
>>> On Jul 1, 10:59 am, Fred Marshall <fmarshallx@remove_the_xacm.org>
>>> wrote:
>>>> r b-j,
>>>> Well, maybe I've had it wrong all these years but I'd say that the
>>>> windowing method starts with N frequency samples where N is the length
>>>> of the filter you want.
>>> so then, since the DFT and iDFT are bijective (i love using fancy-
>>> pants words), why window?  if h[n] has N samples and N degrees of
>>> freedom, so does H[k].  you specify your N frequency samples and you
>>> can hit it perfectly with no windowing.
>>> so, that seems curious to me.
>>> r b-j
>> You window because the frequency response between those points can be
>> nasty.  Either you don't window and convolve those samples with a
>> matching sinc er.. Dirichlet or you window and convolve those samples
>> with something else.
> 
> you're always windowing it (unless your impulse response is an
> infinitely-repeating periodic function.  you draw your N (or more)
> discrete points in the frequency domain, you iDFT it and what you have
> there is still a periodic sequence.  if you use just one cycle of that
> (and zero the rest, as you would for an FIR), then you've applied the
> rectangular window.
> 
>> The trade is that the Dirichlet matches all the samples exactly.  Other
>> windows don't.  Some match them all except the adjacent two as in the
>> 1/2  1  1/2 sum of adjacent Dirichlets.  It still has zeros in the sum
>> at all the other sample points - so the convolution doesn't perturb
>> their values.  And, the values between points are better behaved.
> 
> the values between points exhibit less "high frequency" behavior or
> are less smooth.  that's normally thought of as "better".  and the
> reason is that in the other domain (where the window is applied),
> you've reduced the amplitude of the points further from h[0] (the
> "high frequency" points, but they're really the high time-displacement
> points that fuel the "high frequency" wiggling in the frequency
> domain.
> 
>> And making N larger to begin with doesn't affect any of this in my way
>> of looking at it.
> 
> it allows you to draw the frequency response more densely.  so you can
> explicitly specify how "wiggly" you want it between the sparser points
> you have with your smaller "N".  but to do so, the impulse response is
> very large, too large.  so you have to shorten it and that is
> windowing.  then, starting with the ideal (but too long) impulse
> response, you can try shortening it in a variety of different ways
> (using different windows) and see how well you do.  you can also try
> shortening it to different lengths (using whatever window that looks
> best) to make a tradeoff decision on FIR length and its performance.
> 
> it's very similar in philosophy to what we do with the windowed-sinc
> LPF design.  start with the ideal (which is too long) and see what
> happens when you settle for something shorter.
> 
> r b-j
> 
> 

r b-j,

OK.  That interesting - the approaches have the same underlying 
ingredients and are only different by virtue of how the starting points 
are viewed.

In the "infinite-N" case you would be seeking N based on responses and, 
I guess would decide on a non-rectangular window (or windows) to use in 
the process.  As soon as you do the first test case it's the same as the 
"fixed-N" approach.  It rather assumes that you don't know N very well, 
right?

In the "fixed-N" case you would be seeking acceptable response by 
choosing non-rectangular window functions of length N.  And then, if you 
aren't satisfied, iterate N.  Just as the "infinite-N" case does in the end.

I think of both as the "windowing method".  How you approach it depends 
on how closely you figure you know N and your own preference.  Either 
way, just like using P-M, you have to choose an N at some point.

My way of roughly choosing N is to look at the minimum transition band 
width W and spec N so that the temporal length of the filter is approx 
1.5/W.  And, for even rougher estimates just 1/W which is what I usually 
use when discussing the "physics".  You get to pick based on where you 
put the break points in the transition - so that choice usually gets 
mixed up with band ripple - particularly with minimax criterion.  The 
lower the ripple, the wider the transition and the longer the filter. 
But, it's sure that the filter can't be any *shorter*!!  :-)

And, then, if you like you can use one of those length-estimating formulas.

Fred
0
Reply fmarshallx (24) 7/5/2010 8:30:29 PM

48 Replies
213 Views

(page loaded in 1.056 seconds)


Reply: