View Raw SPL
/*****************************************************************************
* *
* TF2CAS.SPL Copyright (C) 2015 DSP Development Corporation *
* All Rights Reserved *
* *
* Author: Randy Race *
* *
* Synopsis: Converts direct form filter coefficients to cascade form *
* *
* Revisions: 5 Aug 2015 RRR Creation *
* *
*****************************************************************************/
#if @HELP_TF2CAS
TF2CAS
Purpose:
Converts direct form filter coefficients to cascade form.
Format:
TF2CAS(b, a)
(num, den, gain) = TF2CAS(b, a)
b - a series, numerator (i.e. zero) coefficients
a - a series, denominator (i.e. pole) coefficients
Returns:
A series, the cascade form of the filter coefficients.
(num, den, gain) = TF2CAS(b, a) returns the numerator,
denominator cascade coefficients as separate Nx3 arrays where
each row represents a 2nd order stage. The gain is returned
as a separate scalar.
Example:
b = {1};
a = {1, -0.5, 0.2};
c = tf2cas(b, a);
c == {1, 1, 0, 0, -0.5, 0.2}
The direct form filter coefficients represent the following
Z transform:
1
H(z) = -------------------
-1 -2
1 - 0.5z + 0.2z
The resulting 2nd order cascade coefficients represent the
equivalent Z transform:
-1 -2
1 + 0.0z + 0.0z
H(z) = 1 * -------------------
-1 -2
1 - 0.5z + 0.2z
Since the original coefficients represent a single 2nd order
system, the cascade denominator coefficients are identical.
Example:
b = {1};
a = {1, -0.5, 0.2, 0.1};
c = tf2cas(b, a);
c == {1, 1, 0, 0, -0.754856, 0.392379, 1, 0, 0, 0.254856, 0}
The filter coefficients represent the following Z transform:
1
H(z) = -----------------------------
-1 -2 -3
1 - 0.5z + 0.2z + 0.1z
The resulting 2nd order cascade coefficients represent the
equivalent Z transform:
-1 -2
1 + 0.0z + 0.0z
H(z) = 1 * ----------------------------
-1 -2
1 - 0.754856z + 0.392379z
-1 -2
1 + 0.0z + 0.0z
* ------------------------
-1 -2
1 - 0.254856z + 0.0z
The direct form coefficients are represented by two stages
of 2nd order cascade coefficients.
Example:
b = {1};
a = {1, -0.5, 0.2, 0.1};
(num, den, gain) = tf2cas(b, a);
num = {{1, 0, 0},
{1, 0, 0}}
den = {{1, -0.754856, 0.392379},
{1, 0.254856, 0.0}}
gain == 1
Same as the previous example except the cascade coefficients
are returned as separate 2x3 arrays where each row represents
the coefficients of a 2nd order stage.
Remarks:
TF2CAS converts direct form filter coefficients to cascaded
2nd order bi-quad form. The direct form coefficients represent
the following Z transform:
-1 -Q
Y(z) b[1] + b[2]*z + ... + b[Q+1]*z
H(z) = ---- = -------------------------------------
-1 -P
X(z) 1 + a[2]*z + ... + a[P+1]*z
jw
where z = e unit circle (frequency response)
Q order of zeros (numerator)
P order of poles (denominator)
The cascade form of the coefficients represent the equivalent
Z transform:
-1 -2 -1 -2
b10 + b11 z + b12 z b20 + b21 z + b22 z
H(z) = G * ________________________ * ________________________ ...
-1 -2 -1 -2
1 + a11 z + a12 z 1 + a21 z + a22 z
The cascade filter coefficients are returned as a single column
series with the coefficients in the following order:
G, b10, b11, b12, a11, a12, b20, b21, b22, a21, a22, ...
The cascade form of the coefficients can be processed by
the CASCADE functions. Filtering with cascade form
coefficients can reduce numerical roundoff errors compared
to the same filtering operation in direct form. However,
because numeric instabilities may be present in the direct
form coefficients, it is best to design the filter using
the cascade form directly, without conversion. This is
the format produced by the IIR filter design routines
provided by DADiSP/Filters.
See CAS2TF to convert cascade form coefficients to direct form.
See Also:
Cas2tf
Cas2zp
Zp2cas
Cascade
Filteq
#endif
/* convert direct form B(z)/(A(z) to cascaded biquad coefficients */
tf2cas(b, a, Fs)
{
local c, g, b1, a1, z, p, k, len;
if (outargc < 2)
{
echo("Converting Direct Form to Cascaded Biquad Form ...");
}
/* check input args */
(b, a, Fs) = tf2cas_args(argc, b, a, Fs);
if (length(a) == 1)
{
/* fir case */
(z, p, k) = tf2zpk(b, a);
}
else
{
/* same size */
len = max(length(b), length(a));
b = extract(b, 1, len);
a = extract(a, 1, len);
/* zeros, poles and gain */
(z, p, k) = tf2zp(b, a);
}
/* cascade form */
c = zp2cas(z, p, k);
if (outargc > 1)
{
g = c[1];
c = extract(c, 2, -1);
c = ravel(c, 5);
b1 = unravel(extract(c, 1, 3));
a1 = unravel(extract(c, 4, 2));
setdeltax(b1, 1/Fs);
setdeltax(a1, 1/Fs);
setxoffset(b1, xoffset(b));
setxoffset(a1, xoffset(a));
return(b1, a1, g);
}
else
{
/*
* standard DADiSP cascade format:
*
* G, b10, b11, b12, a11, a12, b20, b21, b22, a21, a22, ...
*/
setdeltax(c, 1/Fs);
setxoffset(c, xoffset(b));
return(c);
}
}
/* get sample rate from coefficients */
tf2cas_rate(b, a)
{
local cas_rate;
if (rate(b) == rate(a))
{
cas_rate = rate(b);
}
else if (length(b) == 1)
{
cas_rate = rate(a);
}
else if (length(a) == 1)
{
cas_rate = rate(b);
}
else
{
cas_rate = 1.0;
}
return(cas_rate);
}
/* parse tf2cas input arguments */
tf2cas_args(cnt, b, a, Fs)
{
/* parse args */
if (cnt < 3)
{
if (cnt < 2)
{
if (cnt < 1)
{
tf2cas_args_err();
}
if (not(isarray(b)))
{
tf2cas_args_err();
}
a = {1};
}
if (not(isarray(a)))
{
Fs = a;
a = {1};
}
if (numcols(b) == 2)
{
/* B(z), A(z) as columns */
a = col(b, 2);
b = col(b, 1);
}
Fs = tf2cas_rate(b, a);
}
if (not(isarray(b)) || not(isarray(a)) || not(isnumeric(Fs)))
{
tf2cas_args_err();
}
if (length(b) == 0 || length(a) == 0)
{
tf2cas_args_err();
}
return(b, a, Fs);
}
/* error function for tf2cas */
tf2cas_args_err()
{
error("tf2cas - B and A Coefficient Series Required");
}