clear all; close all;
load ab1209_db
pressure=1:800';
%load ab0409_db.mat

data=input(['Which data would you like to plot?:\n'...
            '0 - raw\n'...
            '1 - calibrated\n'...
            'Select:']);
        if data == 0,
            cnv_nwpc = cnv_dlp;
        else
            cnv_nwpc = cnv_cal;
        end;

sts=5;
pot_temp_prov=zeros(800,sts)*NaN;
sal_prov=zeros(800,sts)*NaN;
sigma0_prov=zeros(800,sts)*NaN;
oxygen_prov=zeros(800,sts)*NaN;
prs=zeros(800,sts)*NaN;

% NW Providence ab1202 stations 34 to 38
% NW Providence ab1202 stations 1 to 5
sensor=input('Do you want to use primary or secondary?(1 or 2)','s')
for i=1:5
   sta=i;
   lon(i)=nanmean(cnv_nwpc(sta).longitude);
   lat(i)=nanmean(cnv_nwpc(sta).latitude);
   max_depth=find(cnv_nwpc(sta).pressure==max(cnv_nwpc(sta).pressure))
   prof_prov(i)=max_depth;
   for j=1:max_depth
       if data == 0,
           if sensor == '1',
            pot_temp_prov(j,i)=cnv_nwpc(sta).thetapri(j);
            sal_prov(j,i)=cnv_nwpc(sta).salpri(j);
            oxygen_prov(j,i)=cnv_nwpc(sta).oxmmpri(j);
            sigma0_prov(j,i)=cnv_nwpc(sta).sigmapri(j);
            prs(j,i)=cnv_nwpc(sta).pressure(j);            
           else
            pot_temp_prov(j,i)=cnv_nwpc(sta).thetasec(j);
            sal_prov(j,i)=cnv_nwpc(sta).salsec(j);
            oxygen_prov(j,i)=cnv_nwpc(sta).oxmmsec(j);
            sigma0_prov(j,i)=cnv_nwpc(sta).sigmasec(j);
            prs(j,i)=cnv_nwpc(sta).pressure(j);    
           end;
       else
        pot_temp_prov(j,i)=cnv_nwpc(sta).pottemp(j);
        sal_prov(j,i)=cnv_nwpc(sta).salt(j);
        oxygen_prov(j,i)=cnv_nwpc(sta).oxyumol(j);
        sigma0_prov(j,i)=cnv_nwpc(sta).sigma0(j); 
       end;   
   end
end


%%%%%%%%%%%%%%%%%% Ploting the potential temperature %%%%%%%%%%%%%%%%%%%%%%
f1
%xi=[lat];
%zi=[1:1:800]';
%[XI,ZI]=meshgrid(xi,zi);
%TI=interp2(lat,pressure,pot_temp_prov,XI,ZI,'linear');
%[cs,hc]=contourf(xi,zi,nghbrfilt(TI),[9:1.5:27]);
%[cs,hc]=contourf(XI,ZI,pot_temp_prov,[0:1.5:30],'-k');%nghbrfilt

contourf(lat,pressure,nghbrfilt(pot_temp_prov),[0:1.5:30],'-k')
%[new_ptem]=make_edge_contour(pot_temp_prov,pressure,prof_prov,'prov');
%contourf(lat,pressure,new_ptem,[0:1.5:30],'-k')

colormap(cejet)
colorbar
hold on
revyaxis

ax1=gca;
set(ax1,'FontSize',14,'ydir','reverse','XAxisLocation','top','Xtick', ...
    [26.1 26.2 26.3 26.4],'XtickLabel',[]);
hold on
set(0,'DefaultTextFontSize',14);
text(26.39,-25,'26.4^oN');
text(26.29,-20,'26.3^oN');
text(26.19,-20,'26.2^oN');
text(26.09,-20,'26.1^oN');
[C,h]=contour(lat,pressure,nghbrfilt(pot_temp_prov),[10 15 20],'k');
clabel(C,h);
set(ax1,'FontSize',14,'ydir','reverse','Xtick',[26.1 26.2 26.3 26.4], ...
    'XtickLabel',[]);
set(h,'LineWidth',2);
%ylabel('Pressure (dbar)');
for t=1:length(prof_prov)
    plot([lat(t) lat(t)],[0 prof_prov(t)+50],'k--')
end
xlabel('NW Providence Section')
ax2=gca;
set(ax2,'XAxisLocation','bottom');set(0,'DefaultTextFontSize',14);
set(0,'DefaultTextFontSize',14);
%xlim(26,26.6)
patch([lat 26.433 26.433 26.071 26.071 ], ...
    [prof_prov 760 800 800 300],'w','LineWidth',3);
ylabel('Pressure (dbar)');
text(26.5,550,'Potential Temperature [^{o}C]','rotation',90)
noaalogo([0.13 0.13 0.15 0.15])

print -depsc figures/nwpc_section_pottemp

%%%%%%%%%%%%%%%%%%%%%%%%%% Ploting the salinity %%%%%%%%%%%%%%%%%%%%%%%%%%%

f2
%xi=[lat];
%zi=[1:1:800]';
%[XI,ZI]=meshgrid(xi,zi);
%TI=interp2(lat,pressure,sal_prov,XI,ZI,'linear');
%[cs,hc]=contourf(xi,zi,nghbrfilt(TI),[35:0.1:37]);
%[cs,hc]=contourf(XI,ZI,sal_prov,[0:1.5:30],'-k');


contourf(lat,pressure,nghbrfilt(sal_prov),'-k');
%[new_sal]=make_edge_contour(sal_prov,pressure,prof_prov,'prov');
%contourf(lat,pressure,nghbrfilt(new_sal),'-k')

colormap(cejet)
colorbar
hold on
revyaxis

ax1=gca;
set(ax1,'FontSize',14,'ydir','reverse','XAxisLocation','top','Xtick', ...
    [26.1 26.2 26.3 26.4],'XtickLabel',[]);
hold on
set(0,'DefaultTextFontSize',14);
text(26.39,-20,'26.4^oN');
text(26.29,-20,'26.3^oN');
text(26.19,-20,'26.2^oN');
text(26.09,-20,'26.1^oN');
[C,h]=contour(lat,pressure,nghbrfilt(sal_prov),[35 35.5 36 36.5],'k');
clabel(C,h);
set(ax1,'FontSize',14,'ydir','reverse','Xtick',[26.1 26.2 26.3 26.4], ...
    'XtickLabel',[]);
set(h,'LineWidth',2);
ylabel('Pressure (dbar)');
for t=1:length(prof_prov)
    plot([lat(t) lat(t)],[0 prof_prov(t)+50],'k--')
end
xlabel('NW Providence Section')
ax2=gca;
set(ax2,'XAxisLocation','bottom');set(0,'DefaultTextFontSize',14);
set(0,'DefaultTextFontSize',14);
patch([lat 26.433 26.433 26.071 26.071 ], ...
    [prof_prov 760 800 800 300],'w','LineWidth',3);
ylabel('Pressure (dbar)');
text(26.51,450,'Salinity (psu)','rotation',90)
noaalogo([0.13 0.13 0.15 0.15])

print -depsc figures/nwpc_section_salinity
%%%%%%%%%%%%%%%%%%%%%%% plotting potential density %%%%%%%%%%%%%%%%%%%%%%%%

f3
%xi=[lat];
%zi=[1:1:800]';
%[XI,ZI]=meshgrid(xi,zi);
%TI=interp2(lat,pressure,sigma0_prov,XI,ZI,'linear');
%[cs,hc]=contourf(xi,zi,nghbrfilt(TI),[23.5:0.5:27.5]);

contourf(lat,pressure,nghbrfilt(sigma0_prov),[22:0.5:27],'-k');
%[new_sigma0]=make_edge_contour(sigma0_prov,pressure,prof_prov,'prov');
%contourf(lat,pressure,nghbrfilt(new_sigma0),[23:0.5:27],'-k');

colormap(cejet)
colorbar
hold on
revyaxis

ax1=gca;
set(ax1,'FontSize',14,'ydir','reverse','XAxisLocation','top','Xtick', ...
    [26.1 26.2 26.3 26.4],'XtickLabel',[]);
hold on
set(0,'DefaultTextFontSize',14);
text(26.39,-20,'26.4^oN');
text(26.29,-20,'26.3^oN');
text(26.19,-20,'26.2^oN');
text(26.09,-20,'26.1^oN');
[C,h]=contour(lat,pressure,nghbrfilt(sigma0_prov),[26 27],'k');
clabel(C,h);
set(ax1,'FontSize',14,'ydir','reverse','Xtick',[26.1 26.2 26.3 26.4], ...
    'XtickLabel',[]);
set(h,'LineWidth',2);
ylabel('Pressure (dbar)');
for t=1:length(prof_prov)
    plot([lat(t) lat(t)],[0 prof_prov(t)+50],'k--')
end
xlabel('NW Providence Section')
ax2=gca;
set(ax2,'XAxisLocation','bottom');set(0,'DefaultTextFontSize',14);
set(0,'DefaultTextFontSize',14);
patch([lat 26.433 26.433 26.071 26.071 ], ...
    [prof_prov 760 800 800 300],'w','LineWidth',3);
ylabel('Pressure (dbar)');
text(26.51,450,'Sigma_{\theta}','rotation',90)
noaalogo([0.13 0.13 0.15 0.15])

print -depsc figures/nwpc_section_sigma0
%%%%%%%%%%%%%%%%%%%%%%%%%%%% plotting oxygen %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

f4
%xi=[lat];
%zi=[1:1:800]';
%[XI,ZI]=meshgrid(xi,zi);
%TI=interp2(lat,pressure,oxygen_prov,XI,ZI,'linear');
%[cs,hc]=contourf(xi,zi,nghbrfilt(TI));

contourf(lat,pressure,nghbrfilt(oxygen_prov),10,'-k');
%[new_oxy]=make_edge_contour(oxygen_prov,pressure,prof_prov,'prov');
%contourf(lat,pressure,nghbrfilt(new_oxy),'-k');

colormap(cejet)
colorbar
hold on
revyaxis

ax1=gca;
set(ax1,'FontSize',14,'ydir','reverse','XAxisLocation','top','Xtick', ...
    [26.1 26.2 26.3 26.4],'XtickLabel',[]);
hold on
set(0,'DefaultTextFontSize',14);
text(26.39,-20,'26.4^oN');
text(26.29,-20,'26.3^oN');
text(26.19,-20,'26.2^oN');
text(26.09,-20,'26.1^oN');
%[C,h]=contour(xi,zi,nghbrfilt(TI),[26 27],'k');
%clabel(C,h);
set(ax1,'FontSize',14,'ydir','reverse','Xtick',[26.1 26.2 26.3 26.4], ...
    'XtickLabel',[]);
set(h,'LineWidth',2);
ylabel('Pressure (dbar)');
for t=1:length(prof_prov)
    plot([lat(t) lat(t)],[0 prof_prov(t)+50],'k--')
end
xlabel('NW Providence Section')
ax2=gca;
set(ax2,'XAxisLocation','bottom');set(0,'DefaultTextFontSize',14);
set(0,'DefaultTextFontSize',14);
patch([lat 26.433 26.433 26.071 26.071 ], ...
    [prof_prov 760 800 800 300],'w','LineWidth',3);
ylabel('Pressure (dbar)');
text(26.5,600,'Dissolved Oxygen [{\mu}mol/kg]','rotation',90)
noaalogo([0.13 0.13 0.15 0.15])

print -depsc figures/nwpc_section_oxygen