More Preparation Work
Many of the GMT tools require the area that they operate on to be defined in the command line. We could use the area that is shown by our Stats program, or we could use the GMT tool minmax:
$ minmax smoothed.xyzt smoothed.xyzt: N = 34299 <-61.2469/-61.236> <12.9999/13.0127> <0.850357/28.01> <1.14098e+09/1.14434e+09>
We are interested in the first two groups: longitude (X) and latitude (Y). It isn't intuitive how big an area this is because we are looking a section on ellipsoid. We also have to make a decision on how big a grid size we want to use when we use the GMT program nearneighbor.
Since our data spans less than a nautical mile, we can think in terms of a grid that is in the just a few meters. But GMT wants the grid size in angles of geographic arc, so let's say, 15 meters which is approximately .008 nm or .008 minutes latitude. We can express this as .48 seconds of latitude. Since we are near the equator, for convenience we will also use 0.48 seconds of longitude.
$ blockmedian smoothed.bxyz -I0.48c/0.48c -R-61.2469/-61.236/12.9999/13.0127 -G smoothed.grd -S0.005k
But we have errors:
nearneighbor: Encountered first invalid record near/at line # 34302 nearneighbor: Likely causes: nearneighbor: (1) Invalid x and/or y values, i.e. NaNs or garbage in text strings. nearneighbor: (2) Incorrect data type assumed if -J, -f are not set or set incorrectly. nearneighbor: (3) The -: switch is implied but not set. nearneighbor: (4) Input file in multiple segment format but the -M switch is not set. nearneighbor: Mismatch between actual (0) and expected (3) fields near line 34299
It turns out that GMT is not so simple to use for quick mapping. Since I use GMT so seldom, I thought I should study my particular use of GMT and automate the process with yet another Perl script: bp which takes the xyzt files all the way through some preprocessing steps and on to surface plots with several options.
bp - BathyPlot Perl Script
You can see the steps, arguments, etc. required to bring the xyz data all the way through, if you tell bp to spit that information out by using the --verbose option. The lines starting with an exclamation (!) are GMT commands that are issued by the script. The "Values:" are parameters that the script has set for itself, either from data, command line options, or a configuration/history file that the script maintains. The rest is information from the GMT programs. Notice how involved the arguments are for the GMT programs (hyperlinked to their man pages).
! minmax smoothed.bxyz -C minmax: -61.2469 -61.236 12.9999 13.0127 -28.01 -0.850357 Values: Region => -61.2469/-61.236/12.9999/13.0127 BACKGROUNDClr => 128/128/255 clean => 1 contour => 1 contourAnnt => 5 contourIntv => 1 cptFile => shallowWater.cpt cptTable => shallowWater filter1Dtype => median filter1Dwidth => 3 gradAz => 300 gradNorm => t.65 gridFilter => median gridMethod => surface gridSize => 0.484444444c/0.485052632c gridSizeM => 15 inputFields => 1,2,3, landClr => white maxLat => 13.0127 maxLong => -61.236 maxZ => -0.850357 minLat => 12.9999 minLong => -61.2469 minZ => -28.01 negateZ => 1 over => 1.25 pageX => 8.5 pageY => 11 portrait => 1 projection => M7i psFile => smoothed.ps surfConvergance => 0.25 surfTension => 0.8 symcptTable => bright useBinary => 1 verbose => 1 xMeters => 1172.7 yMeters => 1424.3 ! blockmedian smoothed.bxyz -I0.484444444c/0.485052632c \ -R-61.2469/-61.236/12.9999/13.0127 -fg -C -bo3 >smoothed.blk ! surface smoothed.blk -T0.8 -C0.5 -N50 -fg -Gsmoothed.osgrid \ -R-61.2469/-61.236/12.9999/13.0127 -I0.484444444c/0.485052632c -bi3 ! grdcut smoothed.osgrid -R-61.2469/-61.236/12.9999/13.0127 \ -Gsmoothed.grid --ANNOT_FONT_SIZE_PRIMARY=10p --ANNOT_FONT_SIZE_SECONDARY=12p \ --BASEMAP_TYPE=fancy+ --COLOR_NAN=0/0/0 --D_FORMAT=%12lg --ELLIPSOID=WGS-84 \ --GRID_CROSS_SIZE_PRIMARY=4p --GRID_CROSS_SIZE_SECONDARY=4p --GRID_PEN_PRIMARY=thin,white \ --LABEL_OFFSET=0p --PAPER_MEDIA=A3 --PLOT_DEGREE_FORMAT=ddd.mmF --TICK_LENGTH=0p ! grdgradient smoothed.grid -G1 -A280 -A300 -Nt.65 -V >>1 GMT_grd_is_global: no! grdgradient: Min Mean Max sigma intensities: -51907.4 1170.4 63297.77261.54 ! grdview smoothed.grid -N-1/white -R-61.2469/-61.236/12.9999/13.0127 \ -JM7i -CshallowWater.cpt -E180/90 -Qs -I1 -S -K -P -V --ANNOT_FONT_SIZE_PRIMARY=10p \ --ANNOT_FONT_SIZE_SECONDARY=12p --BASEMAP_TYPE=fancy+ --COLOR_NAN=0/0/0 --D_FORMAT=%12lg \ --ELLIPSOID=WGS-84 --GRID_CROSS_SIZE_PRIMARY=4p --GRID_CROSS_SIZE_SECONDARY=4p \ --GRID_PEN_PRIMARY=thin,white --LABEL_OFFSET=0p --PAPER_MEDIA=A3 --PLOT_DEGREE_FORMAT=ddd.mmF \ --TICK_LENGTH=0p >>smoothed.ps GMT_grd_is_global: no! grdview: Processing shape file GMT_grd_is_global: no! grdview: Find contours grdview: Trace and bin contours... GMT_grd_is_global: no!ontour interval 0 grdview: GMT WARNING: (w-x_min) must equal (NX + eps) * x_inc), where NX is an integer and |eps| <= 0.0001. grdview: GMT WARNING: w reset to 298.753 grdview: Processing illumination file GMT_grd_is_global: no! grdview: GMT WARNING: (w-x_min) must equal (NX + eps) * x_inc), where NX is an integer and |eps| <= 0.0001. grdview: GMT WARNING: w reset to 298.753 grdview: Start creating PostScript plot grdview: Done! ! grdcontour smoothed.grid -A5+k128/50/50+c5+jCB -C1+k255/255/255+c5+jCB \ -S2 -R-61.2469/-61.236/12.9999/13.0127 -JM7i -Wa1p255/100/100 -Wc0.25p0/128/128 \ -O -K -P -V --ANNOT_FONT_SIZE_PRIMARY=10p --ANNOT_FONT_SIZE_SECONDARY=12p \ --BASEMAP_TYPE=fancy+ --COLOR_NAN=0/0/0 --D_FORMAT=%12lg --ELLIPSOID=WGS-84 \ --GRID_CROSS_SIZE_PRIMARY=4p --GRID_CROSS_SIZE_SECONDARY=4p --GRID_PEN_PRIMARY=thin,white \ --LABEL_OFFSET=0p --PAPER_MEDIA=A3 --PLOT_DEGREE_FORMAT=ddd.mmF \ --TICK_LENGTH=0p >>smoothed.ps grdcontour: Allocate memory and read data file GMT_grd_is_global: no! GMT_grd_is_global: no! grdcontour: Tracing the -28 contour grdcontour: Tracing the -27 contour grdcontour: Tracing the -26 contour grdcontour: Tracing the -25 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = -25 grdcontour: Tracing the -24 contour grdcontour: Tracing the -23 contour grdcontour: Tracing the -22 contour grdcontour: Tracing the -21 contour grdcontour: Tracing the -20 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = -20 grdcontour: Tracing the -19 contour grdcontour: Tracing the -18 contour grdcontour: Tracing the -17 contour grdcontour: Tracing the -16 contour grdcontour: Tracing the -15 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = -15 grdcontour: Tracing the -14 contour grdcontour: Tracing the -13 contour grdcontour: Tracing the -12 contour grdcontour: Tracing the -11 contour grdcontour: Tracing the -10 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = -10 grdcontour: Tracing the -9 contour grdcontour: Tracing the -8 contour grdcontour: Tracing the -7 contour grdcontour: Tracing the -6 contour grdcontour: Tracing the -5 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = -5 grdcontour: Tracing the -4 contour grdcontour: Tracing the -3 contour grdcontour: Tracing the -2 contour grdcontour: Tracing the -1 contour grdcontour: Tracing the 0 contour grdcontour: Warning: Your -Gd|D option produced no contour labels for z = 0 grdcontour: Done! ! psbasemap -R-61.2469/-61.236/12.9999/13.0127 -JM7i -Bpf1ca10cWSne \ -Bsf10c -O -P -V --ANNOT_FONT_SIZE_PRIMARY=10p --ANNOT_FONT_SIZE_SECONDARY=12p \ --BASEMAP_TYPE=fancy+ --COLOR_NAN=0/0/0 --D_FORMAT=%12lg --ELLIPSOID=WGS-84\ --GRID_CROSS_SIZE_PRIMARY=4p --GRID_CROSS_SIZE_SECONDARY=4p --GRID_PEN_PRIMARY=thin,white \ --LABEL_OFFSET=0p --PAPER_MEDIA=A3 --PLOT_DEGREE_FORMAT=ddd.mmF \ --TICK_LENGTH=0p >> smoothed.ps psbasemap: Constructing basemap
Perl Scripts
- bp BandPass - a Perl script for quickly generating plots of bathymetric data from xyzt tables