15- calcul couches minces

Code python:  illustration

from numpy import *

from matplotlib.pyplot import *

import tkinter as tk

from tkinter import ttk

def calcul_empilement_couches_minces(materiau_H, materiau_L, substrat, lambda_0, empilement, mode, intervalle, pas, face_arriere=False):

    """

    Calcule la réflectance et la transmittance d'un empilement de couches minces 

    pour les polarisations s et p. Affiche également l'empilement des couches.

    Args:

        materiau_H: Indice de réfraction du matériau à indice élevé (H).

        materiau_L: Indice de réfraction du matériau à indice faible (L).

        substrat: Indice de réfraction du substrat.

        lambda_0: Longueur d'onde de centrage (nm).

        empilement: Liste des épaisseurs optiques des quarts d'onde (QWOT) pour chaque couche.

        mode: 'spectral' ou 'angulaire'.

        intervalle: Tuple (debut, fin) de l'intervalle de longueur d'onde (nm) ou d'angle d'incidence (degrés).

        pas: Pas de l'intervalle.

        face_arriere: Booléen indiquant si la face arrière est prise en compte.

    Returns:

        Un dictionnaire contenant les résultats du calcul (longueurs d'onde, 

        angles d'incidence, réflectance s, réflectance p, transmittance s, 

        transmittance p).

    """

    # 1. Collecte des informations

    N_couches = len(empilement)

    # 2. Construction de tables d'indices

    if mode == 'spectral':

        longueurs_onde = arange(intervalle[0], intervalle[1] + pas, pas)

        NH = materiau_H  # Indice constant pour le mode spectral

        NL = materiau_L  # Indice constant pour le mode spectral

        Nex = substrat  # Indice constant pour le mode spectral

        theta_inc = 0  # Angle d'incidence fixe à 0 degré (incidence normale) 

    elif mode == 'angulaire':

        longueurs_onde = array([lambda_0])  # Longueur d'onde fixe

        NH = materiau_H

        NL = materiau_L

        Nex = substrat

        theta_inc = radians(arange(intervalle[0], intervalle[1] + pas, pas))

    else:

        raise ValueError("Mode invalide. Choisissez 'spectral' ou 'angulaire'.")

    # 3. Calcul des épaisseurs mécaniques

    epaisseurs = []

    for i in range(N_couches):

        indice = NH if i % 2 == 0 else NL

        epaisseur_courante = empilement[i] * lambda_0 / (4 * real(indice))

        epaisseurs.append(epaisseur_courante)

    # 4. Calcul sans face arrière

    if mode == 'spectral':

        Rs = zeros((len(longueurs_onde), 1))

        Rp = zeros((len(longueurs_onde), 1))

        Ts = zeros((len(longueurs_onde), 1))

        Tp = zeros((len(longueurs_onde), 1))

    elif mode == 'angulaire':

        Rs = zeros((len(longueurs_onde), len(theta_inc)))

        Rp = zeros((len(longueurs_onde), len(theta_inc)))

        Ts = zeros((len(longueurs_onde), len(theta_inc)))

        Tp = zeros((len(longueurs_onde), len(theta_inc)))

    for i_lambda, longueur_onde in enumerate(longueurs_onde):

        if mode == 'spectral':

            theta_loop = [theta_inc]

        else:

            theta_loop = theta_inc

        for i_theta, theta in enumerate(theta_loop):

            Ninc = 1  # Air                                         

            Nq1 = Nex

            alpha = Ninc * sin(theta)

            # Calcul pour la polarisation s

            Ms = eye(2, dtype=complex)

            etainc = sqrt(Ninc**2 - alpha**2)

            etasub = sqrt(Nex**2 - alpha**2)

            for i in range(N_couches):

                Ni = NH if i % 2 == 0 else NL

                phi = (2 * pi / longueur_onde) * (sqrt(Ni**2 - alpha**2) * epaisseurs[i])

                eta_s = sqrt(Ni**2 - alpha**2) 

                Mi = array([[cos(phi), (1j / eta_s) * sin(phi)],

                            [1j * eta_s * sin(phi), cos(phi)]])

                Ms = Mi @ Ms

            rs = (etainc * Ms[0, 0] - etasub * Ms[1, 1] + etainc * etasub * Ms[0, 1] - Ms[1, 0]) / (etainc * Ms[0, 0] + etasub * Ms[1, 1] + etainc * etasub * Ms[0, 1] + Ms[1, 0])

            ts = 2 * etainc / (etainc * Ms[0, 0] + etasub * Ms[1, 1] + etainc * etasub * Ms[0, 1] + Ms[1, 0])

            Rs[i_lambda, i_theta] = abs(rs)**2

            Ts[i_lambda, i_theta] = (real(abs(etasub)) / real(abs(etainc))) * abs(ts)**2

            # Calcul pour la polarisation p

            Mp = eye(2, dtype=complex)

            etainc = Ninc**2/sqrt(Ninc**2 - alpha**2)

            etasub = Nex**2/sqrt(Nex**2 - alpha**2)

            for i in range(N_couches):

                Ni = NH if i % 2 == 0 else NL

                phi = (2 * pi / longueur_onde) * (sqrt(Ni**2 - alpha**2) * epaisseurs[i])

                eta_p = Ni**2 / sqrt(Ni**2 - alpha**2)

                Mi = array([[cos(phi), (1j / eta_p) * sin(phi)],

                            [1j * eta_p * sin(phi), cos(phi)]])

                Mp = Mi @ Mp

            rp = (etainc * Mp[0, 0] - etasub * Mp[1, 1] + etainc * etasub * Mp[0, 1] - Mp[1, 0]) / (etainc * Mp[0, 0] + etasub * Mp[1, 1] + etainc * etasub * Mp[0, 1] + Mp[1, 0])

            tp = 2 * etainc / (etainc * Mp[0, 0] + etasub * Mp[1, 1] + etainc * etasub * Mp[0, 1] + Mp[1, 0])

            Rp[i_lambda, i_theta] = abs(rp)**2

            Tp[i_lambda, i_theta] = (real(abs(etasub)) / real(abs(etainc))) * abs(tp)**2

    # 5. Calcul avec face arrière (optionnel)

    if face_arriere:

        # Implémentez les calculs de Rbackside, Rwith Backside et Twith Backside ici

        pass

    # 6. Vérification

    assert all(Rs + Ts <= 1.0001), "Erreur: Rs + Ts > 1."  

    assert all(Rp + Tp <= 1.0001), "Erreur: Rp + Tp > 1."  

    # Retourner les résultats

    resultats = {

        'longueurs_onde': longueurs_onde,

        'angles_incidence': degrees(theta_inc),

        'Rs': Rs,

        'Rp': Rp,

        'Ts': Ts,

        'Tp': Tp        

    }

    return resultats, epaisseurs

def afficher_parametres(materiau_H, materiau_L, substrat, lambda_0, empilement, mode, intervalle, pas, epaisseurs):

    """Crée une nouvelle fenêtre pour afficher les paramètres utilisés."""

    fenetre_parametres = tk.Toplevel()

    fenetre_parametres.title("Paramètres utilisés")

    # Créer un Treeview pour afficher les paramètres

    tree = ttk.Treeview(fenetre_parametres)

    # Définir les colonnes

    tree["columns"] = ("parametre", "valeur")

    tree.column("#0", width=0, stretch=tk.NO)  # Colonne fantôme

    tree.column("parametre", anchor=tk.W, width=200)

    tree.column("valeur", anchor=tk.W, width=200)

    # Définir les en-têtes

    tree.heading("#0", text="", anchor=tk.W)

    tree.heading("parametre", text="Paramètre", anchor=tk.W)

    tree.heading("valeur", text="Valeur", anchor=tk.W)

    # Ajouter les données

    tree.insert(parent='', index='end', iid=0, text="", values=("Indice matériau H", f"{materiau_H:.2f}"))

    tree.insert(parent='', index='end', iid=1, text="", values=("Indice matériau L", f"{materiau_L:.2f}"))

    tree.insert(parent='', index='end', iid=2, text="", values=("Indice substrat", f"{substrat:.2f}"))

    tree.insert(parent='', index='end', iid=3, text="", values=("Longueur d'onde de centrage (nm)", f"{lambda_0}"))

    tree.insert(parent='', index='end', iid=4, text="", values=("Mode", mode))

    tree.insert(parent='', index='end', iid=5, text="", values=("Intervalle", f"{intervalle}"))

    tree.insert(parent='', index='end', iid=6, text="", values=("Pas", f"{pas}"))

    # Empilement

    tree.insert(parent='', index='end', iid=7, text="", values=("Empilement (QWOT)", ""))

    for i, epaisseur_optique in enumerate(empilement):

      tree.insert(parent=7, index='end', iid=f'7-{i}', text="", values=(f"  Couche {i+1}", f"{epaisseur_optique}"))

    # Epaisseurs

    tree.insert(parent='', index='end', iid=8, text="", values=("Epaisseurs (nm)", ""))

    for i, epaisseur in enumerate(epaisseurs):

        tree.insert(parent=8, index='end', iid=f'8-{i}', text="", values=(f"  Couche {i+1}", f"{epaisseur:.2f}"))

    tree.pack()

def recuperer_parametres():

    """Récupère les paramètres depuis l'interface utilisateur."""

    mat_H = materiau_H_real.get() - 1j * materiau_H_imag.get()

    mat_L = materiau_L_real.get() - 1j * materiau_L_imag.get()

    sub = substrat.get()

    l0 = lambda_0.get()

    emp = [float(x) for x in empilement_str.get().split(',')]

    m = mode.get()

    if m == "spectral":

        inter = (intervalle_spectral_debut.get(), intervalle_spectral_fin.get())

        p = pas_spectral.get()

    else:

        inter = (intervalle_angulaire_debut.get(), intervalle_angulaire_fin.get())

        p = pas_angulaire.get()

    return mat_H, mat_L, sub, l0, emp, m, inter, p

def tracer_graphiques(resultats_spectral, resultats_angulaire, epaisseurs, materiau_H_real, materiau_H_imag, materiau_L_real, materiau_L_imag, substrat, empilement_str):

    """Trace les graphiques de réflectance, transmittance et l'empilement des couches."""

    fig, axes = subplots(2, 2, figsize=(10, 8))  # 2x2 grid for 4 plots

    # Tracé spectral

    axes[0, 0].plot(resultats_spectral['longueurs_onde'], resultats_spectral['Rs'][:, 0], label='Rs')

    axes[0, 0].plot(resultats_spectral['longueurs_onde'], resultats_spectral['Ts'][:, 0], label='Ts')

    axes[0, 0].set_xlabel('Longueur d onde (nm)')

    axes[0, 0].set_ylabel('Reflectance / Transmittance')

    axes[0, 0].set_title('Tracé spectral (incidence normale)')

    axes[0, 0].legend()

    # Tracé angulaire

    axes[0, 1].plot(resultats_angulaire['angles_incidence'], resultats_angulaire['Rs'][0, :], label='Rs', linestyle='--')

    axes[0, 1].plot(resultats_angulaire['angles_incidence'], resultats_angulaire['Rp'][0, :], label='Rp', linestyle='--')

    axes[0, 1].plot(resultats_angulaire['angles_incidence'], resultats_angulaire['Ts'][0, :], label='Ts')

    axes[0, 1].plot(resultats_angulaire['angles_incidence'], resultats_angulaire['Tp'][0, :], label='Tp')

    axes[0, 1].set_xlabel("Angle d incidence (degrés)")

    axes[0, 1].set_ylabel('Reflectance / Transmittance')

    axes[0, 1].set_title('Tracé angulaire')

    axes[0, 1].legend()

    # --- Graphique : Partie réelle de l'indice en fonction de l'épaisseur ---

    indices = [materiau_H_real - 1j* materiau_H_imag if i % 2 == 0 else materiau_L_real - 1j* materiau_L_imag for i in range(len(empilement_str.split(',')))]

    epaisseur_cumulee = cumsum(epaisseurs)

    partie_reelle_indices = [real(n) for n in indices[::-1]]

    # Ajout des points fictifs

    x_coords = concatenate(([-50, 0], epaisseur_cumulee,[epaisseur_cumulee[-1]+1,epaisseur_cumulee[-1]+51]))

    y_coords = [real(substrat), real(substrat)] + partie_reelle_indices + [1,1]

    axes[1, 0].plot(x_coords, y_coords, drawstyle='steps-pre')

    axes[1, 0].set_xlabel('Epaisseur cumulée (nm)')

    axes[1, 0].set_ylabel('Partie réelle de l\'indice')

    axes[1, 0].set_title('Partie réelle de l\'indice en fonction de l\'épaisseur')

    # Graphique d'empilement

    afficher_empilement(epaisseurs, indices, axes[1, 1])  # Pass axes[1, 1] to the function

    tight_layout()

    show()

def afficher_empilement(epaisseurs, indices, ax_empilement):

    """Affiche l'empilement des couches."""

    epaisseurs_nm = array(epaisseurs)

    epaisseurs_nm = epaisseurs_nm[::-1]

    couleurs = ['blue' if i % 2 == 0 else 'red' for i in range(len(indices))]

    ax_empilement.barh(range(len(indices)-1,-1,-1), epaisseurs_nm, align='center', color=couleurs)

    ax_empilement.set_yticks(range(len(indices)))

    ax_empilement.set_yticklabels([f"Couche {i+1} (n={n:.2f})" for i, n in enumerate(indices)])

    ax_empilement.set_xlabel('Epaisseur (nm)')

    ax_empilement.set_title('Empilement des couches (substrat en bas)')

    for i, (epaisseur, indice) in enumerate(zip(epaisseurs_nm, indices)):

        ax_empilement.text(epaisseur / 2, len(indices)-i-1, f"{epaisseur:.1f} nm", va='center', ha='center', color='white')

    print("Liste des couches (depuis le substrat):")

    for i, (epaisseur, indice) in enumerate(zip(epaisseurs_nm[::-1], indices[::-1])):

        print(f"Couche {i+1}: epaisseur = {epaisseur:.1f} nm, indice = {indice:.2f}")

def lancer_calcul():

    """Fonction principale pour lancer le calcul et l'affichage."""

    mat_H, mat_L, sub, l0, emp, m, inter, p = recuperer_parametres()

    resultats_spectral, epaisseurs = calcul_empilement_couches_minces(mat_H, mat_L, sub, l0, emp, 'spectral', (intervalle_spectral_debut.get(), intervalle_spectral_fin.get()), pas_spectral.get())

    resultats_angulaire, epaisseurs = calcul_empilement_couches_minces(mat_H, mat_L, sub, l0, emp, 'angulaire', (intervalle_angulaire_debut.get(), intervalle_angulaire_fin.get()), pas_angulaire.get())

    afficher_parametres(mat_H, mat_L, sub, l0, emp, m, inter, p, epaisseurs)

    tracer_graphiques(resultats_spectral, resultats_angulaire, epaisseurs, float(materiau_H_real.get()), float(materiau_H_imag.get()), float(materiau_L_real.get()), float(materiau_L_imag.get()), float(substrat.get()), empilement_str.get())

    fenetre_principale.destroy()

# --- Interface graphique ---

fenetre_principale = tk.Tk()

fenetre_principale.title("Paramètres de simulation")

# Variables

materiau_H_real = tk.DoubleVar(value=2.1)

materiau_H_imag = tk.DoubleVar(value=0.001)

materiau_L_real = tk.DoubleVar(value=1.5)

materiau_L_imag = tk.DoubleVar(value=0.0001)

substrat = tk.DoubleVar(value=1.45)

lambda_0 = tk.DoubleVar(value=500)

empilement_str = tk.StringVar(value="3,3,3,3,3,3,3")

mode = tk.StringVar(value="spectral")

intervalle_spectral_debut = tk.DoubleVar(value=150)

intervalle_spectral_fin = tk.DoubleVar(value=700)

intervalle_angulaire_debut = tk.DoubleVar(value=0)

intervalle_angulaire_fin = tk.DoubleVar(value=89)

pas_spectral = tk.DoubleVar(value=1)

pas_angulaire = tk.DoubleVar(value=1)

# Création des widgets

# Matériau H

tk.Label(fenetre_principale, text="Matériau H (H)").grid(row=0, column=0, sticky="w")

tk.Label(fenetre_principale, text="Partie réelle:").grid(row=1, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=materiau_H_real).grid(row=1, column=1)

tk.Label(fenetre_principale, text="Partie imaginaire:").grid(row=2, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=materiau_H_imag).grid(row=2, column=1)

# Matériau L

tk.Label(fenetre_principale, text="Matériau L (L)").grid(row=3, column=0, sticky="w")

tk.Label(fenetre_principale, text="Partie réelle:").grid(row=4, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=materiau_L_real).grid(row=4, column=1)

tk.Label(fenetre_principale, text="Partie imaginaire:").grid(row=5, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=materiau_L_imag).grid(row=5, column=1)

# Substrat

tk.Label(fenetre_principale, text="Substrat").grid(row=6, column=0, sticky="w")

tk.Label(fenetre_principale, text="Indice:").grid(row=7, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=substrat).grid(row=7, column=1)

# Longueur d'onde de centrage

tk.Label(fenetre_principale, text="Longueur d'onde de centrage (nm)").grid(row=8, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=lambda_0).grid(row=8, column=1)

# Empilement

tk.Label(fenetre_principale, text="Empilement (QWOT, séparés par des virgules)").grid(row=9, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=empilement_str).grid(row=9, column=1)

# Mode

tk.Label(fenetre_principale, text="Mode").grid(row=10, column=0, sticky="w")

tk.Radiobutton(fenetre_principale, text="Spectral", variable=mode, value="spectral").grid(row=11, column=0, sticky="w")

tk.Radiobutton(fenetre_principale, text="Angulaire", variable=mode, value="angulaire").grid(row=12, column=0, sticky="w")

# Intervalle spectral

tk.Label(fenetre_principale, text="Intervalle spectral (nm)").grid(row=13, column=0, sticky="w")

tk.Label(fenetre_principale, text="Début:").grid(row=14, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=intervalle_spectral_debut).grid(row=14, column=1)

tk.Label(fenetre_principale, text="Fin:").grid(row=15, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=intervalle_spectral_fin).grid(row=15, column=1)

# Intervalle angulaire

tk.Label(fenetre_principale, text="Intervalle angulaire (degrés)").grid(row=16, column=0, sticky="w")

tk.Label(fenetre_principale, text="Début:").grid(row=17, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=intervalle_angulaire_debut).grid(row=17, column=1)

tk.Label(fenetre_principale, text="Fin:").grid(row=18, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=intervalle_angulaire_fin).grid(row=18, column=1)

# Pas

tk.Label(fenetre_principale, text="Pas").grid(row=19, column=0, sticky="w")

tk.Label(fenetre_principale, text="Spectral (nm):").grid(row=20, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=pas_spectral).grid(row=20, column=1)

tk.Label(fenetre_principale, text="Angulaire (degrés):").grid(row=21, column=0, sticky="w")

tk.Entry(fenetre_principale, textvariable=pas_angulaire).grid(row=21, column=1)

# Bouton de calcul

tk.Button(fenetre_principale, text="Lancer le calcul", command=lancer_calcul).grid(row=22, column=0, columnspan=2)

fenetre_principale.mainloop()