s_align.m

地震、测井方面matlab代码,解释的比较详细

      [lagi,cci]=max_corr(aa,refi,-2*maxshift,2*maxshift);
    end
%   [cci,lagi]=cornor(aa(:,ii),refi,'max');
%    [lagi,cci]=max_corr(aa(:,ii),refi,-2*maxshift,2*maxshift);
    lag(ii)=lagi-ia+iae;
    cc(ii)=cci; 
    temp=zeros(nref,1);
    ie=min([lagi+nref,nsamp]);
    ia2=max([1,lagi+1]);
    ia1=max([1,nref-ie+1]);
    temp(ia1:ia1+ie-ia2)=aa(ia2:ie,ii);
%    [lagix,ccix]=max_corr(temp,refi,-4,4)
    refi=refi*(1-update)+update*temp;
    
%    refi(1:ntemp-lagi)=refi(1:ntemp-lagi)+update*aa(lagi+1:ntemp,ii); % Update reference signal
 
%   refi=refi*(1-update)+update*aa(lagi+1:lagi+nref,ii);    
  end

  refi=ref;
  for ii=itr+1:ntr
    if method == 1 
      [cci,lagi]=cornor(aa,refi,'max');
    else
      [lagi,cci]=max_corr(aa,refi,-2*maxshift,2*maxshift);
    end
%    [cci,lagi]=cornor(aa(:,ii),refi,'max');
%    [lagi,cci]=max_corr(aa(:,ii),refi,-2*maxshift,2*maxshift);
    lag(ii)=lagi-ia+iae;
    cc(ii)=cci; 
    temp=zeros(nref,1);
    ie=min([lagi+nref,nsamp]);
    ia2=max([1,lagi+1]);
    ia1=max([1,nref-ie+1]);
    temp(ia1:ia1+ie-ia2)=aa(ia2:ie,ii);
%    [lagix,ccix]=max_corr(temp,refi,-4,4)
    refi=refi*(1-update)+update*temp;
  end 

               case 'yes'         % Adapt signal and update window
  iai=ia;                         % Initialize start of search window
  iei=ie;
  lag=zeros(ntr,1);
  cc=ones(ntr,1);
  refi=ref;
  for ii=itr-1:-1:1
     minshift=-min([maxshift,iai-1]);
     iae=iai+minshift;
%    iae=max([iai-maxshift,1]);
     aa=a(iae:min([iei+maxshift,nsamp]),ii);
     if method == 1 
        [cci,lagi]=cornor(aa,refi,'max');
     else
%      [lagi,cci]=max_corr(aa,refi,-maxshift,maxshift);
      [ lagi,cci]=max_corr(aa,refi,0,-minshift+maxshift);
     end 
%    [cci,lagi]=cornor(aa,refi,'max');
%    [lagi,cci]=max_corr(refi,aa,-maxshift,maxshift);
%    lag(ii)=lagi-ia+iae;
     lag(ii)=lag(ii+1)-lagi-minshift;
     cc(ii)=cci;
%    if ii <=20 & ii >= 21
%      keyboard
%    end
%    iai=ia+lagi;     % Update start of search window
%    iei=ie+lagi;
     iai=ia-lag(ii);     % Update start of search window
     iei=ie-lag(ii);

% 	Update reference signal                    
%		try
%    mplot(refi,aa(lagi+1:lagi+nref));
     if lagi >= 0
        refi=refi*(1-update)+update*aa(lagi+1:lagi+nref);
     else
       nref=nref+lagi;
       refi=refi(1-lagi:end)+update*aa(1:nref);
     end
%                catch
%    keyboard    
%                end
  end

  iai=ia;             % Initialize start of search window
  iei=ie;
  refi=ref;
  nref=length(ref);

  for ii=itr+1:ntr
     minshift=-min([maxshift,iai-1]);
     iae=iai+minshift;
%    iae=max([iai-maxshift,1]);
     aa=a(iae:min([iei+maxshift,nsamp]),ii);
     if method == 1 
      [cci,lagi]=cornor(aa,refi,'max');
     else
%      [lagi,cci]=max_corr(aa,refi,-maxshift,maxshift);
        [lagi,cci]=max_corr(aa,refi,0,-minshift+maxshift);
     end 
%    [cci,lagi]=cornor(aa,refi,'max');
%    [lagi,cci]=max_corr(aa,refi,-2*maxshift,2*maxshift);
%    lag(ii)=lagi-ia+iae;
     lag(ii)=lag(ii-1)-lagi-minshift;
     cc(ii)=cci;
%    keyboard
%    iai=ia+lagi;     % Update start of search window
%    iei=ie+lagi;
     iai=ia-lag(ii);     % Update start of search window
     iei=ie-lag(ii);

% 	Update reference signal                    
%		try
%    mplot(refi,aa(lagi+1:lagi+nref));
     if lagi >= 0
        refi=refi*(1-update)+update*aa(lagi+1:lagi+nref);
     else
        nref=nref+lagi;
        refi=refi(1-lagi:end)+update*aa(1:nref);
     end
%                catch
%    keyboard    
%                end
  end

               case 'otherwise'
  error([' Unknown option for adapt: ',adapt])
               end

else
   error([' Update parameter not within [0,1]: ',num2str(update)]) 
end      




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [lag,cc]=trace_align_adapt(a,itr,ia,ie,maxshift,update,method)
% Function aligns a "signal" in columns of matrix "a" with the signal in a 
% reference column. 
% It outputs the shifts required for alignment with this reference signal 
% and the coefficient of correlation.
%                [lag,cc]=align(a,itr,ia,ie,maxshift,update,method)
% INPUT
% a         matrix whose columns are to be aligned
% itr       index of column to use as reference for alignment
% ia        index of first row to use for the signal in the reference column
% ie        index of last row to use for the signal in the reference column
% maxshift  maximum shift (up or down) allowed
% update    numerical factor specifying if the reference signal is to be updated to account
%           for systematic column-to-column changes in the signal (e.g. NMO stretch).
%           "update" must satisfy 0 <= update <= 1.
%           new_reference_signal=(1-update)*old_reference_signal + update*signal_on_last_trace
% method    option regarding measurement of similarity
% OUTPUT
% lag       shifts found (of course, lag(itr)=0). Positive values of the lag 
%           mean a shift down, negative values a shift up.
% cc        crosscorrelation between the reference signal and the signal
%            on each trace. (Of course, cc(itr)=1).


[nsamp,ntr]=size(a);
ref=a(ia:1:ie,itr);
nref=length(ref);
%index=zeros(ntr,1);
cc=ones(ntr,1);

iai=ia;                         % Initialize start of search window
iei=ie;
lag=zeros(ntr,1);
cc=ones(ntr,1);
refi=ref;
for ii=itr-1:-1:1
  minshift=-min([maxshift,iai-1]);
  iae=iai+minshift;
  aa=a(iae:min([iei+maxshift,nsamp]),ii);
  if method == 1 
    [cci,lagi]=cornor(aa,refi,'max');
  else
    [lagi,cci]=max_corr(aa,refi,0,-minshift+maxshift);
  end 
  lag(ii)=lag(ii+1)-lagi-minshift;
  cc(ii)=cci;

  iai=ia-lag(ii);   % Update start of search window
  iei=ie-lag(ii);

% 	Update reference signal          
  if lagi >= 0
  tie=iae+lagi+nref-1;
    if tie > nsamp
      refi=refi*(1-update)+[update*a(iae+lagi:nsamp,ii);zeros(tie-nsamp,1)];
    else
      refi=refi*(1-update)+update*a(iae+lagi:tie,ii);
    end
  else
    refi=refi+[zeros(1-lagi,1);update*aa(1:nref+lagi-1)];
  end

end

iai=ia;             % Initialize start of search window
iei=ie;
refi=ref;
for ii=itr+1:ntr
  minshift=-min([maxshift,iai-1]);
  iae=iai+minshift;
  aa=a(iae:min([iei+maxshift,nsamp]),ii);
  if method == 1 
    [cci,lagi]=cornor(aa,refi,'max');
  else
    [lagi,cci]=max_corr(aa,refi,0,-minshift+maxshift);
  end 
  lag(ii)=lag(ii-1)-lagi-minshift;
  cc(ii)=cci;
  iai=ia-lag(ii);     % Update start of search window
  iei=ie-lag(ii);

% 	Update reference signal         
  if lagi >= 0
    tie=iae+lagi+nref-1;
    if tie > nsamp
      refi=refi*(1-update)+[update*a(iae+lagi:nsamp,ii);zeros(tie-nsamp,1)];
    else
      refi=refi*(1-update)+update*a(iae+lagi:tie,ii);
    end 
  else
    nref=nref+lagi;
%    refi=refi(1-lagi:end)+update*aa(1:nref);
    refi=refi+[zeros(1-lagi,1);update*aa(1:nref+lagi-1)];
  end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [lag,cc]=trace_align_noadapt(a,itr,ia,ie,maxshift,update,method)
% Function aligns a "signal" in columns of matrix "a" with the signal in a 
% reference column. 
% It outputs the shifts required for alignment with this reference signal 
% and the coefficient of correlation.
% INPUT
% a         matrix whose columns are to be aligned
% itr       index of column to use as reference for alignment
% ia        index of first row to use for the signal in the reference column
% ie        index of last row to use for the signal in the reference column
% maxshift  maximum shift (up or down) allowed
% update    numerical factor specifying if the reference signal is to be updated to 
%           accountfor systematic column-to-column changes in the signal .
%           (e.g. NMO stretch) "update" must satisfy 0 <= update <= 1.
%           new_reference_signal=(1-update)*old_reference_signal 
%                               +update*signal_on_last_trace
% method    option regarding measurement of similarity
% OUTPUT
% lag       shifts found (of course, lag(itr)=0). Positive values of the lag mean a shift
%           down, negative values a shift up.
% cc        crosscorrelation between the reference signal and the signal on each trace.
%           (Of course, cc(itr)=1).
%                [lag,cc]=align(a,itr,ia,ie,maxshift,adapt,update)

[nsamp,ntr]=size(a);
ref=a(ia:1:ie,itr);
nref=length(ref);
cc=ones(ntr,1);

minshift=-min([maxshift,ia-1]);
%iae=max([ia-maxshift,1]);
iae=ia+minshift;
aa=a(iae:min([ie+maxshift,nsamp]),:);
lag=zeros(ntr,1);
cc=ones(ntr,1);
refi=ref;
for ii=itr-1:-1:1
  if method == 1 
    [cci,lagi]=cornor(aa(:,ii),refi,'max');
  else
    [lagi,cci]=max_corr(aa(:,ii),refi,0,maxshift-minshift);
  end

  lag(ii)=-lagi-minshift;
  cc(ii)=cci;
  ia1=max(ia-lag(ii),1);
  ia1=min(nsamp-nref+1,ia1);
  refi=(1-update)*refi+update*a(ia1:ia1+nref-1,ii);
end

refi=ref;
for ii=itr+1:ntr
  if method == 1 
    [cci,lagi]=cornor(aa(:,ii),refi,'max');
  else
    [lagi,cci]=max_corr(aa(:,ii),refi,-2*maxshift,2*maxshift);
  end
  lag(ii)=-lagi-minshift;
  cc(ii)=cci;

  ia1=max(ia-lag(ii),1);
  ia1=min(nsamp-nref+1,ia1);
  refi=(1-update)*refi+update*a(ia1:ia1+nref-1,ii);
end