#!/usr/bin/env python3

# © Concours Centrale-Supélec 2017
#
# Ce document est diffusé sous la licence Créative Commons Paternité-Pas
# d'Utilisation Commerciale-Partage des Conditions Initiales à l'Identique.
# http://creativecommons.org/licenses/by-nc-sa/2.0/fr/
#
# Cela signifie que vous pouvez librement le reproduire, distribuer,
# communiquer ou modifier à condition de :
# - mentionner qu'il s'agit d'un sujet du Concours Centrale-Supélec,
# - ne pas en faire d'usage commercial,
# - distribuer sous une licence similaire toute version transformée,
#   modifiée ou adaptée.

## Importations

import matplotlib.pyplot as plt
import numpy as np

## Valeurs numériques

N = 1.5
c = 3.0E8
omega0 = 3.14E15
E0 = 1   # éclairement incident en unités arbitraires

## Calcul des éclairements E1 et E2 des deux ondes réfléchies

def E1(n):
    valeur  = (((1-n)/(1+n))**2)*E0
    return valeur

def E2(n):
    valeur  = ( (2/(1+n)*(n-N)/(n+N)*2*n/(1+N))**2 )*E0
    return valeur

n_liste = np.linspace(1,N,1000)
E1_liste = [ E1(n) for n in n_liste ]
E2_liste = [ E2(n) for n in n_liste ]
Er_liste = [ (np.sqrt(E1(n)) - np.sqrt(E2(n)))**2 for n in n_liste ]

plt.close()  # ferme toutes les figures précédemment ouvertes
plt.figure(1)
plt.plot(n_liste,E1_liste,color = "blue")
plt.plot(n_liste,E2_liste, color = "red")
plt.plot(n_liste,Er_liste, color = "green")
plt.grid()
plt.xlabel("Indice n")
plt.ylabel("E1,E2 et Er")
plt.show()

## À adapter ou à justifier

n = 1.21  # obtenu avec le programme précédent
J0 = 1    # densité spectrale en unité arbitraire
Delta_omega = 2.35E15

alpha1 = ( (1-n)/(n+1) )**2
alpha2 = ( 2/(n+1)*(n-N)/(n+N)*2*n/(n+1) )**2


def delta(e):
    return(2*n*e)

def sinc(x) :
    if x == 0 :
        return 1
    else :
       return np.sin(x)/x

def R(e):
    valeur = alpha1+alpha2+2*np.sqrt(alpha1*alpha2)*sinc(Delta_omega*delta(e)/(2*c))*np.cos(omega0*delta(e)/c)
    return valeur

plt.figure(2)
e_liste = np.linspace(0,500E-9,100)
R_liste = [ R(e) for e in e_liste ]
plt.plot(e_liste,R_liste)
plt.xlabel("Epaisseur")
plt.ylabel("Coefficient de réflexion")
plt.grid()
plt.show()

## Etude de R en fonction de i pour e et n fixés

e = 119E-9

def delta_oblique(i):
    resultat = 2*n*e*np.sqrt(1-np.sin(i)*np.sin(i)/(n**2))
    return(resultat)

def R_oblique(i):
    valeur = alpha1+alpha2+2*np.sqrt(alpha1*alpha2)*sinc(Delta_omega*delta_oblique(i)/(2*c))*np.cos(omega0*delta_oblique(i)/c)
    return valeur


plt.figure(3)
i_liste = np.linspace(0,np.pi/2,100)
R_oblique_liste = [ R_oblique(i) for i in i_liste ]
plt.plot(i_liste,R_oblique_liste)
plt.xlabel("Angle d'incidence")
plt.ylabel("Coefficient de réflexion")
plt.grid()
plt.show()