Next: Internals of the film-loop Up: WAVEMOVIE PROGRAM Previous: Earth surface boundary condition

## Frames changing with time

For a film loop to make sense to a viewer, the subject of the movie must be periodic, and organized so that the last frame leads naturally into the first. In the movie created by wavemovie() there is a parameter lambda that controls the basic repetition rate of wave pulses fired onto the screen from the top. When a wavelet travels one-quarter of the way down the frame, another is sent in. This is defined by the line

Take any point in (x,z)-space. The signal there will be a superposition of sinusoids of various frequencies, .We can choose what frequencies we will use in the calculation and what amplitudes and phases we will attach to the initial conditions at those frequencies. Here we will simply take uniformly spaced sinusoids of unit amplitude and no phase. The nw frequencies are ,,..., nw.The lowest frequency dw = must be inversely proportional to the wavelength lambda

Finally, the time duration of the film loop must equal the period of the lowest-frequency sinusoid

This latter equation defines the time interval on the line

If you use more frequencies, you might like the result better because the wave pulses will be shorter, and the number of wavelengths between the pulses will increase. Thus the quiet zones between the pulses will get quieter. The frequency components can be weighted differently--but this becomes a digression into simple Fourier analysis.
# from par:  integer n3:nt=12, n2:nx=48, n1:nz=96,  nw=2, nlam=4
# from par:  real dx=2, dz=1, v=1
#
subroutine wavemovie( nt, nx, nz, nw, nlam, dx,dz,v, p, cd, q)
integer it,nt,ix,nx,iz,nz,iw,nw, nlam
real dx,dz,v, phase,pi2,z0,x0,dt,dw,lambda,w,wov,x, p(nz,nx,nt)
complex aa,a,b,c,cshift, cd(nx),q(nx)
lambda=nz*dz/nlam;  pi2=2.*3.141592; dw=v*pi2/lambda;  dt=pi2/(nt*dw)
x0 = nx*dx/3; z0 = nz*dz/3
call null( p, nz*nx*nt)
do iw = 1,nw {                                  # superimpose nw frequencies
w =  iw*dw;   wov = w/v                 # frequency / velocity
do ix = 1,nx {                          # initial conditions for a
x = ix*dx-x0;                   # collapsing spherical wave
phase = -wov*sqrt( z0**2+x**2)
q(ix) = cexp( cmplx( 0.,phase))
}
aa = (0.,1.)*dz/(4.*dx**2*wov)          # tridiagonal matrix coefficients
a = -aa;      b = 1.+2.*aa;      c = -aa
do iz = 1,nz {                          # extrapolation in depth
do ix = 2,nx-1                  # diffraction term
cd(ix) = aa*q(ix+1) + (1.-2.*aa)*q(ix) + aa*q(ix-1)
cd(1) = 0.;      cd(nx) = 0.
call ctris( nx,-a,a,b,c,-c,cd,q)
# Solves complex tridiagonal equations
cshift = cexp( cmplx( 0.,wov*dz))
do ix = 1,nx                    # shifting term
q(ix) = q(ix) * cshift
do it=1,nt {                    # evolution in time
cshift = cexp( cmplx( 0.,-w*it*dt))
do ix = 1,nx
p(iz,ix,it) = p(iz,ix,it) + q(ix)*cshift
}
}
}
return;         end


Next: Internals of the film-loop Up: WAVEMOVIE PROGRAM Previous: Earth surface boundary condition
Stanford Exploration Project
12/26/2000