% Read in balloon data RAW files.
% DEW 10/30/2000

%Read in wind data file
clear

ymax = 5;
tmin = -30;
wmax = 40;
%read_parameters_rhb; %Reading parameter files flux conf from D:/data/lb/flux_eval_rhb_06
way_raw_sonde = 'Y:\Knorr\sondes\';
sfn=input('Input sonde  directory.... 080321.002 ','s');
hh=input('Input hour of launch ','s');
YY = sfn(1:2);
MM = sfn(3:4);
DD = sfn(5:6);
ptu=[way_raw_sonde, sfn '\' sfn(1:6) hh,'.ptu'];
wnd=[way_raw_sonde, sfn '\' sfn(1:6) hh,'.wnd'];

% Test for 23 launch  File name needs 1 day subtracted
if hh == '23'
    dd = str2num(DD);
    dd = dd -1;
    if dd < 10
        dds(1) = '0';
        dds(2) = num2str(dd);
        wnd(32:33) = dds(1:2)
        ptu(32:33) = dds(1:2)
    else
        wnd(32:33) = num2str(dd);
        ptu(32:33) = num2str(dd);
    end
end
%%dew
%wnd = '08032105.wnd'
%ptu = '08032105.wnd'

fid2=fopen(wnd, 'r');
display('Opened '); wnd
if fid2 > 0 % file found

    %
    %
    %Complete listing of Fine Structure file
    %File: C:\Metgraph\archive\001007.001\00100716.06e
    %
    %Started at       7 October 2000 16:06 UTC
    %
    %Station  :      313
    %Location :      33.24 N 128.89 W        10 m/MSL
    %
    %Sounding type : PTU GPS
    %RS-number     : 805103110
    %Sounding processor serial number: 0
    %
    %
    %Ground check omitted
    %
    %Clouds    : /////
    %Special 1 : /////
    %Special 2 : /////
    %
    %Telemetry noise level  : 80 %
    %
    %                 P      T      U
    %Accepted (%) :  35.3   64.2   60.7
    %Replaced (%) :  52.1   23.2   26.7
    %Rejected (%) :  12.6   12.6   12.6
    %
    %GPS data quality
    %    Valid raw wind ratio        :   73.0 %
    %    Valid raw wind on 5 ratio   :   73.0 %
    %    Invalid raw wind on 4 ratio :    1.0 %
    %    Poor PDOP ratio             :    1.0 %
    %    Unidentified signals ratio  :    0.0 %
    %    Unused                      :    0.0 %
    %    Mean satellite track count  :    4.7
    %    Mean track count valid      :    6.4
    %    Gaps over 1 min             :    0.0 min
    %
    %
    %Reason for termination : Manual stop
    %
    %   Time Pressure      Height   Speed  Direction   Significance flags
    %
    % min  s      hPa     gpm     m   m/s    deg  Automatic  Operator
    %
    %
    % Read in all header lines
    for n = 1:46
        h = fgetl(fid2);
        % print out Launch time
        if n == 6
            baltm = h;
            HR = baltm(31:32);
            if str2num(HR) <= 9
                hr(1) = '0';
                hr(2) = HR(2);
            else
                hr = HR;
            end
            mn = baltm(34:35);
            tit(1:6) = wnd(28:33);
            tit(7) = '.';
            tit(8:9) = hr;
            tit(10:11) = mn;
        end
    end
    % Read in data lines
    %   Time Pressure      Height   Speed  Direction   Significance flags
    % min  s      hPa     gpm     m   m/s    deg  Automatic  Operator
    %   0  0   1017.0      10 /////   3.6     12       dfv         dfv
    %   0  2   1015.5      23 /////   3.6     11
    %   0  4   1014.7      29 /////   3.7     11
    wt = 0;
    %for i=1:1500
    while feof(fid2)==0,
        clear data
        data = fgetl(fid2);
        data(1,41:43);
        data(1,14:16);
        if feof(fid2) == 1
            break
        end

        if (data(1,41:43) ~= '///') &  (data(1,14:16) ~= '///' )     % Test for flagged data
            wt = wt + 1; % count good wind pts
            balw(wt,1) = str2num(data(1,1:4)); %min
            balw(wt,3) = str2num(data(1,8:16));%p
            balw(wt,2) = str2num(data(1,17:24))/1000.;%ht msl km
            balw(wt,6) = str2num(data(1,37:43));%WD
            balw(wt,7) = str2num(data(1,31:36));%WS
            balw(wt,8) = str2num(data(1,43:49));%dz
        end
    end
    lb = size(balw);
    if lb(1) > 10
        %De-Bob balloon sounding
        bct = 1;
        % save 1st point
        bal1(bct,1) = balw(1,1); %min
        bal1(bct,4) = balw(1,2); %ht
        bal1(bct,3) = balw(1,3); %P
        bal1(bct,6) = balw(1,6); %wd
        bal1(bct,7) = balw(1,7); %ws
        bal1(bct,8) = balw(1,8); %dz
        blast = balw(1,2);
        for bb = 2:(length(balw))
            balw(bb,2);
            if balw(bb,2) > blast
                bct = bct + 1;
                bal1(bct,1) = balw(bb,1); %min
                bal1(bct,2) = balw(bb,2); %ht
                bal1(bct,3) = balw(bb,3); %P
                bal1(bct,6) = balw(bb,6); %wd
                bal1(bct,7) = balw(bb,7); %ws
                bal1(bct,8) = balw(bb,8); %dz
                blast = balw(bb,2);  % save last good pt assuming lower pt is good
            end
        end
        % put de-bobbed data back into array for later comparisons
        %clear bal

        baltmpw = balw; % save original data for testing
        balw = bal1; % copy de-bobbed file into original name
    end % test > 10 pts to de-bob
    if 1
        %Wind plot
        figure(2)

        subplot(1,2,1),plot(balw(:,6),balw(:,2),'k.');
        axis([0,360,0,ymax])
        xlabel('Dir True')
        ylabel('Ht km-MSL')
        %hold

        subplot(1,2,2),plot(balw(:,7),balw(:,2),'k.');
        title('')
        axis([0,wmax,0,ymax])
        title(tit)
        xlabel('Speed m/s')
        ylabel('Ht km-MSL')

        figure(4)
        subplot(1,2,1),plot(balw(:,6),balw(:,2),'k.');
        axis([0,360,0,max(balw(:,2))])
        xlabel('Dir True')
        ylabel('Ht km-MSL')
        %hold

        subplot(1,2,2),plot(balw(:,7),balw(:,2),'k.');
        axis([0,max(balw(:,7)),0,max(balw(:,2))])
        title(tit)
        xlabel('Speed m/s')
        ylabel('Ht km-MSL')
    end  % skip wind plot
    pause(1)

    fclose(fid2); % wnd file
end %file found

%%%%%%%%%%%%%%%%
%Read in Temp data file
fid3=fopen(ptu, 'r');
display('Opened '); ptu

if fid3 > 0 % PTU file found
    %
    %
    %Complete listing of Fine Structure file
    %File: C:\Metgraph\archive\001007.001\00100716.06e
    %
    %Started at       7 October 2000 16:06 UTC
    %
    %Station  :      313
    %Location :      33.24 N 128.89 W        10 m/MSL
    %
    %Sounding type : PTU GPS
    %RS-number     : 805103110
    %Sounding processor serial number: 0
    %
    %
    %Ground check omitted
    %
    %Clouds    : /////
    %Special 1 : /////
    %Special 2 : /////
    %
    %Telemetry noise level  : 80 %
    %
    %                 P      T      U
    %Accepted (%) :  35.3   64.2   60.7
    %Replaced (%) :  52.1   23.2   26.7
    %Rejected (%) :  12.6   12.6   12.6
    %
    %GPS data quality
    %    Valid raw wind ratio        :   73.0 %
    %    Valid raw wind on 5 ratio   :   73.0 %
    %    Invalid raw wind on 4 ratio :    1.0 %
    %    Poor PDOP ratio             :    1.0 %
    %    Unidentified signals ratio  :    0.0 %
    %    Unused                      :    0.0 %
    %    Mean satellite track count  :    4.7
    %    Mean track count valid      :    6.4
    %    Gaps over 1 min             :    0.0 min
    %
    %
    %Reason for termination : Manual stop
    %
    %   Time Pressure   Height  Temperature  RH   Dewp   Significance flags
    %
    % min  s      hPa      gpm     deg C      %  deg C  Automatic  Operator
    %
    % Read in all header lines
    for n = 1:46
        h = fgetl(fid3);
        % print out Launch time
        if n == 6
            h;
        end
    end
    % Read in data lines
    %   Time Pressure   Height  Temperature  RH   Dewp   Significance flags
    % min  s      hPa      gpm     deg C      %  deg C  Automatic  Operator
    %   0  0   1017.0       10      18.3     58    9.9 TU          TU
    %   0  2   1015.5       23      18.3     55    9.1
    %   0  4   1014.7       29      18.2     54    8.8
    tt = 0;
    %for i=1:length(balt)
    while feof(fid2)==0,
        clear data
        data = fgetl(fid3);
        if feof(fid3) == 1
            break
        end


        if (data(1,33:35) ~= '///') &  (data(1,14:16) ~= '///' )   % Test for flagged data
            tt = tt + 1; %count good data pts
            balt(tt,1) = str2num(data(1,1:4)); %min
            balt(tt,3) = str2num(data(1,8:16));%p
            balt(tt,4) = str2num(data(1,26:35));%Temp
            if (data(1,40:42)) ~= '///'  % good temp, but flagged RH
            balt(tt,5) = str2num(data(1,36:42));%RH
            else
            balt(tt,5) = -99;%RH
            end
            balt(tt,2) = str2num(data(1,17:25))/1000.;%ht msl km
        end
    end
    %De-Bob balloon sounding
    clear bal1 bb
    bb = 1;
    % Test for > 10 pts to de-bob
    lb = size(balw);
    if lb(1) > 10
        bct = 1;
        %save 1st point
        bal1(bct,1) = balt(bb,1); %min
        bal1(bct,2) = balt(1,2); %ht
        bal1(bct,3) = balt(1,3); %P
        bal1(bct,4) = balt(1,4); %T
        bal1(bct,5) = balt(1,5); %RH
        blast = balt(1,2);
        for bb = 2:(length(balt))
            if balt(bb,2) > blast
                bct = bct + 1;
                bal1(bct,1) = balt(bb,1); %min
                bal1(bct,2) = balt(bb,2); %ht
                bal1(bct,3) = balt(bb,3); %P
                bal1(bct,4) = balt(bb,4); %T
                bal1(bct,5) = balt(bb,5); %RH
                blast = balt(bb,2);  % save last good pt assuming lower pt is good
            end
        end
        % put de-bobbed data back into ba array for later comparisons
        %clear bal

        baltmpt = balt; % save original data for testing
        clear balt
        balt = bal1; % copy de-bobbed file into original name
        fnb = wnd(28:33);
    end % > 10 pts to de-bob
    %%%%%%%%%%%%%%%%%

    % plotting data
    %temp plot
    if 1
        figure(1)
        %subplot(1,1,1),plot(balt(:,4),balt(:,2),'b')
        subplot(1,2,1),plot(balt(:,4),balt(:,2),'k.');
        axis([tmin,5,0,ymax])
        xlabel('Temp C')
        ylabel('Ht km-MSL')
        %hold

        subplot(1,2,2),plot(balt(:,5),balt(:,2),'k.');
        title(tit)
        axis([0,100,0,ymax])
        xlabel('RH %')
        ylabel('Ht km-MSL')

        figure(3)
        %subplot(1,1,1),plot(balt(:,4),balt(:,2),'b')
        subplot(1,2,1),plot(balt(:,4),balt(:,2),'k.');
        axis([min(balt(:,4)),max(balt(:,4)),0,max(balt(:,2))])
        xlabel('Temp C')
        ylabel('Ht km-MSL')
        %hold

        subplot(1,2,2),plot(balt(:,5),balt(:,2),'k.');
        title(tit)
        axis([0,100,0,max(balt(:,2))])
        xlabel('RH %')
        ylabel('Ht km-MSL')

    end
    %fclose(fid3); % ptu file
end % ptu file found

if 1
    %      cd 'G:\Cruise_Data\EPIC 2001'
    figure(2)
    print_buffer = ['c:\data\icealot\data\sonde\Raw_Images\ICEALOT_Sonde_WND_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn '.png'];
    print('-dpng ', print_buffer);

    figure(1)
    print_buffer = ['c:\data\icealot\data\sonde\Raw_Images\ICEALOT_Sonde_TRH_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn '.png'];
    print('-dpng ', print_buffer);

    figure(3)
    print_buffer = ['c:\data\icealot\data\sonde\Raw_Images\ICEALOT_Sonde_TRH_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn 'a.png'];
    print('-dpng ', print_buffer);

    figure(4)
    print_buffer = ['c:\data\icealot\data\sonde\Raw_Images\ICEALOT_Sonde_WND_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn 'a.png'];
    print('-dpng ', print_buffer);
end

if 1
    % 4 panel plot for PMEL WEB
    figure(5)
    %subplot(1,1,1),plot(balt(:,4),balt(:,2),'b')
    subplot(2,2,1),plot(balt(:,4),balt(:,2),'k.');
    axis([min(balt(:,4)),max(balt(:,4)),0,max(balt(:,2))])
    xlabel('Temp C')
    ylabel('Ht km-MSL')
    %hold

    subplot(2,2,2),plot(balt(:,5),balt(:,2),'k.');
    title(tit)
    axis([0,100,0,max(balt(:,2))])
    xlabel('RH %')
    ylabel('Ht km-MSL')

    subplot(2,2,3),plot(balw(:,6),balw(:,2),'k.');
    axis([0,360,0,max(balw(:,2))])
    xlabel('Dir True')
    ylabel('Ht km-MSL')
    %hold

    subplot(2,2,4),plot(balw(:,7),balw(:,2),'k.');
    axis([0,max(balw(:,7)),0,max(balw(:,2))])
    title(tit)
    xlabel('Speed m/s')
    ylabel('Ht km-MSL')

    print_buffer = ['c:\data\icealot\data\sonde\Raw_Images\ICEALOT_Sonde_4panel_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn 'a.png'];
    print('-dpng ', print_buffer);
end

%Save temperature profile to a file for calculating pwv
if 1
    f=['c:\data\icealot\data\sonde\ICEALOT_sonde_' fnb(1:2) '_' fnb(3:4) '_' fnb(5:6) '_' hr mn '.txt'];
    flist=fopen(f,'w');
    save( f, 'balt', '-ascii','-tabs')
    fclose('all');
end;
