next up previous [pdf]

Next: Compiling Up: Fomel: RSF API Previous: Compiling

Fortran-77 interface

The Fortran-77 clip function is listed below.

	program Clipit
	implicit none
	integer n1, n2, i1, i2
	integer*8 in, out
	integer*8 sf_input, sf_output, sf_leftsize, sf_gettype
	logical sf_getfloat, sf_histint
	real clip, trace(1000)

	call sf_init()
	in = sf_input("in")
	out = sf_output("out")

	if (3 .ne. sf_gettype(in)) 
     &  call sf_error("Need float input")

	if (.not. sf_histint(in,"n1",n1)) then
	   call sf_error("No n1= in input")
	else if (n1 > 1000) then
	   call sf_error("n1 is too long")
	end if
	n2 = sf_leftsize(in,1)

	if (.not. sf_getfloat("clip",clip)) 
     &  call sf_error("Need clip=")

	do 10 i2=1, n2
	   call sf_floatread(trace,n1,in)

	   do 20 i1=1, n1
	      if (trace(i1) >  clip) then
		 trace(i1)=clip
	      else if (trace(i1) < -clip) then
		 trace(i1)=-clip
	      end if
 20	   continue

	   call sf_floatwrite(trace,n1,out)
 10	continue

	stop
	end
Let us examine it in detail.

The program starts with a call to sf_init, which initializes the command-line interface.

	call sf_init()
	in = sf_input("in")
The input and output files are created with calls to sf_input and sf_output. Because of the absence of derived types in Fortran-77, we use simple integer pointers to represent RSF files. Both sf_input and sf_output accept a character string, which may refer to a file name or a file tag. For example, if the command line contains vel=velocity.rsf, then both sf_input("velocity.rsf") and sf_input("vel") are acceptable. Two tags are special: "in" refers to the file in the standard input and "out" refers to the file in the standard output.


	if (3 .ne. sf_gettype(in)) 
RSF files can store data of different types (character, integer, floating point, complex). The function sf_gettype checks the type of data stored in the RSF file. We make sure that the type corresponds to floating-point numbers. If not, the program is aborted with an error message, using the sf_error function. It is generally a good idea to check the input for user errors and, if they cannot be corrected, to take a safe exit.


	if (.not. sf_histint(in,"n1",n1)) then
	   call sf_error("No n1= in input")
	else if (n1 > 1000) then
	   call sf_error("n1 is too long")
	end if
Conceptually, the RSF data model is a multidimensional hypercube. By convention, the dimensions of the cube are stored in n1=, n2=, etc. parameters. The n1 parameter refers to the fastest axis. If the input dataset is a collection of traces, n1 refers to the trace length. We extract it using the sf_histint function (integer parameter from history) and abort if no value for n1 is found. Since Fortran-77 cannot easily handle dynamic allocation, we also need to check that n1 is not larger than the size of the statically allocated array. We could proceed in a similar fashion, extracting n2, n3, etc. If we are interested in the total number of traces, like in the clip example, a shortcut is to use the sf_leftsize function. Calling sf_leftsize(in,0) returns the total number of elements in the hypercube (the product of n1, n2, etc.), calling sf_leftsize(in,1) returns the number of traces (the product of n2, n3, etc.), calling sf_leftsize(in,2) returns the product of n3, n4, etc. By calling sf_leftsize, we avoid the need to extract additional parameters for the hypercube dimensions that we are not interested in.


	if (.not. sf_getfloat("clip",clip)) 
The clip parameter is read from the command line, where it can be specified, for example, as clip=10. The parameter has the float type, therefore we read it with the sf_getfloat function. If no clip= parameter is found among the command line arguments, the program is aborted with an error message using the sf_error function.


	do 10 i2=1, n2
	   call sf_floatread(trace,n1,in)

	   do 20 i1=1, n1
	      if (trace(i1) >  clip) then
		 trace(i1)=clip
	      else if (trace(i1) < -clip) then
		 trace(i1)=-clip
	      end if
 20	   continue

	   call sf_floatwrite(trace,n1,out)
Finally, we do the actual work: loop over input traces, reading, clipping, and writing out each trace.



Subsections
next up previous [pdf]

Next: Compiling Up: Fomel: RSF API Previous: Compiling

2013-04-08