%% temppaw33.m is the third part of the temppaw3#.m series. It 
%% determines the corrections due to the current Kp. See temppaw31.m
%% for more details.
%% It is designed to be called by plasimod.m and should probably not
%% be called otherwise as it requires a specific input format.
%% 

function [t0d,t0dc,t0n,t0nc,g0d,g0dc,g0n,g0nc]=temppaw33(t0d,t0dc,t0n,t0nc,g0d,g0dc,g0n,g0nc,kp,sl)

% t0d = day base temperature (K)
% t0dc = conjugate day base temperature (K)
% t0n = night base temperature (K)
% t0nc = conjugate night base temperature (K)
% g0d = day gradient temperature (K/km)
% g0dc = conjugate day gradient temperature (K/km)
% g0n = night gradient temperature (K/km)
% g0nc = conjuagte night gradient temperature (K/km)
% kp = 3-hourly Kp (scalar)
% sl = sl parameter

% Modify the values of G0 and T0 to reflect changes in magnetic activity
latm3=asin(sqrt(sl));                                  %magnetic latitude at 300km altitude on L shell
l3i=abs(latm3)>deg2rad(46);                            %index of latm3 values above 46 degrees
if any(l3i)
 l3=abs(latm3(l3i));                                   %latm3 values above 46 degrees
 l3=rad2deg(l3);                                       %convert l3 to degrees
 dk=37+1.33*(l3-46)-37*cos(0.135*(l3-46));             %Dk parameter (eqn 17)
 t0dif=dk.*(kp-1.5);                                   %T0 correction
 g0ddif=(t0d(l3i)./(t0d(l3i)+t0dif)).^(2.5);           %day G0 correction term
 g0ddifc=(t0dc(l3i)./(t0dc(l3i)+t0dif)).^(2.5);        %conjugate day G0 correction term
 g0ndif=(t0n(l3i)./(t0n(l3i)+t0dif)).^(2.5);           %night G0 correction term 
 g0ndifc=(t0nc(l3i)./(t0nc(l3i)+t0dif)).^(2.5);        %conjugate night G0 correction term 
 g0d(l3i)=g0d(l3i).*g0ddif;                            %modifed G0 day (eqn 18b)
 g0dc(l3i)=g0dc(l3i).*g0ddifc;                         %modifed conjugate G0 day (eqn 18b)
 g0n(l3i)=g0n(l3i).*g0ndif;                            %modifed G0 night (eqn 18b)
 g0nc(l3i)=g0nc(l3i).*g0ndifc;                         %modifed conjugate G0 night (eqn 18b)
 t0d(l3i)=t0d(l3i)+t0dif;                              %modifed T0 day (eqn 18a)
 t0dc(l3i)=t0dc(l3i)+t0dif;                            %modifed conjugate T0 day (eqn 18a)
 t0n(l3i)=t0n(l3i)+t0dif;                              %modifed T0 night (eqn 18a)
 t0nc(l3i)=t0nc(l3i)+t0dif;                            %modifed conjugate T0 night (eqn 18a)
end 



