Vitesses radiales et astrométrie : Page pour l'impression

Présentation

L'observation directe d'exoplanètes est délicate à cause du grand contraste entre l'étoile et la planète. La plupart des détections sont faites en utilisant un effet de la planète sur l'étoile. Par exemple, la technique des transits consiste à mesurer une baisse de luminosité reçue, éventuellement due à une planète passant devant l'étoile observée.

Les techniques d'observations dont il est question dans ce cours reposent sur un principe simple: si des planètes orbitent autour d'une étoile, alors celle-ci aura un mouvement dépendant des planètes. En l'analysant, on infère la présence ou non de planètes et si oui, leurs orbites.

En effet, en astrométrie on mesure les variations de position angulaire de l'étoile et la détection par vitesse radiale consiste à mesurer la vitesse de l'étoile sur la ligne de visée. Le phénomène physique exploité et la réduction des données est similaire, c'est pourquoi ces deux méthodes seront présentées conjointement. Ce cours a deux objectifs:

Donner une vue d'ensemble de la détection, pour en comprendre les enjeux scientifiques et techniques
Donner les outils pour réaliser un code informatique simple d'analyse de données.

On présentera un modèle d'observations "sans bruit" relativement réaliste, les principales sources d'incertitude, le principe des instruments utilisés et des éléments de traitement de données.

Prérequis

Bases d'algèbre linéaire (opérations sur les matrices, produit scalaire, produit vectoriel) pour la modélisation et l'analyse des données
Notions de base d'optique ondulatoire (intensité lumineuse, spectre lumineux, diffraction)
Notions de probabilité (variable aléatoire, espérance, variance), pour l'analyse des données.

Décrire

Vitesses radiales

Lorsque la distance entre une source d'onde et un observateur varie dans le temps, la longueur d'onde vue par l'observateur varie aussi. Ce phénomène, appelé effet doppler, est ce qui permet la détection du mouvement radial de l'étoile. Cet effet est géométrique, et présent pour tous les types d'ondes. Un exemple bien connu est celui de la sirène d'ambulance: le son qu'elle émet a l'air plus grave une fois qu'elle est passée devant nous. Pour le comprendre, considérons une personne qui se baigne sur une plage (suffisament loin pour que les vagues se brisent derrière elle). Si elle reste statique les vagues vont l'atteindre avec une certaine période. Selon qu'elle s'avance vers le large ou revient vers la plage elle rencontrera les vagues plus fréquemment ou moins fréquemment. De la même manière, plus vite l'étoile s'avance vers la Terre, plus la lumière reçue est décalée vers les hautes fréquences.

L'énergie par longueur d'onde d'une source lumineuse, ou de manière équivalente, le nombre de photons reçus par longueur d'onde, est appelé le spectre de cette source. Cette définition fait intervenir une mesure, le spectre dépend donc de l'observateur. La lumière provenant d'une étoile est composée à priori d'une infinité de longueurs d'ondes. Certaines d'entre elles sont absorbées par l'atmosphère de l'étoile, de sorte qu'elles sont absentes du spectre mesuré. On parle de raies d'absorption. Le spectre d'une étoile a un aspect proche de la figure 1. Si la source lumineuse et l'observateur sont immobiles l'un par rapport à l'autre, l'observateur verra un certain profil spectral $I_0(\lambda)$ (intensité par longueur d'onde $\lambda)$ . Si maintenant la source et l'observateur ont une vitesse relative $\delta V(t)$ au temps , l'observateur verra un spectre décalé $I_t(\lambda) = I_0(\lambda-\delta \lambda(t))$ avec $\frac{\delta \lambda(t)} \lambda \approx \frac{\delta V(t)}{c}$ où est la vitesse de la lumière, dans l'hypothèse où $\delta V$ est très petit par rapport à . En mesurant $\delta \lambda$ on peut remonter à $\delta V$ . La présence d'une planète autour d'une étoile engendre un mouvement périodique de celle-ci, donc un décalage périodique du spectre (voir figure 2).

La mesure du décalage du spectre est utilisée depuis la fin du XIXe siècle pour détecter des étoiles binaires. La vitesse de l'étoile dont on observe le spectre (étoile cible) diminue avec le rapport de masse du compagnon et de l'étoile cible. Plus faible est la masse du compagnon plus il est difficile de le détecter. Les techniques ont été perfectionnées pour détecter des compagnons de masses de plus en plus faible.

Remarque: Une des limitations des mesures de vitesse radiale est qu'on ne mesure qu'une masse minimale car on mesure la projection du mouvement selon une direction. En notant l'angle entre la ligne de visée et le plan orbital d'un objet de masse , on mesure sin(i) .

En 1989, Latham trouve un objet d'une masse minimale de 11 masses de Jupiter (notée M_J ) en orbite autour d'une étoile de type solaire. Etant donné que la limite admise à la masse d'une planète (au delà de laquelle le deutérium entre en fusion) est de 13 M_J , la conclusion de l'article est encore valide: il est possible que ce compagnon soit une Naine brune ou une super Jupiter, selon la valeur de sin(i) . La première planète confirmée comme telle est découverte en 1995 par Michel Mayor et Didier Queloz à l'observatoire de Haute Provence.

Pour donner une idée des ordres de grandeur considérées, la Terre engendre un déplacement du Soleil d'environ 9 cm/s, Jupiter de 12.5 m/s. Un déplacement d'1 m par seconde correspond à un décalage relatif de longueurs d'onde $\frac{\Delta \lambda}{\lambda} = 3 . 10^{-9}$ . La précision attendue doit être maintenue sur plusieurs mois, voire plusieurs années.

Spectre du Soleil dans le visible

Représentation de la lumière du Soleil diffractée dans le visible. Pour une exploration plus précise du spectre du Soleil, voir cette page.

Crédit : Observatoire de Paris

Décalage du spectre

L'étoile et la planète orbitent autour de leur centre de gravité. L'observateur (représenté par une lunette) voit le spectre reçu se décaler et en infère la vitesse V de l'étoile au cours du temps.

Crédit : Observatoire de Paris, ASM, E. Pécontal

Astrométrie

L'astrométrie désigne la mesure de la position des astres. Chronologiquement, l'astrométrie peut être considérée comme la première discipline de l'astronomie. Les anciens s'étaient déjà aperçus que certains astres semblent mobiles: les planètes, mais l'observation à l'oeil nu et avec les premiers instruments astronomiques sont trop imprécises pour détecter le mouvement des étoiles. Ils pensaient donc qu'elles sont fixes, mais ce n'est pas du tout le cas.

Plus précisément, on mesure la position des astres sur la Sphère céleste: Il s'agit d'une manière d'appeler l'ensemble des coordonnées angulaires d'un système de coordonnées sphériques. Selon le contexte, le centre de centre de ce repère est le barycentre du système solaire ou celui de la Terre, l'observateur, ou un autre point. Les directions du repère fixe peuvent être définies par rapport à des étoiles très lointaines (quasars), ou l'intersection du plan de l'orbite de la Terre (l'ecliptique) et de l'équateur.

Les observations doivent être exprimées dans un même système de référence. En l'occurrence, le repère fixe choisi est le référentiel barycentrique du système solaire qui est un référentiel galiléen. De ce fait, on peut modéliser simplement la trajectoire de l'étoile dans ce référentiel. Si l'étoile a des compagnons planétaires, le mouvement dû aux planètes projeté aura l'aspect d'ellipses imbriquées les unes dans les autres (voire figure). Le mouvement rectiligne uniforme a une amplitude bien supérieure au mouvement dû aux planètes, il n'est pas représenté sur la figure pour cette raison.

Une mesure astrométrique comporte toujours plusieurs étoiles dans le champ, de sorte que l'on peut mesurer le déplacement de l'étoile cible par rapport aux étoiles du champ. Ce qui nous intéresse est en effet un mouvement différentiel (on compare les mesures les unes aux autres, on ne cherche pas une position absolue). Un mouvement global des étoiles du champ provient des turbulences atmosphériques ou de bruits instrumentaux, qui peuvent être corrigés.

En mars 2015, deux planètes on été découvertes par astrométrie (DE0823-49 b et HD 176051 b). La résolution nécessaire à la détection de planètes est difficilement atteignable à cause des perturbations atmosphériques. Pour les corriger, il est possible d'utiliser une technique appelée astrométrie différentielle. Sa présentation sort du cadre de ce cours. Il est probable que la mission astrométrique spatiale GAIA, dont la précision est de l'ordre de 20 μas permettra de faire de nombreuses découvertes.

Les mouvements de l'étoile autour du centre de gravité du système planétaire sont de l'ordre du rayon de l'étoile en général. Par exemple, le mouvement de Jupiter entraine un mouvement quasi circulaire du Soleil dont le rayon est environ 1.06 rayon solaire.

L'ordre de grandeur du déplacement observé dépend de la distance de l'étoile, des masses de l'étoile et de la planète et du demi-grand axe de l'orbite.

ordres de grandeur
Type de planète	Distance à l'étoile	Amplitude du mouvement
Jupiter	10 parsec	500 µas
	100 parsec	50 µas
Terre	10 parsec	0.3 µas
	100 parsec	0.03 µas

Projection du mouvement de deux planètes sur la sphère céleste

Observations idéales (sans bruit, mouvement de l'étoile uniquement dû aux planètes) d'un système à deux planètes sur cinq ans. En bleu on représente tout le parcours de la planète et en rouge les quarante-cinq mesures de position.

Exercices

Auteur: Nathan Hara

Questions générales

Difficulté : ☆

Question 1)

Pourquoi peut-on dire que la vitesse radiale et l'astrométrie se basent sur des effets dynamiques ?
Le système étoile + planètes tourne autour de son centre de masse (ou centre de gravité). Si une seule planète de masse orbite autour d'une étoile de masse à une distance de que vaut la distance de l'étoile à , notée $r^\star$ ? Que pensez vous de l'affirmation: "à orbite fixée, plus une planète est massive plus elle est facile à détecter" ?

Auteur: Nathan Hara

Vitesses radiales

Difficulté : ☆

Question 1)

L'effet Doppler existe-t-il seulement pour la lumière ou pour tous les types d'ondes ?
Est-il dû au déplacement de la source, du receveur, ou à leur vitesse relative ?

Auteur: Nathan Hara

Astrométrie

Difficulté : ☆

Question 1)

Une seconde d'arc est égale à un degré divisé par 3600. Une micro seconde d'arc (µas) est un millionième d'une seconde d'arc.

Exprimer une micro seconde d'arc en radians
On considère un bâton de jonglage à extrémités enflammées de 50cm. A quelle distance doit se trouver le jongleur pour que les deux extrémités apparaissent à une une microseconde d'arc l'une de l'autre dans notre champ de vision ? Comparer avec les ordres de grandeurs des mouvements donnés page précédente.

Notez que si on observait uniquement l'effet de la planète si petit soit-il, on aurait toujours accès à toute l'information. Malheureusement, de nombreuses autre sources perturbent le signal, comme on va le voir dans les pages suivantes.

Que mesure-t-on ?

Le principe de ces techniques de détection est simple, mais pour obtenir la précision désirée il faut prendre en compte une grande quantité d'effets affectant les observations. Les mesures sont issues in fine de capteurs CCD, situés au plan focal d'un télescope ou en sortie d'un spectrographe, respectivement dans les cas de l'astrométrie et de la mesure de vitesses radiales. Ces capteurs fonctionnent par effet photoélectrique: lorsqu'un photon les percute, un électron est émis (si le photon a une énergie supérieure à une certaine limite). Ils sont exposés à la lumière pendant un certain temps appelé temps d'intégration. Le nombre d'électrons reçus pendant ce temps est ensuite compté pixel par pixel. Finalement, on obtient un tableau de nombre: le nombre de photons reçus par pixel de la caméra CCD. A chaque mesure est attachée une erreur, calculée selon une méthode explicitée par l'observateur. On effectue ensuite une série d'opérations mathématiques sur les mesures obtenues pour en extraire l'information souhaitée.

Avant d'arriver sur ces capteurs, les photons passent par les instruments, par l'atmosphère, par le milieu interstellaire. De plus, le mécanisme d'émission des photons par les étoiles est complexe: même si l'étoile était fixe par rapport à l'observateur, son spectre et sa position sembleraient variables. Comme les instruments ne permettent pas de résoudre angulairement l'étoile (elle n'apparait que sur un pixel), on mesure la lumière moyenne de sa photosphère, c'est à dire la fine partie de son enveloppe dont la lumière nous parvient. L'astrométrie est sensible à son photocentre. Enfin, il est possible que la lumière reçue provienne partiellement d'une autre source céleste située à proximité de l'étoile, ou du Soleil.

Chacune de ces étapes affecte le signal reçu, de sorte que l'information recherchée n'en représente qu'une petite partie. Par exemple, la vitese radiale apparente d'une étoile est de l'ordre de quelques dizaines de km/s, les techniques actuelles permettent de réduire le bruit instrumental à un peu moins de 0.3-0.5 m/s, qui est aussi l'ordre de grandeur de l'erreur due à une étoile de type solaire. Le signal d'une planète tellurique est de l'ordre de 0.5 m/s, c'est à dire moins que l'ordre de grandeur du bruit. Pour chercher des signaux de plus en plus faibles, il faut améliorer à la fois la modélisation du signal, des instruments, et des méthodes de traitement de données.

Représentation du trajet de la lumière

En mauve, les sources de lumière et en bleu les milieux traversés.

Signal et bruit

Idéalement, on voudrait directement avoir accès aux orbites et aux masses exactes des planètes autour d'une étoile donnée. C'est jusqu'à présent impossible, et probablement pour longtemps ! Les techniques dont il est question ici visent à mesurer la position ou la vitesse de l'étoile, ce sont lessignaux recherchés. Cependant, comme on l'a vu page précédente, la quantité effectivement mesurée contient non seulement le signal mais aussi de nombreux processus aléatoires qu'on appelle bruits. Notre problématique est d'extraire les paramètres du signal (les paramètres des orbites) dans des mesures entachées de bruits.

Pour savoir si un signal sera détectable à priori, on calcule le rapport signal sur bruit qui est le rapport de "la puissance du signal" sur "la puissance du bruit". Ce rapport est fondamental. Le sens précis de cette expression varie selon le contexte mais typiquement, si un signal $x \in \mathbb{R}$ et un bruit $b \in R$ forment une mesure y = x + b , le rapport signal sur bruit est $\frac{x}{b}$ . Plus ce rapport est grand, plus la mesure est précise.

Supposons que l'on veuille estimer une variation de position angulaire ou de vitesse (resp. $\Delta x$ et $\Delta V$ ) d'une étoile de type solaire à 10 parsec autour de laquelle une planète de type jupiter orbite. On a $\Delta x \approx 50 \; \mu as$ et $\Delta V \approx 12 \;m/s$ Les mesures sont contaminées par des variations aléatoires d'amplitudes $x_b \approx 20 \mu$ as (précision de la mission Gaia, mission astrométrique la plus précise) et $V_b \approx 0.5$ m/s (ordre de grandeur pour les meilleurs spectrographes actuels) . Le rapport signal sur bruit est de l'ordre de $\frac{\Delta x}{x_b} = \frac{50}{20} = 2.5$ et $\frac{\Delta V}{V_b} = \frac{12}{0.5} = 24$ . La technique par vitesse radiale est pour l'instant plus précise.

Lorsqu'on mesure une quantité modélisée par une variable alétoire (par exemple le nombre de photons reçus pendant une seconde), le rapport signal sur bruit peut être défini comme le rapport de la moyenne et de l'écart-type.

Attention: il ne faut pas confondre l'amplitude d'un effet indésirable et le bruit qui lui est associé. Par exemple, le mouvement de la Terre dans le système solaire induit une vitesse apparente de plusieurs dizaines de km par secondes. Cependant, sa position est connue avec une très bonne précision, de sorte que l'incertitude liée à la soustraction du mouvement de la Terre est de l'ordre d'1 m/s. Remarque: On distingue en général la précision et l'exactitude. On peut disposer d'un instrument très précis mais comme d'autres effets perturbent la mesure, la mesure donne une valeur inexacte de ce que l'on cherche à mesurer. Dans l'exemple précédent, si on ne soustrait pas la vitesse de la Terre, l'instrument peut être aussi précis qu'on veut on aura un signal parasite entre 1000 et 60000 fois plus gros que celui qu'on cherche à mesurer.

Un problème d'estimation

La problématique de la détection par vitesse radiale et astrométrie est la suivante: quelles sont les orbites et les masses des planètes compagnons de l'étoile cible ? Comme les mesures sont entachées de bruit, on ne peut pas donner des paramètres orbitaux exacts. On estime aussi la "confiance" dans les valeurs des paramètres, souvent donnée sous forme de probabilités. On est aussi amené à se demander si toutes les planètes ont été détectées, et si les détections annoncées sont vraisemblables. En particulier, on étudie la stabilité du système trouvé.

La détection de planètes extrasolaires, comme d'autres problématiques de détection peut être subdivisé comme suit:

Le problème direct, ou modélisation,Notre problème est bien plus complexe que celui de la page précédente. On ne cherche pas à estimer un paramère à partir d'une mesure , mais un vecteur de paramètre $\theta$ à partir de plusieurs mesures aux instants $T=(t_k)_{k=1..m}$ , et les mesures dépendent de ces paramètres et des instants d'observations par une fonction : on a $y = f(\theta,T) + \epsilon$ , où $\epsilon$ est une variable aléatoire suivant une certaine loi de probabilité qui modélise l'incertitude. Le problème direct est la construction de la fonction . En l'occurrence, on utilise les équations bien connues de la mécanique céleste. "Si des planètes avec des paramètres orbitaux connus orbitent autour d'une étoile de masse connue observée à des instants connus, alors les observations seront ".Une première approche est présentée dans la section 1 de la partie "comprendre"
Le problème inverse: étant donné des observations aux instants , =1.. et un modèle (éventuellement un bruit ε), quels sont les paramètres ? C'est le problème de la réduction des données expérimentales, qui est l'objet de la section 3 de la partie "comprendre".
L'instrumentation: Comment construire des instruments performants et modéliser leurs mesures ? La prise en compte des bruits instrumentaux fait partie du modèle.
L'interprétation des résultats: Une fois que les paramètres orbitaux pour plusieurs systèmes planétaires ont été estimés (ainsi que l'incertitude associée), on peut se demander ce que ces résultats signifient Il existe dans le cas des exoplanètes une littérature abondante donnant des scénarios de formation planétaire expliquant les valeurs prises par un système donnée, ou la distribution des éléments orbitaux observés.

Ce cours se focalisera sur les trois premières problématiques, car la quatrième est transverse à toutes les techniques d'observation et fait l'objet d'un cours à part entière.

Exercices

Modélisation

Difficulté : ☆

Question 1)

Quelles sont les précisions annoncées des spectrographes européens ELODIE, HARPS et de celui à venir ESPRESSO ?
Les observateurs attendent que le Soleil soit couché depuis plusieurs dizaines de minutes avant de commencer leurs observations, pourquoi ?
Quel désavantages a-t-on à observer une étoile proche de la Lune ?
Pourquoi observer en altitude ?

Effets à prendre en compte

Les observations traitées sont des mesures de position sur le plan focal d'un télescope et des mesures de spectres lumineux; respectivement pour l'astrométrie et les vitesses radiales. Faire l'hypothèse que les mesures ne sont que le résultat du mouvement d'une étoile parfaitement homogène donnerait des résultats invraisemblables. La fonction de la page précédente prend en général en compte plusieurs effets présentés ici. Compte tenu du format du cours, seuls certains d'entre eux seront développés dans la suite.

Bruit de photon: il s'agit d'un bruit dû à la quantification de la lumière. Les photons sont émis aléatoirement selon ce qu'on appelle un "flux poissonien". Si un détecteur observe une source émettant en moyenne photons sur un certain intervalle de temps $\Delta t$ , le nombre de détections pendant $\Delta t$ sera une variable aléatoire de moyenne et de variance . Cette notion est expliquée plus en détail page bruit de photon
Bruit instrumental: les instruments se comportent en pratique différemment de leurs idéalisation de manière imprévisible. Le miroir du télescope est sensible aux variations de température, l'utilisation d'équipement électronique induit du bruit grenaille (similaire au bruit de photon) et du bruit thermique. Il s'agit de toutes les erreurs dues au trajet de la lumière depuis l'entrée du télescope jusqu'à la lecture du détecteur CCD. Ce bruit peut être déterminé par étalonnage.
Les bruits stellaires: A chaque longueur d'onde , on observe le flux lumineux provenant de toute la partie visible de l'étoile. D'un point de vue mathématique, est intégré sur toute le disque visible de la photosphère de l'étoile: . Au rayon et à l'angle , dépend de la température, de la composition et de la vitesse de la matière par rapport à l'observateur en . est un facteur variant selon les modèles dits de "limb darkening". On peut modéliser des propriétés des variations du flux à une longueur d'onde donnée avec un modèle physique de l'étoile. En pratique, on prend en compte trois types d'effets en compte:
- les oscillations acoustiques (p-modes): la propagation d'ondes mécaniques dans la couche extérieure de l'étoile provoque une contraction/dilatation radiale, qui a une signature de quelques centimètres par seconde dans les données de vitesses radiales (cette vidéo montre une simulation des ondes mécaniques au sein d'une étoile).
- les tâches dues à l'activité magnétique: Le champ magnétique à la surface de l'étoile peut "piéger" localement du plasma, qui se refroidit et devient donc moins lumineux. Ce phénomène a des effets à court terme: une tache présente sur un des côtés de l'étoile brise la symétrie entre la partie "rouge" et la partie "bleue" (voir figure). D'autre part, sur une période d'observation le nombre de taches peut varier de zéro à plusieurs centaines. Lorsque la convection est inhibée par le champ magnétique l'étoile apparaît plus froide, globalement décalée vers le rouge.
- La granulation: la surface de l'étoile étant plus froide que son centre, du gaz chaud remonte à la surface, se refroidit en émettant de la lumière puis redescend. Ce phénomène de convection crée des inhomogénéités de luminosité sur la suface de l'étoile dont la somme donne une vitesse radiale apparente, ou un déplacement du photocentre.
Selon le type d'étoile, ces bruits sont d'importances variables. Certaines étoiles sont trop actives pour avoir une chance de détecter une planète de masse tellurique.

Intensité globale

L'intensité dépend de la composition, de la température

et de la vitesse $\overrightarrow{v}$ en chaque point du disque apparent de l'étoile

Effet de la rotation de l'étoile sur le spectre

La lumière provenant de la moitié de l'étoile ayant un mouvement vers l'observateur est décalé vers les hautes fréquences (vers le bleu). L'autre moitié est décalée vers les basses fréquences (vers le rouge). Dans l'hypothèse où l'étoile est sphérique, et a une luminosité identique partout sur sa surface, le décalage vers le rouge et celui vers le bleu ne fait qu'élargir les raies spectrales. Si une tache est présente, ici sur la partie bleue, la symétrie est brisée et le déficit de lumière entraine un décalage du spectre vers le rouge.

Effets à prendre en compte (2)

Contamination du fond: à cause de la diffraction, il est impossible d'observer un point du ciel, observe toujours le flux moyen d'une certaine région du ciel. Ce flux ne provient pas nécessairement uniquement de l'objet observé. Cette prise en compte est particulièrement importante lorsque l'étoile est dans un système binaire. En effet, Comme un compagnon planétaire, l'autre étoile engendre un mouvement elliptique de l'étoile cible de plus grande amplitude que celui provoqué par une planète. De plus, si la séparation angulaire entre les deux étoiles est plus petite que le champs de vue, la lumière provenant de la seconde étoile contamine les observations.
Système multiple: l'étoile cible peut avoir un mouvement dû à la présence d'autres étoiles proches.
Les perturbations atmosphériques: les turbulences de l'atmosphère aux hautes altitudes font que le signal est diffracté aléatoirement.
La parallaxe: Le mouvement de l'étoile est modélisé dans un référentiel galiléen dont le référentiel barycentrique du système solaire est une approximation. Cependant, les observations ne sont évidemment pas disponibles depuis le centre du système solaire. Il faut prendre en compte le mouvement de l'observateur, qui se manifeste par un mouvement apparent de l'astre.

Quantifier l'incertitude

Pour pouvoir annoncer qu'une planète a été découverte, il faut avoir une certaine confiance dans le résultat. Cependant, comme on vient de le voir, les mesures sont contaminées par de nombreuses sources qui peuvent être aléatoires. La détection d'exoplanète, comme toutes les problématiques de détection en astronomie, donne l'occasion de présenter des modélisations de phénomènes aléatoires. Introduire ces outils permettra en particulier de définir proprement ce que sont les barres d'erreurs. On donnera un sens à une phrase comme: la période d'une orbite estimée est 100 jours, plus ou moins 10 jours à trois sigmas.

Le cadre mathématique utilisé pour les prendre en compte est la théorie des statistiques. Les quantités physiques dont le comportement est imprévisible sont modélisées par des variables aléatoires qui prennent une certaine valeur à chaque mesure. Etant donné que cette notion intervient dans:

La modélisation des variations de luminosité de l'étoile
La modélisation des perturbations atmosphériques
Le calcul de la précision de l'instrument
Le calcul des incertitudes sur les paramètres des orbites estimées

On donnera une description fonctionnelle, permettant de comprendre le cours sans avoir besoin de rentrer dans les détails des statistiques. Le lecteur intéressé pourra se référer à un ouvrage spécialisé, par exemple le cours de Didier Pelat (en libre accès), cité dans la bibliographie.

Une variable alétoire peut être vue comme un programme informatique donnant une valeur réelle à chaque fois qu'on lui demande. La valeur retournée dépendra de la distribution de probabilité de cette variable aléatoire, notée f(x) . Par exemple, si on veut modéliser le lancer d'une pièce de monnaie, chaque lancé correspond à une requête au programme qui retournera "face" ou "pile" avec une probabilité 1/2. Dans le cas de distributions continues, prendra une valeur entre x_1 et x_2 avec une probabilité $\int_{x_1}^{x_2} f(x) dx$ . La distribution la plus utilisée est la gaussienne (parfois appelée courbe en cloche), elle vérifie $f(x) = \frac{1}{\sqrt{2 \pi} \sigma} e^{-\frac{(x-\mu)^2}{\sigma^2}}$ où $\mu$ est sa moyenne et $\sigma^2$ est sa variance. L'analyse des données repose sur l'hypothèse que nous observations des réalisations d'une variable aléatoire. Comme si réaliser une expérience consistait à demander à un ordinateur de sortir une valeur (la mesure) avec une certaine loi de probabilité.

Les variables aléatoires ont deux caractéristiques particulièrement importantes: leur moyenne et leur variance, définie respectivement comme $\mu = \int_{-\infty}^\infty xf(x)dx$ et $\sigma^2 = \int_{-\infty}^\infty (x-\mu)^2f(x)dx$ . La racine carére de la variance est appelée l'écart-type. Cette quantité donne "l'écart-typique" à la moyenne des réalisations. Dans la plupart des cas et en particulier dans le cas gaussien, l'écart-type a une interprétation univoque: il quantifie les déviations à la moyenne. Plus il est petit, plus les valeurs éloignées de la moyenne seront improbables. Dans la limite ou l'écart-type tend vers 0, on aura une variable aléatoire certaine, qui prend la valeur $\mu$ avec une probabilité 1. En physique, on manipule souvent des incertitudes, qui sont implicitement modélisées comme des écart-types de variables aléatoires. Si un observateur dit "j'ai mesuré une vitesse radiale de y = 30.28781 km/s avec une incertitude de 0.00115 km/s", implicitement on suppose qu'il existe une variable aléatoire dont la moyenne correspond à la "vraie" valeur de la vitesse radiale et dont l'écart-type vaut 0.00115 km/s.

Remarque: nous utilisons une terminologie vague, pour des définitions précises voir les références.

Exercices

Auteur: Nathan Hara

Modélisation des vitesses radiales

Difficulté : ☆

Question 1)

Citer des éléments à prendre en compte dans la modélisation des mesures de vitesses radiales ou astrométriques (on ne demande pas une liste très précise, encore moins exhaustive).

Modélisations

Difficulté : ☆☆

Question 1)

La démarche consistant à écrire un problème direct comportant une partie déterministe et une partie aléatoire n'est pas propre à la détection d'exoplanète. En majorité - si ce n'est en totalité - les détections reposent sur la même démarche. Dans chaque item de la liste suivante, on donne un phénomène à observer et un moyen d'observation. Dans chaque cas, lister des effets déterministes et/ou aléatoires à prendre en compte dans le modèle. Par exemple: "Mesurer une tension avec un voltmètre" peut se modéliser par u = u_v + b , où est la tension mesurée, u_v la tension vraie et un bruit gaussien dû aux fluctuations intrinsèques de la tension et aux erreurs de mesure. On peut aussi supposer que le voltmètre réalise en réalité une moyenne de tensions sur un certain intervalle de temps.

Un lancé de dé
Le point d'impact d'un tir ballistique (projectile sans système de propulsion propre) mesuré au gps
Des tendances de ventes de vêtement dans un magasin à partir du livre des comptes tenu par un employé
La mesure de position des corps du système solaire (attention ils ne peuvent pas être modélisés par des sources ponctuelles).

Comprendre

Modélisation - Problème direct

Présentation du modèle

Dans cette section on répond à la question suivante: pour un système planétaire donné, quelles seront les observations ? Dans la section précédente on a introduit la fonction $f: \theta \in \mathbf{R}^n,t \in \mathbf{R}^m , \epsilon \in X^p \rightarrow y \in \mathbf{R}^m$ , où désigne un espace de variables aléatoires, qui donne des observations en fonction d'instants d'observation , d'erreurs alétoires $\epsilon$ et de paramètre du modèle $\theta$ incluant masse de la planète, de l'étoile, paramètres de l'orbites, etc. qui seront précisés. En l'occurrence, l'erreur est bien représentée par un "bruit additif", c'est à dire les observations sont de la forme $y=f(\theta,t) + \epsilon$ , où désigne une fonction modélisant la physique du phénomène observée, celle qui donnerait les observations "parfaites", sans aucune source d'erreur aléatoire. L'erreur $\epsilon$ modélise tous les phénomènes aléatoires intervenant dans les mesures. Selon la distinction adoptée ici, les paramètres du modèle peuvent inclure une source indésirable mais non aléatoire.

L'objet de ce chapitre est de construire cette fonction f_p avec une modélisation physique. La modélisation la plus précise, utilisant tout ce que l'on sait de la physique serait ici inutilement complexe, car compte tenu des erreurs de mesures la différence avec certains modèles plus simples serait si faible qu'il serait impossible de la détecter avec un niveau de confiance acceptable. Le niveau de précision du modèle présenté ici est couramment utilisé par les observateurs.

Remarque: En pratique, les algorithmes d'estimation de paramètres sont des algorithmes d'optimisation, qui nécessitent tous de donner un ou plusieurs points de départs pour la recherche. Lorsque des modèles plus complexes sont utilisés, les estimations précises peuvent être utilisées comme tels.

On présentera d'abord une modélisation physique de l'objet étudié, puis d'autres effets physiques à prendre en compte. Les bruits instrumentaux seront traités dans le chapitre suivant. Les effets physiques "indépendants de l'observateur" qui seront pris en compte sont:

Le mouvement de l'étoile dû aux planètes.
Le mouvement propre: mouvement de translation du barycentre du système {Etoile+Planètes}
Des effets stellaires variables (Ondes de pression, granularité, activité)

Ensuite, on modélise "comment" la lumière émise par l'étoile nous parvient, ce qui amènera à présenter:

La prise en compte du mouvement de la Terre (vitesse de celle-ci dans le référentiel barycentrique du système solaire pour les vitesses radiales, parallaxe et réduction des observations en astrométrie.
Les perturbations atmosphériques
Les variations de parallaxe et de mouvement propre.

Problème à deux corps

Problème à deux corps newtonien

Equation de Newton

Les planètes et l'étoile ont un certain volume. Cependant, on peut montrer que lorsque la distance entre deux corps en intéraction gravitationnelle augmente, leur comportement se rapproche de plus en plus de celui de deux points matériels. On néglige aussi les effets relativistes, de sorte que l'on est ramené à la modélisation de l'intéraction de deux points matériels, dont la résolution va suivre. Considérons deux points matériels (planète) et (pour "star" de sorte à éviter un conflit de notation avec l'anomalie excentrique ) de masses respectives et et l'origine d'un repère Galiléen. Le principe fondamental de la dynamique s'écrit: $m \frac{d^2\overrightarrow{OP}}{dt^2} = -GmM \frac{\overrightarrow{S P}}{S P^3}$ où est la constante universelle de gravitation.. En posant $\overrightarrow{r} = \overrightarrow{SP}$ et $\mu = G(M+m)$ on obtient:

$\frac{d^2 \overrightarrow{r}}{dt^2} = -G(m+M) \frac{ \overrightarrow{r}}{r^3} = -\mu \frac{\overrightarrow{r}}{r^3}$

On va montrer que la solution de cette équation décrit une conique plane, une ellipse dans le cas des planètes liées à une étoile. On suppose que le système est isolé, on peut donc choisir comme origine du repère le barycentre du système, ce qui permet d'obtenir facilement $\overrightarrow{OP}$ et $\overrightarrow{OS}$ par la relation $m\overrightarrow{OP} + M\overrightarrow{OS} = \overrightarrow{0}$ .

On suppose que le système est isolé, donc le mouvement du barycentre du système {Etoile+planètes} est rectiligne uniforme dans un référentiel galiléen (par exemple le référentiel barycentrique du système solaire).

Quantités conservées

L'équation de Newton est une équation différentielle de degré deux sur des vecteurs de $\mathbf{R}^3$ . Pour la résoudre il faut trouver six quantités conservées au cours du mouvement (ou intégrales premières du mouvement) indépendantes. En l'occurrence on peut facilement montrer que deux vecteurs sont conservés au cours du mouvement (ce qui fait deux fois trois composantes, on a bien six scalaires conservés). Notons $\overrightarrow{r}} = r \overrightarrow{u}$ , où $\overrightarrow{u}$ est un vecteur unitaire (et donc est la norme de $\overrightarrow{r}$ ) et $\dot{\overrightarrow{a}}$ la dérivée par rapport au temps d'un vecteur $\overrightarrow{a}$ . On a la conservation du moment cinétique par unité de masse $\overrightarrow{G} = \overrightarrow{r} \wedge \dot{\overrightarrow{r}}$ et du vecteur excentricité $\overrightarrow{P} = \frac{ \dot{\overrightarrow{r}} \wedge \overrightarrow{G} }{\mu} - \overrightarrow{u}$ au cours du mouvement c'est à dire leur dérivée temporelle est nulle

En effet, soit $\overrightarrow{r}(t)$ une solution de l'équation de Newton. Alors $\frac{d \overrightarrow{G}}{dt}(t) = \frac{d( \overrightarrow{r}(t) \wedge \dot{\overrightarrow{r}(t)})}{dt} = \dot{\overightarrow{r}}(t) \wedge \dot{\overightarrow{r}}(t) + \overightarrow{r}(t) \wedge \ddot{\overightarrow{r}(t)} =0$ car d'après l'équation de Newton, $\ddot{\overrightarrow{r}$ est colinéaire à $\overrightarrow{r}$ . D'autre part, $\frac{d\overrightarrow{P}}{dt}(t) = \frac{\ddot{\overrightarrow{r}}(t) \wedge \overrightarrow{G}(t) }{\mu} - \dot{\overrightarrow{u}}(t)$ , or

$\ddot{\overrightarrow{r}} \wedge \overrightarrow{G} = - \frac{\mu}{r^3} \overrightarrow{r} \wedge (\overrightarrow{r} \wedge \dot{\overrightarrow{r}}) = - \frac{\mu}{r^3} \overrightarrow{r} \wedge (\overrightarrow{r} \wedge r \dot{\overrightarrow{u}}) = - \frac{\mu}{r^3} r \overrightarrow{u} \wedge (r \overrightarrow{u} \wedge r \dot{\overrightarrow{u}}) = \mu \dot{\overrightarrow{u}}$

Donc $\frac{d\overrightarrow{P}}{dt} = \overrightarrow{0}$ .

Géométrie de l'orbite

Les quantités conservées définies page précédente permettent de donner une description géométrique de l'évolution de $\overrightarrow{r}$ .

On déduit de la conservation du moment cinétique que le mouvement est plan. En effet, $\overrightarrow{r}\cdot\overrightarrow{G} = \overrightarrow{r} \cdot (\overrightarrow{r} \wedge \dot{ \overrightarrow{r}) }$ et comme $\overrightarrow{r} \cdot (\overrightarrow{r} \wedge \dot{ \overrightarrow{r}) }$ est un produit mixte, il est invariant par permutation circulaire: $\overrightarrow{r} \cdot (\overrightarrow{r} \wedge \dot{ \overrightarrow{r}) } = \dot{\overrightarrow{r}} \cdot (\overrightarrow{r} \wedge \overrightarrow{r}) } = \overrightarrow{0}$ car le produit vectoriel de deux vecteurs colinéaires est nul. Le vecteur $\overrightarrow{r}$ est orthogonal à $\overrightarrow{G}$ à tout instant, autrement dit le mouvement est dans un plan orthogonal à $\overrightarrow{G}$ .

Notons l'angle entre $\overrightarrow{r}$ et $\overrightarrow{P}$ et posons $e = \|\overrightarrow{P} \|$ . Alors comme $\overrightarrow{u}$ est unitaire, $\overrightarrow{P} \cdot \overrightarrow{u} = e \cos v$ d'autre part en remplaçant $\overrightarrow{P}$ par sa définition, on a $\overrightarrow{P} \cdot \overrightarrow{u} = \frac{\overrightrarrow{r} \cdot ( \dot{\overrightarrow{r}} \wedge \overrightarrow{G} )}{r \mu} -1$ et comme $\overrightarrow{r} \cdot ( \dot{\overrightarrow{r}} \wedge \overrightarrow{G} )$ est un produit mixte, $\overrightarrow{r} \cdot ( \dot{\overrightarrow{r}} \wedge \overrightarrow{G} ) = \overrightarrow{G} \cdot ( \overrightarrow{r} \wedge \dot{\overrightarrow{r}} ) = \overrightarrow{G} \cdot \overrightarrow{G} = G^2$ . où est la norme de $\overrightarrow{G}$ . On obtient alors en fonction de :

$r(v) = \frac{p}{1+e \cos v}$

Où $p = \frac{G^2}{\mu}$ , qui est l'équation polaire d'une conique du plan. Cette équation donne une paramétrisation de la solution en fonction de , appelée anomalie vraie. Cependant, nous voulons exprimer la solution en fonction du temps, l'objet de la page suivante est d'exhiber une relation entre et le temps. On sait que lorsque e < 1 , il s'agit de l'équation d'une ellipse. On peut montrer géométriquement que p = a(1-e^2) où est le demi-grand axe de l'ellipse. Si e = 1 ou e>1 , la trajectoire est respectivement parabolique ou hyperbolique. Dans ces deux cas le mouvement n'est pas borné, il concernerait une planète en phase d'éjection, événement dont l'observation est très improbable et indiscernable d'une planète à très longue période ou du mouvement propre sur les données actuelles.

L'orbite est dans un plan perpendiculaire à $\overrightarrow{G}$ , et le vecteur $\overrightarrow{P}$ est parallèle à $\overrightarrow{FA}$ . Notons $\overrightarrow{I}= \frac{\overrightarrow{P}}{\|\overrightarrow{P} \|}$ et $\overrightarrow{K}= \frac{\overrightarrow{G}}{\|\overrightarrow{G} \|}$ . On introduit un vecteur $\overrightarrow{J}$ de sorte que $(\overrightarrow{I}, \overrightarrow{J}, \overrightarrow{K})$ forme une base orthonormale.

Géométrie du mouvement elliptique

Le mouvement est dans le plan $(\overrightarrow{P}, \overrightarrow{Y})$

Loi des aires

Nous avons montré que la solution au problème des deux corps est plane, et peut s'exprimer en fonction de l'anomalie vraie , angle entre le vecteur excentricité et le vecteur $\overrightarrow{r}$ .

$r(v) = \frac{a(1-e^2)}{1+e \cos v}$

Afin d'exprimer en fonction du temps on va introduire successivement deux variables, et appelées respectivement 'anomalie moyenne et 'anomalie excentrique.

Loi des Aires

Avec les notations précédentes, $\overrightarrow{u}(v) = \cos v \overrightarrow{I} + \sin v \overrightarrow{J}$ , donc $\dot{\overrightarrow{u}}(v) = -\sin v \overrightarrow{I} + \cos v \overrightarrow{J}$ en remplaçant dans, $\overrightarrow{G} = r \overrightarrow{u} \wedge (\dot{r}\overrightarrow{u}} + r \dot{\overrightarrow{u}})$ on obtient:

$r^2 \frac{dv}{dt}=G$

Cette équation s'appelle la Loi des aires et signifie que $\overrightarrow{r}$ balaie les aires à vitesse constant . En particulier l'aire totale de l'ellipse est parcourue en un certain temps fixe: le mouvement est périodique. L'aire de l'ellipse vaut:

$\pi a b = \pi a^2 \sqrt{1-e^2} = \int_{0}^{2\pi} \frac{1}{2} r^2 dv = \int_{0}^T \frac{G}{2}dt = \frac{T}{2} \sqrt{\mu a (1-e^2)}$

On retrouve bien la troisième loi de Kepler : $n^2 a^3 = \mu$ avec $n = \frac{2 \pi}{T}$ . On définit alors l'anomalie moyenne par M = n(t-t_0) où t_0 est le temps de passage au périastre ( sur la figure), qui est proportionnelle à l'aire FAC (voir figure).

Exercices

Conservation de l'énergie

Difficulté : ☆☆

Question 1)

Pour établir l'équation du mouvement, les conservations du moment cinétique et du vecteur eccentricité suffisent. On a une autre intégrale du mouvement: l'énergie. Considérons deux points matériels et de masses respectives et . On pose $\overrightarrow{r} = \overrightarrow{SP}$ , $r = \|\overrightarrow{r}\|$ et $\mu = \mathcal{G}(m+M)$ où $\mathcal{G}$ est la constante universelle de gravitation.

Montrer que $\frac{d^2\overrightarrow{r}}{dt^2} = -\mu \frac{\overrightarrow{r}}{r^3}$
Calculer l'énergie potentielle associée à la force de gravitation
On note $v = \left\| \frac{d \overrightarrow{r}}{dt} \right\|$ . Déduire du théorème de l'énergie mécanique que $\frac{1}{2} v^2 - \frac{\mu}{r}$ est une quantité constante au cours du mouvement
Quel est le point de l'orbite où la vitesse de l'étoile est la plus élevée ? la moins élevée ?

Equation de Kepler

L'équation de Kepler est une relation entre l'anomalie excentrique et l'anomalie moyenne, cette page présente un moyen de l'établir. L'aire $\mathcal{A}(FAC)$ est proportionnelle à l'anomalie moyenne .

$\mathcal{A} (AFC) = \frac{M}{2 \pi} \pi a^2 \sqrt{1-e^2} = \frac{1}{2} a^2 \sqrt{1-e^2} M$

L'ellipse de la trajectoire est obtenue par une affinité sur l'axe $\overrightarrow{Y}$ de rapport $\frac{b}{a} = \sqrt{1-e^2}$ . Donc

$\mathcal{A}(AFC') = \frac{\mathcal{A}(AFC)}{\sqrt{1-e^2}} = \frac{a^2 M}{2}$

Par ailleurs en notant l'angle $\widehat{AOC'}$

$\mathcal{A}(AOC') = \frac{1}{2} a^2 E = \mathcal{A}(FOC') + \mathcal{A}(AFC')$

L'aire du triangle FOC' s'obtient facilement car $HC' = HC / \sqrt{1-e^2}$

$\mathcal{A}(FOC') = \frac{1}{2} a^2e \sin E$

On a finalement l'équation de Kepler

$E - e\sin E = M$

Cette équation est "transcendante", en conqéquence il n'existe pas d'expression analytique de en fonction de . Cependant, on peut développer en puissances de .

Anomalies

Représentation de l'anomalie vraie

, excentrique

et moyenne

Solutions en fonction du temps

Il reste à trouver une relation entre et , ce qui s'obtient aisément en exprimant la position (X,Y) de dans le plan $(\overrightarrow{I}, \overrightarrow{J})$ :

$\left\{ \begin{array}{l l l} X =& r \cos v =& a(\cos E - e) \\ Y =& r \sin v =& a \sqrt{1-e^2} \sin E\end{array}$

On en déduit les relations utiles: $\begin{array}{l l} \cos v = \frac{\cos E - e}{1 - e \cos E}, \quad & \sin v = \frac{\sqrt{1-e^2} \sin E}{1 - e \cos E} \end{array}$

A ce point, rappelons que l'objectif est d'exprimer la position et la vitesse de l'étoile comme des observables. En exprimant en fonction de puis en fonction de on a la position du corps à un instant quelconque. En ce qui concerne l'astrométrie, on peut s'arrêter aux équations ci-dessus et simplement faire un changement de référentiel, c'ést à dire exprimer la projection de $X(t)\overrightarrow{I} + Y(t)\overrightarrow{J}$ sur la sphère céleste et choisir un jeu de paramètres $\theta$ à ajuster. Pour les vitesses radiales nous avons besoin de $\frac{dX}{dt} = \frac{dX}{dE} \frac{dE}{dM} \frac{dM}{dt}$ . En dérivant l'équation de Kepler par rapport à on obtient:

$\frac{dE}{dM} = \frac{a}{r}$

D'où:

$\left\{ \begin{array}{l l l} \dot{X} =& -na\frac{\sin E}{1-e \cos E} &= -\frac{na}{\sqrt{1-e^2}} \sin v \\ \dot{Y} =& na\sqrt{1-e^2}\frac{\cos E}{1-e \cos E} &= \frac{na}{\sqrt{1-e^2}} (e+\cos v) \end{array}$

En pratique, il n'est pas nécessaire de calculer en fonction de . Les expressions ci-dessus suffisent. Pour mémoire, en remplaçant $\cos v$ dans la formule trigonométrique $\tan^2 \frac{v}{2} = \frac{1-\cos v}{1+\cos v}$ , on obtient:

$\tan \frac{v}{2} = \sqrt{\frac{1+e}{1-e}} \tan \frac{E}{2}$

Changement de référentiel

Les observations sont disponibles dans un référentiel $(\overrightarrow{i}, \overrightarrow{j}, \overrightarrow{k})$ où $\overrightarrow{k}$ est la direction d'observation et $(\overrightarrow{i}, \overrightarrow{j})$ sont choisis de sorte que le repère est orthonormé direct. Le repère $(\overrightarrow{I}, \overrightarrow{J}, \overrightarrow{K})$ est lui aussi orthormé direct. La matrice de passage du repère orbital au repère d'observation est donc une rotation, que l'on décompose en trois rotations dont les angles ont des noms usuels.

Une rotation d'angle $\omega$ d'axe $\overrightarrow{K}$ (parallèle au moment cinétique), cet angle est appelé argument du périastre. Le repère obtenu par cette rotation est noté $(\overrightarrow{I'}, \overrightarrow{J'}, \overrightarrow{K})$ .
Une rotation d'angle d'axe $\overrightarrow{I'}$ , c'est l' inclinaison
Une rotation d'angle $\Omega$ d'axe $\overrightarrow{k}$ , c'est l'ascension droite au noeud ascendant

$\left( \begin{array}{cc} x & \dot{x} \\ y & \dot{y} \\ z & \dot{z} \end{array} \right) = \mathcal{R}_3(\Omega) \mathcal{R}_1(i) \mathcal{R}_3(\omega) \left( \begin{array}{cc} X & \dot{X} \\ Y & \dot{Y} \\ Z & \dot{Z} \end{array} \right)$

Où:

$\mathcal{R}_1(\theta) = \left( \begin{array}{ccc} 1 & 0 & 0 \\ 0 & \cos \theta & -\sin \theta \\ 0 & \sin \theta & \cos \theta \end{array} \right) \quad ; \quad \mathcal{R}_3(\theta) = \left( \begin{array}{ccc} \cos \theta & -\sin \theta & 0 \\ \sin \theta & \cos \theta & 0 \\ 0 & 0 & 1 \end{array} \right)$

Soit le barycentre du système {planète+étoile}, on note $\mathcal{R}_3(\Omega) \mathcal{R}_1(i) \mathcal{R}_3(\omega) = \left( \begin{array}{ccc} B & G & \star \\ A & F & \star \\ C & D & \star \end{array} \right)$ où A, B, F, G sont appelées les constantes de Thiele-Innes. On conserve la notation classique pour ces constantes, qui sautent quelques lettres de l'alphabet pour une raison inconnue des auteurs. La notation C,D n'est en revanche qu'une convention pour ce cours et ne se trouve pas spécialement dans la littérature. On ne donne pas de nom particulier aux éléments de la dernière colonne de la matrice car étant donné que Z=0 , ils n'apparaissent jamais dans les calculs

Paramètres d'un mouvement à deux corps newtonien

Le plan de référence est ici le plan d'observation.

Paramètres orbitaux classiques

On peut caractériser l'orbite par les éléments suivants:

: le demi-grand axe
: l'excentricité
: la période
$\omega$ : l'argument du périastre
: l'inclinaison
$\Omega$ : l'ascension droite au noeud ascendant

Ces éléments donnent la géométrie de l'orbite. Pour déterminer la position de la planète à un instant donné, il faut de plus connaître l'instant de son passage au périastre t_p . On peut alors calculer $M = 2\pi \frac{t-t_p}{P}$ connaisant on peut calculer l'anomalie excentrique par l'équation de Kepler, puis l'anomalie vraie . On en déduit la position sur l'ellipse par l'équation donnée page "Loi des aires", $r(v) = \frac{a(1-e^2)}{1+e \cos v}$ . Enfin, la position sur l'orbite est donnée par les rotations explicitées ci-dessus.

Conclusion pour le problème à deux corps

Dans un référentiel galiléen quelconque de centre , le mouvement de l'étoile en fonction du temps peut s'écrire:

$\begin{array}{ccc}\overrightarrow{OS}(t) &= \left( \begin{array}{c} x(t) \\ y(t) \\ z(t) \end{array}\right) &= \overrightarrow{OB}_0 + \alpha t \overrightarrow{u} + \overrightarrow{BS}(t) \\ \dot{\overrightarrow{OS}}(t) &= \left( \begin{array}{c} \dot{x}(t) \\ \dot{y}(t) \\ \dot{z} (t) \end{array}\right) &= \alpha \overrightarrow{u} + \dot{\overrightarrow{BS}}(t) \end{array}$

Où $\overrightarrow{OB}_0$ est la position du barycentre du système {Etoile, Planète} à t=0 , $\overrightarrow{u}$ est la direction du mouvement de et $\alpha$ son module. $\overrightarrow{BS}(t) = - \frac{m}{m+M} \overrightarrow{r}$ où $\overrightarrow{r}$ vérifiel'équation de Newton $\frac{d^2 \overrightarrow{r}}{dt^2} = -\mu \frac{\overrightarrow{r}}{r^3}$ , qui est une équation différentielle de degré deux sur l'espace. Lorsque $\mu$ est fixé et la position et la vitesse à t=0 , $(\overrightarrow{r}_0, \dot{\overrightarrow{r}}_0)$ sont connus, la position et la vitesse sont données par le flot: $\begin{array}{cccc } \Phi_{\mu}:& \mathbb{R} \times \mathbb{R}^3 \times \mathbb{R}^3 & \rightarrow & \mathbb{R}^3 \times \mathbb{R}^3 \\ & (t,\overrightarrow{r}_0, \dot{\overrightarrow{r}}_0)& \rightarrow &(\overrightarrow{r}(t), \dot{\overrightarrow{r}}(t))\end{array}$ . En d'autres termes, les sept paramètres $(\mu, X_0, Y_0, Z_0, \dot{X}_0, \dot{Y}_0,\dot{Z}_0)$ définissent une orbite de manière univoque. On peut faire un changement de variables pour décrire l'orbite par un autre jeu de sept paramètres, par exemple $a, e, P, \omega, i, \Omega$ et t_p , qui sont les paramètres classiques présentés dans la section précédente. La plupart des auteurs les utilisent pour ajuster le mouvement des planètes, mais ils ont l'inconvénient d'être très sensibles aux erreurs pour de faibles excentricités et inclinaisons. Pour palier à ce problème on définit:

$\varpi = \omega + \Omega$ : la longitude du périhélie
$k = e \cos \varpi$
$h = e\sin \varpi$
: le demi-grand axe
: l'inclinaison
$\Omega$ : l'ascension droite au noeud ascendant
$\lambda = M + \varpi$ : la longitude

Ces éléments sont des fonctions continues de l'inclinaison et de l'excentricité.

Exercices

Auteur: Nathan Hara

Bonnes et mauvaises configurations d'observation

Difficulté : ☆☆

Question 1)

On définit un repère (x,y,z) où la direction est selon la droite observateur - étoile. Avec les notations des pages précédentes, pour quelles valeurs de l'inclinaison les observations de vitesses radiales ne donnent aucune information sur le système ? lequelles rendent la configuration optimale ? Même question pour les angles et $\Omega$ en astrométrie.

Auteur: Nathan Hara

Qu'est-ce que les vitesses radiales permettent de mesurer ?

Difficulté : ☆☆☆

Question 1)

On définit un repère (x,y,z) où la direction est selon la droite observateur - étoile. On note $\left( \begin{array}{c} \dot{X} \\ \dot{Y} \\ 0 \end{array} \right)$ les composantes du mouvement dans le plan orbital. Calculer la composante du mouvement de l'étoile selon la direction d'observation. On mettra les observations sous la forme $K ( \cos(\omega+\nu) + e\cos(\omega))$ , où est une constante à déterminer. Quels paramètres de l'orbite peut-on observer par vitesse radiale ? Montrer en particulier qu'on ne peut pas déterminer la masse exacte de la planète.

Modèle de la trajectoire

Les expressions établies jusqu'ici permettent de paramétrer le mouvement d'une étoile autour de laquelle une planète orbite, en faisant l'hypothèse que les deux corps sont ponctuels, forment un système isolé. Nous allons utiliser ces expressions pour donner un modèle général, lorsque n_p planètes orbitent autour de l'étoile.

En notant $\overrightarrow{r}_0$ et respectivement la position et la masse de l'étoile et $\overrightarrow{r_k}, m_k$ , $k = \mathbf{1}..n_p$ les positions et masses des n_p planètes , les équations de la mécanique (classique) dans le référentiel barycentrique de ce système sont:

$\begin{array}{cccccc} \ddot{\overrightarrow{r}_0} & = &-G\sum\limits_{k=1}^{n_p} m_k \frac{\overrightarrow{r_0}-\overrightarrow{r_k}}{\|\overrightarrow{r_0}-\overrightarrow{r_k}\|^3} & &=& -GM\left(\sum\limits_{k=1}^{n_p} \frac{m_k}{M} \frac{\overrightarrow{r_0}-\overrightarrow{r_k}}{\|\overrightarrow{r_0}-\overrightarrow{r_k}\|^3 \right) \\ \vdots & & & & & \\ \ddot{\overrightarrow{r}_j} & = &-GM\frac{\overrightarrow{r_j}-\overrightarrow{r_0}}{\|\overrightarrow{r_j}-\overrightarrow{r_0}\|^3}& -G\sum\limits_{k=1, k\not=j}^{n_p} m_k \frac{\overrightarrow{r_j}-\overrightarrow{r_k}}{\|\overrightarrow{r_j}-\overrightarrow{r_k}\|^3}& = &-GM \left(\frac{\overrightarrow{r_j}-\overrightarrow{r_0}}{\|\overrightarrow{r_j}-\overrightarrow{r_0}\|^3} - \sum\limits_{k=1,k\not=j}^{n_p} \frac{m_k}{M} \frac{\overrightarrow{r_j}-\overrightarrow{r_k}}{\|\overrightarrow{r_j}-\overrightarrow{r_k}\|^3} \right) \\ \vdots & & & & & \\ \ddot{\overrightarrow{r}_{n_p}} & = &-GM\frac{\overrightarrow{r_{n_p}}-\overrightarrow{r_0}}{\|\overrightarrow{r_{n_p}}-\overrightarrow{r_0}\|^3}& -G\sum\limits_{k=1}^{n_p-1} m_k \frac{\overrightarrow{r_{n_p}}-\overrightarrow{r_k}}{\|\overrightarrow{r_{n_p}}-\overrightarrow{r_k}\|^3} &=& -GM \left(\frac{\overrightarrow{r_{n_p}}-\overrightarrow{r_0}}{\|\overrightarrow{r_{n_p}}-\overrightarrow{r_0}\|^3} - \sum\limits_{k=1}^{n_p-1} \frac{m_k}{M} \frac{\overrightarrow{r_{n_p}}-\overrightarrow{r_k}}{\|\overrightarrow{r_{n_p}}-\overrightarrow{r_k}\|^3} \right) \end{array}$

En négligeant tous les termes du type $\frac{m_k}{M} \frac{\overrightarrow{r_j} - \overrightarrow{r_k}}{\|\overrightarrow{r_j} - \overrightarrow{r_k} \|^3}$ , c'est à dire l'intéraction entre les planètes, on obtient:

$\begin{array}{cccc} \ddot{\overrightarrow{r}_0} & =& -GM\left(\sum\limits_{k=1}^{n_p} \frac{m_k}{M} \frac{\overrightarrow{r_0}-\overrightarrow{r_k}}{\|\overrightarrow{r_0}-\overrightarrow{r_k}\|^3 \right) = -\sum\limits_{k=1}^{n_p} \frac{m_k}{M} \ddot{\overrightarrow{r_k}} \\ \vdots & & & \\ \ddot{\overrightarrow{r}_j} & = &-GM\frac{\overrightarrow{r_j}-\overrightarrow{r_0}}{\|\overrightarrow{r_j}-\overrightarrow{r_0}\|^3}& \\ \vdots & & & \\ \ddot{\overrightarrow{r}_{n_p}} & = &-GM\frac{\overrightarrow{r_{n_p}}-\overrightarrow{r_0}}{\|\overrightarrow{r_{n_p}}-\overrightarrow{r_0}\|^3}& \end{array}$

En résolvant les n_p problèmes à deux corps associés à chacune des planètes, par $M \overrightarrow{r}_0 + \sum\limits_{k=1}^{n_p} m_k \overrightarrow{r_k} = \overrightarrow{0}$ (dont la première équation du système ci-dessus est la dérivée seconde) on obtient le mouvement de l'étoile. Le modèle de trajectoire complet dans un référentiel galiléen est:

$\begin{array}{ccccccc}\overrightarrow{OS}(t) &=& \left( \begin{array}{c} x(t) \\ y(t) \\ z(t) \end{array}\right) &=& \overrightarrow{OB}_0 + \alpha t \overrightarrow{u} + \overrightarrow{BS}(t) & =& \overrightarrow{OB}_0 + \alpha t \overrightarrow{u} - \frac{1}{M}\sum\limits_{k=1}^{n_p} m_k \overrightarrow{r_k} \\ \dot{\overrightarrow{OS}}(t) &= &\left( \begin{array}{c} \dot{x}(t) \\ \dot{y}(t) \\ \dot{z} (t) \end{array}\right) &=& \alpha \overrightarrow{u} + \dot{\overrightarrow{BS}}(t)&=& \alpha \overrightarrow{u} - \frac{1}{M}\sum\limits_{k=1}^{n_p} m_k \dot{\overrightarrow{r_k}} \end{array}$

Trajectoire observée depuis le barycentre du système solaire

Le modèle précédent donne le mouvement de l'étoile observée dans un référentiel galiléen, en particulier le référentiel barycentrique du système solaire $(O,\overrightarrow{I},\overrightarrow{J},\overrightarrow{K})$ , que l'on munit d'un repère de coordonnées sphériques $(O, \rho, \alpha, \delta)$ repérant l'étoile . Supposons que l'observateur est situé en . A un instant il mesure:

soit la vitesse radiale $\dot{\|\overrightarrow{OS}\|}(t) = \dot{\rho}(t)$
soit la position sur la sphère céleste $(\alpha(t), \delta(t))$ dans le cas de mesures astrométriques

On peut facilement montrer avec un développement limité qu'au premier ordre en $\frac{\| \overrightarrow{BS} \|}{\| \overrightarrow{OS}\|}$ , $\rho(t) = \rho_B(t)+ \overrightarrow{e_{\rho}} \cdot \overrightarrow{BS}(t)$ , $\alpha(t) \cos \delta(t) = \alpha_B(t)\cos \delta(t)+ \frac{1}{\rho_B} \overrightarrow{e_{\alpha}} \cdot \overrightarrow{BS}(t)$ , de même $\delta(t) = \delta_B(t)+ \frac{1}{\rho_B} \overrightarrow{e_{\delta}} \cdot \overrightarrow{BS}(t)$ (à faire en exercice). Supposons que n_p planètes indéxées par gravitent autour de l'étoile en et $(\overrightarrow{I_k}, \overrightarrow{J_k}, \overrightarrow{K_k})$ le repère orbital de la planète . On projette le mouvement $\left( \begin{array}{c} X_k(t) \\ Y_k(t) \\ 0\end{array} \right)$ dans le repère orthonormé direct $(\overrightarrow{e_{\rho}}, \overrightarrow{e_{\alpha}}, \overrightarrow{e_{\delta}})$ associé aux coordonnées sphériques $(O, \rho_B, \alpha_B, \delta_B)$ (voir figure). Avec les notations de la page Changement de référentiel:

$\Large \left( \begin{array}{c} \alpha(t) \cos \delta_B \\ \delta(t) \\ \rho(t) \end{array} \right) = \left( \begin{array}{c} \alpha_B(t)\cos \delta_B \\ \delta_B(t) \\ \rho_B(t) \end{array} \right) - \frac{1}{M} \sum\limits_{k=1}^{n_p} m_k \left( \begin{array}{ccc} \frac{A_k}{\rho} & \frac{F_k}{\rho} & \star \\ \frac{B_k}{\rho} & \frac{G_k}{\rho} & \star \\ C_k & D_k & \star \end{array} \right) \left( \begin{array}{c} X_k(t) \\ Y_k(t) \\ 0 \end{array} \right)$

Observation depuis le barycentre du système solaire

Position de l'étoile en coordonnées sphériques dans le référérentiel barycentrique du système solaire et dans le référentiel translaté au centre de masse de la Terre. Les observations sont disponibles dans le référentiel terrestre, représenté en rouge, qui dépend de l'instant de mesure

Trajectoire observée depuis la Terre

Le modèle précédent donne l'évolution de l'étoile dans le repère barycentrique du système solaire (RBSS). Or les observations sont disponibles depuis la Terre. Si est la position de l'observateur, l'étoile cible on a $\overrightarrow{TS} = \overrightarrow{TO} + \overrightarrow{OS}$

La mesure de vitesse radiale est celle de $\dot{\overrightarrow{TS}}$ , il faut donc soustraire des observations $\dot{\overrightarrow{TO}}$
En astrométrie on observe la projection du mouvement de l'étoile sur la sphère céleste à un certain instant. Non seulement il faut estimer $\overrightarrow{TO}$ comme dans le cas des vitesses radiales, mais aussi s'assurer que toutes les observations sont exprimées dans le même référentiel (par exemple le référenciel barycentrique translaté, sur la figure) alors que chaque observation est dans un référentiel dépendant de l'observation (en rouge sur la figure). L'expression des mesures prises à des instants différents selon une même échelle se pose aussi pour les vitesses radiales, mais passe par l'étalonnage de l'instrument de mesure.

La dérivée temporelle des vecteur est toujours définie par rapport à un référentiel. Ici la notation $\dot{a}$ désigne la dérivée temporelle par rapport au RBSS.

La détermination de la position du centre de masse de la Terre par rapport au barycentre du système solaire est un sujet à part entière. La trajectoire d'un corps céleste au cours du temps dans un référentiel donné est appelée une éphéméride. Les principaux laboratoires de calcul des éphémérides sont le JPL (NASA) et l'IMCCE (Observatoire de Paris). Les liens envoient sur les générateurs en lignes d'éphémérides respectifs des deux laboratoires. .

Vitesses radiales

La vitesse de l'observateur par rapport au barycentre du système solaire peut se décomposer en $\dot{\overrightarrow{OT}} = \dot{\overrightarrow{OC}} + \dot{\overrightarrow{CT}}$ . La partie $\overrightarrow{OC}$ est donné par les éphémérides, comme expliqué plus haut. La correction de la deuxième partie est essentielle: l'observateur parcourt deux fois le rayon terrestre en une nuit, ce qui donne une vitesse d'environ 300 m/s. Si M(t) désigne la matrice de changement de repère entre le référentiel lié à la Terre et le référentiel barycentrique du système solaire, on a $\dot{\overrightarrow{CT}} = \dot{M}(t) \overrightarrow{CT}_{R_t}$ où $\overrightarrow{CT}_{R_t}$ désigne la position de l'observateur dans le référentiel terrestre.

La détermination de est aussi un sujet à part entière, appelé "rotation de la Terre". Les paramètres de rotations officiels sont donnés par le Service de rotation de la Terre au SYRTE (Observatoire de Paris).

Parallaxe (astrométrie)

En ce qui concerne l'astrométrie, il faut exprimer la relation entre les positions mesurées et la position dans le RBSS, ce qui peut se décomposer en deux étapes: passer du référentiel terrestre au RBS translaté au centre de masse de la Terre (passer du référentiel rouge au référentiel noir à droite sur la figure), puis passer du référentiel translaté au RBS. Comme on le verra plus tard, il est aussi possible de se passer de cette étape en prenant un champ contenant des étoiles de référence.

Notons (x,y,z) la position de la Terre dans le RBS. En exprimant $\overrightarrow{OS}$ dans le RBSS de deux manières on obtient une relation entre les coordonnées $(\rho_B, \alpha_B, \delta_B)$ et $(\rho_T, \alpha_T, \delta_T)$ : $\begin{array}{lll} \rho_T \cos \delta_T \cos \alpha_T & = & \rho_T \cos \delta_S \cos \alpha_S -x\\ \rho_T \cos \delta_T \sin \alpha_T & = & \rho_T \cos \delta_S \sin \alpha_S - y \\ \rho_T \sin \delta_T & = & \rho_B \sin \delta_B -z\end{array}$

Lorsque est suffisamment loin, ces expressions différenciées au voisinage de $(\rho_B, \alpha_B, \delta_B)$ donnent au premier ordre en $\Delta \alpha = \alpha_T - \alpha_B$ , $\Delta \delta = \delta_T - \delta_B$ et $\Delta \rho = \rho_T - \rho_B$

$\begin{array}{lll} \Delta \alpha \cos \delta &=& \frac{x}{\rho_B} \sin \alpha - \frac{y}{\rho_B} \cos \alpha \\ \Delta \delta &=& \left( \frac{x}{\rho_B} \cos \alpha + \frac{y}{\rho_B} \sin \alpha \left) \sin \delta - \frac{z}{\rho_B} \cos \delta \\ \Delta \rho &= &-x \cos \alpha \cos \delta - y\sin \alpha \cos \delta + z \sin \delta \end{array}$

On appelle la quantité $\varphi = \frac{1}{\rho_B} \approx \frac{1}{\rho}$ la parallaxe de l'étoile (en général notée $\pi$ ou $\varpi$ , qui sont des symboles déjà utilisés dans le cours). Elle considérée comme constante au cours des observations et est ajustée aux observations en astrométrie.

Remarque: la procédure de changement de référentiel passe par des changements d'échelle de temps (UTC, UT1, TDB...) qui ne seront pas détaillés ici.

Changement de coordonnées

Accélération de perspective et changement de parallaxe (astrométrie)

Accélération de perspective

L'accélération de perspective est un effet purement géométrique, qui tient à la définition du mouvement propre d'une étoile. Le mouvement propre est en effet défini comme $\mu = \frac{V\sin \theta}{\rho}$ où est la vitesse du système observé, $\rho$ la distance entre l'observateur et le système, $\theta$ est l'angle entre la ligne de visée et la vitesse du système observé. Nous supposons que les mesures astrométriques sont sur un plan. Cependant, on mesure la projection du mouvement sur une sphère, donc $\theta$ n'est pas constant au cours du temps. En dérivant par rapport au temps, comme est constant:

$\frac{d \mu }{dt} = - \frac{V}{\rho^2} \sin \theta \frac{d\rho}{dt}+\frac{V}{\rho} \cos \theta \frac{d\theta}{dt}$

On voit sur le triangle OSP que $\frac{d\theta}{dt} = - \mu$ et $V \cos\theta =V_r=\frac{d\rho}{dt}$ . D'où:

$\frac{d\mu}{dt} = - \frac{\mu}{\rho} \frac{d\rho}{dt} - \frac{\mu}{\rho} V_r = -\frac{2\mu}{\rho} V_r$

En pratique on ajuste un terme quadratique $\eta$ , où le mouvement angulaire sur la sphère céleste est donné par $\mu (t) = \mu_0 t + \eta t^2$ .

Variation de parallaxe

La variation de parallaxe vaut: $\dot{\varphi}= -\frac{1}{\rho^2} \frac{d\rho}{dt} = -\frac{1}{\rho^2} V_r= - \varphi^2 V_r$ . Comme elle est du deuxième ordre en $\varphi$ , elle n'est en général pas prise en compte.

Changement de coordonnées

Définition de $\theta$

Etoile

L'ajustement des paramètres orbitaux permet de connaître $\mu = G(m+M)$ ou bien $\mu \sin i$ dans le cas des vitesses radiales. Il est impossible de distinguer les masses et séparément a priori. Cependant, on peut mesurer la masse de l'étoile de sorte à lever l'indetermination.

La modélisation des étoiles permet de distinguer les effets sur le spectre de la variation de leur flux lumineux de la présence de compagnons planétaires. Ces variations sont modélisées par des variables aléatoires suivant une certaine loi de probabilité, dont l'amplitude varie de plusieurs ordres de grandeur selon le type d'étoile. Plus l'étoile est active, plus l'amplitude du mouvement dû à la planète doit être grande pour distinguer la planète du bruit. En particulier certains types d'étoiles sont trop actives pour pouvoir détecter des potentielles super-Terres compagnon.

La description physique des étoiles est complexe car de nombreux phénomènes, tous interdépendants, ont lieu: convection, radiation, magnétisme... Les modèles utilisés en détection par vitesses radiales comprennent trois phénomènes

Les ondes de pression se propageant dans les couches extrérieures de l'étoile entrainent des compressions-dilatations de celles-ci, ce qui change la vitesse radiale moyenne observée par rapport à une photosphère statique.Cet effet affecte principalement les observations de vitesses radiales.
La granulation: phénomène convectif (c'est à dire lié au déplacement de matière dans l'étoile)
L'activité de l'étoile: la présence de régions sombres ou lumineuses à la période de rotation de l'étoile

L'impact de ces effets sur les mesures n'est pas simple à quantifier. Selon le type d'étoile et les instants d'observations, l'effet peut être très variable. Certains auteurs estiment les incertitudes par des simulations numérique, d'autres par des modèles théoriques. Nous donnons une description qualitative de ces phénomènes et une modélisation possible, mais il faut garder à l'esprit que c'est un sujet de recherche ouvert. Dans le cas des vitesses radiales, ces bruits sont classiquement modélisées par un processus stochastique d'une certaine densité spectrale de puissance, notion importante en statistique, qui est définie dans cette page.

Masse de l'étoile

La mesure de la masse de l'étoile peut se faire de deux manières

Par ajustement aux observations de modèles théoriques d'intérieurs stellaires.
Par mesure du mouvement dans un système binaire

Dans le premier cas, si la distance à l'étoile est connue (par exemple par mesure de parallaxe), on peut mesurer sa luminosité intrinsèque (si on ne connaît pas la distance on ne mesure évidemment que la luminosité apparente). Par son spectre, on peut mesurer sa température effective $T_{eff}$ . Des modèles d'intérieurs stellaires permettent ensuite d'évaluer la masse. Cette estimation peut être rafinée avec un modèle d'atmosphère stellaire. On peut alors avoir la gravité à la surface de l'étoile . Comme $L \approx T_{eff}^4 R^2$ et $g \approx \frac{M}{R^2}$ , on peut avoir une estimation de la masse.

Dans le cas des étoiles binaires (systèmes de deux étoiles), le spectre présente des raies des deux étoiles. Le mouvement de ces raies se fait à la même fréquence, mais dans des directions opposées (lorsqu'une étoile approche l'autre s'éloigne). La période de ces mouvements est liée à la masse du système par l'équation de Kepler. L'amplitude relative de ces mouvements permet de déterminer la masse des deux étoiles séparément. Comme dans le cas de la détermination des orbites des planètes, la masse n'est connue qu'à un faceteur $\sin i$ près. Pour lever cette indetermination, il faut déterminer l'inclinaison de l'orbite par rapport à l'observateur . Si on observe des eclipses (une binaire passe devant l'autre), $\approx 0$ . Si le système n'est pas dans cette configuration, on peut séparer angulairement les deux étoiles par des techniques d'interférométrie.

Processus stochastique et Densité spectrale de puissance

Processus stochastique

La notion de densité spectrale de puissance (DSP) n'est pas simple à définir, cependant très utilisée dans la littérature de traitement du signal. Nous donnons une définition mathématique pour qu'il n'y ait pas d'ambiguités mais compte tenu de la sophistication des notions introduites, le lecteur pourra se référer à la description qualitative suivante.

La densité spectrale de puissance est une propriété relative à plusieurs variables aléatoires. Les familles de variables aléatoires peuvent par exemple représenter des mesures sur lesquelles on a une incertitude. A chaque instant de mesure on associe une variable alétoire qui a une certaine densité de probabilité. En physique théorique ou en économie, on rencontre des processus stochastiques continus - typiquement le mouvement brownien, qui représente des mouvements d'atomes ou des fluctuations de prix. Formellement, un processus stochastique est une famille de variables aléatoires indexées par un ensemble totalement ordonné , toutes définies sur le même espace de probabilité ( $\forall t \in T, X_t \in (\Omega,\mathcal{F},\mathbb{P}) )$ . Dans ce cours on aura seulement besoin de $T=\mathbb{R}$ ou $\mathbb{N}$ . On note $\mathbb{E}$ l'espérance mathématique.

Dans le cas général, la densité de probabilité de la variable aléatoire $X_{t}$ (pour $t \in \mathbb{R}$ ) dépend des valeurs prises à d'autres "instants" par les autres variables aléatoires. En particulier on peut s'intéresser à une éventuelle probabilité de périodicité. Par exemple si on modélise un nombre de ventes de vêtement par jour, on verra des ventes plus importantes au moment des soldes (à peu près tous les six mois). La densité spectrale de puissance est un outil qui permet de visualiser ce genre de périodicité. Dans la section suivante, on voit que si on prend une famille de variables aléatoires certaines, c'est à dire que X_t vaut une certaine valeur réelle x(t) avec la probabilité 1, la DSP en une fréquence $\omega$ est égale à $|\hat{x}(\omega)|^2$ , où $\hat{x}$ est la transformée de Fourier de . Si maintenant X_t est une variable alétoire, la DSP sera la "transformée de Fourier typique" d'une réalisation de X_t .

Pour définir cette notion mathématiquement, on doit d'abord introduire les notions de convergences et intégrales en moyenne quadratique. Pour plus de précision le lecteur peut se référer au cours de Timo Koski à KTH.

Intégrale en moyenne quadratique (Mean Square Integral)

Rappelons d'abord que si Y_1, ..., Y_n sont des variables aléatoires et $\phi: \mathbb{R}^n\rightarrow \mathbb{R}^m$ une fonction mesurable alors $\phi(Y_1,...,Y_n)$ est une variable aléatoire. En particulier, si $\lambda$ est un scalaire, $\lambda Y_1$ et Y_1 + Y_2 +...+Y_n sont des variables aléatoires. Pour un rappel sur les variables aléatoires, voir le cours de Didier Pelat.

Soit $(\Omega,\mathcal{F},\mathbb{P})$ un espace de probabilités, on dit que la suite de variables aléatoires $(Y_n)_{n \in \mathbb{N}}$ telle que $E\{Y_n\} < + \infty$ , définies sur cet espace converge en moyenne quadratique si et seulement si:

$\lim\limits_{n \rightarrow \infty} \mathbb{E}\{|Y_n-Y |^2\}=0$

Soit $(X_t)_{\mathbf{t \in \mathbb{R}}}$ un processus stochastique continu ( $T=\mathbb{R}$ ) tel que chacune des variables aléatoires X_t a une espérance finie ( $E\{X_t\} < + \infty$ ). L'intégrale en moyenne quadratique du processus $(X_t)_{t \in T}$ sur l'intervalle [a,b] est définie comme la limite en moyenne quadratique (lorsqu'elle existe) de:

$\sum\limits_{k=1}^n X_{t_k} (t_k-t_{k-1}})$

Pour a=t_0 < t_1 <... <t_n=b et $\lim\limts_{n\rightarrow\infty}\max (t_k-t_{k-1}) = 0$ . On la note alors $\int_b^b X_t dt$ .

On définit alors la densité spectrale de puissance comme:

$P_X(\omega) = \lim\limits_{T\rightarrow \infty} \mathbb{E} \{|\hat{x}_T(\omega)|^2\}$

Où $\hat{x}_T(\omega) = \frac{1}{\sqrt{T}} \int_0^T X_t e^{-i \omega t} dt$

Cette définition un peu complexe peut être vue comme une généralisation de la transformée de Fourier à des processus stochastiques. En effet, lorsque le processus X_t est telle que X_t = x(t) avec une probabilité 1, l'intégrale en moyenne quadratique se comporte comme l'intérale de Riemann, alors $P(\omega)$ est le carré du module de la transformée de Fourier de la fonction réelle d'une variable réelle x(t) . Dans le cas où les X_t sont aléatoire, est le carré de la transformée de Fourier "en moyenne" des réalisations de X_t. Par exemple si X_t modélise une tension mesurée au cours du temps dans une expérience d'électronique réalisée un grand nombre de fois, donnant profils de tension u_k(t) à l'expérience (des réalisations du processus stochastique U_t ), la moyenne des carrés du module des transformées de Fourier des u_k(t) notée $\tilde{P}_n(\omega)$ sera approximativement égal à $P_U(\omega)$ . Si le nombre d'expérience tend vers l'infini $\tilde{P_n} \rightarrow P_U$ en norme 2.

Dans le cas d'un processus stationnaire discret ( $T= \mathbb{Z}$ ), on peut directement définir $\hat{x}_n(\omega)=\frac{1}{\sqrt{2n+1}} \sum\limits_{k=-n}^n X_ke^{-i\omega t_k}$ et $P_X(\omega) = \lim\limits_{n\rightarrow \infty} \mathbb{E} \{|\hat{x}_n(\omega)|^2\}$ .

Processus stationnaires au sens large

La densité spectrale de puissance a une définition plus simple lorsque le processus est stationnaire, c'est à dire lorsque le processus $(X_t)_{t\in T}$ vérifié:

Il existe $\mu \in \mathbb{R}$ tel que $\forall t \in T, \mathbb{E}\{X_t\} =\mu$
La fonction d'autocorrélation $\begin{array}{cccc}C: &\mathbb{R} \times \mathbb{R}: & \rightarrow & \mathbb{R} \\ &t,s & \rightarrow& \mathbb{E}\{X^{\star}_s X_t \} \end{array}$ ne dépend que de , i. e. $\forall t,s, t',s'$ tels que , on a . On note appelle alors habituellement la fonction d'autocorrélation , telle que

Les processus stationnaires modélisent des phénomènes qui ont une certaine invariance dans le temps, en particulier la covariance ne dépend pas de de manière absolue, mais de manière relative à un autre instant. Dans ce cas, la densité spectrale de puissance est égale au carré du module de la transformée de Fourier de la fonction . L'équivalence avec la définition de la densité spectrale de puissance donnée plus haut est établie par le théorème de Wiener-Khinchin.

Granulation

Le phénomène de granulation est lié à la convection du gaz dans l'étoile. La lumière rayonnée par le gaz chaud remontant à la surface va vers l'observateur, la longueur d'onde reçue est donc décalée vers le bleu. En rayonnant, le gaz se refroidit, puis repart vers le centre de l'étoile. Etant moins chaud, il émet moins de lumière, si bien que la lumière est globalement décalée vers le bleu. Ce phénomène est variable dans le temps, donc le décalage vers le bleu aussi. Cette variation peut apparaître dans le spectre et créer des fréquences parasites. La nature aléatoire du phénomène fait que même après ajustement, il reste un bruit résiduel. Pour une étoile de type solaire, il est de l'ordre de 0.5 - 1 m/s sur une observation.

La granulation est en général modélisée par un processus stochastique dits de bruits en créneaui (popcorn noise ou burst noise en anglais). Il s'agit de processus stochastiques pouvant prendre deux valeurs, par exemple -1 ou 1 avec une probabilité de changement suivant une loi de Poisson (loi exponentielle). Si à la valeur passe de 1 à -1, la densité de probabilité pour que la valeur passe à 1 à $t+ \Delta t$ est $e^{-\frac{t}{\tau}}$ où $\tau$ est un réel positif. Pour les vitesses radiales, la densité spectrale de puissance de ces bruits peut être modélisée par::

$P(\nu) = \frac{4 \sigma ^2 \tau}{1+(2\pi\nu\tau)^2}}$

Cette modélisation, due à Harvey (1985) a depuis été revue et d'autres densités spectrales de puissances ont été proposée à partir de simulations 3d de convection au sein d'une étoile. En pratique, le bruit dû à la granulation apparaitra comme un signal périodique de l'ordre de cinq minutes. Cependant, on observe aussi des phénomènes appelés meso-granulation et super-granulation sur des échelles de temps plus longues. La contribution totale de ces bruits est:

$P(\nu) = \frac{4 \sigma_g ^2 \tau_g}{1+(2\pi\nu\tau_g)^2}}+\frac{4 \sigma_{mg} ^2 \tau_{mg}}{1+(2\pi\nu\tau_{mg})^2}}+\frac{4 \sigma_{sg} ^2 \tau_{sg}}{1+(2\pi\nu\tau_{sg})^2}}$

Où les indices g, mg et sg se réfèrent respectivement à la granulation, la méso granulation et la super-granulation. En anticipant un peu sur le troisième chapitre, lorsque ces bruits sont pris en compte, les valeurs des $\sigma$ et $\tau$ sont ajustés sur le spectre de puissance du signal.

Activité magnétique

La formation d'arcs de champ magnétique à la surface de l'étoile inhibe le mouvement des particules, donc réduit la température et provoque donc des tâches sombres. Cet effet à des effets à court termes (à la fréquence de rotation de l'étoile $\approx 1$ mois), et à plus long terme à travers des cycles d'activité magnétique.

A court trerme, la tache introduit une dissymétrie entre la partie de l'étoile tournant vers l'observateur, et la partie s'en éloignant, ce qui engendre un décalage du spectre mesuré. D'autre part, la tache engendre un déplacement du photocentre de l'étoile périodique, pouvant être confondu avec la présence d'une planète. Pour éviter ces confusions, on estime la période de rotation de l'étoile par spectroscopie, et on ajuste des sinusoïdes à cette période et ses premières harmoniques.

L'activité magnétique peut se mesurer à travers divers indicateurs, dont on analyse les corrélations.

Toujours selon le modèle de Harvey (1985), la densité spectrale de puissance du bruit de vitesse radiale induit par une tache solaire est:

$P(\nu) = \frac{4 \sigma_{am} ^2 \tau_{am}}{1+(2\pi\nu\tau_{am})^2}}$

Effet de la rotation de l'étoile sur le spectre

Oscillations radiales (vitesses radiales)

Des inhomogénéités de densité dues aux mouvements convectifs font que des ondes mécaniques se propagent au sein des étoiles. Certaines de ces ondes sont radiales, ce qui provoque un mouvement d'ensemble de la photosphère qui a une signature sur le décalage du spectre mesuré. L'étude de ces ondes est un domaine de la physique stellaire appelé "astérosismologie". Etant donné que la théorie est accessible au niveau licence, nous en donnons des principes généraux.

La théorie procède comme suit: on écrit localement 1) l'équation du mouvement linéarisée au premier ordre au voisinage d'un état d'équilibre, 2) la conservation de la masse ou équation de continuité, 3) l'équation de Poisson, liant le potentiel gravitationnel et la densité, 4) Le premier principe de la thermodynamique. On néglige l'effet du champ magnétique. En général, on fait l'hypothèse que le terme de transfert thermique dans le 1er principe est nul.

Perturbations atmosphériques

Vitesses radiales

L'atmosphère peut influencer la mesure de vitesses radiales de deux manières:

Par atténuation de l'intensité, ce qui augmente le bruit de photon
Par la modification du contenu spectral du signal

La présence de lignes spectrales d'émission ou d'absorption est difficile à corriger. C'est pourquoi on ne considère que des plages de fréquences où l'intensité des raies atmosphériques est inférieure à 1/10000 de l'intensité de la cible. De plus, la réponse de l'atmosphère dépend de la longueur d'onde. Le spectre obtenu est pondéré de sorte à corriger ces inhomogénéités.

Astrométrie

Avant de parvenir au télescope, la lumière issue de l'étoile traverse l'atmosphère. Les turbulences aux hautes altitudes créent des inhomogénéités de densité, donc d'indice optique, qui distordent le front d'onde de la lumière pénétrant dans l'atmosphère. La modélisation classique de ce phénomène passe par l'introduction de "Structure functions" qui modélisent la densité de probabilité des champs de vitesses, densité etc.

Afin de ne pas surcharger le cours, nous ne donnerons qu'un modèle très simple de ce phénomène. Les mouvements de turbulence forment des "cellules" de composition homogène et d'une taille typique d_0 . Vue du dessus, l'atmosphère se comporte approximativement comme un tableau de pixels de taille d_0 dont chaque cellule introduit une déviation angulaire du faisceau incident. Cette déviation est de l'ordre de $\approx \frac{\lambda}{d_0}$ . Dans le cas limite où une cellule a un indice 0 et les cellules avoisinantes ont un indice 1 (sont opaques), on retrouve l'ordre de grandeur de la diffraction. Les fronts d'ondes issus de ces différents pixels arrivent au télescope avec des angles différents, donc donnent chacun une image différente appelé "speckle". L'union de ces speckle forme une tache lumineuse élargie de taille environ égale à $1.22\frac{\lambda}{d_0}$ .

Récapitulation

Les effets expliqués précédemment sont pour la plupart périodiques, donc peuvent être confondus avec le signal d'une planète. Le tableau suivant récapitule ces effets, leurs amplitudes et périodes typiques pour des étoiles de type solaire.

Effets physiques
Effet	Description	Amplitude (vitesse radiale)	Amplitude (astrométrie)	Echelle de temps
Mouvement des planètes	Le mouvement des planètes autour de l'étoile engendre un mouvement périodique	1 cm/s - 200 m/s	Terre: 0.3, 0.03 $\mu$ as Jupiter: 500, 20 $\mu$ as à 1 et 10 parsec resp.	de un jour à plusieurs centaines d'années
Mouvement propre	Mouvement rectiligne uniforme	Dizaines de kilomètres par seconde	10 - 1000 mas pour les étoiles observables	Mouvement non périodique
Taches	La présence de taches sombres ou lumineuses dues à l'activité magnétique provoque une disymétris entre la luminosité de la partie bleue et la partie rouge de l'étoile, engendrant un décalage du spectre ou du photocentre.	1 m/s	Avec ces trois effets, 0.5 - 10 $\mu$ UA	Période de rotation de l'étoile (quelques jours)
Activité magnétique (long terme)	Le nombre de tâche à la surface de l'étoile peut varier de quelques unes à plusieurs centaines. Par exemple le soleil a une périodicité de 11 ans. L'effet précédent est modulé par ces variations à long terme.	10 m/s		1 - 10 ans
Granulation	Le mouvement convectif à l'intérieur de l'étoile provoque un mouvement du gaz dans la photosphère	0.5 - 1 m/s		Granulation: quelques minutes, mesogranulation:
Oscillations radiales (p-modes)	La propagations d'ondes acoustiques dans le manteau de l'étoile entraîne une oscillation de celui-ci	0.5 - 1 m/s	inconnu	5 - 10 minutes
Système multiple	La présence d'autres étoiles orbitant autour de l'étoile cible engendre un mouvement décrit par les mêmes équations que celles utilisées pour les planètes. Les autres étoiles étant beaucoup plus massives que des planètes, l'effet est plus important.	1 - 30 km/s	0.1 - 1 as/an	10 - 100 ans

effets observationnels
Effet	Description	Amplitude (vitesse radiale)	Amplitude (astrométrie)	Echelle de temps
Mouvement de l'observateur dans le RBSS	Comme l'observateur est en mouvement dans un référentiel galiléen,	30 km/s (après soustraction, on a une erreur $\approx$ 0.5 m/s à cause des incertitudes sur la rotation de la Terre et les éphémérides)	0.01 - 0.1 as/an	Une année
Perturbations atmosphériques	La présence de raies spectrales atmosphériques peut perturber les observations par vitesses radiale, et les turbulences dévient les fronts d'ondes, engendrant un déplacement apparent des sources.	0.5 m/s	1 mas	Quelques minutes

Instrumentation & Observations

Objectifs

Jusqu'à présent, on a modélisé la lumière qui parvient à un observateur théorique situé au voisinage de la Terre. Dans cette section, on présentera l'aspect concret de l'observation, en particulier les instruments. En astrométrie on mesure la position sur le detecteur CCD sur le plan focal du télescope, pour les vitesses radiales on utilise un spectrographe dont l'entrée est située au point focal du télescope . Il est essentiel de bien comprendre le fonctionnement des instruments pour des raisons évidentes: ils constituent notre seule source d'information et leur étude permet de mieux modéliser les mesures, donc d'exploiter au mieux les données.

Dans le cas de l'astrométrie comme celui des vitesses radiales, deux conditions sont nécessaires pour obtenir des mesures suffisamment significatives:

La quantité de lumière collectée doit être suffisament grande pour obtenir le signal sur bruit désiré
Les erreurs instrumentales doivent être suffisamment petites pour avoir un rapport signal sur bruit acceptable.

En astrométrie comme en détection par vitesses radiales on place un dispositif sur le plan focal d'un télescope, respectivement des capteurs CCD et l'entrée d'un spectromètre (ou spectrographe). C'est pourquoi les principes généraux des télescopes seront présentés. On évoquera le fonctionnement des spectrographes utilisés pour les détections par mesures de vitesses radiales.

Spectrographe ELODIE

Cet appareil a permis la détection de 51 Pegasi b en 1995 par Michel Mayor et Didier Queloz à l'Observatoire de Haute Provence

Crédit : Observatoire de Haute Provence

Télescope

Un télescope est un appareil permettant de recueillir un rayonnement eléctromagnétique. Pour observer un rayonneme nt dans le visible ou dans l'infrarouge proche, on utilise des télescopes à miroir parabolique. Pour éviter les ambiguités, on définira le plan focal comme le plan perpendiculaire à l'axe optique (ici l'axe de symétrie du télescope) passant par le foyer, et on fait l'hypothèse que les rayons reçus font un angle faible avec l'axe optique. Dans ces conditions, La relation donnant la distance au point focal de l'image d'un rayon arrivant avec un angle $\alpha$ sur le plan focal est en première approximation $FF_1 \approx f\tan \alpha \approx f \alpha$ , où est la distance focale. Pour les angles faibles, on peut travailler avec une lentille équivalente au télescope, de même diamètre et distance focale. On va introduire trois notions de bases sur les télescopes: le champ, la résolution angulaire et la vitesse d'acquisition.

Le champ est la portion du ciel observée par le détecteur du télescope. Comme le détecteur est au plan focal, il est égal à $\tan \frac{s}{f} \approx \frac{s}{f}$

La résolution angulaireest l'angle minimal entre deux sources permettant de les séparer par le système de détection. Cette définition est vague, et Supposons qu'une source ponctuelle émettant à une longueur d'onde $\lambda$ soit placée en un point du plan focal . A cause du phénomène de diffraction, la lumière ne sera pas émise selon une direction unique, mais son énergie sera répartie sur certains angles centrés sur $\alpha = \frac{FP}{f}$ . Comme la lumière suit le même trajet dans les deux sens, ce détecteur reçoit de la lumière provenant de ces angles.

Le rapport du diamètre et de la distance focale du télescope donne la "vitesse" de l'instrument. En effet, considérons une source circulaire de taille angulaire $\alpha$ , son image sur le plan focal est un cercle d'aire $4\pi(f\alpha)^2$ . L'intensité observée est proportionnel à l'aire du télescope, donc l'énergie par unité de temps reçue est proportionnelle à $\left(\frac{d}{f}\right)^2$ . On définit l'ouverture du télescope par $R = \frac{f}{d}$ Le rapport signal sur bruit des mesures de CCD est égal à $\frac{1}{\mu \Delta t}$ où $\mu$ est le nombre moyen d'électron par unité de temps et $\Delta t$ est le temps d'intégration (voir Bruit de photon). On peut calculer le temps d'intégration minimal pour obtenir un certain signal sur bruit compte tenu de l'ouverture du télescope, de la taille de la source et de son intensité. Plus l'aire et grande, plus la distance focale est petite, et plus l'instrument collecte rapidement le nombre de photons requis.

Cependant, lorsque l'ouverture augmente, la résolution angulaire diminue.

Miroir parabolique

Les rayons perpendiculaires à l'axe de symétrie d'une parabole parviennent au point focal

.Ceux qui arrivent avec un certain angle arrivent légèrement décalés.

CCD

Un détecteur CCD est composé de trois partie: des électrodes en silicium, un isolant en SiO_2 et une jonction de semiconducteurs NP. Les électrodes en polysilicium sont reliées périodoquement par un fil conducteur (toutes les trois électrodes sur la figure) de sorte à créer un profil de potentiel électrique alternant puits et régions plates. Le principe d'un capteur CCD est le suivant (voir figure)

Lorsqu'un photon percute un atome du slicium dopé N, un électron est émis par effet photo-électrique.
Une jonction de semi-conducteurs dopés crée un champ électrique (exactement comme dans une diode). Comme la force qui s'exerce sur l'électron est $\overrightarrow{F} = - e \overrightarrow{E}$ où $\overrightarrow{E}$ est le champ électrique, est la charge de l'électron, et que le champ est orienté vers le bas l'électron est piégé sous l'isolant.
L'électron est confiné dans un puit de potentiel électrique, selon l'endroit où il a été émis
Ces trois premières étapes ont lieu pendant tout le temps d'intégration. Ensuite, le nombre d'électron par puit de potentiel est compté, voir (b) sur la figure. Les puits de potentiels sont progressivement décalés et le nombre d'électrons s'y trouvant est compté par puit.

Le contrôle de l'erreur induite par la CCD est primordial et peut s'avérer très complexe. Nous ne rentrerons pas dans ces considérations.

Schéma de principe d'un capteur CCD

Bruit de photon

Le bruit de photon est une limitation fondammentale en astronomie car il est dû à la nature de l'émission de la lumière. En effet, les instants d'émission de deux photons sont indépendants les uns des autres.On peut montrer mathématiquement qu'une suite d'évènement indépendants, sans mémoire et stationnaire est nécessairement un flux de Poisson. Dans notre cas, si on note le temps d'attente entre deux émission de photon, la probabilité que soit supérieure à t+s sachant que T>s est $e^{-\mu t}$ où $\mu$ est une constante.

Cette constante a une interprétation physique: plus elle est grande, plus $e^{-\mu t}$ diminue rapidement à donné, c'est à dire plus la probabilité que soit long diminue. On peut montrer $\mu$ est le nombre d'évènement par unité de temps moyen, proportionnelle à l'intensité. Rappelons que le détecteur capte des photons pendant le temps d'intégration $\Delta t$ . Le nombre de photons reçus pendant ce temps est une variable aléatoire suivant la loi:

$\text{Pr}\{N=n\} = \frac{\mu \Delta t}{n!} e^{-\mu \Delta t}$

Qui est appelée loi de Poisson. La moyenne et la variance d'une telle loi sont égales et valent $\mu \Delta t$ . En conséquence, l'écart-type vaut $\sqrt{\mu \Delta t}$ donc le rapport signal sur bruit (Signal-to-Noise Ratio) est:

$SNR = \frac{\text{ecart-type}}{moyenne} = \frac{\sqrt{\mu \Delta t}}{\mu \Delta t} = \frac{1}{\sqrt{\mu \Delta t}}$

Exercices

Auteur: Nathan Hara

Magnitude et temps d'observation

Difficulté : ☆☆

Question 1)

Au deuxième siècle avant J.-C., l'astronome Hipparque a classé cinquante étoile par ordre de brillance en 6 catégories, les plus brillantes occupant la première. L'échelle de magniture apparente moderne est sous la forme $m = -2.5 \log_{10} \frac{F}{F_0}$ où est le flux lumineux reçu et un flux de référence. Le choix de 2.5 et $F_0 = 3.52 \times 10^{-23}$ W $m^{-2}$ $Hz^{-1}$ fait que les étoiles de classe 1 d'Hipparque aient une magnitude comprise entre 0 et 1, les classes 2 on une magnitude entre 1 et 2 etc. De ce fait, il semble naturel que les objets les moins lumineux visibles à l'oeil nu soient environ de magnitude 6. On considère un télescope de diamètre , on suppose qu'on observe dans le visible $\Delta \lambda = [400,800]$ nm, et que le flux lumineux reçu ne dépend pas de la longueur d'onde sur cette plage.

Sachant qu'un photon a une énergie $h \frac{c}{\lambda}$ où $\lambda$ est la longueur d'onde et la vitesse de la lumière, montrer que le flux de photon reçus $\mu$ (photons par seconde) est lié au flux par $\mu = \pi \frac{D^2}{4} \frac{\Delta \lambda}{hc} F$
Donner l'expression du temps d'intégration pour obtenir un rapport signal sur bruit en fonction de la magnitude observée, de et $\Delta \lambda$
Le bruit de photon est il plus grand pour les étoiles peu lumineuses ou très lumineuses ?

Effet Doppler

On considère une source d'ondes et un observateur. Selon leur vitesse relative, la fréquence reçue par l'observateur varie: c'est ce qu'on appelle l'effet Doppler. Dans le cadre de la physique classique et d'une onde harmonique (purement sinusoidale), on peut le calculer simplement. En effet, notons la distance entre la source et l'observateur, la vitesse de l'onde et. Par définition d'une onde harmonique de pulsation $\omega$ et d'amplitude , l'amplitude mesurée à et vaut

$A\cos \left( \omega\left(t-\frac{x}{c}\right) + \phi\right)$

Où $\phi$ donne en particulier la phase en t=x=0 . Si la distance varie avec le temps selon x=x_0+vt , alors en et on mesure

$A\cos \left( \omega\left(1-\frac{v}{c}\right)t - \omega\frac{x_0}{c} + \phi\right) = A\cos \left( \omega't + \phi'\right)$

Donc pour l'observateur, tout se passe comme s'il recevait une onde de pulsation $\omega' = \omega \left(1-\frac{v}{c} \right)$ ou de longueur d'onde $\lambda' = \frac{1}{1-\frac{v}{c} } \lambda \approx \left(1+\frac{v}{c} \right)\lambda$ et de phase $\phi' =- \omega\frac{x_0}{c} + \phi$ .

Cependant, lorsque est la vitesse de la lumière, la relativité générale donne une description beaucoup plus précise du phénomène. Nous ne présenterons pas le calcul menant à l'expression de l'effet Doppler dans ce cadre mais en donnons l'expression:

$\lambda' = \lambda \frac{1 + \frac{v}{c} \cos \theta}{\sqrt{1-\frac{v^2}{c^2}}}$

Où $\theta$ est l'angle entre la direction de la vitesse relative et la ligne de visée observateur-source. Remarquons que lorsque $\theta$ tend vers $\pi/2$ et v/c tend vers 0, l'expression tend vers l'expression classique.

Spectrographe

L'observation des vitesses radiales nécessite de mesurer des longueurs d'onde très précisément, pour cela on utilise des spectrographe. Les équipes américaines et européennes utilisent des appareils différents, mais dans les deux cas ce sont des spectrographes d'échelle. Le principe d'un tel instrument est d'observer simultanément plusieurs ordres élevés de diffraction à l'aide de deux diffractions successives. La lumière est d'abord diffractée par un premier réseau. Un dispositif, appelé "echelle grating" est placé à un certain angle (blazing angle) du premier réseau de sorte à recevoir des ordres élevés de la première diffraction, qui sont diffractés à nouveau.

Ce dispositif permet "d'étaler" le spectre de sorte qu'une rangée de détecteurs CCD reçoit des longueurs d'ondes très proches, ce qui permet une haute résolution spectrale. En contrepartie, l'énergie est elle aussi répartie, ce qui augmente le temps d'intégration nécessaire pour recevoir suffisamment de lumière pour obtenir un certain rapport signal sur bruit.

D'une mesure à l'autre, à cause de variations internes à l'instrument (température, pression), une longueur d'onde donnée peut se décaler. Comme les mesures doivent pouvoir être comparées entre elles; ce problème doit être résolu efficacement: il faut étalonner l'instrument. Sur ce point, les instruments européens et américains diffèrent. Pour les premiers: ELODIE, CORALIE, HARPS, HARPS-N, l'étalonnage se fait en observant simultanément l'étoile cible est une source dont le spectre est connu. ELODIE observe le ciel, HARPS une lampe thorium-argon calibrée et HARPS-N utilise deux calibrations: un spectre de Fabry-Perot et un "laser frequency comb" (un laser dont le spectre est constitué de raies régulièrement espacées). Les instruments américains font passer la lumière par une cavité contenant de l'iode, dont la position des raies d'absorption est connue. Dans les deux cas on peut comparer Les raies du spectre de référence et celles de l'étoile observée. Si leur déplacement est corrélé (elles se décalent simultanément), il est dû à l'instrument.

La résolution spectrale de ces spectrographes, c'est à dire le rapport $R=\frac{\lambda}{\delta \lambda}$ d'une longeur d'onde $\lambda$ et de la sa variation détectable par le dispositi f $\delta \lambda$ est de l'ordre de 100000. Des simulations numériques (Hatzes & Cochran 1992) ont montré que l'écart type sur la mesure finale de vitesse radiale $\sigma_{rv}$ vérifie:

$\sigma_{rv} = k I^{-\frac{1}{2}} {\Lambda}^{-\frac{1}{2}} R^{-1}$

Où est l'intensité reçue, $\Lambda$ est la plage de fréquences considérées et est une constante de proportionnalité. Comme certaines longueurs d'ondes jugées contaminées peuvent être exclues de certaines mesures, cette valeur varie d'une mesure à l'autre.

Principe du spectrographe d'échelle

Le rayon incident est d'abord diffracté sur une grille standard (std. grating), puis à nouveau diffracté par le réseau d'échelle. Les trois spectres finalement obtenus sont reçus par des capteurs CCD.

Crédit : "Echelle Principle" by Boris Považay (Cardiff University) - Own work. Licensed under CC BY-SA 2.5 via Wikimedia Commons

Cross-correlation Function

Nous avons expliqué que la mesure des vitesses radiales se fait par mesure du déplacement du spectre. Cependant, nous n'avons pas précisé comment ce déplacement était mesuré. Il y a a priori beaucoup d'estimateurs classiques. La méthode majoritairement adoptée repose sur une "cross correlation function" (CCF) ou fonction de corrélation.

Les étoiles ont certains types de spectres que l'on sait reconnaître. Pour une étoile donnée, on définit un "masque" $M(\lambda)$ valant 1 pour $\lambda$ correspondant à une raie d'absorption de l'étoile et 0 ailleurs. Le spectre de l'étoile observée $I(\lambda)$ est multiplié par ce masque décalé d'une valeur $\delta \lambda$ et on mesure $CCF(\delta \lambda) = \int_{\lambda_0+\delta \lambda}^{\lambda_1+\delta \lambda} M(\lambda + \delta \lambda) I(\lambda) d \lambda$ . Où $\lambda_0$ et $\lambda_1$ désignent les bornes inférieures et supérieures du spectre observé. Ensuite, une fonction gaussienne $Ae^{\frac{(\lambda- \mu_{\lambda})^2}{\sigma^2}}$ est ajustée sur la CCF. La valeur $\mu_\lambda$ correspondant au minimum de la fonction ajustée est prise comme valeur moyenne du déplacement. L'analyse consiste ensuite à comparer les $\mu_\lambda$ issus d'observations différentes.

La CCF n'est elle même pas symétrique. Ses propriétés d'assymétrie sont analysées, car elles sont significatives d'effets physiques.

Réduction - Problème inverse

Objectif

La modélisation physique du phénomène nous a permis d'obtenir un modèle $f(\theta,T)$ où $\theta \in \mathbb{R}^n$ sont les paramètres du modèle et $T \in \mathbb{R}^m$ sont les instants d'observation. De nombreuses sources d'erreurs ont aussi été listées. On se pose maintenant la question suivante: comment estimer les paramètres du modèle compte tenu des observations ?

Après avoir listé diverses sources de signal et de bruits de l'émission de la lumière à la valeur donnée par le déteteur, nous allons maintenant donner une expression finale au modèle. Cette expression n'est pas universelle, et selon l'étoile observée, la précision recherchée, d'autres formulations peuvent être préférables.

Après avoir établi ce modèle, on donnera quelques principes d'analyse statistique et des moyens algorithmiques pour estimer les paramètres du modèle. On verra en particulier que l'analyse "dans le domaine fréquenciel" est particulièrement importante.

Modèle final

Le modèle comportera une partie déterministe et une partie aléatoire. Pour l'astrométrie, la position $x(t) = \alpha(t) \cos \delta(t)$ $y(t)) = \delta(t)$ sur la sphère céleste est

$\begin{array}{ccc} x(t) = x_0 + \mu_x t + \eta_x t^2 + \sum\limits_{k=1}^{n_p} B_j X_j(t) + G_j Y_j(t) + \varphi \Pi_x(t) + S_x(t) + \Delta x_{atm}+ \epsilon_{S_x} + \epsilon_{x,atm} + \epsilon_x\\ y(t) = y_0 + \mu_y t + \eta_y t^2 + \sum\limits_{k=1}^{n_p} A_j X_j(t) + F_j Y_j(t) + \varphi \Pi_y(t) + S_y(t) + \Delta y_{atm} + \epsilon_{S_y}+ \epsilon_{y,atm} + \epsilon_y \end{array}$

Où (x_0,y_0) est la position initiale de l'étoile, $\mu_x, \mu_y$ sont les composantes du mouvement propre $\eta_x, \eta_y$ sont des termes d'accélération de perspective, B_j, G_j, A_j, F_j sont les constantes de Thiele-Innes de la planète , X_j et Y_j sont ses coordonnées sur son plan orbital, $\varphi$ est la parallaxe, $\Pi_x, \Pi_y$ sont les coefficients parametrant le mouvement de la Terre, Les $\epsilon$ sont les bruits résiduels modélisés par des bruits gaussiens.: $\epsilon_{S_x}, \epsilon_{S_y}$ sont les bruits stellaires, $\epsilon_{x,atm}, \epsilon_{y,atm}$ sont les bruits atmosphériques et $\epsilon_x, \epsilon_y$ représentent des bruits instrumentaux.

$V(t) = \dot{\overrightarrow{OT}} + \mu_z + K\left( \cos(v+\omega) + \cos \omega \right) + S_z(t) + z_{atm}(t)+ \epsilon_{S} + \epsilon_{mes}$

Où $\dot{\overrightarrow{OT}$ est la vitesse de l'observateur dans le référentiel barycentrique du système solaire, $\mu_z$ est la composante du mouvement propre dans la direction radiale, S_z(t) est le signal stellaire dû à la granulation, aux oscillations et à l'activité, $\epsilon_S$ est le bruit stellaire résiduel et $\epsilon_{mes}$ le bruit associé à la mesure.

Dans les deux cas, les techniques de réduction de données visent à trouver des paramètres qui sont "plausibles", en l'occurrence, qui reproduisent les observations.

Traitement statistique

Rappelons qu'une variable aléatoire peut être vue comme un programme informatique qui délivre des valeurs suivant une certaine distribution de probabilité lorsqu'on lui demande. Le problème que nous posons maintenant est équivalent au suivant. Supposons qu'un ordinateur ait en mémoire des paramètres $\theta \in \mathbb{R}^p$ (nombre de planètes, leurs caractéristiques orbitales, la période de rotation de l'étoile etc.). S'il n'y avait pas de bruit, une mesure à l'instant t_k reviendrait à demander à l'ordinateur d'évaluer une fonction $f(\theta,t_k)$ . Nous connaissons t_k et (c'est l'un des modèles de la page précédente), mais nous ne connaissons pas $\theta$ . Notre but est de le déterminer à partir des mesures $y_k = f (\theta,t_k)$ . Ce principe est similaire à la résolution d'une énigme: quelqu'un connaît une information $\theta$ et nous essayons de la deviner en posant une question. Dans le cas sans bruit, celui qui pose l'énigme ne nous induit pas en erreur, mais cela ne veut pas dire que la résolution est facile !

Exemple: on veut trouver les paramètres d'une fonction affine du temps f(a,b,t) = at+c (ici $\theta = \left( \begin{array}{c} a \\c \end{array} \right)$ ). On évalue la valeur de en t_1 , on obtient y_1 = at_1+c : on a deux inconnues pour une équation, on ne peut pas résoudre. Si on a y_2 = at_2 + c , avec $t_2 \neq t_1$ , alors on a $a = \frac{y_2-y_1}{t_2-t_1}$ et c = y_1 - at_1 = y_2 - at_2 .

Notre cas est plus compliqué. Les valeurs que nous obtenons sont $y_k = f(\theta, t_k) + b_k$ où b_k est la réalisation d'une variable aléatoire $\epsilon_k$ . A l'appel numéro du programme, l'ordinateur fait appel à un autre programme qui délivre une variable alétoire selon une certaine loi. En l'occurrence, nous supposons cette loi gaussienne. Si nous faisons mesures, tout se passe comme si un programme principal évaluait la fonction $f(\theta, t_k)$ et programmes secondaires $\epsilon_k$ , k=1..m retournent chacun une valeur b_k . Si nous pouvions remonter le temps et faire les mesures plusieurs fois aux mêmes instants, on aurait des des vecteurs de mesures $y_1 = \left( \begin{array}{c} y_1 \\y_2 \\ \\y_m \end{array} \right) = f(\theta,T) + b_1$ , puis $y_2 = f(\theta,T) + b_2$ etc. (notez qu'ici b_1 et b_2 sont des vecteurs.

Comme nous ne connaissons que la loi suivie par les variables $\epsilon_k$ , ils nous est impossible de connaître les paramètres avec certitude. On leur attache une "erreur", qui quantifie l'incertitude que l'on a sur eux. En reprenant l'exemple précédent on mesure y_1 = at_1 + c +b_1 et y_2 = at_2 + c +b_2 . Si on estime et avec les mêmes formules, on fera une erreur $\frac{\sqrt{\sigma_1^2 + \sigma_2^2}}{t_1-t_2}$ sur et une erreur $\sqrt{\left(\frac{t_1}{t_2-t_1} \sigma_2\left)^2 + \right((1-\frac{t_1}{t_2-t_1}) \sigma_1 \left)^2}$ sur (admis).

Bruits

Dans le modèle $f(\theta,T) + \epsilon$ , le symbole $\epsilon$ désigne un bruit gaussien. Comme ils apparaissent constamment en détection de planètes extrasolaires et ailleurs, nous allons en donner quelques propriétés.

A une expérience donnée, $\epsilon$ prendra une valeur imprévisible. La probabilité que la valeur de soit comprise entre b_1 et b_2 est $\int_{b_1}^{b_2} f(x) dx$ où f(x) est la densité de probabilité de $\epsilon$ . Dire que $\epsilon$ est un bruit gaussien veut dire que sa densité est de la forme $f(x) = \frac{1}{\sqrt{2 \pi} \sigma} e^{-\frac{(x-\mu)^2}{\sigma^2}}$ où $\mu$ et $\sigma$ sont des réels, qui sont égaux respectivement à la moyenne et à l'écart-type de $\epsilon$ . On note souvent $\epsilon \sim g(\mu, \sigma^2)$ , qui signifie " $\epsilon$ suit une loi gaussienne de moyenne $\mu$ et de variance $\sigma^2$ . Dans la plupart des cas, le bruit est de moyenne nulle (c'est le cas ici).

Dans le modèle, des bruits d'origines différentes s'additionnent. Sachant que le résidu de l'activité stellaire que nous n'avons pas ajusté $\epsilon_S$ et le bruit de mesure $\epsilon_{mes}$ suiven une certaine loi, quelle loi suivra $\epsilon = \epsilon_{S} + \epsilon_{mes}$ ? Nous pouvons déjà dire que la moyenne de $\epsilon$ sera égale à la somme des moyennes de $\epsilon_S$ et $\epsilon_{mes}$ car l'espérance est un opérateur linéaire. Peut-on dire plus ? Si ces bruits dépendaient l'un de l'autre, la réponse pourrait être complexe. En l'occurrence, la physique de l'étoile cible et les erreurs instrumentales sont totalement indépendantes. On peut montrer que dans ces conditions, la variance de $\epsilon$ est égale à la somme des variances de $\epsilon_S$ et $\epsilon_{mes}$ . Nous pouvons même aller plus loin car la somme de deux variables gaussiennes indépendante est une variable gaussienne. En résumé, $\epsilon \sim g\left(\mu_{S} +\mu_{mes}, \sigma_{S}^2+ \sigma_{mes}^2\right)$ en l'occurrence $\mu_S$ et $\mu_{mes}$ sont nulles.

Lorsqu'on dispose de plusieurs mesures, à l'expérience numéro on a un certain bruit b_k réalisation d'une variable $\epsilon_k$ de densité f_k . La plupart du temps, on fait l'hypothèse que les brutis $\epsilon_k$ sont indépendants, c'est à dire que la probabilité d'obtenir le bruit b_k à l'expérience ne dépend pas des valeurs prises aux expériences précédentes et suivantes. Lorsque ce n'est pas le cas on parle de bruits corrélés. Pour les caractériser, on utilise souvent leur densité spectrale de puissance. Un certain profil de densité spectrale correspond à une "couleur" du bruit.

A retenir: la somme de variables gaussienne indépendantes $\epsilon = \sum\limits_{k=1}^n \epsilon_k$ où $\epsilon_k \sim g\left( \mu_k, \sigma_k^2)$ est une variable gaussienne suivant la loi $g \left( \sum\limits_{k=1}^n \mu_k , \frac{1}{n} \sum\limits_{k=1}^n \sigma_k^2 \right)$ .

Signification statistique

Il esiste plusieurs outils pour s'assurer qu'une détection possible n'est pas due au bruit. L'un des plus utilisé est le test de signification (significance en anglais), qui consiste à calculer la probabilité d'avoir le signal observé "au moins aussi grand" s'il n'y avait en réalité que du bruit. Par exemple, qupposons que l'on veuille mesurer une quantité qui est perturbée par un bruit gaussien additif de moyenne nulle et d'écart-type $\sigma=1$ , donnant une mesure y = a+b .

Si il n'y avait en réalité pas de signal ( a=0 ), les mesures seraient uniquement dues au bruit. On imagine deux cas de figures:

1er cas: on mesure . La probabilité d'avoir obtenue une valeur au moins aussi grande par hasard, c'est à dire la probabilité que $|y| \geq 0.1$ sachant que est $\int_{-\infty}^{-0.1} e^{-x^2}dx + \int_{0.1}^{\infty} e^{-x^2}dx = 2 \int_{0.1}^{\infty} e^{-x^2}dx \approx 0.818$
2ème cas: on mesure . La probabilité pour que alors que est $\int_{-\infty}^{-10} e^{-x^2}dx + \int_{10}^{\infty} e^{-x^2}dx = 2 \int_{10}^{\infty} e^{-x^2}dx \approx 2.06 \times 10^{-9}$

Ces valeurs ont une interprétation: si on réalisait exactement le même type de mesure alors qu'il n'y a pas de signal, on observerait $|y| \geq 0.1$ dans 81.8 % des cas et $|y| \geq 10$ dans $2.06 \times 10^{-7}$ % des cas.Dans le premier cas, la probabilité d'avoir un signal aussi grand que celui que l'on a mesuré est grande. On ne peut pas assurer qu'un signal a été détecté. Par contre, dans le deuxième cas on serait dans un des deux cas sur un milliard où le signal serait dû au bruit. On peut alors dire qu'on a détecté un signal à $10 \;\sigma$ , car la valeur de l'écart-type du bruit est , sa moyenne est , donc on a $\frac{y}{\sigma-0} = \frac{10}{1} = 10$ .

Détecter un signal "à 10 sigmas" est un luxe que l'on peut rarement se payer. Les détections sont annoncées plutôt pour des valeurs de sigmas.

Remarque importante: on calcule la probabilité d'avoir les observations sachant qu'il n'y a pas de signal et non la probabilité d'avoir un signal sachant les observations qui est une quantité qui a davantage de sens. Le calcul de cette dernière quantité se fait dans le cadre du calcul bayésien, outil très puissant qui ne sera pas développé dans ce cours.

Exercices

Auteur: Nathan Hara

Espérance et variance de la moyenne empirique

Difficulté : ☆☆☆

Question 1)

On rappelle que l'espérance et la variance d'une variable alétoire de densité de probabilité sont données par $\mathbb{E}(X) = \int_{-\infty}^{+\infty} xf(x)dx$ et $\text{Var}(X) = \mathbb{E}\left((X-\mathbb{E}(X))^2\right)= \int_{-\infty}^{+\infty} (x-\mathbb{E}(X))^2f(x)dx$

Soit une variable aléatoire et un réel. Montrer que $\mathbb{E}(aX) = a \mathbb{E}(X)$ et $\text{Var}(aX) = a^2 \text{Var}(X)$
Soient et deux variables aléatoires indépendantes. Montrer que $\mathbb{E}(X+Y) = \mathbb{E}(X) + \mathbb{E}(Y)$ .
On définit la covariance par $\text{Cov}(X,Y) = \mathbb{E}\{ (X-\mathbb{E}(X)) (Y-\mathbb{E}(Y))\}$ . Montrer que $\text{Var}(X+Y) = \text{Var}(X) + \text{Var}(Y) + 2 \text{Cov}(X,Y)$
Nous allons maintenant voir un cas où la précision d'un estimateur est facilement calculable. On considère mesures entachées de bruits gaussiens d'une quantité fixe . Plus précisément, la mesure numéro est modélisée par une variable aléatoire $X_j = K + \epsilon_j$ où $\epsilon_j$ est un bruit gaussien de moyenne nulle et de variance $\sigma^2$ . On suppose que les mesures sont indépendantes, ce qui implique en particulier que $\text{Cov}(\epsilon_j,\epsilon_i) = 0$ pour $i \neq j$ . Pour obtenir une estimation de , on fait la moyenne empirique des expériences, c'est à dire $M = \frac{X_1 + X_2 +... X_m}{m}$ . Montrer que $\mathbb{E}(M) = K$ et $\text{Var}(M) = \frac{\sigma^2}{m}$
La précision de l'estimateur augmente-t-elle avec le nombre de mesures ?

Loi du chi 2

Nous nous sommes toujours ramenés à des modèles du type: modèle déterministe + bruit gaussien. Afin de vérifier que les observations sont compatibles avec le modèle, on étudie les résidus, définis comme "les observation - le modèle ajusté".

La loi du $\chi^2$ est un outil commode pour étudier le comportement de plusieurs variables gaussiennes. Considérons d'abord une famille de variables aléatoires gaussiennes indépendantes $(\epsilon_k)_{k=1..m}$ , de moyenne nulle et de variance unité. On forme la quantité $\chi^2_m= \epsilon_1^2 + \epsilon_2^2 + ... + \epsilon_m^2$ . Comme les $\epsilon_k$ sont des variables aléatoires, les $\epsilon_k^2$ le sont aussi. La somme de variables alétoires étant toujours une variable alétoire, $\chi^2_m$ suit une certaine loi de probabilité. Dans l'analogie avec un programme informatique, la variable alétoire $\chi^2_m$ se comporte comme un programme qui appelle programmes générant une variables gaussiennes, puis additionne leurs carrés. Elle est appelée loi du $\chi^2$ à m degrés de liberté. On peut montrer qu'en moyenne une variable gaussienne au carré $\epsilon^2$ a une moyenne de 1. En conséquence, $\chi^2_m$ vaudra typiquement .

Pourquoi cette loi serait utile pour notre cas ? Si le modèle est bien ajusté, les résidus doivent se comporter comme un bruit gaussien. En supposant que les erreurs sont toutes indépendantes, de moyenne nulle et de variance unité, les résidus $r_k = y-f(\theta,t_k})$ en sont une réalisation. Donc $R^2=\sum\limits_{k=1}^m r_k^2$ est une réalisation d'une loi du $\chi^2$ à degrés de liberté. Si R^2 est de l'ordre de , le modèle est cohérent. Sinon, le modèle ou les paramètres ajustés sont à revoir.

En pratique, les erreurs $\epsilon__k$ ne sont évidemment pas de variance unité et parfois pas indépendantes. Par contre on peut à bon droit supposer qu'elles sont de moyenne nulle. Pour se ramener au cas précédent, on calcule non pas une réalisation de $\epsilon_1^2+...+\epsilon_m^2$ mais de $\chi^2_m'=\epsilon^T \Sigma^{-1} \epsilon$ où $\epsilon = \left( \begin{array}{c} \epsilon_1 \\ \epsilon_2 \\ \\ \epsilon_m \end{array} \right)$ , $\epsilon^T$ est sa transposée et $\Sigma$ est la matrice des variances-covariances de $\epsilon$ . Dans le cas où les $\epsilon_k$ sont indépendantes, la matrice des variances-covariances est diagonale, son -ème terme diagonal étant $\frac{1}{\sigma_k^2}$ , soit l'inverse de la variance de $\epsilon_k$ .

La vraisemblance

Etant donné des paramètres $\theta$ , le modèle global $f(\theta,T)+\epsilon$ est une variable aléatoire: il est somme d'une variable aléatoire valant $f(\theta,T)$ avec une probabilité 1 et d'un vecteur de variables aléatoires gaussienne $\epsilon$ . A ce titre, il a une certaine densité de probabilité que l'on note $L(y|\theta)$ . Le symbole | se lisant "sachant". La lettre L vient de Likelihood, qui veut dire vraisemblance en anglais. Il s'agit dans l'idée de la probabilité d'obtenir $y = f(\theta,T)+\epsilon$ pour une valeur de $\theta$ donnée..

La fonction $\theta \rightarrow L(y|\theta)$ est souvent appelée "fonction de vraisemblance". La valeur de $\theta$ maximisant $L(y|\theta)$ est appelé l'estimateur du maximum de vraisemblance. Il a de bonnes propritétés statistiques. En effet, on peut montrer que c'est un estimateur:

Non biaisé: lorsque le nombre de mesure tend vers l'infini, l'estimateur du maximum de vraisemblance tend vers la bonne valeur $\theta$
De variance minimale: parmi les estimateurs non biaisés, il a la variance la plus faible que l'on puisse espérer.

Dans notre cas, si les $\epsilon_k$ sont des variables indépendantes, leur densité de probabilité jointe est égale au produit de leurs densité de probabilité. $L(y|\theta) = \prod\limits_{k=1}^m g_k(y-f(\theta,t_k))$ où g_k est la densité de probabilité de la variable $\epsilon_k$ . De plus, si ces varibles sont gaussiennes et indépendantes, on a: $L(y|\theta) = \prod\limits_{k=1}^m \frac{1}{\sqrt{2\pi}\sigma_k} e^{-\frac{(y-f(\theta,t_k))^2}{\sigma_k^2}} = \frac{1}{(\sqrt{2\pi})^m \prod\limits_{k=1}^m \sigma_k} e^{-\sum\limits_{k=1}^m \frac{(y-f(\theta,t_k))^2}{\sigma_k^2}}}$

Méthode des moindres carrés

Pour des bruits gaussiens indépendants, maximiser la vraisemblance, revient à minimiser $\theta \rightarrow F(y,\theta)=-ln(L(y|\theta)) =\sum\limits_{k=1}^m \frac{(y-f(\theta,t_k))^2}{\sigma_k^2}}$ puisque $x \rightarrow e^{-x}$ est une fonction décroissante. La méthode consistant à minimiser $F(y,\theta)$ s'appelle la méthode des moindres carrés. C'est la méthode d'estimation de loin la plus utilisée dans tous les domaines. Elle est parfois utilisée quand les bruits ne sont pas gaussiens, mais il faut garder à l'esprit qu'elle n'a alors plus de propriétés statistiques sympathiques (sauf quand le modèle est linéaire en $\theta)$ .

Dans notre cas, les paramètres $\theta$ sont les éléments des orbites, les paramètres du bruit stellaire, du mouvement propre, etc. La fonction $F(y,\theta)$ a donc de nombreux paramètres, et trouver son minimum global est une tâche ardue qui fait l'objet d'une littérature très vaste.

Lorsque le modèle est linéaire en $\theta$ i. e. $f(T,\theta) = A \theta$ où est une certaine matrice dépendant des instants d'observation $T = (t_k)_{k=1..m}$ , l'ajustement est beaucoup plus simple car il a une solution explicite (voir mini-projet). On essaye de se ramener autant que possible à des ajustements linéaires. La plupart du temps, on estime les paramètre les uns après les autres, puis un ajustement global est réalisé. Une démarche classique consiste à:

Ajuster le mouvement propre et la parallaxe, puis soustraire ces ajustements du signal initial
Chercher les périodes éventuelles dans le signal à l'aide du périogogramme (voir page suivante)
Si on trouve une période, on ajuste une orbite, que l'on soustrait du signal
On recommence les étapes 2 et 3 jusqu'à ce que le périodogramme ne présente plus de pic significatif

Sur la figure, on représente les étapes d'un ajustement d'un signal astrométrique simulé (2 planètes, 45 observations). De gauche à droite et de haut en bas:

Signal brut en rouge En vert signal après correction des bruits instrumentaux et atmosphériques . Le mouvement propre (tendance linéaire) et la parallaxe (oscillation) dominent très nettement le signal.
En vert: le signal planétaire réel, sans aucun bruit $(s^{\star})$ . En rouge: signal estimé en soustrayant le mouvement propre, l'accélération de perspective et la parallaxe au signal .
En rouge: signal . En vert : orbite d'une planète ajustée sur .
En rouge : signal dont on a soustrait . En vert :orbite d'une planète ajustée sur .
En rouge: signal . En vert: , ajustement de deux planètes en partant des estimations et .
En vert, signal $(s^{\star})$ . On ajuste tous les paramètres simultanément sur en partant des estimations précédentes. Le signal en rouge représente le signal planétaire de cet ajustement.

Exemple d'ajustement en astrométrie (données simulées)

Six étapes d'ajustement aux données. Pour représenter la chronologie des mesures, on trace des segments reliant deux points de mesure consécutifs.

Périodogramme

Evaluer les résidus sur une grille d'un modèle à paramètres, où chacun d'eux peut prendre valeurs recquiert n^p évaluations, ce qui devient rapidement ingérable numériquement. Les planètes ont un mouvement périodique, donc il est raisonnable de checher des signaux périodiques dans le signal en ne faisant varier que la période du signal recherché. Pour des signaux échantillonnés à intervalles réguliers, on utilise la transformée de Fourier. Le périodogramme est un moyen de checher des signaux périodiques dans des données échantillonnées irrégulièrement. On les notera . Le périodogramme de Lomb-Scargle d'un signal $(y_k)_{k=1..m}$ échantillonné aux instants $(t_j)_{j=1..m}$ est défini comme suit pour une fréquence quelconque $\omega$ :

$P_y(\omega) = \frac{1}{2} \left( \frac{ \left( \sum\limits_{j=1}^m y_j \cos \omega(t_j-\tau) \right)^2}{ \sum\limits_{j=1}^m \cos^2 \omega(t_j-\tau)}} + \frac{ \left( \sum\limits_{j=1}^m y_j \sin \omega(t_j-\tau) \right)^2}{ \sum\limits_{j=1}^m \sin^2 \omega(t_j-\tau)}} \right)$ où $\tau$ vérifie:

$\tan 2 \omega \tau = \frac{\sum\limits_{j=1}^m \sin 2\omega t_j}{\sum\limits_{j=1}^m \cos 2\omega t_j}$

Cette expression est équivalente à $P_y(\omega) = A^2 +B^2$ où et sont les paramètres minimisant $\sum\limits_{j=1}^m (y_j-A\cos(\omega t_j) - B\sin(\omega t_j))^2$ . Le modèle $A,B, t \rightarrow A \cos \omega t + B \sin \omega t$ est linéaire en A, B , on a donc une solution explicite à la minimisation.

Le périodogramme a une propriété très intéressante: si le signal d'entrée est un bruit gaussien $\epsilon$ de variance unité, p_0 une valeur réelle fixée et $\omega$ une fréquence quelconque, la probabilité que $P_{\epsilon}(\omega)$ dépasse p_0 est $Pr\{P_\epsilon(\omega)>p_0\} = e^{-p_0}$ . En d'autres termes, la probabilité qu'une valeur du périodogramme à $\omega$ fixée soit "au moins aussi grand que p_0 " par hasard décroît exponentiellement. Supposons que l'on ait un signal y(t) = x(t) + b(t) où b(t) est un bruit gaussien de variance unité et nous trouvons un pic de taille p_0 , on calcule la probabilité de trouver un pic au moins aussi grand si le signal n'est composé que de bruit: $e^{-p_0}$ . Si cette valeur est petite, on pourra confirmer la détection d'un signal avec une erreur de fausse alarme de $e^{-p_0}$ . Ce procédé n'est autre qu'un test de signification statistique.

La figure montre un exemple de périodogramme. Il s'agit d'un périodogramme d'une des coordonnées d'un signal astrométrique simulé dont on a soustrait le mouvement propre et la parallaxe. En bleu, on représente un périodogramme idéal, sans bruit, avec 10000 observations. Le périodogramme représenté en rouge est lui calculé pour 45 observations. Le pic le plus haut correspond bien à une fréquence réelle. Par contre, le deuxième pic le plus important (à 0.37 rad/s) ne correspond pas à une sinusoïde. C'est ce qu'on appelle un alias de la fréquence principale.

En pratique, la variance du bruit n'est pas unitaire et dépend de l'instant de mesure. On peut corriger ce problème en minimisant un critère pondéré $\sum\limits_{j=1}^m \frac{(y_j-A\cos(\omega t_j) - B\sin(\omega t_j))^2}{\sigma_j^2}$ où $\sigma_j$ est la variance du bruit à la mesure .

Exemple de périodogramme

En rouge: périodogramme d'un système simulé à deux planètes après soustraction du mouvement propre et de la parallaxe (signal très peu bruité), pour 45 mesures. En bleu: périodogramme du même système sans bruit et avec 10000 observations.

Se tester

Questions qualitatives

Expliquer avec la notion de signal sur bruit qu'avoir plus de mesures rend la détection plus précise.
Quels avantages et inconvénients a-t-on à utiliser un télescope au sol plus large dans le cas des vitesses radiales ? De l'astrométrie ?
La vitesse radiale par rapport à l'observateur ne dépend pas de la distance de l'étoile. Expliquer pourquoi l'observation d'un système planétaire donné est plus précise quand celui-ci est plus proche de l'observateur.
Pourquoi, à orbite fixée, les planètes plus massives sont plus faciles à détecter ?
Pourquoi, à masse fixée, les orbites plus proches de l'étoile sont plus faciles à détecter ?
Pourquoi, pour une planète donnée et un type d'étoile donné, plus l'étoile est massive plus la détection est difficile ?
Vous observez le spectre d'une étoile dont les oscillations radiales ont une période de 5 minutes, pourrait-on réduire l'influence de ce phénomène en choisissant un certain temps d'intégration ?
Voyez vous des biais détection possibles ?

Mini projet

But

L'objectif de ce mini projet est de réaliser un code permettant de calculer un périodogramme. On commence par établir la formule des moindres carrés dans le cas d'un modèle linéaire.

Moindres carrés linéaires

On se donne des observations et on considère une matrice réelle de taille m,n avec m>n et de rang . On veut calculer $\hat{\theta}$ , minimisant $F(\theta)=\sum\limits_{k=1}^m \left(\frac{y_k - A \theta}{\sigma_k}\right)^2$ . On pose matrice diagonale dont les éléments diagonaux valent $\frac{1}{\sigma_k}$ et $\Sigma = W^2$ . En termes vectoriels, $F(\theta) = \|W(AX-y) \|^2$ . Comme on veut trouver le minimum de $F(\theta)$ et que cette fonction est différentiable, ce minimum vérifie $dF(\theta)$ où est la différentielle de . Montrer que la solution de $dF(\theta) = 0$ est $\hat{\theta} = (A^T \Sigma A)^{-1} A^T \Sigma$ où désigne la transposition. Montrer que c'est un minimum global. (On pourra utiliser le fait que $\|W(AX-y) \|^2 = (y-AX)^T\Sigma (y-AX)$ .

Périodogramme

Nous voulons trouver la sinusoïde qui a la distance minimale aux données au sens des moindres carrés. En d'autres termes, on veut ajuster $\theta_{\omega} = \left( \begin{array}{c} a_{\omega} \\ b_{\omega} \end{array} \right)$ $F_\omega(\theta)=\sum\limits_{k=1}^m \left(\frac{y_k - a_{\omega} \cos \omega t_k - b_{\omega} \sin \omega t_k}{\sigma_k}\right)^2$ où les t_k sont des instants d'observation. Calculer avec la formule précédente $\theta$ minimisant $F_\omega(\theta)$ . On notera $P(\omega) = a_{\omega}^2 + b_{\omega}^2$

Implémentation

Dans le langage informatique de votre choix, écrire un programme permettant de:

Lire des données provenant d'un fichier texte et les mettre dans un tableau, en particulier extraire les instants d'observation, vitesses radiales et erreurs associées.
Calculer $P(\omega)$ pour une fréquence $\omega$ donnée
Tracer $P(\omega)$ pour $\omega$ allant de 0 à $4 \pi$
Sur ce lien ou ce lien , copier coller le tableau de données dans un fichier texte et en faire l'analyse.