% Plot UM 100m logger data
% 150m sonic 1/2, 100m sonic 3/4, 50m sonic 5/6 even = N
% UMass sonic 2D Model 81000
% u is E/W + fm E
% v is N/S + fm N
% w - is up/dn + is up
% on N boom +v = 334
% on S boom +v = 154

clear
close all

read_parameters_UM; %Reading parameter files flux conf

pltit = 1;

for x = 19:19 % loop through days
    close all
    xs = sprintf('%02d',x)
    f150 = ['z:\BAO\Tower\UMASS\Processed\sonic_150s_200708' xs '.dat'];
    f100 = ['z:\BAO\Tower\UMASS\Processed\sonic_100s_200708' xs '.dat'];
    f50 =  ['z:\BAO\Tower\UMASS\Processed\sonic_050s_200708' xs '.dat'];
  
    f150
    f100
    f50
    
    st150 = load(f150);
    st100 = load(f100);
    st50 = load(f50);
    %load w:\data\umass\data\processed\SONIC_100m\sonic_1_20070825.dat  %1 day 1 sec averaged data
    %st150 = sonic_1_20070825;
    %load w:\data\umass\data\processed\SONIC_100m\sonic_3_20070825.dat  %1 day 1 sec averaged data
    %st100 = sonic_3_20070825;
    %load w:\data\umass\data\processed\SONIC_100m\sonic_5_20070825.dat  %1 day 1 sec averaged data
    %st50 = sonic_5_20070825;
    tit = 'BAO UMass 2007 S';
    % yy  mm dd hh mn ss    u       v     w       t     npt bat
    %2007 08 20 00 00 00    3.89   0.50  -0.40   27.10  32 11.1
    %2007 08 20 00 00 01    3.51   0.85  -0.17   27.09  32 11.1

    %%%%%%% PSD 10m data
    %1    2   3   4     5     6      7      8       9      10        11     12     13    14
    %ID  YR  YD  hhmm  T      RH     WSv   WSs     Dir    WDt       Bat    press  precip WDT
    %10,2007,165,1405,15.095,50.858,1.810,-999.000,238.61,-999.0000,13.515,843.
    %01,0.000,212.610

    st150(:,13) = datenum(st150(:,1),st150(:,2),st150(:,3),st150(:,4),st150(:,5),st150(:,6))-datenum(st150(:,1)-1,12,31);
    st100(:,13) = datenum(st100(:,1),st100(:,2),st100(:,3),st100(:,4),st100(:,5),st100(:,6))-datenum(st100(:,1)-1,12,31);
    st50(:,13) = datenum(st50(:,1),st50(:,2),st50(:,3),st50(:,4),st50(:,5),st50(:,6))-datenum(st50(:,1)-1,12,31);

    u50 = st50(:,7);
    v50 = st50(:,8);
    w50 = st50(:,9);
    t50 = st50(:,10);
    yys50 = num2str(st50(1,1));
    mms50 = num2str(st50(1,2));
    dds50 = num2str(st50(1,3));
    hhs50 = num2str(st50(1,4));
    mns50 = num2str(st50(1,5));
    sss50 = num2str(st50(1,6));

    u100 = st100(:,7);
    v100 = st100(:,8);
    w100 = st100(:,9);
    t100 = st100(:,10);
    yys100 = num2str(st100(1,1));
    mms100 = num2str(st100(1,2));
    dds100 = num2str(st100(1,3));
    hhs100 = num2str(st100(1,4));
    mns100 = num2str(st100(1,5));
    sss100 = num2str(st100(1,6));

    u150 = st150(:,7);
    v150 = st150(:,8);
    w150 = st150(:,9);
    t150 = st150(:,10);
    yys150 = num2str(st150(1,1));
    mms150 = num2str(st150(1,2));
    dds150 = num2str(st150(1,3));
    hhs150 = num2str(st150(1,4));
    mns150 = num2str(st150(1,5));
    sss150 = num2str(st150(1,6));


    %  Plot Temperatures
    figure(1)
    plot(st50(:,13),st50(:,10),'k')
    hold
    plot(st100(:,13),st100(:,10),'r')
    plot(st150(:,13),st150(:,10),'b')
    legend('50','100','150')
    ylabel('T (C)')
    xlabel('Year Day')
    mx = max(st50(:,13));
    xlim([st50(1,13) mx])
    ylim([min(st50(:,10))-5 max(st50(:,10)+5)])
    title(tit);
    if pltit
        print_buffer = [way_images_sonic 'UM_T_' yys150 mms150 dds150 '_50_100_150S.png'];
        print('-dpng ', print_buffer);
    end

    S50=sqrt(u50.*u50+v50.*v50);
    %dir=(atan2(-v,u)*180/pi) - 26;
    dir50=(atan2(u50,v50)*180/pi);
    for j = 1:length(dir50)
        if dir50(j) < 0
            dir50(j) = 360 + dir50(j);
        end
                 dir50(j) = dir50(j) - 206;
        if dir50(j) < 0
            dir50(j) = 360 + dir50(j);
        end
    end

    S100=sqrt(u100.*u100+v100.*v100);
    %dir=(atan2(-v,u)*180/pi) - 26;
    dir100=(atan2(u100,v100)*180/pi);
    for j = 1:length(dir100)
        if dir100(j) < 0  
            dir100(j) = 360 + dir100(j);
        end
                 dir100(j) = dir100(j) - 206;
        if dir100(j) < 0
            dir100(j) = 360 + dir100(j);
        end
    end


    S150=sqrt(u150.*u150+v150.*v150);
    %dir=(atan2(-v,u)*180/pi) - 26;
    dir150=(atan2(u150,v150)*180/pi);
    for j = 1:length(dir150)
        if dir150(j) < 0
            dir150(j) = 360 + dir150(j);
        end
         dir150(j) = dir150(j) - 206;
        if dir150(j) < 0
            dir150(j) = 360 + dir150(j);
        end
    end

    % on N boom +v = 334
    % on S boom +v = 154
    % 10m is on N boom so we need to add 26 degrees to directions

    figure(2)
    subplot(2,1,1),plot(st50(:,13),S50,'k')
    hold
    subplot(2,1,1),plot(st100(:,13),S100,'r')
    subplot(2,1,1),plot(st150(:,13),S150,'b')

    legend('50','100','150')
    ylabel('Spd (m/s)')
    xlabel('Year Day')
    xlim([st50(1,13) mx])
    ylim([0 20])
    title(tit);

    subplot(2,1,2),plot(st50(:,13),dir50,'k')
    hold
    subplot(2,1,2),plot(st100(:,13),dir100,'r')
    subplot(2,1,2),plot(st150(:,13),dir150,'b')
    ylabel('Direction')
    xlabel('Year Day')
    xlim([st50(1,13) mx])
    ylim([0 360])
    if pltit
        print_buffer = [way_images_sonic 'UM_SPDDIR_' yys150 mms150 dds150 '_50_100_150S.png'];
        print('-dpng ', print_buffer);
    end


    figure(3)
    plot(st50(:,13),st50(:,7:9))
    legend('u','v','w')
    ylabel('Spd (m/s)')
    xlabel('Year Day')
    xlim([st50(1,13) mx])
    ylim([-10 10])
    title('50');
    if pltit
        print_buffer = [way_images_sonic 'UM_uvw_' yys150 mms150 dds150 '_50S.png'];
        print('-dpng ', print_buffer);
    end
    figure(4)
    plot(st100(:,13),st100(:,7:9))
    legend('u','v','w')
    ylabel('Spd (m/s)')
    xlabel('Year Day')
    xlim([st100(1,13) mx])
    ylim([-10 10])
    title('100');
    if pltit
        print_buffer = [way_images_sonic 'UM_uvw_' yys150 mms150 dds150 '_100S.png'];
        print('-dpng ', print_buffer);
    end
    figure(5)
    plot(st150(:,13),st150(:,7:9))
    legend('u','v','w')
    ylabel('Spd (m/s)')
    xlabel('Year Day')
    xlim([st150(1,13) mx])
    ylim([-10 10])
    title('150m');
    if pltit
        print_buffer = [way_images_sonic 'UM_uvw_' yys150 mms150 dds150 '_150S.png'];
        print('-dpng ', print_buffer);
    end

    figure(6)
    plot(st50(:,13),st50(:,11),'k.')
    hold
    plot(st100(:,13),st100(:,11),'r.')
    plot(st150(:,13),st150(:,11),'bo')
    legend('50','100','150','Location','Southwest')
    title('BAO 2007 UMass sonic npts');
    ylim([0 34]);
    if pltit
        print_buffer = [way_images_sonic 'UM_npts_' yys150 mms150 dds150 '_50_100_150S.png'];
        print('-dpng ', print_buffer);
    end

    
    clear st150 st100 st50
end   % day loop