Представление периодических сигналов рядом Фурье

Содержание

Вводные замечания

Периодический сигнал. Тригонометрический ряд Фурье

Ряд Фурье в комплексной форме

Некоторые пояснения к ряду Фурье в комплексной форме

Спектр периодических сигналов

Выводы

Программная реализация в библиотеке DSPL

Список литературы

DSPL-2.0 — свободная библиотека алгоритмов цифровой обработки сигналов

Распространяется под лицензией LGPL v3

Страница проекта на SourceForge

Вводные замечания

Жан-Батист Жозеф Фурье
1768–1830

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

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

Периодический сигнал. Тригонометрический ряд Фурье

Пусть имеется периодический сигнал s(t) непрерывного времени , который повторяется с периодом с, т.е. s(t+nT) = s(t) , где — произвольное целое число.

В качестве примера на рисунке 1 показана последовательность прямоугольных импульсов длительности $\tau$ c, повторяющиеся с периодом с.

Рисунок 1. Периодическая последовательность
прямоугольных импульсов

Из курса математического анализа известно [1, стр. 158], что система тригонометрических функций

(1)

с кратными частотами $\omega_n = \Delta \omega n$ , где $\Delta \omega = \frac{2\pi}{T}$ рад/с,

— целое число, образует ортонормированный базис для разложения периодических сигналов с периодом

, удовлетворяющих условиям Дирихле [1, стр. 165].

Иоганн Петер Густав Лежён-Дирихле
1805–1859

Условия Дирихле сходимости ряда Фурье требуют, чтобы периодический сигнал s(t)

был задан на сегменте $\left[ -\frac{T}{2}, \frac{T}{2} \right]$ , при этом удовлетворял следующим условиям:

Сигнал не должен иметь бесконечных значений.
Сигнал должен быть кусочно-непрерывным, т.е. иметь конечное число точек разрыва первого рода (скачки и устранимые разрывы).
Сигнал должен быть кусочно-монотонным и иметь конечное число экстремумов.

Например, периодическая функция $s(t) = \tg(t)$ не удовлетворяет условиям Дирихле, потому что функция $s(t) = \tg(t)$ имеет разрывы второго рода и принимает бесконечные значения при $t = \frac{\pi}{2} + k\pi$ , где — произвольное целое. Таким образом, функция $s(t) = \tg(t)$ не может быть представлена рядом Фурье. Также можно привести пример функции $s(t) = \sin\left(\frac{1}{\sin(t)}\right)$ , которая является ограниченной, но также не удовлетворяет условиям Дирихле, поскольку имеет бесконечное число точек экстремума при приближении к нулю. График функции $s(t) = \sin\left(\frac{1}{\sin(t)}\right)$ показан на рисунке 2.

График функции : а — два периода повторения; б — в окрестности

Рисунок 2. График функции $s(t) = \sin\left(\frac{1}{\sin(t)}\right)$ :
а — два периода повторения; б — в окрестности t = 0

На рисунке 2а показано два периода повторения функции $s(t) = \sin\left(\frac{1}{\sin(t)}\right)$ , а на рисунке 2б — область в окрестности t = 0 . Можно видеть, что при приближении к нулю, частота колебаний бесконечно возрастает, и такая функция не может быть представлена рядом Фурье, потому что она не является кусочно-монотонной.

Необходимо заметить, что на практике не бывает сигналов с бесконечными значениями тока или напряжения. Функции с бесконечным числом экстремумов типа $s(t) = \sin\left(\frac{1}{\sin(t)}\right)$ также в прикладных задачах не встречаются. Все реальные периодические сигналы удовлетворяют условиям Дирихле и могут быть представлены бесконечным тригонометрическим рядом Фурье вида:

(2)

В выражении (2) коэффициент a_0

задает постоянную составляющую периодического сигнала s(t)

Во всех точках, где сигнал непрерывен, ряд Фурье (2) сходится к значениям данного сигнала, а в точках разрыва первого рода — к среднему значению $\frac{1}{2}\big(s(t-0)+s(t+0)\big)$ , где s(t-0) и s(t+0) — пределы слева и справа от точки разрыва соответственно.

Также из курса математического анализа известно [2, стр. 500], что использование усеченного ряда Фурье, содержащего только первых членов вместо бесконечной суммы, приводит к приближенному представлению сигнала s(t) :

(3)

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

Рисунок 3. Приближение сигналов усеченным рядом Фурье:
а — прямоугольных импульсов; б — пилообразного сигнала

Ряд Фурье в комплексной форме

В предыдущем параграфе мы рассмотрели тригонометрический ряд Фурье для разложения произвольного периодического сигнала s(t) , удовлетворяющего условиям Дирихле. Применив формулу Эйлера, можно показать:

(4)

Тогда тригонометрический ряд Фурье (2) с учетом (4):

(5)

Таким образом, периодический сигнал s(t)

может быть представлен суммой постоянной составляющей $\frac{a_0}{2}$ и комплексных экспонент, вращающихся с частотами $\pm\omega_n$ с коэффициентами $\frac{a_n-jb_n}{2}$ для положительных частот $\omega_n$ , и $\frac{a_n+jb_n}{2}$ для комплексных экспонент, вращающихся с отрицательными частотами $-\omega_n$ .

Рассмотрим коэффициенты для комплексных экспонент, вращающихся с положительными частотами $\omega_n$ :

(6)

Аналогично, коэффициенты для комплексных экспонент, вращающихся с отрицательными частотами $-\omega_n$ :

(7)

Выражения (6) и (7) совпадают, кроме того постоянную составляющую $\frac{a_0}{2}$ также можно записать через комплексную экспоненту на нулевой частоте:

(8)

Таким образом, (5) с учетом (6)–(8) можно представить как единую сумму при индексации

от минус бесконечности до бесконечности:

(9)

Леонард Эйлер
1707–1783

Выражение (9) представляет собой ряд Фурье в комплексной форме. Коэффициенты ряда Фурье в комплексной форме $S(\omega_n)$ связаны с коэффициентами a_n

ряда в тригонометрической форме, и определяются как для положительных, так и для отрицательных частот $\omega_n$ . Индекс

в обозначении частоты $\omega_n$ указывает номер дискретной гармоники, причем отрицательные индексы соответствуют отрицательным частотам $\omega_n$ .

Из выражения (2) следует, что для вещественного сигнала s(t) коэффициенты a_n и b_n ряда (2) также являются вещественными. Однако (9) ставит в соответствие вещественному сигналу s(t) , набор комплексно-сопряженных коэффициентов $S(\omega_n)$ , относящихся как положительным, так и к отрицательным частотам $\omega_n$ .

Некоторые пояснения к ряду Фурье в комплексной форме

В предыдущем параграфе мы осуществили переход от тригонометрического ряда Фурье (2) к ряду Фурье в комплексной форме (9). В результате, вместо разложения периодических сигналов в базисе вещественных тригонометрических функций, мы получили разложение в базисе комплексных экспонент, с комплексными коэффициентами $S(\omega_n)$ , да еще и появились отрицательные частоты в разложении! Поскольку данный вопрос часто встречает непонимание, то необходимо дать некоторые пояснения.

Во-первых, работать с комплексными экспонентами в большинстве случаев проще, чем с тригонометрическими функциями. Например, при умножении и делении комплексных экспонент достаточно лишь сложить (вычесть) показатели, в то время как формулы умножения и деления тригонометрических функций более громоздкие.

Дифференцировать и интегрировать экспоненты, пусть даже комплексные, также проще, чем тригонометрические функции, которые постоянно меняются при дифференцировании и интегрировании (синус превращается в косинус и наоборот).

Если сигнал периодический и вещественный, то тригонометрический ряд Фурье (2) кажется более наглядным, потому что все коэффициенты разложения a_0 , a_n и b_n остаются вещественными. Однако, часто приходится иметь дело с комплексными периодическими сигналами (например, при модуляции и демодуляции используют квадратурное представление комплексной огибающей). В этом случае при использовании тригонометрического ряда Фурье все коэффициенты a_0 , a_n и b_n разложения (2) станут комплексными, в то время как при использовании ряда Фурье в комплексной форме (9) будет использованы одни и те же коэффициенты разложения $S(\omega_n)$ как для вещественных, так и для комплексных входных сигналов.

Ну и наконец, необходимо остановиться на пояснении отрицательных частот, которые появились в (9). Этот вопрос часто вызывает непонимание. В повседневной жизни мы не сталкиваемся с отрицательными частотами. Например, мы никогда не настраиваем свой радиоприемник на отрицательную частоту. Давайте рассмотрим следующую аналогию из механики. Пусть имеется механический пружинный маятник, который совершает свободные колебания с некоторой частотой $\omega_0$ . Может ли маятник колебаться с отрицательной частотой $\omega_0$ ? Конечно нет. Как не бывает радиостанций, выходящих в эфир на отрицательных частотах, так и частота колебаний маятника не может быть отрицательной. Но пружинный маятник — одномерный объект (маятник совершает колебания вдоль одной прямой).

Мы можем также привести еще одну аналогию из механики: колесо, вращающееся с частотой $\omega_0$ . Колесо, в отличие от маятника вращается, т.е. точка на поверхности колеса перемещается в плоскости, а не просто совершает колебания вдоль одной прямой. Поэтому для однозначного задания вращения колеса, задать частоту вращения недостаточно, потому что необходимо задать также направление вращения. Вот именно для этого мы и можем использовать знак частоты.

Так, если колесо вращается с угловой частотой $\omega_0$ рад/с против часовой стрелки, то считаем, что колесо вращается с положительной частотой, а если по направлению часовой стрелки, то частота вращения будет отрицательной. Таким образом, для задания вращения отрицательная частота перестает быть бессмыслицей и указывает направление вращения.

А теперь самое главное, что мы должны понять. Колебание одномерного объекта (например, пружинного маятника) может быть представлено как сумма вращений двух векторов, показанных на рисунке 4.

Рисунок 4. Колебание пружинного маятника
как сумма вращений двух векторов
на комплексной плоскости

Маятник совершает колебания вдоль вещественной оси комплексной плоскости с частотой $\omega_0$ по гармоническому закону $\cos(\omega_0 t)$ . Движение маятника показано горизонтальным вектором. Верхний вектор совершает вращения на комплексной плоскости с положительной частотой $\omega_0$ (против часовой стрелки), а нижний вектор вращается с отрицательной частотой $-\omega_0$ (по направлению часовой стрелки). Рисунок 4 наглядно иллюстрирует хорошо известное из курса тригонометрии соотношение:

(10)

Таким образом, ряд Фурье в комплексной форме (9) представляет периодические одномерные сигналы s(t)

как сумму векторов на комплексной плоскости, вращающихся с положительными и отрицательными частотами. При этом обратим внимание, что в случае вещественного сигнала согласно (9) коэффициенты разложения $S(\omega_n)$ для отрицательных частот $\omega_n < 0$ являются комплексно-сопряженными соответствующим коэффициентам для положительных частот $\omega_n$ . В случае комплексного сигнала это свойство коэффициентов не выполняется ввиду того, что a_n

также являются комплексными.

Спектр периодических сигналов

Ряд Фурье в комплексной форме представляет собой разложение периодического сигнала s(t) в сумму комплексных экспонент, вращающихся с положительными и отрицательными частотами кратными $\Delta\omega=\frac{2\pi}{T}$ рад/c с соответствующими комплексными коэффициентами $S(\omega_n)$ , которые определяют спектр сигнала s(t) . Комплексные коэффициенты $S(\omega_n)$ могут быть представлены по формуле Эйлера как $S(\omega_n) = |S(\omega_n)|\exp\big(j\Phi(\omega_n)\big)$ , где $|S(\omega_n)|$ — амплитудный спектр, a $\Phi(\omega_n)$ — фазовый спектр.

Поскольку периодические сигналы раскладываются в ряд только на фиксированной сетке частот $\omega_n = \Delta\omega n$ , то спектр $S(\omega_n)$ периодических сигналов является линейчатым (дискретным).

Рисунок 5. Спектр периодической последовательности
прямоугольных импульсов:
а — амплитудный спектр; б — фазовый спектр

На рисунке 5 приведен пример амплитудного $|S(\omega_n)|$ и фазового $\Phi(\omega_n)$ спектра периодической последовательности прямоугольных импульсов (см. рисунок 1) при T = 4 с, длительности импульса $\tau = 2$ c и амплитуде импульсов A = 2 В.

Амплитудный спектр исходного вещественного сигнала является симметричным относительно нулевой частоты, а фазовый спектр — антисимметричным. При этом заметим, что значения фазового спектра $\Phi(\omega_n) = \pi$ и $\Phi(\omega_n) = -\pi$ соответствуют одной и той же точке комплексной плоскости $\exp(\pm j\pi) = -1$ .

Можно сделать вывод, что все коэффициенты разложения приведенного сигнала являются чисто вещественными, и фазовый спектр $\Phi(\omega_n) = \pm \pi$ соответствует отрицательным коэффициентам $S(\omega_n)$ .

Обратим внимание, что размерность амплитудного спектра $|S(\omega_n)|$ совпадает с размерностью сигнала s(t) . Если s(t) описывает изменение напряжения во времени, измеряемое в вольт, то амплитуды гармоник спектра $|S(\omega_n)|$ также будут иметь размерность вольт.

Выводы

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

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

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

Программная реализация в библиотеке DSPL

Данные для построения рисунков данного раздела были просчитаны при использовании библиотеки DSPL-2.0

Ниже приведён исходный код программы расчета данных для построения рисунка 2:

fourier_series_dirichlet_ex.c


/*******************************************************************************
 * Данная программа производит расчет данных для построения графика функции
 * 
 *                          s(t) = sin(1/sin(t))
 * 
 * Расчет функции производится в диапазоне t от -pi до pi и от -0.2 до 0.2.
 * 
 * Сохраненные данные:
 * 
 * dat/fourier_series_signal0.txt - данные функции при t от -2*pi до 2*pi; 
 * dat/fourier_series_signal1.txt - данные функции при t от -0.2  до 0.2.
 *  
 ******************************************************************************/ 

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "dspl.h"


/* размер массивов */
#define N  10000

int main(int argc, char* argv[])
{
  double t[N]; /* массив аргументов функции */
  double s[N]; /* массив значений функции s(t) = sin(1/sin(t)) */
  void* hdspl; /* DSPL handle        */
  void* hplot; /* GNUPLOT handle     */
  int k;
  
  /* Загружаем libdspl.dll */
  hdspl = dspl_load();
  if(!hdspl)
  {
    printf("Cannot to load libdspl!\n");
    return 0;
  }

  /* заполняем массив t от -2*pi до 2*pi */
  linspace(-M_2PI, M_2PI, N, DSPL_PERIODIC, t);
  
  /* рассчитываем функцию s(t) = sin(1/sin(t)) 
   * для каждого значения переменной t[k]      */
  for(k = 0; k < N; k++)
    s[k] = sin(1.0/sin(t[k]));
    
  /* Сохраняем данные в файл dat/fourier_series_signal0.txt */  
  writetxt(t, s, N, "dat/fourier_series_signal0.txt");


  /* заполняем массив t от -0.2 до 0.2 */
  linspace(-0.2, 0.2, N, DSPL_PERIODIC, t);
  
  /* рассчитываем функцию s(t) = sin(1/sin(t)) 
   * для каждого значения переменной t[k]      */
  for(k = 0; k < N; k++)
    s[k] = sin(1.0/sin(t[k]));

  /* Сохраняем данные в файл dat/fourier_series_signal1.txt */  
  writetxt(t, s, N, "dat/fourier_series_signal1.txt");
  
  /* plotting by GNUPLOT */
  gnuplot_create(argc, argv, 560, 380, "img/dirichlet_ex.png", &hplot);
  gnuplot_cmd(hplot, "set grid");
  gnuplot_cmd(hplot, "unset key");
  gnuplot_cmd(hplot, "set xlabel 'x'");
  gnuplot_cmd(hplot, "set ylabel 'sin(1/sin(x))'");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_signal1.txt' with lines");
  gnuplot_close(hplot);
  
  /* Выгружаем libdspl.dll из памяти */
  dspl_free(hdspl);
  

  return 0;
}

Исходный код программы расчета данных для построения графиков приближения периодической последовательности прямоугольных импульсов и пилообразного сигнала усеченным рядом Фурье (рисунок 3):

fourier_series_rec.c


/*******************************************************************************
 * Данная программа производит расчет приближенного представления 
 * периодической последовательности прямоугольных импульсов и 
 * пилообразного сигнала усеченным рядом Фурье. 
 * 
 * Сохраненные данные:
 * 
 * dat/fourier_series_pimp0.txt - Последовательность прямоугольных импульсов
 *  
 * dat/fourier_series_pimp_rec_5.txt -  восстановленная последовательность 
 *                                      прямоугольных импульсов при 
 *                                      использовании 2 гармоник
 * 
 * dat/fourier_series_pimp_rec_9.txt -  восстановленная последовательность 
 *                                      прямоугольных импульсов при 
 *                                      использовании 4 гармоник
 *
 * dat/fourier_series_pimp_rec_21.txt - восстановленная последовательность 
 *                                      прямоугольных импульсов при 
 *                                      использовании 10 гармоник
 * 
 * dat/fourier_series_pimp_rec_61.txt - восстановленная последовательность 
 *                                      прямоугольных импульсов при 
 *                                      использовании 10 гармоник
 * 
 * dat/fourier_series_saw0.txt      -   Исходный пилообразный сигнал
 * 
 * dat/fourier_series_saw_rec_5.txt -   восстановленная пилообразный сигнал 
 *                                      при использовании 2 гармоник
 *                                      
 * dat/fourier_series_saw_rec_9.txt -   восстановленная пилообразный сигнал 
 *                                      при использовании 4 гармоник
 * 
 * dat/fourier_series_saw_rec_21.txt -  восстановленная пилообразный сигнал 
 *                                      при использовании 10 гармоник
 * 
 * dat/fourier_series_saw_rec_61.txt -  восстановленная пилообразный сигнал 
 *                                      при использовании 30 гармоник
 * 
 *  
 ******************************************************************************/ 

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "dspl.h"

/* размер массивов исходных и восстановленных сигналов */
#define N  1000

/* Период */
#define T  2

/* Длительность импульса */
#define TAU 1

/* Максимальное количество коэффициентов ряда Фурье */
#define MAX_M 61

/*  Количество периодов */
#define P 3.0

int main(int argc, char* argv[])
{
  double t[N];         /* время                      */
  double s[N];         /* исходный сигнал            */
  double x[N];         /* восстановленный сигнал     */
  double w[MAX_M];     /* массив частоты             */
  complex_t S[MAX_M];  /* Спектр                     */
  
  complex_t xc[N];     /* восстановленный по спектру */
                       /* комплексный сигнал         */
   
  /*
   * Количество спектральных гармоник усеченного ряда
   * Заметим, что 5 гармоник усеченного ряда содержат
   * две пары комплексно-сопряженных спектральных компонент
   * и одну постоянную составляющую.
   * Это означает, что ряд из 5 гармоник в комплексной форме соответствует
   * двум значениям a_n и b_n ряда в тригонометрической форме.
   * Аналогично M = 21 соответствует 10 значениям a_n и b_n ряда
   * в тригонометрической форме.
   */ 
  int M[4] = {5, 9, 21, MAX_M};

  void* hdspl; /* DSPL handle        */
  void* hplot; /* GNUPLOT handle     */
  
  
  int k, n;
  char fname[64];  /* имя файла для сохранения рассчитанных данных */


  /* Загружаем libdspl.dll */
  hdspl = dspl_load();
  if(!hdspl)
  {
    printf("Cannot to load libdspl!\n");
    return 0;
  }
  
  /* Заполняю вектор времени для P периодов повторения
   * периодических сигналов
   */
  linspace(-T*0.5*P, T*0.5*P, N, DSPL_PERIODIC, t);
  
  /*--------------------------------------------------------------------------*/
  /* последовательность прямоугольных импульсов */
  signal_pimp(t, N, 1.0, TAU, 0.0, T, s);
  
  /* сохраняем в файл dat/fourier_series_pimp0.txt */
  writetxt(t, s, N, "dat/fourier_series_pimp0.txt");
  
  /* цикл по разному количеству гармоник усеченного ряда */
  for(k = 0; k < 4; k++)
  {
    /* расчет спектра для текущего усеченного ряда */
    fourier_series_dec(t, s, N, T, M[k], w, S);
    
    /* нормировка потому что P периодов в исходном сигнале */
    for(n = 0; n < M[k]; n++)
    {
      RE(S[n]) /= P;
      IM(S[n]) /= P;
    }
  
    /* восстанавливаю сигнал по усеченному ряду */
    fourier_series_rec(w, S, M[k], t, N, xc);
  
    /* Комплексный восстановленный сигнал имеет 
     * мнимую часть на уровне машинной точности */
    cmplx2re(xc, N, x, NULL);
  
    /* сохраняю в файл для последующего построения графика */
    sprintf(fname, "dat/fourier_series_pimp_rec_%d.txt", M[k]);
    writetxt(t, x, N, fname);
  }
  
  /*--------------------------------------------------------------------------*/
  /* Пилообразный сигнал */
  signal_saw(t, N, 0.5, 0.0, T, s);
  for(n = 0; n < N; n++)
    s[n] += 0.5;
  
  /* сохраняем в файл dat/fourier_series_saw0.txt */
  writetxt(t, s, N, "dat/fourier_series_saw0.txt");
  
  /* цикл по разному количеству гармоник усеченного ряда */
  for(k = 0; k < 4; k++)
  {
    /* расчет спектра для текущего усеченного ряда */
    fourier_series_dec(t, s, N, T, M[k], w, S);
    
    /* нормировка потому что P периодов в исходном сигнале */
    for(n = 0; n < M[k]; n++)
    {
      RE(S[n]) /= P;
      IM(S[n]) /= P;
    }
    
    /* восстанавливаю сигнал по усеченному ряду */
    fourier_series_rec(w, S, M[k], t, N, xc);
  
    /* Комплексный восстановленный сигнал имеет 
     * мнимую часть на уровне машинной точности */
    cmplx2re(xc, N, x, NULL);
  
    /* сохраняю в файл для последующего построения графика */
    sprintf(fname, "dat/fourier_series_saw_rec_%d.txt", M[k]);
    writetxt(t, x, N, fname);
  }
  
  
  /* Построение графиков пакетом GNUPLOT */
  gnuplot_create(argc, argv, 800, 640,  "img/fourier_series_rec.png", &hplot);
  gnuplot_cmd(hplot, "unset key");
  gnuplot_cmd(hplot, "set multiplot layout 4,2 rowsfirst");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_pimp0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_pimp_rec_5.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_saw0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_saw_rec_5.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_pimp0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_pimp_rec_9.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_saw0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_saw_rec_9.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_pimp0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_pimp_rec_21.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_saw0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_saw_rec_21.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_pimp0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_pimp_rec_61.txt' with lines");
  gnuplot_cmd(hplot, "plot 'dat/fourier_series_saw0.txt' with lines,\\");
  gnuplot_cmd(hplot, "     'dat/fourier_series_saw_rec_61.txt' with lines");
  gnuplot_cmd(hplot, "unset multiplot");
  gnuplot_close(hplot);


  /* Выгружаем libdspl.dll из памяти */
  dspl_free(hdspl);

  return 0;
}

Исходный код программы расчета амплитудного и фазового спектра (рисунок 5):

fourier_series_pimp_spectrum.c


/*******************************************************************************
 * Данная программа производит расчет амплитудного и фазвого спектра 
 * последовательности прямоугольных импульсов
 * 
 * Сохраненные данные:
 * 
 * dat/fourier_series_pimp_mag.txt - амплитудный спектр 
 * dat/fourier_series_pimp_phi.txt - фазовый спектр
 *  
 ******************************************************************************/ 

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "dspl.h"

/* количество точек прямоугольного 
 * импульса на одном периоде повторения */
#define N  1000

/* Период повторения сигнала */
#define T  4

/* Длительность импульса */
#define TAU 2

/* Количество гармоник спектра */
#define K 41

int main(int argc, char* argv[])
{
  double t[N];    /* время              */
  double s[N];    /* исходный сигнал    */
  double w[K];    /* массив частоты     */
  complex_t S[K]; /* спектр             */
  double mag[K];  /* амплитудный спектр */
  double phi[K];  /* фазовый спектр     */
  void* hdspl;    /* DSPL handle        */
  void* hplot;    /* GNUPLOT handle     */
  int n;


  /* Загружаем libdspl.dll */
  hdspl = dspl_load();
  if(!hdspl)
  {
    printf("cannot to load libdspl!\n");
    return 0;
  }

  /* заполняем массив времени для одного периода сигнала */
  linspace(-T*0.5, T*0.5, N, DSPL_SYMMETRIC, t);

  /* сигнал */
  signal_pimp(t, N, 2.0, TAU, 0.0, T, s);

  /* расчет комплексного спектра */
  fourier_series_dec(t, s, N, T, K, w, S);

  /* Амплитудный и фазовый спектры */
  for(n = 0; n < K; n++)
  {
    mag[n] = ABS(S[n]);
    phi[n] = atan2(IM(S[n]), RE(S[n]));
  }

  /* Сохранение амплитудного спектра в файл */
  writetxt(w, mag, K, "dat/fourier_series_pimp_mag.txt");

  /* Сохранение фазового спектра в файл */
  writetxt(w, phi, K, "dat/fourier_series_pimp_phi.txt");

  /* Построение графиков пакетом GNUPLOT */
  gnuplot_create(argc, argv, 560, 380, "img/spectrum.png", &hplot);
  gnuplot_cmd(hplot, "unset key");
  gnuplot_cmd(hplot, "set grid");
  gnuplot_cmd(hplot, "set lmargin 8");
  gnuplot_cmd(hplot, "set multiplot layout 2,1 rowsfirst");
  gnuplot_cmd(hplot, "set xlabel 'Частота, рад/с'");
  gnuplot_cmd(hplot, "set ylabel 'Амплитудный спектр'");
  gnuplot_cmd(hplot, "plot[-10*pi:10*pi] 'dat/fourier_series_pimp_mag.txt' with  impulses lt 1 ,\\");
  gnuplot_cmd(hplot, "'dat/fourier_series_pimp_mag.txt' with points pt 7 ps 0.5 lt 1");
  gnuplot_cmd(hplot, "set ylabel 'Фазовый спектр'");
  gnuplot_cmd(hplot, "plot[-10*pi:10*pi] 'dat/fourier_series_pimp_phi.txt' with  impulses lt 1 ,\\");
  gnuplot_cmd(hplot, "'dat/fourier_series_pimp_phi.txt' with points pt 7 ps 0.5 lt 1");
  gnuplot_cmd(hplot, "unset multiplot");
  gnuplot_close(hplot);
  
  /* Выгружаем libdspl.dll из памяти */
  dspl_free(hdspl);

  return 0;

}

Смотри также

Некоторые свойства разложения периодических сигналов в ряд Фурье
Спектр периодической последовательности прямоугольных импульсов
Преобразование Фурье непериодических сигналов

Список литературы

[1] Воробьев Н.Н. Теория рядов. Москва, Наука, Главная редакция физико-математической литературы, 1979, 408 с.

[2] Каханер Д., Моулер К., Нэш С. Численные методы и математическое обеспечение. Москва, Мир, 1998.

[3] Баскаков, С.И. Радиотехнические цепи и сигналы. Москва, ЛЕНАНД, 2016, 528 c. ISBN 978-5-9710-2464-4

[4] Гоноровский И.С. Радиотехнические цепи и сигналы Москва, Советское радио, 1977, 608 c.

[5] Дёч, Г. Руководство по практическому применению преобразования Лапласа. Москва, Наука, 1965, 288 c.

[6] Bracewell R. The Fourier Transform and Its Applications McGraw-Hills, 1986, 474 c. ISBN 0-07-007-015-6

Последнее изменение страницы: 12.05.2022 (19:42:33)

Страница создана Latex to HTML translator ver. 5.20.11.14