Skip to content

Latest commit

 

History

History
174 lines (127 loc) · 8.21 KB

File metadata and controls

174 lines (127 loc) · 8.21 KB

RFDFWIM — User Manual

This manual describes the workflow, file formats, main parameters, and module roles for the RFDFWIM (Radar Frequency-Domain Full-Waveform Inversion using MATLAB) package.


1. Workflow overview

1.1 Standard three-step workflow

  1. Create models and acquisition
    Run create_models_mkl.m to generate:

    • Permittivity model: input/MYMOD_eps.dat
    • Conductivity model: input/MYMOD_sig.dat
    • Source positions: input/MYMOD_source.dat
    • Receiver positions: input/MYMOD_receiver.dat
  2. Forward modelling (generate “observed” data)
    Run GPRFM.m. It:

    • Reads the above input files
    • Runs 2D FDFD forward modelling with CFS-PML
    • Saves synthetic data as obs/obs_data.mat (variable precrnew)
    • Optionally produces wavefield and model figures in output/
  3. Full-waveform inversion
    Run RFDFWI.m. It:

    • Loads obs/obs_data.mat as observed data
    • Builds a smooth initial model from the true model (Gaussian filter)
    • Runs FWI (gradient + optional Hessian, L-BFGS, Tikhonov, Wolfe line search)
    • Saves models, gradients, wavefields, and figures in outputmyINV/

This is an inverse crime setup: the same forward code generates the data and is used inside the inversion.


2. Directory and file roles

2.1 Root scripts

Script Role
create_models_mkl.m Builds two-cross (2C) synthetic models and 4-sided acquisition geometry; writes MYMOD_*.dat into input/.
GPRFM.m Forward modelling driver: reads input, runs FDFD for all shots/frequencies, saves obs/obs_data.mat and optional outputs in output/.
RFDFWI.m Main FWI driver: reads observed data and input models, runs iterative FWI, writes all inversion outputs to outputmyINV/.
inp_GPRmodel.m Defines model parameters used mainly for conductivity (grid, PML, file names, plotting).
inp_GPRmodel1.m Defines model and inversion parameters (grid, frequencies, PML, CFS-PML, L-BFGS, Wolfe, Tikhonov, paths).

2.2 Forward package (forward/)

File Role
forward_shotNEW.m Core forward routine: for each frequency, builds system matrix with CFS-PML, solves for wavefield, samples at receivers. Returns wavefield te and data precrnew.
RHS_TE1.m Builds right-hand side for TE mode for a given shot and frequency.
PML_9p_CFS_stag_para.m CFS-PML parameters for 9-point stencil.
imp_A_TE_9p_cfs_PML_stag1_para.m Assembles impedance matrix for TE, 9-point, CFS-PML, staggered grid (variant 1).
imp_A_TE_9p_cfs_PML_stag2_para.m Same as above, staggered variant 2 (MKL).
read_acq_jan.m Reads input/MYMOD_receiver.dat and input/MYMOD_source.dat; returns acquisition struct (e.g. acq.xrec, acq.yrec, acq.xshot, acq.yshot, acq.nxrec, acq.nyrec, acq.nxshot, acq.nyshot, acq.ntr, acq.nsrc).
shift_acq.m Shifts acquisition coordinates when PML is used (extended grid).
readeps1.m Loads permittivity model from input/MYMOD_eps.dat.
readsig1.m Loads conductivity model from input/MYMOD_sig.dat.
extend_eps1JAN.m Extends ε (and σ) model with PML padding.
extract_model_ext.m Extracts physical domain from extended (PML) model.
dotp.m Dot product utility (e.g. for inner products).

2.3 Inversion package (inv/)

File Role
grad_obj_hessMKLnew.m Computes forward wavefield, adjoint wavefield, data residual, L2 objective, gradients (σ, ε), and pseudo-Hessian (σ, ε) for TE mode.
cal_radargMKLnew.m Simulated data at receivers from forward wavefield.
cal_resjMKLnew.m Residual and L2 contribution per frequency/shot; optional Wiener weighting.
RHS_adj_TEm.m Right-hand side for adjoint TE solve.
RHS_TE_hess.m Right-hand side for Hessian-related solve.
ass_grad_TEMKLnew.m Assembles gradient contributions from wavefields.
copy_grad_frame_TE.m Copies gradient into frame (e.g. for PML handling).
norm_matrix_TE.m Norm of gradient matrix (for scaling).
scale_grad_TE.m Scales material parameters (σ, ε) for inversion.
Tikhonov_cost_TE.m Tikhonov regularisation term in cost function.
Tikhonov_grad_TE.m Gradient of Tikhonov term.
steep_descent.m Steepest-descent direction (negative gradient).
LBFGS_TEmNEW.m L-BFGS quasi-Newton update; returns preconditioned descent directions for σ and ε.
pseudo_hess_shin_ass_TEmNEW.m Pseudo-Hessian (Shin et al.) assembly for TE.
check_descent_TE.m Checks if search direction is a descent direction; can trigger L-BFGS reset.
wolfe_TENEW.m Wolfe line search; updates σ and ε and returns step length.
calc_mat_change_wolfe_multi_para_TE.m Material parameter update for multi-parameter Wolfe.
plot_wavefld.m Plots wavefield for figures.
dotp_matrix.m Matrix inner product.

2.4 Input file formats

  • MYMOD_eps.dat, MYMOD_sig.dat
    ASCII matrices: one value per grid point (row/column order as in create_models_mkl.m and readeps1/readsig1). Size must match model1.nx, model1.ny (e.g. 180×180).

  • MYMOD_receiver.dat, MYMOD_source.dat
    Two columns: x (distance), y (depth) in metres. One row per receiver or source. Read by read_acq_jan.m.

2.5 Output folders

  • output/ — Forward modelling figures and exports (e.g. true model, wavefield).
  • outputmyINV/ — Inversion: figures (e.g. figs/), saved models (models/), gradients (grad/, grad_tikh/, gradm/), wavefields (wavefield/), Hessian-related (Hgrad/), and .mat files (e.g. FinalL2.mat).
  • obs/ — Contains obs_data.mat (variable precrnew) produced by GPRFM.

3. Main parameters (inp_GPRmodel1.m and inp_GPRmodel.m)

3.1 Grid and domain

  • model1.nx, model1.ny — Number of grid points in x and y (e.g. 180).
  • model1.dh — Grid spacing [m] (e.g. 0.05).
  • model1.npml — PML thickness in grid points (e.g. 10).
  • model1.free — Free surface flag (0 = no free surface in tests).

3.2 Physics and constants

  • model1.eps0 — Vacuum permittivity (8.85e-12).
  • model1.sig0 — Reference conductivity for scaling (e.g. 5.6e-3).
  • model1.mu0 — Magnetic constant.

3.3 Frequencies (FWI)

In RFDFWI.m the frequency vector is set as:

f = [50e6 60e6 70e6 80e6 90e6 100e6 125e6 150e6 175e6 200e6];  % Hz

Simultaneous multi-frequency inversion (e.g. Lavoué et al., Layek & Sengupta).

3.4 Inversion and optimisation

  • itermax — Max FWI iterations (e.g. 1500).
  • model1.wolfe — 1 = use Wolfe line search; otherwise fixed step.
  • model1.LAMBDA_1, model1.LAMBDA_2 — Tikhonov weights.
  • model1.nlbfgs — L-BFGS history length.
  • model1.C1, model1.C2 — Wolfe constants.
  • model1.SRTRADIUS — Source taper radius [m].
  • model1.sigr_low, model1.sigr_high, model1.epsr_low, model1.epsr_high — Bounds for σ and ε.
  • model1.STAG — 1 = staggered variant 1; 2 = staggered variant 2.

3.5 Acquisition

  • model1.ACQMY — 1 = 4-sided acquisition (as in create_models_mkl).
  • Source/receiver files: input/MYMOD_source.dat, input/MYMOD_receiver.dat (see read_acq_jan.m).

4. Convergence and stopping (RFDFWI.m)

The main loop stops when any of the following holds:

  • ratio = L2(iter) / L2(1) ≤ 5e-5
  • iter = itermax
  • diff ≤ PRO (relative change of L2, with PRO=0 in the default script)

Final models are written and figures exported when the loop exits.


5. Colormap

  • seismic.map — MATLAB colormap file in the root directory; used for model and wavefield plots (e.g. load('seismic.map'); colormap(flipud(seismic))).

6. References (short)

  • Lavoué et al. (2014) — FWI of GPR data.
  • Layek & Sengupta (2019, 2021, 2023) — Multi-parameter GPR FDFD FWI.
  • Meles et al. (2011) — GPR FDTD adjoint, two-cross style setup.
  • Shin et al. (2001) — Pseudo-Hessian.
  • Nocedal & Wright (2006) — L-BFGS and Wolfe conditions.

For more detail on algorithms and parameters, see the header comments in the corresponding .m files and the README/INSTALLATION/COMMANDS documents.