%PressurePlot.m
%
%Program which plots low-passed wind stres, hourly averaged pressure and 
%low-passed pressure from a given station. (Normalized by the mean)
%
%Calls plfilt.m,ridfill_nan.m 
%
%Soupy Alexander
%12/7/2001
%
%Modified 3/7/02 to add wind stress

%Set Limits for plots for consistency
pmin = -200;
pmax = 200;

stationID = input('Choose a station: ');

%Select the proper date files for the selected station.
if stationID == 'A';
    presfile = netcdf('5952tcp-a1h.nc');
elseif stationID == 'B';
    presfile = netcdf('5972p-a1h.nc');
elseif stationID == 'C';
    presfile = netcdf('5992tcp-a1h.nc');
elseif stationID == 'D';
    presfile = netcdf('6012tcp-a1h.nc');
elseif stationID =='E';
    presfile = netcdf('6032tcp-a1h.nc');
elseif stationID == 'F';
    presfile = netcdf('6042pt-a1h.nc');
else
    error('Not a proper station ID.  Be sure to use single quotes and capital letters');
end

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

%Extract useful information from the file
t = presfile{'time'}(:);
t2 = presfile{'time2'}(:);

time_h = singleJD(t,t2);

%Station F data is from a SBE-39; measures in dbar which need to be
%converted to mbar to match the other stations.
if stationID == 'F';
    pres_h = presfile{'P_1'}(:);
    pres_h = 100.*pres_h;
else
    pres_h = presfile{'P_4023'}(:);
end

[pres_h] = ridfill_nan(pres_h);

[pres_lp,time_lp] = plfilt(pres_h,time_h,6);


%Wind Stress

%Wind Stress from Ambrose Lighthouse
ambrose = netcdf('alsn6-a.nc','nowrite');
wu_ambrose = ambrose{'WU_422'}(:);
wv_ambrose = ambrose{'WV_423'}(:);
time_ambrose = ambrose{'time'}(:);
time2_ambrose = ambrose{'time2'}(:);
jdtime_ambrose = singleJD(time_ambrose,time2_ambrose);

data_scale_wind = 3;
plot_scale_wind = 1;

subplot(3,1,1)
title(['Low-passed Wind Stress and Pressure Data at Station ' stationID])
whiskWindNY(jdtime_ambrose,wu_ambrose,wv_ambrose,data_scale_wind,plot_scale_wind);
ylabel('dynes/cm^2')
text(julian(2000,1,1,00)*plot_scale_wind,0,'Estimated 10 m Wind Stress, Ambrose Light');


%Hourly Averaged Data
subplot(3,1,2)
plot(time_h,pres_h-nanmean(pres_h))
axis([startjd endjd pmin pmax])
ylabel('mb')
set(gca,'ytick',[-200:50:200])
gregaxd(time_h,5)

%Low-Pass Filtered Data
subplot(3,1,3)
plot(time_lp,pres_lp-nanmean(pres_lp))
axis([startjd endjd pmin pmax])
ylabel('mb')
set(gca,'ytick',[-200:50:200])
gregaxd(time_lp,5)

if stationID == 'F';
    title('DATA IS FROM SBE-39, NOT MIDAS');
end

orient landscape