%Read WHOI and SCS Event and PSD files file
% Each Event is stopped and started to creat daily files (UTC)
% DEW LB 7/10/2013
disp('Need to run read_epicz to create 5 min all file');
disp('Read_SCS_WHOTS_2013');

clear
close all

cruise='WHOTS_2013';     % Cruise name, for example 'STRATUS_2009', 'TOTO_1969', etc...
ship='Hiialakai';                % Research vessel
data_drive='D:\';
Lacie_drive='G:\';
plotit = 1;

ev = input('Input SCS event number to read (011) ','s');

way_proc_data_flux = [data_drive cruise  '\' ship '\flux\Processed\'];     % Processed data path
way_images_flux = [data_drive cruise  '\' ship '\flux\Raw_Images\' ];

%Ludo's read
clear dir
dr=dir(['D:\WHOTS_2013\Hiialakai\SCS\Raw\WHOTS-OBS_' ev '.elg']);
e=dr(1).name;
disp(['Reading scs file ',e]);
flist=fopen(['D:\WHOTS_2013\Hiialakai\SCS\Raw\WHOTS-OBS_' ev '.elg']);
temp=[];scsdata=[];
if flist>0,                                 %if the file exists,
    while ~feof(flist)
        %             temp=textscan(flist,'%f %f %f %*c %f %f %f %*c %f %f %*c %f %f %*c %f %f %*c %f %f %*c %f %f %f %f %f  %*[^\n]','delimiter', ' ', 'treatAsEmpty', {'SBE-21','****','FWD','Temp'}, 'headerlines', 1,  'emptyvalue', NaN);
        temp=textscan(flist,'%2f/%2f/%4f %2f:%2f:%2f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %*4c %f %f %f %*6c %f %*c %f %*c %*4c %f %f %f  %*[^\n]','delimiter', ', ', 'headerlines', 1,  'emptyvalue', NaN);
        scsdata=[scsdata cell2mat(temp)'];
    end;  %end while
end;      %     end if list
scs = scsdata;

% a=['D:\WHOTS_2013\Hiialakai\SCS\Raw\WHOTS-OBS_' bds '.dat'];
% scs = load(a);
% 1-3 Date
% 4-6 Time UTC
% 7 SAMOS-Lat-VALUE
% 8 SAMOS-Lon-VALUE
% 9 SAMOS-Gyro-VALUE
% 10 SAMOS-COG-VALUE
% 11 SAMOS-SOG-VALUE
% 12 SAMOS-AirTemp-VALUE                 **** RMYoung deg C
% 13 SAMOS-Barometer-VALUE             Bad RMYoung
% 14 SAMOS-RelHumidity-VALUE            **** RMY %
% 15 SAMOS-RelWind-Dir-VALUE
% 16 SAMOS-RelWind-Spd-VALUE
% 17 SAMOS-TrueWind-Dir-VALUE           *** Prop-Vane deg T
% 18 SAMOS-TrueWind-Spd-VALUE              *** Prop-Vane kts
% 19 SAMOS-TSG-Temp-VALUE              **** SST deg C
% 20 SAMOS-TSG-Cond-VALUE
% 21 SAMOS-TSG-Sal-VALUE
% 22 SAMOS-TSG-SV-VALUE
% 23 WINCH-PAYOUT
% 24 WINCH-SPD
% 25 WINCH-TENSION
% 26 WIND-85004-Ultrasonic-Anemometer   Rel dir deg
% 27 WIND-85004-Ultrasonic-Anemometer   Rel spd kts
% 28 Depth-Meters
% 29 Vaisala-Barometer--               Good Station pressure mb
% 30 Baro-Pressure  **Bad RMY

% Calculate Decimal year day
yd_scs=datenum(scsdata(3,:),scsdata(1,:),scsdata(2,:),scsdata(4,:),scsdata(5,:),scsdata(6,:))-datenum(scsdata(3,1),1,0);   %Day-Of-Year
bd = fix(yd_scs(1));
ed = bd + 1;

bds = sprintf('%03d',bd);
dds = sprintf('%02d',scsdata(2,1));
mms = sprintf('%02d',scsdata(1,1));


g=['D:\WHOTS_2013\Hiialakai\flux\Processed\WHOTS_2013_proc_met_1min_2013_' mms '_' dds '_' bds '.txt'];
met=load(g)';
% 1  jd_pc;  %decimal JD
% 2  pirm; %
% 3  pspm; %
% 4  Tc1;  %
% 5  Td1;  %
% 6  Tsea;
% 7  Tvais;    %
% 8  Rhvais;   %
% 9  org;  %
% 10  org_carrier; %
% 11  aspir_on;    %
% 12  press;    %

s=['D:\WHOTS_2013\Hiialakai\flux\Processed\WHOTS_2013_proc_son_1min_2013_' mms '_' dds '_' bds '.txt'];
son=load(s)';
% jd_son;	%decimal JD
% U;
% V;
% w;
% Tsonic;
% dir;

%  Load 5 min all file
b=['D:\WHOTS_2013\Hiialakai\flux\Processed\WHOTS_2013_PSD_flux_5min_all.txt'];
met5=load(b)';
% jdy  1  julian day to Matlab format
% U   2   true wind,PSD sonic (m/s)
% dir   3   true wind direction, PSD sonic (deg)
% tsnk   4   sea snake temperature,PSD, 0.05 m depth (C)
% ta=x(:,7);%air temperature, PSD (C) (z=15 m)
% press   31   PSD atmospheric pressure, mb
% RH  32    PSD RHvais, %

% Calculating Sea Level Pressure using Barometer height
% ESRL ht = 14.38m
% ship ht = 12.02m
% NOAA Tech training

%ESRL
Pz1 = met(12,:); % Barometer pressure
Ra = 287.05;  % Constant J/kgK
T1 = met(7,:);  % Temperature K
g = 9.81;  % gravity m/s
Z1 = 14.38;  % height of presure sensor in m
%Ps = ssealevel pressure
s1 = (9.81*Z1);
s2 = (Ra*(T1+274.15));
Ps1 = exp(s1./s2).*Pz1;
clear s1 s2
%Ship
Pz2 = scs(29,:); % Barometer pressure
Ra = 287.05;  % Constant J/kgK
T2 = scs(12,:);  % Temperature K
g = 9.81;  % gravity m/s
Z2 = 12.02;  % height of presure sensor in m
%Ps = ssealevel pressure

s1 = (9.81*Z2);
s2 = (Ra*(T2+274.15));
Ps2 = exp(s1./s2).*Pz2;

mPs2 = mean(Ps2);
minp = fix(min(Ps2));
maxp = fix(max(Ps2));
pp = find(Ps1 > mPs2-3 );
mPs1 = mean(Ps1(pp));

% figure
% plot(yd_scs,scs(29,:),'r')  %Vaisala-Barometer-
% hold
% plot(met(1,:),met(12,:),'k.')
% plot(yd_scs,Ps2,'g')  %Vaisala-Barometer-
% plot(met(1,pp),Ps1(pp),'b')
% ylabel('Pressure (mb)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Sea Level Pressure');
% grid
% xlim([ bd ed ]);
% ylim([ minp-2 maxp+2]);
% legend('SCS station','PSD station','SCS sea level','PSD sea level','location','best');
%
% if plotit
%     print_buffer = [way_images_flux cruise '_SCS_PSD_Psealvl_' bds '.png'];
%     print('-dpng ', print_buffer);
% end
% %  Plot RMY and other pressures
% if 0
%     figure
%     plot(yd_scs,scs(:,35),'b') %Baro-Pressure
%     hold
%     plot(yd_scs,scs(:,34),'r') %Vaisala-Barometer-
%     plot(yd_scs,scs(:,13),'g')  %Samos-Barometer-Value
%     plot(met(1,:),met(12,:),'k.')  % PSD
%     ylabel('Pressure (mb)');
%     xlabel('Year Day (UTC)')
%     title('WHOTS 2013 Pressure');
%     grid
%     xlim([ bd ed ]);
%     legend('RMY','SCS','RMY','PSD','location','best');
%
%     if plotit
%         print_buffer = [way_images_flux cruise '_SCS_PSD_Pall_' bds '.png'];
%         print('-dpng ', print_buffer);
%     end
%
% end
%
% figure
% plot(yd_scs,scs(19,:),'b') % Samos-TSG-Temp
% hold
% plot(met(1,:),met(6,:),'k')  % Seasnake
% ylabel('SST (C)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 SST');
% grid
% xlim([ bd ed ]);
% legend('SCS','PSD','location','best');
%
% if plotit
%     print_buffer = [way_images_flux cruise '_SCS_PSD_SST_' bds '.png'];
%     print('-dpng ', print_buffer);
% end
%
% figure
% subplot(2,1,1),plot(yd_scs,scs(12,:),'b') % Samos-AirTemp
% hold
% subplot(2,1,1),plot(met(1,:),met(7,:),'r')  %PSD Vai Temp
% subplot(2,1,1),plot(son(1,:),son(5,:),'k')  %PSD Son Temp
% ylabel('Air T (C)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Air Temp');
% grid
% xlim([ bd ed ]);
% legend('SCS','PSD Vai','PSD Son','location','best','orientation','horizontal');
%
% subplot(2,1,2),plot(yd_scs,scs(14,:),'b') % Samos-RH-Temp
% hold
% subplot(2,1,2),plot(met(1,:),met(8,:),'k') % PSD RH
% ylabel('RH(%)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 RH');
% ylim([ 0 100]);
% grid
% xlim([ bd ed ]);
% legend('SCS','PSD Vai','location','best');
%
% if plotit
%     print_buffer = [way_images_flux cruise '_SCS_PSD_TRH_' bds '.png'];
%     print('-dpng ', print_buffer);
% end
%
%
% sonspd = sqrt((son(2,:).^2)+(son(3,:).^2));
%
% figure
% subplot(2,1,1),plot(yd_scs,scs(17,:),'b.') % Samos-TrueWind-Dir
% hold
% subplot(2,1,1),plot(met5(1,:),met5(3,:),'k.') % PSD Son
% ylabel('Dir (T)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Dir True');
% ylim([ 0 360]);
% grid
% xlim([ bd ed ]);
% legend('SCS P-V','PSD Son','location','best');
%
% subplot(2,1,2),plot(yd_scs,scs(18,:)*.55,'b') % Samos TruWind Speed
% hold
% subplot(2,1,2),plot(met5(1,:),met5(2,:),'k.') % PSD TrueWind Speed
% ylabel('Spd (m/s)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Spd True');
% ylim([ 0 25]);
% xlim([ bd ed ]);
% grid
%
% if plotit
%     print_buffer = [way_images_flux cruise '_SCS_PSD_DirSpdT_' bds '.png'];
%     print('-dpng ', print_buffer);
% end
%
%
% figure
% subplot(2,1,1),plot(yd_scs,scs(15,:),'b.') % Samos WindDir-Rel
% hold
% plot(yd_scs,scs(26,:),'r.')  % Ultra Sonic
% plot(son(1,:),son(6,:),'k.')  % PSD Son
% ylabel('Dir (Rel)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Dir Rel');
% ylim([ 0 360]);
% grid
% xlim([ bd ed ]);
% legend('SCS P-V','SCS Son','PSD Son','location','best');
%
% subplot(2,1,2),plot(yd_scs,scs(16,:)*.55,'b') % Samos-WindRel Rel
% hold
% plot(yd_scs,scs(27,:)*.55,'r') % Ultra Sonic
% plot(son(1,:),sonspd,'k')  % PSD SOn
% ylabel('Spd (m/s)');
% xlabel('Year Day (UTC)')
% title('WHOTS 2013 Spd Rel');
% ylim([ 0 25]);
% grid
% xlim([ bd ed ]);
%
% if plotit
%     print_buffer = [way_images_flux cruise '_SCS_PSD_DirSpdR_' bds '.png'];
%     print('-dpng ', print_buffer);
% end

% load and plot WHOI data Times are UTC 1 hr data...wind dir are
% oceanographic 180 out from met
%New WHOI files
% whots 10 SN 15 13-Jul-2013 20:18
% winds are reported as wind TO east and wind TO north
% magnetic correction applied: 9.81
% Please see http://uop.whoi.edu/UOP_realtime_files.html
%year month  day yday hr    BPR   ATMP  SST    WNDE   WNDN     HRH     SWR   LWR    RN   CND
%                    UTC    mb    DegC  DegC    m/s    m/s      %%     W/m^2 W/m^2   mm   S/m
%   2013    7   11  192  6   1021  25.55  25.90   -9.3   -3.2     75      2    388      0   5.43

%
%WHOTS 9  is 2012 Logger 9 and 10 are 1 and 2  Yd day is 2nd year subtract
%365
%WHOTS 10  is 2013 Logger 15 and 8 are 1 and 2

flistw=fopen(['D:\WHOTS_2013\Hiialakai\WHOI\whots_10_lgr_15.txt']);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
ic = 0;
while ~feof(flistw)
    ic = ic + 1;
    tt = fgetl(flistw);%if the file exists,
    [data] = sscanf(tt,'%f');
    whoi1015(ic,1) = data(1);  % YR
    whoi1015(ic,2) = data(2);  % Mon
    whoi1015(ic,3) = data(3);  % Day
    whoi1015(ic,4) = data(4);  % YD
    whoi1015(ic,5) = data(5);  % HR
    whoi1015(ic,6) = data(6); % P mb
    whoi1015(ic,7) = data(7);  % Ta C
    whoi1015(ic,8) = data(8);  % SST C
    whoi1015(ic,9) = data(9);  %WNDE
    whoi1015(ic,10) = data(10);% WNDN
    whoi1015(ic,11) = data(11);  % RH
    whoi1015(ic,12) = data(12);  % SW
    whoi1015(ic,13) = data(13);  % LW
    whoi1015(ic,14) = data(14);  % Rain
    whoi1015(ic,15) = data(15);  % CND
    whoi1015(ic,16) = (whoi1015(ic,4) + whoi1015(ic,5)./24);  % spd
    whoi1015(ic,17) = sqrt((whoi1015(ic,9)^2) + (whoi1015(ic,10)^2));  % Decimal day
    drad = atan2(whoi1015(ic,9),whoi1015(ic,10));
    whoi1015(ic,18) = drad*(180/pi); % dir
    if whoi1015(ic,18) < 0
        whoi1015(ic,18) = whoi1015(ic,18) + 360;
    end
end;  %end while

flistw=fopen(['D:\WHOTS_2013\Hiialakai\WHOI\whots_10_lgr_8.txt']);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
ic = 0;
while ~feof(flistw)
    ic = ic + 1;
    tt = fgetl(flistw);%if the file exists,
    [data] = sscanf(tt,'%f');
    whoi108(ic,1) = data(1);  % YR
    whoi108(ic,2) = data(2);  % Mon
    whoi108(ic,3) = data(3);  % Day
    whoi108(ic,4) = data(4);  % YD
    whoi108(ic,5) = data(5);  % HR
    whoi108(ic,6) = data(6); % P mb
    whoi108(ic,7) = data(7);  % Ta C
    whoi108(ic,8) = data(8);  % SST C
    whoi108(ic,9) = data(9);  %WNDE
    whoi108(ic,10) = data(10);% WNDN
    whoi108(ic,11) = data(11);  % RH
    whoi108(ic,12) = data(12);  % SW
    whoi108(ic,13) = data(13);  % LW
    whoi108(ic,14) = data(14);  % Rain
    whoi108(ic,15) = data(15);  % CND
    whoi108(ic,16) = whoi108(ic,4) + whoi108(ic,5)./24;  % spd
    whoi108(ic,17) = sqrt((whoi108(ic,9)^2) + (whoi108(ic,10)^2));  % Decimal day
    drad = atan2(whoi108(ic,9),whoi108(ic,10));
    whoi108(ic,18) = drad*(180/pi); % dir
    if whoi108(ic,18) < 0
        whoi108(ic,18) = whoi108(ic,18) + 360;
    end
end;  %end while

flistw=fopen(['D:\WHOTS_2013\Hiialakai\WHOI\whots_9_lgr_10.txt']);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
ic = 0;
while ~feof(flistw)
    ic = ic + 1;
    tt = fgetl(flistw);%if the file exists,
    [data] = sscanf(tt,'%f');
    whoi910(ic,1) = data(1);  % YR
    whoi910(ic,2) = data(2);  % Mon
    whoi910(ic,3) = data(3);  % Day
    whoi910(ic,4) = data(4);  % YD
    whoi910(ic,5) = data(5);  % HR
    whoi910(ic,6) = data(6); % P mb
    whoi910(ic,7) = data(7);  % Ta C
    whoi910(ic,8) = data(8);  % SST C
    whoi910(ic,9) = data(9);  %WNDE
    whoi910(ic,10) = data(10);% WNDN
    whoi910(ic,11) = data(11);  % RH
    whoi910(ic,12) = data(12);  % SW
    whoi910(ic,13) = data(13);  % LW
    whoi910(ic,14) = data(14);  % Rain
    whoi910(ic,15) = data(15);  % CND
    whoi910(ic,16) = (whoi910(ic,4) + whoi910(ic,5)./24)-365;  % YD
    whoi910(ic,17) = sqrt((whoi910(ic,9)^2) + (whoi910(ic,10)^2));  % Decimal day
    drad = atan2(whoi910(ic,9),whoi910(ic,10));
    whoi910(ic,18) = drad*(180/pi); % dir
    if whoi910(ic,18) < 0
        whoi910(ic,18) = whoi910(ic,18) + 360;
    end
end;  %end while

flistw=fopen(['D:\WHOTS_2013\Hiialakai\WHOI\whots_9_lgr_9.txt']);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
fgets(flistw);
ic = 0;
while ~feof(flistw)
    ic = ic + 1;
    tt = fgetl(flistw);%if the file exists,
    [data] = sscanf(tt,'%f');
    whoi99(ic,1) = data(1);  % YR
    whoi99(ic,2) = data(2);  % Mon
    whoi99(ic,3) = data(3);  % Day
    whoi99(ic,4) = data(4);  % YD
    whoi99(ic,5) = data(5);  % HR
    whoi99(ic,6) = data(6); % P mb
    whoi99(ic,7) = data(7);  % Ta C
    whoi99(ic,8) = data(8);  % SST C
    whoi99(ic,9) = data(9);  %WNDE
    whoi99(ic,10) = data(10);% WNDN
    whoi99(ic,11) = data(11);  % RH
    whoi99(ic,12) = data(12);  % SW
    whoi99(ic,13) = data(13);  % LW
    whoi99(ic,14) = data(14);  % Rain
    whoi99(ic,15) = data(15);  % CND
    whoi99(ic,16) = (whoi99(ic,4) + whoi99(ic,5)./24)-365;  % yd
    whoi99(ic,17) = sqrt((whoi99(ic,9)^2) + (whoi99(ic,10)^2));  % Decimal day
    drad = atan2(whoi99(ic,9),whoi99(ic,10));
    whoi99(ic,18) = drad*(180/pi); % dir
    if whoi99(ic,18) < 0
        whoi99(ic,18) = whoi99(ic,18) + 360;
    end
end;  %end while
if plotit
    print_buffer = [way_images_flux cruise '_WHOI_DirT_' bds '.png'];
    print('-dpng ', print_buffer);
end

figure
plot(met5(1,:),met5(2,:),'k.') % PSD Son
hold
plot(whoi99(:,16),whoi99(:,17),'go') %
plot(whoi910(:,16),whoi910(:,17),'c*') %
plot(whoi1015(:,16),whoi1015(:,17),'rd') %
plot(whoi108(:,16),whoi108(:,17),'b+') %

ylabel('Spd (m/s)');
xlabel('Year Day (UTC)')
title('WHOTS 2013 Spd True');
ylim([ 0 20]);
grid
xlim([ bd ed ]);
legend('PSD','W9-1','W9-2','W10-1','W10-2','location','best','orientation','horizontal');
if plotit
    print_buffer = [way_images_flux cruise '_WHOI_SpdT_' bds '.png'];
    print('-dpng ', print_buffer);
end

figure
plot(met5(1,:),met5(3,:),'k') % PSD Son
hold
plot(whoi99(:,16),whoi99(:,18)-180,'go') %
plot(whoi910(:,16),whoi910(:,18)-180,'c*') %
plot(whoi1015(:,16),whoi1015(:,18)-180,'rd') %
plot(whoi108(:,16),whoi108(:,18)-180,'b+') %

ylabel('Dir (T)');
xlabel('Year Day (UTC)')
title('WHOTS 2013 Dir True');
ylim([ 0 180]);
xlim([ bd ed ]);
grid
legend('PSD','W9-1','W9-2','W10-1','W10-2','location','best','orientation','horizontal');