readsegy.m

这是用matlab对segy数据进行处理

% 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=[];

endian='ieee-be'; % Big Endian is default


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
    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],0)
    segyid = fopen(filename,'r',endian);
else
    endian='ieee-be';
    SegymatVerbose([mfilename,' : ENDIAN SET TO ',endian],0)
    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
    SegymatVerbose([mfilename,' - Using supplied SegyHeader'])
end

% APPLY CHANGES TO SEGY HEADER IF NEEDE
if isempty(revision)==0,
    SegyHeader.SegyFormatRevisionNumber=revision;
    SegymatVerbose([mfilename,' - Manually set SEG-Y revision to ',num2str(revision)])
end
if isempty(dsf)==0,
    SegyHeader.DataSampleFormat=dsf;
end



% JUST SOME INFORMATION TO WRITE TO SCREEN :

Revision=SegyHeader.SegyFormatRevisionNumber;
if Revision>0, Revision=1; end

if (SegyHeader.DataSampleFormat>length(SegyHeader.Rev(Revision+1).DataSampleFormat));
    SegymatVerbose([mfilename,' : WARNING : YOU HAVE SELECTED (OR THE FILE IS FORMATTED SUCH THAT) A DATASAMPLE FORMAT THAT IS NOT DEFINED. \nREMEBER IEEE IS NOT SPECIFIED IN THE SEGY REV0 STANDARD !'])

    if (Revision==0)
        SegymatVerbose([mfilename,' : TRYING TO USE REVISION 1 AS OPPOSED TO REVISION 0'])