UNIVERZITA JANA EVANGELISTY PURKYNĚ V ÚSTÍ NAD LABEM FAKULTA ŽIVOTNÍHO PROSTŘEDÍ
___________________________________________________________________________
"
Furierova transformace"
Seminární práce z předmětu Dálkový průzkum Země
Marcela Bartošová, Veronika Bláhová OŽP, 3.ročník Datum: 20.2. 2006
Most 2006
1
Furierova transformace Furierova transformace je matematická metoda, která dovoluje analyzovat průběh libovolného signálu a převést jej na součet sinusových signálů vhodných frekvencí a amplitud. V obrazovém ‚signálu‘ pak nejvyšší nalezené frekvence odpovídají čárové frekvenci, která musí být zaznamenána. Fourierova transformace je modifikací Fourierovy řady a je užitečná pro řešení mnoha různých problémů. Používá se např. pro převedení řešení diferenciálních rovnic na řešení algebraických rovnic nebo pro frekvenční analýzu časově proměnných signálů. V oblasti zpracování obrazů je možné Fourierovu transformaci uplatnit pro úpravy kvality obrazů, ale také pro vyhodnocování prostorových frekvencí, což lze s výhodou použít pro vyhodnocování interferenčních řádů v obrazech interferogramů. Dvojrozměrná Fourierova transformace umožní převést rozložení obrazových intenzit I(x, y) vyhodnocovaného obrazu na obraz prostorových frekvencí F(fx, fy). Definicí dvojrozměrné Fourierovy transformace je vztah
Podobně jako v oblasti signálů spojitých, je možné i v oblasti signálů diskrétních definovat transformaci, která bude diskrétní obdobou Fourierovy transformace ve spojité oblasti. Tuto transformaci nazýváme Discrete Fourier Transform - DFT - Diskrétní Fourierova transformace. Ale vzhledem k tomu, že výpočet DFT vyžaduje značný počet násobení ( ), což je časově nejnáročnější operace, byl vyvinut algoritmus, umožňující značné urychlení výpočtu. Tento algoritmus je označován Fast Fourier Transform - FFT - Rychlá Fourierova Transformace.
2
Někteří fyzikové považují Fourierovu transformaci za fyzikální fenomén, nejen za nástroj pro matematické výpočty. Využívá ji například kvantová mechanika. Konvolucí vektorů a = (a0, a1, a2, . . . , an-1), b = (b0, b1, b2, . . . , bn-1) rozumíme vektor a×b = c = (c0, c1, c2, . . . , c2n-2), pro který platí: cj = Σk=0..j akbj-k Příklad konvoluce: Součin polynomů stupně n, reprezentovaných vektorem koeficientů. Definice poskytuje algoritmus s asymptotickou složitostí θ(n2). Strategie výpočtu v čase θ(n. log n) převodem na hodnotovou reprezentaci: 0. Zdvojení délky vektorů, protože součin (konvoluce) má dvojnásobnou délku 1. Převod na reprezentaci hodnotami - přímá FFT v čase θ(n. log n) 2. Vynásobení obrazů člen po členu - θ(n) 3. Převod zpět na koeficienty - inverzní FFT v čase θ(n. log n)
Definiční vztahy pro Fourierovu transformaci Fourierův integrál:
Jedná se o komplexní integrál s parametrem, který definuje transformaci (obecně komplexní) funkce f(t) na její Fourierův obraz F(ω). Zpětná transformace je dána vztahem:
Diskrétní Fourierova transformace a její inverze:
, kde k = 0, 1, 2, ... n - 1
, kde j = 0, 1, 2, ... n - 1 Fourierova transformace je vzájemně jednoznačné lineární zobrazení! Diskrétní verze je jednoznačně popsaná čtvercovou maticí. Algoritmus pro inverzní transformaci je jen drobnou modifikací algoritmu pro přímou transformaci. Pro komplexní exponencielu se často se používá značení:
3
Veličina ωn se nazývá "twiddle factor" (otáčecí činitel), v jiné literatuře "complex n-th root of unity" (n-tá komplexní odmocnina jedničky). DFT tedy můžeme psát ve tvaru:
Při počítačovém zpracování digitalizovaných obrazů se ale používá tzv. diskrétní dvojrozměrná Fourierova transformace. Výsledkem této transformace je pak obraz četností všech prostorových frekvencí ve vyhodnocovaném obraze, a to v souřadném systému fx, fy. Fourierova transformace je vhodná především pro vyhodnocování tvarově složitých interferogramů získaných při seřízení interferometru na konečnou šířku interferenčních proužků, kde můžeme obvykle lépe rozlišit, zda interferenční řád roste či klesá, viz obr. 16-9.
Obr. Interferogram neizotermního vzduchového proudu z vyústky obtékajícího překážku
Z obrazu četností frekvencí lze vyvozovat různé závěry o prostorových frekvencích ve vyhodnocovaném obraze. Při vyhodnocování interferogramů nás obvykle zajímají složky nejnižší prostorové frekvence. Některé postupy zpracování vizualizačních experimentů však provádějí i různé úpravy obrazu četností frekvencí a pak aplikují zpětnou Fourierovu transformaci, která umožní získat opět rozložení intenzit I(x, y) ve zkoumaném obraze, ale s upravenými prostorovými frekvencemi. Např. po odstranění vysokých prostorových frekvencí v obraze frekvencí umožní zpětná Fourierova transformace získat původní obraz intenzit bez zrnitosti. Ponecháme-li v obraze frekvencí jen určité oblasti frekvencí, můžeme po aplikaci zpětné Fourierovy transformace získat obraz se zvýrazněnými strukturami, např. v leteckém snímku města lze zvýraznit ulice vedoucí stejným směrem apod.
Využití Fourierovy transformace: a) Analýza signálu
4
Analyzovaný diskrétní signál je rozdělen na krátké časové úseky. Na každém z nich se spočte Fourierův obraz. Jeho absolutní hodnota v bodě f je amplituda odpovídající frekvenci f. Na vodorovnou osu promítneme čas, na svislou osu frekvenci a amplitudě dáme význam barvy. Získáme spektrogram, který přehledně zobrazuje časový vývoj signálu.
Na obrázku je spektrogam zvuku nízko letícího letadla. Je na něm dobře patrný Dopplerův efekt - změna frekvence zvuku v závislosti na pohybu zdroje. b) Spektra Spektra generovaných obrázků jsou počítána funkcí FFT v Matlabu. Spektra jsou počítána ze světelnosti (kanál Y). Kvůli vyniknutí detailů jsou ve spektrech dělány úpravy ( snížení jasu ). Jsou zde vykreslena "jednostraná" amplitudová a fázová spektra. Z estetického hlediska je zde také přidán přímo výsledek FFT z Matlabu. FFT počítá i záporné frekvenční pásmo, ale to je z důvodu symetrie spekter je stejné jako kladné frekvenční pásmo. Z tohoto důvodu není potřeba uchovávat tyto informace. c) Frekvenční obraz obrázků
Obrázek
Obrázek
5
Fázové spektrum
Fázové spektrum
Amplitudové spektrum
Amplitudové spektrum
Celé amplitudové spektrum
Celé amplitudové spektrum
6
d) Komprese obrázků Pro kompresi těchto obrázků je kvůli objektivnímu posouzení z hlediska transformací využit stejný algoritmus jako u komprese s využitím DCT. Pro všechny obrázky je opět stejné nastavení. U této transformace by se dal algoritmus trochu modifikovat k dosažení lepších kompresních poměrů díky symetrii frekvenčního spektra OBR1 Obr1 před kompresí
Obr1 po kompresi
Obrázek sestavený ze zahozených koeficientů
Statistická data - Matlab Velikost obrázku před kompresí: 1843200 B Velikost obrázku po kompresi: 46640 B Kompresní poměr: 1 : 39.5197 Ušetřené místo v procentech 97.46962%
OBR2 Obr2 před kompresí
Obr2 po kompresi
7
Obrázek sestavený ze zahozených koeficientů
Statistická data - Matlab Velikost obrázku před kompresí: 1843200 B Velikost obrázku po kompresi: 258960 B Kompresní poměr: 1 : 7.117702 Ušetřené místo v procentech 85.95052%
OBR3 Obr3 před kompresí
Obr3 po kompresi
Obrázek sestavený ze zahozených koeficientů
Statistická data - Matlab
8
Velikost obrázku před kompresí: 1843200 B Velikost obrázku po kompresi: 26928 B Kompresní poměr: 1 : 68.4492 Ušetřené místo v procentech 98.53906%
e) Frekvenční analýza rotačních strojů
Co to je? Pravidelné sledování strojů, jejich vibrací, vyhodnocování vibračních spekter a určování závad. K čemu je to dobré? Na základě vyhodnocení vibračních spekter rotačních strojů lze získat přesný přehled o jeho technickém stavu, predikovat pravděpodobné závady a určit zbývající životnost zařízení. Na základě těchto podkladů je možné naplánovat opravu stroje a předejít tak krizovým stavům a neočekávaným haváriím. Opačným případem je zbytečně uspěchaná oprava, kdy se z preventivních důvodů výměňují i ty díly, které by bylo možné ještě bezpečně provozovat a častějším opravám, než by bylo potřeba. Oba extrémy, jak krizový stav, tak uspěchaná oprava, zvyšují náklady na údržbu stroje. Pomocí vibračních spekter lze tedy určit nejvhodnější okamžik opravy, což se projeví předevšim na nákladech, ale i zvýšené životnosti zařízení. Není to moc drahé? Není. Samozřejmě záleží na množství sledovaných strojů. Pokud někdo má na starosti dvě čerpadla, asi se mu to příliš nevyplatí. Pokud ale má strojů více a náklady na jejich údržbu tvoří nezanedbatelnou položku, potom mu vibrační diagnostika může tyto náklady výrazně snížit, a to tak, že dokáže určit vhodný okamžik opravy. Ten může nastat při dosažení mezních
9
vibrací určitého elementu, anebo při výskytu většího počtu závad najednou. Nedochází tak k tomu, že se opraví jedna součást, po čase další, a tak náklady rostou a rostou. Jak to funguje? Technik přijde k danému stroji se speciálním zařízením — vibrometrem a pomocí snímačů naměří v předem definovaných bodech průběhy vibrací. Ty jsou obvykle ještě tím samým zařízením převedeny použitím tzv. rychlé furierovy transformace (FFT –– fast furier transformation) na frekvenční spektrum, ze kterého se dají vyčíst případné závady. Přiklad: stroj má poškozené kuličkové ložisko –– na vnějším kroužku je vada. Pokud má lož isko např. 12 kuliček, pak každá kulička, která projde přes toto poškození, způsobí chvění, tedy 12-krát za jednu otáčku. Jsou-li provozní otáčky stroje 1 500 min-1, potom se tato vada projeví ve vibračním spektru jako dvanáctinásobek základních otáček, tj. 12×1500=18 000 min-1 atp.
Používané programy: • •
fft.pas spectrum.pas
fft.pas a spectrum.pas jsou programy, na kterých jsem zmíněné unity testoval. První pracuje v textovém módu a převádí vstupní posloupnost komplexních čísel na její Fourierův obraz a zpět. Celý výpočet probíhá v jednom poli, kde se vstupní hodnoty přepíšou výstupními. Délka posloupnosti musí být mocnina dvojky. Druhý program (spectrum.pas) pracuje také s jediným vstupním vektorem. Spočítá jeho spektrum a zobrazí graficky jeho absolutní hodnotu. Ve výsledném grafu je tedy na x-ové ose frekvence a na y-ové amplituda.
Prameny: http://hyperkrychle.cz/fftpas.html http://hyperkrychle.cz/praktika/01.html http://www.inforum.cz/inforum2001/prispevky/psohlavec.htm http://www-dt.fme.vutbr.cz/users/pavelek/optika/1607.htm http://hyperkrychle.cz/praktika/05.html http://www.euromont.cz/cz/sp_cz.html
10
