libdspl-2.0
Digital Signal Processing Algorithm Library
IIR filters design.

Functions

int DSPL_API butter_ap (double rp, int ord, double *b, double *a)
 Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Butterworth filter. More...
 
int DSPL_API butter_ap_zp (int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
 Function calculates arrays of zeros and poles for analog normlized lowpass Batterworth filter transfer function \( H(s) \) order ord . More...
 
int DSPL_API cheby1_ap (double rp, int ord, double *b, double *a)
 Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Chebyshev type 1 filter. More...
 
int DSPL_API cheby1_ap_zp (int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
 Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 1 filter transfer function \( H(s) \) order ord . More...
 
int DSPL_API cheby2_ap (double rs, int ord, double *b, double *a)
 Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Chebyshev type 2 filter. More...
 
int DSPL_API cheby2_ap_zp (int ord, double rs, complex_t *z, int *nz, complex_t *p, int *np)
 Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 2 filter transfer function \( H(s) \) order ord . More...
 
int DSPL_API ellip_ap (double rp, double rs, int ord, double *b, double *a)
 Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass elliptic filter order ord with passband ripple rp dB and stopband suppression equals rs dB. More...
 
int DSPL_API ellip_ap_zp (int ord, double rp, double rs, complex_t *z, int *nz, complex_t *p, int *np)
 Function calculates arrays of zeros and poles for analog normlized lowpass elliptic filter transfer function \( H(s) \) order ord . More...
 
int DSPL_API low2high (double *b, double *a, int ord, double w0, double w1, double *beta, double *alpha)
 Lowpass to highpass filter frequency transform. More...
 
int DSPL_API low2low (double *b, double *a, int ord, double w0, double w1, double *beta, double *alpha)
 
int DSPL_API ratcompos (double *b, double *a, int n, double *c, double *d, int p, double *beta, double *alpha)
 Rational composition. More...
 
int DSPL_API bilinear (double *bs, double *as, int ord, double *bz, double *az)
 Transform a s-plane analog filter transfer function \(H(s)\) to the digital filter transfer function \(H(z)\). More...
 
int DSPL_API iir (double rp, double rs, int ord, double w0, double w1, int type, double *b, double *a)
 Digital IIR filter design. More...
 

Detailed Description

Function Documentation

◆ bilinear()

int bilinear ( double *  bs,
double *  as,
int  ord,
double *  bz,
double *  az 
)

Transform a s-plane analog filter transfer function \(H(s)\) to the digital filter transfer function \(H(z)\).


Bilinear transform is rational composition:

\[ s \leftarrow \frac{1 - z^{-1}}{1 - z^{-1}}. \]

Digital filter order, passband magnitude ripple and stopband suppression still the same after bilinear transform as analog filter.

Frequency \(\Omega\) of analog filter and frequency
\(\omega\) of digital filter relations:

\[ \Omega = \tan(\omega / 2). \]

Parameters
[in]bsPointer to the vector of analog filter \(H(s)\) numerator coefficients. Vector size is [ord+1 x 1].

[in]asPointer to the vector of analog filter \(H(s)\) denominator coefficients vector. Vector size is [ord+1 x 1].

[in]ordAnalog and digital filters order.

[out]bzPointer to the vector of digital filter \(H(z)\) numerator coefficients after bilinear transform. Vector size is [ord+1 x 1].
Memory must be allocated.

[out]azPointer to the vector of digital filter \(H(z)\) denominator coefficients after bilinear transform. Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if bilinear transform is calculated successfully.
Else code error.

Example:

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
#define N 1000
#define ORD 4
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
double w[N], h[N];
complex_t hz[N];
double bs[ORD+1], as[ORD+1];
double bz[ORD+1], az[ORD+1];
int err, k;
hdspl = dspl_load(); /* Load DSPL function */
/* normalized analog lowpass filter Chebyshev type 1 */
err = cheby1_ap(1.0, ORD, bs, as);
if(err != RES_OK)
{
printf("cheby1_ap error code %d\n", err);
return err;
}
/* Bilinear transform */
err = bilinear(bs, as, ORD, bz, az);
if(err != RES_OK)
{
printf("bilinear error code %d\n", err);
return err;
}
/* Print coefficients */
for(k = 0; k < ORD+1; k++)
printf("bz[%d] = %7.3f az[%d] = %7.3f\n", k, bz[k], k, az[k]);
/* Digital filter magnitude */
linspace(0, M_PI, N, DSPL_PERIODIC, w);
freqz(bz, az, ORD, w, N, hz);
for(k = 0; k < N; k++)
{
h[k] = 10.0 * log10 (ABSSQR(hz[k])); /* Logarithmic scale */
w[k] /= M_PI; /* frequency from 0 to 1 */
}
writetxt(w,h,N,"dat/bilinear.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 560, 380, "img/bilinear.png", &hplot);
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set xlabel 'normalized frequency'");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-80:5]");
gnuplot_cmd(hplot, "plot 'dat/bilinear.txt' with lines");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return err;
}

This program calculates the transfer function \(H(s)\) of analog Chebyshev filter of the first kind, with a cutoff frequency of 1 rad/s, and produces bilinear trandform to digital filter, with a normilized cutoff frequency equals 0.5.

Result:

bz[0] =     0.246        az[0] =     4.425
bz[1] =     0.983        az[1] =    -3.318
bz[2] =     1.474        az[2] =     4.746
bz[3] =     0.983        az[3] =    -2.477
bz[4] =     0.246        az[4] =     1.034
err = 0

In addition, the frequency response of the resulting digital filter is calculated and plotted by GNUPLOT package.

Author
Sergey Bakhurin www.dsplib.org

Definition at line 210 of file filter_iir.c.

Referenced by iir().

◆ butter_ap()

int butter_ap ( double  Rp,
int  ord,
double *  b,
double *  a 
)

Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Butterworth filter.


Analog normalized lowpass filter magnitude ripple equals \( -R_p \) dB for angular frequency \( \omega \) from 0 to 1 rad/s.

Parameters
[in]RpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[in]ordFilter order.
Filter coefficients number equals ord+1 for numerator and denominator of transfer function \( H(s) \)

[out]bPointer to the vector of transfer function \(H(s)\) numerator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

[out]aPointer to the vector of transfer function \(H(s)\) denominator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Example:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 3
/* Frequency response vector size */
#define N 1000
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
/* Load DSPL functions */
hdspl = dspl_load();
double a[ORD+1]; /* H(s) numerator coefficients vector */
double b[ORD+1]; /* H(s) denominator coefficients vector */
double Rp = 1.0; /* Magnitude ripple from 0 to 1 rad/s */
double w[N]; /* Angular frequency (rad/s) */
double mag[N]; /* Filter Magnitude (dB) */
double phi[N]; /* Phase response */
double tau[N]; /* Group delay */
int k;
/* H(s) coefficients calculation */
int res = butter_ap(Rp, ORD, b, a);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* Print H(s) coefficients */
for(k = 0; k < ORD+1; k++)
printf("b[%2d] = %9.3f a[%2d] = %9.3f\n", k, b[k], k, a[k]);
/* Frequency in logarithmic scale from 0.01 to 100 rad/s */
logspace(-2.0, 2.0, N , DSPL_SYMMETRIC, w);
/* Filter frequency parameter calculation */
filter_freq_resp(b, a, ORD, w, N,
DSPL_FLAG_LOGMAG|DSPL_FLAG_UNWRAP|DSPL_FLAG_ANALOG,
mag, phi, tau);
/* Write Magnitude, phase response and group delay to the files */
writetxt(w, mag, N, "dat/butter_ap_test_mag.txt");
writetxt(w, phi, N, "dat/butter_ap_test_phi.txt");
writetxt(w, tau, N, "dat/butter_ap_test_tau.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 920, 260, "img/butter_ap_test.png", &hplot);
gnuplot_cmd(hplot, "set logscale x");
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "set xlabel 'frequency, rad/s'");
gnuplot_cmd(hplot, "set multiplot layout 1,3 rowsfirst");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-100:5]");
gnuplot_cmd(hplot, "plot 'dat/butter_ap_test_mag.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Phase response, rad'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/butter_ap_test_phi.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Groupdelay, sec'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/butter_ap_test_tau.txt' with lines");
gnuplot_cmd(hplot, "unset multiplot");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

b[ 0] =     1.965     a[ 0] =     1.965
b[ 1] =     0.000     a[ 1] =     3.138
b[ 2] =     0.000     a[ 2] =     2.505
b[ 3] =     0.000     a[ 3] =     1.000 


In dat folder will be created 3 files:

butter_ap_test_mag.txt    magnitude
butter_ap_test_phi.txt    phase response
butter_ap_test_tau.txt    group delay

In addition, GNUPLOT will build the following graphs from data stored in files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 175 of file filter_ap.c.

◆ butter_ap_zp()

int butter_ap_zp ( int  ord,
double  rp,
complex_t z,
int *  nz,
complex_t p,
int *  np 
)

Function calculates arrays of zeros and poles for analog normlized lowpass Batterworth filter transfer function \( H(s) \) order ord .


Analog normalized lowpass filter magnitude ripple equals \( -R_p \) dB for angular frequency \( \omega \) from 0 to 1 rad/s.

Parameters
[in]ordFilter order.
Number of zeros and poles of filter can be less or equal ord.

[in]rpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[out]zPointer to the \( H(s) \) zeros array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]nzPointer to the variable which keep number of finite zeros \( H(s) \).
Number of finite zeros which was calculated and saved in vector z.
Pointer cannot be NULL.

[out]pPointer to the \( H(s) \) poles array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]npPointer to the variable which keep number of calculated poles of \( H(s) \).
Pointer cannot be NULL.

Returns
RES_OK if zeros and poles is calculated successfully.
Else code error.
Note
Normalized Butterworth lowpass filter has no finite zeros. So z vector will not changed and in pointer nz will write 0 value.

Example of normalized Butterworth lowpass filter zeros and poles calculation:

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 7
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
/* Load DSPL functions */
hdspl = dspl_load();
complex_t z[ORD+1]; /* H(s) zeros vector */
complex_t p[ORD+1]; /* H(s) poles vector */
double Rp = 1.0; /* Magnitude ripple from 0 to 1 rad/s */
int res, k, nz, np;
/* Zeros and poles vectors calculation */
res = butter_ap_zp(ORD, Rp, z, &nz, p, &np);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* print H(s) zeros values */
printf("Butterworth filter zeros: %d\n", nz);
for(k = 0; k < nz; k++)
printf("z[%2d] = %9.3f %+9.3f j\n", k, RE(z[k]), IM(z[k]));
/* print H(s) poles values */
printf("Butterworth filter poles: %d\n", np);
for(k = 0; k < np; k++)
printf("p[%2d] = %9.3f %+9.3f j\n", k, RE(p[k]), IM(p[k]));
/* Write complex poles to the file */
writetxt_cmplx(p, ORD, "dat/butter_ap_zp.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 440, 360, "img/butter_ap_zp_test.png", &hplot);
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set xzeroaxis");
gnuplot_cmd(hplot, "set yzeroaxis");
gnuplot_cmd(hplot, "set xlabel 'sigma'");
gnuplot_cmd(hplot, "set ylabel 'jw'");
gnuplot_cmd(hplot, "set size square");
gnuplot_cmd(hplot, "set xrange [-1.5:1.5]");
gnuplot_cmd(hplot, "set yrange [-1.5:1.5]");
gnuplot_cmd(hplot, "plot 'dat/butter_ap_zp.txt' with points pt 2");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

Butterworth filter zeros: 0
Butterworth filter poles: 7
p[ 0] =    -1.101    +0.000 j
p[ 1] =    -0.245    +1.074 j
p[ 2] =    -0.245    -1.074 j
p[ 3] =    -0.687    +0.861 j
p[ 4] =    -0.687    -0.861 j
p[ 5] =    -0.992    +0.478 j
p[ 6] =    -0.992    -0.478 j  


In dat folder will be created butter_ap_zp.txt file.
In addition, GNUPLOT will build the following graphs from data stored in dat/butter_ap_zp.txt file:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 402 of file filter_ap.c.

Referenced by butter_ap().

◆ cheby1_ap()

int cheby1_ap ( double  Rp,
int  ord,
double *  b,
double *  a 
)

Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Chebyshev type 1 filter.


Analog normalized lowpass filter magnitude ripple equals \( -R_p \) dB for angular frequency \( \omega \) from 0 to 1 rad/s.

Parameters
[in]RpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[in]ordFilter order.
Filter coefficients number equals ord+1 for numerator and denominator of transfer function \( H(s) \)

[out]bPointer to the vector of transfer function \(H(s)\) numerator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

[out]aPointer to the vector of transfer function \(H(s)\) denominator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Example:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 4
/* Frequency response vector size */
#define N 1000
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
hdspl = dspl_load(); /* Load DSPL function */
double a[ORD+1]; /* H(s) numerator coefficients vector */
double b[ORD+1]; /* H(s) denominator coefficients vector */
double Rp = 1.0; /* Magnitude ripple from 0 to 1 rad/s */
double w[N]; /* Angular frequency (rad/s) */
double mag[N]; /* Filter Magnitude (dB) */
double phi[N]; /* Phase response */
double tau[N]; /* Group delay */
int k;
/* H(s) coefficients calculation */
int res = cheby1_ap(Rp, ORD, b, a);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* Print H(s) coefficients */
for(k = 0; k < ORD+1; k++)
printf("b[%2d] = %9.3f a[%2d] = %9.3f\n", k, b[k], k, a[k]);
/* Frequency in logarithmic scale from 0.01 to 100 rad/s */
logspace(-2.0, 2.0, N , DSPL_SYMMETRIC, w);
/* Filter frequency parameter calculation */
filter_freq_resp(b, a, ORD, w, N,
DSPL_FLAG_LOGMAG|DSPL_FLAG_UNWRAP | DSPL_FLAG_ANALOG,
mag, phi, tau);
/* Write Magnitude, phase response and group delay to the files */
writetxt(w, mag, N, "dat/cheby1_ap_test_mag.txt");
writetxt(w, phi, N, "dat/cheby1_ap_test_phi.txt");
writetxt(w, tau, N, "dat/cheby1_ap_test_tau.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 920, 260, "img/cheby1_ap_test.png", &hplot);
gnuplot_cmd(hplot, "set logscale x");
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "set xlabel 'frequency, rad/s'");
gnuplot_cmd(hplot, "set multiplot layout 1,3 rowsfirst");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-100:5]");
gnuplot_cmd(hplot, "plot 'dat/cheby1_ap_test_mag.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Phase response, rad'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/cheby1_ap_test_phi.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Groupdelay, sec'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/cheby1_ap_test_tau.txt' with lines");
gnuplot_cmd(hplot, "unset multiplot");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

b[ 0] =     0.125     a[ 0] =     0.177
b[ 1] =     0.000     a[ 1] =     0.405
b[ 2] =     0.000     a[ 2] =     1.169
b[ 3] =     0.000     a[ 3] =     0.582
b[ 4] =     0.000     a[ 4] =     1.000


In dat folder will be created 3 files:

cheby1_ap_test_mag.txt    magnitude
cheby1_ap_test_phi.txt    phase response
cheby1_ap_test_tau.txt    group delay

In addition, GNUPLOT will build the following graphs from data stored in files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 601 of file filter_ap.c.

◆ cheby1_ap_zp()

int cheby1_ap_zp ( int  ord,
double  rp,
complex_t z,
int *  nz,
complex_t p,
int *  np 
)

Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 1 filter transfer function \( H(s) \) order ord .


Analog normalized lowpass filter magnitude ripple equals \( -R_p \) dB for angular frequency \( \omega \) from 0 to 1 rad/s.

Parameters
[in]ordFilter order.
Number of zeros and poles of filter can be less or equal ord.

[in]rpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[out]zPointer to the \( H(s) \) zeros array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]nzPointer to the variable which keep number of finite zeros \( H(s) \).
Number of finite zeros which was calculated and saved in vector z.
Pointer cannot be NULL.

[out]pPointer to the \( H(s) \) poles array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]npPointer to the variable which keep number of calculated poles of \( H(s) \).
Pointer cannot be NULL.

Returns
RES_OK if zeros and poles is calculated successfully.
Else code error.
Note
Normalized Chebyshev type 1 lowpass filter has no finite zeros. So z vector will not changed and in pointer nz will write 0 value.
Example of normalized Chebyshev type 1 lowpass filter zeros and poles calculation:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 7
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
/* Load DSPL functions */
hdspl = dspl_load();
complex_t z[ORD+1]; /* H(s) zeros vector */
complex_t p[ORD+1]; /* H(s) poles vector */
double Rp = 0.5; /* Magnitude ripple from 0 to 1 rad/s */
int res, k, nz, np;
/* Zeros and poles vectors calculation */
res = cheby1_ap_zp(ORD, Rp, z, &nz, p, &np);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* print H(s) zeros values */
printf("Chebyshev type 1 filter zeros: %d\n", nz);
for(k = 0; k < nz; k++)
printf("z[%2d] = %9.3f %+9.3f j\n", k, RE(z[k]), IM(z[k]));
/* print H(s) poles values */
printf("Chebyshev type 1 filter poles: %d\n", np);
for(k = 0; k < np; k++)
printf("p[%2d] = %9.3f %+9.3f j\n", k, RE(p[k]), IM(p[k]));
/* Write complex poles to the file */
writetxt_cmplx(p, ORD, "dat/cheby1_ap_zp.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 440, 360, "img/cheby1_ap_zp_test.png", &hplot);
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set xzeroaxis");
gnuplot_cmd(hplot, "set yzeroaxis");
gnuplot_cmd(hplot, "set xlabel 'sigma'");
gnuplot_cmd(hplot, "set ylabel 'jw'");
gnuplot_cmd(hplot, "set size square");
gnuplot_cmd(hplot, "set xrange [-1.5:1.5]");
gnuplot_cmd(hplot, "set yrange [-1.5:1.5]");
gnuplot_cmd(hplot, "plot 'dat/cheby1_ap_zp.txt' with points pt 2");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

Chebyshev type 1 filter zeros: 0
Chebyshev type 1 filter poles: 7
p[ 0] =    -0.256    +0.000 j
p[ 1] =    -0.057    +1.006 j
p[ 2] =    -0.057    -1.006 j
p[ 3] =    -0.160    +0.807 j
p[ 4] =    -0.160    -0.807 j
p[ 5] =    -0.231    +0.448 j
p[ 6] =    -0.231    -0.448 j 


In dat folder will be created cheby1_ap_zp.txt file.
In addition, GNUPLOT will build the following graphs from data stored in dat/cheby1_ap_zp.txt file:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 832 of file filter_ap.c.

Referenced by cheby1_ap().

◆ cheby2_ap()

int cheby2_ap ( double  Rs,
int  ord,
double *  b,
double *  a 
)

Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass Chebyshev type 2 filter.


Analog normalized Chebyshev type 2 filter lowpass filter has \(Rs\) dB suppression in stopband. Also analog normalized Chebyshev type 2 filter magnitude equals \(-Rs\) dB for angular frequency \(\omega = 1\) rad/s.

Parameters
[in]RsSuppression level in stopband (dB).
This parameter sets filter supression for \(\omega \geq 1\) rad/s frequency.
Parameter must be positive.

[in]ordFilter order.
Filter coefficients number equals ord+1 for numerator and denominator of transfer function \( H(s) \)

[out]bPointer to the vector of transfer function \(H(s)\) numerator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

[out]aPointer to the vector of transfer function \(H(s)\) denominator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Example:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 4
/* Frequency response vector size */
#define N 1000
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
double a[ORD+1]; /* H(s) numerator coefficients vector */
double b[ORD+1]; /* H(s) denominator coefficients vector */
double Rs = 60.0; /* Stopband suppression (dB) */
double w[N]; /* Angular frequency (rad/s) */
double mag[N]; /* Filter Magnitude (dB) */
double phi[N]; /* Phase response */
double tau[N]; /* Group delay */
int k;
/* Load DSPL function */
hdspl = dspl_load();
if(!hdspl)
{
printf("libdspl loading error!\n");
return -1;
}
/* H(s) coefficients calculation */
int res = cheby2_ap(Rs, ORD, b, a);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* Print H(s) coefficients */
for(k = 0; k < ORD+1; k++)
printf("b[%2d] = %9.3f a[%2d] = %9.3f\n", k, b[k], k, a[k]);
/* Frequency in logarithmic scale from 0.01 to 100 rad/s */
logspace(-2.0, 2.0, N , DSPL_SYMMETRIC, w);
/* Filter frequency parameter calculation */
filter_freq_resp(b, a, ORD, w, N,
DSPL_FLAG_LOGMAG|DSPL_FLAG_UNWRAP | DSPL_FLAG_ANALOG,
mag, phi, tau);
/* Write Magnitude, phase response and group delay to the files */
writetxt(w, mag, N, "dat/cheby2_ap_test_mag.txt");
writetxt(w, phi, N, "dat/cheby2_ap_test_phi.txt");
writetxt(w, tau, N, "dat/cheby2_ap_test_tau.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 920, 260, "img/cheby2_ap_test.png", &hplot);
gnuplot_cmd(hplot, "set logscale x");
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "set xlabel 'frequency, rad/s'");
gnuplot_cmd(hplot, "set multiplot layout 1,3 rowsfirst");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-100:5]");
gnuplot_cmd(hplot, "plot 'dat/cheby2_ap_test_mag.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Phase response, rad'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/cheby2_ap_test_phi.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Groupdelay, sec'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/cheby2_ap_test_tau.txt' with lines");
gnuplot_cmd(hplot, "unset multiplot");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

b[ 0] =     0.008     a[ 0] =     0.008
b[ 1] =     0.000     a[ 1] =     0.068
b[ 2] =     0.008     a[ 2] =     0.300
b[ 3] =     0.000     a[ 3] =     0.774
b[ 4] =     0.001     a[ 4] =     1.000 


In dat folder will be created 3 files:

cheby2_ap_test_mag.txt    magnitude
cheby2_ap_test_phi.txt    phase response
cheby2_ap_test_tau.txt    group delay

In addition, GNUPLOT will build the following graphs from data stored in files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 1037 of file filter_ap.c.

◆ cheby2_ap_zp()

int cheby2_ap_zp ( int  ord,
double  rs,
complex_t z,
int *  nz,
complex_t p,
int *  np 
)

Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 2 filter transfer function \( H(s) \) order ord .


Analog normalized Chebyshev type 2 filter lowpass filter has \(Rs\) dB suppression in stopband. Also analog normalized Chebyshev type 2 filter magnitude equals \(-Rs\) dB for angular frequency \(\omega = 1\) rad/s.

Parameters
[in]ordFilter order.
Number of zeros and poles of filter can be less or equal ord.

[in]rsSuppression level in stopband (dB).
This parameter sets filter supression for \(\omega \geq 1\) rad/s frequency.
Parameter must be positive.

[out]zPointer to the \( H(s) \) zeros array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]nzPointer to the variable which keep number of finite zeros \( H(s) \).
Number of finite zeros which was calculated and saved in vector z.
Pointer cannot be NULL.

[out]pPointer to the \( H(s) \) poles array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]npPointer to the variable which keep number of calculated poles of \( H(s) \).
Pointer cannot be NULL.

Returns
RES_OK if zeros and poles is calculated successfully.
Else code error.
Example of normalized Chebyshev type 2 lowpass filter zeros and poles calculation:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 7
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
/* Load DSPL functions */
hdspl = dspl_load();
complex_t z[ORD+1]; /* H(s) zeros vector */
complex_t p[ORD+1]; /* H(s) poles vector */
double Rs = 40; /* Stopband suppression, dB */
int res, k, nz, np;
/* Zeros and poles vectors calculation */
res = cheby2_ap_zp(ORD, Rs, z, &nz, p, &np);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* print H(s) zeros values */
printf("Chebyshev type 2 filter zeros: %d\n", nz);
for(k = 0; k < nz; k++)
printf("z[%2d] = %9.3f %+9.3f j\n", k, RE(z[k]), IM(z[k]));
/* print H(s) poles values */
printf("Chebyshev type 2 filter poles: %d\n", np);
for(k = 0; k < np; k++)
printf("p[%2d] = %9.3f %+9.3f j\n", k, RE(p[k]), IM(p[k]));
/* Write complex poles to the file */
writetxt_cmplx(z, nz, "dat/cheby2_ap_z.txt");
writetxt_cmplx(p, np, "dat/cheby2_ap_p.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 440, 360, "img/cheby2_ap_zp_test.png", &hplot);
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set xzeroaxis");
gnuplot_cmd(hplot, "set yzeroaxis");
gnuplot_cmd(hplot, "set xlabel 'sigma'");
gnuplot_cmd(hplot, "set ylabel 'jw'");
gnuplot_cmd(hplot, "set size square");
gnuplot_cmd(hplot, "set xrange [-3:3]");
gnuplot_cmd(hplot, "set yrange [-3:3]");
gnuplot_cmd(hplot, "plot 'dat/cheby2_ap_p.txt' with points pt 2, \\");
gnuplot_cmd(hplot, " 'dat/cheby2_ap_z.txt' with points pt 6");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

Chebyshev type 2 filter zeros: 6
z[ 0] =     0.000    +1.026 j
z[ 1] =     0.000    -1.026 j
z[ 2] =     0.000    +1.279 j
z[ 3] =     0.000    -1.279 j
z[ 4] =     0.000    +2.305 j
z[ 5] =     0.000    -2.305 j
Chebyshev type 2 filter poles: 7
p[ 0] =    -1.203    +0.000 j
p[ 1] =    -0.113    +0.772 j
p[ 2] =    -0.113    -0.772 j
p[ 3] =    -0.398    +0.781 j
p[ 4] =    -0.398    -0.781 j
p[ 5] =    -0.852    +0.642 j
p[ 6] =    -0.852    -0.642 j 


In dat folder will be created cheby2_ap_z.txt and cheby2_ap_z.txt files which keeps zeros and poles vectors.
In addition, GNUPLOT will build the following graphs from data stored in the files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 1304 of file filter_ap.c.

Referenced by cheby2_ap().

◆ ellip_ap()

int ellip_ap ( double  rp,
double  rs,
int  ord,
double *  b,
double *  a 
)

Function calculates the transfer function \( H(s) \) coefficients of analog normalized lowpass elliptic filter order ord with passband ripple rp dB and stopband suppression equals rs dB.


Parameters
[in]rpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[in]rsSuppression level in stopband (dB).
This parameter sets filter supression for \(\omega \geq 1\) rad/s frequency.
Parameter must be positive.

[in]ordFilter order.
Filter coefficients number equals ord+1 for numerator and denominator of transfer function \( H(s) \)

[out]bPointer to the vector of transfer function \(H(s)\) numerator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

[out]aPointer to the vector of transfer function \(H(s)\) denominator coefficient.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Example:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 4
/* Frequency response vector size */
#define N 1000
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
hdspl = dspl_load(); /* Load DSPL function */
double a[ORD+1]; /* H(s) numerator coefficients vector */
double b[ORD+1]; /* H(s) denominator coefficients vector */
double Rp = 1.0; /* Magnitude ripple from 0 to 1 rad/s */
double Rs = 60.0;/* Stopband suppression (dB) */
double w[N]; /* Angular frequency (rad/s) */
double mag[N]; /* Filter Magnitude (dB) */
double phi[N]; /* Phase response */
double tau[N]; /* Group delay */
int k;
/* H(s) coefficients calculation */
int res = ellip_ap(Rp, Rs, ORD, b, a);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* Print H(s) coefficients */
for(k = 0; k < ORD+1; k++)
printf("b[%2d] = %9.3f a[%2d] = %9.3f\n", k, b[k], k, a[k]);
/* Frequency in logarithmic scale from 0.01 to 100 rad/s */
logspace(-2.0, 2.0, N , DSPL_SYMMETRIC, w);
/* Filter frequency parameter calculation */
filter_freq_resp(b, a, ORD, w, N,
DSPL_FLAG_LOGMAG|DSPL_FLAG_UNWRAP | DSPL_FLAG_ANALOG,
mag, phi, tau);
/* Write Magnitude, phase response and group delay to the files */
writetxt(w, mag, N, "dat/ellip_ap_test_mag.txt");
writetxt(w, phi, N, "dat/ellip_ap_test_phi.txt");
writetxt(w, tau, N, "dat/ellip_ap_test_tau.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 920, 260, "img/ellip_ap_test.png", &hplot);
gnuplot_cmd(hplot, "set logscale x");
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "set xlabel 'frequency, rad/s'");
gnuplot_cmd(hplot, "set multiplot layout 1,3 rowsfirst");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-100:5]");
gnuplot_cmd(hplot, "plot 'dat/ellip_ap_test_mag.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Phase response, rad'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/ellip_ap_test_phi.txt' with lines");
gnuplot_cmd(hplot, "set ylabel 'Groupdelay, sec'");
gnuplot_cmd(hplot, "unset yrange");
gnuplot_cmd(hplot, "plot 'dat/ellip_ap_test_tau.txt' with lines");
gnuplot_cmd(hplot, "unset multiplot");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

b[ 0] =     0.268     a[ 0] =     0.301
b[ 1] =     0.000     a[ 1] =     0.764
b[ 2] =     0.045     a[ 2] =     1.472
b[ 3] =     0.000     a[ 3] =     0.948
b[ 4] =     0.001     a[ 4] =     1.000


In dat folder will be created 3 files:

ellip_ap_test_mag.txt    magnitude
ellip_ap_test_phi.txt    phase response
ellip_ap_test_tau.txt    group delay

In addition, GNUPLOT will build the following graphs from data stored in files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 1538 of file filter_ap.c.

◆ ellip_ap_zp()

int ellip_ap_zp ( int  ord,
double  rp,
double  rs,
complex_t z,
int *  nz,
complex_t p,
int *  np 
)

Function calculates arrays of zeros and poles for analog normlized lowpass elliptic filter transfer function \( H(s) \) order ord .


Parameters
[in]ordFilter order.
Number of zeros and poles of filter can be less or equal ord.

[in]rpMagnitude ripple in passband (dB).
This parameter sets maximum filter distortion from 0 to 1 rad/s frequency.
Parameter must be positive.

[in]rsSuppression level in stopband (dB).
This parameter sets filter suppression for \(\omega \geq 1\) rad/s frequency.
Parameter must be positive.

[out]zPointer to the \( H(s) \) zeros array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]nzPointer to the variable which keep number of finite zeros \( H(s) \).
Number of finite zeros which was calculated and saved in vector z.
Pointer cannot be NULL.

[out]pPointer to the \( H(s) \) poles array.
Maximum vector size is [ord x 1].
Memory must be allocated for maximum vector size.

[out]npPointer to the variable which keep number of calculated poles of \( H(s) \).
Pointer cannot be NULL.

Returns
RES_OK if zeros and poles is calculated successfully.
Else code error.
Example of normalized elliptic lowpass filter zeros and poles calculation:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Filter order */
#define ORD 7
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
/* Load DSPL functions */
hdspl = dspl_load();
complex_t z[ORD+1]; /* H(s) zeros vector */
complex_t p[ORD+1]; /* H(s) poles vector */
double Rp = 1.0; /* Magnitude ripple from 0 to 1 rad/s */
double Rs = 40; /* Stopband suppression, dB */
int res, k, nz, np;
/* Zeros and poles vectors calculation */
res = ellip_ap_zp(ORD, Rp, Rs, z, &nz, p, &np);
if(res != RES_OK)
printf("error code = 0x%8x\n", res);
/* print H(s) zeros values */
printf("Elliptic filter zeros: %d\n", nz);
for(k = 0; k < nz; k++)
printf("z[%2d] = %9.3f %+9.3f j\n", k, RE(z[k]), IM(z[k]));
/* print H(s) poles values */
printf("Elliptic filter poles: %d\n", np);
for(k = 0; k < np; k++)
printf("p[%2d] = %9.3f %+9.3f j\n", k, RE(p[k]), IM(p[k]));
/* Write complex poles to the file */
writetxt_cmplx(z, nz, "dat/ellip_ap_z.txt");
writetxt_cmplx(p, np, "dat/ellip_ap_p.txt");
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 440, 360, "img/ellip_ap_zp_test.png", &hplot);
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set xzeroaxis");
gnuplot_cmd(hplot, "set yzeroaxis");
gnuplot_cmd(hplot, "set xlabel 'sigma'");
gnuplot_cmd(hplot, "set ylabel 'jw'");
gnuplot_cmd(hplot, "set size square");
gnuplot_cmd(hplot, "set xrange [-2:2]");
gnuplot_cmd(hplot, "set yrange [-2:2]");
gnuplot_cmd(hplot, "plot 'dat/ellip_ap_p.txt' with points pt 2, \\");
gnuplot_cmd(hplot, " 'dat/ellip_ap_z.txt' with points pt 6");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return res;
}

Result:

Elliptic filter zeros: 6
z[ 0] =     0.000    +1.053 j
z[ 1] =     0.000    -1.053 j
z[ 2] =     0.000    +1.136 j
z[ 3] =     0.000    -1.136 j
z[ 4] =     0.000    +1.626 j
z[ 5] =     0.000    -1.626 j
Elliptic filter poles: 7
p[ 0] =    -0.358    +0.000 j
p[ 1] =    -0.011    +1.000 j
p[ 2] =    -0.011    -1.000 j
p[ 3] =    -0.060    +0.940 j
p[ 4] =    -0.060    -0.940 j
p[ 5] =    -0.206    +0.689 j
p[ 6] =    -0.206    -0.689 j


In dat folder will be created ellip_ap_z.txt and ellip_ap_z.txt files which keeps zeros and poles vectors.
In addition, GNUPLOT will build the following graphs from data stored in the files:

Author
Sergey Bakhurin www.dsplib.org

Definition at line 1786 of file filter_ap.c.

Referenced by ellip_ap().

◆ iir()

int iir ( double  rp,
double  rs,
int  ord,
double  w0,
double  w1,
int  type,
double *  b,
double *  a 
)

Digital IIR filter design.


The function calculates the coefficients of the digital IIR filter transfer fucntion \( H(z) \). Filter coeffitients can be used in filter_iir function

Parameters
[in]rpMagnitude ripple in passband (dB).

[in]rsSuppression level in stopband (dB).

[in]ordFilter order.
Number of \(H(z)\) numerator and denominator coefficients is ord+1.
For bandpass and bandstop filters ord must be even.

[in]w0Normalized cutoff frequency (from 0 to 1) for lowpass or highpass filter.
Or left normalized cutoff frequency (from 0 to 1) for bandpass and bandstop filter.

[in]w1Right normalized cutoff frequency (from 0 to 1) for bandpass and bandstop filter.
This parameter is ingnored for lowpass and highpass filters.
[in]typeFilter type.
This patameter sets combination of filter type (one of follow):
DSPL_FILTER_LPF   - lowpass    filter;
DSPL_FILTER_HPF   - highpass filter;
DSPL_FILTER_BPASS - bandpass filter;
DSPL_FILTER_BSTOP - bandstop filter,
and of filter approximation type (one of follow):
DSPL_FILTER_BUTTER - Butterworth filter;
DSPL_FILTER_CHEBY1 - Chebyshev of the first kind filter;
DSPL_FILTER_CHEBY2 - Chebyshev of the second kind filter;
DSPL_FILTER_ELLIP  - Elliptic filter.


[out]bPointer to the transfer function \(H(z)\) numerator coefficients vector.
Vector size is ord+1.
Memory must be allocated.

[out]aPointer to the transfer function \(H(z)\) denominator coefficients vector.
Vector size is ord+1.

Returns
RES_OK if filter is calculated successfully.
Else code error.

Example:

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
/* Low-pass filter order */
#define LPF_ORD 6
/* High-pass filter order */
#define HPF_ORD 6
/* band-pass filter order */
#define BPF_ORD 12
/* Band-stop filter order */
#define BSF_ORD 12
/* Maximum filter order */
#define MAX_ORD BPF_ORD
/* Stopband suppression (dB) */
#define RS 60.0
/* Pass-band maximum distortion (dB) */
#define RP 2.0
/* Frequency response vector size */
#define N 1024
/*******************************************************************************
* function calculates filter frequency response and save magnitude to
* the text file.
* params: b - pointer to the transfer fuction H(z) numerator vector
* a - pointer to the transfer fuction H(z) denominator vector
* ord - filter order
* n - number of magnitude vector size
* fn - file name
******************************************************************************/
void freq_resp_write2txt(double* b, double* a, int ord, int n, char* fn)
{
double *w = NULL, *mag = NULL;
int k;
w = (double*)malloc(n*sizeof(double));
mag = (double*)malloc(n*sizeof(double));
/* Normalized frequency from 0 to pi */
linspace(0, M_PI, n , DSPL_PERIODIC, w);
/* Magnitude (dB) calculation */
filter_freq_resp(b, a, ord, w, n, DSPL_FLAG_LOGMAG, mag, NULL, NULL);
/* Frequency normalization from 0 to 1 */
for(k = 0; k < N; k++)
w[k] /= M_PI;
/* Save magnitude to the txt file */
writetxt(w, mag, n, fn);
free(w);
free(mag);
}
/*******************************************************************************
* Main program
******************************************************************************/
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
void* hplot; /* GNUPLOT handle */
hdspl = dspl_load(); /* Load DSPL functions */
/* Transfer function H(z) coeff. vectors */
double a[MAX_ORD+1], b[MAX_ORD+1];
/*------------------------------------------------------------------------*/
/* LPF Batterworth */
iir(RP, RS, LPF_ORD, 0.3, 0.0, DSPL_FILTER_BUTTER | DSPL_FILTER_LPF, b, a);
freq_resp_write2txt(b, a, LPF_ORD, N, "dat/iir_butter_lpf.txt");
/* HPF Batterworth */
iir(RP, RS, HPF_ORD, 0.3, 0.0, DSPL_FILTER_BUTTER | DSPL_FILTER_HPF, b, a);
freq_resp_write2txt(b, a, HPF_ORD, N, "dat/iir_butter_hpf.txt");
/* Band-pass Batterworth */
iir(RP, RS, BPF_ORD, 0.3, 0.7, DSPL_FILTER_BUTTER | DSPL_FILTER_BPASS, b, a);
freq_resp_write2txt(b, a, BPF_ORD, N, "dat/iir_butter_bpf.txt");
/* Band-stop Batterworth */
iir(RP, RS, BSF_ORD, 0.3, 0.7, DSPL_FILTER_BUTTER | DSPL_FILTER_BSTOP, b, a);
freq_resp_write2txt(b, a, BSF_ORD, N, "dat/iir_butter_bsf.txt");
/*------------------------------------------------------------------------*/
/* LPF Chebyshev type 1 */
iir(RP, RS, LPF_ORD, 0.3, 0.0, DSPL_FILTER_CHEBY1 | DSPL_FILTER_LPF, b, a);
freq_resp_write2txt(b, a, LPF_ORD, N, "dat/iir_cheby1_lpf.txt");
/* HPF Chebyshev type 1 */
iir(RP, RS, HPF_ORD, 0.3, 0.0, DSPL_FILTER_CHEBY1 | DSPL_FILTER_HPF, b, a);
freq_resp_write2txt(b, a, HPF_ORD, N, "dat/iir_cheby1_hpf.txt");
/* Bnad-pass Chebyshev type 1 */
iir(RP, RS, BPF_ORD, 0.3, 0.7, DSPL_FILTER_CHEBY1 | DSPL_FILTER_BPASS, b, a);
freq_resp_write2txt(b, a, BPF_ORD, N, "dat/iir_cheby1_bpf.txt");
/* Bnad-stop Chebyshev type 1 */
iir(RP, RS, BSF_ORD, 0.3, 0.7, DSPL_FILTER_CHEBY1 | DSPL_FILTER_BSTOP, b, a);
freq_resp_write2txt(b, a, BSF_ORD, N, "dat/iir_cheby1_bsf.txt");
/*------------------------------------------------------------------------*/
/* LPF Chebyshev type 2 */
iir(RP, RS, LPF_ORD, 0.3, 0.0, DSPL_FILTER_CHEBY2 | DSPL_FILTER_LPF, b, a);
freq_resp_write2txt(b, a, LPF_ORD, N, "dat/iir_cheby2_lpf.txt");
/* HPF Chebyshev type 2 */
iir(RP, RS, HPF_ORD, 0.3, 0.0, DSPL_FILTER_CHEBY2 | DSPL_FILTER_HPF, b, a);
freq_resp_write2txt(b, a, HPF_ORD, N, "dat/iir_cheby2_hpf.txt");
/* Band-pass Chebyshev type 2 */
iir(RP, RS, BPF_ORD, 0.3, 0.7, DSPL_FILTER_CHEBY2 | DSPL_FILTER_BPASS, b, a);
freq_resp_write2txt(b, a, BPF_ORD, N, "dat/iir_cheby2_bpf.txt");
/* Band-stop Chebyshev type 2 */
iir(RP, RS, BSF_ORD, 0.3, 0.7, DSPL_FILTER_CHEBY2 | DSPL_FILTER_BSTOP, b, a);
freq_resp_write2txt(b, a, BSF_ORD, N, "dat/iir_cheby2_bsf.txt");
/*------------------------------------------------------------------------*/
/* LPF Elliptic */
iir(RP, RS, LPF_ORD, 0.3, 0.0, DSPL_FILTER_ELLIP | DSPL_FILTER_LPF, b, a);
freq_resp_write2txt(b, a, LPF_ORD, N, "dat/iir_ellip_lpf.txt");
/* HPF Elliptic */
iir(RP, RS, HPF_ORD, 0.3, 0.0, DSPL_FILTER_ELLIP | DSPL_FILTER_HPF, b, a);
freq_resp_write2txt(b, a, HPF_ORD, N, "dat/iir_ellip_hpf.txt");
/* Band-pass Elliptic */
iir(RP, RS, BPF_ORD, 0.3, 0.7, DSPL_FILTER_ELLIP | DSPL_FILTER_BPASS, b, a);
freq_resp_write2txt(b, a, BPF_ORD, N, "dat/iir_ellip_bpf.txt");
/* Band-stop Elliptic */
iir(RP, RS, BSF_ORD, 0.3, 0.7, DSPL_FILTER_ELLIP | DSPL_FILTER_BSTOP, b, a);
freq_resp_write2txt(b, a, BSF_ORD, N, "dat/iir_ellip_bsf.txt");
/*------------------------------------------------------------------------*/
/* plotting by GNUPLOT */
gnuplot_create(argc, argv, 920, 840, "img/iir_test.png", &hplot);
gnuplot_cmd(hplot, "unset key");
gnuplot_cmd(hplot, "set grid");
gnuplot_cmd(hplot, "set xlabel 'normalized frequency'");
gnuplot_cmd(hplot, "set ylabel 'Magnitude, dB'");
gnuplot_cmd(hplot, "set yrange [-100:5]");
gnuplot_cmd(hplot, "set xtics 0,1");
gnuplot_cmd(hplot, "set xtics add ('0.3' 0.3)");
gnuplot_cmd(hplot, "set xtics add ('0.7' 0.7)");
gnuplot_cmd(hplot, "set xtics add ('1' 1)");
gnuplot_cmd(hplot, "set multiplot layout 4,4 rowsfirst");
gnuplot_cmd(hplot, "plot 'dat/iir_butter_lpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_butter_hpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_butter_bpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_butter_bsf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby1_lpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby1_hpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby1_bpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby1_bsf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby2_lpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby2_hpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby2_bpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_cheby2_bsf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_ellip_lpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_ellip_hpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_ellip_bpf.txt' with lines");
gnuplot_cmd(hplot, "plot 'dat/iir_ellip_bsf.txt' with lines");
gnuplot_cmd(hplot, "unset multiplot");
gnuplot_close(hplot);
/* free dspl handle */
dspl_free(hdspl);
return 0;
}

This program calcultes filter coefficients for different flags type.

In addition, the filters magnitudes is calculated and plotted by GNUPLOT package.

Author
Sergey Bakhurin www.dsplib.org

Definition at line 404 of file filter_iir.c.

◆ low2high()

int low2high ( double *  b,
double *  a,
int  ord,
double  w0,
double  w1,
double *  beta,
double *  alpha 
)

Lowpass to highpass filter frequency transform.


Function transforms lowpass filter transfer function \( H(s) \) to the highpass filter transfer function \( F(s) \).

Filter order, magnitude ripple in passband and stopband supression still the same.

Parameters
[in]bPointer to the lowpass filter transfer function \(H(s)\) numerator coefficients vector.
Vector size is [ord+1 x 1].

[in]aPointer to the lowpass filter transfer function \(H(s)\) denominator coefficients vector.
Vector size is [ord+1 x 1].

[in]ordFilter order.

[in]w0Lowpass filter cutoff frequency.

[in]w1Highpass filter cutoff frequency after transformation.

[in,out]betaPointer to the highwpass filter transfer function \(F(s)\) numerator coefficients vector after transformation.
Vector size is [ord+1 x 1].
Memory must be allocated.

[in,out]alphaPointer to the highwpass filter transfer function \(F(s)\) denominator coefficients vector after transformation.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Author
Sergey Bakhurin www.dsplib.org

Definition at line 280 of file filter_ft.c.

Referenced by iir().

◆ low2low()

int low2low ( double *  b,
double *  a,
int  ord,
double  w0,
double  w1,
double *  beta,
double *  alpha 
)

Lowpass to lowpass filter frequency transform

Function transforms lowpass filter transfer function \( H(s) \) to the lowpass filter transfer function \( F(s) \) with other cutoff frequency.

Filter order, magnitude ripple in passband and stopband supression still the same.

Parameters
[in]bPointer to the input lowpass filter transfer function \(H(s)\) numerator coefficients vector.
Vector size is [ord+1 x 1].

[in]aPointer to the input lowpass filter transfer function \(H(s)\) denominator coefficients vector.
Vector size is [ord+1 x 1].

[in]ordFilter order.

[in]w0Input lowpass filter cutoff frequency.

[in]w1Lowpass filter cutoff frequency after transformation.

[in,out]betaPointer to the lowpass filter transfer function \(F(s)\) numerator coefficients vector after transformation.
Vector size is [ord+1 x 1].
Memory must be allocated.

[in,out]alphaPointer to the lowpass filter transfer function \(F(s)\) denominator coefficients vector after transformation.
Vector size is [ord+1 x 1].
Memory must be allocated.

Returns
RES_OK if filter coefficients is calculated successfully.
Else code error.
Author
Sergey Bakhurin www.dsplib.org

Definition at line 426 of file filter_ft.c.

Referenced by iir().

◆ ratcompos()

int ratcompos ( double *  b,
double *  a,
int  n,
double *  c,
double *  d,
int  p,
double *  beta,
double *  alpha 
)

Rational composition.



Function calcultes composition \(Y(s) = (H \circ F)(s) = H(F(s))\), here

\[ H(s) = \frac{\sum\limits_{m = 0}^{n} b_m s^m} {\sum\limits_{k = 0}^{n} a_k s^k}, \quad F(s) = \frac{\sum\limits_{m = 0}^{p} d_m s^m} {\sum\limits_{k = 0}^{p} c_k s^k}, \quad Y(s) = \frac{\sum\limits_{m = 0}^{n p} \beta_m s^m} {\sum\limits_{k = 0}^{n p} \alpha_k s^k} \]

This function is using for filter frequency transform.

Parameters
[in]bPointer to the \(H(s)\) polynomial function numerator coefficients vector.
Vector size is [n+1 x 1].

[in]aPointer to the \(H(s)\) polynomial function denominator coefficients vector.
Vector size is [n+1 x 1].

[in]nOrder of \(H(s)\) numerator and denominator polynomials.

[in]cPointer to the \(F(s)\) polynomial function numerator coefficients vector.
Vector size is [p+1 x 1].

[in]dPointer to the \(F(s)\) polynomial function denominator coefficients vector.
Vector size is [p+1 x 1].

[in]pOrder of \(F(s)\) numerator and denominator polynomials.

[in,out]betaPointer to the numerator coefficients vector of \(Y(s) = (H \circ F)(s)\).
Vector size is [n*p+1 x 1].
Memory must be allocated.

[in,out]alphaPointer to the denominator coefficients vector of \(Y(s) = (H \circ F)(s)\).
Vector size is [n*p+1 x 1].
Memory must be allocated.

Returns
RES_OK if rational composition is calculated successfully.
Else code error.
Author
Sergey Bakhurin www.dsplib.org

Definition at line 602 of file filter_ft.c.

Referenced by bilinear(), low2high(), and low2low().

int DSPL_API ellip_ap_zp(int ord, double rp, double rs, complex_t *z, int *nz, complex_t *p, int *np)
Function calculates arrays of zeros and poles for analog normlized lowpass elliptic filter transfer f...
Definition: filter_ap.c:1786
int DSPL_API cheby1_ap(double rp, int ord, double *b, double *a)
Function calculates the transfer function coefficients of analog normalized lowpass Chebyshev type 1...
Definition: filter_ap.c:601
int DSPL_API butter_ap_zp(int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
Function calculates arrays of zeros and poles for analog normlized lowpass Batterworth filter transfe...
Definition: filter_ap.c:402
int DSPL_API bilinear(double *bs, double *as, int ord, double *bz, double *az)
Transform a s-plane analog filter transfer function to the digital filter transfer function .
Definition: filter_iir.c:210
#define RE(x)
Macro sets real part of the complex number.
Definition: dspl.h:359
void DSPL_API gnuplot_close(void *h)
Close GNUPLOT handle.
Definition: gnuplot.c:322
int DSPL_API iir(double rp, double rs, int ord, double w0, double w1, int type, double *b, double *a)
Digital IIR filter design.
Definition: filter_iir.c:404
int DSPL_API logspace(double x0, double x1, int n, int type, double *x)
Function fills a vector with n logarithmically spaced elements between and .
Definition: array.c:1192
void DSPL_API gnuplot_cmd(void *h, char *cmd)
Function sends cmd command to GNUPLOT corresponds to h handle.
Definition: gnuplot.c:387
int DSPL_API butter_ap(double rp, int ord, double *b, double *a)
Function calculates the transfer function coefficients of analog normalized lowpass Butterworth filt...
Definition: filter_ap.c:175
int DSPL_API writetxt(double *x, double *y, int n, char *fn)
Save real data to the text file fn. .
Definition: inout.c:464
int DSPL_API freqz(double *b, double *a, int ord, double *w, int n, complex_t *h)
Function calculates the digital filter frequency response corresponds to transfer function .
Definition: filter_an.c:788
double complex_t[2]
Complex data type.
Definition: dspl.h:86
int DSPL_API cheby1_ap_zp(int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 1 filter tr...
Definition: filter_ap.c:832
int DSPL_API filter_freq_resp(double *b, double *a, int ord, double *w, int n, int flag, double *mag, double *phi, double *tau)
Magnitude, phase response and group delay vectors calculation for digital or analog filter correspond...
Definition: filter_an.c:238
int DSPL_API ellip_ap(double rp, double rs, int ord, double *b, double *a)
Function calculates the transfer function coefficients of analog normalized lowpass elliptic filter ...
Definition: filter_ap.c:1538
#define RES_OK
The function completed correctly. No errors.
Definition: dspl.h:497
#define ABSSQR(x)
The macro returns the square of the modulus of a complex number x.
Definition: dspl.h:473
int DSPL_API linspace(double x0, double x1, int n, int type, double *x)
Function fills a vector with n linearly spaced elements between x0 and x1.
Definition: array.c:1009
#define IM(x)
Macro sets imaginary part of the complex number.
Definition: dspl.h:417
int DSPL_API cheby2_ap_zp(int ord, double rs, complex_t *z, int *nz, complex_t *p, int *np)
Function calculates arrays of zeros and poles for analog normlized lowpass Chebyshev type 2 filter tr...
Definition: filter_ap.c:1304
int DSPL_API gnuplot_create(int argc, char *argv[], int w, int h, char *fn_png, void **hplot)
Create GNUPLOT chart.
Definition: gnuplot.c:200
int DSPL_API cheby2_ap(double rs, int ord, double *b, double *a)
Function calculates the transfer function coefficients of analog normalized lowpass Chebyshev type 2...
Definition: filter_ap.c:1037