logo elektroda
logo elektroda
X
logo elektroda
Adblock/uBlockOrigin/AdGuard mogą powodować znikanie niektórych postów z powodu nowej reguły.

Do doświadczonych-Algorytm Goertzela, pytania.

lukashb 15 Lis 2011 22:22 3765 9
  • #1 10143451
    lukashb
    Poziom 39  
    Witam! Niedawno miałem na forum temat n/t próbkowania i zapisu do tablicy sygnału sinusoidalnego. Gdy się z tym uporałem, wziąłem się za pisanie algorytmu Goertzela. Szukalem w sieci dużo. Między innymi ta biblia:
    http://www.mstarlabs.com/dsp/goertzel/goertzel.html
    a także ta:
    http://www.eetimes.com/design/embedded/4024443/The-Goertzel-Algorithm
    Koniec końców znalazłem jakąś implementację na PIC, przerobiłem ją na swoje potrzeby (przeliczenie stałych coeff i zapisanie ich w programie - po co liczyć to "na żywo" skoro można je stablicować? chyba, że w tym jest jakiś cel? ).
    W zasadzie co do próbkowania niemam zastrzeżeń. Dla podawanego np. 1kHz mamy 8 punktów na okres sygnału czyli 8000/1000 = 8 (fsamplerate/freqsignal=x pkt). Dla potomnych, jak wygląda owa sprawa w Excelu po przesłaniu tablicy UARTEM:
    Do doświadczonych-Algorytm Goertzela, pytania.

    Tablica docelowo dla prób posiada 100 elementów, w obliczeniach jest tak podstawiane do wzorów na K i na Coeff.
    Jednak nurtuje mnie inne pytanie. Jak powinien się zachowywać wynik obliczeń (magnituda)w przypadku:
    a)sygnału sinusoidalnego o poprawnej częstotliwości
    b)sygnału sinusoidalnego ale o innej częstotliwości (poza pasmem bin width)
    c)mowy, szumu i ogólnych zakłuceń

    Co do pierwszego to zapewne magnituda odpowie gwałtownym wzrostem gdy zbliżamy się do częstotliwości zadeklarowanej, ale co w przypadku B i C?
    Pytam ponieważ kod który napisałem, nie bardzo działa. Tzn, wiadomo, ADC bez sygnału podaje liczbę 512 - 520, magnituda wylicza się na 2052 i tak w sumie pozostaje. Gdy podaję teraz sygnał gwałtownie wzrasta do ok. 64 tysięcy, ale spada do wartości nie kiedy zerowych, czyli
    początek sygnału 2052/2052/2052/11432/45322/55634/64533/34453/23453/1100/0/2300/2052/2052052 koniec sygnału
    .Nie dziwiło by mnie to zachowanie gdyby nie fakt, że tak samo zachowuje się dla częstotliwości innej od zadeklarowanej, a z tego co wiem nie powinno być odpowiedzi tak znacznej (jakaś może być po to się deklaruje próg zadziałania by uzyskać odpowiednią energię). Załączam kod:
    Kod: C / C++
    Zaloguj się, aby zobaczyć kod

    W zasadzie, nie wiem czy dobry dział, ale to podchodzi pod przetwarzanie sygnałów, bo jako takich problemów z programowaniem tegoż nie mam (tak mi się wydaje). Tu algorytm przerobiony dla jednego elementu w tablicy (element 0) dla częstotliwości 1124Hz. Proszę o jakąś podpowiedź, dlaczego tak się to zachowuje. Dziękuję za zainteresowanie. Pozdrawiam.
  • Pomocny post
    #2 10144397
    nibbit
    Poziom 20  
    Patrząc jak to jest FFT to składowa stała jest zawsze na zerowym prążku. Dziwne, że podawanie 0Hz ma taki wpływ na detekcję wyższych częstotliwości. Dla testu możesz podać na stałe wartość, np 512 i zobaczyć jak się zachowa algorytm ale jak na moje powinieneś dostać 0.

    Inna sprawa to używanie double w embedded to porażka jeśli chodzi o czas wykonywania. Spójrz lepiej jak wygląda implementacja algorytmu na liczbach całkowitych (http://www.keil.com/download/docs/325.asp). Btw u Ciebie prędkość nie ma znaczenia bo i tak zatrzymujesz timer.

    Teraz nie mam dostępu do matlaba ale wieczorem może sprawdzę jaka jest odpowiedź algorytmu na skladową stałą o różnej wartości (np.0 i 512).
  • #3 10145719
    lukashb
    Poziom 39  
    nibbit napisał:
    Patrząc jak to jest FFT to składowa stała jest zawsze na zerowym prążku. Dziwne, że podawanie 0Hz ma taki wpływ na detekcję wyższych częstotliwości. Dla testu możesz podać na stałe wartość, np 512 i zobaczyć jak się zachowa algorytm ale jak na moje powinieneś dostać 0.

    Inna sprawa to używanie double w embedded to porażka jeśli chodzi o czas wykonywania. Spójrz lepiej jak wygląda implementacja algorytmu na liczbach całkowitych (http://www.keil.com/download/docs/325.asp). Btw u Ciebie prędkość nie ma znaczenia bo i tak zatrzymujesz timer.

    Teraz nie mam dostępu do matlaba ale wieczorem może sprawdzę jaka jest odpowiedź algorytmu na skladową stałą o różnej wartości (np.0 i 512).

    Po kolei. ADC bez sygnału sam podaje 512, bo tak jest zrobiony dzielnik do przesunięcia poziomu masy do pomiaru sinusa - by obie połówki się załapały.
    Gdy jest 512 to wynik dostaję 2052. Double używam bo double to zmienna liczb rzeczywistych. Problemem jest taka sama odpowiedź algorytmu na różne częstotliwości. W zasadzie taka sama, bo tak czy inaczej skacze do 64 tysięcy wynik. Zasymuluj jak możesz w MatLabie algorytm bo sam jestem ciekaw co z tego wyjdzie. Pozdrawiam
  • Pomocny post
    #4 10145863
    nibbit
    Poziom 20  
    No więc tak jak mówiłem, gdy podaję wartość stałą to wszędzie jest 0 nie licząc prążka na 0Hz. A Ty jak sam pisałeś masz około 512 bo coś tam skacze na wejściu i dlatego chciałem żebyś sprawdził co się stanie gdy wpiszesz na sztywno 512.

    Druga sprawa, wynik ładujesz do tablicy unsigned char. Na jakim procku to robisz? Bo w DSP char jest 16-bit ale jeśli nie to 8-bit więc to nie ma prawa działać skoro przetwornik masz 10-bit.
    I jak już z tym dojdziesz do ładu to chyba, żeby to w tej wersji algorytmu działało to będziesz musiał dane z przetwornika wpisywać do tablicy zmiennych ze znakiem i odejmować 512, tak by sygnał miał -512..511.
  • #5 10146812
    lukashb
    Poziom 39  
    nibbit napisał:
    No więc tak jak mówiłem, gdy podaję wartość stałą to wszędzie jest 0 nie licząc prążka na 0Hz. A Ty jak sam pisałeś masz około 512 bo coś tam skacze na wejściu i dlatego chciałem żebyś sprawdził co się stanie gdy wpiszesz na sztywno 512.

    Druga sprawa, wynik ładujesz do tablicy unsigned char. Na jakim procku to robisz? Bo w DSP char jest 16-bit ale jeśli nie to 8-bit więc to nie ma prawa działać skoro przetwornik masz 10-bit.
    I jak już z tym dojdziesz do ładu to chyba, żeby to w tej wersji algorytmu działało to będziesz musiał dane z przetwornika wpisywać do tablicy zmiennych ze znakiem i odejmować 512, tak by sygnał miał -512..511.

    Tak jest :D Coś mi się zaczyna rozjaśniać. Znaczy to tyle, że Goertzel wymaga (no jasne) sinusa zapisanego dyskretnie, ale nie w formie jak ja mam czyli powiecmy od 480 (ujemny szczyt sinusa) do 530 (analogicznie, dodatni szczyt sinusa) tylko względem poziomu 0 np. -120,+95 itp. Czyli od wartości odczytanej z rejestru ADC odejmuję 512 i wpisuje do tablicy. Trzeba zmienić będzie także z unsigned na signed by zmieścić się od -128 do 127. Wtedy przy ewidentnym braku sygnału mam zawsze 0 co by nie było. Oczywista sprawa to zmiana ADC na 8 bit i odczyt z ADCH, czyli będzie wtedy wynik nie 0..1023 tylko 0..255, ale to nie będzie problemem, powinno dać sobie radę. Dobrze, sprobuję to zastosować na żywo, najpierw z tablicą i -512 i ADC 8 bit a potem cały Goertzel. To musi działać, jest tak proste, że nie może być inaczej :D Pozdrawiam
  • #6 11545743
    sigmalko
    Poziom 10  
    Napisałem program, który pobiera z mikrofonu po 205 próbek (1 bajt na sampla, 8000 Hz próbkowania, mono) i liczy dla nich:
    - algorytmem goertzela wielkości dla poszczególnych częstotliwości z listy LOW i HIGH dla DTMF
    - max amplitudę dla 205 próbek
    - całkowitą energię dla 205 próbek, czyli dodaję do siebie kwadraty poszczególnych próbek

    Na ekranie prezentuję tylko tylko te sytuacje w której zarówno z grupy LOW jak i HIGH, największa wielkość z goertzela jest przynajmniej 2.5 raza większa od drugiej najwiekszej z grupy.

    Wyniki jednak mnie nie zadawalają bo oprócz DTMF'ów wykrywają mi się też inne dźwięki.


    Gdy mówię sobie "raz dwa trzy":
    =====================================================================
    totalEnergy=17348, maxAmplitude=17, low(159, 51), high(68, 24),
    totalEnergy=13584, maxAmplitude=19, low(192, 42), high(90, 32),
    totalEnergy=3404, maxAmplitude=12, low(161, 12), high(23, 7),
    totalEnergy=19274, maxAmplitude=23, low(111, 36), high(65, 24),
    totalEnergy=8999, maxAmplitude=15, low(139, 51), high(106, 39),

    Gdy sobie gwiżdżę:
    =====================================================================
    totalEnergy=1775, maxAmplitude=5, low(8, 3), high(414, 11),
    totalEnergy=1088, maxAmplitude=4, low(13, 5), high(314, 21),
    totalEnergy=1761, maxAmplitude=5, low(13, 4), high(412, 20),
    totalEnergy=1821, maxAmplitude=5, low(11, 4), high(427, 6),

    Gdy puszczę jakąś muzykę instrumentalną
    =====================================================================
    totalEnergy=17488, maxAmplitude=25, low(242, 68), high(618, 224),
    totalEnergy=15711, maxAmplitude=23, low(250, 95), high(585, 231),
    totalEnergy=13164, maxAmplitude=19, low(459, 40), high(170, 25),
    totalEnergy=12989, maxAmplitude=24, low(290, 50), high(363, 62),

    Gdy wcisnę "1" blisko mikrofonu
    =====================================================================
    totalEnergy=160569, maxAmplitude=46, low(592, 28), high(4003, 61),
    totalEnergy=186092, maxAmplitude=52, low(767, 57), high(4291, 15),
    totalEnergy=195950, maxAmplitude=52, low(803, 41), high(4402, 12),
    totalEnergy=209440, maxAmplitude=54, low(820, 44), high(4553, 37),
    totalEnergy=218549, maxAmplitude=54, low(812, 44), high(4657, 25),

    Gdy wcisnę "1" daleko mikrofonu
    =====================================================================
    totalEnergy=2982, maxAmplitude=8, low(114, 18), high(489, 47),
    totalEnergy=1648, maxAmplitude=6, low(124, 16), high(380, 14),
    totalEnergy=786, maxAmplitude=4, low(143, 8), high(223, 8),
    totalEnergy=1344, maxAmplitude=5, low(158, 7), high(328, 11),
    totalEnergy=949, maxAmplitude=5, low(158, 11), high(255, 11),

    Gdy wcisnę "2"
    =====================================================================
    totalEnergy=71423, maxAmplitude=34, low(482, 84), high(2407, 182),
    totalEnergy=61362, maxAmplitude=29, low(589, 52), high(2196, 153),
    totalEnergy=54813, maxAmplitude=29, low(731, 68), high(2057, 152),
    totalEnergy=48103, maxAmplitude=30, low(765, 66), high(1883, 148),
    totalEnergy=45178, maxAmplitude=27, low(742, 25), high(1839, 126),

    Gdy wcisnę "3"
    =====================================================================
    totalEnergy=147564, maxAmplitude=45, low(839, 81), high(3659, 157),
    totalEnergy=146385, maxAmplitude=46, low(837, 76), high(3641, 147),
    totalEnergy=146507, maxAmplitude=45, low(810, 46), high(3649, 134),
    totalEnergy=144311, maxAmplitude=45, low(798, 40), high(3620, 149),
    totalEnergy=144930, maxAmplitude=44, low(802, 44), high(3627, 150),
    totalEnergy=146732, maxAmplitude=44, low(804, 22), high(3653, 138),
    totalEnergy=145719, maxAmplitude=44, low(785, 27), high(3643, 139),

    Gdy wcisnę "4"
    =====================================================================
    totalEnergy=241023, maxAmplitude=58, low(984, 140), high(4842, 55),
    totalEnergy=241295, maxAmplitude=59, low(981, 150), high(4845, 32),
    totalEnergy=243241, maxAmplitude=58, low(974, 155), high(4867, 51),
    totalEnergy=244121, maxAmplitude=58, low(965, 164), high(4878, 40),
    totalEnergy=241596, maxAmplitude=59, low(985, 151), high(4848, 52),
    totalEnergy=241108, maxAmplitude=58, low(997, 151), high(4840, 36),


    Gdy wcisnę "5"
    =====================================================================
    totalEnergy=87655, maxAmplitude=37, low(920, 201), high(2548, 188),
    totalEnergy=88526, maxAmplitude=38, low(1000, 149), high(2543, 199),
    totalEnergy=86908, maxAmplitude=38, low(921, 203), high(2530, 184),
    totalEnergy=88619, maxAmplitude=38, low(992, 146), high(2549, 203),
    totalEnergy=87700, maxAmplitude=38, low(937, 191), high(2545, 191),
    totalEnergy=87417, maxAmplitude=36, low(989, 141), high(2529, 191),
    totalEnergy=85893, maxAmplitude=38, low(959, 183), high(2505, 188),


    Gdy wcisnę "6"
    =====================================================================
    totalEnergy=58355, maxAmplitude=32, low(923, 143), high(2131, 84),
    totalEnergy=58382, maxAmplitude=31, low(914, 146), high(2133, 84),
    totalEnergy=58572, maxAmplitude=32, low(899, 164), high(2145, 83),
    totalEnergy=56925, maxAmplitude=31, low(926, 147), high(2101, 96),
    totalEnergy=54076, maxAmplitude=30, low(920, 164), high(2034, 75),
    totalEnergy=51812, maxAmplitude=30, low(952, 137), high(1969, 77),
    totalEnergy=46964, maxAmplitude=29, low(949, 136), high(1853, 88),


    Gdy wcisnę "7"
    =====================================================================
    totalEnergy=271262, maxAmplitude=68, low(2059, 141), high(4825, 57),
    totalEnergy=270163, maxAmplitude=68, low(2058, 146), high(4815, 44),
    totalEnergy=267584, maxAmplitude=68, low(2062, 149), high(4786, 30),
    totalEnergy=267603, maxAmplitude=68, low(2084, 139), high(4776, 18),
    totalEnergy=266303, maxAmplitude=68, low(2086, 141), high(4761, 12),
    totalEnergy=266785, maxAmplitude=68, low(2075, 136), high(4772, 16),
    totalEnergy=268688, maxAmplitude=68, low(2065, 146), high(4795, 31),
    totalEnergy=270235, maxAmplitude=68, low(2066, 144), high(4812, 47),


    Gdy wcisnę "8"
    =====================================================================
    totalEnergy=116487, maxAmplitude=47, low(2240, 180), high(2351, 201),
    totalEnergy=115251, maxAmplitude=46, low(2229, 153), high(2329, 166),
    totalEnergy=112543, maxAmplitude=47, low(2190, 201), high(2307, 159),
    totalEnergy=117447, maxAmplitude=47, low(2272, 181), high(2344, 196),
    totalEnergy=114474, maxAmplitude=47, low(2211, 183), high(2332, 172),
    totalEnergy=114110, maxAmplitude=48, low(2214, 192), high(2312, 160),
    totalEnergy=118413, maxAmplitude=46, low(2286, 186), high(2353, 207),
    totalEnergy=113543, maxAmplitude=48, low(2216, 175), high(2298, 173),


    Gdy wcisnę "9"
    =====================================================================
    totalEnergy=217479, maxAmplitude=63, low(2359, 159), high(3901, 166),
    totalEnergy=218476, maxAmplitude=63, low(2347, 153), high(3922, 155),
    totalEnergy=215700, maxAmplitude=63, low(2321, 183), high(3899, 144),
    totalEnergy=215663, maxAmplitude=63, low(2322, 188), high(3895, 160),
    totalEnergy=218480, maxAmplitude=63, low(2341, 165), high(3926, 152),
    totalEnergy=218853, maxAmplitude=62, low(2323, 184), high(3943, 164),
    totalEnergy=216500, maxAmplitude=63, low(2299, 202), high(3923, 147),


    Gdy wcisnę "0"
    =====================================================================
    totalEnergy=180965, maxAmplitude=59, low(2761, 239), high(2927, 202),
    totalEnergy=178536, maxAmplitude=58, low(2746, 224), high(2885, 234),
    totalEnergy=182626, maxAmplitude=59, low(2791, 260), high(2920, 198),
    totalEnergy=181429, maxAmplitude=59, low(2762, 234), high(2921, 223),
    totalEnergy=180788, maxAmplitude=58, low(2719, 207), high(2949, 214),
    totalEnergy=181534, maxAmplitude=60, low(2715, 209), high(2970, 237),






    Co powinienem sobie jeszcze policzyć i co do czego porównać żeby filtrować tylko DTMFy i pozbyć się reagowania algorytmu na moje piszczenie, burczenie, muzykę czy też zwykłe odliczanie ?
  • #7 11545885
    lukashb
    Poziom 39  
    Sam jestem ciekaw. Zakończyłem co prawda próby, bo mi algo nie ruszył poprawnie. Ale także się obawiałem już wtedy tych właśnie problemów o których piszesz. Na moje oko, tu może chodzi o nie przesterowanie wejścia, zniekształcenia tego co podajesz itp. Do Goertzela najlepiej sinus, a gwizd to wiadomo co może mieć tym bardziej muzyka. Ale nie powinno to źle działać, bo jakoś się to buduje i chodzi. Ale patrzę i w czym problem :D masz wynik magnitudy w postaci low i high i teraz tak, gwizd: low(13, 5), high(314, 21), , a teraz DTMF low(592, 28), high(4003, 61). Magnituda jest liczona dla współczynników przeliczonych dla konkretnej częstotliwości badanej, czyli dla numeru 1 masz ewidentnie dużo większe wskazania niż przy mowie itp. Dałbyś może zmienną poziom i powiedzmy sprawdzałbyś czy magnituda większa od poziomu jak tak to ton ok, jak nie to liczymy dalej. Poza tym, podzieliłbyś się kodem? dalej mnie temat interesuje i chcę dojść do poznania jak się powinno napisać. Na co kolega pisze? tzn jaka rodzina?
  • #8 11546474
    sigmalko
    Poziom 10  
    lukashb napisał:
    Poza tym, podzieliłbyś się kodem?




    Kod: Java
    Zaloguj się, aby zobaczyć kod





    Kod: Java
    Zaloguj się, aby zobaczyć kod






    Kod: Java
    Zaloguj się, aby zobaczyć kod






    Kod: Java
    Zaloguj się, aby zobaczyć kod






    Wycinek pozostałości


    Kod: Java
    Zaloguj się, aby zobaczyć kod

    [/code]

    Dodano po 1 [minuty]:

    Ale nie bardzo podoba mi się koncepcja w której trzeba zaszyć w kodzie jakieś współczynniki w rodzaju "jeżeli coś tam jest większe niż 15dB to uznaj, że ok". Fajnie rozwiązać problem bez stałych.

    Dodano po 1 [godziny] 1 [minuty]:

    Obecnie zawęziłęm krąg poprzez:

    [1]
    uwzględnienie forward twist / reverse twist. Czyli
    jeśli sumaNiskichCzestotliwosci>sumaWysokichCZestotliwosci to jesli
    Math.abs(10 * Math.log10(sumLow / sumHigh) >= 8
    czyli jesli jest > 8 dB to odrzucam sygnał

    jeśli sumaNiskichCzestotliwosci<sumaWysokichCZestotliwosci to jesli
    Math.abs(10 * Math.log10(sumLow / sumHigh) >= 4
    czyli jesli jest > 4 dB to odrzucam sygnał

    [2]
    Siła sygnału będąca sumą najwyżej energii z LOW i najwyższej z HIGH ma być większa niż 75% sumy wszystkich 8 energii.


    Mniej dźwięków wykrywa się jako DTMF. Nadal wykrywa się muzyka.

    Jeśli chodzi o niższe wartości low/high to wystarczy że podgłośnię muzykę i przebiję progi/wskaźniki zaszywane w kodzie.





    Gdzieś czytałem, że można liczyć jeszcze:
    - drugą harmoniczną dla każdej z 8 częstotliwości; no ale jeśli policzę, to co dalej. Pisali, że 1-sza harmoniczna ma być przynajmniej 8 razy większa niż 2-ga. Ale wtedy nie wykrywa mi np. dtmfa "8". Czytałem też że rozpatrywanie drugiej harmonicznej wcale nie jest wiarygodne.
  • #9 11546950
    sigmalko
    Poziom 10  
    lukashb napisał:
    Magnituda jest liczona dla współczynników przeliczonych dla konkretnej częstotliwości badanej





    Interesuje mnie jeszcze relacja "totalEnergy" energii całkowitej sygnału (suma kwadratów 205 próbek) do energii uzyskanych z algorytmu Goertzela.
    Gdy wciskam "4" tuż obok mikrofonu, to uzyskuję:
    - 36 904 całkowitej energii sygnału oraz
    - 1 164 to suma energii z Goertzela dla 4 częstotliwości LOW
    - 1 784 to suma energii z Goertzela dla 4 częstotliwości HIGH
    - LOW + HIGH = 1164+1784=2 536

    Pytanie: gdzie się podziała reszta energii, tj. 36 904 - 2 536 = 34 368 ? Energia z palety częstotliwości dtmf to tylko 6,9 % całkowitej energii sygnału.

    Czy:
    - źle policzyłem energię?
    - pozostała energia to szum albo inne dźwięki (raczej nic wtedy nie mówiłem i w tle też niczego nie było słychać)?
    - pozostała energia sobie wyciekła po wszystkich częstotliwościach z zakresu 0 - 8000 Hz, bo np. nie użyłem sobie funkcji okna?
REKLAMA