sgywritv.m

这是matlab在地球物理数据处理方面的源码

function fid = sgywritv(sgyfile,traces,t,trahead)
% Program to automatically write  a SEGY file when the 101 element header is available.  FOR WRITING VISTA FORMAT SEGY FILES (hence 'v')
% Program returns matrix traces with amplitudes
% and time base t.
% Can easily add more variables if required
% Refer to SEG digital tape standards booklet.
% Written by D. Schmitt, January 1, 1999
% The input file must have the full filename including
% the standard .sgy designation at the end.
% To run [fileid] = sgywritv('filename',tracematrix,timebase,traceheadermatrix);
% Trace and Binary Headers are described. No Ascii Header is 
% written but just filled.


fclose('all');
if length(sgyfile) > 12
	error('Input filename is too long, please rewrite')
end
if sgyfile((length(sgyfile)-3):length(sgyfile)) ~= ['.sgy'] %This doesn't work
	error('Input filename must have .sgy appended')
end	

% fid = fopen(sgyfile,'w','b')  % for workstation mode
fid = fopen(sgyfile,'w','l');  % for standard PC Vista mode

[m,n] = size(traces);

l1 = ['C  1   University of Alberta 1998 Field School                                  '];
  count = fwrite(fid,l1,'char'); 
l1 = ['C  2   Allerston Data Set, Stations run W to E                                  '];
  count = fwrite(fid,l1,'char'); 
l1 = ['C  3   Receiver and Station Spacing = 4 m                                       '];
  count = fwrite(fid,l1,'char'); 
l1 = ['C  4   Station 1 at position 1 mile W of Hall                                   '];
  count = fwrite(fid,l1,'char'); 
l1 = ['C  5   Obs:  Beaty, Saeed, Whelehan, Reynard, Forde, Meglis, Klein, Osbourne    '];
  count = fwrite(fid,l1,'char'); 
l1 = ['C  6   NOTE NOTE - Files written by sgywritv - Vista headers, offset/position   '];
  count = fwrite(fid,l1,'char'); 
filler = ['0'];
for i = 1:(3200-480)  % 3200 byte ascii header
	count = fwrite(fid,filler,'char');  % 3200
end
 

% Binary header information, 400 bytes

binhead = [ 1  	% jobid
         1  	% line number *
	 1  	% reel number *
	n	% Number of Data Traces *
 	0 	% Number of Auxiliary Traces *
	(t(2)-t(1))*1e6 	% Sample interval in microseconds *
 	(t(2)-t(1))*1e6 	% Original sample interval in microseconds
 	m 	% Number of samples per trace *
	m 	% Original number of saples per trace
	1 	% Data Sample Format 1=Float(4), 2=Fixed(4),3 = Fixed(2)*
	24	% CDP Fold (expected per CDP ensemble)
 	1	% Trace sorting code, 1 = as recorded
 	1	% Vertical Sum Code (Stack?), N = number of vertical sums. 
 	0	% Sweep Frequency at Start
	0	% Sweep Frequency at End
	0	% Sweep Length (seconds)
	0	% Sweep Type Code
	0	% Trace Number of Sweep Channel
	0	% Sweep Taper Length at start (in ms)
	0	% Sweep Taper Length at End (in ms)
	0	% Taper Type
	0	% Correlated Data Traces 1 = No 2 = Yes
	0	% Binary Gain Recovered 1 = Yes  2 = No
	0	% Amplitude Recovery Method 1=none, 2=spherical, 3=AGC, 4=other
	1	% Measurement system 1= meters, 2 = feet
	1	% Polarity 1 = increase P or upward case motion gives negative
	0];	% Vibratory polarity code - see SEG specs.
binhead(28:197) = zeros(170,1);  % Fill out remainer of binary header

count = fwrite(fid,binhead(1:3),'int32');  %  The first three are 4 byte
count = fwrite(fid,binhead(4:197),'int16');  % Remainder are in 2 byte form 

% Write trace headers from trahead and data from traces.

for i = 1:n

	count = fwrite(fid,trahead(1:7,i),'int32');  % Bytes 1-28
	count = fwrite(fid,trahead(8:11,i),'int16'); % Bytes 29-36
	count = fwrite(fid,trahead(12:19,i),'int32'); % Bytes 37-68
	count = fwrite(fid,trahead(20:21,i),'int16'); % Bytes 69-72
	count = fwrite(fid,trahead(22:25,i),'int32'); % Bytes 73-88
	count = fwrite(fid,trahead(26:71,i),'int16'); % Bytes 89-180
	count = fwrite(fid,trahead(72:101,i),'int16'); % Fill bytes, unassigned 181-240 

	fwrite(fid,traces(:,i),'float'); %  Write trace data
end
fclose(fid);