Fixdepth ac9 pnb.m
Jump to navigation
Jump to search
% FIXDEPTH_AC9.M
% 29 SEP 97, 04 JUN 2000
% John Ubante, SRW
% Input: stripped mer, ac9 files [time depth]
% Output: returns scale and offset of regression of ac9 on mer
% steps_ac9$cruisename:
% [startindex, X, Y, endindex, X, Y, avedepth, avetime]
% depthoffset_$cruisename: [scale offset]
% Assumptions: you are in a $cruisename/acmpc directory
% Purpose: to find values that will align stupid ac9
% values to less stupid mer values
% Calls: {shrink,xnr,elimexp,jexcise,combine}
% Maintain: {filename, columnarraylist, disp_text}
% Example: fixdepth_ac9
% Notes:
% 06/2000 Doesn't account for instrument position on cage.
% Changed to function, pass arrays retrieve from mat file.
%
function [scale,offset,matrix]=fixdepth_ac9(mer,ac9)
% declare local variables
mat_in_flag=-1.0e+40;
flag=-9.9e+35;
matrixnull=[];
matrix=matrixnull;
yesno=2;
ctr=0;
matrixoffset=0;
matrixbrand='mer';
numsteps=[];
xmax=0; % get cruisename
pwdd = pwd;
indslash = find(pwdd(1,:)=='/');
cruisename = pwdd(1,indslash(1,size(indslash,2)-1)+1:indslash(1,size(indslash,2))-1);
clear indslash pwdd
% load mer and ac9 data
%load mat_in.depth;
% $$$ % separate mer and ac9 data
% $$$ indflag=find(mat_in(:,1)==mat_in_flag);
% $$$ mer=mat_in(1:indflag-1,:);
% $$$ ac9=mat_in(indflag+1:size(mat_in,1),:);
% $$$ clear flag mat_in indflag
% check for sign on depth and
% make positive
if mean(mer(:,2))<0
mer(:,2)=-mer(:,2);
end % end if
if mean(ac9(:,2))<0
ac9(:,2)=-ac9(:,2);
end % end if
% check for validity of data
if size(mer,2)~=2
disp('Your mer data is not in 2 column format');
end
if size(ac9,2)~=2
disp('Your ac9 data is not in 2 column format')
end
% begin step interval loop
dataset=mer;
while(yesno)
ctr=ctr+1;
% if first step of data ask for number
% of step intervals
if ctr==1
% plot data real big
fig=figure;
set (gcf,'units','inches', 'pos',[2.5 2 8 5])
plot(dataset(:,1),dataset(:,2))
set(gca,'ydir','rev')
titlestr=['Depth vs. Index in ',matrixbrand,' data of ',cruisename];
title(titlestr);
grid on;
disp(' ')
disp('A whole step is defined as a step preceeded and followed by sloped data.')
disp('Count the bottom step only if both plots display a leading slope.')
numsteps(yesno)=input('Count the whole steps enter the value: ');
disp('Select the boundaries of your first interval.')
else
disp('Select the boundaries of your next interval.')
end % end if
% get zoom points
[x,y]=ginput(2);
% zoom in
if y(1)<y(2)
if x(1)<x(2)
axis([x(1),x(2),y(1),y(2)])
else
axis([x(2),x(1),y(1),y(2)])
end
else
if x(1)<x(2)
axis([x(1),x(2),y(2),y(1)])
else
axis([x(2),x(1),y(2),y(1)])
end % end if
end % end if
% if y(1) > y(2)
% axis([x(1),x(2),y(2),y(1)]);
% else
% axis([x,y]);
% end
hold on;
grid off;
plot(dataset(:,1),dataset(:,2),'ro');
titlestr=['Zoom of step interval #',num2str(ctr),' of ',cruisename];
title(titlestr);
plot([x(1), x(2)], [0, 0], 'g-.');
hold off;
% get start and end points of step interval
disp('click on start and finish of step interval')
[x,y]=ginput(2);
% convert x(:) to indices
for i=1:2
find(x(i)<=dataset(:,1));
xind(i)=ans(1);
end
% find interval's averages
aved=mean(dataset(xind(1):xind(2),2));
avet=mean(dataset(xind(1):xind(2),1));
% add step interval data to matrix
%matrix=[matrix;x(1),dataset(xind(1),:),x(2),dataset(xind(2),:),aved,avet];
matrix=[matrix;xind(1),dataset(xind(1),:),xind(2),dataset(xind(2),:),aved,avet];
% plot dataset real big with
% already checked step intervals
plot(dataset(:,1),dataset(:,2),'g')
set(gca,'ydir','rev')
titlestr=['Depth vs. Index in ',matrixbrand,' data of ',cruisename];
title(titlestr)
hold on
% loop to plot done step intervals
for i=matrixoffset+1:size(matrix,1),
plot(dataset(matrix(i,1):matrix(i,4),1),dataset(matrix(i,1):matrix(i,4),2),'m')
end % end for loop
clear i
plot(matrix(matrixoffset+1:size(matrix,1),2),matrix(matrixoffset+1:size(matrix,1),3),'ro')
plot(matrix(matrixoffset+1:size(matrix,1),5),matrix(matrixoffset+1:size(matrix,1),6),'ro')
grid on;
% after numsteps intervals, ask if points ok
if ctr==numsteps(yesno)
titlestr=['Selected steps in ',matrixbrand,' of ',cruisename];
title(titlestr);
done=input('Are these points ok? ','s');
% if points ok, exit loop
if (done=='y'|done=='Y')
% eval(['save steps_ac9' cruisename '.dat matrix -ascii']);
disp('You RULE!!');
yesno=yesno-1;
dataset=ac9;
matrixnull=matrix;
ctr=0;
matrixoffset=size(matrix,1);
matrixbrand='ac9';
grid on
set (gcf, 'pos',[0.9 5 5.5 5])
else
matrix=matrixnull;
ctr=0;
end % end if
end % end if
end % end while loop
clear yesno titlestr aved avet ctr done fig matrixnull matrixoffset x y dataset
clear matrixbrand
temp=numsteps(1);
numsteps(1)=numsteps(2);
numsteps(2)=temp;
clear temp
grid on
set (gcf, 'pos',[6.5 5 5.5 5])
hold off
% DIAGNOSTIC PLOT: superimpose all step intervals from both datasets on to one plot.
% The x axis will represent index_end - index_start. The y axis will be depth.
% Because of the different sample rates between the two instruments, rescale the
% x axis for the ac9 to match its corresponding mer step interval. DAMN!!!
% FINISH LATER
figb=figure;
set(gca,'ydir','rev')
% set up plot
% begin looping through each line of matrix
dataset=mer;
for i=1:numsteps(1)
hold on
x=dataset(matrix(i,1):matrix(i,4),1)-dataset(matrix(i,1),1);
if max(x)> xmax
xmax=max(x);
end % end if
y=dataset(matrix(i,1):matrix(i,4),2)-dataset(matrix(i,1),2);
plot(x,y,'ro')
plot(x,y,'g')
hold off
end % end for
hold on
plot([0 xmax],[0 0],'w--')
hold off
clear xmax
done=input('Are these steps ok? ','s');
if (done=='y'|done=='Y')
close(figb)
clear figb
else
disp('re-run by typing "fixdepth_ac9" at the prompt');
return%break;
end % end if
% DO LATER: if numsteps(1)~=numsteps(2) then remove the steps in the cruise
% with more steps until the number of steps are same, ensuring that only steps
% without corresponding steps are eliminated.
% DO NOW: Find offset and scale for ac9 data [depthave timeave]
dt_mer=matrix(1:numsteps(1),7:8);
dt_ac9=matrix(numsteps(1)+1:sum(numsteps),7:8);
[p,S]=polyfit(dt_mer(:,1),dt_ac9(:,1),1);
% eval(['save scaleoffset' cruisename '.dat p -ascii']);
scale=p(1);
offset=p(2);
r=corrcoef(dt_mer(:,1)',dt_ac9(:,1)')
sprintf('YOUR SCALE _IS_: %3.3e',scale)
sprintf('YOUR OFFSET _IS_: %3.3e',offset)
plot(dt_mer(:,1),dt_ac9(:,1),'o');
hold on
plot(dt_mer(:,1),polyval(p,dt_mer(:,1)),'k')
ax=axis;
axis([min(ax(1),ax(3)) max(ax(2),ax(4)) min(ax(1),ax(3)) max(ax(2),ax(4))]);
return
