Adaptive Filter | ADPCM | ADSP | ADSP-2181 | Aliasing | AMR | Anti-Aliasing | ARMA | Autocorrelation | AutoCovariance | Beamforming | Bessel | Blackfin | Butterworth | C6713 | CCS | Chebyshev | CIC Filter | Circular Convolution | Code Composer Studio | Comb Filter | Compression | Convolution | Cross Correlation | DCT | Decimation | Deconvolution | Demodulation | DM642 | DSP Boards | DSP/BIOS | DTMF | Echo Cancellation | Equalization | Equalizer | ETSI | EZLITE (Ez-kit Lite) | FFT | FFTW | FIR Filter | Fixed Point | FSK | G.711 | G.723 | G.729 | Gaussian Noise | Goertzel | GPIO | Hilbert Transform | IFFT | IIR Filter | Interpolation | Invariance | JTAG | Kalman | Laplace Transform | Levinson | LPC | McBSP | MIPS | Modulation | MPEG | Multirate | Notch Filter | Nyquist | OFDM | Oversampling | Pink Noise | Pitch | PLL | Polyphase | QAM | QDMA | Quantization | Quantizer | Radar | Random Noise | Reed Solomon | Remez | Resampling | RTDX | Sampling | Sharc | TI C6711 | Undersampling | Viterbi | Wavelets | White Noise | Wiener Filter | Windowing | XDS510PP | Z Transform

in article 1...@c13g2000cwb.googlegroups.com, Jaime A.
Aranguren C., SanJaaC Electronics at j...@ieee.org wrote on
01/10/2005 10:51:

> My sample rate is 240 Hz. The filters is some sort of band pass filter,
> bell shaped. Frequencies of interest are from 0.05 Hz up to 40Hz, being
> the maximum at 8.8 Hz.

how is it that 8.8 Hz is "maximum" instead of 40 Hz?

 
> Could you please give me some example on how to get a IIR filter
> designed with more terms than what's proided by MATLAB's bilinear
> function? This one gives me a 6th order IIR filter

then your analog filter is 6th order also.

at http://hyperphysics.phy-astr.gsu.edu/hbase/math/lnseries.html there is a
list of series for ln(x).


 ln(x) = 2*[(x-1)/(x+1) + (1/3)*((x-1)/(x+1))^3 + (1/5)*((x-1)/(x+1))^5 ...]


the exact mapping from s to z is:

    z = e^(s*T)

and

    s = (1/T) * ln(z) = Fs * ln(z)

the BLT take the first term of the series above for an approximation of s:

    s = Fs * ln(z) ~= 2*Fs*[(z-1)/(z+1)]

with another term:

    s  ~= 2*Fs*[(z-1)/(z+1) + (1/3)*((z-1)/(z+1))^3]

now you're using an 18th order digital filter to approximate your analog
filter.  but it should be a better approximation.  now i am not sure, but i
think i heard that adding that 2nd term made an unstable digital filter even
though the poles in the s-plane were stable.  i do not know if this is right
or not.  but it doesn't seem to me to be so.  we can find out ourselves by
trying it.


-- 

r b-j                  r...@audioimagination.com

"Imagination is more important than knowledge."

______________________________
What are the most bookmarked web sites about "fpga"? Find out here.

"Jon Harris" <g...@hotmail.com> wrote in message 
news:3...@individual.net...
>
> Another question for r b-j and others, I've heard of using an alternate
> definition for the bilinear transform in order to alter the mapping 
> between the
> s and z domains but don't know much about this topic.  For example, the 
> standard
> BLT yields a good match up until maybe about half Nyquist and then it 
> diverges,
> reaching the analog value at infinity at Nyquist.  Could one use an 
> alternate
> mapping that for example made the digital filter match well up to say
> .95*Nyquist and then diverge _really_ fast?  Or maybe even more unusual to 
> have
> a mapping where the bulk of the error was in the low frequency terms?  How 
> does
> one go about this?
>
>

Hello Jon,

With the BLT, there is only one parameter to adjust. We had a thread a while 
back on the BLT, where I described how to match this up. Now certainly one 
may use a mapping other than the BLT to map between the s and z domains. And 
R-BJ gave one extension to the BLT. However by the time this much effort is 
done, I'd just go and design the best fit (given whatever constraints are 
important to the problem) IIR via a Fletcher-Powell method.

Clay

______________________________
Are you a member of DSPRelated.com? If you are, you can now access a powerful tool to discover new web sites. Details here.

On 10 Jan 2005 07:51:15 -0800, "Jaime A. Aranguren C., SanJaaC
Electronics" <j...@ieee.org> wrote:

>>My sample rate is 240 Hz. The filters is some sort of band pass filter,
>>bell shaped. Frequencies of interest are from 0.05 Hz up to 40Hz, being
>>the maximum at 8.8 Hz.

Can you post the Laplace-domain transfer function?  Or at least a
set of frequency/magnitude/phase measurements?

Greg

______________________________
Discover some of the greatest websites! Details here.

Hi Robert,

robert bristow-johnson wrote:
> > My sample rate is 240 Hz. The filters is some sort of band pass
filter,
> > bell shaped. Frequencies of interest are from 0.05 Hz up to 40Hz,
being
> > the maximum at 8.8 Hz.
>
> how is it that 8.8 Hz is "maximum" instead of 40 Hz?

I meant the maximun of the bell shaped magnitude response is located at
8.8 Hz.

> > Could you please give me some example on how to get a IIR filter
> > designed with more terms than what's proided by MATLAB's bilinear
> > function? This one gives me a 6th order IIR filter
>
> then your analog filter is 6th order also.
>
> at http://hyperphysics.phy-astr.gsu.edu/hbase/math/lnseries.html
there is a
> list of series for ln(x).
>
>
>  ln(x) = 2*[(x-1)/(x+1) + (1/3)*((x-1)/(x+1))^3 +
(1/5)*((x-1)/(x+1))^5 ...]
>
>
> the exact mapping from s to z is:
>
>     z = e^(s*T)
>
> and
>
>     s = (1/T) * ln(z) = Fs * ln(z)
>
> the BLT take the first term of the series above for an approximation
of s:
>
>     s = Fs * ln(z) ~= 2*Fs*[(z-1)/(z+1)]
>
> with another term:
>
>     s  ~= 2*Fs*[(z-1)/(z+1) + (1/3)*((z-1)/(z+1))^3]
>
> now you're using an 18th order digital filter to approximate your
analog
> filter.  but it should be a better approximation.  now i am not sure,
but i
> think i heard that adding that 2nd term made an unstable digital
filter even
> though the poles in the s-plane were stable.  i do not know if this
is right
> or not.  but it doesn't seem to me to be so.  we can find out
ourselves by
> trying it.

I'll try two things: Google for "Frequency domain filtering", which
gave a bunch of hits. Also, I'll take my basic DSP books, and try to
design my filters with the second aprox for the s.

Regards,

JaaC

______________________________
Discover some of the greatest websites! Details here.

Hello Greg,

Thanks for your interest and support. Following your post is a MATLAB
script, on which I design the analog prototype filter (which actually
is the cascade of two sections), see its response, design two
"matching" IIR filters, with BLT and Impulse Invariant methods, and
then I compare them to the full response of the analog prototype. I
decided to go for the BLT. Then I convert into SOS (Second Order
Sections) for analyzing with Filter Design and Analysis tool, basicly
to study the effects of single precision floating point implementation,
which is what I have in my SHARC (21065L)

In my script totH and totW represent the magnitud and phase response of
the prototype filter. I've thought of this approach for FFT based,
frequency domain filtering:

- Take FFT of input signal.
- Multiply it by totH and totW (multiply the FFTs, in other words,
either in Magnitude-Phase representation or Real-Imaginary
representation)
- Tahe IFFT of the above multiplication.

As I understand it, this could be enough for frequency domain
filtering. Questions:

- FFT length requeriments?
- Magnitude-Phase or Real-Imaginary representation?
- Windowing?

I hope to have your comments.

TIA,

JaaC


Greg Berchin wrote:
> On 10 Jan 2005 07:51:15 -0800, "Jaime A. Aranguren C., SanJaaC
> Electronics" <j...@ieee.org> wrote:
>
> >>My sample rate is 240 Hz. The filters is some sort of band pass
filter,
> >>bell shaped. Frequencies of interest are from 0.05 Hz up to 40Hz,
being
> >>the maximum at 8.8 Hz.
>
> Can you post the Laplace-domain transfer function?  Or at least a
> set of frequency/magnitude/phase measurements?
>
> Greg

% Sampling frequency
Fs = 240;
% N points for frequency analysis
N = 8192;

% Analog prototype, first section
firstB = [1.6357*2*pi*9.077169 0];
firstA = [1 2*2*pi*4.167375 (2*pi*9.077169)^2];
[firstH, firstW] = freqs(firstB,firstA, N);
firstMag = 20*log10(abs(firstH));
firstF = firstW/(2*pi);
figure
semilogx(firstF, firstMag,'r');
hold on

% Analog prototype, second section
secondB = [1/(2*pi*2.939902) 1];
secondA = [1/((2*pi*1.394468)*(2*pi*17.31512))
(1/(2*pi*1.394468))+(1/(2*pi*17.3152)) 1];
[secondH, secondW] = freqs(secondB,secondA, N);
secondMag = 20*log10(abs(secondH));
secondF = secondW/(2*pi);
semilogx(secondF, secondMag,'m');

% Full filter
totB = conv(firstB, secondB);
totA = conv(firstA, secondA);
[totH, totW] = freqs(totB,totA, N);
totMag = 20*log10(abs(totH));
totPha = unwrap(angle(totH));
totF = totW/(2*pi);
semilogx(totF, totMag,'b');
title('Lamp model 120V Partial & Full Responses, log scale')
grid on

% Plot analog sections, full analog filter
figure
plot(firstF, firstMag,'r');
hold on
plot(secondF, secondMag,'m');
plot(totF, totMag,'b');
title('Lamp model 120V Partial & Full Responses, lin scale')
grid on

% Design digital filter, bilinear transform method
[bilinB, bilinA] = bilinear(totB, totA, Fs, 8);
bilinH = freqz(bilinB, bilinA, totF/(4*pi^2));
bilinMag = 20*log10(abs(bilinH));
bilinPha = unwrap(angle(bilinH));

% Design digital filter, impulse invariant method
[impinvB, impinvA] = impinvar(totB, totA, Fs);
impinvH = freqz(impinvB, impinvA, totF/(4*pi^2));
impinvMag = 20*log10(abs(impinvH));
impinvPha = unwrap(angle(impinvH));

% Plot, logarithmic axis
figure;
semilogx(totF, totMag, 'b');
hold on;
semilogx(totF, bilinMag, 'k');
hold on;
semilogx(totF, impinvMag, 'r');
title('Log Plot. Analog Filter: Blue. Bilinear Transf: Black. Impulse
Inv: Red')
grid on, zoom on

% Plot, linear axis
figure;
plot(totF, totMag, 'b');
hold on;
plot(totF, bilinMag, 'k');
hold on;
plot(totF, impinvMag, 'r');
title('Lin Plot. Analog Filter: Blue. Bilinear Transf: Black. Impulse
Inv: Red')
grid on, zoom on

% Plot, logarithmic axis
figure;
semilogx(totF, totPha, 'b');
hold on;
semilogx(totF, bilinPha, 'k');
hold on;
semilogx(totF, impinvPha, 'r');
title('Log Plot. Analog Filter: Blue. Bilinear Transf: Black. Impulse
Inv: Red')
grid on, zoom on

% Plot, linear axis
figure;
plot(totF, totPha, 'b');
hold on;
plot(totF, bilinPha, 'k');
hold on;
plot(totF, impinvPha, 'r');
title('Lin Plot. Analog Filter: Blue. Bilinear Transf: Black. Impulse
Inv: Red')
grid on, zoom on


% Error. Bilinear transform
figure;
semilogx(totF, totMag-bilinMag)
title('Log Plot. Error: Bilinear')
grid on, zoom on
figure;
plot(totF, totMag-bilinMag)
grid on, zoom on
title('Lin Plot. Error: Bilinear')

% Error. Impulse invariant
figure;
semilogx(totF, totMag-impinvMag)
title('Log Plot. Error: Impulse Inv.')
grid on, zoom on
figure;
plot(totF, totMag-impinvMag)
grid on, zoom on
title('Lin Plot. Error: Impulse Inv.')

% Error. Bilinear Impulse invariant
figure;
semilogx(totF, totMag-bilinMag, totF, totMag-impinvMag)
title('Log Plot. Error')
grid on, zoom on
figure;
plot(totF, totMag-bilinMag, totF, totMag-impinvMag)
grid on, zoom on
title('Lin Plot. Error')

% Selected filter si bilinear transform. Convert to Second Order
Sections
[modSOS, modG] = tf2sos(bilinB, bilinA, 'DOWN', 2)

______________________________
Are you a member of DSPRelated.com? If you are, you can now access a powerful tool to discover new web sites. Details here.

Hello Clay,

Clay S. Turner wrote:

> With the BLT, there is only one parameter to adjust. We had a thread
a while
> back on the BLT, where I described how to match this up.

Could you please post the URL of the thread?

> However by the time this much effort is
> done, I'd just go and design the best fit (given whatever constraints
are
> important to the problem) IIR via a Fletcher-Powell method.

May I know about information sources for this design method?

Finally, let me tell you, if I already have not done before, that I
find your article "Recursive Discrete-Time Sinusoidal Oscillators" IEEE
Sig Proc Mag, May 2003, to be a real trasure about the topic for the
practicing engineer. Thanks for that jewel!

> Clay

JaaC

______________________________
Discover some of the greatest websites! Details here.

Hello Clay,

Clay S. Turner wrote:

> With the BLT, there is only one parameter to adjust. We had a thread
a while
> back on the BLT, where I described how to match this up.

Could you please post the URL of the thread?

> However by the time this much effort is
> done, I'd just go and design the best fit (given whatever constraints
are
> important to the problem) IIR via a Fletcher-Powell method.

May I know about information sources for this design method?

Finally, let me tell you, if I already have not done before, that I
find your article "Recursive Discrete-Time Sinusoidal Oscillators" IEEE
Sig Proc Mag, May 2003, to be a real trasure about the topic for the
practicing engineer. Thanks for that jewel!

> Clay

JaaC

______________________________
What are the most bookmarked web sites about "linux tutorial"? Find out here.

Hello Clay,

Clay S. Turner wrote:

> With the BLT, there is only one parameter to adjust. We had a thread
a while
> back on the BLT, where I described how to match this up.

Could you please post the URL of the thread?

> However by the time this much effort is
> done, I'd just go and design the best fit (given whatever constraints
are
> important to the problem) IIR via a Fletcher-Powell method.

May I know about information sources for this design method?

Finally, let me tell you, if I already have not done before, that I
find your article "Recursive Discrete-Time Sinusoidal Oscillators" IEEE
Sig Proc Mag, May 2003, to be a real trasure about the topic for the
practicing engineer. Thanks for that jewel!

> Clay

JaaC

______________________________
Are you a member of DSPRelated.com? If you are, you can now access a powerful tool to discover new web sites. Details here.

Hello Clay,

Clay S. Turner wrote:

> With the BLT, there is only one parameter to adjust. We had a thread
a while
> back on the BLT, where I described how to match this up.

Could you please post the URL of the thread?

> However by the time this much effort is
> done, I'd just go and design the best fit (given whatever constraints
are
> important to the problem) IIR via a Fletcher-Powell method.

May I know about information sources for this design method?

Finally, let me tell you, if I already have not done before, that I
find your article "Recursive Discrete-Time Sinusoidal Oscillators" IEEE
Sig Proc Mag, May 2003, to be a real trasure about the topic for the
practicing engineer. Thanks for that jewel!

> Clay

JaaC

______________________________
Discover some of the greatest websites! Details here.

Hello Jaime,
comments below:

"Jaime A. Aranguren C., SanJaaC Electronics" <j...@ieee.org> 
wrote in message 
news:1...@c13g2000cwb.googlegroups.com...
> Hello Clay,
>
> Clay S. Turner wrote:
>
>> With the BLT, there is only one parameter to adjust. We had a thread
> a while
>> back on the BLT, where I described how to match this up.
>
> Could you please post the URL of the thread?

http://www.dsprelated.com/showmessage/1049/1.php


>
>> However by the time this much effort is
>> done, I'd just go and design the best fit (given whatever constraints
> are
>> important to the problem) IIR via a Fletcher-Powell method.
>
> May I know about information sources for this design method?


There is a program by K. Steiglitz of Princeton that uses a F-P algo. It is 
written in Fortran and is described in the book "Digital Signal Processing: 
Theory, Design, and Implementation" by Peled and Liu, (c) 1976, John Wiley & 
Sons,  TK5102.5.P36

This is an old book but it does contain a few gems.


>
> Finally, let me tell you, if I already have not done before, that I
> find your article "Recursive Discrete-Time Sinusoidal Oscillators" IEEE
> Sig Proc Mag, May 2003, to be a real trasure about the topic for the
> practicing engineer. Thanks for that jewel!

Thanks,

I'm glad that you find it useful.


Clay

______________________________
What are the most bookmarked web sites about "fpga"? Find out here.