Widma częstotliwościowe prostych sygnałów
Tworzenie aplikacji analizatora widma najlepiej rozpocząć od wykonania czegoś prostego. Ponieważ interesuje nas widmo częstotliwościowe to warto sprawdzić jak wygląda ono dla kilku prostych sygnałów. W tym celu możemy wygenerować cztery sygnały, które zawierają od jednej do czterech składowych harmonicznych:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 |
if(mode == 1){ for(i=0; i<512; ++i){ buffer_input[i] = (float32_t) 30*sin(2*PI*freq*i*dt); } } else if(mode == 2){ for(i=0; i<512; ++i){ buffer_input[i] = (float32_t) 20*sin(2*PI*freq*i*dt) + 10*sin(2*PI*2*freq*i*dt); } } else if(mode == 3){ for(i=0; i<512; ++i){ buffer_input[i] = (float32_t) 15*sin(2*PI*freq*i*dt) + 10*sin(2*PI*2*freq*i*dt) + 5*sin(2*PI*3*freq*i*dt); } } else if(mode == 4){ for(i=0; i<512; ++i){ buffer_input[i] = (float32_t) 15*sin(2*PI*freq*i*dt) + 10*sin(2*PI*2*freq*i*dt) + 3*sin(2*PI*3*freq*i*dt) + 2*sin(2*PI*4*freq*i*dt); } } |
Generowanych jest 512 próbek sygnału i liczba ta wynika z wymagań narzuconych przez funkcję wyznaczającą widmo i chęci przedstawienia jak najdłuższego przebiegu na wyświetlaczu (wykorzystane będą tylko pierwsze 256). Wartość zmiennej freq w toku działania aplikacji demonstracyjnej jest zmieniana od 1 do 1024, a zmienna dt jest wirtualnym okresem próbkowania, który wynosi 0.001. Dalej można wyświetlić wygenerowane sygnały:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 |
// Usunięcie poprzedniego przebiegu sygnału LCD_SetForegroundColor(Black); for(i=0; i<255; ++i){ x_LCD_DrawLine_2(i + 32, (uint16_t)(buffer_input_copy[i] + 50), i + 33, (uint16_t)(buffer_input_copy[i+1] + 50)); } // Rysowanie nowego przebiegu sygnału LCD_SetForegroundColor(Yellow); for(i=0; i<255; ++i){ x_LCD_DrawLine_2(i + 32, (uint16_t)(buffer_input[i] + 50), i + 33, (uint16_t)(buffer_input[i+1] + 50)); } // Wyświetlenie informacji o częstotliwościach składowych harmonicznych sygnału if(mode == 1){ sprintf(text,"F=%4d [1]",freq); x_LCD_DrawText(32, 100, text, LeftJustify, NormalOrientation); } else if(mode == 2){ sprintf(text,"F1=%4d [1], F2=%4d [1]",freq,freq*2); x_LCD_DrawText(32, 100, text, LeftJustify, NormalOrientation); } else if(mode == 3){ sprintf(text,"F1=%4d [1], F2=%4d [1]",freq,freq*2); x_LCD_DrawText(32, 100, text, LeftJustify, NormalOrientation); sprintf(text,"F3=%4d [1]",freq*3); x_LCD_DrawText(32, 110, text, LeftJustify, NormalOrientation); } else if(mode == 4){ sprintf(text,"F1=%4d [1], F2=%4d [1]",freq,freq*2); x_LCD_DrawText(32, 100, text, LeftJustify, NormalOrientation); sprintf(text,"F3=%4d [1], F4=%4d [1]",freq*3,freq*4); x_LCD_DrawText(32, 110, text, LeftJustify, NormalOrientation); } // Zachowanie kopii pierwszej połowy próbek aktualnego sygnału for(i=0; i<256; ++i){ buffer_input_copy[i] = buffer_input[i]; buffer_output_mag_copy[i] = buffer_output_mag[i]; } |
Mając przygotowany sygnał wykorzystujemy trzy funkcje z biblioteki DSP i wyznaczamy najpierw transformatę Fouriera danego sygnału, a następnie obliczamy wartości modułów dla każdej pary liczb rzeczywistej i urojonej. Na koniec możemy wyszukać składową sygnału (częstotliwość), która ma największą amplitudę i przeskalować tablicę modułów do zakresu nadającego się do wyświetlenia:
|
1 2 3 4 5 6 7 8 9 10 |
// Wyznaczenie transformaty Fouriera arm_rfft_f32(&S, buffer_input, buffer_output); // Obliczenie modułów arm_cmplx_mag_f32(buffer_output, buffer_output_mag, 512); // Znalezienie składowej harmonicznej sygnału o największej amplitudzie arm_max_f32(buffer_output_mag, 512, &maxvalue, &maxvalueindex); // Skalowanie wartości modułów for(i=0; i<512; ++i){ buffer_output_mag[i] = 100*buffer_output_mag[i]/maxvalue; } |
Pierwsza funkcja przyjmuje trzy argumenty: wskaźnik do struktury typu arm_rfft_instance_f32 opisującej sposób działania transformaty FFT, wskaźnik do buforu wejściowego z wartościami próbek sygnału oraz wskaźnik do buforu wyjściowego w którym będą zapisane zespolone wartości transformaty sygnału. Ogólną ideę funkcjonowania tej funkcji pokazano na rysunku 9. Ważne jest aby w buforze wejściowym były zapisane rzeczywiste wartości (UWAGA: zostaną one zmodyfikowane po wykonaniu się funkcji), czyli:
|
1 |
{real[0], real[1], real[2], real[3], ..} |
Natomiast w buforze wyjściowym, który powinien być dwa razy większy od wejściowego, oprócz wartości rzeczywistych będą znajdować się również wartości urojone:
|
1 |
{real[0], imag[0], Real[1], imag[1], ...} |
W przypadku struktury przekazywanej w pierwszym argumencie nie ma potrzeby ręcznego uzupełniania jej, gdyż można do tego celu wykorzystać funkcję arm_rfft_init_f32():
|
1 2 3 4 5 6 7 8 9 |
// Wystarczą deklaracje dwóch struktur arm_rfft_instance_f32 S; arm_cfft_radix4_instance_f32 S_CFFT; // oraz inicjalizacja, gdzie: // arg#3 – liczba próbek w buforze wejściowym; możliwe wartości to: 128, 512, 2048 // arg#4 – kierunek transformacji: 0 – prosta, 1 – odwrotna // arg#5 – uporządkowanie wartości w buforze wyjściowym: 0 – odwrócone, 1 - normalne arm_rfft_init_f32(&S, &S_CFFT, 512, 0, 1); |
Rys. 9. Idea działania funkcji arm_rfft_f32()
Druga funkcja przyjmuje również trzy argumenty: wskaźnik do buforu wejściowego z wartościami par liczb rzeczywistych i urojonych, wskaźnik do buforu wyjściowego do którego zostaną zapisane wartości rzeczywiste oraz ostatni parametr – ilość liczb zespolonych. Funkcja oblicza moduł dla każdej pary wartości liczby rzeczywistej i urojonej zgodnie ze wzorem:
Trzecia funkcja znajduje maksymalną wartość w tablicy buffer_output_mag o długości 512 elementów i zapisuje ją do zmiennej maxvalue, a indeks komórki w jakiej się ona znajduje jest zapisywany w zmiennej maxvalueindex.
Fot. 10. Sygnał o jednej składowej harmonicznej i jego widmo
Fot. 11. Sygnał o dwóch składowych harmonicznych i jego widmo
Fot. 12. Sygnał o trzech składowych harmonicznych i jego widmo
Fot. 13. Sygnał o czterech składowych harmonicznych i jego widmo
Na koniec, po przeskalowaniu wartości modułów, można wyświetlić otrzymane wyniki (fotografie 10…13). Wystarczy przedstawić pierwsze 256 elementów tablicy buffer_output_mag, gdyż druga połowa stanowi lustrzane odbicie pierwszej:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 |
// Usunięcie poprzedniego widma LCD_SetForegroundColor(Black); for(i=0; i<256; ++i){ x_LCD_DrawLine_2(i + 32, 130 + 100, i + 32, 130 + (uint16_t)(100 - (uint16_t)buffer_output_mag_copy[i+1])); } // Rysowanie nowego widma LCD_SetForegroundColor(Green); for(i=0; i<256; ++i){ x_LCD_DrawLine_2(i + 32, 130 + 100, i + 32, 130 + (uint16_t)(100 - (uint16_t)buffer_output_mag[i+1])); } |
