Fourierova transformace
Transformace
Zpracovánív
transform. oblasti
Zpracovánív
časové oblasti
x(n) X(n)
Základní idea transformace
InverzníTransformace
transform. oblastičasové oblasti
x(n)‘ X(n)‘
Typy Fourierových transformací
Diskrétní Fourierova transformaceN-bodová DFT signálu s N vzorky:
1,,1,0)(][1
0
/2 −==∑−
=
− NkenxkXN
n
NnkjDFT L
π
[ ]∑−
=
−=1
)/2sin()/2cos()(][N
DFT NnkjNnknxkX ππ
zpětná transformace IDFT :
1,,1,0][1
][1
0
/2 −== ∑−
=
NnekXN
nxN
k
NnkjDFT L
π
Jelikož ej2πnk/N je periodická, je periodická i DFT a IDFT ⇒ počítáme vzorky pouze přes jednu periodu.
∑=0n
Ximag
Ф
Xmag
Xm(k)=Xreal(k)+jXimag(k)
Polární tvar DFT
Xreal
22 )()( kXkXXX imagrealmagmag +==
= −Φ )(
)(tan 1
kX
kXX
real
imag
222 )()()()( kXkXkXkX imagrealmagPS +==
• Při použití polární reprezentace DFT – pozor na následující možné problémy :
– správnou konverzi fáze - sw většinou vrací fázový úhel v radiánech a to v rozsahu
– při výpočtu fáze pozor na nulovou reálnou část ( přetečení) (fáze je v tomto případě ±90º
– pozor na správnou konverzi úhlu z intervalu na interval – fáze u velmi nízkých amplitud, které se ztrácí v šumů může chaoticky kmitat
okolo nulové hodnoty
– fázová charakteristika se opakuje s periodou 2π – to způsobuje v některých případech nespojitost ve fázové charakteristice
– Amplitudová frekvenční charakteristika je vždy kladná – problémy – Amplitudová frekvenční charakteristika je vždy kladná – problémy mohou nastat pokud imaginární část transformovaného signálu je celá nulová -> může docházet ke prudké změně ve fázové charakteristice mezi –π a π
Př.: Výpočet DFT z definice: x[n]={1,2,1,0}
jeeenxXk
enxXk
jj
n
nj
DFT
nDFT
20121][]1[:1
40121][]0[:0
3
23
0
2
3
0
0
=+++===
−=+++===
=+++===
−−−
−−
=
−
=
∑
∑
∑
πππ
πππ
jeeenxXk
eeenxXk
jj
n
nj
DFT
jj
n
jnDFT
20121][]3[:3
00121][]2[:2
32
33
0
2
3
23
0
=+++===
=+++===
−−
=
−
−−
=
−
∑
∑
πππ
πππ
XDFT[k]={4,-2j,0,2j}
Symetrie DFT:
DFT reálného signálu vykazuje komplexně sdruženou symetrii okolo počátku tj.:
XDFT[-k]=XDFT* [k]
Protože je DFT periodické XDFT[-k]=XDFT [N-k]
Výpočet hodnot pro k=0 a k=N/2 (pro sudé N)
∑ ∑
∑∑−
=
−
=
−
=
−
=
−==
−==
1
0
1
0
1
0
1
0
][)1(1
]2
[][1
]0[
][)1(]2
[][]0[
N
k
N
kDFT
kDFT
N
n
nDFT
N
nDFT
kXN
NxkX
Nx
nxN
XnxX
Vlastnosti DFT
1. Linearita k1x1(n)+ k2x2(n) ↔ k1X1(n)+ k2X2(n)2. Periodičnost - funkce x(n) a X(n) jsou periodické s periodou P=N3. Kruhový časový posun :
celénkXennxkn
Nj
−⋅↔−−
0
)2
(
0 ],[][0
π
⇒ posun v čase způsobí změnu ve fázi4. Kruhový frekvenční posun
celénkXennx −⋅↔− 00 ],[][
celéknxekkXkn
Nj
−⋅↔−−
0
)2
(
0 ],[][0
π
4. Periodická konvoluce v časové oblasti
Periodickou konvoluci dvou sekvencí délky N určíme jako součin N-bodových transformací.
5. Periodická korelace
][][][][ 2121 kXkXnxnx ⋅↔∗
][][][][ *2121 kXkXnxnx ⋅↔∗∗
6. Obraz obrácené posloupnosti
7. Vlastnosti spektra reálné posloupnosti
][][ kXnx −↔−
][][][ kNXkXkX −=−= ∗∗
)()(
)()(
)](Im[)](Im[
)](Re[)](Re[
kNk
kNXkX
kNXkX
kNXkX
−−=−=
−−=−=
φφ
5. Periodická konvoluce ve frekvenční oblasti
8. Vlastnosti spektra reálné a sudé posloupnosti• je-li x[n] reálná a sudá je i X[k] reálná sudá
9. Vlastnosti spektra reálné a liché posloupnosti• je-li x[n] reálná a lichá, pak je X[k] imaginární, li chá
10. Alternativní vzorec pro výpočet IDFT
][][1
][][ 2121 kXkXNnxnx ∗↔⋅
*1
0
2* )(
1)(
= ∑
−
=
−N
k
knN
jekX
Nnx
π
K výpočtu inverzní transformace je možné použít algoritmů pro výpočet DFT:
• nejprve obrátíme znaménka hodnot imaginární části X(k),• vypočteme DFT• obrátíme znaménka imaginárních částí vypočtených hodnot• výsledek vydělíme N
Vlastnosti fázové charakteristiky
Velikost DFT vzorků
• Je-li x(n) reálný vstupní signál složený ze sinusovek s odpovídající amplitudou A0 a celočíselným počtem cyklů přes N vzorků, je velikost DFT vzorku odpovídající sinusovky Mr daná vztahem
0
NAM r = 20AM r =
• Pro komplexní vstupní signál s velikostí A0 (tj. A0ej2πft) je výstupnívelikost DFT vzorků
NAM r 0=
Rozlišení DFT
Amplitudové, fázové, výkonové spektrum
Amplitudové spektrum
Replikace signálu a nulová interpolace spektra
Existuje vztah mezi replikací v jedné oblasti a nulovou interpolací v druhé oblasti:
Je-li x[n] M-krát replikován, je DFT(x[n]) interpolován (nulová interpolace) a násoben M:
x[n] → DFT → XDFTk]
{ x[n], x[n], … ,x[n] } ↔ M XDFT↑[k/M]
M-krát
M-násobná interpolace:
x↑[n/M] ↔ {XDFTk], XDFTk], … , XDFTk] }
M-krát
DFT dvojice
Impulz : {1,0,0, …, 0} ↔ {1, 1, 1, …, 1} konstanta
Konstanta: {1, 1, 1, …, 1} ↔ {N, 0,0, … , 0, 0} impulz
Exponenciála:
N
kj
nn
eπ
α
αα 21
1−
−
−↔
Sinusoida:
Neα1−
)]([5.0][5.02cos 000 kNkNkkN
N
kn −−+−↔
δδπ⇓
N-bodová sinusoida s periodou N a F0=k/N má pouze dva nenulové vzorky
Inverzní DFT
1,,1,0][1
][1
0
/2 −== ∑−
=
NnekXN
nxN
k
NnkjDFT L
π
DFT a IDFT mají podobný vztah, až na znaménko u exponenciály a normalizační faktor 1/N. Je-li DFT symetrická (komplexně sdružená symetrie) ⇒ x[n] může být vyjádřena jako součet sinusoid, protože dvojici AejΘ v k0 a v N-k0 odpovídá sinusoida
Θ+N
kn
N
A 02cos2 π
XDFT[0] odpovídá dc složce XDFT[0]/N a DFT vzorek XDFT[N/2] (pro sudé N) odpovídá frekvenci F= 0.5 a odpovídá signálu Kcos(nπ), kde K=(1/N)XDFT[N/2]
Př.: XDFT[k] = {4, -2j, 0, 2j}
2)2024(25.0][25.0]1[:1
1)0121(25.0][25.0]0[:0
3
2
3
23
0
2
3
0
0
=++−⋅=⋅==
=++−⋅=⋅==
=+++⋅=⋅==
∑
∑
∑
=
=
πππ
πππjj
k
kj
DFT
kDFT
jejeekXxn
ekXxn
0)2024(25.0][25.0]3[:3
1)2024(25.0][25.0]2[:2
2
3
2
33
0
2
3
33
0
=++−⋅=⋅==
=++−⋅=⋅==
∑
∑
=
=πππ
πππ
jj
k
kj
DFT
jj
k
jkDFT
jejeekXxn
jejeekXxn
IDFT → x[n]={1,2,1,0}
*1
0
2* )(
1)(
= ∑
−
=
−N
k
knN
jekX
Nnx
π
K výpočtu inverzní transformace je možné použít algoritmů pro výpočet DFT:
Alternativní vzorec pro výpočet IDFT
DFT: • nejprve obrátíme znaménka hodnot imaginární části X(k),• vypočteme DFT• obrátíme znaménka imaginárních částí vypočtených hodnot• výsledek vydělíme N
DFT- volba frekvenční osy
Prosakování ve spektru (Spectrum leakage)
• Prosakování se vyskytuje tehdy, pokud u vzorkovaného analogového signálu počítáme DFT z N vzorků a v těchto vzorcích není obsažen celočíselný počet period sinusoid obsažených ve vstupním signálu.
• Při prosakování se v DFT vyskytnou spektrální čáry i jinde než v ±f0 (f0 je frekvence vstupního signálu
• Prosakování je možné eliminovat, ale nelze jej odstranit úplně.
• K vysvětlení prosakování ve spektru se používá Fourierův obraz reálné kosinové funkce
]/)(sin[)](sin[
21
)(
),/2cos()(
]/)()([
Nmk
mkemX
NnKnx
Nmkmknjr
r
+−
−⋅=
=
−−−
ππ
ππ
]/)(sin[)](sin[
21
]/)(sin[2
]/)()([
Nmk
mke
Nmk
Nmkmknj
++⋅+
−
−−+
ππ
ππ
NnmjN
nw enxnwmX
/21
0
)()()( π==
=
⋅=∑
K omezení vlivu prosakování ve spektru se před provedením DFT násobí vstupní signál tzv. vyhlazovacím oknem.
Maticový zápis DFT a IDFT
N
j
N eWπ2−
= [ ] [ ]
[ ] [ ][ ] 1,...,1,011,...,1,0
1
0
*
1
0
−==
−==
∑
∑−
=
−
=
NnWkXN
nx
NkWnxkX
N
k
nkNDFT
N
n
nkNDFT
Označíme
Nx1 matice
[ ][ ][ ]
[ ]
[ ][ ][ ]
[ ]
−
=
− −−−−
−
−
1
2
1
0
1
2
1
0
)1)(1()1(210
0
)1(2420
1210
00000
Nx
x
x
x
WWWW
W
WWWW
WWWW
WWWWW
NX
X
X
X
NNN
NN
NNN
N
NNNNN
NNNNN
NNNNN
M
L
MOMM
L
L
M
N DFT rovnic lze zapsat maticově X=WNX
Nx1 matice
NxN matice
Maticový zápis IDFT
[ ]
[ ]TNN
T
NN
WN
W
XWN
XXWX
∗−
∗−
=
=⇒=
1
1
1
1
IDFT maticeN
• výběr N vzorků vzorkovaného signálu = násobení vzorkovaného signálu pravoúhlým oknem
• kromě pravoúhlého okna lze použít i jiné typy oken• každé okno má určitou frekvenční charakteristiku a podle toho je vhodné
pro zpracování určitého typu signálu
Spektrální vyhlazení časovými okny
Základní charakteristiky okenZákladní charakteristiky oken
• Amplitudová charakteristika
3
6
3707.0
log20
65.0
log20
WdBP
P
WdBP
P
⇒=
⇒=
Další charakteristiky oken:
koherentní zesílení (coherent gain)
ENBW (equivalent noise bandwith)
[ ]∑−
=
=1
0
1 N
k
kwN
CG
[ ]∑− 1
2N
kwN
SL (scallop loss)
[ ]
[ ]∑
∑−
=
==1
0
0
2
N
k
k
kw
kwNENBW
[ ]
[ ]∑
∑−
=
−
=
−
= 1
0
1
0log20 N
k
N
k
N
kj
kw
ekw
SL
π
Používaná DFT okna a jejich charakteristiky
Používaná DFT okna a jejich charakteristiky
Bartlett
von Hann
Hamming
Blackmann
Kaiser
• Zvětšování délky okna (u všech oken) klesá šířka hlavního laloku, velikost postranních laloků se téměř nemění
• Ideální případ amplitudové charakteristiky (pro okno dané délky) - úzký (a co nejvyšší) hlavní lalok a nízké postranní lalok y -> protichůdnépožadavky - většinou s klesající šířkou hlavního laloku roste amplituda postranních laloků
• Dynamické rozlišení - schopnost odlišit velké změny v amplitudě signálu.
• Frekven ční rozlišení okna - schopnost od sebe odlišit sinusovky s podobnou amplitudou a blízkou frekvencí. Platí: při použití vyhlazovacích oken od sebe nelze odlišit dvě sinusovky s Platí: při použití vyhlazovacích oken od sebe nelze odlišit dvě sinusovky s frekvencemi nižšími než je šířka hlavního laloku okna.
N
kW
S
fF M ==
∆=∆
ke zmenšení ∆F musíme zmenšit WM → zvětšit N (ale ne přidáním nul !!!) → více vzorků signálu
Příklad: x(t)=A1cos(2πf0t)+ A2cos(2π(f0+∆f)t)A1 = A2 = 1 f0=30Hz S=128HzMáme N vzorků, doplnit počet nulami na NFFTJaké je nejmenší ∆f pro: pravoúhlé okno von Hannovo okno
pro následující N a NFFT: a) N=256, NFFT=2048
b) N=512, NFFT=2048
c) N=256, NFFT=4096
⇒ nelze zvětšit rozlišení p řidáním nul, ale pouze zv ětšením délky signálu
A2 = 0.05 (26dB pod A1 tj. 20log(1/0.05) )
• Výkonová spektrální hustota (PSD – power spectral density) Rxx(f) analogového výkonového nebo náhodného signálu x(t) je Fourierova transformace autokorelační funkce rxx(t).
• Rxx(f) je reálná nezáporná sudá funkce s hodnotou Rxx(0) rovnou průměrnému výkonu signálu x(t).
• PSD se používá k odhadu spektra vzorkovaného signálu konečné délky
Dva způsoby odhadu spektra :
Odhad spektra
– neparametrický odhad – o signálu nic nevíme a odhadujeme spektrum– parametrický odhad – odhad spektra modelujeme jako koeficienty filtru
buzeného šumem
Neparametrické odhady: • Periodogram – periodogram P[k] je založena na DFT (FFT) N-vzorkované řady x[n] a je definován jako:
[ ] [ ] 21 kXN
kP DFT=
• dobrý odhad pro deterministické, pásmově omezené, výkonové signály, vzorkované vyšší frekvencí než Nyquistova
• Špatné pro signály poškozené šumem
Welchova metoda• je založená na průměrování překrývajících se periodogramových odhadů• Princip:
– signál je rozdělen na k překrývajících se M vzorkovaných úseků (s překrytím D vzorků, D=50-75%). Každý úsek je násoben M vzorkovým oknem w(n) kvůli omezení prosakování.
– PSD každého úseku je odhadnuto periodogramem – PSD každého úseku je odhadnuto periodogramem – k periodogramů je zprůměrováno – M-bodový průměrný periodogram
• Metoda je vhodná pro detekci frekvencí, které jsou blízko sebe (dobré frekvenční rozlišení)
Bartlettova metoda • Jako Welchova metoda, ale segmenty se nepřekrývají a nejsou násobeny
oknem
Blackman-Tukey metoda: • metoda určuje PSD z autokorelace násobené oknem• Princip:
– N vzorkový signál x[n] je doplněn nulami na 2N vzorků výsledkem je signál y[n]– určí se periodická autokorelace y[n] nalezením YFFT[k] a určením IFFT součinu
YFFT[k] YFFT[k]* a dostaneme 2N bodový autokorelační odhad rxx[n].– autokorelační odhad rxx[n] je násoben M-bodovým oknem (kvůli vyhlazení
spektra a redukci vlivu špatného autokorelačního odhadu signálu konečné délky)– určíme FFT M-bodové autokorelace násobené oknem a dostaneme vyhlazený
periodogram– menší M (užší okno) lepší vyhlazení, ale může dojít k maskování některých – menší M (užší okno) lepší vyhlazení, ale může dojít k maskování některých
špiček (peaků) nebo k překrytí ostrých detailů. Obvykle se používá M v rozsahu M=0.1N až M=0.5N
– jako okno pro vyhlazení se nejčastěji používá Bartlettovo okno.
• metoda je vhodná pro detekci dobře oddělených peaků s rozdílnou velikostí (dobré dynamické rozlišení)
Rychlá Fourierova transformaceFFTFFT
2)2
(2
22
22)2
(2
)(22
113
2
1
Njj
NkN
NNkj
kj
nN
Nn
NN
njN
Nnj
nN
NnN
NNnj
Nnj
Wee
WWee
WWee
===
−==
==
−−
−−
+−+−
++−−
ππ
ππ
ππ
ππ
Algoritmus využívá symetrie a periodicity exponenciály WN=e-j2πn/N
2
222)
2(2
4 NNN
jN
jWWee ==
−− ππ
Základní výsledky využité ve FFT výpočtech:
1-bodová transformace : DFT jediného čísla A je číslo A2-bodová transformace : DFT 2-bodové řady
XDFT[0]=x[0] + x[1] XDFT[1]=x[0] - x[1]
• Radix 2-FFT transformace využívá toho, že N-bodová DFT může být zapsána jako součet dvou N/2 bodových transformací vytvořených ze sudých a lichých vzorků.
[ ] [ ] [ ] [ ]
[ ] [ ]
[ ] [ ]∑∑
∑∑
∑∑∑
−−
−
=
−
=
−
=
+−
=
−
=
=++=
=++==
12
12
12
0
2
12
0
2
12
0
)12(
12
0
21
0
122
122
NN
N
n
nkN
kN
N
n
nkN
N
n
knN
N
n
nkN
N
n
nkNDFT
WnxWWnx
WnxWnxWnxkX
Označíme Xe[k] – sudé vzorky (sudý index) a Xo[k] - liché vzorky (lichý index), pak
[ ] [ ]∑∑==
++=2
0 2
2
0 2
122n
nkN
kN
n
nkN WnxWWnx
[ ] [ ] [ ] 1,,2,1,0 −=+= NkkXWkXkX okNeDFT L
[ ] [ ] [ ]
[ ] [ ] 12
1,,2,1,0
222
12
1,,2,1,0
2
−=−=
=
++
+=
+
−=+=
+
NkkXWkX
NkXW
NkX
NkX
NkkXWkXkX
okN
e
oN
k
Ne
DFT
okN
eDFT
L
L
Na základě periodicity lze odvodit:
tzv. Danielson-Lanzos lemma
A=Xe[k]
B=Xo[k]
Výpočet FFT – decimace ve frekvenci DIF
Výpočet FFT – decimace v čase DIT
Výpočet FFT - porovnání