% SORTOMP2.M: sorts water type results and residuals,
% latitude and longitude into increasing pressure and
% writes the data as sections into files.
%
% NOTES: 1. This routine reads output from a file. It does
% not require output from omp2.m in the workspace
% but takes output from omp2.m from the file which was
% generated by that program.
% 3. This routine assumes that you have a directory/folder
% "output" in the current work directory and places all
% files into that directory/folder.
%
%---------------------------------------------
% This program is part of the OMP package from:
% Institut fuer Meereskunde
% J. Karstensen
% Troplowitzstr. 7
% 22529 Hamburg
% Germany
%
% It was written by Matthias Tomczak.
%
% BUGS: matthias.tomczak@flinders.edu.au
%--------------------------------------------
%
% This routine produces ASCII output in the following format:
%
% <number of depth levels> <number of stations>
% <maxvar> <minvar>
% <depth level 1> <depth level 2> <depth level 3> ... etc.
% <station location 1> <station location 2> <station location 3> ... etc.
% <data station 1 depth 1> <data station 2 depth 1> <data station 3 depth 1> ... etc.
% <data station 1 depth 2> <data station 2 depth 2> <data station 3 depth 2> ... etc.
% <data station 1 depth 3> <data station 2 depth 3> <data station 3 depth 3> ... etc.
% ... etc.
%
% This format is read automatically by Transform (Spyglass/Fortner) and Noesys (Fortner).
%
% Notes: 1. "Station location" means "Distance from first station" in km.
% 2. maxvar and minvar are the maximum and minimum value for the data. It is set
% to zero, which forces Transform to calculate these values. This is done because
% in most cases max(max(data)) and min(min(data)) turns out to be NaN.
% If you are sure that there are no NaNs in your data you can include the
% calculation of maxvar and minvar from the data in the code.
%
%name variables
varname = [''];
varname(1,1:11) = [' latitude'];
varname(2,1:11) = [' longitude'];
varname(3,1:11) = [' pressure'];
varname(4,1:11) = [' salinity'];
varname(5,1:11) = [' pot. temp.'];
varname(6,1:11) = [' oxygen'];
varname(7,1:11) = [' phosphate'];
varname(8,1:11) = [' nitrate'];
varname(9,1:11) = [' silicate'];
varname(10,1:11) = [' pot. vort.'];
varname(11,1:11) = ['mass consn.'];
filname = [''];
filname(1,1:8) = [' '];
filname(2,1:8) = [' '];
filname(3,1:8) = [' '];
filname(4,1:8) = ['filsalin'];
filname(5,1:8) = ['filptemp'];
filname(6,1:8) = ['filoxyge'];
filname(7,1:8) = ['filphosp'];
filname(8,1:8) = ['filnitra'];
filname(9,1:8) = ['filsilic'];
filname(10,1:8) = ['filpvort'];
filname(11,1:8) = ['filtempe'];
errname(1,1:8) = [' '];
errname(2,1:8) = [' '];
errname(3,1:8) = [' '];
errname(4,1:8) = ['errsalin'];
errname(5,1:8) = ['errptemp'];
errname(6,1:8) = ['erroxyge'];
errname(7,1:8) = ['errphosp'];
errname(8,1:8) = ['errnitra'];
errname(9,1:8) = ['errsilic'];
errname(10,1:8) = ['errpvort'];
errname(11,1:8) = ['errmcons'];
%switch to output directory/folder for input and output
drswitch('This program reads data from');
%define input file
incontrol = input('Which file for data input? [result] ','s');
if length(incontrol) > 0
dataset = incontrol;
else
dataset = 'result';
end
eval(['load ' dataset]);
for i = 1:nr_of_wm
wmname(i,1:7) = ['filewm' num2str(i)];
end
clear ptsort
maxvar = 0;
minvar = 0;
%warning to check existing files
m = 1;
disp(' ')
disp('WARNING! This program will create the following files:')
for i = 4:10
if esx(i) == 1 disp([filname(i,:) ' (for ' varname(i,1:11) ')']); m = m+1; end
end
disp(' ')
for i = 4:10
if esx(i) == 1 disp([errname(i,:) ' (for ' varname(i,1:11) ' residual)']); end
end
disp([errname(11,:) ' (for ' varname(11,1:11) ' residual)']);
disp(' ')
for i = 1:nr_of_wm
disp([wmname(i,1:7) ' (for ' tit_index(5*(i-1)+1:5*(i-1)+5) ')']);
end
want = 'n';
drswitch('It will place the files into');
incontrol = input('Do you want to continue (y/n)? [n] ','s');
if incontrol== 'y' want = 'y'; end
if want == 'n' return; end
%main block, only run if want == 'y'
%find start and end of stations
statind=[find(diff(press)<0) length(press)];
%calculate distance between stations
[dist,phaseangle] = sw_dist(lat,long,'km');
cumdist=[0 cumsum(dist)];
cumdist = cumdist(statind);
statind = [0 statind];
lcd = length(cumdist);
%define interpolation levels
incontrol = input('Which interpolation file do you want to use? [testlevels] ','s');
disp(' ')
if length(incontrol) > 0
dataset = incontrol;
else
dataset = 'testlevels';
end
eval(['load ' dataset]);
llv = length(levels);
%check for duplicates and separate them by 0.5 m.
for i=2:length(press)
if press(i-1) == press(i) press(i-1) = press(i-1) - 0.5; end
end
%initialize data and result matrices
ptsort(length(statind)-1,length(levels)) = NaN;
for i=1:length(statind)-1
for j=1:length(levels)
ptsort(i,j) = NaN;
end
end
%sort data into station format
for nn = 4:10
if esx(nn) == 1
disp(['working on ' varname(nn,1:11) ' extraction'])
for i=2:length(statind)
j = 1;
while press(statind(i-1)+1)>levels(j) & j<length(levels)
j = j+1;
end
k = length(levels);
while press(statind(i))<levels(k) & k>1
k = k-1;
end
sortwrt %stores data in ptsort; see sortwrt.m
end
ptsort = ptsort';
fout = filname(nn,1:8);
sout = sprintf('save %s llv lcd maxvar minvar levels cumdist ptsort -ascii',fout);
eval(sout);
disp(['file ' filname(nn,1:8) ' created and saved'])
disp(' ')
ptsort = ptsort';
end %esx
end %nn
for nn = 1:7
if esx(nn+3) == 1
disp(['working on ' varname(nn+3,1:11) ' residual extraction'])
for i=2:length(statind)
j = 1;
while press(statind(i-1)+1)>levels(j) & j<length(levels)
j = j+1;
end
k = length(levels);
while press(statind(i))<levels(k) & k>1
k = k-1;
end
if j < k
ptsort(i-1,j:k) = interp1(press(statind(i-1)+1:statind(i)), ...
err(nn,statind(i-1)+1:statind(i)),levels(j:k))';
end
end
ptsort = ptsort';
fout = errname(nn+3,1:8);
sout = sprintf('save %s llv lcd maxvar minvar levels cumdist ptsort -ascii',fout);
eval(sout);
disp(['file ' errname(nn+3,1:8) ' created and saved'])
disp(' ')
ptsort = ptsort';
end %esx
end %nn
disp(['working on ' varname(11,1:11) ' residual extraction'])
for i=2:length(statind)
j = 1;
while press(statind(i-1)+1)>levels(j) & j<length(levels)
j = j+1;
end
k = length(levels);
while press(statind(i))<levels(k) & k>1
k = k-1;
end
if j < k
ptsort(i-1,j:k) = interp1(press(statind(i-1)+1:statind(i)), ...
err(m,statind(i-1)+1:statind(i)),levels(j:k))';
end
end
ptsort = 100 * ptsort';
fout = errname(11,1:8);
sout = sprintf('save %s llv lcd maxvar minvar levels cumdist ptsort -ascii',fout);
eval(sout);
disp(['file ' errname(11,1:8) ' created and saved'])
disp(' ')
ptsort = ptsort';
for nn = 1:nr_of_wm
disp(['working on ' tit_index(5*(nn-1)+1:5*(nn-1)+5) ' extraction'])
for i=2:length(statind)
j = 1;
while press(statind(i-1)+1)>levels(j) & j<length(levels)
j = j+1;
end
k = length(levels);
while press(statind(i))<levels(k) & k>1
k = k-1;
end
if j < k
ptsort(i-1,j:k) = interp1(press(statind(i-1)+1:statind(i)), ...
A(nn,statind(i-1)+1:statind(i)),levels(j:k))';
end
end
ptsort = 100 * ptsort';
fout = wmname(nn,1:7);
sout = sprintf('save %s llv lcd maxvar minvar levels cumdist ptsort -ascii',fout);
eval(sout);
disp(['file ' wmname(nn,1:7) ' created and saved'])
disp(' ')
ptsort = ptsort';
end %esx
end %nn
clear ptsort levels cumdist dist phaseangle llv lcd maxvar minvar nn fout sout want
clear ans lat long oxy ph ni si press ptemp pvort temp sal tit_index statind
