%% Analyse d'erreur sur la mesure d'un diviseur de Tension
% Table of Contents
% Schéma du diviseur
%     Paramètres
% Résolution du système d'équations
%         Détermination des estimateurs des résistances internes
% Calcul d'erreur
%     Simulation de la mesure
%         Données du système réel
%         Calcul de la tension de sortie
%         Calcul des résistances internes
%     Calcul des incertitudes sur les résistances internes
%% Schéma du diviseur
%% Paramètres
% Tension d'entrée : U_{in}
% Tension de sortie: U_{out}
% Impédance interne: R_1 et R_2
% Impédance externe: Z_{ex}
%
%% Résolution du système d'équations
%% Détermination des estimateurs des résistances internes
% On détermine A, le rapport des tensions d'entrée et de sortie à partir des lois de Kirchhof. L'équation

% Impédances
syms R_1 R_2 R_ex 

% Tensions
syms U_in U_out

% Lois de Kirchhoff
eqn_Kirchhof=(U_in-U_out)/R_1==U_out*(1/R_2+1/R_ex);
disp(eqn_Kirchhof)

% permet de déterminer

A=solve(eqn_Kirchhof,U_out)/U_in

%
% Le calcul de R_1 et R_2 se fait à partir de deux mesures de A, le rapport des tensions, pour deux impédances externes différentes:
% Pour la première mesure A_1,  R_{ex}=1 \Omega
% Pour la deuxième mesure A_2, R_{ex}=1k\Omega par exemple et le rapport des tensions est simplifié

A_1=A
A_2=R_2/(R_1+R_2)

% Cela établit un système d'équations qui permet de déterminer les estimateurs de R_1 et R_2

% Création des variables correspondant aux rapports des mesures de tension
syms U_in1 U_out1 U_in2 U_out2
% Système d'équations
eqn=[A_1==U_out1/U_in1, A_2==U_out2/U_in2];
disp(eqn)

%


% Résolution du système d'équation
Estimateurs=solve(eqn,[R_1 R_2])
% Les solutions sont données dans la structure 'Estimateur'

%
% On place les solutions dans les variables r_1 et r_2

r_1=simplify(collect(Estimateurs.R_1,[U_in1,U_in2]),10);

r_1=collect(r_1,R_ex)
r_2=collect(Estimateurs.R_2,R_ex)

% La variable R_{ex} apparaissant dans l'équation est larésistance externe pour le premier essay à basse résistance. La haute résistance du second essai a été simplifié et n'apparait plus dans le calcul.
%
%% Calcul d'erreur
% On commence par calculer systématiquement les dérivées individuelles \frac{\partial R_i}{\partial x_j}. On néglige l'incertitude sur R_{ex}.

% dérivées de R1
d_11=diff(r_1,U_in1);
d_12=diff(r_1,U_out1);
d_13=diff(r_1,U_in2);
d_14=diff(r_1,U_out2);

% dérivées de R2
d_21=diff(r_2,U_in1);
d_22=simplify(diff(r_2,U_out1));
d_23=simplify(diff(r_2,U_in2));
d_24=simplify(diff(r_2,U_out2));

%
% On définit des variables pour les incertitudes de mesure.

syms deltaU_in1 deltaU_out1 deltaU_in2 deltaU_out2 

% On calcule l'incertitude sur R_1et R_2 en utilisant les dérivées individuelles:

delta_R1=d_11*deltaU_in1...
    -d_12*deltaU_out1...
    -d_13*deltaU_in2...
    +d_14*deltaU_out2

delta_R2=collect(d_21*deltaU_in1...
    -d_22*deltaU_out1...
    -d_23*deltaU_in2...
    +d_24*deltaU_out2,R_ex)

%
%% Simulation de la mesure
% On simule le système en définissant les résistances R_1et R_2 et on calcule les tensions U_{in} et U_{out} que l'on mesurerait dans la réalité.
%% Données du système réel
% On fixe les résistances interne (en Ohms):

R1=100;
R2=10;

% On fixe les résistances externes pour la première et la deuxième mesure (en Ohms):

Z1=1;
Z2=1000;

% On fixe la tension d'entrée (en Volt):

Uin=100;

% On fixe l'erreur relative des mesures de tension:

erreur_relative=0.01

%
%% Calcul de la tension de sortie
% On calcul la tension de sortie pour différent cas:
%
% 1) Le cas de basse impédance

% Tension de sortie pour la basse résistance
Uout1=Uin*R2*Z1/((R1+R2)*Z1+R1*R2);
Uout1=Uout1*(1+erreur_relative*randn(1,1))

% Cela donne le rapport de tension entre sortie et entrée

A1=Uout1/Uin

% On calcule l'incertitude des deux mesures de tension en tenant compte des données du fabricant

% incertitude de la tension d'entrée
[Pin,Din]=incertitude_AM500(Uin);
delta_Uin1=Uin*Pin+Din

% incertitude de la tension de sortie
[Pout1,Dout1]=incertitude_AM500(Uout1);
delta_Uout1=Uout1*Pout1+Dout1

%
%
% 2) Le cas de haute impédance, calculé avec simplification

Uout2=Uin*R2/(R1+R2);
Uout2=Uout2*(1+erreur_relative*randn(1,1))

% Cela donne le rapport de tension entre sortie et entrée

A2=Uout2/Uin

% On calcule l'incertitude de la mesure de tension en fonction des données du fabricant

delta_Uin2=Uin*Pin+Din
[Pout2,Dout2]=incertitude_AM500(Uout2);
delta_Uout2=Uout2*0.01+0.3
errorbar([1 2 3],[Uin,Uout1,Uout2],[Uin*erreur_relative,delta_Uout1,delta_Uout2],'.r',"LineWidth",2)
grid on
xlim([0.5,3.5])
set(gca,"YLim",[0 110],"XTick",[1,2,3],"XTickLabel",["U_{in}","U_{out1}","U_{out2}"])
title("Mesures de tension")
ylabel("Volts")

% 3) Le cas de haute impédance, calculé sans simplification

Uout=Uin*R2*Z2/((R1+R2)*Z2+R1*R2);
Uout=Uout*(1+erreur_relative*randn(1,1))

% On calcul la différence relative entre les résultats avec et sans simlification

delta=(Uout-Uout2)/Uout

% On vérifie si la simplification est praticable
%
%% Calcul des résistances internes
% On calcule les deux résistances

R1_est=Z1*(A2-A1)/A2/A1
R2_est=Z1*(A2-A1)/(1-A2)/A1

% On peut observer que la mesure simule correctement le système défini.
%
%% Calcul des incertitudes sur les résistances internes
% On calcule les erreurs de mesure

dR1=double(subs(delta_R1,...
{R_ex,U_in1,U_out1,U_in2,U_out2,deltaU_in1,deltaU_out1,deltaU_in2,deltaU_out2},...
{Z1,Uin,Uout1,Uin,Uout2,delta_Uin1,delta_Uout1,delta_Uin2,delta_Uout2}))

dR2=double(subs(delta_R2,...
{R_ex,U_in1,U_out1,U_in2,U_out2,deltaU_in1,deltaU_out1,deltaU_in2,deltaU_out2},...
{Z1,Uin,Uout1,Uin,Uout2,delta_Uin1,delta_Uout1,delta_Uin2,delta_Uout2}))

% Présentation graphique des résultats

errorbar([1 2 ],[R1 R2],[dR1 dR2],'or',"LineWidth",2)
grid on
xlim([0.5,2.5])
set(gca,"XTick",[1,2],"XTickLabel",["R_1","R_2"])
title("Résistances internes")
ylabel("Ohms")

%
%
%

function [pourcent,digit]=incertitude_AM500(V)

if V < 0.4
    disp("AM500: Plage de 400.0 mV")
    pourcent=0.012;
    digit=0.0003;
elseif V < 4
    disp("AM500: Plage de 4.000 V")
    pourcent=0.01;
    digit=0.003;
elseif V < 40
    disp("AM500: Plage de 40.00 V")
    pourcent=0.01;
    digit=0.03;
elseif V < 400
    disp("AM500: Plage de 400.0 V")
    pourcent=0.01;
    digit=0.3;
else
    disp("AM500: Plage de 600 V")
    pourcent=0.012;
    digit=3;
end


end





