% Matlab script to demonstrate overturn detection and length scale
% Calculation from CTD data.

% 2005 Mar 24
%% 0. Test program setup input
clear all; close all; paper;

% From requirements for Galbraith's reorangeV1.1 package
%
%   Finding overturns in CTD profiles requires good data to start with.
%   There should be no depth reversals in the downcast, otherwise it would
%   be meaningless to interpret density inversions associated with them.
%   The profiles are reorder using density first, then the record number
%   (or index) as a second key instead of pressure.  Therefore bit-noise
%   reversals in the pressure column should bot affect the results.
% 
%   Density spikes will cause large reordering ranges which will mask any
%   other features in the associated depth span.  To find such spikes, you 
%   could plot the Thorpe fluctuation (measured minus reordered density)
%   profile v. depth and look for single large values among noise-related
%   values. The `thorpe' command will calculate the Thorpe fluctuation
%   profile for you.



% Get the data, can have them specify a filename with no header
% or send a matrix of data with the proper columns. Can select with
% isstr if the input var is a str, then load the data from a file as below

% Load his example CTD Profiles, they have no header, the header should be:
% Pressure Temperature Salinity Density
% NOTE: he is passing z as (+) values, so his z-points downward.
%				He is also using sigma_T or sigma_theta (paper says sigma_theta)
p1 = './reorange-V1.1/example/01102605.dat';
p2 = './reorange-V1.1/example/01126339.dat';
rho0 = 1024;  % define a constant reference density to pass in
P1 = load(p1);
P2 = load(p2);

% Calculate the runlengths in the profiles
run1 = runlength(P1(:,4)-sort(P1(:,4)));
run2 = runlength(P2(:,4)-sort(P2(:,4)));

% Calculate the cutoff runlength value for each profile
c1 = ot_cutoff(run1)
c2 = ot_cutoff(run2)
% The cutoff based on the combined profiles
cb = ot_cutoff([run1;run2])


overturns(p1)