function fit_diode_curve_sem_drain_body()
    global data
    data = readmatrix("DRAIN_BODY2023_02_10_150122.csv", "Delimiter",";");

    global I_L
    I_L = 10e-11;%find_light_generated_current(data(:, 1), data(:, 2));
    fprintf("Light generated current I_L = %g\n", I_L);
    shifted_current = data(:, 2) - I_L;
    [I_S, n] = fit_diode_forward(data(:, 1), shifted_current, 0.0);
    fprintf("I_S = %g\n", I_S);
    fprintf("n = %f\n", n);
    [a, b] = fit_diode_reverse(data(:, 1), shifted_current, I_S);
    fprintf("a = %g b = %g\n", a, b);
    semilogy(data(:,1), abs(data(:,2)), ".r")
    hold on
    plot_domain = linspace(min(data(:,1)), 10.0, length(data(:,1)));
    %semilogy(data(:,1), abs(shifted_current), ".g")
    semilogy(plot_domain, abs(shockley(plot_domain, I_S, n, a, b)+I_L), "-b")
    %semilogy(plot_domain, abs(shockley(plot_domain, I_S, n, a, b)), "-k")
    %legend('Experimental Data', 'Exp. Data I_L Shifted','Fitted Model', 'Fitted Model I_L Shifted','Location','northwest')
    legend('Experimental Data', 'Fitted Model', 'Location','northwest')
    %semilogy(plot_domain, abs(shockley_fit(plot_domain)), "-k")
    hold off
    xlabel("Voltage / V")
    ylabel("Current / A")
    title("Characteristic Curve of p-MOSFET")
    grid()

    return;

function current = shockley(voltage, I_S, n, a, b)
    %current = zeros(size(voltage))
    current = I_S.*(exp(1.60217e-19*voltage./(n.*293.15.*1.3806e-23))-1.0);
    current(voltage<a) = current(voltage<a); %- (b.*(a-voltage(voltage<a))).^(1/2);



function [I_sat, n] = fit_diode_forward(voltage, current, resistance)
    threshold_voltage = 4.2;
    upper_voltage_threshold = 10.0;
    k_B = physconst('Boltzmann');
    global I_L
    voltage = voltage - (current+I_L).*resistance;
    current = current(voltage>threshold_voltage);
    voltage = voltage(voltage>threshold_voltage);
    current = current(voltage<upper_voltage_threshold);
    voltage = voltage(voltage<upper_voltage_threshold);
    p = polyfit(voltage, log(current), 1);
    I_sat = exp(p(2));
    n = 1./(293.15*k_B*p(1)./(1.602e-19));


function [a, b] = fit_diode_reverse(voltage, current, I_S)
    k_B = physconst('Boltzmann');
    T = 293.15;
    threshold_voltage = -0.1;
    current = current(voltage<threshold_voltage);
    voltage = voltage(voltage<threshold_voltage);
    current_hat = current + I_S;
    p12 = polyfit(voltage, current_hat.^2, 1);
    p13 = polyfit(voltage, current_hat.^3, 1);
    %I_R = I_S + sqrt(b)*sqrt(a-V) for abrupt junction
    b = -p12(1);
    a = p12(2)./b;

    %{
function current = simulate_diode_reverse_AlGaAs(voltage, x)
    q = 1.602176e-19;
    k = physconst("Boltzmann");
    m_0 =
    %source = ioffe.ru/SVA/NSM/AlGaAs/index.html
    eps_s = 12.90-2.84.*x;
    mu_p = 200-550.*x250.*x.^2;
    mu_n = 8e3-2.2e4.*x+1e4.*x^2;
    m_p = (0.063+0.083.*x).*m_0;
    m_n = ()
    tau_p = mu_p*m_p./q; % have diffusion coeff and mobility
    tau_n = mu_n*m_n./q;
    D_p = 9.2-24.0*x+18.5.*x.^2;
    D_n = 200.0-550.0.*x+250.0.*x.^2;
    E_g = 1.424+1.155.*x+0.37.*x.^2; %eV
    N_C = 2.5e19.*(.063+0.083.*x).^(3/2);% eff. DOS in conduction band
    N_V = 2.5e19.*(.51+.25.*x).^(3/2)% eff. DOS in valence band
    N = %=N_A or N_D which ever is larger
    T = 273.15 + 23;
    W_D = sqrt(2*eps_s./(q*N).*(psi_bi - voltage - 2.*k*T/q)); %one sided abrupt junction
    n_i = sqrt(N_C.*N_V).*exp(-E_g./(k.*T));
    % n + N_A = p + N_D
    % p*n = n_i^2
    tau_g = %
    I_ge = q*n_i.*W_D./tau_g
    I_sat = q.*sqrt(D_p./tau_p).*n_i.^2/N_D;% for one sided abrupt junction with p_n0 >> _p0
%}




function light_current = find_light_generated_current(voltage, current)
    zero_volt_index = find(min(abs(voltage))==abs(voltage));
    if voltage(zero_volt_index)<0.0
        voltage1 = voltage(zero_volt_index);
        current1 = current(zero_volt_index);
        voltage2 = voltage(zero_volt_index+1);
        current2 = current(zero_volt_index+1);
    else
        voltage1 = voltage(zero_volt_index-1);
        current1 = current(zero_volt_index-1);
        voltage2 = voltage(zero_volt_index);
        current2 = current(zero_volt_index);
    end
    light_current = (current1.*(voltage2-0.0)-current2.*(voltage1))...
        ./(voltage2-voltage1);

Light generated current I_L = 1e-10
I_S = 4.61933e-12
n = 29.503181
a = -1.62781 b = 6.1705e-20

Published with MATLAB® R2023a