МІНІСТЕРСТВО ОСВІТИ І НАУКИ УКРАЇНИ
НАЦІОНАЛЬНИЙ УНІВЕРСИТЕТ “ЛЬВІВСЬКА ПОЛІТЕХНІКА”
Кафедра ЕОМ
Розрахунково-графічна робота
з дисципліни «Цифрова обробка сигналів»
на тему: «ДОСЛІДЖЕННЯ ДИСКРЕТНОГО ПЕРЕТВОРЕННЯ ФУР’Є ТА ЙОГО ВЛАСТИВОСТЕЙ
»
Мета роботи: освоїти методику обчислення спектральних характеристик сигналу за допомогою перетворення Фур’є.
Спектральний аналіз періодичних сигналів
Нехай - періодичний сигнал з періодом Т. Якщо він описується неперервною або кусочно-неперервною функцією, то його можна подати у вигляді ряду Фур'є
, (1)
, (2)
де - основна гармоніка.
Коефіцієнти називають частотним спектром,
- амплітудним спектром,
- фазовим спектром.
Нехай , . Тобто, пряме ДПФ наближає (за формулою чисельного інтегрування прямокутників) коефіцієнти розкладу сигналу в ряд Фур'є: , .
Рівність Парсеваля тут має вигляд:
.
При збільшенні , зменшенні , похибка (методу) такого представлення зменшується. Подібне має місце і при оберненому перетворенні, тобто наближенні сигналу відрізком ряду Фур'є, коли обмежуються границі підсумування в формулі (1).
Поняття нормованої частоти
Нормованою частотою називається відношення поточної частоти до частоти дискретизації. Тобто:
Таким чином, залежно від обраної шкали частот, основна смуга відповідає областям: ; ; ; .
Інверсія спектру дійсного сигналу
Мета інверсії спектральних складових в тому, що в основній смузі частот довільна складова спектру (тут мається на увазі нормована частота) повинна виявитися на «протилежній» частоті - (відбитися). При цьому у аргументу спектра додатково змінюється знак. Це ілюструється рисунком, наведеним далі.
Операція інверсії спектру реалізовується за допомогою зсуву вправо на частоту .
ЗАВДАННЯ
Отримати аналітичні вирази для знаходження частотного спектру заданого варіантом сигналу (рівність (2)).
Обрахувати 256 спектральних коефіцієнтів та привести таблицю їх значень.
Привести графік зміни спектру для заданої кількості коефіцієнтів.
Порівняти (в табличному та графічному вигляді) отримані значення з результатами, обрахованими в лабораторній роботі №3.
Пояснити різницю між коефіцієнтами ряду Фур’є та дискретного перетворення Фур’є.
Відтворити вхідну послідовність за допомогою наближення рядом Фур’є (рівність (1)) для 256 коефіцієнтів та привести її графік.
Порівняти вхідну та наближену послідовності.
Обрахувати амплітудну, фазову характеристики та привести їх графіки.
Знайти енергію сигналу за формулою Парсеваля.
Здійснити операцію інверсії спектру.
За допомогою оберненого ПФ знайти сигнал, що відповідає оберненому спектру та порівняти його із заданим.
Зробити висновки про спектральні характеристики заданого сигналу.
Завдання
Варіант 26
Сигнал задано наступним чином:
Аналітичні вирази для знаходження частотного спектру заданого варіантом сигналу:
Виходячи із заданого графіка, сигнал описується фінітною функцією, тривалістю 15с. Перші 9с сигнал має значення -4, наступні 6с сигнал має значення 8. Таким чином, маємо функцію, що описує заданий сигнал:
;
Для знаходження частотного спектру слід підставити аналітичний вираз вхідного сигналу у формулу:
,
, де - основна гармоніка, де Т=15.
Виходячи з завдання, маємо: T=15; . Тому :
Нехай , тоді
256 спектральних коефіцієнтів, таблиця значень:
Для знаходження спектральних характеристик в дискретних точках слід створити програму для виконання у SCILAB, використавши отримані аналітичні вирази. Текст програми приведено нижче.
Таблиця значень спектральних коефіцієнтів.
Графік зміни спектру для заданої кількості коефіцієнтів:
На рисунку приведено графіки дійсної та уявної частини комплексного частотного спектру заданого сигналу.
Отримані результати, обраховані в лабораторній роботі №3:
Таблиця значень спектральних коефіцієнтів.
На рисунках приведено графіки дійсних та уявних частин спектру заданого сигналу.
Рис.2. Порівняння реальних частин
Рис.3. Порівняння уявних частин
Різниця між коефіцієнтами ряду Фур’є та дискретного перетворення Фур’є:
Різниця між коефіцієнтами ряду Фур’є та дискретного перетворення Фур’є,
зумовлена наступним:
Нехай , . Тобто, пряме ДПФ наближає (за формулою чисельного інтегрування прямокутників) коефіцієнти розкладу сигналу в ряд Фур'є: , .
При збільшенні , зменшенні , похибка (методу) такого представлення зменшується. Подібне має місце і при оберненому перетворенні, тобто наближенні сигналу відрізком ряду Фур'є, коли обмежуються границі підсумування в формулі (1).
Вхідна послідовність за допомогою наближення рядом Фур’є для 256 коефіцієнтів та її графік:
Таблиця вхідної послідовності отриманої за допомогою наближення рядом Фур’є.
Порівняння вхідної та наближеної послідовностей.
Рис.4. Порівняння вхідної та наближеної послідовностей
Амплітудна та фазова характеристики та їх графіки:
Рис.5. Амплітудна та фазова характеристики
Енергія сигналу за формулою Парсеваля:
Під теоремою Парсеваля звичайно розуміють унітарність перетворення Фур'є. Тобто сума (або інтеграл) квадрата функції дорівнює (дорівнює) сумі (або інтегралу) квадрата результату перетворення. Слід зауважити, що загальний вигляд теореми Парсеваля часто називають теоремою Планшереля або Узагальненої формулою Релея. Теорема була доведена для рядів Марком-Антуаном Парсевалем в 1799 і була пізніше застосована до ряду Фур'є.
Запис теореми має вигляд:
де позначає неперервне перетворення Фур'є, яке пов'язує часовий або просторовий сигнал з його поданням в частотній області .
У дискретному вигляді теорему записують таким чином:
,
де являє собою дискретне перетворення Фур'є сигналу , що має N відліків. Теорема Парсеваля встановлює рівність між енергією сигналу і енергією його спектру.
Знайдемо енергію сигналу за формулою Парсеваля.
Операція інверсії спектру:
Рис.6. Інвертований спектр.
ОПФ шуканого сигналу, що відповідає оберненому спектру та порівняня його із заданим:
Рис.7. Порівняння заданого та інвертованого сигналів.
Висновки:
При виконанні даної розрахункової роботи було вивчено методи обчислення спектральних характеристик сигналу за допомогою перетворення Фур’є, а також вдосконалено навики роботи з середовищем SciLab. Дослідження спектральних характеристик сигналу проводилися двома шляхами - аналітичним та практичним, з подальшим порівнянням результатів обох методів.
Код програми:
clc
clear
close,close,close,close,close,close
//Завдання 1,2
//Аналітичні вирази для знаходження частотного спектру
N = 256;clc
clear
close,close,close,close,close,close
//Завдання 1,2
//Аналітичні вирази для знаходження частотного спектру
N = 256;
w=1:1:N;
for k=1:1:255
A =-2*(%pi*%i*(k)/(15));
c(k) = (1/15)*( - (12/A)*exp(6*A) + (12/A)*exp(12*A) - (4/A)*exp(15*A) + 4/A );
end;
c0 = 0.5333;
c=[ c0; c]
//Завдання 6
//Спектральний аналіз, обчислення ряду Фурє
T = 15;
t0 = -3;
dt=T/N;
t=t0:dt:12-dt;
for n=1:1:256
s = 0;
for k=1:1:255
A =2*(%pi*%i*(k-1)/(15));;
s = s + (c(k) * exp(A*t(n)));
end;
x(n) = s;//відтворена послідовність
end;
//Завдання 7
//Вигляд вхідного сигналу, та обчислення його ШПФ
T=15;m=8;N=2^m;
dt=T/N;
t11=-3:dt:6;
t12=6:dt:12-dt;
t=[t11,t12];
x11=0*t11-4;
x12=0*t12+8;
xx=[x11, x12];
w=1:1:N;
Sx = fft(xx)/(N);//ШПФ
//Завдання 7
//Знаходження абсолютного значення(порівняння)
for k = 1:1:256
b(k) = abs((x(k))-(xx(k)));
end;
//Завдання 8
//Знаходження амплітудного спектру
for k = 1:1:256
Ak(k) = sqrt(imag(c(k))*imag(c(k)) + real(c(k))*real(c(k))); //амплітудний спектр
end;
//Знаходження фазового спектру
for k = 1:1:256
Fk(k) = -atan(imag(c(k))/real(c(k)));
end;
for k = 1:1:256
Ee(k) = real(c(k))*real(c(k));
end;
//Завдання 9
//Знаходження енергії
E = sum(Ak) * sum(Ak);//Енергія
//Завдання 3
figure(1,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, real(c)), title('real(c)','position',[N/2,max(real(c))/2]);
subplot(2,1,2), plot2d(w, imag(c)), title('imag(c)','position',[N/2,max(imag(c))/2]);
//Завдання 4
figure(2,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, real(Sx),);
title('real(Sx)')
subplot(2,1,2), plot2d(w, real(c),);
title('real(c)')
figure(7,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, imag(Sx),);
title('imag(Sx)')
subplot(2,1,2), plot2d(w, imag(c),);
title('imag(c)')
//Завдання 6,7
figure(3,'BackgroundColor',[1,1,1]);
subplot(3,1,1), plot2d(t, (xx),3), title('Input data - x(n)','position',[T/2,max(real(x))/2]);
subplot(3,1,2), plot2d(t, (x),4), title('Reprodused data - x(n)','position',[T/2,max(real(x))/2]);
subplot(3,1,3), plot2d(t, b,), title('Comparison data - x(n)','position',[T/2,max(real(x))/2]);
//Завдання 8
figure(4,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, Ak,), title('Ak','position',[T/2,max(real(x))/2]);
subplot(2,1,2), plot2d(w, Fk,), title('Fk','position',[T/2,max(real(x))/2]);
//Завдання 10
for k = 1:1:255
newC(k) = c(256-k);
end
newC=[newC;c0]
qwerty=[newC;c]
k = 1:1:512
figure(5,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(k, real(qwerty)), title('real(newC)','position',[N/2,max(real(qwerty))/2]);
subplot(2,1,2), plot2d(k, imag(qwerty)), title('imag(newC)','position',[N/2,max(imag(qwerty))/2]);
T = 15;
t0 = -3;
dt=T/N;
qw = N*ifft(qwerty);//OШПФ
//Завдання 11
k=1:1:512
figure(6,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(t,(xx),3), title('Input data - x(n)','position',[N/2,max(real(c))/2]);
subplot(2,1,2), plot2d(k, (qw),), title('Reprodused data - x(n)','position',[T/2,max(real(qw))/2]);
w=1:1:N;
for k=1:1:255
A =-2*(%pi*%i*(k)/(15));
c(k) = (1/15)*( - (12/A)*exp(6*A) + (12/A)*exp(12*A) - (4/A)*exp(15*A) + 4/A );
end;
c0 = 0.5333;
c=[ c0; c]
//Завдання 6
//Спектральний аналіз, обчислення ряду Фурє
T = 15;
t0 = -3;
dt=T/N;
t=t0:dt:12-dt;
for n=1:1:256
s = 0;
for k=1:1:255
A =2*(%pi*%i*(k-1)/(15));;
s = s + (c(k) * exp(A*t(n)));
end;
x(n) = s;//відтворена послідовність
end;
//Завдання 7
//Вигляд вхідного сигналу, та обчислення його ШПФ
T=15;m=8;N=2^m;
dt=T/N;
t11=-3:dt:6;
t12=6:dt:12-dt;
t=[t11,t12];
x11=0*t11-4;
x12=0*t12+8;
xx=[x11, x12];
w=1:1:N;
Sx = fft(xx)/(N);//ШПФ
//Завдання 7
//Знаходження абсолютного значення(порівняння)
for k = 1:1:256
b(k) = abs((x(k))-(xx(k)));
end;
//Завдання 8
//Знаходження амплітудного спектру
for k = 1:1:256
Ak(k) = sqrt(imag(c(k))*imag(c(k)) + real(c(k))*real(c(k))); //амплітудний спектр
end;
//Знаходження фазового спектру
for k = 1:1:256
Fk(k) = -atan(imag(c(k))/real(c(k)));
end;
for k = 1:1:256
Ee(k) = real(c(k))*real(c(k));
end;
//Завдання 9
//Знаходження енергії
E = sum(Ak) * sum(Ak);//Енергія
//Завдання 3
figure(1,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, real(c)), title('real(c)','position',[N/2,max(real(c))/2]);
subplot(2,1,2), plot2d(w, imag(c)), title('imag(c)','position',[N/2,max(imag(c))/2]);
//Завдання 4
figure(2,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, real(Sx),);
title('real(Sx)')
subplot(2,1,2), plot2d(w, real(c),);
title('real(c)')
figure(7,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, imag(Sx),);
title('imag(Sx)')
subplot(2,1,2), plot2d(w, imag(c),);
title('imag(c)')
//Завдання 6,7
figure(3,'BackgroundColor',[1,1,1]);
subplot(3,1,1), plot2d(t, (xx),3), title('Input data - x(n)','position',[T/2,max(real(x))/2]);
subplot(3,1,2), plot2d(t, (x),4), title('Reprodused data - x(n)','position',[T/2,max(real(x))/2]);
subplot(3,1,3), plot2d(t, b,), title('Comparison data - x(n)','position',[T/2,max(real(x))/2]);
//Завдання 8
figure(4,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(w, Ak,), title('Ak','position',[T/2,max(real(x))/2]);
subplot(2,1,2), plot2d(w, Fk,), title('Fk','position',[T/2,max(real(x))/2]);
//Завдання 10
for k = 1:1:255
newC(k) = c(256-k);
end
newC=[newC;c0]
qwerty=[newC;c]
k = 1:1:512
figure(5,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(k, real(qwerty)), title('real(newC)','position',[N/2,max(real(qwerty))/2]);
subplot(2,1,2), plot2d(k, imag(qwerty)), title('imag(newC)','position',[N/2,max(imag(qwerty))/2]);
T = 15;
t0 = -3;
dt=T/N;
qw = N*ifft(qwerty);//OШПФ
//Завдання 11
k=1:1:512
figure(6,'BackgroundColor',[1,1,1]);
subplot(2,1,1), plot2d(t,(xx),3), title('Input data - x(n)','position',[N/2,max(real(c))/2]);
subplot(2,1,2), plot2d(k, (qw),), title('Reprodused data - x(n)','position',[T/2,max(real(qw))/2]);