# Gaussian packet migration for a given source-receiver configuration
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# see  Zacek, K. (2005): Gaussian packet prestack depth migration, 
#      In: Seismic Waves in Complex 3-D Structures, Report 15, pp.29-48.
#
# Calculating source-receiver coordinates, ray parameters
  NX=116  NT=750  DX=25.0  DT=0.004
  VSURF=1500.0
  FMIN=5.0  FMAX=55.0  FREF=35.0
  PARK=0.8  PARN=0.5   KAPPA=1.253314137
# Models mar-mod.dat and mar-modc.dat differ only by coordinates order,
# analogously for models mar-mo3.dat and mar-mo3c.dat
  MODEL='mar-modc.dat' # should be used for Fig. 7a in Zacek (2005)
 #MODEL='mar-mo3c.dat' # should be used for Fig. 7c in Zacek (2005)
 #MODEL='mar-modx.dat' # used by Zacek (2005), Figs. 7a and 7c
  RPAR='gpm-rp.out'
  REC='gpm-rec.out'
  SRC='gpm-src.out'
  gpstep:

# Ray tracing from a given source
  MODEL='mar-mod.dat'  # used by Zacek (2005), Fig. 7a
 #MODEL='mar-mo3.dat'  # used by Zacek (2005), Fig. 7c
  SRC='gpm-src.out'    # File with the coordinates of the initial point
  REC=' '              # No receivers in this example
  DCRT='mgp-crt1.dat'  # Numerical parameters for complete ray tracing
  CODE='codep.dat'     # Ray codes
  INIPAR=2             # Geographic-like spherical ray parameters
  RPAR='mgp-crt1.dat'  # Ray-parameter domain and control
  WRIT='mgp-writ.dat'  # Names and extent of output files
  CRTLOG='mgp-log.out' # Output log file
  crt:                 # Running program crt.for with these data

# Calculating travel times and amplitudes of the incident wave field
# N1=751  N2=2301  D1=4   D2=4
  N1=151  N2=461   D1=20  D2=20
  NUM='mtt-num.out'  MTT='mtt-tt.out'
  AMP='mtt-apr.out'  AMP2D3=1.0
  mtt:

# Back-propagation of the wavefield from receivers 
# and computation of migrated sections
  MODEL='mar-modc.dat' # used by Zacek (2005), Fig. 7a
 #MODEL='mar-mo3c.dat' # used by Zacek (2005), Fig. 7c
  SRC='gpm-rec.out'
  DCRT='mgp-crt2.dat'
  INIPAR=3
  RPAR='gpm-rp.out'
  crt:

  IXTZU=1  IXTZD=151  IXTZL=1  IXTZR=461
  DT=0.05  KAPPA=1.0
  GPMSEC='gpms.out'
  gpmig:

# Plotting migrated sections
  N1=2301  N2=751
  VREF=0   CREF=0.166667  VMIN=-999999000000.0  # VCIRC=1
  GRD='gpms.out'
  grdps: