🥇Přijímání celého Bluetooth najednou na SDR s CUDA? Snadné

V poslední době se mi moji kolegové z "oboru" nezávisle na sobě začali ptát: jak získat z jednoho SDR přijímače současně všechny kanály Bluetooth? Šířka pásma to umožňuje, existují SDR s výstupní šířkou 80 MHz a více. Samozřejmě to lze udělat na FPGA, ale čas na takový vývoj bude značný. Už dlouho jsem věděl, že to lze na GPU udělat docela snadno, ale aby to bylo takhle!

Standard Bluetooth definuje fyzickou vrstvu ve dvou verzích: Classic a Low Energy. Specifikace je zde. Dokument je ohromně velký, číst ho celý je pro mozek nebezpečné. Naštěstí velké společnosti, které vyrábějí měřicí techniku, mají prostředky na vytváření přehledných dokumentů na toto téma. Tektronix a National Instruments, například. Vůbec nemám šanci soutěžit s nimi co do kvality prezentace materiálu. Zájemce prosím, aby se studiem pověnovali podle odkazů.

Vše, co potřebuji vědět o fyzické vrstvě pro vytvoření víc kanálového filtru, je krok frekvenční mřížky a rychlost modulace. Tyto informace jsou shrnuty v tabulce v jednom z uvedených dokumentů:

Takže potřebujeme rozdělit šířku pásma 80 MHz na 79 filtrů se krokem nastavení 1 MHz a současně na 40 filtrů se krokem nastavení 2 MHz. Frekvence sledování výstupů filtrů by měly být 1 MHz a 2 MHz, příslušně.

Potřebujeme tedy dvě hřebínkové struktury filtrů.

Nejprve vybereme parametry těchto filtrů na základě šířek signálů Bluetooth Classic a Bluetooth Low Energy. Potřebujeme jejich impulzní charakteristiky, abychom spočítali zatížení výpočetního zařízení filtru. Je třeba podotknout, že délky impulzních charakteristik jsme vybrali na základě požadavků "rychlého" filtračního algoritmu. Podstata se od toho nemění. A počet koeficientů impulzní charakteristiky by neměl být příliš velký, aby byl filtr realizovatelný na rozumném výpočetním zařízení.

Pro filtry s krokem 1 MHz zvolíme pásmo propustnosti NF (polovina pásma propustnosti pásmového filtru) 500 kHz, délku impulzní charakteristiky přizpůsobíme na 480 odvodů. Pro filtry s krokem 2 MHz zvolíme tyto parametry 1 MHz a 240 odvodů, příslušně. Typ okna volíme Kaiser. Vypočítáme impulzní charakteristiky v programu filterDesigner a exportujeme je ve formátu C-header:

Snímky z filterDesigner

Úlohu lze řešit přímo: vytvořit odpovídající množství filtrů DDC (Digitální dolní převodník). Tento přístup je dobrý pro FPGA, kde je možná úspora díky snížení bitové šířky výpočetních bloků prvních stupňů. Také FPGA je nejvíce energeticky efektivní způsob realizace. Avšak náklady na práci jsou v tomto případě nejvyšší.

Při aplikaci filtrů formou hrebenatky na dnes populárních GPU se objevuje možnost realizace složitějšího algoritmu: polifázové hrebenatky filtrů na základě FFT, která je na CUDA dostupná z knihovny. V zahraniční literatuře se algoritmus nazývá Polyphase nebo WOLA (Weight, Overlap and Add) FFT Filterbank. Lenost kreslit mi nedovoluje provést názorné vysvětlení vlastními silami. Na internetu je mnoho materiálů na toto téma, zejména názorný graf je vyhotoven zde na stránce 11 (velké díky váženým autorům), tady je:

Pokud dovolíte, pokusím se schéma zpracování vysvětlit svými slovy. Slabší povahy prosím, ať nečtou.

Pokud dovolíte, pokusím se schéma zpracování vysvětlit svými slovy v rámci svých metodických možností. FFT je konvoluce vstupního signálu s celým spektrem komplexních ortogonálních harmonických, které se vejdu do intervalu impulsní charakteristiky. Impulsní charakteristika filtru, na kterou je signál násoben před vstupem FFT, je modulována tímto spektrem harmonických. Jinými slovy, obálka impulsních charakteristik výsledných filtrů hrebenatky je vyňata z závorky. Dále, spektrum harmonických je zredukováno na určité množství, vzhledem k rozšíření pásma propusti filtru ve srovnání s filtrem s obdélníkovým oknem. Na obrázku vidíme redukci na čtyři. Jinými slovy, po rozšíření pásma oknem Kaiser (s současným zvýšením útlumu v pásmu zpoždění) již potřebujeme pouze čtvrtinu filtrů. Ostatní jsou nadbytečné, jejich frekvenční charakteristiky se překrývají. Z čtyř bodů FFT jdoucích za sebou vybíráme pouze nulový, jehož výpočtem je sečtení čtyř
vstupních bodů, pořízených v čase, rovné čtvrtině délky původního FFT.

Vybereme hardware, který máme po ruce. Je to vstupní deska od společnosti „Instrumentální Systémy“ FMC126P. O ní jsem již psal v jednom z předchozích článků. článku. Do slotu FMC desky byl vložen submodul téže společnosti s převodníkem AD9371 s šířkou pásma 100 MHz. Celý tok z převodníku může být neustále přenášen do počítače pro zpracování.

Vybereme grafickou kartu s GPU GTX 1050. (Lhal jsem, to ona si nás vybrala: bylo to jediné, co bylo po ruce, byla vytržena z počítače pro výpočty antén, ale o to více bylo překvapující vidět fungující hřebínek). Přejdeme k softwarové části.

Bohužel, kvůli licencím nemůžeme zveřejnit celý kód. Můžeme ukázat pouze jádra GPU. Nicméně, zbytek kódu není tak zajímavý.

Toto je jádro, které provádí násobení signálu oknem a sčítání, a obal pro jeho volání:

__global__ void cuComplexMultiplyWindowKernel(const cuComplex *data, const float *window, size_t windowSize, cuComplex *result) {
    __shared__ cuComplex multiplicationResult[480];
    multiplicationResult[threadIdx.x] = cuComplexMultiplyFloat(data[threadIdx.x + windowSize / 4 * blockIdx.x], window[threadIdx.x]);
    __syncthreads();
    cuComplex sum;
    sum.x = sum.y = 0;
    if (threadIdx.x < windowSize / 4) {
        for(int i = 0; i < 4; i++) {
            sum = cuComplexAdd(sum, multiplicationResult[threadIdx.x + i * windowSize / 4]);
        }
        result[threadIdx.x + windowSize / 4 * blockIdx.x] = sum;
    }
}

cudaError_t cuComplexMultiplyWindow(const cuComplex *data, const float *window, size_t windowSize, cuComplex *result, size_t dataSize, cudaStream_t stream) {
    size_t windowStep = windowSize / 4;
    cuComplexMultiplyWindowKernel<<>>(data, window, windowSize, result);
    return cudaGetLastError();
}

Kód zpracování signálu, který vyvolává toto jádro, přesně opakuje schéma algoritmu zobrazeného na obrázku výše, proto nevidím důvod ho zde uvádět.

Hřebínky byly ověřeny na výstupním spektru kanálů v reálném čase. Na vstup AD9371 byl podáván signál generátoru signálů 2450 MHz, selektivita filtrů odpovídala výpočtům.