Principe du traitement GPS pour l’estimation du délai troposphérique par GPS

Les erreurs 

Il existe plusieurs types d’erreurs qui affectent les mesures GPS :

Les erreurs d’orbites

Les erreurs d’orbites sont dues aux imprécisions de positionnement des satellites, en effet il existe des prévisions d’orbites mais celles-ci ne sont pas assez précises et cela influe directement sur la précision des résultats dans la méthode du positionnement absolu. La précision peut être améliorée en utilisant des orbites recalculées à 10cm à postériori fournies par l’IGS. Par la suite, nous utiliserons ces orbites recalculées, elles seront récupérées sur le serveur ftp: ftp://cddis.gsfc.nasa.gov/pub/gps/products/. Au final, ces erreurs sont millmétriques.

Les erreurs d’horloges satellites

Les erreurs de synchronisation des horloges des satellites entraînent des erreurs sur les résultats du calcul de la position. Par exemple, un décalage d’un millionième de seconde de 300 mètres d’erreur sur la position. En réalité, les corrections d’horloges sont contenues dans le message de navigation. D’autres corrections plus précises existent, comme les erreurs d’horloges recalculées par l’IGS disponibles a postériori. Dans la suite, nous utiliserons ces fichiers horloges recalculés récupérés sur le serveur ftp de l’IGS.

Les erreurs d’horloge du récepteur

Les erreurs de synchronisation des horloges du récepteur sont plus importantes que celles des satellites. En effet, les récepteurs possèdent des horloges à quartz qui sont moins précises que les horloges atomiques des satellites. C’est pour cela, que l’on ajoute les décalages d’horloges du récepteur dans les inconnues et qu’il est nécessaire d’observer au moins quatre satellites.

Les perturbations du signal par la ionosphère

La ionosphère est la partie de l’atmosphère située au-dessus de 60 kilomètres d’altitude.
Ces effets sont variables dans le temps et ils dépendent de l’activité solaire. C’est un milieu dispersif, c’est-à-dire que la vitesse de l’onde dépend de la fréquence de l’onde. La traversée de l’onde GPS dans ce milieu provoque une erreur sur la mesure. Dans le cas du positionnement GPS les effets de l’ionosphère se traduisent par un retard de l’onde de code et une avance sur les phases qui occasionne une erreur due à la mesure de la pseudo-distance de 0 à 50m. Dans le cas du positionnement absolu, la combinaison des deux fréquences émises par les satellites permet de réduire jusqu’à quelques millimètres les effets de l’ionosphère. En positionnement relatif sur de petites lignes de base, ces perturbations sont éliminées par doubles différences.

Les perturbations du signal par la troposphère

Lors du passage de l’onde dans l’atmosphère, celle-ci subit deux perturbations, l’une lors de son passage dans l’ionosphère et l’autre durant la traversée dans la troposphère.
Dans la troposphère, partie de l’atmosphère comprise entre le sol et 17 kilomètres d’altitude, la perturbation subit par l’onde dans cette partie de l’atmosphère n’est pas la même que dans la ionosphère : elle est fonction de l’humidité, de la pression, de la température et de l’élévation du satellite. Cette erreur due à la troposphère se traduit elle aussi par un retard de l’onde que nous détaillerons au chapitre 1.4. C’est à partir de ce retard que l’estimation de la vapeur d’eau GPS est possible.

Les variations du centre de phase des antennes

Les antennes GPS effectuent leurs mesures au centre de phase électrique de l’antenne.
Celui-ci n’étant pas confondu avec le centre physique de l’antenne, il est nécessaire de quantifier cet écart. De plus, les signaux GPS sont sur deux fréquences différentes et le centre de phase diffère en fonction de la fréquence. L’excentrement du centre de phase dépend aussi de l’élévation et de l’azimut du satellite observé. Le fabricant ou des centres de recherche fournissent ces écarts du centre de phase avec le centre physique de l’antenne. Dans notre cas, tous ces paramètres seront répertoriés par référence d’antenne dans un fichier. [Doerflinger, 1997]

Les modèles de surcharges océaniques

La Terre subit des déformations dues à la surcharge des marées océaniques sur les plateaux continentaux. Plus précisément, le phénomène des marées provoque des déplacements d’avant en arrière de grande quantité d’eau sur les côtes continentales. Ces mouvements d’eau entraînent des redistributions périodiques des centres de masses du fond de l’océan ce qui implique des mouvements verticaux d’un point terrestre à chaque marée. Il existe plusieurs modèles mathématiques qui prennent en compte ce phénomène. Le modèle le plus utilisé est le FES 2004, dans la suite, les effets de la surcharge océanique seront récupérés sur le site internet : http://holt.oso.chalmers.se/loading/index.htmlqui calcule cette surcharge à partir de la position la station et du modèle de surcharges océaniques choisi. Les surcharges ne sont pas à proprement parlé une erreur mais un déplacement de la position du GPS. Afin de pouvoir comparer deux mesures à des temps différents, il est néanmoins nécessaire de les corriger.

Les paramètres de rotation de la Terre

L’axe de rotation de la terre bouge sur une durée de quatre ans dans un carré de 15 mètres de côté et l’incertitude de son évaluation est d’environ un centimètre. La variation la plus importante est radiale et donc sur la coordonnée verticale. Le service international de la rotation de la terre (IERS) détermine les paramètres décrivant l’orientation de la terre comme la rotation des pôles, correction des modèles de précession nutation. Cette détermination sert dans de nombreux domaines comme l’astronomie, la navigation et la géodésie spatiale. Les résultats de ces calculs sont mis à la disposition des utilisateurs sur le serveur ftp : ftp://cddis.gsfc.nasa.gov/pub/gps/products/. S’ils ne sont pas pris en compte, la répétabilité des mesures se dégradent de l’ordre d’un centimètreenviron sur un an.

Validation du traitement GPSPPP

L’enjeu de cette partie est de valider le logiciel GPSPPP qui repose sur la technique de Positionnement Ponctuel Précis. Cette méthode de calcul présente des avantages par rapport aux méthodes classiques, entre autre un temps de calcul largement réduit et l’utilisation d’un seul récepteur par rapport au positionnement relatif. Pour réaliser cette validation, une campagne de mesure déjà validée par d’autres logiciels sera utilisée.

Présentation du logiciel 

Le logiciel GPSPPP a été réalisé par le laboratoire des Ressources Naturelles du Canada (RNCan), il peut être utilisable directement et gratuitement sur le site internet : http://www.geod.nrcan.gc.ca/online_data_f.php. Ce logiciel permet de calculer les positions de stations fixes, mais aussi les positions de stations cinématiques. Il travaille en deux systèmes de coordonnées le système Canadien NAD 83 et l’ITRS que nous utiliserons par la suite.[KOUBA et al.]
Les avantages principaux de ce logiciel sont :
• de permettre la réalisation de calculs de manière bien plus rapide que des logiciels travaillant en positionnement relatif.
• d’utiliser qu’une seule station réceptrice.
Le positionnement absolu reste a priori moins précis que le positionnement relatif, nous allons quantifier les écarts entre les résultats obtenus par GPSPPP et d’autres résultats de logiciel de positionnement relatif et absolu pour mesurer l’écart de précision qui subsiste entre les deux techniques.
Le laboratoire de Géosciences Montpellier possédant sa propre version de GPSPPP, nous l’utiliserons donc pour tous nos calculs. Dans notre cas, cette version ne possède pas d’interface c’est-à-dire qu’il faudra l’utiliser depuis un environnement LINUX.

Influence de la fonction de projection

Le logiciel GPSPPP possède deux modèles mathématiques de fonctions de projection permettant de ramener les délais troposphériques à des délais zénithaux totaux (cf page 8). La fonction GMF est utilisée par défaut quand les fichiers contenant les mesures météorologiques sont absents et la fonction VMF quand les fichiers sont présents.

Evolution à long terme des coordonnées de la station de Marseille 

Dans le but de connaître les variations Nord, Est et altimétrique de la station MARS de Marseille, on a tracé la série temporelle sur 10 ans (référentiel ITRF 2005). Ces séries temporelles sont issues du site internet du RENAG. La période de mesures sur laquelle les calculs de GPSPPP sont effectués est définie par l’encadré sur la figure 6 ci-dessous. Ce test a pour but de vérifier que GPSPPP ne donne pas de résultats aberrants au premier ordre.
Les tendances de vitesses de déplacements de la station de Marseille calculées à l’aide du logiciel GPSPPP montrent de fortes corrélations avec la série temporelle calculée par le RENAG. Pour la station de Marseille, les vitesses de déplacement sont pour la composante Nord de 16 mm.an -1 , sur la composante Est de 28mm.an -1 et sur la composante altimétrique de -5mm.an -1 . Cela représente un écart respectif avec les résultats du RENAG de 0,1mm.an -1 , de 8,5mm.an -1 et de 6,5mm.an -1 . Au premier ordre, on retrouve les tendances observées sur la figure 6 avec le traitement du RENAG. Cette conclusion est à modérer car la période de mesure de la campagne VAPIMED est très courte il faudra plusieurs années d’observations pour voir les vitesses obtenues par GPSPPP converger vers celles calculées par le RENAG.
Quant à la station AJA2 d’Ajaccio elle a le même type de mouvement que celle de Marseille en théorie, car ces deux sites sont sur la même plaque tectonique. Il est vrai que les résultats de GPSPPP montrent des similitudes entre la stationde Marseille et d’Ajaccio. Les vitesses de déplacement de cette station sont pour la composante Nord de 21 mm.an -1 , sur la composante Est de 29 mm.an -1 et sur la composante altimétrique de 3 mm.an -1 . Cela représente un écart respectif avec les résultats du RENAG de 5,1mm.an -1 , de 9,5mm.an -1 et de 1,6 mm.an -1 . Les écarts des vitesses de déplacement sont donc du même ordre de grandeur entre la station de Marseille calculées par le RENAG et les vitesses dedéplacement de la station AJA2 calculées avec GPSPPP. Donc, nous pouvons dire que les résultats fournis par GPSPPP ne sont pas aberrants à première vue.

Test externes statiques

Le mode statique nécessite d’avoir des coordonnées approchées du point à calculer.
Comme cette station fait partie des réseaux RENAG et RGP, on peut connaître précisément les coordonnées du point. Dans la suite du rapport, le réseau GNSS permanent (RGP) sera noté IGN.
Nous allons réaliser le calcul de la station de Marseille à l’aide du logiciel GPSPPP en mode statique sur la période de mesure de la campagne VAPIMED. Ces résultats seront comparés avec ceux obtenus par les traitements du RENAG et de l’IGN. Afin de pouvoir comparer les résultats de la station d’Ajaccio, un autre calcul a été réalisé avec le logiciel Gipsy à partir des fichiers RINEX issus de cette campagne. Le site internet de l’APPS (Automatic Precise Positioning Service) du service global differential GPS system accessible à l’adresse suivante: http://apps.gdgps.net/apps_file_upload.php permet d’utiliser le logiciel Gipsy gratuitement. Il permet d’effectuer des calculs à partir de fichiers RINEX de taille maximale de 10 Mo pour des utilisateurs enregistrés sur le site et le calcul est réalisé en quelques secondes pour une journée.

Marseille: Mars

Les comparaisons des coordonnées se feront à partirdes variations Nord, Est et hauteurellipsoïdale. Dans le but de comparer plus facilement les résultats de GPSPPP avec ceuxd’autres logiciels, les courbes tracées seront les différences respectives entre les dn, de et dh calculées par GPSPPP et les autres logiciels. La comparaison des variations de GPSPPP par rapport aux autres logiciels nous intéresse, c’est pour cela que les différences de dn, de et dh seront ramenées à leurs moyennes respectives.
Pour visualiser les différences des positions issues des différents traitements, on a retiré respectivement sur chaque composante le résultat du calcul de GPSPPP. Les notations des composantes sur les graphiques seront dn, de, dh. Ensuite, on a calculé la moyenne et l’écart-type pour chaque composante.
Au premier ordre, les différentes courbes suivent les mêmes tendances, ce qui signifieque GPSPPP fait un calcul cohérent par rapport aux autres logiciels. Au deuxième ordre, les calculs du RENAG, d’APPS et de l’IGN sont à moins de 3 millimètres pour la composante Nord, et Est entre 4 et 6 mm pour la composante altimétrique. Ce test a été aussi réalisé entre deux logiciels qui utilisent la méthode du positionnement relatif et les résultats sont du même ordre de grandeur à savoir que pour la composante Nord l’écart est de 1,7 mm, pour la composante Est l’écart est de 2,1 mm et pour la composante hauteur l’écart est de 4,5 mm. Les écarts IGN-RENAG et IGN-GPSPPP et RENAG-GPSPPP sontdu même ordre de grandeur pour ce test.

Optimisation des résultats

Le logiciel GPSPPP en mode cinématique cause des décrochages en début de sessions durant la série temporelle des ZTD lorsque la station est en mouvement (figure 16). Nous faisons l’hypothèse que ces décrochages proviennentde mauvaises coordonnées approchées en entrée du calcul de GPSPPP.
Nous allons donc tenter de vérifier cette hypothèseen rentrant manuellement les dernières coordonnées calculées par le logiciel de la journée précédente en coordonnées approchées pour l’initialisation du calcul de la nouvelle journée. Cette opération sera répétée sur une période de 5 jours. Comme le montre la figure en annexe 1a, il n’y a strictement aucune différence entre les deux calculs. Cela signifie que les décrochages ne sont pas dus à de mauvaises coordonnées approchées.
Ce test ne montrant aucune amélioration, nous allons regarder si ce type d’erreur se produit sur la station fixe de Marseille. Pour cela, on va comparer les calculs en mode statique et en mode cinématique. Pour ce dernier, il est inutile de changer les coordonnées approchées à chaque journée, la station étant fixe, celles-ci ne se modifies pas d’un jour sur l’autre.
Comme le montre la figure en annexe 1b, les décrochages des valeurs de GPSPPP sont du même ordre de grandeur entre les deux calculs soit quelques dizaines de centimètres.
Cependant, le calcul statique montre une valeur aberrante à chaque début de journée. Tandis que dans l’autre cas pour le calcul cinématique, plusieurs valeurs aberrantes apparaissent au cours des deux premières heures et lors de la dernière heure de la journée.
Pour tenter de remédier à ce problème nous allons construire des fichiers RINEX de 48 heures en concaténant deux fichiers de 24 heures dedeux journées successives et cela pour chaque jour. Nous allons regarder ces décrochages sur les coordonnées calculées par le logiciel GPSPPP dans le cas de la station sur le bateau. Au niveau planimétrique, aussi bien les fichiers de données concaténés que les autres ne montrent pas de décrochages significatifs, mis en part la première valeur calculée au changement de journée comme le montre l’annexe 1c. Quant à la composante altimétrique, elle montre plus de différences en début de journée. La figure 17 montre la série temporelle des résultats de la composante altimétrique avec les fichiers concaténés et ceux non concaténés (Série sur 1j). Les fichiers concaténés ont été scindés en deux séries, l’une où le premier jour du fichier est pair (Série paire) et l’autre où le premier jour est impair (Série impaire).
La série temporelle de la hauteur ellipsoïdale figure 17 montre que les trois courbes sont quasiment superposées tout au long de la journée sauf en début de calcul de session où les résultats présentant des décrochages. Dans le cas des fichiers concaténés, ces valeurs aberrantes sont inexistantes lors de la transition entre la première et la deuxième journée.

Traitement de la campagne HyMex

Près de 400 scientifiques d’une vingtaine de pays travaillent sur cette campagne. Il est prévu qu’elle dure 10 ans entre 2010 et 2020 avec un budget estimé à huit millions d’euros. L’objectif de cette campagne est de mieux comprendre le cycle de l’eau en Méditerranée. En effet, le GIEC (Groupe intergouvernemental d’experts sur l’évolution du climat) prévoit une forte augmentation des évènements météorologiques violents mais aussi les sècheresses dans les décennies à venir. La méditerranée étant l’une des zones du monde souffrant le plus de la pénurie d’eau, ces changements climatiques pourraient avoir des conséquences désastreuses. Plusieurs types de mesures seront réalisées dans les airs, la mer et sur terre notamment dans la région méditerranée Nord Occidentale par l’intermédiaire de radiosondages, de ballons dérivants, de mesures océaniques à partir de bateaux, de bouées dérivantes.
La partie de la campagne HyMex à laquelle nous allons nous intéresser sont les mesures GPS réalisées sur un Ferry : le Marfret Niolon transportant des marchandises entre plusieurs villes du pourtour de la Méditerranée tel que Marseille, Barcelone, Alger, Mostaganem (Algérie)…
Depuis le 13 janvier 2012 un GPS a été installé surle bateau, il se situe sur l’avant du navire. La période des mesures que nous allons étudier se situe entre le 13 janvier 2012 et le 14 avril 2013. Voici sur la figure 20 les principaux trajets du bateau.

Traitement des données avec les paramètres de la campagne VAPIMED 

Pour traiter les données avec le logiciel GPSPPP, il est nécessaire de convertir les fichiers T02 en fichier RINEX. Pour cela, nous allons utiliser la fonction TEQC (Translation, Editing, and Quality Check) qui permet de réaliser cette étape. Ensuite, il faut importer les données en respectant l’architecture du logiciel, à savoir en rangeant les fichiers RINEX dans des dossiers par année et des sous-dossiers par jour.
D’autre part, il faut aussi récupérer sur leur siteinternet respectif les fichiers des horloges des satellites (CLK), les fichiers des orbites (SP3), les fichiers des paramètres de rotation de la Terre (ERP), les fichiers VMF ainsi que les fichiers de surcharges océaniques. Le calcul de la campagne HyMex est réalisé en concaténant les fichiers RINEX et en conservant seulement la partie de traitement de midi de la première journée à midi de la deuxième journée comme pour la campagne VAPIMED.

Visualisation des erreurs sur les composantes Nord, Est et altitude

La solution venant immédiatement à l’idée est de tracer directement les trois composantes de déplacement du bateau, mais cette solution ne permettra pas de visualiser les erreurs si elles sont présentes. En effet, le bateau effectue des déplacements sur de longues distances et les erreurs sont très petites par rapport à ses déplacements. C’est pour cela qu’à partir des fichiers concaténés, nous allons créer deux séries de données, la première sera la série impaire dont le premier jour sera impair et la deuxième la série paire dont le premier jour sera pair. Ensuite, sur le même principe que pour la campagne VAPIMED, nous allons soustraire les jours en communs de ces deux séries pour obtenir la figure 22.

Conclusion 

Nous avons réalisé différents tests dans l’objectif de valider le logiciel GPSPPP pour l’estimation du contenu intégré en vapeur d’eau troposphérique par GPS embarqué.
La première étape consistait à réaliser des calculs avec le logiciel GPSPPP sur une campagne déjà validée. Pour cela, nous avons effectués des calculs en mode statique puis en mode cinématique sur deux stations fixes MARS et AJA2. Les résultats réalisés avec ces deux stations fixes nous ont permis de montrer que les résultats obtenus par le logiciel GPSPPP étaient en accord avec d’autres logiciels comme (Gipsy, Bernese, Rocken) c’est-à-dire un écart de 4mm en planimétrique et de 9mm en vertical. De plus, les écarts des délais zénithaux totaux en statique sont en moyenne compris entre 5 et 10 mm en comparaison avec les autres logiciels. Le logiciel GPSPPP est donc validé pour le traitement de station statique.
La deuxième étape consistait à refaire les calculs de la campagne VAPIMED mais cette fois-ci les calculs ont été réalisés avec les données issues d’un ferry qui faisait des liaisons entre Marseille et Ajaccio. Dans ce cas, nous n’avons pu comparer que les valeurs des délais zénithaux totaux et non les valeurs des positions obtenues par GPSPPP par faute de données des autres logiciels. Les résultats obtenus avec GPSPPP en comparaison avec les autres logiciels (Bernese, Rocken) montrent également des résultats concordants. Cette comparaison ne met pas en évidence de dégradations de résultats pour les données issues du bateau par rapport aux autres calculs. Toutefois, il subsiste des décrochages en début de calcul qui n’étaient pas présents pour le calcul de stations statiques.
C’est pour cela qu’il a été nécessaire de réaliser une concaténation de fichiers RINEX sur 48 heures pour conserver seulement les 24 heures au centre de la session.
La troisième étape consistait à traiter la nouvelle campagne HyMex pour laquelle il n’y avait aucun point de comparaison avec d’autres logiciels comme pour la campagne VAPIMED. Le premier calcul a donc été réalisé avec la configuration optimale de GPSPPP obtenue lors de la campagne VAPIMED. Les premiers résultats des délais zénithaux totaux ne se sont pas avérés concluants par rapport à l’ordrede grandeur physique qui est compris entre 2m et 2.5 environ, les valeurs que nous avons obtenues étaient pour certaines négatives jusqu’à -8m par endroit. Nous avons voulu vérifier si le logiciel parvenait à converger vers des solutions stables durant la journée mais sans résultats satisfaisants, celles-ci ne convergeaient pas dans certains cas. Nous nous sommes ensuite intéressés aux valeurs des composantes planimétriques et verticales, restant elles mêmes non concluantes. Nous avons choisi une zone d’étude dans laquelle le logiciel ne convergeait pas. La solution retenue la plus appropriée fut l’actualisation des coordonnéesapprochées avec les dernières coordonnées obtenues de la journée précédente par le logiciel. De plus, il fut nécessaire de mettre les paramètres de surcharges océaniques à zéro pour que le logiciel puisse converger vers unesolution stable. Chose étonnante puisque l’ordre degrandeur des surcharges océaniques est dequelques dizaines de centimètres et ne devrait donc pas empêcher le logiciel de converger.
Après modification du code GPSPPP, les résultats obtenus avec des fichiers non concaténés sont satisfaisants aussi bien en position qu’en délais zénithaux totaux. Toutefois, comme les calculs ont été réalisés avec des fichiers non concaténés, il subsiste des décrochages en début de calcul.
Les améliorations futures qui pourront être apportées au logiciel GPSPPP seraient l’automatisation de l’actualisation des coordonnéesapprochées pour des fichiers concaténés, ou encore faire adapter les paramètres de surcharges océaniques à la zone géographique pour être réintroduits dans le logiciel. Dans notre cas, cette dernière modification n’aurait pas vraiment d’importance car les mesures GPS sont effectuées sur l’eau.
C’est au travers de mesures sur un bateau en mouvement que la méthode PPP montre ses avantages, en effet avec les méthodes de positionnement relatif, il n’est pas possible de réaliser des mesures en cinématique avec des lignesde bases supérieures à 50km. De plus, cette méthode est beaucoup plus rapide en temps de traitement par rapport aux méthodes conventionnelles.

Le rapport de stage ou le pfe est un document d’analyse, de synthèse et d’évaluation de votre apprentissage, c’est pour cela rapport-gratuit.com propose le téléchargement des modèles complet de projet de fin d’étude, rapport de stage, mémoire, pfe, thèse, pour connaître la méthodologie à avoir et savoir comment construire les parties d’un projet de fin d’étude.

Table des matières
Introduction
1. Principe du traitement GPS pour l’estimation du délai troposphérique par GPS
1.1. Le positionnement par GPS
1.1.1. La mesure de pseudo-distance
1.1.2. La mesure de la phase
1.1.3. Le positionnement absolu traditionnel
1.1.4. Le positionnement relatif
1.1.5. Le positionnement absolu avec la méthode PPP
1.2. Le référentiel utilisé
1.3. Les erreurs
1.4. L’estimation du retard troposphérique
2. Validation du traitement GPSPPP
2.1. Présentation du logiciel
2.2. La campagne VAPIMED
2.3. Tests internes du logiciel
2.3.1. Influence des paramètres de rotation de la Terre (ERP : Earth Rotation parameter)
2.3.2. Influence de la fonction de projection
2.3.3. Evolution à long terme des coordonnées de lastation de Marseille
2.4. Positionnement d’un point avec GPSPPP
2.4.1. Paramétrages des calculs
2.4.2. Test externes statiques
2.4.3. Tests externes cinématiques
2.5. Les délais zénithaux totaux
2.5.1. Tests externes statiques
2.5.2. Tests externes cinématiques
2.5.3. Les données issues du bateau
2.5.4. Optimisation des résultats
2.6. Conclusion des différents tests
3. Traitement de la campagne HyMex
3.1. Installation du récepteur et récupération des données sur le bateau
3.2. Traitement des données avec les paramètres de la campagne VAPIMED
3.3. Recherche des erreurs
3.3.1. Visualisation des erreurs sur les composantes Nord, Est et altitude
3.3.2. Visualisation des erreurs des ZTD des fichiers non concaténés sur une journée
3.3.3. Sélection d’une zone d’étude
3.4. Solution proposée
3.5. Calcul de la campagne avec la mise à jour de GPSPPP
Conclusion 
Bibliographie
Annexe 1a
Annexe 1b
Annexe 1c
Annexe 2
Index des figures 
Index des tables 
Index des illustrations

Lire le rapport complet

Télécharger aussi :

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *