clear all;
close all;

% Dec 2012  Program to calculate the ave of CHAMP PLP Te and Ne 
%    from the C/NOFS estimates over various seasons

lt24 = 0.5:1:23.5;

% CHAMP PLP Ne and Te
%chme07=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2007_290')
 chme07=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2007_15')
 chte07=dlmread(chme07,'',0,0);
%chme08=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2008_290')
 chme08=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2006_15')
 chte08=dlmread(chme08,'',0,0);
%chme09=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2009_290')
 chme09=sprintf('/hao/aim3/emery/cetic_climatology/champ_plp_cnofs_vs_lt_2005_15')
 chte09=dlmread(chme09,'',0,0);
% outcnofs = sprintf('champ_plp_cnofs_vs_lt_290_2007-9')
  outcnofs = sprintf('champ_plp_cnofs_vs_lt_015_2005-7')
  OUTC = fopen(outcnofs,'w');

htn = 999.;
htx = 0.;
for nh=1:24
 numt = chte07(nh,2)+chte08(nh,2)+chte09(nh,2);
 if (numt>0)
  chte(nh) = (chte07(nh,2)*chte07(nh,3) + chte08(nh,2)*chte08(nh,3) + chte09(nh,2)*chte09(nh,3))/numt;
  chne(nh) = (chte07(nh,2)*chte07(nh,4) + chte08(nh,2)*chte08(nh,4) + chte09(nh,2)*chte09(nh,4))/numt;
  chht1 = (chte07(nh,2)*chte07(nh,5) + chte08(nh,2)*chte08(nh,5) + chte09(nh,2)*chte09(nh,5))/numt;
  chht2 = (chte07(nh,2)*chte07(nh,6) + chte08(nh,2)*chte08(nh,6) + chte09(nh,2)*chte09(nh,6))/numt;
  chhtav = 0.5*(chht1+chht2);	 
  chhtsig = 0.;
  numstd = 0;
  if (chte07(nh,2)>0)
   numstd = numstd + 2;
   chhtsig = chhtsig + (chhtav-chte07(nh,5))^2 + (chhtav-chte07(nh,6))^2;
  end
  if (chte08(nh,2)>0)
   numstd = numstd + 2;
   chhtsig = chhtsig + (chhtav-chte08(nh,5))^2 + (chhtav-chte08(nh,6))^2;
  end
  if (chte09(nh,2)>0)
   numstd = numstd + 2;
   chhtsig = chhtsig + (chhtav-chte09(nh,5))^2 + (chhtav-chte09(nh,6))^2;
  end
  chhtstd(nh) = sqrt(chhtsig/numstd);
  chht(nh) = chhtav - 7.;  % Subtract 7km since Re eq 7km larger
  htn = min(chht(nh),htn);
  htx = max(chht(nh),htx);
 else
   chte(nh)=0.;
   chne(nh)=0.;
   chht(nh)=0.;
   chhtstd(nh)=0.;
 end
   fprintf(OUTC,'%5.1f %d %8.2f %12.4e %6.2f %6.2f \n',lt24(nh),numt,chte(nh),chne(nh),chht(nh),chhtstd(nh));
end	% for nh=1:24
stdx = max(chhtstd(nh))
htn,htx

% Put NaN where numt=0
%chte(3)=NaN;
%chne(3)=NaN;
%chte(15)=NaN;
%chne(15)=NaN;

% make plots
zero=[0,0];
diur=[0,24];
fsize=12;
fsize2=15;

tex=3000.;
ten=500.;
tn=600;
tx=3200;
 nden = 3;
 xden = 6;
cn=3.0;
cx=7.0;

%titne = sprintf('CHAMP PLP lg10Ne cm-3 290-055 2007-9 alt~%d-%d',htn-stdx,htx+stdx);
%titte = sprintf('CHAMP PLP Te (K) 290-055 2007-9 alt~%d-%d',htn-stdx,htx+stdx);
titne = sprintf('CHAMP PLP lg10Ne cm-3 015-145 2005-7 alt~%d-%d',htn-stdx,htx+stdx);
titte = sprintf('CHAMP PLP Te (K) 015-145 2005-7 alt~%d-%d',htn-stdx,htx+stdx);

dark_green = [0 0.5 0];
green_blue = [0 0.75 0.50];
purple = [0.48 0.06 0.89];
light_orange = [1 0.69 0.39];
yellow_green = [0.75 0.75 0];
light_purple = [0.7 0.78 1];
aqua = [0.04 0.52 0.78];

figure;
set (gcf,'position',[142 99 647 822])
H=gcf;
clf
hold on;
title(titte,'FontSize',fsize2,'FontWeight','bold');
plot (diur,zero,'--k');
plot (lt24,chte,'Color',yellow_green,'LineWidth',2);
text (8.5,0.96*(tx-tn)+tn,'Te CHAMP PLP','Color',yellow_green,'FontSize',fsize,'FontWeight','bold'); 
ylabel('Ti,e K','FontSize',fsize,'FontWeight','bold');
xlabel('LT(h)','FontSize',fsize,'FontWeight','bold');
set (gca,'xlim',[0,24],'ylim',[tn tx],'xtick',([0:2:24]),'FontSize',fsize,'box','on');

hold off;

figure;
set (gcf,'position',[142 99 647 822])
H=gcf;
clf
hold on;
title(titne,'FontSize',fsize2,'FontWeight','bold');
plot (diur,zero,'--k');
plot (lt24,log10(chne),'Color',yellow_green,'LineWidth',2);
text (9.0,0.95*(cx-cn)+cn,'Ne CHAMP PLP','Color',yellow_green,'FontSize',fsize,'FontWeight','bold');
ylabel('Ne cm^3','FontSize',fsize,'FontWeight','bold');
xlabel('LT(h)','FontSize',fsize,'FontWeight','bold');
set (gca,'xlim',[0,24],'ylim',[cn cx],'xtick',([0:2:24]),'FontSize',fsize,'box','on');

hold off;
