disp('openophir_csp')
% Inputs data from ASCII data files
%  Note:  hhh=hour

A0_ophir=1.051;
A1_ophir=-0.0017;
k_ophir=0.051897;
b_ophir=0.734214;

fclose('all'); 
if hhh<10,
   dfx2=['0',num2str(hhh)];
else,
   dfx2=num2str(hhh);
end;

if ddd<10,
   dfx1=['00',num2str(ddd)];
elseif ddd<100,
   dfx1=['0',num2str(ddd)];
else,
   dfx1=num2str(ddd);
end;

dfx=[dfx1 dfx2];

%   Open OPHIR file

dfl=['/cdrom/P4_',dfx,'.DAT'];%disp(dfl)
flist=fopen(dfl,'r');
headr=fgetl(flist);%disp(headr)
headr=fgetl(flist);
%disp('     ###############     Ophir Evaluation     ###############     ')
z=fscanf(flist,'%g,%g,%g,%g,%g',[5,inf]);
if ~isempty(z);
   Tairlong=interp(Tairv,fix(length(z)/length(Tairv)));
   while length(Tairlong)<length(z),
      Tairlong=[Tairlong Tairlong(length(Tairlong))];
   end;   
   Tcan=Tairlong+12;
   a_ophir=A0_ophir*(1+A1_ophir.*Tcan);
   z=z(:,1:length(z)-1);
   t1=z(1,:)-z(1,1);t2=z(2,:)/10000-z(2,1)/10000;
   Oclear=z(4,:);Oabsorb=z(3,:);
   tophir=t1+t2;clear t1 t2;

   firstq=fix(length(Oclear)*1/6);
   seconq=fix(length(Oclear)*2/6);
   thirdq=fix(length(Oclear)*3/6);
   fourtq=fix(length(Oclear)*4/6);
   fifthq=fix(length(Oclear)*5/6);
   sixthq=fix(length(Oclear)*6/6);

   for pp=1:6;
      if pp==1,zz=1:firstq;end;
      if pp==2,zz=firstq+1:seconq;end;
      if pp==3,zz=seconq+1:thirdq;end;
      if pp==4,zz=thirdq+1:fourtq;end;
      if pp==5,zz=fourtq+1:fifthq;end;
      if pp==6,zz=fifthq+1:sixthq;end;
      
      Oclear_10=detrend(Oclear(zz));Oabsorb_10=detrend(Oabsorb(zz));
      med_Oabsorb=median(Oabsorb_10);med_Oclear=median(Oclear_10);
      sort_Oabsorb=sort(Oabsorb_10);sort_Oclear=sort(Oclear_10);
      sig_Oabsorb=0.5*(sort_Oabsorb(fix(2/3*length(sort_Oabsorb)))-sort_Oabsorb(fix(1/3*length(sort_Oabsorb))));
      sig_Oclear=0.5*(sort_Oclear(fix(2/3*length(sort_Oclear)))-sort_Oclear(fix(1/3*length(sort_Oclear))));
      if abs(Oclear_10(1)-med_Oclear) > 7*sig_Oclear,
         Oclear_10(1)=med_Oclear;
      end;
      if abs(Oabsorb_10(1)-med_Oabsorb) > 7*sig_Oabsorb,
         Oabsorb_10(1)=med_Oabsorb;
      end;
      for kk=2:length(Oclear_10);
         if abs(Oclear_10(kk)-med_Oclear) > 7*sig_Oclear,
            Oclear_10(kk)=Oclear_10(kk-1);
         end;   
         if abs(Oabsorb_10(kk)-med_Oabsorb) > 7*sig_Oabsorb,
            Oabsorb_10(kk)=Oabsorb_10(kk-1);
         end;
      end;

      rc=(Oabsorb_10+mean(Oabsorb(zz)))./(Oclear_10+mean(Oclear(zz)));
      Ophir_stat_clear=std(Oclear_10);
      Ophir_stat_absorb=std(Oabsorb_10);

      if pp==1,
         stdC_1=Ophir_stat_clear;
         stdA_1=Ophir_stat_absorb;
         medC_1=median(Oclear(zz));
         medA_1=median(Oabsorb(zz));
      elseif pp==2,
         stdC_2=Ophir_stat_clear;
         stdA_2=Ophir_stat_absorb;
         medC_2=median(Oclear(zz));
         medA_2=median(Oabsorb(zz));
      elseif pp==3,
         stdC_3=Ophir_stat_clear;
         stdA_3=Ophir_stat_absorb;
         medC_3=median(Oclear(zz));
         medA_3=median(Oabsorb(zz));
      elseif pp==4,
         stdC_4=Ophir_stat_clear;
         stdA_4=Ophir_stat_absorb;
         medC_4=median(Oclear(zz));
         medA_4=median(Oabsorb(zz));
      elseif pp==5,
         stdC_5=Ophir_stat_clear;
         stdA_5=Ophir_stat_absorb;
         medC_5=median(Oclear(zz));
         medA_5=median(Oabsorb(zz));
      elseif pp==6,
         stdC_6=Ophir_stat_clear;
         stdA_6=Ophir_stat_absorb;
         medC_6=median(Oclear(zz));
         medA_6=median(Oabsorb(zz));
      end;
      abshum=(((-log(rc./a_ophir(zz)))./k_ophir).^(1/b_ophir)).*((Tairlong(zz)+273.15)/273.15).^1.54;
      abshum_=detrend(abshum);
      llsx=length(abshum);
      if llsx > 65535,
         specsize=65536;
      elseif llsx > 32767,
         specsize=32768;   
      elseif llsx > 16383,
         specsize=16384;   
      elseif llsx > 8191,
         specsize=8192;   
      elseif llsx > 4095,
         specsize=4096;   
      elseif llsx > 2047,
         specsize=2048;   
      elseif llsx > 1023,
         specsize=1024;   
      elseif llsx > 511,
         specsize=512;   
      elseif llsx > 255,
         specsize=256;   
      elseif llsx > 127,
         specsize=128;   
      elseif llsx > 63,
         specsize=64;   
      elseif llsx > 31,
         specsize=32;   
      elseif llsx > 15,
         specsize=16;   
      end;

      if pp==1,
         qq_1=cov(abshum_);
         [Ch_1,Fo]=psd(abshum_,specsize,20);[Chs_1,Foph]=specsmooth(Ch_1,20);%figure;loglog(Foph,Chs_1.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 0']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_1.*Foph,Foph(20:length(Foph)),(Chs_1(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 00'])
      elseif pp==2,
         qq_2=cov(abshum_);
         [Ch_2,Fo]=psd(abshum_,specsize,20);[Chs_2,Foph]=specsmooth(Ch_2,20);%figure;loglog(Foph,Chs_2.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 10']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_2.*Foph,Foph(20:length(Foph)),(Chs_2(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 10'])
      elseif pp==3,
         qq_3=cov(abshum_);
         [Ch_3,Fo]=psd(abshum_,specsize,20);[Chs_3,Foph]=specsmooth(Ch_3,20);%figure;loglog(Foph,Chs_3.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 20']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_3.*Foph,Foph(20:length(Foph)),(Chs_3(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 20'])
      elseif pp==4,
         qq_4=cov(abshum_);
         [Ch_4,Fo]=psd(abshum_,specsize,20);[Chs_4,Foph]=specsmooth(Ch_4,20);%figure;loglog(Foph,Chs_4.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 30']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_4.*Foph,Foph(20:length(Foph)),(Chs_4(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 30'])
      elseif pp==5,
         qq_5=cov(abshum_);
         [Ch_5,Fo]=psd(abshum_,specsize,20);[Chs_5,Foph]=specsmooth(Ch_5,20);%figure;loglog(Foph,Chs_5.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 40']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_5.*Foph,Foph(20:length(Foph)),(Chs_5(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 40'])
      elseif pp==6,
         qq_6=cov(abshum_);
         [Ch_6,Fo]=psd(abshum_,specsize,20);[Chs_6,Foph]=specsmooth(Ch_6,20);%figure;loglog(Foph,Chs_6.*Foph);title(['JD',dfx1,' Hr',dfx2,' UTC  Spectrum Ophir abshum 50']);xlabel('freq in Hz');ylabel('F*S(abshum)');
         figure;loglog(Foph,Chs_6.*Foph,Foph(20:length(Foph)),(Chs_6(20)/Foph(20)^(-2/3))*Foph(20:length(Foph)).^(-2/3));title(['humidity  day ',num2str(ddd),'    hour ',num2str(hhh),'    min 50'])
      end;
      clear abshum_;
   end;

   %if length(abshum) > length(sx),
   %   abshum=abshum(1:length(sx));
   %end;
else,
   abshum=(1:72000).*0;
   [Chs_1,Foph]=psd(abshum,1024,20);
   [Chs_2,Foph]=psd(abshum,1024,20);
   [Chs_3,Foph]=psd(abshum,1024,20);
   [Chs_4,Foph]=psd(abshum,1024,20);
   [Chs_5,Foph]=psd(abshum,1024,20);
   [Chs_6,Foph]=psd(abshum,1024,20);
end;

clear z x sort_Oclear sort_Oabsorb sig_Oclear sig_Oabsorb rc pp mm med_Oclear med_Oabsorb kk zz headr
clear flist dfx dfl ans abshum_mean a Ophir_stat_clear Ophir_stat_absorb Oclear_10 Oabsorb_10
clear Ch llsx specsize a_ophir A0_ophir A1_ophir k_ophir b_ophir rc abshum
clear Oabsorb Oclear sixthq fifthq fourtq thirdq seconq firstq tophir tsons
clear Tairlong Tcan a_ophir

fclose('all');