%
% can't use station 6 - just asfter the shoft because the pump was bad.
% use after station 7
tdiff = cnv_dlp(7).temppri - cnv_dlp(7).tempsec

II=find(cnv_dlp(7).pressure>=3300);
JJ=find(cnv_dlp(7).pressure<3300);


[m_deep, b_deep, r, sm, sb, rms_err, mean_per_err, mean_per_err2] = ...
    do_mbari_fits2 (cnv_dlp(7).pressure(II), tdiff(II), 1);
[m_shal, b_shal, r, sm1, sb, rms_err, mean_per_err, mean_per_err2] = ...
    do_mbari_fits2 (cnv_dlp(7).pressure(JJ), tdiff(JJ), 1);

figure;
    plot (cnv_dlp(7).pressure(II), tdiff(II))
    hold on
    title(['Station ' num2str(cnv_dlp(7).station)])
    plot(cnv_dlp(7).pressure(II), m_deep(5) * cnv_dlp(7).pressure(II) + b_deep(5),  'k', 'linewidth', 2)

    plot (cnv_dlp(7).pressure(JJ), tdiff(JJ),'r')
    plot(cnv_dlp(7).pressure(JJ), m_shal(5) * cnv_dlp(7).pressure(JJ) + b_shal(5),  'k', 'linewidth', 2)
    set(gca,'ylim',[-0.002 +0.002]) 

    
ista_idx = 24;    
II=find(cnv_dlp(ista_idx).pressure>=3300);
JJ=find(cnv_dlp(ista_idx).pressure<3300);
tdiff = cnv_dlp(ista_idx).temppri - cnv_dlp(ista_idx).tempsec;
[m_deep24, b_deep24, r2, sm2, sb2, rms_err2, mean_per_err2, mean_per_err22] = ...
    do_mbari_fits2 (cnv_dlp(ista_idx).pressure(II), tdiff(II), 1);
[m_shal24, b_shal24, r3, sm3, sb3, rms_err3, mean_per_err3, mean_per_err23] = ...
    do_mbari_fits2 (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ), 1);

figure;
    plot (cnv_dlp(ista_idx).pressure(II), tdiff(II))
    hold on
    title(['Station ' num2str(cnv_dlp(ista_idx).station)])
    plot(cnv_dlp(ista_idx).pressure(II), m_deep(5) * cnv_dlp(ista_idx).pressure(II) + b_deep(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(II), m_deep24(5) * cnv_dlp(ista_idx).pressure(II) + b_deep24(5),  'k--', 'linewidth', 2)

    plot (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ),'r')
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal24(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal24(5),  'k--', 'linewidth', 2)
    set(gca,'ylim',[-0.002 +0.002]) 
        

%
% check which slope/offset works best
%
ista_idx = 31;    
II=find(cnv_dlp(ista_idx).pressure>=3300);
JJ=find(cnv_dlp(ista_idx).pressure<3300);
tdiff = cnv_dlp(ista_idx).temppri - cnv_dlp(ista_idx).tempsec;
figure;
    plot (cnv_dlp(ista_idx).pressure(II), tdiff(II))
    hold on
    title(['Station ' num2str(cnv_dlp(ista_idx).station)])
    plot(cnv_dlp(ista_idx).pressure(II), m_deep(5) * cnv_dlp(ista_idx).pressure(II) + b_deep(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(II), m_deep24(5) * cnv_dlp(ista_idx).pressure(II) + b_deep24(5),  'k--', 'linewidth', 2)

    plot (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ),'r')
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal24(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal24(5),  'k--', 'linewidth', 2)
    set(gca,'ylim',[-0.002 +0.002]) 


ista_idx = 8;    
II=find(cnv_dlp(ista_idx).pressure>=3300);
JJ=find(cnv_dlp(ista_idx).pressure<3300);
tdiff = cnv_dlp(ista_idx).temppri - cnv_dlp(ista_idx).tempsec;
figure;
    plot (cnv_dlp(ista_idx).pressure(II), tdiff(II))
    hold on
    title(['Station ' num2str(cnv_dlp(ista_idx).station)])
    plot(cnv_dlp(ista_idx).pressure(II), m_deep(5) * cnv_dlp(ista_idx).pressure(II) + b_deep(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(II), m_deep24(5) * cnv_dlp(ista_idx).pressure(II) + b_deep24(5),  'k--', 'linewidth', 2)

    plot (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ),'r')
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal(5),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal24(5) * cnv_dlp(ista_idx).pressure(JJ) + b_shal24(5),  'k--', 'linewidth', 2)
    set(gca,'ylim',[-0.002 +0.002]) 

    
%
% Now use several stations all at once
%
depth_match = 3300;
shallow_depth_match = 1200;
ista_start = 8;
ista_stop = 27;
temp_pressure = vertcat(cnv_dlp(ista_start:ista_stop).pressure);
temp_pri = vertcat(cnv_dlp(ista_start:ista_stop).temppri);
temp_sec = vertcat(cnv_dlp(ista_start:ista_stop).tempsec);
temp_tdiff = temp_pri-temp_sec;

II=find(temp_pressure>=depth_match);
JJ=find(temp_pressure<depth_match & temp_pressure>=shallow_depth_match);
[m_all_deep, b_all_deep, r, sm, sb, rms_err, mean_per_err, mean_per_err2] = ...
    do_mbari_fits2 (temp_pressure(II), temp_tdiff(II), 1);
[m_all_shal, b_all_shal, r, sm, sb, rms_err, mean_per_err, mean_per_err2] = ...
    do_mbari_fits2 (temp_pressure(JJ), temp_tdiff(JJ), 1);
    
ista_idx = 9; 
ifit_to_use = 1;
for ista_idx = ista_start: ista_stop
    II=find(cnv_dlp(ista_idx).pressure>=depth_match);
    JJ=find(cnv_dlp(ista_idx).pressure<depth_match);
    tdiff = cnv_dlp(ista_idx).temppri - cnv_dlp(ista_idx).tempsec;
    figure;
    plot (cnv_dlp(ista_idx).pressure(II), tdiff(II))
    hold on
    title(['Station ' num2str(cnv_dlp(ista_idx).station)])
    plot(cnv_dlp(ista_idx).pressure(II), m_deep(ifit_to_use) * cnv_dlp(ista_idx).pressure(II) + b_deep(ifit_to_use),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(II), m_deep24(ifit_to_use) * cnv_dlp(ista_idx).pressure(II) + b_deep24(ifit_to_use),  'k--', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(II), m_all_deep(ifit_to_use) * cnv_dlp(ista_idx).pressure(II) + b_all_deep(ifit_to_use),  'g--', 'linewidth', 2)

    plot (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ),'r')
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal(ifit_to_use) * cnv_dlp(ista_idx).pressure(JJ) + b_shal(ifit_to_use),  'k', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(JJ), m_shal24(ifit_to_use) * cnv_dlp(ista_idx).pressure(JJ) + b_shal24(ifit_to_use),  'k--', 'linewidth', 2)
    plot(cnv_dlp(ista_idx).pressure(JJ), m_all_shal(ifit_to_use) * cnv_dlp(ista_idx).pressure(JJ) + b_all_shal(ifit_to_use),  'g--', 'linewidth', 2)
    set(gca,'ylim',[-0.002 +0.002]) 
end

close all

depth_match = 3300;
m_use_deep = m_all_deep(1);
b_use_deep = b_all_deep(1);
%
% make sure T correction is continuous.
% to do that figure out correction at 3000 m
%
tdiff_3000 = m_use_deep * depth_match + b_use_deep;
%
% now use an upper ocean fit that starts at:
% (3000., tdiff_3000)
% Equation for new line is:
% m_use_shall = (tdiff - tdiff_3000) / (pr - 3000.)
% rearranging
% tdiff = m_use_shall*pr + m_use_shall*(-3000)+tdiff_3000
% tdiff = m_use_shall*pr + new_b_shall
m_use_shal = m_all_shal(1);
b_use_shal = m_use_shal*(-depth_match)+tdiff_3000;

ifit_to_use = 1;
for ista_idx = ista_start: ista_stop
    II=find(cnv_dlp(ista_idx).pressure>=depth_match);
    JJ=find(cnv_dlp(ista_idx).pressure<depth_match);
    tdiff = cnv_dlp(ista_idx).temppri - cnv_dlp(ista_idx).tempsec;
    figure;
        plot (cnv_dlp(ista_idx).pressure(II), tdiff(II))
        hold on
        title(['Station ' num2str(cnv_dlp(ista_idx).station)])
        plot(cnv_dlp(ista_idx).pressure(II), m_use_deep * cnv_dlp(ista_idx).pressure(II) + b_use_deep,  'k', 'linewidth', 2)

        plot(depth_match, tdiff_3000, 'ok','markersize',20)
        
        plot (cnv_dlp(ista_idx).pressure(JJ), tdiff(JJ),'r')
        plot(cnv_dlp(ista_idx).pressure(JJ), m_use_shal * cnv_dlp(ista_idx).pressure(JJ) + b_use_shal,  'k', 'linewidth', 2)
        set(gca,'ylim',[-0.002 +0.002]) 
end


