// **************
// **************
// ** Donnéees **
// **************
// **************

chemin = ""; // depends on your config
nomfic = "verre_HOQ.csv"

donnees = csvRead(chemin+nomfic, ";");
lambda = 1e-3*donnees(:, 1); // en µm
nexp = donnees(:, 2);

KL0 = [1, 0.2, 1, 0.006, 0.02, 100];

// ***************
// ***************
// ** Fonctions **
// ***************
// ***************

function [n] = sellmeier(KL, lambda)
    K1 = KL(1); K2 = KL(2); K3 = KL(3);
    L1 = KL(4); L2 = KL(5); L3 = KL(6);
    lambdacarre = lambda.*lambda;
    foo1 = lambdacarre - L1;
    foo2 = lambdacarre - L2;
    foo3 = lambdacarre - L3;
    foo = 1./(foo1.*foo2.*foo3);
    n = sqrt(1 +...
    lambdacarre.*(K1*foo.*foo2.*foo3 + K2*foo.*foo1.*foo3 + K3*foo.*foo1.*foo2));
endfunction

function [err] = ecart(KL, lambda, n)
    err = sellmeier(KL, lambda) - n;
endfunction

// ****************
// ****************
// ** Régression **
// ****************
// ****************

[S, KLopt] = leastsq(list(ecart, lambda, nexp), KL0);

// ***********
// ***********
// ** Tracé **
// ***********
// ***********

scf(1);
clf();
plot(lambda, sellmeier(KLopt, lambda));
plot(lambda, nexp, "+");