Функции анализа аналоговых и цифровых фильтров

Функции

int freqs (double *b, double *a, int ord, double *w, int n, complex_t *h)
 Расчет комплексного коэффициента передачи $ H(j \omega)$ аналогового фильтра. Подробнее...
 
int freqs_resp (double *b, double *a, int ord, double *w, int n, int flag, double *h, double *phi, double *tau)
 Расчет амплитудно-частотной (АЧХ), фазочастотной характеристик (ФЧХ), а также группового времени запаздывания (ГВЗ) аналогового фильтра. Подробнее...
 
int freqz (double *b, double *a, int ord, double *w, int n, complex_t *h)
 Расчет комплексного коэффициента передачи $ H \left(e^{j \omega} \right)$ цифрового фильтра. Подробнее...
 

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

Функции анализа аналоговых и цифровых фильтров

Функции

int freqs ( double *  b,
double *  a,
int  ord,
double *  w,
int  n,
complex_t h 
)

Расчет комплексного коэффициента передачи $ H(j \omega)$ аналогового фильтра.


Функция рассчитывает значения комплексного коэффициента передачи $H(j \omega)$ аналогового фильтра, заданного коэффициентами передаточной функции $H(s)$:

\[ H(s) = \frac{\sum_{k = 0}^{N} b_k \cdot s^k}{\sum_{m = 0}^{N} a_m \cdot s^m}, \]

где $N$ - порядок фильтра (параметр ord).

Комплексный коэффициент передачи рассчитывается путем подстановки $s = j \omega $.

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

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

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

[in]wУказатель на вектор значений циклической частоты $\omega$ (рад/с), для которого будет рассчитан комплексный коэффициент передачи $H(j \omega)$.
Размер вектора [n x 1].

[in]nРазмер вектора циклической частоты w.

[out]hУказатель на вектор комплексного коэффициента передачи $H(j \omega)$, рассчитанного для циклической частоты w.
Размер вектора [n x 1].
Память должна быть выделена.

Возвращает
RES_OK Комплексноый коэффициент передачи рассчитан успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org
int freqs_resp ( double *  b,
double *  a,
int  ord,
double *  w,
int  n,
int  flag,
double *  h,
double *  phi,
double *  tau 
)

Расчет амплитудно-частотной (АЧХ), фазочастотной характеристик (ФЧХ), а также группового времени запаздывания (ГВЗ) аналогового фильтра.


Функция рассчитывает АЧХ, ФЧХ и ГВЗ аналоговго фильтра, заданного передаточной характеристикой $H(s)$

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

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

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

[in]wУказатель на вектор значений циклической частоты $\omega$ (рад/с), для которого будет рассчитаны АЧХ, ФЧХ и ГВЗ.
Размер вектора [n x 1].

[in]nРазмер вектора циклической частоты w.

[in]flagКомбинация флагов, которые задают расчет параметров:
DSPL_FLAG_LOG        АЧХ рассчитывать в логарифмическом масштабе
DSPL_FLAG_UNWRAP     раскрывать периодичность ФЧХ`
[out]hУказатель на вектор АЧХ.
Размер вектора [n x 1].
Память должна быть выделена.
Если указатель NULL, то рассчет АЧХ не производится.

[out]phiУказатель на вектор ФЧХ.
Размер вектора [n x 1].
Память должна быть выделена.
Если указатель NULL, то рассчет ФЧХ не производится.

[out]tauУказатель на вектор ГВЗ.
Размер вектора [n x 1].
Память должна быть выделена.
Если указатель NULL, то рассчет ГВЗ не производится.

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

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
#define ORD 4
#define N   1000
int main()
{
    void* handle;           // DSPL handle
    handle = dspl_load();   // Load DSPL function
    double a[ORD+1], b[ORD+1];
    double Rp = 3.0;
    double w[N], mag[N], phi[N], tau[N];
   
    int k;
    int res = butter_ap(Rp, ORD, b, a);
    if(res != RES_OK)
        printf("error code = 0x%8x\n", res);
    
    for(k = 0; k < ORD+1; k++)
        printf("b[%2d] = %9.3f     a[%2d] = %9.3f\n", k, b[k], k, a[k]);
    logspace(-2.0, 2.0, N , DSPL_SYMMETRIC, w);
    freqs_resp(b, a, ORD, w, N, DSPL_FLAG_LOG|DSPL_FLAG_UNWRAP, mag, phi, tau);
    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");
    dspl_free(handle);      // free dspl handle
    res = system("gnuplot gnuplot/butter_ap_test.plt");
    return 0;
}

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

    b[ 0] =     1.002     a[ 0] =     1.002
    b[ 1] =     0.000     a[ 1] =     2.618
    b[ 2] =     0.000     a[ 2] =     3.418
    b[ 3] =     0.000     a[ 3] =     2.615
    b[ 4] =     0.000     a[ 4] =     1.000



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

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

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

АЧХ фильтра: butter_ap_test_mag.png

butter_ap_test_mag.png

ФЧХ фильтра: butter_ap_test_phi.png

butter_ap_test_phi.png

ГВЗ фильтра: butter_ap_test_tau.png

butter_ap_test_tau.png
Возвращает
RES_OK Параметры фильтра рассчитаны успешно.
В противном случае код ошибки.
Автор
Бахурин Сергей www.dsplib.org
int freqz ( double *  b,
double *  a,
int  ord,
double *  w,
int  n,
complex_t h 
)

Расчет комплексного коэффициента передачи $ H \left(e^{j \omega} \right)$ цифрового фильтра.


Функция рассчитывает значения комплексного коэффициента передачи $ H \left(e^{j \omega} \right)$ цифрового фильтра, заданного коэффициентами передаточной функции $H(z)$:

\[ H(z) = \frac{\sum_{k = 0}^{N} b_k \cdot z^{-k}}{\sum_{m = 0}^{N} a_m \cdot z^{-m}}, \]

где $N$ - порядок фильтра (параметр ord).

Комплексный коэффициент передачи рассчитывается путем подстановки $z = e^{j \omega} $.

Аргументы
[in]bУказатель на вектор коэффициентов числителя передаточной функции $H(z)$.
Размер вектора [ord+1 x 1].

[in]aУказатель на вектор коэффициентов знаменателя передаточной функции $H(z)$.
Размер вектора [ord+1 x 1].

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

[in]wУказатель на вектор значений нормированной циклической частоты $\omega$, для которого будет рассчитан комплексный коэффициент передачи $ H \left(e^{j \omega} \right)$.
Размер вектора [n x 1].

[in]nРазмер вектора нормированной циклической частоты w.

[out]hУказатель на вектор комплексного коэффициента передачи $ H \left(e^{j \omega} \right)$, рассчитанного для циклической частоты w.
Размер вектора [n x 1].
Память должна быть выделена.

Возвращает
RES_OK Комплексный коэффициент передачи расчитан успешно.
В противном случае код ошибки.
Заметки
Комплексный коэффициент передачи $ H \left(e^{j \omega} \right)$ цифрового фильтра представляет собой $ 2 \pi-$периодическую функцию нормированной циклической частоты $\omega$. Поэтому анализ цифровых фильтров целесообразно вести на одном периоде повторения $ H \left(e^{j \omega} \right)$, т.е. в интервале $\omega$ от 0 до $2 \pi$, или от $-\pi$ до $ \pi$.
Кроме того известно, что для фильтра с вещественными коэффициентами $ H \left(e^{j \omega} \right) = H^* \left(e^{-j \omega} \right)$, а значит, анализ цифрового фильтра с вещественными коэффициентами достаточно вести для нормированной частоты $\omega$ от 0 до $\pi$.
Автор
Бахурин Сергей www.dsplib.org