Міністерство освіти і науки України
Національний університет «Львівська політехніка»
Кафедра ЕОМ
Лабораторна робота №3
«Швидкі алгоритми обчислення дискретних тригонометричних перетворень»
Варіант: 2
1.Мета роботи: дослідити швидкі алгоритми дискретних тригонометричних перетворень і порівняти їх з алгоритмами безпосереднього обчислення тригонометричних перетворень.
2.Загальні відомості:
Алгоритми швидкого перетворення Фур’є (ШПФ) можуть бути отримані за допомогою послідовного застосування операції розкладу одномірного масиву вхідних відліків сигналу на двохмірний. Ця операція здійснена тільки у випадку, коли N (довжина послідовності) є складним числом (N = N1 ( N2 ( ... ( Nj). Якщо N просте, його неможливо розкласти на прості співмножники; для такої послідовності алгоритмів ШПФ не існує. В більшості практичних випадків вхідну послідовність штучно продовжують додаванням нульових відліків до отримання N як складного числа.
Для характеристики розкладу використовують поняття “основа”, якщо всі співмножники однакові (N1 = N2 = ... = Nj) і “рощеплена основа”, якщо співмножники неоднакові.
Приклад. N = 64 = 2 ( 2 ( 2 ( 2 ( 2 ( 2 - основа 2.
N = 64 = 4 ( 4 ( 4 - основа 4.
N = 128 = 2 ( 4 ( 4 ( 4 - рощеплена основа 2-4.
Дискретне перетворення Фур’є
Дискретне перетворення Фур’є (ДПФ) кінцевої послідовності {x(n)}, 0 ( n ( N-1 визначається як
, k = 0, 1, ..., N-1 (1)
де .
Суть алгоритмів ШПФ в тому, що якщо N складне число і є степенем двійки (N=2m), то вихідна послідовність розбивається на дві коротші послідовності, ДПФ яких можуть бути скомбіновані таким чином, щоб утворилось ДПФ вихідної послідовності.
Методика побудови алгоритмів ШПФ наступна. Введемо дві (N/2) - точкові послідовності {x1(n)} і {x2(n)} з парних і непарних членів x(n) відповідно, x1(n) = x(2n) і x2(n) = x(2n + 1), n = 0, 1, ..., N/2-1. Тоді формулу розкладу можна записати так:
де X1(k) і X2(k) є N/2 - точкові ДПФ парних і непарних відліків вихідної послідовності.
Застосовуючи цю формулу розкладу X1(k) і X2(k) до пори, поки X1(k) і X2(k) не стануть двохточковим ДПФ, отримаємо алгоритм ШПФ з прорідженням в часі.
Для отримання іншої розповсюдженої форми алгоритмів ШПФ вихідну послідовність розбивають на дві (N/2) - точкові послідовності таким чином: перша складається з перших N/2 відліків, а друга з решти, x1(n) = x(n) і x2(n) = x(n + N/2), n = 0, 1, ..., N/2-1.
При такому розбитті N - точкове ДПФ можна записати так
де X1(k) і X2(k) є N/2 - точкові ДПФ першої та другої половини відліків вихідної послідовності.
Алгоритми ШПФ, що отримані з застосуванням цієї методики називаються алгоритмами з прорідженням по частоті.
3.Порядок виконання роботи:
1. Застосовуючи методики розкладу, отримати алгоритм обчислення заданого перетворення за певними “основою” і розмірністю.
2. Скласти процедуру на мові високого рівня для безпосереднього (прямого) обчислення тригонометричного перетворення.
3. Скласти процедуру на мові високого рівня для обчислення швидкого перетворення по алгоритму отриманому в п.1.
4. Виміряти часи виконання процедур п.2 і п.3.
5. Порівняти часи виконання процедур п.2 і п.3, і пояснити отримані результати.
4.Завдання:
Швидке перетворення Фур`є за основою 4 над 256-точками
4.Моделювання функції:
Параметри функції:
Форма сигнала – синус; Частота – 20; Амплітуда – 5; Фаза - 0; Інтервал виборок 0.001
Текст програми на мові С:
Main.cpp
//-------------------------------------------------------------------
#include <vcl.h>
#pragma hdrstop
#include "Main.h"
#include <fstream.h>
#include <math.h>
//#include <SysUtils.h>
//---------------------------------------------------------------------------
#pragma package(smart_init)
#pragma resource "*.dfm"
TForm1 *Form1;
double a[512];
//---------------------------------------------------------------------------
__fastcall TForm1::TForm1(TComponent* Owner)
: TForm(Owner)
{ FFT();
initall();
fft4();
}
//---------------------------------------------------------------------------void TForm1::fastfouriertransform(real& a, int nn, bool inversefft)
void TForm1::FFT()
{ double temp;
// Відкриваю вхідний файл
ifstream file1;
file1.open("sin3.xy");
// Читаю дані в масив
for(int kk=0;kk<512;kk++)
{ file1>>temp;
file1>>a[kk];
}
file1.close();
// Виводжу дані у таблицю
ValueListEditor1->ColCount=2;
String str1,str2;
for(int tt=0;tt<512;tt=tt+2)
{ str1=FloatToStrF(a[tt],ffGeneral,7,10);
str2=FloatToStrF(a[tt+1],ffGeneral,7,10);
ValueListEditor1->InsertRow(str1,str2,true);
}
}
//---------------------------------------------------------------------------
void TForm1::initall()
{
int i,k,l;
//m_x init
for (i = 0; i < N; i++)
{ m_x[i].re = a[i*2];
m_x[i].im = a[i*2+1];
}
//m_om - order massive init
k = (int)(log(BASE)/log(2));
l = (int)(log(N)/log(2)) - k;
for (i = 0; i < N; i++)
m_om[i] = ((i % BASE) << l) + i/BASE;
//m_w init
for (i = 0; i < N/2; i++)
{
m_w[i].re = cos((double)(2*i*PI)/(double)N);
m_w[i].im = sin((double)(2*i*PI)/(double)N);
}
}
void TForm1::baseoperation(complex x0, complex x1, complex x2, complex x3, complex *y0, complex *y1, complex *y2, complex *y3)
{
(*y0).re = x0.re + x1.re + x2.re + x3.re;
(*y0).im = x0.im + x1.im + x2.im + x3.im;
(*y1).re = x0.re + x1.im - x2.re - x3.im;
(*y1).im = x0.im - x1.re - x2.im + x3.re;
(*y2).re = x0.re - x1.re + x2.re - x3.re;
(*y2).im = x0.im - x1.im + x2.im - x3.im;
(*y3).re = x0.re - x1.im - x2.re + x3.im;
(*y3).im = x0.im + x1.re - x2.im - x3.re;
}
void TForm1::fft4()
{
int i,imax,j,x1,x2;
complex temp1,temp2;
for (i = 0; i < N; i = i + 4)
{
baseoperation(m_x[m_om[i]],m_x[m_om[i+1]],m_x[m_om[i+2]],m_x[m_om[i+3]],
&m_fft4[i],&m_fft4[i+1],&m_fft4[i+2],&m_fft4[i+3]);
}
for(imax = BASE-1; imax < N/2; imax = (imax+1)*2 - 1)
{
for(j = 0; j < N; j = j + (imax+1)*2)
{
x1 = j; x2 = j + (imax+1);
for (i = 0; i <= imax; i++)
{
temp1.re = m_fft4[x1+i].re + m_fft4[x2+i].re*m_w[i].re - m_fft4[x2+i].im*m_w[i].im;
temp1.im = m_fft4[x1+i].im + m_fft4[x2+i].re*m_w[i].im + m_w[i].re*m_fft4[x2+i].im;
temp2.re = m_fft4[x1+i].re - m_fft4[x2+i].re*m_w[i].re + m_fft4[x2+i].im*m_w[i].im;
temp2.im = m_fft4[x1+i].im - m_fft4[x2+i].re*m_w[i].im - m_w[i].re*m_fft4[x2+i].im;
m_fft4[x1+i].re = temp1.re;
m_fft4[x1+i].im = temp1.im;
m_fft4[x2+i].re = temp2.re;
m_fft4[x2+i].im = temp2.im;
}
}
}
// Виводжу дані у таблицю 2
ValueListEditor2->ColCount=2;
String str1,str2;
for(int tt=0;tt<512;tt=tt+2)
{ str1=FloatToStrF(m_fft4[tt].re,ffGeneral,7,10);
str2=FloatToStrF(m_fft4[tt+1].im,ffGeneral,7,10);
ValueListEditor2->InsertRow(str1,str2,true);
}
// Відкриваю вихідний файл
ofstream file1;
file1.open("fft.xy");
// Записую дані в файл
for(int kk=0;kk<512;kk++)
file1<<m_fft4[kk].re<<" "<<m_fft4[kk].im<<"\n";;
file1.close();
}
Main.h
//-------------------------------------------------------------------
#ifndef mainH
#define mainH
//---------------------------------------------------------------------------
#include <Classes.hpp>
#include <Controls.hpp>
#include <StdCtrls.hpp>
#include <Forms.hpp>
#include <Grids.hpp>
#include <ValEdit.hpp>
#include <Buttons.hpp>
//---------------------------------------------------------------------------
#define N 256
#define BASE 4
#define PI 3.14159265358979323846
class TForm1 : public TForm
{
__published: // IDE-managed Components
TValueListEditor *ValueListEditor1;
TLabel *Label1;
TValueListEditor *ValueListEditor2;
TLabel *Label2;
private: // User declarations
public: // User declarations
struct complex
{ double re;
double im;
};
__fastcall TForm1(TComponent* Owner);
void FFT();
void fft4();
void baseoperation(complex x0,complex x1, complex x2, complex x3, complex *y0, complex *y1, complex *y2, complex *y3);
complex m_x[N];
complex m_w[N/2];
complex m_fft4[N];
int m_om[N];
void initall();
};
//---------------------------------------------------------------------------
extern PACKAGE TForm1 *Form1;
//---------------------------------------------------------------------------
#endif
Приклад виконання:
5.Висновок: на лабораторній роботі досліджувала швидкі алгоритми дискретних тригонометричних перетворень.