function I = Cell_UJT(V,G,TaC)
%Code produced by the team from the CubeSat SwissCube
%modified by Philippe Ledent (OUFTI-1)
%modified by Francois Grosjean (OUFTI-2)
%
%Model of 28% efficiency solar cell from Spectrolab : I=f(V,T)
%Use of function : I = Cell_GaAs(V,G,TaC)
%V = Voltage on cell terminals [V]
%G = relative insolation [-] (G=1 => 1367 W/m^2)
%TaC = temperature of the cell in operation [Celsius]

%Boltzman constant
k = 1.38e-23;
%Electric charge
q = 1.60e-19;
%Quality factor of the diode (1<n<2)
n = 1.5;
%Band gap voltage (1.12eV < Vg < 1.757eV)
Vg = 1.75;
%Reference values
A=26.62;%Area of the cell [cm2]
Tref = 273 + 28; %temperature
Voc_Tref = 2.660; %open circuit voltage (G=1 et T=Tref)
Isc_Tref = 17.05*A/1000; %short circuit current (G=1 et T=Tref)[A]
%Temperature of the cell in operation
TaK = 273 + TaC;
%Photo-current thermal coefficient
K0 = 2.72e-4/Isc_Tref;
%Photo-current (G=1 et T=Tref)
Iph = Isc_Tref * G * (1 + K0*(TaK - Tref));
%Diode saturation current (T=Tref)
Id_Tref = Isc_Tref / (exp(q*Voc_Tref/(n*k*Tref))-1);
%Diode saturation current (T=Tak)
Id = Id_Tref * (TaK/Tref)^(3/n) * exp(-(Vg * q/(n*k))*(1/TaK - 1/Tref));
%Calculation of serie resistance Rs
Vpmax_Tref = 2.350; %voltage at maximum power for T=Tref
Ipmax_Tref = 16.3*A/1000; %current at maximum power for T=Tref
Rs = (Voc_Tref-Vpmax_Tref)/Ipmax_Tref;
%Iterative calculation of Isc
I = zeros(size(V));
for i=1:10 %number of iteration = 10
    for j=1:length(V)
        I(j) = I(j) - (Iph - I(j) - Id*(exp(q*(V(j)+I(j)*Rs)/(n*k*Tref)) -1))/...
            (-1 - (Id*(exp(q*(V(j)+I(j)*Rs)/(n*k*Tref)) -1))*q*Rs/(n*k*Tref));
    end
end