readsegy.m

matlab编写的地震数据格式转换

% ReadSegy : Reads a SEG Y rev 1 formatted file
%
% Call :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename);
%
% To read time slice 0.5<t<5 :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'trange',.5,3);
% Skip every 5th trace :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'jump',5);
% Read data in a CDP header range : 5000<cdp<5800 :
% (change cdp to any other valid TraceHeader value)
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'minmax','cdp'5000,5800);
% Read only the header values (Data will return empty)
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'SkipData',1);
%
% SEG-Y format revision number can be '0' (1975) or 
% '100' (similar to '1') (2002).
% By default the SEG-Y format revision number is read in the 
% binary header, but this can be overruled using :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'revision',0);
%
% Read using a specific Data Sample Format :
% Rev 0, IBM FLOATING POINT
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'revision',0,'dsf',1);
% Rev 1, IEEE FLOATING POINT
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'revision',1,'dsf',5);
%
% A SegyHeader can be forced on the SEG-Y file using :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'SegyHeader',SegyHeader);
% The SegyHeader can be obtain by GetSegyHeader(segyfilename), and
% then edited.
%
% To read using little endian :
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'endian','l');
%
% Combine any combination of the above
% [Data,SegyTraceHeaders,SegyHeader]=ReadSegy(filename,'jump',1,'minmax','cdp',5300,5400);
%
%
% Plot the data using e.g. 
% imagesc([SegyTraceHeaders.cdp],SegyHeader.time,Data);
% wiggle([SegyTraceHeaders.TraceNumber],SegyHeader.time,Data);
%
% (C) 2003-2004, Thomas Mejer Hansen, [e-mail]
%

% Implemented using the syntax of the SEG-Y revised format :
% SEGY-Y rev 0, SEG-Y rev 1 as described in 
% http://seg.org/publications/tech-stand/
%
% Extended Textual Header is not yet tested
% If you would like it implemented, please send me an SEGY file with
% that sort of information, as well as a description of the segy file
%
%
% (C) 2001-2004, Thomas Mejer Hansen, [e-mail]/[e-mail]
% 
%    This program is free software; you can redistribute it and/or modify
%    it under the terms of the GNU General Public License as published by
%    the Free Software Foundation; either version 2 of the License, or
%    (at your option) any later version.
%
%    This program is distributed in the hope that it will be useful,
%    but WITHOUT ANY WARRANTY; without even the implied warranty of
%    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%    GNU General Public License for more details.
%
%    You should have received a copy of the GNU General Public License
%    along with this program; if not, write to the Free Software
%    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
%
%
% 
%
% 0.1 : INitial Release
% 0.2 : Added SkipData var, to skip reading of data.
% 0.3 : May 01, 2002 
%       Added ability to read in ever 'jump' traces.
%       Added ability to read in time range.
%       Added abiliy to read in header range (ex. mincdp to maxcdp).
%

% TODO : WHEN READING ONLY PART OF DATASET MAKE SURE TO ADJUST THE SEGY
% HEADER ACCORDINGLY !!!!!!

function [Data,SegyTraceHeaders,SegyHeader,HeaderInfo]=ReadSegy(filename,varargin);


  if isoctave
    doWaitBar=0; % [1] show progress bar gui
    mfilename='ReadSegy';
  else
    doWaitBar=1;
    mfilename='ReadSegy';
  end
    
  dsf=[];
  revision=[];
  endian_tight=[];
  tmin=[];tmax=[];
  headermin=[];headermax=[];header=[];
  jump=[];
  SkipData=[];
  
  
  SegymatVerbose([mfilename,' : reading ',filename])
  
  if ~(exist(filename)==2'),
    SegymatVerbose([mfilename,' : ', filename,' does not exist !'])
    Data=[];SegyTraceHeaders=[];SegyHeader=[];HeaderInfo=[];
    return
  end
  
  % IF ONLY 'filename', AND one outpuet HAS BEEN
  % SPECIFIED AS IN/OUTPUT, THEN USE THE FAST
  % ALGORITHM FOR READING.
  if (nargin==1)&(nargout==1)
    [Data]=ReadSegyFast(filename);
    return
  end
  
  SegymatVerbose([mfilename,' - Checking Varargin'],90)
  
  ninput=nargin;
  % NEXT TWO LINES TO ENUSRE THAT VARARGIN CAN BE PASSED TO FUNCTION
  if ninput==2
    % CALL USING VARARGIN
    ninput=1+length(varargin{1});
    varargin=varargin{1};
  else
    % DIRECT CALL
    ninput=length(varargin);
  end
  
  
  
  % TRANSFORM VARARGING INTO PARAMETERS
  cargin=1;
  while (cargin<ninput)
    
    SegymatVerbose([mfilename,' - Converting varargin, ',num2str(cargin)],90)
    
    if strcmp(varargin{cargin},'jump')
      cargin=cargin+1;
      eval(['jump=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['JUMP : Read only every ',num2str(jump),'th trace'])
    end
    
    if strcmp(varargin{cargin},'minmax')
      cargin=cargin+1;
      eval(['header=''',varargin{cargin},''';']);
      cargin=cargin+1;
      eval(['headermin=',num2str(varargin{cargin}),';']);
      cargin=cargin+1;
      eval(['headermax=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['MIN MAX : Using header ',header,' from ',num2str(headermin),' to ',num2str(headermax)])
    end    
    
    if strcmp(varargin{cargin},'trange')
      cargin=cargin+1;
      eval(['tmin=',num2str(varargin{cargin}),';']);
      cargin=cargin+1;
      eval(['tmax=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['TRANGE : tmin=',num2str(tmin),' tmax=',num2str(tmax)])
    end    
    
    % ENDIAN FORMAT 
    endian='ieee-be'; % Big Endian is default
    if strcmp(varargin{cargin},'endian')
      cargin=cargin+1;
      eval(['endian_tight=varargin{cargin};'])
      if endian_tight=='l',
        SegymatVerbose(['USING LITTLE ENDIAN TYPE'])
        endian='ieee-le';
      else
        SegymatVerbose(['USING BIG ENDIAN TYPE'])
      end
    end    
    
    
    if strcmp(varargin{cargin},'revision')
      cargin=cargin+1;
      eval(['revision=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['USING REVISION : rev=',num2str(revision)])
    end    
    
    if strcmp(varargin{cargin},'dsf')
      cargin=cargin+1;
      eval(['dsf=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['USING Data Sample Format : dsf=',num2str(dsf)])
    end    
    
    if strcmp(varargin{cargin},'SkipData')
      cargin=cargin+1;
      eval(['SkipData=',num2str(varargin{cargin}),';']);
      SegymatVerbose(['SKIPPING DATA - READING ONLY HEADERS'])
    end    
    
    if strcmp(varargin{cargin},'SegyHeader')
      cargin=cargin+1;
      SegyHeader=varargin{cargin};
      SegymatVerbose(['USING LOADED SEGYHEADER'])
    end    
    
    cargin=cargin+1;
    
  end
  
  
  if isempty(SkipData)==1,
    SegymatVerbose([mfilename,' : Skip data is not set (defautls to 0)'],90)
    SkipData=0; % [0] READ ONLY HEADER VALUES, [1] READ IN ALL DATA
  end
  
  
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % OPEN FILE HANDLE
  
  if exist('endian')==1,
    SegymatVerbose([mfilename,' : ENDIAN : ',endian],90)
    segyid = fopen(filename,'r',endian);   
  else
    endian='ieee-be';
    SegymatVerbose([mfilename,' : ENDIAN SET TO ',endian],90)
    segyid = fopen(filename,'r','ieee-be');  % ALL DISK FILES ARE IN BIG
  end                                        % ENDIAN FORMAT, ACCORDING TO 
                                             % SEGY Y rev 1
  
  
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % BINARY HEADERS
  if exist('SegyHeader')==0
    SegyHeader=GetSegyHeader(segyid);
  else