%% apmsis2.m produces the required data arrays containing AP
%% infomration over the past 2 1/2 days as required by 
%% MSIS-86 model. Each row is for a 3 hour period, with the
%% columns being day number, UT, required AP values for 
%% MSIS-86. The required Ap values are daily AP, and then the
%% following current 3 hour Ap values - current, 3 hour previously,
%% 6 hours previously, 9 hours previously, average of eight 3 hour
%% values from 12 to 33 hours previously and average of eight 3 hour
%% values from 36 to 59 hours previously. Note that the daily Ap is
%% the average of the previous 24 hours.
%% The files are saved every year, with the naming format being
%% apm# where # = year.
%% apmsis2.m is the same as apmsis.m, except that it saves the data
%% to yearly files, not monthly. 

clear

years=input('Starting File Year? ');
yearf=input('Finishing File Year? ');

nop=yearf-years+1;                           %number of files to be converted
ap21=zeros(1,21);                            %current daily Ap plus latest 20 3-hour Ap values

for i=1:nop
 if i==1
  yeari=years;                          %year index
 else
  yeari=yeari+1;                        %increment year index
 end 
 
 eval(['load a',int2str(yeari)])        %load in i'th year data file
 ap3h=ap3h';                            %take transpose of 3-hour Ap matrix data
 
 dn=0;                                  %current day number
 h3i=0;                                 %current 3-hour index
 jmax=8*length(apv);                    %number of 3-hour periods in current year
 apmat=zeros(jmax,9);                   %results matrix for i'th year

 for j=1:jmax                           %loop through each day of i'th year
  if rem(j-1,8)==0                      %check if first 3-hour period of new day
   dn=dn+1;                             %update day number
   timut=0;                             %reset UT
  else 
   timut=timut+3;                       %update new UT
  end
   
  ap21(3:21)=ap21(2:20);                %move pass 3-hour Ap values back one time interval
  h3i=h3i+1;                            %increment h3i by one
  ap21(2)=ap3h(h3i);                    %j'th 3-hour Ap increment
  ap21(1)=mean(ap21(2:9));              %set daily Ap for 3hours equal to mean of latest eight 3-hour increments
  
  %Record the j'th Ap values
  apmat(j,:)=[dn timut ap21(1:5) mean(ap21(6:13)) mean(ap21(14:21))];   
 end
 
 %Save the results for the i'th year, j'th month
 eval(['apm',int2str(yeari),'=apmat;'])
 eval(['save apm',int2str(yeari),' apm',int2str(yeari)])  
 eval(['clear apm',int2str(yeari)]) 
end 

 