%P_CS_OBS_Plotmaker
%m-file to create plots of the pressure standard deviation, current speed, and 
%optical backscatter data from each station (A-E) at the NY Bight.  
%Labels/Location of text in general was modified for appearance using 
%Adobe Illustrator.  Station F not included since none of the 5 plotted
%variables are available from this station.
%
%Soupy Alexander, 11/18/2001
%
%Requires NetCDF toolbox, Soupy's tools "singleJD.m" and "gregaxdNM.m", and Rich's
%tool "gregaxd.m"

%Determine the file(s) corresponding to the station of interest

station = input('Enter the station of interest, use single quotes ');

if station == 'A';
    Pres_file = netcdf('5952var-a1h_d1.nc','nowrite');
    CS_file_1 = netcdf('5952v-a1h_d1.nc','nowrite');
    CS_file_2 = netcdf('5952v-a1h_d2.nc','nowrite');
    OBS_file_1 = netcdf('5952tcp-a1h.nc','nowrite');
    OBS_file_2 = netcdf('5952obs-a1h_d2.nc','nowrite');
elseif station == 'B';
    Pres_file = netcdf('5972var-a1h_d1.nc','nowrite');
    CS_file_1 = netcdf('5972v-a1h_d1.nc','nowrite');
    CS_file_2 = netcdf('5972v-a1h_d2.nc','nowrite');
    OBS_file_1 = netcdf('5972tc-a1h.nc','nowrite');
    OBS_file_2 = netcdf('5972obs-a1h_d2.nc','nowrite');
elseif station == 'C';
    Pres_file = netcdf('5992var-a1h.nc','nowrite');
    CS_file_1 = netcdf('5992v-a1h.nc','nowrite');
    OBS_file_1 = netcdf('5992tcp-a1h.nc','nowrite');
    OBS_file_2 = netcdf('5992obs-a1h_d2.nc','nowrite');
elseif station == 'D';
    Pres_file = netcdf('6012var-a1h_d1.nc','nowrite');
    CS_file_1 = netcdf('6012v-a1h_d1.nc','nowrite');
    CS_file_2 = netcdf('6012v-a1h_d2.nc','nowrite');
    OBS_file_1 = netcdf('6012tcp-a1h.nc','nowrite');
    OBS_file_2 = netcdf('6012obs-a1h_d2.nc','nowrite');
elseif station == 'E';
    Pres_file = netcdf('6032var-a1h.nc','nowrite');
    CS_file_1 = netcdf('6032v-a1h.nc','nowrite');
    OBS_file_1 = netcdf('6032tcp-a1h.nc','nowrite');
    OBS_file_2 = netcdf('6032obs-a1h_d2.nc','nowrite');
end

%Set the time limits for the plot
startjd = julian(1999,12,01,00);
endjd = julian(2000,4,20,00);

%Pull out the relavent data
Pres_t = Pres_file{'time'}(:);
Pres_t2 = Pres_file{'time2'}(:);
Pres_time = singleJD(Pres_t,Pres_t2);
Pres_data = Pres_file{'SDP_850'}(:);
Pres_depth = Pres_file{'depth'}(:);

CS_1_t = CS_file_1{'time'}(:);
CS_1_t2 = CS_file_1{'time2'}(:);
CS_1_time = singleJD(CS_1_t,CS_1_t2);
CS_1_data = CS_file_1{'CS_300'}(:);
CS_1_depth = CS_file_1{'depth'}(:);

%2nd Current Speed must be handled different since 2 stations have no 2nd depth
if station == 'A' | station == 'B' | station == 'D';
    CS_2_t = CS_file_2{'time'}(:);
    CS_2_t2 = CS_file_2{'time2'}(:);
    CS_2_time = singleJD(CS_2_t,CS_2_t2);
    CS_2_data = CS_file_2{'CS_300'}(:);
    CS_2_depth = CS_file_2{'depth'}(:);
end

OBS_1_t = OBS_file_1{'time'}(:);
OBS_1_t2 = OBS_file_1{'time2'}(:);
OBS_1_time = singleJD(OBS_1_t,OBS_1_t2);
OBS_1_data = OBS_file_1{'NEP_56'}(:);
OBS_1_depth = OBS_file_1{'depth'}(:);

OBS_2_t = OBS_file_2{'time'}(:);
OBS_2_t2 = OBS_file_2{'time2'}(:);
OBS_2_time = singleJD(OBS_2_t,OBS_2_t2);
OBS_2_data = OBS_file_2{'NEP_56'}(:);
OBS_2_depth = OBS_file_2{'depth'}(:);

%Plot the st. deviation of pressure in the top plot
subplot(5,1,1)
badpoints = find(Pres_data > 50);
Pres_data(badpoints) = NaN;
plot(Pres_time,Pres_data)
axis([startjd endjd 0 40])
ylabel('mb')
text(julian(2000,1,1,00),30,['St. Deviation of Pressure at ' num2str(Pres_depth) ' m'])
gregaxdNM(Pres_time,5)
set(gca,'ytick',[0 10 20 30 40],'yticklabel',[0 10 20 30 40])
title(['Station ' station])

%Plot the current speeds in subplots 2 & 3
subplot(5,1,2)
plot(CS_1_time,CS_1_data)
axis([startjd endjd 0 55])
ylabel('cm/s')
text(julian(2000,1,1,00),45,['Current Speed at ' num2str(CS_1_depth) ' m'])
gregaxdNM(CS_1_time,5)
set(gca,'ytick',[0 10 20 30 40 50],'yticklabel',[0 10 20 30 40 50])

if station == 'C' | station == 'E';
    subplot(5,1,3)
    axis([startjd endjd 0 55])
    ylabel('cm/s')
    text(julian(2000,1,1,00),45,'No second measure of near bottom current speed')
    gregaxdNM(CS_1_time,5)
    set(gca,'ytick',[0 10 20 30 40 50],'yticklabel',[0 10 20 30 40 50])
    set(gca,'box','on')
else    
    subplot(5,1,3)
    plot(CS_2_time,CS_2_data)
    axis([startjd endjd 0 55])
    ylabel('cm/s')
    text(julian(2000,1,1,00),45,['Current Speed at ' num2str(CS_2_depth) ' m'])
    gregaxdNM(CS_2_time,5)
    set(gca,'ytick',[0 10 20 30 40 50],'yticklabel',[0 10 20 30 40 50])
    set(gca,'box','on')
end

%Plot the optical backscatter in subplots 4 & 5
subplot(5,1,4)
plot(OBS_1_time,OBS_1_data)
axis([startjd endjd 0 10])
ylabel('v')
text(julian(2000,1,1,00),3,['Optical Backscatter at ' num2str(OBS_1_depth) ' m'])
gregaxdNM(OBS_1_time,5)
set(gca,'ytick',[0 2 4 6 8 10],'yticklabel',[0 2 4 6 8 10])

subplot(5,1,5)
plot(OBS_2_time,OBS_2_data)
axis([startjd endjd 0 10])
ylabel('v')
text(julian(2000,1,1,00),3,['Optical Backscatter at ' num2str(OBS_2_depth) ' m'])
gregaxd(OBS_2_time,5)
set(gca,'ytick',[0 2 4 6 8 10],'yticklabel',[0 2 4 6 8 10])

orient landscape