Ground-Penetrating Radar

Analysis of Tiwanaku’s Data

GPR Survey: Akapana, West Side (2002)

Step 1: Conversion from Waveforms to Intensity Maps

Done with proprietary software “GPR Process”

Step 2: Intensity Normalization

Necessary to eliminate “tiling” effects between patches surveyed under distinct environmental conditions / equipment settings


It is also necessary to eliminate intensity fluctuations across depth layers (caused by absorption / scattering of GPR waves by ground)

  1. Find tile boundaries (automatically)
  2. Compute median intensities
  3. Re-scale each tile
\

Step 3: Detection of Elongated Structures

Elongated structures are difficult to segment from “background”. They are strongly affected by GPR noise, broken into little pieces

Step 3a: Apply steerable filters to image

Each streerable filter is a detector for features with a particular orientation and thickness (outputs of 12 out of 360 shown below)

Step 3b: Line-voting (Hough transform)

Each pixel on filter outputs “votes” for one line that contains it Votes are accumulated on space of all lines (polar coordinates)

Step 3c: Line clipping into segments

  1. Convert lines back to image coordinates (pixels)
  2. Within each line, find segments with intense filter responses

Step 4: Multi-layer analysis

  1. Stack linear structures found at each depth layer
  2. Filter out structures that are “too dim” compared to others
  3. Mark features that are under other features with same orientation as “unreliable”

Click here to see video (Alt+Enter for full screen):

  • Reliable features appear as green “bricks”
  • Unreliable features appear as gray “bricks”
  • Almost all structures shallow (6-8ns)

Let’s now focus on an 80 x 80 tile with some deeper structures

  • Reliable features appear as green “bricks”
  • Unreliable features appear as gray “bricks”
  • Thin, linear structures recovered well
  • Thick, blob-like / massive structures missed

Click here to see video (Alt+Enter for full screen)

Step 5: Segmentation

Goal: recover thick structures in addition to thin ones Approach: use Normalized Cuts framework of J. Shi et al

  • Choose metric of (dis)similarity between pixels
  • Choose number of image regions (classes), K
  • Relax resulting graph cut problem (NP-complete) into generalized eigensystem (K-1 eigenvectors needed)
  • Discretize eigenvectors so that their elements become class labels (1 … K)
  • Compute mean intensity within each segmented region
  • Choose regions with mean intensity above threshold

Step 5a: Choice of similarity metric

Option 1: dissimilarity exponential on intensity difference

Good: requires less eigenvectors than other alternative (efficient)

Bad: results are too noisy (hard to visualize relevant structures)

Option 2: apply edge detector; dissimilarity increases with strength of edges crosses in the path between two pixels

Less efficient (more eigeinvectors), but less noisy (easier to visualize)

Step 5b: Number of regions

Step 5-final: Thresholding

Step 6: Visualization

Goal: convert segmented slices into smooth 3D surfaces Approach: use software developed by U. Penn. alumnus Marcelo Siqueira in his Ph.D. thesis (2005)

  • Stack binary slices intro binary volume (3D)
  • “Flip” small subset of voxels in 3D space to make volume well composed (i.e., singularity-free)
  • Approximate volume boundaries by smooth surfaces
  • Generate vertices on surfaces and triangulate them using algorithm with topologic & geometric guarantees

Final mesh, under Gouraud shading. Click here for movie

Compare this figure with previous one: thin wall not reconstructed

Final step: Use of Priors in thresholding

Goal:

Add to the 3D reconstruction the thin, approximately linear structures that are missed due to inaccurate segmentation

Approach:

Use results of Step 3 (Dectection of Elongated Structures) to build maps of prior constraints for binarization Steerable filters identical to those of Step 3a are slid along each detected linear structure:

How priors are used:

Binarization threshold is selectively lowered for points with high values (reddish) on map of priors

Priors in thresholding: Comparative results

Drawback after 3D volume is generated:

Resulting reconstruction is more cluttered, harder to visualize, due to the added set of thin, linear features

However:

By examining resulting volume carefully, relevant linear structures that were missing can now be found!

Status of Project

Stable set of tools for analysis and visualization has been produced

  • Results similar to those shown here have been obtained for entire W. Akapana set (42 patches, each 1600 sq m, 15 ns deep)

To do list:

  • Perform comparative measurements between reconstructions and ground-truth structures (on sites that have been excavated)
  • Develop and test methodology to solve inverse wave-propagation problem in order to (a) be able to infer electromagnetic properties of soil accurately (b) be able to account for effects of inter-reflections, absorption and scattering of GPR signals
  • Use estimated electromagnetic properties of soil in order to perform depth correction, which is essential to fuse 3D shapes of multiple patches into single reconstruction of entire area