IIR Hilbert - Matlab help anyone

Hi,
I'm (trying to) implement an IIR Hilbert Transformer.  I found a
reference to a design which provides the z domain poles and zeros as
listed at the bottom of this post.  The function is 12th order and I
figured that I could implement this conveniently with 6 cascaded
biquads (which I happen to have available in my FPGA DSP design).

I calculated out the IIR filter coeffs by taking the conjugate pole
pairs and zero pairs and multiplying them out to provide a function
like z^2 + aZ + b to provide the coefficients.  These are listed below
also:

I then simply transplanted these coefficients into my biquad DSP (b
coeffs do not need to be negated as is often the case).

The response I get is not what I expected.  It is flat in magnitude
with 45dB gain (should be 0dB) but the phase constantly varies instead
of being roughly 90 degrees for most of the band.

I also computed the frequency response in Excel and this does not
appear to give the response expected either.  I assume it's finger
trouble but would someone out there mind having a look at this problem
in Matlab?  I don't have Matlab or Scilab (and not trained up) and I
think it would be relatively straightforward for an experienced person
to check these calculations.

Many thanks if you can spare some time to look at this and please feel
free to call me stupid if I have made an obviously bad assumption.

Zeros (complex)
=========================
    a                	jb
-0.834946		-1.459074
-0.834946		1.459074
0.834946		-1.459074
0.834946		1.459074
0		-1.707483
0		1.707483
-1.362397		-0.79655
-1.362397		0.79655
1.362397		-0.79655
1.362397		0.79655
-1.118507		0
1.118507		0

Poles (complex)
==============
-0.894071		0
0.894071		0
-0.547048		-0.319837
-0.547048		0.319837
0.547048		-0.319837
0.547048		0.319837
-0.295469		-0.516331
-0.295469		0.516331
0.295469		-0.516331
0.295469		0.516331
0		-0.585696
0		0.585696

Coefficients (numerator coeffs first then denominator):
======================================
a0	a1	a2
1	-1.669892	2.82603176
1	-1.669892	2.82603176
1	1.669892	2.82603176
1	0	2.915498195
1	-2.724794	2.490617488
1	2.724794	2.490617488
1	0	-1.251057909

(1)	b1	b2
1	0	-0.799362953
1	-1.094096	0.401557221
1	1.094096	0.401557221
1	-0.590938	0.353899632
1	0.590938	0.353899632
1	0	0.343039804

I made a mistake in my post - the first line of a coefficients was
repeated so it should read:

Coefficients (numerator coeffs first then denominator):
======================================
a0      a1      a2
1       -1.669892       2.82603176
1       1.669892        2.82603176
1       0       2.915498195
1       -2.724794       2.490617488
1       2.724794        2.490617488
1       0       -1.251057909

In the mean time, I'm currently installing Scilab to see if I can
solve this myself.

OK, I managed to install Scilab and plotted the frequency response of
my H(z) and it matches what I see in my Excel model.  So either my
coefficient calculations are wrong or the original poles and zeros
provided are wrong.  I'm more inclined to think that it is me at fault
here!

Does anyone know if it is possible to define a syslin in Scilab (not
sure of the equivalent Matlab language) in terms of complex poles and
zeros rather than real coefficients of z?

If I could do this then I could test the poles and zero's as published
- this would then allow me to narrow the problem down to my
coefficient calculations.

davew wrote:

> In the mean time, I'm currently installing Scilab to see if I can
> solve this myself.

And I am glad to tell the World that I am currently installing the 
PL-2303 USB to serial driver.

VLV

I'm just trying to point out that I'm not asking for help without
trying to help myself!
Good luck with your driver installation.

blah blah

Here is the original paper:

http://wwwtw.vub.ac.be/elec/Papers%20on%20web/Papers/RikPintelon/IM90Kollar-Optimal%20FIR%20an%20IIR.pdf

I created a H(z) in Scilab using the poles and zeros as specified on
page 5 of this document.  I made a pole zero plot which is identical
to that in the document, however the magnitude and phase plots are
nothing like those given.  Am I doing something wrong or is the
document in error?  The author told me that as far as he remembers,
the filter "worked" but he hadn't time to re-visit it.

Can anyone offer an explanation?

Thanks.

On 6 July, 10:46, davew <david.wo...@gmail.com> wrote:
> Here is the original paper:
>
> http://wwwtw.vub.ac.be/elec/Papers%20on%20web/Papers/RikPintelon/IM90...
>
> I created a H(z) in Scilab using the poles and zeros as specified on
> page 5 of this document. &#4294967295;I made a pole zero plot which is identical
> to that in the document, however the magnitude and phase plots are
> nothing like those given. &#4294967295;Am I doing something wrong or is the
> document in error? &#4294967295;The author told me that as far as he remembers,
> the filter "worked" but he hadn't time to re-visit it.
>
> Can anyone offer an explanation?
>
> Thanks.

and here is the Scilab code if anyone wants to try it/investigate
further:

z1 = -0.834946	+ -1.459074 * %i;
z2 = -0.834946	+ 1.459074 * %i;
z3 = 0.834946	+ -1.459074 * %i;
z4 = 0.834946	+ 1.459074 * %i;
z5 = 0	+ -1.707483 * %i;
z6 = 0	+ 1.707483 * %i;
z7 = -1.362397	+ -0.79655 * %i;
z8 = -1.362397	+ 0.79655 * %i;
z9 = 1.362397	+ -0.79655 * %i;
z10 = 1.362397	+ 0.79655 * %i;
z11 = -1.118507	+ 0 * %i;
z12 = 1.118507	+ 0 * %i;

p1 = -0.894071	+ 0 * %i;
p2 = 0.894071	+ 0 * %i;
p3 = -0.547048	+ -0.319837 * %i;
p4 = -0.547048	+ 0.319837 * %i;
p5 = 0.547048	+ -0.319837 * %i;
p6 = 0.547048	+ 0.319837 * %i;
p7 = -0.295469	+ -0.516331 * %i;
p8 = -0.295469	+ 0.516331 * %i;
p9 = 0.295469	+ -0.516331 * %i;
p10 = 0.295469	+ 0.516331 * %i;
p11 = 0	+ -0.585696 * %i;
p12 = 0	+ 0.585696 * %i;

H = (real((%z - z1) * (%z - z2)) * real((%z - z3) * (%z - z4)) * real
((%z - z5) * (%z - z6)) * real((%z - z7) * (%z - z8)) * real((%z - z9)
* (%z - z10)) * real((%z - z11) * (%z - z12))) / (real((%z - p1) * (%z
- p2)) * real((%z - p3) * (%z - p4)) * real((%z - p5) * (%z - p6)) *
real((%z - p7) * (%z - p8)) * real((%z - p9) * (%z - p10)) * real((%z
- p11) * (%z - p12)))

sys=syslin('d', H);

fmin=0.001;
fmax=0.5;

scf(1);clf;
plzr(H);

scf(2);clf;
bode(sys,fmin,fmax); //Plots frequency response in Bode diagram

On 6 July, 11:22, davew <david.wo...@gmail.com> wrote:
> On 6 July, 10:46, davew <david.wo...@gmail.com> wrote:
>
> > Here is the original paper:
>
> >http://wwwtw.vub.ac.be/elec/Papers%20on%20web/Papers/RikPintelon/IM90...
>
> > I created a H(z) in Scilab using the poles and zeros as specified on
> > page 5 of this document. &#4294967295;I made a pole zero plot which is identical
> > to that in the document, however the magnitude and phase plots are
> > nothing like those given. &#4294967295;Am I doing something wrong or is the
> > document in error? &#4294967295;The author told me that as far as he remembers,
> > the filter "worked" but he hadn't time to re-visit it.
>
> > Can anyone offer an explanation?
>
> > Thanks.
>
> and here is the Scilab code if anyone wants to try it/investigate
> further:
>
> z1 = -0.834946 &#4294967295;+ -1.459074 * %i;
> z2 = -0.834946 &#4294967295;+ 1.459074 * %i;
> z3 = 0.834946 &#4294967295; + -1.459074 * %i;
> z4 = 0.834946 &#4294967295; + 1.459074 * %i;
> z5 = 0 &#4294967295;+ -1.707483 * %i;
> z6 = 0 &#4294967295;+ 1.707483 * %i;
> z7 = -1.362397 &#4294967295;+ -0.79655 * %i;
> z8 = -1.362397 &#4294967295;+ 0.79655 * %i;
> z9 = 1.362397 &#4294967295; + -0.79655 * %i;
> z10 = 1.362397 &#4294967295;+ 0.79655 * %i;
> z11 = -1.118507 + 0 * %i;
> z12 = 1.118507 &#4294967295;+ 0 * %i;
>
> p1 = -0.894071 &#4294967295;+ 0 * %i;
> p2 = 0.894071 &#4294967295; + 0 * %i;
> p3 = -0.547048 &#4294967295;+ -0.319837 * %i;
> p4 = -0.547048 &#4294967295;+ 0.319837 * %i;
> p5 = 0.547048 &#4294967295; + -0.319837 * %i;
> p6 = 0.547048 &#4294967295; + 0.319837 * %i;
> p7 = -0.295469 &#4294967295;+ -0.516331 * %i;
> p8 = -0.295469 &#4294967295;+ 0.516331 * %i;
> p9 = 0.295469 &#4294967295; + -0.516331 * %i;
> p10 = 0.295469 &#4294967295;+ 0.516331 * %i;
> p11 = 0 + -0.585696 * %i;
> p12 = 0 + 0.585696 * %i;
>
> H = (real((%z - z1) * (%z - z2)) * real((%z - z3) * (%z - z4)) * real
> ((%z - z5) * (%z - z6)) * real((%z - z7) * (%z - z8)) * real((%z - z9)
> * (%z - z10)) * real((%z - z11) * (%z - z12))) / (real((%z - p1) * (%z
> - p2)) * real((%z - p3) * (%z - p4)) * real((%z - p5) * (%z - p6)) *
> real((%z - p7) * (%z - p8)) * real((%z - p9) * (%z - p10)) * real((%z
> - p11) * (%z - p12)))
>
> sys=syslin('d', H);
>
> fmin=0.001;
> fmax=0.5;
>
> scf(1);clf;
> plzr(H);
>
> scf(2);clf;
> bode(sys,fmin,fmax); //Plots frequency response in Bode diagram

can anyone offer an explanation?