s_append.m

matlab源代码

function seisout=s_append(seis1,seis2,varargin)
% Function combines two seismic data sets. On output the traces of the 
% first input data set are followed by the traces of the second data set.
% The two data sets must satisfy the following two conditions:
%           seis1.step=seis2.step;
%           (seis1.first-seis2.first)/seis1.step is integer
% Unless input parameter "ttype" is {'max','min'} and input parameter "htype" is 'i'
% some seismic samples and/or header values of "seisout" may not be defined by the 
% input data.
% If seis1.null and seis2.null exist and are the same this null value is used 
%     for undefined data
% If only one of the data sets has a null value defined it is used. 
% If neither of the data sets has a null value defined NaN is used.
%
% Written by: E. R.: January 28, 2000
% Last updated: August 20, 2003: fix bug in assignment of null value
%
%           seisout=s_append(seis1,seis2,varargin)
% INPUT 
% seis1   seismic structure
% seis2   seismic structure
% varargin    one or more cell arrays; the first element of each cell array is a 
%             keyword, the other elements are parameters. Presently, keywords are:
%    'header_null' Null value to use for non-existing header values in the output data set
%    'null'    Null value to use for non-existing trace samples in the output data set
%    'times'   This two-cell string defines what should be done if the two data sets have 
%         differing first-sample and/or last-sample times
%         times{1} = 'min': seisout.start=min([seis1.first,seis2.first]);
%         times{1} = 'max': seisout.start=max([seis1.first,seis2.first]);
%         times{2} = 'min': seisout.start=min([seis1.last,seis2.last]);
%         times{2} = 'max': seisout.start=max([seis1.last,seis2.last]);
%         Default: {'times','min','max'}
%    'headers'   This string parameter defines what should be done if the two data sets have
%         different headers
%         {'headers','u'} : Use the union of the headers (default)
%         {'headers','i'} ; Use the intersection of the headers (the headers they have in common)
%
%         See also: s_select

%       Set defaults parameters
param.header_null=NaN;
param.null=NaN;
param.times={'min','max'};
param.headers='u';

%       Decode and assign input arguments
param=assign_input(param,varargin);

ttype=param.times;
htype=param.headers;

% Check compatibility of data sets
if ~strcmpi(seis1.units,seis2.units)
  error(['Input data sets have different depth/time/frequncy units ("',seis1.units,'" vs. "', ...
     seis2.units,'")'])
end

if seis1.step ~= seis2.step
  error([' Input data sets have different sample intervals (',num2str(seis1.step),' vs. ', ...
     num2str(seis2.step),' ',seis1.units,')'])
end

%if ~exist('ttype','var') | isempty(param.ttype)
%   ttype={'min','max'};
%end
%if ~exist('htype','var') | isempty(param.htype)
%  htype='u';
%end

seisout.type=seis1.type;
seisout.tag=seis1.tag;
seisout.name='';
step=seis1.step;
seisout.step=step;
seisout.units=seis1.units;

%       Convert to compact header representation (if headers exist)
%seis1=s_new2compact(seis1);
%seis2=s_new2compact(seis2);
                   
null_value=param.null;
header_null_value=param.header_null;

fdiff=(seis1.first-seis2.first)/seis1.step;
ldiff=(seis1.last-seis2.last)/seis1.step;
[nsamp1,ntr1]=size(seis1.traces);
[nsamp2,ntr2]=size(seis2.traces);


if fdiff == 0 & ldiff == 0       % Simple case where both data sets have the same start and end times
  seisout.traces=[seis1.traces,seis2.traces];
  seisout.first=seis1.first;
  seisout.last=seis1.last;

else
  if strcmpi(ttype(1),'max')
    first=max([seis1.first,seis2.first]);
  elseif strcmpi(ttype(1),'min')
    first=min([seis1.first,seis2.first]);
  else
    error(['Unknown input parameter ttype(1) (',ttype(1),')']) 
  end

  if strcmpi(ttype(2),'max')
     last=max([seis1.last,seis2.last]);
  elseif strcmpi(ttype(2),'min')
     last=min([seis1.last,seis2.last]);
  else
    error(['Unknown input parameter ttype(1) (',ttype(1),')']) 
  end

  seisout.first=first;
  seisout.last=last;
  nsamp=round((last-first)/step)+1;

  seisout.traces=null_value*zeros(nsamp,ntr1+ntr2);

  ia=round((first-seis1.first)/step);
  if ia >= 0
    seisout.traces(1:min([nsamp1-ia,nsamp]),1:ntr1)=seis1.traces(ia+1:min([nsamp1,nsamp+ia]),:);
  else
    seisout.traces(1-ia:min([nsamp1-ia,nsamp]),1:ntr1)=seis1.traces(1:min([nsamp1,nsamp+ia]),:);  
  end

  ia=round((first-seis2.first)/step);
  if ia >= 0
    seisout.traces(1:min([nsamp2-ia,nsamp]),ntr1+1:ntr1+ntr2)= ...
               seis2.traces(ia+1:min([nsamp2,nsamp+ia]),:);
  else
    seisout.traces(1-ia:min([nsamp2-ia,nsamp]),ntr1+1:ntr1+ntr2)= ...
               seis2.traces(1:min([nsamp2,nsamp+ia]),:);  
  end
%  if seisout.first < seis1.first | seisout.first < seis2.first | ...
%     seisout.last  > seis1.last  | seisout.last  > seis2.last
%     seisout.null=null_value;
%  end
end
        
%      Set headers of output data set
if isfield(seis1,'headers') & isfield(seis2,'headers')   % Both data sets have headers
  seisout=merge_headers(seisout,seis1,seis2,htype,header_null_value);

elseif isfield(seis1,'headers')     % Only the first data set has headers
  seisout.header_info=seis1.header_info;
  if header_null_value == 0
    seisout.headers=[seis1.headers,zeros(size(seis1.headers,1),ntr2)];
  else
    seisout.headers=[seis1.headers,header_null_value*zeros(size(seis1.headers,1),ntr2)];
    if isnan(header_null_value)
      seisout.header_null=NaN;
    end
  end

elseif isfield(seis2,'headers')     % Only the second data set has headers
  seisout.header_info=seis2.header_info;
  if header_null_value == 0
    seisout.headers=[zeros(size(seis2.headers,1),ntr1),seis2.headers];
  else
    seisout.headers=[header_null_value*zeros(size(seis2.headers,1),ntr1),seis2.headers];
    if isnan(header_null_value)
      seisout.header_null=NaN;
    end
  end

else                                % Neither input data set has headers
   % Do nothing
end   
 
%     Copy common fields that have not yet been set and that have the same value 
%     in both input data sets (string or number, no arrays)
set_fields={'traces','headers','header_info','first','last','step', ...
            'units','history'};
fields1=fieldnames(seis1);
fields2=fieldnames(seis2);
fields=intersect(fields1,fields2);
fields=fields(~ismember(fields,set_fields));
for ii=1:length(fields);
  fval1=getfield(seis1,fields{ii});
  fval2=getfield(seis2,fields{ii});
        try
  if strcmpi(fval1,fval2)
     seisout=setfield(seisout,fields{ii},fval1);
  end
        catch
        end

        try
  if fval1 == fval2
     seisout=setfield(seisout,fields{ii},fval1);
  end
        catch
        end
end

%	Set null field, if necessary
if isfield(seis1,'null') & isfield(seis2,'null')
  if ~isnan(seis1.null) & ~isnan(seis2.null)
    if seis1.null == seis2.null
      param.null=seis1.null;
    end
  end
end
if ~isempty(find(isnan(seisout.traces)))
  if ~isfield(seisout,'null')
    seisout.null=param.null;
  end
else
  if isfield(seisout,'null')
    seisout=rmfield(seisout,'null');
  end
end


%    Append history field
if isfield(seis1,'history') & isfield(seis2,'history')
  seisout.history=seis1.history;
  seisout=s_history(seisout,'append',['[',ttype{1},' ',ttype{2},'], ',htype]);
  seisout=s_history(seisout,'merge',seis2.history);
end

%if isfield(seisout,'null') & ~isnan(seisout.null)
%   if isfield(seis1,'null') | isfield(seis2,'null')
%      seisout.null=NaN;
%   else
%      seisout=rmfield(seisout,'null');
%   end
%end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function seisout=merge_headers(seisout,seis1,seis2,htype,header_null_value)
% Function merges headers of the two seismic input data sets
% INPUT
% seisout seismic output data set with fields that have already been set
% seis1, seis2    seismic input data sets
% htype   This string parameter defines what should be done if the two data sets have
%         different headers
%         htype='u': Use the union of the headers (default)
%         htype='i'; Use the intersection of the headers (the headers they have in common)
% header_null_value   Null value for missing header values
% OUTPUT
% seisout seismic output data set
%     seisout=merge_headers(seis1,seis2,htype)

ntr1=size(seis1.headers,2);
ntr2=size(seis2.headers,2);

headers1=seis1.header_info(:,1);
headers2=seis2.header_info(:,1);
iheaders=intersect(headers1,headers2);
hidx1=find(ismember(headers1,iheaders));
hi1=seis1.header_info(hidx1,:);   
hidx2=find(ismember(headers2,iheaders));
hi2=seis2.header_info(hidx2,:);

%        Check that units of measurement are the same
lh1=length({hi1{:,2}});
if lh1 ~= sum(ismember({hi1{:,2}}',{hi2{:,2}}))
   spaces=blanks(lh1)';
   commas=char(44*ones(lh1,1));
   units_of_measurement=[char({hi1{:,1}}),spaces,char({hi1{:,2}}),commas, ...
                        spaces,char({hi2{:,2}}),spaces,char({hi2{:,1}})];
   disp(units_of_measurement)
   error('Units of measurement in common headers of the two input data sets differ')
end
seisout.header_info=hi1;
seisout.headers=[seis1.headers(hidx1,:),seis2.headers(hidx2,:)];

if strcmpi(htype,'u')
   if isempty(iheaders)
      idx1=1:length(headers1(:,1));
      idx2=1:length(headers2(:,1));
   else
      [dummy,idx1]=setdiff(headers1(:,1),iheaders(:,1));
      [dummy,idx2]=setdiff(headers2(:,1),iheaders(:,1));
   end
   if ~isempty(idx1) | ~isempty(idx2)
      seisout.header_null=header_null_value;
   end
   seisout.header_info=[seisout.header_info; ...
      seis1.header_info(idx1,:);seis2.header_info(idx2,:)];
   seisout.headers=[seisout.headers;
            seis1.headers(idx1,:),header_null_value*zeros(length(idx1),ntr2);
            header_null_value*zeros(length(idx2),ntr1),seis2.headers(idx2,:)];
elseif ~strcmpi(htype,'i')
   error(['Unknown parameter htype (',htype,')'])
end