libdspl-2.0
Библиотека алгоритмов цифровой обработки сигналов
Расчет БИХ-фильтров

Функции расчета цифровых БИХ-фильтров. Подробнее...

Функции

int bilinear (double *bs, double *as, int ord, double *bz, double *az)
 Билинейное преобразование передаточной характеристики аналогового фильтра \( H(s) \), в передаточную характеристику цифрового фильтра \( H(z) \). Подробнее...
 
int butter_ap (double rp, int ord, double *b, double *a)
 Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Баттерворта. Подробнее...
 
int butter_ap_zp (int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
 Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Баттерворта. Подробнее...
 
int cheby1_ap (double rp, int ord, double *b, double *a)
 Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Чебышёва первого рода. Подробнее...
 
int cheby1_ap_zp (int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
 Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Чебышёва первого рода. Подробнее...
 
int cheby2_ap (double rs, int ord, double *b, double *a)
 Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Чебышёва второго рода. Подробнее...
 
int cheby2_ap_zp (int ord, double rs, complex_t *z, int *nz, complex_t *p, int *np)
 Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Чебышёва второго рода. Подробнее...
 
int ellip_ap (double rp, double rs, int ord, double *b, double *a)
 Расчет передаточной характеристики \( H(s) \) аналогового нормированного эллиптического ФНЧ. Подробнее...
 
int ellip_ap_zp (int ord, double rp, double rs, complex_t *z, int *nz, complex_t *p, int *np)
 Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного эллиптического ФНЧ. Подробнее...
 
int filter_zp2ab (complex_t *z, int nz, complex_t *p, int np, int ord, double *b, double *a)
 Функция пересчета нулей и полюсов аналогового фильтра в коэффициенты передаточной характеристики \( H(s) \). Подробнее...
 
int iir (double rp, double rs, int ord, double w0, double w1, int type, double *b, double *a)
 Функция расчета коэффициентов передаточной характеристики \( H(z) \) цифрового фильтра БИХ. Подробнее...
 
int low2high (double *b, double *a, int ord, double w0, double w1, double *beta, double *alpha)
 Частотное преобразование ФНЧ-ФВЧ Подробнее...
 
int low2low (double *b, double *a, int ord, double w0, double w1, double *beta, double *alpha)
 Частотное преобразование ФНЧ-ФНЧ Подробнее...
 
int ratcompos (double *b, double *a, int n, double *c, double *d, int p, double *beta, double *alpha)
 Рациональная композиця Подробнее...
 

Подробное описание

Функции расчета цифровых БИХ-фильтров.

Функции

◆ bilinear()

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

Билинейное преобразование передаточной характеристики аналогового фильтра \( H(s) \), в передаточную характеристику цифрового фильтра \( H(z) \).


Функция рассчитывает коэффициенты передаточной характеристики \( H(z) \) цифрового фильтра путем дробно-рациональной подстановки вида

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

Порядок цифрового фильтра при этом остается равным порядку аналогового фильтра, а ось частот \( \Omega \) аналогового фильтра связана c осью частот \( \omega \) цифрового фильтра соотношением:

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

Аргументы
[in]bsУказатель на вектор коэффициентов числителя передаточной функции \( H(s) \) исходного аналогового фильтра.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]asУказатель на вектор коэффициентов знаменателя передаточной функции \( H(s) \) исходного аналогового фильтра.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточных функций \( H(s) \) и \(H(z)\) аналогового и цифрового фильтров равно ord+1.

[out]bzУказатель на вектор коэффициентов числителя передаточной функции \( H(z) \) полученного цифрового фильтра.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]azУказатель на вектор коэффициентов знаменателя передаточной функции \(H(z)\) полученного цифрового фильтра.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — фильтр рассчитан успешно.
В противном случае код ошибки.

Пример использования функции bilinear:

#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;
}

Данная программа производит расчет передаточной характеристики аналогового фильтра Чебышева первого рода, с частотой среза равной 1 рад/с, и производит билинейное преобразование в цифровой, с частотой среза равной 0.5.

Результат работы программы:

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

Кроме этого производится расчет АЧХ полученного цифрового фильтра и строится график АЧХ пакетом GNUPLOT

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле bilinear.c строка 209

Используется в iir().

◆ butter_ap()

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

Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Баттерворта.


Функция рассчитывает коэффициенты передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Баттерворта порядка ord с частотой среза 1 рад/с по уровню \( -R_p \) дБ.

Аргументы
[in]rpНеравномерность АЧХ в полосе пропускания (дБ).
Параметр задает уровень искажений в полосе от 0 до 1 рад/с.
Значение должно быть положительным.

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточной функции \( H(s) \) равно ord+1.

[out]bУказатель на вектор коэффициентов числителя передаточной функции \( H(s) \).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \( H(s) \).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — фильтр рассчитан успешно.
В противном случае код ошибки.

Пример использования функции butter_ap:
#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();
/* H(s) numerator coefficients vector */
double a[ORD+1];
/* H(s) denominator coefficients vector */
double b[ORD+1];
/* Magnitude ripple from 0 to 1 rad/s */
double Rp = 1.0;
/* Angular frequency (rad/s) */
double w[N];
/* Filter Magnitude (dB) */
double mag[N];
/* Phase response */
double phi[N];
/* Group delay */
double tau[N];
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;
}

Результат работы программы:

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


В каталоге dat будут созданы три файла:

butter_ap_test_mag.txt АЧХ фильтра
butter_ap_test_phi.txt ФЧХ фильтра
butter_ap_test_tau.txt ГВЗ фильтра

Кроме того программа GNUPLOT произведет построение следующих графиков по сохраненным в файлах данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле butter_ap.c строка 172

◆ butter_ap_zp()

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

Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Баттерворта.


Функция рассчитывает значения нулей и полюсов передаточной функции \( H(s)\) аналогового нормированного ФНЧ Баттерворта порядка ord с частотой среза 1 рад/с по уровню \(-R_p\) дБ.

Аргументы
[in]ordПорядок фильтра.

[in]rpНеравномерность АЧХ в полосе пропускания (дБ).
Параметр задает уровень искажений в полосе от 0 до 1 рад/с.
Значение должно быть положительным.

[out]zУказатель на массив комплексных нулей передаточной характеристики \( H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]nzУказатель на переменную количества нулей передаточной характеристики \( H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор z.
Память должна быть выделена.

[out]pУказатель на массив комплексных полюсов передаточной характеристики \( H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]npУказатель на переменную количества полюсов передаточной характеристики \( H(s)\).
По данному укащзателю будет записано количество нулей фильтра, которые были рассчитны и помещены в вектор p.
Память должна быть выделена.

Возвращает
RES_OK — массивы нулей и полюсов рассчитаны успешно.
В противном случае код ошибки.

Заметки
Нормированный ФНЧ Баттерворта не имеет нулей, поэтому массив нулей z не будет изменен, а по указателю nz будет записан 0.

Пример программы рассчета нулей и полюсов нормированного ФНЧ Баттерворта:

#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;
}

Результат выполнения программы:

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  


В каталоге dat будет создан файл butter_ap_zp.txt.
Пакет GNUPLOT произведет построение карты полюсов по сохранненным в dat/butter_ap_zp.txt данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле butter_ap_zp.c строка 207

Используется в butter_ap().

◆ cheby1_ap()

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

Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Чебышёва первого рода.


Функция рассчитывает коэффициенты передаточной характеристики \( H(s)\) аналогового нормированного ФНЧ Чебышёва первого рода порядка ord с частотой среза 1 рад/с по уровню \(-R_p\) дБ.

Особенностью фильтра Чебышёва первого рода являются равноволновые пульсации АЧХ в полосе пропускания.

Аргументы
[in]rpНеравномерность АЧХ в полосе пропускания (дБ).
Параметр задает уровень искажений в полосе от 0 до 1 рад/с.
Значение должно быть положительным.

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточной функции \( H(s)\) равно ord+1.

[out]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — фильтр рассчитан успешно.
В противном случае код ошибки.

Пример использования функции cheby1_ap:
#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;
}

Результат работы программы:

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


В каталоге dat будут созданы три файла:

cheby1_ap_test_mag.txt    АЧХ фильтра
cheby1_ap_test_phi.txt    ФЧХ фильтра
cheby1_ap_test_tau.txt    ГВЗ фильтра


Кроме того программа GNUPLOT произведет построение следующих графиков по сохраненным в файлах данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле cheby1_ap.c строка 179

◆ cheby1_ap_zp()

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

Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Чебышёва первого рода.


Функция рассчитывает значения нулей и полюсов передаточной функции \( H(s)\) аналогового нормированного ФНЧ Чебышёва первого рода порядка ord с частотой среза 1 рад/с по уровню \(-R_p\) дБ, с неравномерностью в полосе пропускания \( R_p \) дБ.

Аргументы
[in]ordПорядок фильтра.

[in]rpНеравномерность АЧХ в полосе пропускания (дБ).
Параметр задает уровень искажений в полосе от 0 до 1 рад/с.
Значение должно быть положительным.

[out]zУказатель на массив комплексных нулей передаточной характеристики \( H(s)\).
Максимальный размер вектора [ord x 1].
Память должна быть выделена.

[out]nzУказатель на переменную количества нулей передаточной функции \(H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор z.
Память должна быть выделена.

[out]pУказатель на массив комплексных полюсов передаточной характеристики \(H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]npУказатель на переменную количества полюсов передаточной функции \( H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор p.
Память должна быть выделена.

Возвращает
RES_OK — массивы нулей и полюсов рассчитаны успешно.
В противном случае код ошибки.
Заметки
Нормированный ФНЧ Чебышёва первого рода не имеет нулей, поэтому массив нулей z не будет изменен, а по указателю nz будет записан 0.

Пример программы рассчета нулей и полюсов нормированного ФНЧ Чебышева первого рода:

#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;
}

Результат выполнения программы:

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


В каталоге dat будет создан файл cheby1_ap_zp.txt.
Пакет GNUPLOT произведет построение карты полюсов по сохранненным в dat/cheby1_ap_zp.txt данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле cheby1_ap_zp.c строка 204

Используется в cheby1_ap().

◆ cheby2_ap()

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

Расчет передаточной характеристики \( H(s) \) аналогового нормированного ФНЧ Чебышёва второго рода.


Функция рассчитывает коэффициенты передаточной характеристики \(H(s)\) аналогового нормированного ФНЧ Чебышёва второго рода порядка ord с частотой заграждения 1 рад/с по уровню \(-R_s\) дБ.
Особенностью фильтра Чебышёва второго рода являются:
1) равноволновые пульсации АЧХ в полосе заграждения.
2) уровень АЧХ \(H(j\cdot 1) = -R_s\) дБ.

Аргументы
[in]rsУровень подавления в полосе пропускания (дБ).
Значение должно быть положительным.

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточной функции \(H(s)\) равно ord+1.

[out]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Пример использования функции cheby1_ap:
#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;
}

Результат работы программы:

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 


В каталоге dat будут созданы три файла:

cheby2_ap_test_mag.txt    АЧХ фильтра   
cheby2_ap_test_phi.txt    ФЧХ фильтра
cheby2_ap_test_tau.txt    ГВЗ фильтра


Кроме того программа GNUPLOT произведет построение следующих графиков по сохраненным в файлах данным:

Возвращает
RES_OK — фильтр рассчитан успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле cheby2_ap.c строка 182

◆ cheby2_ap_zp()

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

Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного ФНЧ Чебышёва второго рода.


Функция рассчитывает значения нулей и полюсов передаточной функции \(H(s)\) аналогового нормированного ФНЧ Чебышёва второго рода порядка ord с частотой заграждения 1 рад/с по уровню \(-R_s\) дБ.

Аргументы
[in]ordПорядок фильтра.

[in]rsУровень подавления АЧХ в полосе загражения (дБ).
Параметр задает уровень подавления сигнала в полосе частот от 1 рад/с и выше.
Значение должно быть положительным.

[out]zУказатель на массив комплексных нулей передаточной функции \(H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]nzУказатель на переменную количества нулей передаточной функции \(H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор z.
Память должна быть выделена.

[out]pУказатель на массив комплексных полюсов передаточной функции \(H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]npУказатель на переменную количества полюсов передаточной функции \(H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор p.
Память должна быть выделена.

Возвращает
RES_OK — массивы нулей и полюсов рассчитаны успешно.
В противном случае код ошибки.
Пример использования функции cheby2_ap_zp:

Пример программы рассчета нулей и полюсов нормированного ФНЧ Чебышева первого рода:

#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;
}

Результат выполнения программы:

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 


В каталоге dat будет создан файлы cheby2_ap_z.txt и cheby2_ap_z.txt, хранящие наборы нулей и полюсов на комплексной плоскости.
Пакет GNUPLOT произведет построение карты полюсов по сохранненным в dat/cheby2_ap_z.txt и dat/cheby2_ap_p.txt данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле cheby2_ap_zp.c строка 216

Используется в cheby2_ap().

◆ ellip_ap()

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

Расчет передаточной характеристики \( H(s) \) аналогового нормированного эллиптического ФНЧ.


Функция рассчитывает коэффициенты передаточной характеристики \(H(s)\) аналогового нормированного эллиптического ФНЧ порядка ord с частотой среза 1 рад/с по уровню \(-R_p\) дБ.
Особенностью эллиптического фильтра являются равноволновые пульсации АЧХ как в полосе пропускания, так и в полосе заграждения, в результате чего обеспечиваеся минимальная переходная полоса фильтра.

Аргументы
[in]rpУровень пульсаций в полосе пропускания (дБ).
Значение должно быть положительным.

[in]rsУровень подавления в полосе заграждения (дБ).
Значение должно быть положительным.

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточной функции \(H(s)\) равно ord+1.

[out]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Пример использования функции ellip_ap:
#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;
}

Результат работы программы:

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


В каталоге dat будут созданы три файла:

ellip_ap_test_mag.txt    АЧХ фильтра   
ellip_ap_test_phi.txt    ФЧХ фильтра
ellip_ap_test_tau.txt    ГВЗ фильтра


Кроме того программа GNUPLOT произведет построение следующих графиков по сохраненным в файлах данным:

Возвращает
RES_OK — фильтр рассчитан успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле ellip_ap.c строка 191

◆ ellip_ap_zp()

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

Расчет массивов нулей и полюсов передаточной функции \( H(s) \) аналогового нормированного эллиптического ФНЧ.


Аргументы
[in]ordПорядок фильтра.

[in]rpНеравномерность АЧХ в полосе пропускания (дБ).
Параметр задает уровень искажений в полосе от 0 до 1 рад/с.
Значение должно быть положительным.

[in]rsУровень подавления АЧХ в полосе загражения (дБ).
Параметр задает уровень подавления сигнала в полосе частот от 1 рад/с и выше.
Значение должно быть положительным.

[out]zУказатель на массив комплексных нулей передаточной функции \(H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]nzУказатель на переменную количества нулей передаточной функции \(H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор z.
Память должна быть выделена.

[out]pУказатель на массив комплексных полюсов передаточной функции \(H(s)\).
Максимальный размер вектора вектора [ord x 1].
Память должна быть выделена.

[out]npУказатель на переменную количества полюсов передаточной функции \(H(s)\).
По данному указателю будет записано количество нулей фильтра, которые были рассчитаны и помещены в вектор p.
Память должна быть выделена.

Возвращает
RES_OK — массивы нулей и полюсов рассчитаны успешно.
В противном случае код ошибки.
Пример использования функции cheby2_ap_zp:

Пример программы рассчета нулей и полюсов нормированного эллиптического ФНЧ :

#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;
}

Результат выполнения программы:

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


В каталоге dat будет создан файлы ellip_ap_z.txt и ellip_ap_z.txt, хранящие наборы нулей и полюсов на комплексной плоскости.
Пакет GNUPLOT произведет построение карты полюсов по сохранненным в dat/ellip_ap_z.txt и dat/ellip_ap_p.txt данным:

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле ellip_ap_zp.c строка 213

Используется в ellip_ap().

◆ filter_zp2ab()

int filter_zp2ab ( complex_t z,
int  nz,
complex_t p,
int  np,
int  ord,
double *  b,
double *  a 
)

Функция пересчета нулей и полюсов аналогового фильтра в коэффициенты передаточной характеристики \( H(s) \).


\[ H(s) = \frac{\sum_{n = 0}^{N_z} b_n \cdot s^n}{\sum_{m = 0}^{N_p} a_m \cdot s^m} = \frac{\prod_{n = 0}^{N_z}(s-z_n)}{\prod_{m = 0}^{N_p} (s-p_m)} \]

Аргументы
[in]zУказатель на массив нулей передаточной характеристики.
Размер вектора [nz x 1].
Указатель может быть NULL если фильтр не имеет конечных нулей (nz=0).

[in]nzРазмер вектора нулей передаточной характеристики (может быть равен 0).

[in]pУказатель на массив полюсов передаточной характеристики.
Размер вектора [np x 1].
Указатель не может быть NULL.
Память должна быть выделена.

[in]npРазмер вектора полюсов передаточной характеристики (не может быть равен 0).

[in]ordПорядок фильтра для которого рассчитаны нули и полюса.
Количество коэффициентов числителя и знаменателя передаточной функции \(H(s)\) равно ord+1.

[out]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\).
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — пересчет произведен успешно.
В противном случае код ошибки.
Заметки
Функция возвращает вещественные значения коэффициентов b и a передаточной функции. Это означает, что вектора нулей и полюсов должны хранить вещественные значения или комплексно-сопряженные пары нулей и полюсов, потому что мнимая часть коэффициентов b и a игнорируется и не сохраняется.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле filter_zp2ab.c строка 166

Используется в butter_ap(), cheby1_ap(), cheby2_ap() и ellip_ap().

◆ iir()

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

Функция расчета коэффициентов передаточной характеристики \( H(z) \) цифрового фильтра БИХ.



Функция рассчитывает коэффициенты передаточной характеристики \( H(z) \) цифрового фильтра, которые могут быть использованы в функции filter_iir

Аргументы
[in]rpУровень неравномерности квадрата АЧХ в полосе пропускания фильтра (дБ).
Размер вектора [ord+1 x 1].

[in]rsУровень подавления в полосе заграждения фильтра (дБ).

[in]ordПорядок фильтра.
Количество коэффициентов числителя и знаменателя передаточной функции \(H(z)\) цифрового фильтров равно ord+1.
Для полосовых и режекторных фильтров параметр ord должен быть чётным.

[in]w0Нормированная частота среза ФНЧ или ФВЧ, или левая частота среза для полосового и режекторного фильтра.

[in]w1Правая частота среза полосового и режекторного фильтра.
Данный параметр игнорируется для ФНЧ и ФВЧ.

[in]typeТип фильтра.
Данный параметр определяет тип фильтра и образуется набором флагов типа фильтра:
DSPL_FILTER_LPF - фильтр нижних частот;
DSPL_FILTER_HPF - фильтр верхних частот;
DSPL_FILTER_BPASS - полосовой фильтр;
DSPL_FILTER_BSTOP - режекторный фильтр,
а также флагов типа аппроксимации АЧХ фильтра:
DSPL_FILTER_BUTTER - фильтр Баттерворта;
DSPL_FILTER_CHEBY1 - фильтр Чебышева первого рода;
DSPL_FILTER_CHEBY2 - фильтр Чебышева второго рода;
DSPL_FILTER_ELLIP - эллиптический фильтр.


[out]bУказатель на вектор коэффициентов числителя передаточной функции \(H(z)\).
Размер вектора ord+1.
Память должна быть выделена.

[out]aУказатель на вектор коэффициентов знаменателя передаточной функции \( H(z) \).
Размер вектора ord+1.
Память должна быть выделена.

Возвращает
RES_OK — Фильтр рассчитан успешно.
В противном случае код ошибки.

Пример использования функции:
#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 and check magnitude */
for(k = 0; k < N; k++)
{
w[k] /= M_PI;
/* Set magnitude to -400 dB if it is inf. */
if(isinf(mag[k]))
mag[k] = -400.0;
}
/* 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];
/* Batterworth */
/* 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");
/* Chebyshev type 1 */
/* 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");
/* Chebyshev type 2 */
/* 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");
/* Elliptic */
/* 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;
}

Данная программа производит расчет коэффициентов фильтров при различном сочетании флагов параметра type.

Кроме этого производится расчет АЧХ полученных цифровых фильтров и выводится на график АЧХ пакетом GNUPLOT

Автор
Бахурин Сергей www.dsplib.org

См. определение в файле iir.c строка 245

◆ low2high()

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

Частотное преобразование ФНЧ-ФВЧ


Функция производит перобразование передаточной функции \( H(s) \) аналогового ФНЧ с частотой среза w0 рад/c в передаточную функцию \( F(s) \) аналоговго ФВЧ с частотой среза w1 рад/c.

Неравномерность АЧХ в полосе пропускания, уровень подавления в полосе заграждения и порядок фильтра остаются неизменными.

Аргументы
[in]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\) исходного аналогового ФНЧ.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\) исходного аналогового ФНЧ.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]ordПорядок исходного фильтра и фильтра после переобразования.

[in]w0Частота среза исходного ФНЧ.

[in]w1Требуемая частота среза ФВЧ после преобразования.

[in,out]betaУказатель на вектор коэффициентов числителя передаточной функции \(F(s)\) ФВЧ после преобразования.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in,out]alphaУказатель на вектор коэффициентов знаменателя передаточной функции \(F(s)\) аналогового ФВЧ после преобразования.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — преобразование рассчитано успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле low2high.c строка 146

Используется в iir().

◆ low2low()

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

Частотное преобразование ФНЧ-ФНЧ


Функция производит преобразование передаточной функции \( H(s) \) аналогового ФНЧ с частотой среза w0 рад/c в передаточную функцию \( F(s) \) аналоговго ФНЧ с частотой среза w1 рад/c.

Неравномерность АЧХ в полосе пропускания, уровень подавления в полосе заграждения и порядок фильтра остаются неизменными.

Аргументы
[in]bУказатель на вектор коэффициентов числителя передаточной функции \(H(s)\) исходного аналогового ФНЧ.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]aУказатель на вектор коэффициентов знаменателя передаточной функции \(H(s)\) исходного аналогового ФНЧ.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in]ordПорядок исходного фильтра и фильтра после преобразования.

[in]w0Частота среза исходного ФНЧ.

[in]w1Требуемая частота среза ФНЧ после преобразования.

[in,out]betaУказатель на вектор коэффициентов числителя передаточной функции \(F(s)\) ФНЧ после преобразования.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

[in,out]alphaУказатель на вектор коэффициентов знаменателя передаточной функции \(F(s)\) аналогового ФНЧ после преобразования.
Размер вектора [ord+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — Преоборазование расчитано успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле low2low.c строка 147

Используется в iir().

◆ ratcompos()

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

Рациональная композиця


Функция рассчитывает композицию вида \(Y(s) = (H \circ F)(s) = H(F(s))\), где

\[ 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} \]

Функция рациональной композиции необходима для произведения частотных преобразований передаточных характеристик аналоговых и цифровых фильтров, а также для билинейного преобразования передаточных характеристик аналоговых фильтров в соответствующие передаточные характеристики цифровых фильтров.

Аргументы
[in]bУказатель на вектор коэффициентов числителя функции \(H(s)\).
Размер вектора [n+1 x 1].
Память должна быть выделена.

[in]aУказатель на вектор коэффициентов знаменателя функции \(H(s)\).
Размер вектора [n+1 x 1].
Память должна быть выделена.

[in]nПорядок полиномов рациональной функции \(H(s)\).

[in]cУказатель на вектор коэффициентов числителя функции \(F(s)\).
Размер вектора [p+1 x 1].
Память должна быть выделена.

[in]dУказатель на вектор коэффициентов знаменателя функции \(F(s)\).
Размер вектора [p+1 x 1].
Память должна быть выделена.

[in]pПорядок полиномов рациональной функции \(F(s)\).

[in,out]betaУказатель на вектор коэффициентов числителя функции \(Y(s) = (H \circ F)(s)\).
Размер вектора [n*p+1 x 1].
Память должна быть выделена.

[in,out]alphaУказатель на вектор коэффициентов знаменателя функции \(Y(s) = (H \circ F)(s)\).
Размер вектора [n*p+1 x 1].
Память должна быть выделена.

Возвращает
RES_OK — Рациональная композиция рассчитана успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org

См. определение в файле ratcompos.c строка 180

Используется в bilinear(), low2high() и low2low().

int ellip_ap_zp(int ord, double rp, double rs, complex_t *z, int *nz, complex_t *p, int *np)
Расчет массивов нулей и полюсов передаточной функции аналогового нормированного эллиптического ФНЧ.
Definition: ellip_ap_zp.c:213
int cheby1_ap(double rp, int ord, double *b, double *a)
Расчет передаточной характеристики аналогового нормированного ФНЧ Чебышёва первого рода.
Definition: cheby1_ap.c:179
int butter_ap_zp(int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
Расчет массивов нулей и полюсов передаточной функции аналогового нормированного ФНЧ Баттерворта.
Definition: butter_ap_zp.c:207
int bilinear(double *bs, double *as, int ord, double *bz, double *az)
Билинейное преобразование передаточной характеристики аналогового фильтра , в передаточную характерис...
Definition: bilinear.c:209
#define RE(x)
Макрос определяющий реальную часть комплексного числа.
Definition: dspl.h:420
void gnuplot_close(void *h)
Закрыть хэндл GNUPLOT.
Definition: gnuplot_close.c:80
int iir(double rp, double rs, int ord, double w0, double w1, int type, double *b, double *a)
Функция расчета коэффициентов передаточной характеристики цифрового фильтра БИХ.
Definition: iir.c:245
int logspace(double x0, double x1, int n, int type, double *x)
Функция заполняет массив значениями логарифмической шкале
Definition: logspace.c:178
void * dspl_load()
Произвести динамическую линковку и загрузить функции libdspl-2.0.
void gnuplot_cmd(void *h, char *cmd)
Функция посылает команду cmd пакету GNUPLOT, для построения или оформления графика,...
Definition: gnuplot_cmd.c:82
int butter_ap(double rp, int ord, double *b, double *a)
Расчет передаточной характеристики аналогового нормированного ФНЧ Баттерворта.
Definition: butter_ap.c:172
int writetxt(double *x, double *y, int n, char *fn)
Сохранить вещественные данные в текстовый файл
Definition: writetxt.c:122
int freqz(double *b, double *a, int ord, double *w, int n, complex_t *h)
Расчет комплексного коэффициента передачи цифрового фильтра.
Definition: freqz.c:155
double complex_t[2]
Описание комплексного типа данных.
Definition: dspl.h:86
int cheby1_ap_zp(int ord, double rp, complex_t *z, int *nz, complex_t *p, int *np)
Расчет массивов нулей и полюсов передаточной функции аналогового нормированного ФНЧ Чебышёва первого...
Definition: cheby1_ap_zp.c:204
int filter_freq_resp(double *b, double *a, int ord, double *w, int n, int flag, double *mag, double *phi, double *tau)
Расчет амплитудно-частотной (АЧХ), фазочастотной характеристик (ФЧХ), а также группового времени запа...
void dspl_free(void *handle)
Очищает связанную ранее динамическую библиотеку DSPL-2.0.
int ellip_ap(double rp, double rs, int ord, double *b, double *a)
Расчет передаточной характеристики аналогового нормированного эллиптического ФНЧ.
Definition: ellip_ap.c:191
#define RES_OK
Функция завершилась корректно. Ошибки отсутствуют.
Definition: dspl.h:558
#define ABSSQR(x)
Макрос возвращает квадрат модуля комплексного числа x.
Definition: dspl.h:534
int linspace(double x0, double x1, int n, int type, double *x)
Функция заполняет массив линейно-нарастающими, равноотстоящими значениями от x0 до x1
Definition: linspace.c:169
#define IM(x)
Макрос определяющий мнимую часть комплексного числа.
Definition: dspl.h:478
int cheby2_ap_zp(int ord, double rs, complex_t *z, int *nz, complex_t *p, int *np)
Расчет массивов нулей и полюсов передаточной функции аналогового нормированного ФНЧ Чебышёва второго...
Definition: cheby2_ap_zp.c:216
int gnuplot_create(int argc, char *argv[], int w, int h, char *fn_png, void **hplot)
Создать график GNUPLOT.
int cheby2_ap(double rs, int ord, double *b, double *a)
Расчет передаточной характеристики аналогового нормированного ФНЧ Чебышёва второго рода.
Definition: cheby2_ap.c:182