STM32F4 - konstrukce banky filtrů

Zdravím,

tvořím hlukoměr pomocí mikrokontroléru STM32F407VGT6. Vzorky vstupního signálu z ADC převodníku si ukládám do bufferu. Tyto vzorky musím filtrovat třetino-oktávovou bankou filtrů. Tzn. vzorky musí projít 10-ti filtry. Těchto 10 výsledných filtrovaných vzorků se musí ještě usměrnit a poté z nich vypočítat podle časově váženého filtru efektivní hodnotu signálu. Z té poté vypočítat úrověň hluku v dB.

Problém je, že oněch 10 filtrů je moc výpočetně náročných (využil jsem IIR filtry). A během načítání vzorků z ADC nestíhá mikrokontroler filtraci. Tzn. není to real-time (což bych potřeboval…). Zkoušel jsem využít i DSP knihovnu, která přinesla mírné zlepšení, ale pořád to není dostačující.

Mikrokontroler jede na frekvenci 168 MHz a má Hardwarovou FPU, což mi příjde jako solidní výkon, přesto nestíhá. Mám pocit, že na to jdu nesprávně. Hlukoměry, které se dají koupit, toto podporují, a myslím si, že v nich nemůže být tak výkonný procesor/mikrokontroler.

Jinak počítám s floatovými čísly. Vím že počítat s celýmy čísly by bylo rychlejší, ale díky kvantizaci koeficientů na celá čísla se IIR filtry dost často stávají nestabilními.

Nemá někdo v tomto nějaké zkušenosti?

Jakou máš vzorkovací frekvenci ? Možná by pomohlo počítat to v intech a místo IIR použít FIR. Mohl bys zavést decimaci a nepočítat všechny vzorky pro každý filtr ale jen podle jeho šířky pásma.

Nejsem si jist, zda by u tohoto mcu výrazně pomohla změna na FIR (v int). FIR potřebuje vyšší řád než IIR, ale zas méně MAC na stupeň. MAC(‘MLA’) trvá 2 takty, FMA 3. Na 16b výpočtu ovšem lze použít 1-cyklovou SMUAD/SMLALD (2xMAC). Filtr by však možná bylo nutné rozdělit na více částí a meřítkovat mezivýsledky.

K Radiusovým otázkám bych přidal:
jaký překladač (umí HW float)?
má překladač informaci o použitém mcu (aby mohl použít hw float)?

Pro překladač gcc “GNU Tools for ARM Embedded Processors” jsou parametry dle readme.txt:
-mthumb -mcpu=cortex-m4 -mfloat-abi=hard -mfpu=fpv4-sp-d16,
knihovny ze složky “…/armv7e-m/fpu” (jsou ve složce překladače 2).

Radius: vzorkovací frekvenci mam 44100 Hz. S tou decimací to zní jako dobrý nápad, to vyzkouším. Jinak jsem to zkoušel počítat v intech a bylo to pomalejší než s floatem, což mě překvapilo, ta FPU je docela rychlá…

piityy: mám překladač gcc “GNU Tools for ARM Embedded Processors” jak píšeš (používám coocox IDE).

překladač mi vypisuje toto:
[cc] arm-none-eabi-gcc -mcpu=cortex-m4 -mfpu=fpv4-sp-d16 -mfloat-abi=hard -mthumb -Wall -g -nostartfiles -Wl,-Map=SOUND_LEVEL_METER.map -O0 -Wl,–gc-sections -LC:\CooCox\CoIDE\configuration…

knihovny naimportované mám.

Parametry se zdají být vpořádku až na jeden: “-O0” - máš vypnuté optimalizace. Zkus nejdřív -O2, kdyby to nestačilo, tak -O3.

Tak optimalizace značně pomohla. Zrychlila výpočet 4 násobně :slight_smile:. Děkuji za radu. Mimochodem ještě dotaz, má ta optimalizace i nějaké negativní účinky ?

Blbě se ladí. Překladač kód různě přerovná, nepotřebné části vyhází(např. zpožďovací smyčka je z pohledu překladače zbytečná pokud mu není sdělen opak) a podobně. U O3 se navíc využívá nedefinované chování, takže programátor začátečník se pak může divit že co fungovalo v debug bez optimalizací nefunguje v release s nimi. Např. taková jednoduchá věc jako je přetečení int (celé číslo se znaménkem) není definováno a překladač s tím může naložit dle libosti. Různé verze překladače se v takové situaci mohou chovat různě. Teoreticky může dokonce naložit dle libosti s celým programem, kde takové chování je, ale v praxi se tak naštěstí neděje. Dalším příkladem budiž kód: if(x) if(y) ... else ... Ke kterému if patří ono else? Toto není definováno -> závorkovat

Ohledně filtru - máš nějakou optimalizovanou formu nebo jen otrocký přepis rovnice?

To rozdělení signálu několika filtry - není to z důvodu různé citlivosti sluchu na různé frekvence? Nebylo by v takovém případě lepší počítat FFT a tu pak váhovat onou citlivostí?

Využívam knihovny “CMSIS DSP Software Library” disca.upv.es/aperles/arm_cortex_m3/curset/CMSIS/Documentation/DSP/html/

Ty DSP filtrační funkce by měli být optimalizovány. Když jsem porovnával jejich rychlost a rychlost “vlastního kódu” (tj ten přepis rovnice) tak byly zhruba o 30% rychlejší…

Ano, přesně, je to kvůli citlivosti ucha. A ano FFT by určitě taky šla. Dělám tento projekt jako diplomovou práci a můj vedoucí chce abych pracoval v časové oblasti. Protože se to používá u profi hlukoměrů… A jseš si jistý, že by FFT byla rychlejší ?

Těžko říct, neznáme řády filtrů a moje praktické zkušenosti ohledně výpočetního výkonu se týkají lehce odlišné oblasti. Ovšem za zkoušku by mi to stálo.
Nevýhoda oproti časové oblasti je, že FFT spočítá výkon v okně, nikoli vzorek za vzorkem.
Z šířky pásma nejužšího filtru lze odvodit potřebnou velikost okna.
Navíc se nemusí počítat nová FFT s každým novým vzorkem, ale třeba až po příchodu např. 8 nových vzorků. Určitě to máte namodelováno něčím jako matlab, neměl by tedy být problém si tam efekt takovéto operace vyzkoušet.
A kdyby byly výsledky přijatelné, zdůvodnit postup nedostatečným výpočetním výkonem :wink:

Pokud jde jen o úpravu kmitočtové charakteristiky, máš několik možností k úvaze:

FFT, FIR, případně vhodně řešené filtry IIR.

FFT: Zjistit si jeho potřebnou složitost. Tzn, jaké potřebuješ kmitočtové rozlišení. Spolu se vzorkovacím kmitočtem ti z toho vyjde potřebná délka. Nezapomínat i na to, že je třeba následně provést inverzní FFT nebo počítat výkon jednotlivých složek zvlášť. Nezapomenout i na potřebný překryv oken asi 25% délky FFT, jinak budeš mít v signálu mezery. (podle mě toto není vhodná metoda, bude velmi výpočetně náročná a výsledek nejistý).

FIR: Bude potřeba poměrně dlouhý. Čím vyšší vzorkovací kmitočet a větší potřebné rozlišení, tím bude delší a horší přesnost výpočtu.
Při užití 16bitových dat a koeficientů s 32bitovým akumulátorem je však výpočet s DSP instrukcemi velmi rychlý, pod 2 cykly na 1 odbočku FIRu (reálně asi 1.8 clocku na odbočku). Problém je, že na audio (v měřícím přístroji) se s 16b FIRem chodit nedá. Takže toto bych také vynechal, neboť FIRy s vyšším rozlišením koeficientů (32b) a akumulátoru (64b) budou dost pomalé. (jelikož cortex M4 na to nemá vhodné SIMD instrukce) FIR se používá tehdy, je-li zapotřebí lineární fáze. (digitální zpracování v komunikacích, SDR rádio). FIR filtr se dá udělat s konstantní fází, nebo lineárním sklonem fáze.

Soustava vhodných IIR filtrů je podle mě (jediný?) směr, kam se vydat. Vzhledem k tomu, že hlukoměru je docela buřt nějaká fázová nelinearita (zpoždění signálu nemá vliv na měřený akustický výkon), tak je to skvělý způsob řešení. Potřebnou přenosovou charakteristiku si namodeluj v komplexní s rovině (laplaceova transformace) a výsledné koeficienty použij do IIR filtru. Osobně doporučuji využít filtrů typu biquad, jelikož jsou jednoduché výpočetně a jsou přímým ekvivalentem filtrů analogových s operačními zesilovači. (a jsou stabilnější než dlouhý IIR) jsi-li tedy schopný si namodelovat vhodnou přenosovou charku užitím několika analogových filtrů jako lowpass, highpaass, peak, notch, shelf,… každý z těchto typů filtru lze realizovat filtrem typu biquad. Dejme tomu, že budeš nakonec potřebovat nějakých 10 biquadů typu peak/notch (10pásmový ekvalizér). To je ve výsledku akorát 10 * 5 MAC instrukcí a 10*4 přesuny dat, pokud dobře počítám. Ten procesor se u toho na 168MHz bude vyloženě nudit.

Avšak i u biquad filtrů je háček. Jelikož jsou to IIR filtry, zpětnovazební, nastávají problémy se stabilitou v případě, že se číselné velikosti pólů blíží jednotkové kružnici (klesá přesnost výpočtu, roste nestabilita). To se děje právě při nízkých kmitočtech, pak je třeba užít jiného typu filtru, např state variable filter. Opět poměrně jednoduchá záležitost a přímý ekvivalent analogového protějšku. Výhodou je to, že potřebuje ke své činnosti pouze 2 koeficienty F a Q, které lze velmi snadno vypočítat za běhu programu. (to jde i pro biquad filtry). Ve výpočtech koeficientů se nepoužívá nic hnusnějšího, než jeden sinus nebo tangens. Obě dvě tyhle věci jsou na ARMu s FPU rychlé i při použití knihovních funkcí z math.h (tabulky a aproximace).

Pro IIR (biquad) filtry a audio samozřejmě 32 bitů výpočty, float velkou výhodou (nebudou problémy s přetejkáním). O biquad filtrech se toho dá najít spousta na webu. Zkus se na to podívat.

Pro začátek by se ti mohlo hodit tady tohle: earlevel.com/

Omlouvám se za dlouhý příspěvek, snad ti k něčemu bude.

PS: Nejsem odborník na DSP, jen se o to trochu zajímám, zatím jen začátečnicky…

PPS: Jen mě děsí, že nahoře píšeš “použil jsem 10 IIR filtrů a ono to nestíhá”. Cos tam proboha vyrobil? Jaký řád IIR? jakým způsobem počítaný? 168MHz CPU, 48kHz vzorkování… to máš 3500 dostupných klokanů na vzorek. Tzn při 10 filtrech to znamená 350 cyklů na filtr na vzorek. I kdybych to napsal jako prase, tak to 10 biquadů stihne :open_mouth:

Škoda, že už se asi nedozvíme, jak to dotyčný vyřešil. Nechť aspoň mají inspiraci ostatní zájemci o filtry.