CYFROWA ANALIZA WIDMA SYGNAŁÓW OKRESOWYCH
CYFROWA ANALIZA WIDMA SYGNAŁÓW OKRESOWYCH
Cel
Celem ćwiczenia jest poznanie podstawowych własności dyskretnego przekształcenia Fouriera oraz zwrócenie uwagi na podstawowe źródła błędów i ich konsekwencje w wyniku estymacji widma amplitudowego.
Wprowadzenie
Analizę widmową, zwaną również analizą częstotliwościową, stosuje się do rozwiązywania zagadnień z dziedziny fizyki i techniki. Umożliwia ona określenie częstotliwości składowych zawartych w przebiegu czasowym funkcji, nazywanym dalej sygnałem.
Jednym z najważniejszych narzędzi do analizy widmowej jest przekształcenie Fouriera zwane również transformatą Fouriera, które umożliwia transformowanie sygnału w dziedzinie czasu na równoważną postać w dziedzinie częstotliwości (widmo) i na odwrót. Przebieg czasowy i jego postać częstotliwościowa tworzą parę transformat.
Inaczej mówiąc widmem sygnału nazywamy jego obraz po przekształceniu z dziedziny czasowej do dziedziny częstotliwościowej. Obie reprezentacje są równoważne i zawierają te same informacje o sygnale, przy czym każda z nich przedstawia je w różny sposób. Umożliwia to np. zauważenie w dziedzinie częstotliwościowej (widmo sygnału) pewnych cech sygnału nie dających się zaobserwować w dziedzinie czasowej (przebieg sygnału). Dla sygnałów okresowych (periodycznych) u podstaw analizy widmowej leży twierdzenie mówiące, że prawie każdy taki sygnał można rozłożyć na szereg harmonicznych (szereg Fouriera) czyli nieskończoną sumę składowych będących sygnałami sinusoidalnymi o różnych fazach, amplitudach i częstotliwościach. Przy czym kolejne harmoniczne muszą mieć częstotliwości będące całkowitymi wielokrotnościami częstotliwości podstawowej. Ponadto z przekształcenia Fouriera wynika, że widmo elementu bazowego (podstawowego elementu na jakie rozkładamy sygnał), czyli sygnału sinusoidalnego jest pojedynczym prążkiem leżącym na częstotliwości zgodnej z częstotliwością transformowanego sygnału i mającym wysokość równą amplitudzie tegoż sygnału (Rys. 1.).
Przykładowo widmo sygnału prostokątnego o częstotliwości fg= 4 Hz i amplitudzie A będzie składało się z harmonicznych, z których pierwsza będzie miała częstotliwość fg1= fg=4 Hz, druga fg2=8 Hz, trzecia fg3=12 Hz a n-ta fgn=n*fg=n*4 Hz. Ale ponieważ dla tego sygnału amplitudy parzystych harmonicznych wynoszą zero więc w widmie obserwować będziemy tylko nieparzyste (Rys. 5.) o amplitudach ok. A/n.
Podstawowe informacje o przekształceniu Fouriera i jego własności
Pominiemy matematyczne podstawy przekształcenia Fouriera, omówione zostaną tylko własności oraz przypadek dyskretnej transformaty Fouriera (DTF lub DFT).
Z przekształcenia Fouriera wynika, że każdy ciągły sygnał okresowy można przedstawić jako sumę przebiegów sinusoidalnych o odpowiednio dobranych amplitudach i fazach oraz o częstotliwościach będących całkowitymi wielokrotnościami częstotliwości podstawowej (fg).
Rys. 1. Przebieg czasowy (a) sygnału sinusoidalnego przy zastosowaniu okna prostokątnego i jego widmo (b)
Rys. 2. Przebieg czasowy (a) sygnału sinusoidalnego przy zastosowaniu okna prostokątnego i jego widmo (b). Efekt „przecieku” widma.
Rys. 3. Przykładowa funkcja okna innego niż prostokątne (okno Lanczosa). Kształt okna w dziedzinie czasowej (a) oraz jego widmo (b)
Rys. 4. Zastosowanie funkcji okna. Iloczyn sygnału i funkcji okna (a) oraz widmo tego iloczynu (b)
Rys. 5. Przebieg czasowy (a) sygnału prostokątnego o częstotliwości 4 Hz i jego widmo (b)
Kompletny obraz widmowy (częstotliwościowy) takiego sygnału tworzą dwa widma: amplitudowe i fazowe. W postaci wykresu przedstawia się zwykle widmo amplitudowe. Tworzy je zbiór prążków rozmieszczonych na osi częstotliwości, w odstępach równych 1/T, gdzie fg=(1/T) jest częstotliwością podstawową sygnału. Każdy prążek reprezentuje odpowiednią harmoniczną więc jego długość (wysokość) jest proporcjonalna do amplitudy danej harmonicznej.
Z własności szeregu Fouriera wynika, że widmo amplitudowe sygnału sinusoidalnego zawiera tylko jeden współczynnik Fouriera o wartości różnej od zera. W praktyce zjawisko to wykorzystuje się do oceny jakości „rzeczywistego przebiegu sinusoidalnego” czyli jeśli w widmie takiego przebiegu występuje „coś” jeszcze poza jednym prążkiem oznacza to, że sygnał jest np. odkształcony lub zaszumiony czego może nie być widać w dziedzinie czasowej.
Celem analizy widmowej sygnału okresowego jest określenie wartości współczynników Fouriera
dla kolejnych harmonicznych. Obecnie w praktyce zadanie to realizuje się w następujący sposób: sygnał okresowy jest próbkowany i kwantowany (przez przetworniki analogowo-cyfrowe), a uzyskany zbiór próbek wartości chwilowych zostaje poddany cyfrowej obróbce matematycznej, pozwalającej na obliczenie współczynników Fouriera.
Podstawowe wzory (bez wyprowadzeń matematycznych), przedstawiają się następująco.
I tak dla sygnałów okresowych widmo sygnału jest zbiorem współczynników Fouriera, których wartości można obliczyć ze wzoru
Należy zwrócić uwagę na fakt, że zwiększeniu okresu całkowania T w dziedzinie czasu towarzyszy wzrost rozdzielczości w dziedzinie częstotliwości. W granicznym przypadku, tzn. dla T= ∞ widmo staje się widmem ciągłym. Wtedy funkcja
może być nieokresowa, np. impuls pojedynczy. Ten przypadek opisuje całkowe przekształcenie Fouriera określone wzorem
Funkcja
jest nazywana zespoloną transformatą Fouriera lub widmem zespolonym, funkcji
. Moduł
oznacza widmo amplitudowe a argument widmo fazowe.
Jeżeli sygnał
ma ograniczone widmo do częstotliwości fm, to zgodnie z twierdzeniem Kotielnikowa-Shannona może on być jednoznacznie określony na podstawie próbek wartości chwilowych x(i) sygnału
, pobieranych z częstotliwością próbkowania fs spełniającą warunek fs ≥ 2fm. (x(i) jest to spróbkowany sygnał
a i oznacza numer próbki, czyli dyskretny czas).
Jeżeli dysponuje się nieskończoną liczbą próbek sygnału
(czyli sygnałem x(i) dla i z zakresu od -∞ do +∞) o ograniczonym widmie do częstotliwości fm, to - po uprzednim zastąpieniu całkowania sumowaniem i odpowiednich przekształceniach - można z powyższego wzoru wyznaczyć wszystkie współczynniki tego sygnału. W praktyce jednak zawsze dysponujemy ograniczoną liczbą próbek N. Intuicyjnie jest zrozumiałe, że z ograniczonej liczby próbek można wyznaczyć skończoną liczbę współczynników będących estymatorami współczynników Fouriera, które niekoniecznie będą odwzorowywać rzeczywiste wartości współczynników Fouriera. Opisuje to wzór, który nazywa się dyskretnym przekształceniem Fouriera (DPF lub DFT)
gdzie: F(n) - estymator prawdziwego widma
sygnału
,
n - dyskretna częstotliwość,
x(i)=x(iTs), Ts=1/fs,
i=0, ±1, ±2, ... - dyskretny czas
Na podstawie tego wzoru można wyznaczyć estymatory rzeczywistych współczynników Fouriera, ale tylko dla dyskretnych częstotliwości, które są całkowitą wielokrotnością podstawowej częstotliwości analizy f1
f1 =1/Tw=fs /N
gdzie: Tw - okno czasowe, czyli czas pomiaru Tw=N Ts.
Sygnał o nieskończonym czasie trwania od -∞ do +∞, próbkuje się w ograniczonym czasie. Dla matematycznego zapisu tego faktu wprowadza się funkcje w(i) zwaną funkcją okna danych. Ma ona taką własność, że jest równa 1 dla |i| ≤(N/2) i 0 dla |i|>(N/2), czyli przyjmuje wartość 1 w trakcie pomiaru i wartość 0 przed rozpoczęciem i po zakończeniu pomiaru.
w(i)=1 dla |i|≤N/2 i w(i)=0 dla |i|>N/2
Jak łatwo zauważyć dla skończonej liczby wartości x(i) czyli dla N próbek mierzonego sygnału
, a tym właśnie dysponujemy po wykonaniu pomiaru, iloczyn funkcji w(i) z sygnałem x(i) daje sygnał x(i).
w(i)x(i)=x(i) dla N∈(-∞,+∞)
Można teraz zapisać
Zastosowana funkcja okna nosi nazwę okna prostokątnego. Ma ona istotną wadę, która polega na tym, że w dziedzinie częstotliwości wprowadza tzw. efekt „przecieku” widma będący źródłem błędów estymacji widma (Rys. 2.). Można tego uniknąć w szczególnym przypadku tzn. wtedy gdy badanym sygnałem jest sygnał periodyczny którego częstotliwość fg jest całkowitą wielokrotnością częstotliwości analizy f1. Jest to równoznaczne z tym, że w oknie czasowym Tw mieści się całkowita liczba (k) okresów sygnału (N=Tw*fs=k*Tg*fs=k*fs/fg => N fg=k fs ).
W praktyce rzadko spełniony jest ten warunek, więc można przykładowo zwiększyć liczbę próbek N, ale daje to stosunkowo niewielkie korzyści. Innym sposobem, często skuteczniejszym, zmniejszenia efektu „przecieku” widma (Rys. 4.) jest zastosowanie np. różnych okien danych o „złagodzonych” brzegach (Rys. 3.).
Poza opisaną metodą estymacji widma, która polega na zastosowaniu wprost przekształcenia Fouriera stosuje się różne inne algorytmy cyfrowego przetwarzania sygnałów pozwalające na wyestymowanie widma w taki sposób aby w miarę najdokładniej wyeksponować interesującą nas informację zawartą w sygnale np. rozróżnienie częstotliwości dwóch leżących blisko siebie sygnałów sinusoidalnych.
Dokładność pomiaru składowych widma
Na błąd pomiaru, rozumiany jako różnica między rzeczywistą amplitudą składowej widma analizowanego sygnału a jej wyznaczoną amplitudą, ma wpływ wiele czynników. Podstawowe z nich to błąd kwantyzacji przetwornika analogowo-cyfrowego (A/C), periodyczność zawartych w widmie składowych w stosunku do czasu trwania okna czasowego Tw oraz parametry zastosowanego okna danych.
Błąd kwantyzacji jest cechą zastosowanego przetwornika i nie mamy na niego wpływu. Istotny może jeszcze być błąd obliczeń (ilość wykonywanych operacji, zaokrąglenia itd.) ale przyjmuje się że dokładność obliczeń jest wystarczająca, jeżeli błąd z tego tytułu jest mniejszy niż błąd spowodowany kwantyzacją. W praktyce stosuje się tzw. szybkie przekształcenie Fouriera SPF (potocznie określane z języka angielskiego FFT - Fast Fourier Transform), oparte zwykle na liczbie próbek będących całkowitą potęgą liczby 2. Liczba obliczeń jest wtedy mniejsza niż w zwykłym przekształceniu DPF (DFT), znacznie krótszy jest więc czas obliczeń i mniejsze błędy obliczeń.
Pozostałe błędy zależą od znajomości badanego sygnału, doboru parametrów pomiaru i doboru okna.
Literatura
FRANKIEWICZ I. i zespół, Miernictwo elektroniczne i elektryczne - Ćwiczenia laboratoryjne, ćw. 18, Politechnika Wrocławska skrypt, Wrocław 1992
BEAUCHAM K. G., Przetwarzanie sygnałów metodami analogowymi i cyfrowymi, Warszawa, WNT, 1978.
OPPENHEIM A. V., SHAFER R. W., Cyfrowe przetwarzanie sygnałów, Warszawa WKŁ, 1979.
WOJNAR A., Teoria sygnałów, Warszawa, WNT, 1980
BORODZIEWICZ W., JASZCZAK K., Cyfrowe przetwarzanie sygnałów, Warszawa WNT, 1887
Pytania kontrolne
Co oznacza spełnienie warunku Kotielnikowa-Shannona?
Co to jest widmo sygnału?
Z jaką częstotliwością należy próbkować przebieg prostokątny o częstotliwości fs= 1 kHz, aby spełnić warunek Kotielnikowa-Shannna jeśli pomiaru dokonujemy przetwornikiem o 12 bitowej rozdzielczości (4096 poziomów) na zakresie ±5V
Co oznaczają pojęcia: częstotliwość próbkowania - fs = 1/Ts, okno czasowe - Tw ? Jaki jest między nimi związek?
Jaki jest cel stosowania okien innych niż prostokątne?
Na co należy zwracać uwagę wykonując pomiar sygnału przetwornikiem A/C w celu estymacji widma?
Co to jest efekt „przecieku” widma i jak można się go pozbyć?
Program ćwiczenia
Obsługa programu komputerowego
Zapoznać się obsługą programu komputerowego AS (instrukcja obsługi w załączeniu) przeznaczonego do rejestracji i analizy sygnałów. Zwrócić uwagę na obsługę wykresów oraz wpisywanie parametrów w formatkach. W celu szybszego poznania programu można jako trening wybrać opcję SYGNAŁY\Generacja\Deterministyczne\Suma sinusów i na niej ćwiczyć (dla przyspieszenia obliczeń wpisać częstotliwość próbkowania 512 Hz i długość sygnału 1 sek).
UWAGI dotyczące obsługi programu AS:
a) Podczas wykonywania pomiarów w opcji pomiar muszą być zachowane następujące ustawienia:
-Adres bazowy = 0300
-Numer IRQ od karty = 5
-Zakres [V]= -5..+5
-Typ wyzwalania: Klawiatura
b) Większość algorytmów estymujących widmo w wykorzystywanym programie oparta jest na FFT a jak wspomniano we wprowadzeniu algorytm FFT wymaga liczby próbek będącą całkowitą potęgą liczby 2, więc jeśli nie jest spełniony ten warunek wtedy program uzupełnia zestaw próbek do wymaganej liczby próbkami o wartości zerowej. Przykładowo chcąc policzyć widmo dla sygnału o długości 500 próbek, program najpierw uzupełni ten zestaw próbkami o wartości zerowej do liczby 512 próbek. Oczywiście jeśli długość transformaty FFT wynosi 512, bo jeśli wynosi ona 1024 to program zwiększy długość sygnału odpowiednio do 1024 próbek.
c) Przed przystąpieniem do analizy (opcja ANALIZA) nowego sygnału należy zawsze najpierw wczytać taki sygnał (opcja SYGNAŁY\Lista).
d) Liczbę próbek w trakcie analizy można odczytywać z wykresów - jest to liczba n w prawym górnym rogu wykresu. Istotną informacją jest liczba próbek w trakcie wyświetlania widma. Jest ona zawsze całkowitą potęgą liczby 2 i można z niej zorientować się ile próbek zostało poddanych transformacie FFT. Mianowicie liczba próbek po transformacie FFT jest zawsze o połowę mniejsza niż liczba próbek sygnału wejściowego np. jeśli na wykresie przebiegu czasowego sygnału odczytaliśmy liczbę próbek n=500 a na wykresie widma odczytaliśmy liczbę próbek n=256 to świadczy o tym że do analizy wzięto 512 próbek czyli, że program dodał 12 próbek o wartości zerowej.
e) Rozciąganie i zawężanie wyświetlanych wykresów czyli tzw. lupę można wykonywać poprzez klawisze S i D, natomiast przewijanie (przesuw) identycznie tylko trzymając wciśnięty klawisz <CTRL> (patrz instrukcja).
f) Wyzwalanie pomiaru odbywa się klawiszem <SPACE>
g) Jeśli na stanowisku jest drukarka to można drukować wykresy wciskając klawisz F6.
h) Niektóre opcje programu wymagają dużej liczby operacji a co za tym idzie i czasu, więc należy wtedy po wciśnięciu <SPACE> cierpliwie czekać.
i) Program często wychodząc naprzeciw użytkownikowi sam zmienia i dostosowuje niektóre parametry zwłaszcza po zmianie częstotliwości próbkowania. Dlatego po dokonaniu jakichś zmian w formatce zalecane jest, przed wykonaniem obliczeń, wybranie przy pomocy strzałek każdego pola i upewnienie się, że parametry są takie jakich oczekujemy. Najlepiej zmieniać parametry w kolejności zgodnej z występowaniem w formatce (idąc od góry)
Ocena widma pojedynczego sygnału sinusoidalnego
Przy pomocy programu AS wygenerować pojedynczy sygnał sinusoidalny (opcja SYGNAŁY\Generacja\Deterministyczne\Suma sinusów) przy następujących warunkach:
częstotliwość próbkowania - 512 Hz oraz 500 Hz
długość sygnału - 1 sek.
czas t0 - 0
czas t1 - 1 sek. oraz 0.975 sek.
częstotliwość sygnału - 100 Hz oraz 100,5 Hz.
Pamiętać aby wyłączyć (zaznaczyć numer składowej i wcisnąć DEL lub INS) pozostałe składowe (jest ich 30; przesuw góra/dół klawiszami 2/8 lub strzałki na klawiaturze numerycznej). Dla wszystkich 8 przypadków zaobserwować widma (markerami (W, Z) sprawdzić częstotliwość) i opisać wpływ wprowadzanych zmian na ich kształt. Widmo prezentowane jest po wciśnięciu <SPACE> w trakcie wyświetlania przebiegu czasowego. Ocenić wygląd przebiegu czasowego w kontekście twierdzenia Kotielnikowa-Shanona, oraz sprawdzić jaki wpływ na widmo oraz przebieg czasowy będzie miała zmiana częstotliwości sygnału na 10 Hz (np. przy warunkach fs=500 Hz, długość 1 sek., t1=1 sek.) oraz na 250 Hz (np. przy warunkach fs=512 Hz, długość 1 sek., t1=1 sek.).
Zapisać do plików przypadki o widmie najbardziej przystającym do idealnego oraz najbardziej odbiegającym.
Wykonać pomiar sygnału sinusoidalnego z generatora o częstotliwości ok. 100 Hz (opcja POMIAR\Charakterystyki czasowe). Dobrać następujące warunki pomiaru
Zapis: DO PAMIĘCI
Częstotliwość próbkowania: 500 Hz oraz 512 Hz
Liczba próbek: 1000
Pomiar uruchamia się klawiszem SPACE a po jego zakończeniu należy wcisnąć klawisz F5. Pomiary zapisać w dwóch różnych plikach (wpisać swoją nazwę zbioru). Wykonać konwersję obu plików z formatu MBI na SBS (opcja DANE\Konwersja\ MBI na SBS), następnie wczytać przekonwertowane sygnały (opcja SYGNAŁY\Lista) i wykonać analizę widmową (opcja ANALIZA\ Charakterystyki\Widmo amplitudowe). Zmieniać długość sygnału i długość transformaty FFT tak aby raz długość transformaty była większa od długości sygnału, a raz aby były równe. Program pozwala na wybór długości transformaty FFT większej lub równej długości sygnału, dlatego chcąc dobrać jej długość równą długości sygnału czasami konieczne jest zmniejszenie najpierw tej drugiej. Ocenić i uzasadnić wpływ powyższych zmian na widmo. Czy można w widmie wyróżnić jakieś inne harmoniczne sygnału poza podstawową?.
Wpływ funkcji okna na kształt widma
Wczytać (opcja SYGNAŁY\Lista) kolejno sygnały zapisane w punkcie 2.1 (najlepszy -fs=512 Hz, dług.=1 sek., fsygnału=100 Hz, t1=1 sek.; najgorszy - fs=500 Hz, dług.=1 sek., fsygnału=100,5 Hz, t1=0,975 sek.) i wyestymować ich widma (opcja ANALIZA\Estymacja Widma\Cooley-Tukey) przy zastosowaniu następujących okien: Lanczos,Hanning-Poisson i Bartlet (zmiana okien INS/DEL).
Zachować następujące ustawienia
Parametr okna: 2
Długość płaskiej części okna: 0 %
Długość transformaty FFT: równa lub minimalnie większa od długości sygnału
Zwracać uwagę na to ile próbek badanego sygnału brana jest do analizy - powinna być taka aby obejmowała wszystkie próbki danego sygnału (najlepszy - 512, najgorszy 500)
Po wpisaniu właściwych parametrów wciskać <SPACE> a na ekranie pojawiać będą się kolejno: przebieg czasowy sygnału (jak na Rys. 1a), kształt funkcji okna (jak na Rys. 3a), widmo i parametry funkcji okna (jak na Rys. 3b), iloczyn sygnału i funkcji okna (jak na Rys. 4a) oraz widmo sygnału po zastosowaniu funkcji okna (jak na Rys. 4b).
Wczytać kolejno sygnały zapisane w punkcie 2.2 i przeprowadzić estymacje ich widm jak w punkcie 3.1. Co można powiedzieć o sygnale na podstawie jego widma? Przy jakich parametrach otrzymujemy widmo o kształcie najbardziej zbliżonym do oczekiwanego?
Analiza widma sygnału prostokątnego (*dla studentów lepiej przygotowanych)
Wygenerować sygnał prostokątny sumując sygnały sinusoidalne (opcja SYGNAŁY\Generacja\Deterministyczne\Suma sinusów) przy następujących warunkach: częstotliwość próbkowania 512 Hz, długość sygnału 1 sek, częstotliwości kolejnych składowych n*10 Hz, amplitudy kolejnych składowych 1/n, maksymalna częstotliwość harmonicznej 250 Hz. Pamiętać że n musi być nieparzyste.
Wygenerować sygnał prostokątny (opcja SYGNAŁY\Generacja\Deterministyczne\Suma prostokątów) przy następujących warunkach: częstotliwość próbkowania 512 Hz, długość sygnału 1 sek, amplitudy kolejnych składowych ±1, długości kolejnych składowych 0.1 sek, liczba składowych 10 co daje 5 okresów.
Zmierzyć sygnał prostokątny o częstotliwości ok. 10 Hz przy częstotliwości próbkowania 512 Hz i liczbie próbek 512.
Obliczyć widma sygnałów z punktów 4.1, 4.2 i 4.3 różnymi metodami oraz zmieniając wybrane parametry zbadać ich wpływ na kształt widma. Zmierzyć markerami amplitudy kolejnych harmonicznych. Czy zmniejszają się one jak 1/n? Co można powiedzieć o sygnale z punktu 4.3 i jego widmie w kontekście analizy widm sygnałów z punktów 4.1 i 4.2.
UWAGA:
Do analizy widm sygnałów ważne są następujące parametry:
Częstotliwość próbkowania,
Kształt sygnału (przynajmniej przybliżony),
Częstotliwość sygnału,
Liczba próbek poddawanych analizie - zamiennie długość sygnału,
Liczba próbek po analizie (w widmie) - zamiennie długość transformaty FFT,
Częstotliwości występujących w widmie prążków,
Amplitudy sygnału i jego ewentualnych składowych,
Kształt widma (przynajmniej przybliżony)
Informacja o zastosowanym oknie
10
Zb. Świerczyński
Opracował: mgr inż. Zbigniew Świerczyński 99-10-04