Traitement de signaux via Gnuradio
Write-up du challenge Zodiaque de hackropole (FCSC édition 2021)
Le
Dans ce deuxième writeup consacré aux challenges hackropole, je vous propose de résoudre une épreuve d’analyse d'un signal radio en utilisant le logiciel Gnuradio. À noter que je n’ai bien sûr pas lu de solution pour ce challenge afin de ne pas être spoilé. Par ailleurs je n’en ai pas trouvé sur Internet lors de mes recherches (02/01/2026). Par ailleurs, je ne suis pas du tout un expert en traitement du signal, je me contente d'utiliser quelques notions acquises durant mon parcours scolaire. Si j'ai fait des erreurs, n'hésitez pas à m'en faire part.
Description du challenge
Vous devez démoduler une transmission suspecte pour investiguer un nouveau protocole d’échange de messages. Vous avez la note suivante jointe au signal : - Contenu: texte brut - Taux d'échantillonnage: 32kHz - Constellation: QPSK - Temps symbole: ?? - Encodage: ??
Notions RF
Les ondes électromagnétiques
Les ondes électromagnétiques sont de forme sinusoïdale. Il y a 4 manières de placer de l'information dans un signal sinusoïdal :
- la fréquence : nombre d'oscillation dans un certain temps
- la phase : "position" dans la sinusoïde
- l'amplitude : "hauteur" de la sinusoïde
- éventuellement la polarisation (mais là on pourra pas le visualiser qu'avec notre machine)
Un signal sinusoïdal s'écrit sous la forme $s(t) = A \cos(\omega t + \phi)$ avec : - $A$ : l'amplitude maximale. - $\omega$ : la pulsation en rad/s ($\omega = 2\pi f$). - $t$ : le temps. - $\phi$ : la phase à l'origine (l'angle au temps $t=0$).
Notation complexe
Lors de l'étude des signaux, il est très fréquent de représenter ceux-ci sous forme complexe. La raison est qu'ils permettent de passer d'une ligne unidimensionnelle à un plan bidimensionnel. Les nombres imaginaires permettent ainsi de représenter cette deuxième dimension. A noter qu'un signal radio n'est pas complexe dans la nature. On l'exprime sous une forme complexe pour faciliter les calculs sur celui-ci.
Pour utiliser cette notation, on utilise la formule d'Euler qui dit :
$$ e^{j\theta} = \cos(\theta) + j \sin(\theta) $$
Par exemple, si l'on prend un signal de la forme :
$$ s(t) = A \cos(\omega t + \phi) = A \cos(\theta) $$
Alors sa forme complexe est :
$$ \underline{s}(t) = A e^{j\theta} = A \cos(\theta) + j.A \sin(\theta) $$
avec $\theta = \omega t + \phi$. On peut noter les égalités ci-dessous : $$ \begin{aligned} \text{Re}(\underline{s}(t)) = \text{Re}(e^{j\theta}) = \cos(\theta) \newline \text{arg}(\underline{s}(t)) = \text{arg}(e^{j\theta}) = \theta = \omega t + \phi \end{aligned} $$
L'avantage majeur de cette notation se situe lors de la dérivation.
- En réel : dériver $\cos(\omega t)$ donne $-\omega \sin(\omega t)$ (changement de fonction).
- En complexe : dériver $e^{j\omega t}$ donne $j\omega \cdot e^{j\omega t}$.
Dériver revient simplement à multiplier par $j\omega$.
On peut donc :
- Exprimer un signal sous sa forme complexe
- Le dériver / L'intégrer
- Replacer le signal sous forme réel.
Modulations
La modulation peut être définie comme le processus par lequel un signal est transformé de sa forme originale en une forme adaptée au canal de transmission (source : Wikipedia).
L'objectif est de placer une donnée dans une onde électromagnétique. Cette donnée peut être analogique (un signal continu) ou numérique (un signal binaire, c'est à dire à deux états). Il existe plusieurs moyens pour y parvenir, en agissant notamment : - sur la fréquence du signal (FM pour les signaux analogiques, FSK pour les signaux binaires) - sur l'amplitude du signal (AM pour les signaux analogiques, ASK pour les signaux binaires) - sur la phase du signal (PM pour les signaux analogiques, PSK pour les signaux binaires)

Je vous invite à lire mon article sur la démodulation d'un signal FSK avec audacity disponible ici si vous voulez en savoir plus sur ce type de modulation.
Q-PSK
La vidéo ci-dessous présente (bien) le fonctionnement de la modulation Q-PSK :
Perturbations à la réception
Baseband
Idéalement, un émetteur envoie un signal modulé sur une porteuse. Le récepteur paramètre son matériel pour recevoir le signal à la bonne fréquence, avant de le délivrer au logiciel. De cette manière, le logiciel peut traiter le signal en base band, c'est-à-dire à une fréquence centrée sur 0.

Un des intérêts de l'utilisation d'un signal en ramené en baseband est que l'on traite (et on stocke) que la bande passante utile.
Fréquence d'offset
En réalité, le signal que reçoit le logiciel comporte un offset : l'offset de fréquence (Carrier Frequency Offset ou CFO). Il s'agit de la différence entre la fréquence de l'onde porteuse utilisée par l'émetteur et la fréquence de l'oscillateur local du récepteur. Il peut avoir plusieurs causes :
- Imprécision du matériel : Les quartz (oscillateurs) ne sont pas parfaits. Si un émetteur est à 2,4000 GHz et le récepteur à 2,4001 GHz, il y a un offset de 100 kHz.
- Variations thermiques : La fréquence d'un oscillateur dérive quand la température change.
- L'Effet Doppler : Si l'un des deux appareils est en mouvement (l'émetteur ou le récepteur ou les deux), la fréquence perçue se déplace vers le haut ou vers le bas. Dans ce cas l'offset peut être dynamique (il varie dans le temps).
Pour mieux comprendre la cause de cette fréquence d'offset et ses conséquences au niveau mathématique, je vous propose d'analyser le workflow de réception d'un signal via un récepteur radio. Idéalement, lorsqu'on demande au récepteur d'écouter à la fréquence $f_{rf}$ :
- Le matériel écoute le signal à la fréquence $f_{rf}$. Pour cela il utilise un filtre passe-bande;
- Il paramètre son oscillateur local à cette fréquence et génère un signal sinusoïdale (de fréquence $f_{lo} \approx f_{rf}$) ;
- Il ramène le signal reçu à 0 Hz. Pour cela :
- Il multiplie le signal généré par son oscillateur au signal reçu. En sachant que : $$\cos(A)\cos(B) = \frac{1}{2}[\cos(A-B) + \cos(A+B)]$$ On obtient : $$o(t) = \frac{1}{2} \underbrace{\cos(2\pi (f_{rf} - f_{lo})t + \phi)}{\text{Fréquence Différence}} + \frac{1}{2} \underbrace{\cos(2\pi (f{rf} + f_{lo})t + \phi)}{\text{Fréquence Somme}}$$ Or comme $f{rf} \approx f_{lo}$ : $$o(t) = \frac{1}{2} \underbrace{\cos(\phi)}{\text{Fréquence Différence}} + \frac{1}{2} \underbrace{\cos(4\pi f{rf} + \phi)}_{\text{Fréquence Somme}}$$
- Il fitre avec un passe bas pour ne conserver que le terme issue de la difference des fréquence.
- Il numérise les composantes I et Q ;
- Il envoie ces échantillons via USB à GNU Radio.

En réalité, $f_{lo} \neq f_{rf}$, il existe donc un offset de fréquence, issue de la soustraction $\Delta f = f_{rf}-f_{lo}$. Si on en tient compte dans l'équation obtenue précédemment après filtrage des hautes fréquences, on obtient l'expression suivante :
$$ o(t) = \frac{1}{2} \underbrace{\cos(2\pi (f_{rf} - f_{lo})t + \phi)}_{\text{Fréquence Différence}} = \frac{1}{2} \cos(2\pi\Delta ft + \phi) $$
Désynchronisation temporelle
Pour décoder l'information, le récepteur doit échantillonner le signal au moment optimal, généralement au centre exact des symboles, là où le signal est le plus stable. Si le récepteur est décalé, il peut échantillonner le signal sur les phases de transition entre les symboles.

Le challenge
Le challenge nous propose d'étudier un signal, en nous fournissant quelques informations pour orienter notre étude :
- sa fréquence d'échantillonnage
- son type de modulation (Q-PSK)
Le signal est fourni dans un fichier au format IQ. Il s'agit d'un fichier contenant une suite de nombre complexe. Chaque couple de 2 (ou 4) octets code respectivement la partie réelle et la partie imaginaire d'un échantillon du signal. Le I signifie In-Phase et Q Quadrature.
Analyse du signal
Comme le challenge mentionne l'usage d'une modulation PSK, j'ai commencé par m'intéresser à l'allure de la phase du signal (PSK est une modulation de phase). J'ai ainsi affiché l'évolution de la phase dans le temps, et sa constellation.

La constellation du signal a une forme circulaire, ce qui diffère de celle du Q-PSK attendu qui produit une forme carrée, comme le montre la figure ci-dessous :

Comme on sait qu'il s'agit d'une modulation Q-PSK, cette constellation traduit deux problèmes :
- la présence d'une fréquence d'offset : les points forment un cercle ;
- une désynchronisation temporelle : les points ne forment ni un cercle parfait, ni un carré.
En outre, l'allure de la phase dans le temps présente une caractéristique singulaire :

On peut clairement distinguer deux parties (une en rouge et une en orange). La partie en rouge ressemble à du bruit, ce qui laisse penser que le signal utile est le signal orange. C'est d'autant plus flagrant lorsque l'on zoome sur la phase de transition entre les deux signaux.

On peut tracer 4 lignes parallèles sur le signal.

La phase croit linéairement dans le temps modulo $pi$. En effet, dès qu'une ligne atteint la phase maximale, c'est-à-dire $\pi$, elle « repart » de $-\pi$. On parle de courbe en dent de scie.
Suppression de la fréquence d'offset
Une première idée un peu naïve serait de retirer un signal à dent de scie semblable à celui observé (même amplitude, même fréquence) dans l'objectif de rendre les états "parallèles" à l'axe des abscisses (éviter qu'il croie linéairement). J'ai choisi d'utiliser un signal semblable à celui décrit en rouge sur cette capture :

Pour cela je réalise le circuit suivant :
- Je récupère uniquement la portion utile du signal d'étude en retirant les premiers échantillons via le bloc
Skip head. J'ai identifié en tâtonnant que les 370 premiers points ne faisaient pas partie du signal utile. - Je génère un signal cosinus dont la phase vaut $\text{samp_rate}/1000$. Je calcule sa phase, ce qui me donne un signal en dent de scie (La phase d'un signal sinusoïdale est une dent de scie (linéaire et modulo $2 \pi$)). En ce qui concerne la valeur choisie pour la fréquence ($\text{samp_rate}/1000$), je l'ai déterminé en procédant par tâtonnement : j'ai superposé la phase du signal d'origine avec la phase du signal sinusoïdale que je génère (ma dent de scie), puis j'ai modifié la fréquence de ce dernier jusqu'à ce que les deux courbes se superposent à peu près.



- Je soustrais ma dent de scie à la phase du signal d'étude. On obtient alors les courbes ci-dessous :
- en bleu, le signal d'étude
- en vert, la dent de scie générée (phase d'un signal sinusoïdale de fréquence $\text{samp_rate}/1000$)
- en rouge, le résultat de la soustraction.

On observe qu'il y a encore un effet de dent de scie, et que le signal oscille cette fois entre $6 = 2 \pi$ et $-6 = -2 \pi$. C'est parce que notre résultat doit se calculer modulo $\pi$

Il n'existe pas, à ma connaissance, de bloc permettant de calculer le modulo d'un signal sur Gnuradio. Mais la puissance de ce logiciel réside dans le fait qu'on peut développer nos propres blocs très facilement. Pour créer un bloc modulo, il suffit d'ajouter un bloc Python block, puis de cliquer sur le bouton édition de code de celui-ci, et fin de rentrer le code ci-dessous développé à cette occasion.
"""
Embedded Python Blocks:
Each time this file is saved, GRC will instantiate the first class it finds
to get ports and parameters of your block. The arguments to __init__ will
be the parameters. All of them are required to have default values!
"""
import numpy as np
from gnuradio import gr
class blk(gr.sync_block): # other base classes are basic_block, decim_block, interp_block
"""Embedded Python Block example - a simple multiply const"""
def __init__(self): # only default arguments here
"""arguments to this function show up as parameters in GRC"""
gr.sync_block.__init__(
self,
name='Custom block - modulo', # will show up in GRC
in_sig=[np.float32],
out_sig=[np.float32]
)
def work(self, input_items, output_items):
# Applique le modulo n (ex: 2.0) sur tout le vecteur d'entrée
output_items[0][:] = np.mod(input_items[0], np.pi)
return len(output_items[0])
Le schéma bloc ressemble alors à ça :

Il nous permet d'obtenir le signal suivant :

On distingue clairement 4 états :

Synchronisation temporelle
En revanche, certains points sont incohérents, représentés en rouge sur la capture ci-dessous.

Ces points sont très certainement dus au saut induit par le modulo et ne représente pas des symboles, mais la transition entre des symboles.

Il nous faut maintenant associer chacun des points à un symbole (1,2,3 ou 4). Pour cela, je crée un nouveau bloc python sur mesure. Celui-ci associe chaque point à un état en utilisant la règle suivante : Soit $x$ la valeur du point et $y$ la sortie du bloc associé au point $x$, c'est-à-dire $y = f(x)$ :
- Si $x < -\pi/2$, alors $x = 0$;
- Si $-\pi/2 < x < 0$, alors $x = 1$
- Si $0 < x < \pi/2$, alors $x = 2$
- Si $\pi/2 < x$, alors $x=3$
Le code Python on obtient :
"""
Embedded Python Blocks:
Each time this file is saved, GRC will instantiate the first class it finds
to get ports and parameters of your block. The arguments to __init__ will
be the parameters. All of them are required to have default values!
"""
import numpy as np
from gnuradio import gr
class blk(gr.sync_block): # other base classes are basic_block, decim_block, interp_block
"""Embedded Python Block example - a simple multiply const"""
def __init__(self): # only default arguments here
"""arguments to this function show up as parameters in GRC"""
gr.sync_block.__init__(
self,
name='Decode', # will show up in GRC
in_sig=[np.float32],
out_sig=[np.float32]
)
# if an attribute with the same name as a parameter is found,
# a callback is registered (properties work, too).
def work(self, input_items, output_items):
in0 = input_items[0] # Le tableau de samples en entrée
out = output_items[0] # Le tableau de samples en sortie
# On définit nos seuils (bins)
# Les seuils créent les intervalles : (-inf, -pi/2), [-pi/2, 0), [0, pi/2), [pi/2, +inf)
bins = np.array([-np.pi/2, 0, np.pi/2])
# np.digitize renvoie l'indice du "bac" (bin) pour chaque échantillon
# Si x < -pi/2 -> 0
# Si -pi/2 <= x < 0 -> 1
# etc.
out[:] = np.digitize(in0, bins)
return len(out)
En sortie :

Il faut maintenant supprimer les points inutiles et uniformiser le signal. Après avoir compté le nombre d'échantillons composant les premiers symboles, j'ai trouvé qu'un symbole est codé sur un nombre d'échantillons variant entre 6 et 9.
J'ai créé un nouveau bloc Python qui exécute la logique suivante :
- Il itère sur chaque échantillon :
- Il compte le nombre d'échantillons successifs ayant la même valeur.
- Si $n=6$ échantillons successifs sont égaux, alors écrit la valeur de l'échantillon une fois.
- Sinon il ne fait rien.
- Il compte le nombre d'échantillons successifs ayant la même valeur.
Le code du bloc :
import numpy as np
from gnuradio import gr
class burst_filter(gr.basic_block):
def __init__(self, threshold=5):
gr.basic_block.__init__(
self,
name='BCustom block - Pack n points',
in_sig=[np.float32],
out_sig=[np.float32]
)
self.threshold = threshold
self.current_val = None
self.count = 0
self.buffer = []
def general_work(self, input_items, output_items):
in0 = input_items[0]
out = output_items[0]
out_idx = 0
for val in in0:
if val != self.current_val:
# Changement de valeur : on vérifie si la séquence précédente était trop courte
# Si count < threshold, le contenu de self.buffer est perdu
self.current_val = val
self.count = 1
self.buffer = [val]
else:
self.count += 1
# Si on vient d'atteindre le seuil, on libère tout le buffer accumulé
if self.count == self.threshold:
for b_val in self.buffer:
if out_idx < len(out):
out[out_idx] = b_val
out_idx += 1
# On vide le buffer car maintenant on va écrire en direct
self.buffer = []
# Si on est déjà au-dessus du seuil, on écrit directement
elif self.count > self.threshold:
self.current_val = val
self.count = 1
self.buffer = [val]
else:
# On stocke en attendant de savoir si c'est stable
self.buffer.append(val)
# On informe GNU Radio de combien d'échantillons on a consommé et produit
self.consume(0, len(in0))
return out_idx
J'obtiens ainsi une courbe propre, dans laquelle 6 échantillons codent un seul symbole.

Réechantillonage
J'utilise le bloc Keep 1 in N pour ne conserver qu'un seul échantillon sur les 6.

Enfin, après avoir converti les valeurs en char pour les coder sur un seul octet, je place le résultat dans un fichier binaire.

Interprétation mathématique du signal
Essayons d'exprimer l'équation d'un signal à dents de scie. Sur une période, une dent de scie a pour équation une fonction affine. On a donc :
$$ \begin{aligned} \text{Rampe}(t) = A.t + b \newline \text{Or } \text{Rampe}(0) = - \pi \Rightarrow b = - \pi \newline \text{Et } \text{Rampe}(T) = \pi \Rightarrow A.T - \pi = \pi \newline \iff A = \frac{2.\pi}{T} \newline \text{Donc } \text{Rampe}(t) = \frac{2.\pi}{T}.t - \pi \end{aligned} $$
La rampe est modulo $2.\pi$ (entre $]-\pi,\pi]$), donc
$$ \text{Rampe}(t) = \frac{2.\pi}{T}.t - \pi \mod{2\pi} $$
A cette rampe vient s'additionner les données codées sur 4 symboles, le tout modulo $2\pi$ (entre $]-\pi,\pi]$). La phase s'exprime donc
$$ \begin{aligned} \text{Phase}(t) = (\text{Rampe}(t) + \text{Données}(i) \mod{2\pi}) - \pi \newline \text{Phase}(t) = [(\frac{2.\pi}{T}.t - \pi \mod{2\pi}) + {Données}(i)] \mod{2\pi} - \pi \end{aligned} $$
On peut simplifier cette expression, car $((a + d \mod{b}) + c)\mod{b} = (a + d + c) \mod{b}$. On obtient donc, après simplification :
$$
\text{Phase}(t) = [\frac{2.\pi}{T}.t + {Données}(i) - 2.\pi] \mod{2\pi}
$$
Comme $f = 1/T$, on a :
$$ \text{Phase}(t) = [2.\pi.f.t + {Données}(i) - 2.\pi] \mod{2\pi} \text{où} {Données}(i) = i*\pi/2 $$
D'après Wikipedia, l'équation du Q-PSK est :
$$ s_{i}(t)={\sqrt {\frac {2E_{s}}{T_{s}}}}\cos \left(2\pi f_{c}t+(2i-1){\frac {\pi }{4}}\right),i=1,2,3,4 $$
Cette expression est cohérente à celle que l'on a déterminée. On a trouvé expérimentalement $f_{c} = \frac{32_{kHz}}{1000}$ (fréquence d'échantillonnage / 1000).
⚠️ En réalité, la fréquence $f_{c} = \frac{32_{kHz}}{1000}$ trouvé ici empiriquement n'est pas réellement la fréquence de la porteuse, mais la fréquence d'offset $\Delta f$, c'est-à-dire la différence entre la fréquence réelle de la porteuse (celle utilisée idéalement lors de l'émission) et celle utilisée au moment de la mesure.
Pour résumer, le schéma Gnuradio :
- Ne garde que la portion utile du signal en éliminant les 370 premiers échantillons (
Skip head). On a alors $s_{i}(t)={\sqrt {\frac {2E_{s}}{T_{s}}}}\cos \left(2\pi f_{c}t+(2i-1){\frac {\pi }{4}}\right),i=1,2,3,4$ - Calcule la phase du signal (
Complex to arg). On obtient alors $2\pi f_{c}t+(2i-1){\frac {\pi }{4}}$ - Soustrait le terme $2\pi f_{c}t$ (
Signal sourceetSubstract). On obtient alors $(2i-1){\frac {\pi }{4}}$, mais suréchantillonné (on a des échantillons sur les phases de transition). - Supprime les echantillons en les plaçant à l'état le plus proche avec la fonction $y = f(x)$ tel que
- Si $x < -\pi/2$, alors $x = 0$;
- Si $-\pi/2 < x < 0$, alors $x = 1$
- Si $0 < x < \pi/2$, alors $x = 2$
- Si $\pi/2 < x$, alors $x=3$
- Rééchantillone le signal en ne conservant qu'un échantillon si les 5 échantillons suivants sont identiques.
Décodage du signal
On dispose maintenant un fichier contenant une liste de chiffre entre 0 et 3 compris, caractérisant les valeurs des symboles extraites du signal.

En Q-PSK, chaque symbole représente une combinaison de 2 bits : 00, 01, 10 ou 11. Dans notre cas, bien que nous disposions des symboles, le mapping entre ces symboles et les duos de bits nous est parfaitement inconnu. Afin de déterminer cette association, j'ai développé le script ci-dessous, qui teste le décodage du signal en utilisant l'ensemble des permutations possibles de la liste ['00', '01', '10', '11'].
En partant du postulat que le signal original contient du texte, le script devrait afficher au moins un message intelligible après décodage.
from itertools import permutations
symboles = ['00', '01', '10', '11']
with open("gnu-out.bin","rb") as file:
data = file.read()
for perm_symbole in permutations(symboles):
decoded_data = [perm_symbole[i] for i in data]
decoded_data = "".join(decoded_data)
decoded_data = [decoded_data[i:i+8] for i in range(0, len(decoded_data), 8)]
decoded_ascii = ""
for char in decoded_data:
decoded_ascii += chr(int(char,2))
print(decoded_ascii)
print(perm_symbole)

Et c'est un succès 🥳. Parmi toutes les données affichées, un message intelligible se distingue. On peut noter l'association duo de bits / symboles : (0 -> '11' , 1 -> '10', 2 -> '00' , 3 -> '01').
En revanche, la fin du message est illisible. J'ai passé un moment à comprend d'où venait le problème, jusqu'à lire les warnings retournés par Gnuradio.
Warning: This flow graph may not have flow control: no audio or RF hardware blocks found. Add a Misc->Throttle block to your flow graph to avoid CPU congestion.
J'ai ajouté deux blocs throttle dans le circuit, juste après les deux sources de signaux, ce qui a permis de régler le problème. Voici le schéma bloc complet, tenant compte de ces modifications :

Et le message :

Conclusion
Ce challenge est un bon exercice pour appréhender le traitement d'une trame modulé en Q-PSK, ainsi que l'impact des différents phénomènes perturbateurs (fréquence d'offset, désynchronisation temporelle) pouvant compliquer le décodage. Il existe toutefois une approche plus directe permettant d'obtenir le même résultat, détaillée dans la ici. Celle-ci s'appuie sur des blocs standards pour automatiser les corrections :
- le bloc
Symbol Syncpour assurer la récupération de l'horloge (synchronisation temporelle). - le bloc
Costas Looppour corriger le décalage de fréquence et de phase (fréquence d'offset).
