FAST KIRCHHOFF CODE

Subroutine kirchslow() is very slow in the normal case where $x_{\max} \gg v t_{\max}$ because then many values of t are computed beyond $t_{\max}$.To avoid this, we notice that for fixed offset (ix-iy) and variable depth iz, as depth increases, time it eventually goes beyond the bottom of the mesh and, as soon as this happens, it will continue to happen for all larger values of iz. Thus we can break out of the iz loop the first time we go off the mesh to avoid computing anything beyond as shown in subroutine kirchfast(). (Some quality compromises, limiting the aperature or the dip, also yield speedup, but we avoid those.) Another big speedup arises from reusing square roots. Since the square root depends only on offset and depth, once computed it can be used for all ix. Finally, these changes of variables have left us with more complicated side boundaries, but once we work these out, the inner loops can be devoid of tests and in kirchfast() they are in a form that is highly optimizable by many compilers.

# Kirchhoff migration and diffraction.  (greased lightning)
#
subroutine kirchfast(   conj, add,  velocity, t0,dt,dx, nt,nx, modl, data)
integer ix,iz,it,ih,    conj, add,                      nt,nx,     xstart, xend
real	amp,t,z,h,                  velocity, t0,dt,dx,modl(nt,nx), data(nt,nx)
call conjnull(          conj, add,                     modl,nt*nx,  data,nt*nx)
do ih= -nx, nx {	h = dx * ih 			# h = offset
	do iz= 2, nt {	z = t0 + dt * (iz-1)		# z = travel-time depth
			t = sqrt( z**2  +  (h/velocity)**2 )
			it = 1.5 + (t - t0) / dt
	     if( it > nt )  break
		amp = (z / t)  *  sqrt( nt*dt / t )
		xstart =  1-ih;		xstart = max0(  1, xstart)
		xend   = nx-ih;		xend   = min0( nx, xend)
		if( conj == 0 )
			do ix= xstart, xend
				data(it,ix+ih)=data(it,ix+ih)+modl(iz,ix   )*amp
		else
			do ix= xstart, xend
				modl(iz,ix   )=modl(iz,ix   )+data(it,ix+ih)*amp
		}
	}
return; end

Originally the two Kirchhoff programs produced identical output, but finally I could not resist adding an important feature to the fast program, scale factors $z/t=\cos\theta$and $1/\sqrt{t}$ that are described elsewhere. The fast program, with a few changes, allows for velocity variation with depth. (The break should be conditioned on the maximum velocity followed by the if() conditioned on the velocity at the current depth.) When velocity varies laterally we need to go back to the beginning.

Figure 1 shows an example. The model includes dipping beds, syncline, anticline, fault, unconformity, and buried focus. The result is as expected with a ``bow tie'' at the buried focus. On a video screen, I can see hyperbolic events originating from the unconformity and the fault. At the right edge are a few faint edge artifacts. We could have reduced or eliminated these edge artifacts if we had extended the model to the sides with some empty space.

 

kfgood
Figure 1 Left is the model. Right is diffraction to synthetic data.