Fixdepth ac9 pnb.m

From Pnb
Revision as of 15:16, 5 August 2010 by 128.111.101.185 (talk) (Created page with '<pre> % 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…')
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
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

fixdepth_ac9_pnb