%% tmpfit1.m is designed to determine the mean latitudinal IRI temperature 
%% profiles using both the electron temperture models in the IRI2001. It
%% conducts fits in 1995 (IRI model, so choice of year unimportant except 
%% for effect of neutral temperature which does increase with F107) at 21 
%% March, 21 June and 21 Dec for equinox, summer and winter in both 
%% hemispheres at a height of 400km. 
%%

clear all

%Set up required arrays
latmatm=deg2rad(-90:5:90))';
lonlsh=deg2rad(-180:20:160);
nopy=length(latmatm);
nopx=length(lonlsh);
svy=ones(nopy,1);
svx=ones(1,nopx);
svm=ones(nopy,nopx);
[lonmat,latmat]=mag2geo3(svy*lonlsh,latmatm*svx);
re=6371;
dismat=re+400;
year=1995;
day=21;

timld=13;
timmatd=timl2ut(timld*svm,lonmat);

%Conduct fit at 1300LT across globe on March 21, June 21 and December 21
[tmpi3,tmpe3,tmpn3]=tempiri(lonmat,latmat,dismat*svm,timmatd,year*svm,3*svm,day*svm);
[tmpi6,tmpe6,tmpn6]=tempiri(lonmat,latmat,dismat*svm,timmatd,year*svm,6*svm,day*svm);
[tmpi12,tmpe12,tmpn12]=tempiri(lonmat,latmat,dismat*svm,timmatd,year*svm,12*svm,day*svm);

timln=1;
timmatn=timl2ut(timln*svm,lonmat);

%Conduct fit at 0100 across globe on March 21, June 21 and December 21
[tmpi3n,tmpe3n,tmpn3n]=tempiri(lonmat,latmat,dismat*svm,timmatn,year*svm,3*svm,day*svm);
[tmpi6n,tmpe6n,tmpn6n]=tempiri(lonmat,latmat,dismat*svm,timmatn,year*svm,6*svm,day*svm);
[tmpi12n,tmpe12n,tmpn12n]=tempiri(lonmat,latmat,dismat*svm,timmatn,year*svm,12*svm,day*svm);

%Average the day results longitudinally
tmpe3dav=mean(tmpe3,2);
tmpe6dav=mean(tmpe6,2);
tmpe12dav=mean(tmpe12,2);

%Average the night results longitudinally
tmpe3nav=mean(tmpe3n,2);
tmpe6nav=mean(tmpe6n,2);
tmpe12nav=mean(tmpe12n,2);

latd=rad2deg(latmatm);
latd2=latd(46:91);

%Calculate the mean day equinox, summer and winter profiles
tmpdeq=0.5*(flipud(tmpe3dav(1:46))+tmpe3dav(46:91));
tmpdsu=0.5*(flipud(tmpe12dav(1:46))+tmpe6dav(46:91));
tmpdwn=0.5*(flipud(tmpe6dav(1:46))+tmpe12dav(46:91));
tmpdyr=0.25*(2*tmpdeq+tmpdsu+tmpdwn);

%Calculate the mean night equinox, summer and winter profiles
tmpneq=0.5*(flipud(tmpe3nav(1:46))+tmpe3nav(46:91));
tmpnsu=0.5*(flipud(tmpe12nav(1:46))+tmpe6nav(46:91));
tmpnwn=0.5*(flipud(tmpe6nav(1:46))+tmpe12nav(46:91));
tmpnyr=0.25*(2*tmpneq+tmpnsu+tmpnwn);

sl=sin(deg2rad(latd2)).^2;

