(15) Gridding, contouring, and masking of unconstrained areasΒΆ
This example demonstrates some off the different ways one can use to grid data in GMT, and how to deal with unconstrained areas. We first convert a large remote ASCII file to binary with gmtconvert since the binary file will read and process much faster. Our lower left plot illustrates the results of gridding using a nearest neighbor technique (nearneighbor) which is a local method: No output is given where there are no data. Next (lower right), we use a minimum curvature technique (surface) which is a global method. Hence, the contours cover the entire map although the data are only available for portions of the area (indicated by the gray areas plotted using mask). The top left scenario illustrates how we can create a clip path (using mask) based on the data coverage to eliminate contours outside the constrained area. Finally (top right) we simply employ coast to overlay gray land masses to cover up the unwanted contours, and end by plotting a star at the deepest point on the map with plot. This point was extracted from the grid files using grdinfo. Placement of the four panels is simplified via subplot.
#!/usr/bin/env bash
# GMT EXAMPLE 15
#
# Purpose: Gridding and clipping when data are missing
# GMT modules: blockmedian, gmtconvert, grdclip, grdcontour, grdinfo
# GMT modules: info, nearneighbor, coast, mask, surface, plot
# Unix progs: rm
#
gmt begin ex15
gmt convert @ship_15.txt -bo > ship.b
region=`gmt info ship.b -I1 -bi3d`
gmt subplot begin 2x2 -M0.1i/0.05i -Fs2.9i/0 $region -JM2.9i -BWSne -T"Gridding with missing data"
# Raw nearest neighbor contouring
gmt nearneighbor $region -I10m -S40k -Gship.nc ship.b -bi
gmt grdcontour ship.nc -JM -C250 -A1000 -Gd2i -c1,0
# Grid via surface but mask out area with no data using coastlines
gmt blockmedian ship.b -b3d > ship_10m.b
gmt surface ship_10m.b -Gship.nc -bi
gmt mask -I10m ship.b -T -Glightgray -bi3d -c1,1
gmt grdcontour ship.nc -C250 -L-8000/0 -A1000 -Gd2i
# Grid via surface but mask out area with no data
gmt mask -I10m ship_10m.b -bi3d -c0,0
gmt grdcontour ship.nc -C250 -A1000 -L-8000/0 -Gd2i
gmt mask -C
# Clip data above sealevel then overlay land
gmt grdclip ship.nc -Sa-1/NaN -Gship_clipped.nc
gmt grdcontour ship_clipped.nc -C250 -A1000 -L-8000/0 -Gd2i -c0,1
gmt coast -Ggray -Wthinnest
gmt grdinfo -Cn -M ship.nc | gmt plot -Sa0.15i -Wthick -i10,11
gmt subplot end
gmt end show
rm -f ship.b ship_10m.b ship.nc ship_clipped.nc