A Quick Overview on Geospatial Data Visualization using PyGMT

Utpal Kumar, Feb 2022

Install new environment with necessary libraries

First way (only needs python3 installation)

python -m venv geoviz
source geoviz/bin/activate
pip install pygmt scipy #requires GMT installation for execution

Recommended Way

• Install Miniconda [Light-version of Anaconda]
• Create environment (it's always a good idea to keep separate environment for different projects)

conda create --name geoviz --channel conda-forge pandas pygmt scipy jupyter notebook

Contents

1. Plotting a topographic map
2. Plotting data points on a topographic map
3. Plotting focal mechanism on a map
   1. Using Strike, Dip, and Rake (in Harvard CMT convention)
   2. Using Seismic moment tensor
4. Plot tomographic results on a geographical map

Import libraries

import pygmt
import pandas as pd
import numpy as np
import xarray as xr
from scipy.interpolate import griddata

pygmt.show_versions()

PyGMT information:
  version: v0.5.0
System information:
  python: 3.9.10 | packaged by conda-forge | (main, Jan 30 2022, 18:04:04)  [GCC 9.4.0]
  executable: /opt/conda/bin/python
  machine: Linux-5.4.144+-x86_64-with-glibc2.31
Dependency information:
  numpy: 1.21.5
  pandas: 1.3.5
  xarray: 0.19.0
  netCDF4: 1.5.8
  packaging: 21.3
  ghostscript: 9.54.0
  gmt: 6.3.0
GMT library information:
  binary dir: /opt/conda/bin
  cores: 8
  grid layout: rows
  library path: /opt/conda/lib/libgmt.so
  padding: 2
  plugin dir: /opt/conda/lib/gmt/plugins
  share dir: /opt/conda/share/gmt
  version: 6.3.0

1. Plotting a topographic map

#define etopo data file
# topo_data = 'path_to_local_data_file'
topo_data = '@earth_relief_30s' #30 arc second global relief (SRTM15+V2.1 @ 1.0 km)
# topo_data = '@earth_relief_15s' #15 arc second global relief (SRTM15+V2.1)
# topo_data = '@earth_relief_03s' #3 arc second global relief (SRTM3S)

# Where's the data locally stored??

# define plot geographical range
minlon, maxlon = 60, 95
minlat, maxlat = 0, 25

# Visualization
fig = pygmt.Figure()

# make color pallets
pygmt.makecpt(
    cmap='topo',
    series='-8000/8000/1000',
    continuous=True
)

# plot high res topography
fig.grdimage(
    grid=topo_data,
    region=[minlon, maxlon, minlat, maxlat],
    projection='M4i',
    shading=True,
    frame=True
    )

# plot continents, shorelines, rivers, and borders
fig.coast(
    region=[minlon, maxlon, minlat, maxlat],
    projection='M4i',
    shorelines=True,
    frame=True
    )

# plot the topographic contour lines
fig.grdcontour(
    grid=topo_data,
    interval=4000,
    annotation="4000+f6p",
    limit="-8000/0", #to only display it below 
    pen="a0.15p"
    )

# Plot colorbar
fig.colorbar(
    frame='+l"Topography"',
#     position="x11.5c/6.6c+w6c+jTC+v" #for vertical colorbar
    )

# save figure
save_fig = 0
if not save_fig:
    fig.show() 
    #fig.show(method='external') #open with the default pdf reader
else:
    fig.savefig("topo-plot.png", crop=True, dpi=300, transparent=True)
#     fig.savefig("topo-plot.pdf", crop=True, dpi=720)
    print('Figure saved!')