Détection d’évènements impulsionnels en environnement radioélectrique perturbé

CycloDet : Détecteur de brouilleurs cyclostationnaires 

Présentation du concept

La plupart des signaux de télécommunications présente une périodicité cachée due aux caractéristiques périodiques incluses dans la construction des signaux, comme par exemple la fréquence porteuse, la vitesse de transmission et le codage. Ces paramètres sont généralement cachés par le caractère aléatoire du message à transmettre. Toutefois, en utilisant une approche cyclostationnaire, qui fut introduite par [Gardner, 1987][Gardner, 1985], cette périodicité peut être récupérée. Vu que les signaux cosmiques peuvent être considérés comme stationnaires au second ordre, le critère de cyclostationnarité peut les séparer des signaux de télécommunications. Une vue exhaustive sur la théorie de la cyclostationnarité et ses applications peut être trouvée dans [Gardner et al., 2006] [Serpedin et al., 2005].

Principe de cyclostationnarité

Considérons un canal fréquentiel donné, la fonction d’autocorrélation du signal x(t) dans ce canal (on omet la fréquence pour plus de lisibilité) est définie par : Rx,x∗ (t, τ ) = D x(t + τ 2 )x ∗ (t − τ 2 ) E ∞ . (2.3) Où h.iN est la moyenne temporelle sur N échantillons. h.i∞ correspond donc à une moyenne asymptotique. – Si Rx,x∗ (t, τ ) est indépendant de t le processus est dit stationnaire, – Si Rx,x∗ (t + T0, τ ) = Rx,x∗ (t, τ ). On dit alors que x(t) présente une périodicité de second ordre, et le processus est dit cyclostationnaire au second ordre. On définit la fréquence fondamentale correspondante et ses harmoniques comme les fréquences cycliques. Ces fréquences cycliques sont notées α = k/T0, avec k ∈ Z. Si x(t) est un processus cyclostationnaire de période cyclique T0, alors Rx,x∗ (t, τ ) est T0-périodique et peut être développée en série de Fourier : Rx,x∗ (t, τ ) = X +∞ α R α x,x∗ (τ )e j2παt (2.4) Avec α = k T0 et k ∈ Z, et où Rα x,x∗ (τ ) sont les coefficients de Fourier définis par : R α x,x∗ (τ ) = Rx,x∗ (t, τ )e −j2παt T0 . (2.5) Rα x,x∗ (τ ) est appelée la fonction d’autocorrélation cyclique. On peut aussi calculer la fonction d’autocorrélation cyclique conjuguée en remplaçant x ∗ (t − τ 2 ) par x(t − τ 2 ) dans l’équation 2.3. 2. La mise en oeuvre pratique : En pratique, la fonction d’autocorrélation cyclique se calcule directement suivant cette formule : R α x,x∗ (τ ) = D x(t + τ 2 )x ∗ (t − τ 2 )e −j2παtE N . (2.6) Autrement dit, Rα x,x∗ (τ ) peut s’interpréter comme la transformée de Fourier sur N échantillons du produit de corrélation instantannée x(t + τ 2 )x ∗ (t − τ 2 ) et la distinction entre signal pollué et signal non pollué se fait simplement par le critère suivant : – Si Rα x,x(∗) (τ ) est nulle pour tout α 6= 0, alors s(t) est stationnaire et on considère que les N échantillons sont non pollués. – Par contre, si ∃ α 6= 0, de sorte que Rα x,x(∗) (τ ) 6= 0, alors on a une signature due à la présence d’une RFI et on considère que les N échantillons sont pollués. 3. Le cas τ = 0 : Remarquons que si on se limite au cas τ = 0 alors le critère devient simplement : R α x,x∗ (0) = x(t)x ∗ (t)e −j2παt N = S(t, f)e −j2παt N (2.7) Autrement dit, de simples transformées de Fourier sur chaque canal du plan T-F permettent de mettre en oeuvre le détecteur cyclostationnaire. 

CYCLODET : DÉTECTEUR DE BROUILLEURS CYCLOSTATIONNAIRES

Remarquons également, que pour α = 0, Rα=0 x,x∗ (0) n’est que le calcul de la puissance moyenne dans un canal fréquentiel. Avant de présenter les 2 détecteurs cyclostationnaires que nous proposons, nous allons examiner plus en détail le cas d’une RFI spécifique, la modulation de phase à 2 états (BPSK, Binary phase Shift Keying), car elle sera utilisée pour toutes les simulations par la suite. 2.4.2 Exemple de la modulation BPSK : Nous avons considéré un signal BPSK (Binary Phase Shift Keying) comme étant le brouilleur qui nous intéresse pour les simulations car il présente des fréquences cycliques à la fois pour la forme directe et la forme conjuguée de la corrélation cyclique. Notons I(t) le signal issu d’une modulation BPSK dans un canal de fréquence donnée (pour simplifier les notations, nous ne notons plus le lien à un canal de fréquence) : I(t) = X k akh(t − kTsym + t0)e j(2πfct+ϕ) (2.8) où les ak sont des valeurs réelles prenant de manière aléatoire et indépendante l’amplitude ±σa (c’est le message binaire qui doit être transmis), h(t) est la réponse impulsionnelle du filtre d’émission, Br = 1/Tsym est le baud rate avec Tsym la durée d’un bit, t0 est le retard, fc et ϕ0 sont respectivement la fréquence et la phase de la porteuse. La figure 2.4.a donne un exemple de I(t). Les fonctions d’autocorrélation cyclique et cyclique conjuguée pour le signal I(t) sont données respectivement par les équations Rα I,I∗ (τ ) =    σ 2 a Tsym e −j2παt0 e j2πfcτ r l Tsym h,h∗ (τ ) si α = l Tsym , 0 si α 6= l Tsym . (2.9) Rα I,I (τ ) =    σ 2 a Tsym e j2ϕe −j2π(α−2fc)t0 r l Tsym h,h (τ ) si α = 2fc + l Tsym , 0 si non . (2.10) Avec l ∈ Z et r l Tsym h,h(∗) (τ ) = Z +∞ −∞ h(t + τ 2 )h (∗) (t − τ 2 )e −j2π l Tsym t dt (2.11) Le détail du calcul menant à ces deux équations est présenté dans [Gardner, 1991]. Ce qui va nous intéresser pour la suite ce sont les amplitudes de ces corrélations cycliques pour le cas τ = 0. Les figures 2.4.c et2.4.d tracent ces modules en fonction de la fréquence cyclique à partir d’une simulation comprenant du bruit gaussien et un signal BPSK. La figure 2.4.b donne un exemple sur une observation astronomique réelle. Les corrélations cycliques sont simplement calculées à partir de l’équation 2.7, c’est à dire par simple transformée de Fourier de kx(t)k 2 pour la corrélation cyclique et x 2 (t) pour sa version conjuguée, où x(t) est l’amplitude mesurée dans un canal fréquentiel donné. Les positions des pics correspondent bien aux fréquences cycliques attendues, c’est à dire α = l/Tsym pour la fonction d’autocorrélation cyclique et α = 2fc + l/Tsym pour la fonction d’autocorrélation cyclique conjuguée (l ∈ Z). L’amplitude de ces pics est régie par le terme r l/Tsym h,h(∗) (τ ) (Equ. 2.11) qui dépend de la réponse impulsionnelle du filtre d’émission, h(t). Pour certaines modulations et pour certains filtres d’émission, le choix de τ = 0 n’est pas optimum pour maximiser ces pics. Ce point est illustré en annexe B. Dans la suite du mémoire, nous considérons systématiquement le cas τ = 0.

Les autres modulation

Le tableau 2.2 donne les fréquences cycliques des corrélations cyclique et cyclique conjuguée pour différentes modulations. Pour le détail du calcul : – la modulation d’amplitude (AM) est dans [Gardner, 1987] pages 400 et 401. – BPSK [Gardner, 1991] – La modulation de phase à K états et la modulation d’amplitude en quadrature (QAM) [Feliachi, 2010] 29 CHAPITRE 2. TRAITEMENT DES INTERFÉRENCES POUR LES OBSERVATIONS DE PULSARS (a) 70 80 90 100 110 120 130 −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1 temps (sans unité) amplitude (sans unité) signal I(t); σ a =1, Tsym=8, fc=0.3125, h(t) filtre gaussien de roll−off 0.5 I(t) message binaire ak T sym (b) 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 100 102 104 106 108 110 112 114 116 118 120 Fréquence cyclique (MHz) amplitude (dB) ||Rα x,x*(0) ||, GPS RFI 1/T sym =1.023 MHZ (c) 0 0.1 0.2 0.3 0.4 0.5 −60 −40 −20 0 Fréquence cyclique α amplitude (dB) ||Rα x,x*(0) || bruit seul 0 0.1 0.2 0.3 0.4 0.5 −40 −20 0 20 40 Fréquence cyclique α amplitude (dB) bruit + BPSK 2/T sym 1/T sym (d) 0 0.2 0.4 0.6 0.8 1 −40 −20 0 Fréquence cyclique α amplitude (dB) ||Rα x,x (0) || bruit seul 0 0.2 0.4 0.6 0.8 1 −20 0 20 40 Fréquence cyclique α amplitude (dB) 2f bruit + BPSK c 2f c +1/T sym 2f c −1/T sym 2f c −2/T sym 2f c +2/T sym Figure 2.4 – La modulation BPSK et ses caractéristiques cyclostationnaires. (a) Représentation temporelle schématisée d’un signal BPSK. Les paramètres de ce signal BPSK sont ceux utilisés pour les figures c et d. (b)Cas cyclique ,kRα x,x∗ (0)k, d’une observation au radiotélescope de Westerbork (Hollande) contenant un signal BPSK dû au satellite GPS. On retrouve la signature d’un des codes du GPS à 1023 MHz. (c)Cas cyclique kRα x,x∗ (0)k : (haut) bruit blanc gaussien seul, (bas) bruit blanc gaussien seul + BPSK (d)Cas cyclique conjugué kRα x,x(0)k : (haut) bruit blanc gaussien seul, (bas) bruit blanc gaussien seul + BPSK. α est la fréquence cyclique normalisée, c’est à dire que ∆f = 1 et ∆t = 1. – La modulation de fréquence à phase continue CPFSK [Agee et al., 1990] [Ferréol et al., 2004] [Gardner, 1987]. – La modulation à déplacement minimum gaussien (GMSK) [Agee et al., 1990] [Gardner, 1987] [Napolitano and Spooner, – La répartition en fréquences orthogonales (OFDM) basée sur la OQAM (Offset Quadrature Amplitude Modulation) [Bölcskei, 2001] [Ciblat and Serpendin, 2003]. 2.4.4 CycloDet et CycloDet ’aveugle’ : A partir des concepts présentés précédememnt, nous proposons les critères de détection CycloDet, Dα x∗ , et Cyclodet conjugué, Dα x , qui se définissent comme la corrélation cyclique ou cyclique conjuguée normalisée par rapport à la puissance du signal : Dα x(∗) =    Rα x,x(∗)    R0 x,x∗ . (2.12) La multitude de modulations utilisées de nos jours peut rendre compliquée la connaissance a priori des fréquences cycliques pour chaque signal de télécommunication. Mais surtout, les changements de fréquences et la succession des bancs de filtres dans la chaîne de réception rendent la détermination des paramètres de la modulation difficilement prédictible (mais pas impossible). Or ces paramètres sont nécessaires pour déterminer les fréquences cycliques à utiliser dans notre critère Cyclodet. Pour remédier à cette difficulté, nous proposons un critère aveugle ne nécessitant pas la connaissance a priori du type de modulation utilisée par les RFIs potentiels. L’idée est de calculer les critères définis par 30 2.5. ESTIMATEUR SPECTRAL DU KURTOSIS Modulation Fréquence cyclique Fréquence cyclique conjuguée AM 0 2fc BPSK l/Tsym 2fc + l/Tsym MPSK et QAM l/Tsym 0 CPFSK l/Tsym 2fc + (2l + 1)fd 0 ≤ l ≤ m−2 2 GMSK l/Tsym 2fc ± l 2Tsym OFDM et OQAM q/Q ; q = 0, Q − 1 α0 + q/Q ; q = 0, Q − 1 α0 = (∆fcmod1) Table 2.2 – Fréquences cycliques et fréquences cycliques conjuguées de différentes modulations (l ∈ Z). −6 −4 −2 0 2 4 6 0 500 1000 1500 2000 2500 3000 amplitude (unité arbitraire) distribution de probabilité sans RFI INR= 3dB INR= −3dB INR= −10dB INR= 10dB 0.048 0.040 0.032 0.024 0.016 0.008 Figure 2.5 – Distribution de probabilité d’un bruit blanc de puissance unité additionné à un signal BPSK (cf. Equ. 2.8) de différentes valeurs de puissance PBP SK. l’équation 2.12 pour un ensemble de fréquences cycliques via une Transformée de Fourier rapide (FFT). Il suffit ensuite de prendre le maximum de cette FFT pour disposer d’un critère qui fonctionne en aveugle. Plus précisement, pour le détecteur cyclique aveugle : Dx∗ = max |{z} α=k/N   x(t)x ∗ (t)e −j2παt
N   hx(t)x ∗(t)iN ! (2.13) avec k = 0, .., N − 1. Pour le détecteur cyclique conjugué aveugle, la formulation est similaire sauf qu’il faut remplacer x ∗ (t) par x(t) dans la partie supérieure de la fraction de l’équation précédente. Pour ces deux détecteurs aveugles, le fait que les fréquences cycliques calculées ne correspondent pas forcément à la vraie valeur de la fréquence cyclique du RFI, a des répercussions sur les performances. Ce point a été étudié dans [Weber et al., 2007]. 

Conclusions

Dans cette partie, nous avons présenté le concept de cyclostationarité. Nous avons proposé deux détecteurs basés sur ce concept. Le premier suppose la connaissance a priori des caractéristiques de modulation de l’interférence et le second travaille en aveugle sans connaissance a priori. Avant d’entamer l’analyse des performances de ces détecteusr dans le cadre des observations de pulsar, nous allons présenter succinctement un autre critére de détection qui est lui aussi fondé sur une caractéristique intrinsèque aux RFIs : la non-gaussianité. 2.5 Estimateur spectral du Kurtosis Comme l’illustre la figure 2.5, la non-gaussianité des RFIs est un critère potentiellement discriminant pour séparer signaux naturels et RFIs. L’estimateur du Kurtosis spectral (SK) a été initialement proposé par [Nita et al., 2007], pour la détection et l’élimination des RFIs en temps réel sur un plan T-F en puissance (i.e. S(t, f)). Considérant un plan T-F de M canaux fréquentiels et N instants temporels. Pour chaque canal fréquentiel f = 0, . . . , M − 1, on calcule la puissance instantanée S(t, f), avec t = 0, . . . , N − 1. Pour chaque canal de fréquence f, le critère de détection défini par [Nita and Gary, 2010b], s’écrit : SK(f) = N + 1 N − 1  NP2(f) P 2 1 (f) − 1  (2.14) avec P1(f) = N X−1 t=0 S(t, f) P2(f) = N X−1 t=0 S(t, f) 2 (2.15) Lorsque le plan T-F qui fournit les pavés S(t,f) est issu d’un banc de filtres implanté sous la forme d’une simple transformée de Fourier, [Nita et al., 2007] a calculé les statistiques de ce critère. Il a ainsi montré qu’en l’absence de RFIs et en présence de signaux naturels réels ayant une statistique gaussienne, la moyenne de ce critère est égale à l’unité et sa variance est donnée par : V ar(SK(f)) ≈    4 N f = 1, .., M/2 − 1 24 N f = 0, M/2 (2.16) D’autres travaux, qui ont pour but la généralisation de l’estimateur spectral de Kurtosis SK pour n’importe quelle loi de distribution, et notamment pour des signaux complexes sont disponibles dans [Nita and Gary, 2010a]. [Gary et al., 2010] a présenté le résultat de l’implantation de SK dans le cadre d’une expérimentation sur des signaux d’observations du Soleil en présence d’interférences. Ces résultats ont démontré que l’algorithme SK est particulièrement bien adapté pour identifier les RFIs intermittentes, mais qu’il aussi fonctionne pour certains types de RFI continues. Cet algorithme se caractérise en outre par sa facilité d’implémentation. Il est donc intéressant de le comparer à l’approche que nous proposons. C’est l’objet des sections suivantes.