clear all;

% Dec 2011  Program to average UT and UT+1 to make UT+0.5 from TIE model outputs
%        from 07325-08026 for 8 longitudes from output of tec8lonclim.m in ccmc_output
%   tec8lonclim.m has correct lon, 
%   and calc_med_out.m gets lat at -60 to -65 from -61.25 to -63.75 so is off some
%   Will correct calc_med(12)_out.m to get single lat at -62.5!

% 8 begin lons (for 5 deg glon): 25,90,140,175,200,250,285,345
glon=[25,90,140,175,200,250,285,345];
% for dres TIEGCM, for 25E have glon=27.5+/-1.25deg or 26-29
  
%  The model longitudes are 25 deg wide like COSMIC - not like MIT GPS TEC
% lon_LLL_1_IRI(CTIPE,TIE-GCM)_2007.325-2008.020.txt where LLL=glon-10

% day 42 (Jan 1, 2008) is bad for 0 UT dres W05 TIMED lb
ndays = 61+6;

n36 = 72;

year1 = 2007;
ida0 = 324;

for lon=1:8
 glon1 = glon(lon)+1;
 glon2 = glon(lon)+4;
 gloncen = glon(lon)+2.5;	

% Get name for output file of median TEC values for each glon
% Have double the number of lats (36*2=72)
outfile=sprintf('ccmc_output/lon_%3.3d_07325_08026_2_TIE-GCMc.dat',glon(lon))
OUT = fopen(outfile,'w');

% Get name of 61+6-day output file (24h*61d*36glats=52704 lines each)
infile=sprintf('ccmc_output/lon_%3.3d_2_TIE-GCM_2007.325-2008.026c.txt',glon(lon))

% Open file
INP = fopen(infile,'r');

for n=1:ndays
for nh=1:24
yearwr = year1;
ida = ida0 + n;
jda = ida;
if (ida > 365)
 yearwr = year1 + 1;
 jda = jda - 365;
end

uth = nh-1.5;
mh = nh - 2;
if (nh==1)
 mh = 23;
 uth = 23.5;
 jda = jda - 1;
 if jda==0
  jda = 365;
  yearwr = year1;
 end
end

for nlat=1:n36
 glatcen = -90. + (nlat-1)*2.5;
 line = fgetl(INP);
% year daynum beghr begmin begsec beg_glat beg_glon TEC(TECU) NmF2(cm-5?) hmF2(km)
% dres TIE adds wiav(m/s) on end
 input = sscanf(line,'%f');
iav = 1;
if (n==1 & nh==1)
 iav = 0;
end
% Avoid bad input from o UT 1 Jan 2008 - set to 23 UT 31 Dec 2007
if (nh==1 & n==42)
 input(8) = inputp(nlat,8);
 input(9) = inputp(nlat,9);
 input(10) = inputp(nlat,10);
 input(11) = inputp(nlat,11);
end	
if iav==1
% calculate ave and print out if iav=1
tec = 0.5*(input(8)+inputp(nlat,8));
nmax = 0.5*(input(9)+inputp(nlat,9));
hmax = 0.5*(input(10)+inputp(nlat,10));
avwi = 0.5*(input(11)+inputp(nlat,11));
fprintf(OUT,'%4.4d %3.3d %2.2d 30 0 %5.1f %5.1f %5.2f %e %7.2f %7.2f \n',yearwr,jda,mh,glatcen,gloncen,tec,nmax,hmax,avwi);
end	% iav==1

% Set inputp (previous UT time) for all glats
for k=1:11
 inputp(nlat,k) = input(k);
end
end	% nlat=1:36
end	% nh=1:24
end	% n=1:ndays
  
end	% lon=1:8




