Wprowadzenie
Wyobraź sobie, że trzymasz w ręku gitarę, pociągasz za strunę, a następnie słyszysz dźwięk. Dźwięk ten nazywamy tonem i jest falą akustyczną wytworzoną przez drgającą strunę która dostaje się do naszego układu słuchowego. Gdy pociągasz za inną strunę, słyszysz inny ton, oznacza to, że jakieś parametry musiały ulec zmianie. Gdy nagramy taki ton, otrzymamy wykres, który będzie przedstawiał, jak zmienia się amplituda tego tonu w czasie. Po szarpnięciu struny widzimy, jak struna ta odchyla się coraz słabiej oraz dodatkowo słyszymy, że ton ten staje się coraz cichszy. Oznacza to, że amplituda na takim wykresie będzie malała. Wracając do struny, po jej szarpnięciu dodatkowo widzimy, że struna ta, odchyla się do przodu i do tyłu z określoną częstotliwością. Słyszymy też, że ton ten jest cały czas taki sam, co oznacza, że częstotliwość drgania struny nie zmienia się w czasie. Gdy zerkniemy na tabelę częstotliwości tonów [1], możemy zauważyć, że ton a1 posiada częstotliwość równą 440[Hz]. Wizualizację takiego tonu przedstawiliśmy na Rys. 1. W tym miejscu należy zaznaczyć, że jest to ton wygenerowany sztucznie, dlatego między innymi nie widzimy tu tłumienia czy też szumów.
Gdy teraz jeden podzielimy przez częstotliwość tonu a1, otrzymamy okres, który w przybliżeniu jest równy 0,0023 s. Co przedstawia Rys. 2.
Gdybyśmy chcieli teraz przedstawić jak ton ten wygląda w funkcji częstotliwości to dla idealnego przypadku, wykres wyglądałby jak na Rys. 3.
Dla powyższego przypadku sytuacja jest prosta, jednakże wyobraź sobie, że słuchasz teraz piosenki. Taka piosenka składa się z setek pojedynczych tonów. Nagrywając jakiś utwór, przebieg w funkcji czasu dla takiego przypadku będzie wyglądał mniej więcej jak na Rys. 4.
Natomiast sytuacja komplikuje się, jeżeli ten sam przebieg chcielibyśmy przedstawić w funkcji częstotliwości. W tym przypadku przychodzi nam z pomocą DFT, czyli Dyskretna Transformata Fouriera [2]. Natomiast złożoność DFT wynosi N^2 przez to obliczanie w czasie rzeczywistym mogłoby być bardzo trudne, a nawet niemożliwe dla mikrokontrolerów. Dlatego w tym miejscu przychodzi nam FFT, czyli szybka transformata Fourier, która robi to samo, lecz znacznie szybciej. Złożoność tego algorytmu wynosi N logN. Wspomniane algorytmy opierają się na idei, że każdy sygnał można przedstawić jako sumę fal sinusoidalnych o różnych częstotliwościach, czyli szuka tych składowych i jednocześnie mówi, jak są one silne.
Wizualizacja FFT w praktyce — Audacity
Znając już podstawy teorii, sprawdźmy jak to będzie wyglądało w praktyce. W tym celu postanowiliśmy użyć narzędzia jakim jest Audacity. Audacity jest narzędziem opensource i jest udostępniane na licencji GPL [3] więc każdy bez problemu może korzystać z tego oprogramowania. Do naszej analizy wybrałem następujących 6 sygnałów:
- sinus – 1024Hz Rys 5.
- sinus – 8192Hz Rys 5.
- prostokąt – 1024Hz Rys 6.
- prostokąt – 8192Hz Rys 6.
- trójkąt – 1024Hz Rys 7.
- trójkąt – 8192Hz Rys 7.
Mając przygotowane sygnały w dziedzinie czasu, w oprogramowaniu należy wczytać te pliki i zaznaczyć część albo całość przebiegu, a następnie wybrać opcję “Narysuj widmo…”, która znajduje się w zakładce “Analizuj” Rys 8.
Po przeanalizowaniu wybranych sygnałów możemy stwierdzić, że dla przebiegów sinusoidalnych widoczna jest pojedyncza dominująca częstotliwość, która jest zgodna z częstotliwością naszych wybranych sygnałów Rys 9. oraz Rys 10.
W przypadku przebiegów dla sygnałów prostokątnych sytuacja jest inna. Dla każdego sygnału prostokątnego przedstawiliśmy dwa wykresy. Jeden wykres przedstawiony jest w skali logarytmicznej, a drugi w skali liniowej. Na wykresie w skali logarytmicznej łatwo nam zlokalizować prążek z dominującą częstotliwością, czyli 1024Hz na Rys 11. oraz 8192Hz na Rys 13. Za to na wykresach ze skala liniową Rys 12 i Rys 14. możemy zauważyć pewne zachowanie, która pokazuje nam, że przebieg prostokątny zawiera “prawie wszystkie” częstotliwości. Dodatkowo jak się przyjrzymy, możemy zauważyć, że co pewien okres, występują silniejsze i słabsze częstotliwości, które z coraz to wyższymi częstotliwościami można powiedzieć, że osiągają stałą wartość.
W przypadku sygnałów trójkątnych również zamieściłem po dwa wykresy dla jednego sygnału. Na jednym i drugim wykresie możemy mniej więcej wskazać jaka jest dominująca częstotliwość. Dodatkowo jeżeli przyjrzymy się wykresom w skali liniowej to możemy zauważyć, że siła kolejnych częstotliwości liniowo zmierza ku zeru.
FFT na mikrokontrolerze STM32F401VCT
Po wykonaniu FFT w programie Audiacity przejdźmy do implementacji FFT na mikrokontrolerze, na którym to w czasie rzeczywistym będziemy wykonywali FFT. FFT będzie obliczane na sygnale podanym na wejście analogowe z generatora.
Kroki wdrożenia biblioteki CMSIS do projektu:
- Pobieramy plik nagłówka: arm_math.h z oficjalnego repozytorium Github ARM Include.
- Następnie pobieramy plik obiektowy, w moim przypadku: libarm_cortexM4lf_math.a, także z oficjalnego repozytorium Github ARM GCC. Pobrałem plik lf, ponieważ wspiera sprzętowy FPU (Floating Point Unit).
- Wklejamy plik nagłówkowy w folder Core/Inc naszego projektu.
- Dodajemy biblioteki arm_math w sekcji Include oraz wybieramy procesor, w moim przypadku ARM_MATH_CM4, ponieważ używam STM32F4.
- Umieszczamy plik libarm_cortexM4lf_math.a w folderze, gdzie zapisany jest nasz projekt. Stworzyłem nowy folder w projekcie i tam przeniosłem plik libarm.
- Ustawiamy linker, gdzie dodałem ścieżkę do biblioteki oraz samą bibliotekę libarm.
- Ustawienia timera aby wyzwalał konwersję przetwornika z pożądaną częstotliwością, u mnie 80kHz, do czego użyłem następującego wzoru:
8. Ustawienia przetwornika ADC, źródłem wyzwalania konwersji jest timer 2:
Sygnał z timera 2 informujący o częstotliwości próbkowania:
Częstotliwość próbkowania została ustawiona na 80 kHz poprzez wyzwalanie konwersji przetwornika ADC za pomocą timera 2:
Przebieg z generatora funkcyjnego o częstotliwości 10kHz, na którym będzie wykonywane FFT:
Bufor wejściowy do FFT jest to wycięty fragment przebiegu, oknem prostokątnym, składający się z 256 próbek, przesunięty o składową stałą w celu wyznaczenia największej energii z częstotliwości przebiegu.
Wykorzystując okno prostokątne oraz 2^n okresów wyciek widma nie występuje, co można zaobserwować na Rys. 26. Podczas jednego okresu sinusoidy o częstotliwości 10kHz zbieranych jest 8 próbek, czyli w buforze składającym się z 256 próbek będzie 32 okresy.
Najczęstsze problemy i pułapki
Pomimo że algorytm FFT jest bardzo przydatnym narzędziem to jednak łatwo możemy otrzymać błędne wyniki.
- Jeśli sygnał zawiera częstotliwości wyższe niż połowa częstotliwości próbkowania (Nyquist), wówczas pojawią się fałszywe niższe częstotliwości.
- FFT zakłada, że analizowany fragment sygnału zaczyna się i kończy “periodycznie”, gdzie w praktyce tak nie jest. W takim przypadku energia rozlewa się na sąsiednie częstotliwości. Nazywamy to wyciekiem widma (spectral leakage) [4].
okna
- Jeśli długość analizowanego sygnału jest zbyt mała, nie jesteśmy w stanie rozróżnić blisko położonych częstotliwości.
- Obecność szumu powoduje, że widmo staje się „rozmyte”, a identyfikacja rzeczywistych składowych częstotliwościowych jest trudniejsza.
Temat ten jest bardzo szeroki i nie będziemy omawiać go w tym wpisie. Natomiast należy pamiętać, że między innymi z takimi problemami możemy się spoktać podczas stosowani algorytmu FFT.
- Szerzej o FFT
Na naszym blogu, opisywaliśmy także FFT oraz implementację na mikrkontroerze STM32F4, zapraszamy do zapoznania się:
- Obliczenie FFT na mikrokontrolerze STM32 – #1 – Biblioteka CMSIS, wyciek widma, okna wykorzystywane do FFT, wdrożenie w projekcie
- Obliczenie FFT na mikrokontrolerze STM32 – #2 – Program realizujący FFT, prezentacja wyników, wyciek widma oraz jego przyczyny
- Obliczenie FFT na mikrokontrolerze STM32 – #3 – Pomiar trzech przebiegów za pomocą niezależnych przetworników analogowo-cyfrowych z wykorzystaniem DMA i trybu Combined regular simultaneous + alternate trigger mode w celu dokładnego pomiaru przesunięcia fazowego na procesorze STM32F103
Podsumowanie
We wpisie przedstawiliśmy podstawy analizy sygnałów w dziedzinie częstotliwości, zaczynając od intuicyjnego przykładu dźwięku i jego reprezentacji w czasie oraz częstotliwości. Następnie omówiliśmy różnicę między DFT a FFT oraz uzasadniliśmy użycie FFT w systemach wbudowanych ze względu na niższą złożoność obliczeniową.
W części praktycznej zaprezentowaliśmy analizę widmową różnych sygnałów (sinusoidalnych, prostokątnych i trójkątnych) z wykorzystaniem programu Audacity, wskazując charakterystyczne cechy ich widm.
Kolejno opisaliśmy proces implementacji FFT na mikrokontrolerze STM32F401 z użyciem biblioteki CMSIS, obejmujący konfigurację środowiska, timera, przetwornika ADC oraz przygotowanie danych wejściowych do analizy.
Na końcu przedstawiliśmy najczęstsze problemy związane z FFT, takie jak aliasing, wyciek widma, ograniczona rozdzielczość częstotliwościowa oraz wpływ szumu.
W efekcie pokazano, że FFT jest skutecznym narzędziem do analizy sygnałów w czasie rzeczywistym na mikrokontrolerach, pod warunkiem poprawnej konfiguracji systemu oraz uwzględnienia ograniczeń metody.
Źródła:
[1] http://www.fizykon.org/muzyka/muzyka_tabela_czestotliwosci_tonow.htm
[2] https://esezam.okno.pw.edu.pl/mod/book/view.php?id=34&chapterid=700
[3] https://manual.audacityteam.org/man/license.html
[4] https://dmbp.pl/obliczenie-fft-na-stm32f4/
Autor: Mateusz Pluta & Dominik Bednarski