#!/usr/bin/env python3
# -*- coding: utf-8 -*-
# ====================================================================
# Chapitre 13 - Exercice 46                                   STANDARD
# ====================================================================
"""
Programme permettant de déterminer les paramètres orbitaux de l'orbite
de la comète de Halley à partir de données astronomiques.
"""
print('\nParamètres orbitaux de la comète de Halley')

# ====================================================================
# Partie 1: Exploitation des passages de la comète de Halley au périhélie
# ====================================================================
"""
Cette partie permet d'estimer le périhélie q_H, l'excentricité e_H et
la période de révolution T_H de la comète de Halley, ainsi que leurs
incertitudes-types u(q_H), u(e_H) et u(T_H) associées, à partir du
traitement statistique des données du fichier C13_Ex46_Passage_Halley.txt.
"""
# Extraction des données : NE PAS MODIFIER ===========================
delta_t_Halley,q_Halley,e_Halley,T_Halley=[],[],[],[] 
with open('C13_Ex46_Passage_Halley.txt') as fichier1:
    header=[fichier1.readline() for i in range(4)]
    line = fichier1.readline()
    while line:
        data1,data2,data3,data4,data5=line.split(';')
        delta_t_Halley.append(float(data2)) # ∆t (en jour)
        q_Halley.append(float(data3)) # Périhélie (en ua)
        e_Halley.append(float(data4)) # Excentricité (sans unité)
        T_Halley.append(float(data5)) # Période (en année)
        line = fichier1.readline()

# CODE A COMPLETER ===================================================

# Traitement statistique des données
import numpy as np
n=...A compléter...     # Nombre de valeurs de la série de données

q_H=...A compléter...   # q_H (en ua) <- valeur moyenne de la série de données
u_q_H=...A compléter... # u(q_H) (en ua) <- incertitude-type associée

e_H=...A compléter...   # e_H <- valeur moyenne de la série de données
u_e_H=...A compléter... # u(e_H) <- incertitude-type associée

T_H=...A compléter...   # T_H (en année) <- valeur moyenne de la série de données
u_T_H=...A compléter... # u(T_H) (en année) <- incertitude-type associée

# Affichage du résultat
print('\nTraitement statistique des passages de la comète de Halley')

print('\nPérihélie q_H =',...A compléter...) # Affichage de q_H
print('Incertitude-type  u(q_H) =',...A compléter...) # Affichage de u(q_H)

print('\nExcentricité e_H =',...A compléter...) # Affichage de e_H
print('Incertitude-type  u(e_H) =',...A compléter...) # Affichage de u(e_H)

print('\nPériode T_H =',...A compléter...) # Affichage de T_H
print('Incertitude-type  u(T_H) =',...A compléter...) # Affichage de u(T_H)

# ====================================================================
# Partie 2: Test de la 3e loi de Kepler sur les comètes de type Halley
# ====================================================================
"""
Cette partie permet de tester la 3e loi de Kepler sur les comètes de
type Halley et de déterminer la constante k intervenant dans cette loi
et son incertitude-type u(k) à partir des données du fichier
C13_Ex46_Cometes_Type_Halley.txt regroupant les paramètres orbitaux de
ces comètes à leur dernier passage à leur périhélie.
"""
# Extraction des données : NE PAS MODIFIER ===========================
q,e,T=[],[],[] # Paramètres orbitaux des comètes de type Halley
               # à leur dernier passage au périhélie
with open('C13_Ex46_Cometes_Type_Halley.txt') as fichier2:
    header=[fichier2.readline() for i in range(2)]
    line = fichier2.readline()
    while line:
        data0,data1,data2,data3,data4=line.split('\t')
        q.append(float(data2))  # Périhélie q (en ua)
        e.append(float(data3))  # Excentricité e (sans unité) 
        T.append(float(data4))  # Période T (en année)
        line = fichier2.readline()

# CODE A COMPLETER ===================================================
"""
N=...A compléter...       # Nombre de valeurs de la série de données

# Listes des grandeurs à calculer pour tester la 3e loi de Kepler
T_carre=[...A compléter... for i in range(N)] # Carré des périodes (en année^2)
a_cube=...A compléter...  # Cube des demi-grands axes (en ua^3)

# Représentation graphique permettant de tester la 3e loi de Kepler
from matplotlib import pyplot as plt
plt.plot(...A compléter...,'+',label='Données') 

# Regression linéaire et modélisation de la forme y = k*x+b
from scipy.stats import linregress
k,b,r,pval,u_k=...A compléter... # Détermination de k, b et u(k) 
x=np.linspace(...A compléter...) # Domaine des abscisses pour la modélisation
plt.plot(...A compléter...,':',label ='Modélisation') # Tracé de la modélisation

# Mise en forme de la fenêtre graphique
...A compléter... # Titre du graphe
...A compléter... # Label de l'axe des abscisses
...A compléter... # Label de l'axe des ordonnées
...A compléter... # Limites de l'axe des abscisses
...A compléter... # Limites de l'axe des ordonnées
...A compléter... # Affiche une grille en pointillés
...A compléter... # Affiche la légende
...A compléter... # Affiche la figure
"""

# ====================================================================
# Partie 3: Période de révolution T_H de la comète de Halley lors de
# son dernier passage au périhélie le 9 février 1986
# ====================================================================
"""
Cette partie exploite les résultats de la partie 2 et les mesures issues
des observations astronomiques du passage de la comète de Halley à son
périhélie le 9 février 1986. T_H et u(T_H) sont déterminées par
simulation d'un processus aléatoire tenant compte de la composition des
incertitudes-types.
"""

# CODE A COMPLETER ===================================================
"""
# Création des listes L=[Valeur,Incertitude-type] pour e_H, q_H et k
L_e_H=[...A compléter...] # Liste e_H pour le 9 février 1986
L_q_H=[...A compléter...] # Liste q_H (en ua) pour le 9 février 1986
L_k=[...A compléter...]   # Liste k (en année^2/ua^3)

# Tirage aléatoire d'une valeur selon la loi normale
def Alea(L):
    return np.random.normal(L[0],L[1])

# Simulation d'une distribution d pour la période T_H
d=[]
Iteration=100000
for i in range(Iteration):
    Alea_T_H= ...A compléter... # Valeur aléatoire de T_H
    d.append(Alea_T_H)

# Calcul de T_H_Kepler et u(T_H_Kepler)
T_H_Kepler=...A compléter...   # T_H_Kepler <- Valeur moyenne de la distribution d
u_T_H_Kepler=...A compléter...  # u(T_H_Kepler) <- Ecart-type de la distribution d

# Affichage du résultat
print('\nCalcul de la période de révolution T_H de la comète de Halley' \
      ' lors de son dernier passage au périhélie le 9 février 1986' \
      ' par application de la 3e loi de Kepler :')
print('\nPériode T_H =',...A compléter...)            # Affichage de T_H_Kepler 
print('Incertitude-type u(T_H) =', ...A compléter...) # Affichage de u(T_H_Kepler)
"""