Grdimage centering longitude around 0

Q&A
1

I am trying to use grdimage to create a geotiff from a netcdf bathymetry file. The problem I’m having is that grdimage is reprojecting the output geotiff. I set the projection with -J and that has fixed the latitudes but has centered the longitudes around 0.

The command I am using is:

gmt grdimage rose_atoll_40m.nc -Cabyss.cpt -Igradient.nc -JEPSG:4326 -Arose_atoll_40m_image.tif

Is there a way I can preserve the coordinate system from the original netcdf and the gradient file?

2

Maybe not use this? I mean no -J at all.

3

Thanks for your reply, Joaquim. I added the -J flag because without it I was also having problems. Without adding -J option I am looking at the coordinates with gdalinfo:

$ gmt grdimage rose_atoll_40m.nc -Cabyss.cpt -Igradient.nc -Arose_atoll_40m_image.tif

$ gdalinfo rose_atoll_40m.nc
…
Corner Coordinates:
Upper Left (-168.2660964, -14.3795543) (168d15’57.95"W, 14d22’46.40"S)
Lower Left (-168.2660964, -14.6337019) (168d15’57.95"W, 14d38’ 1.33"S)
Upper Right (-168.0027937, -14.3795543) (168d 0’10.06"W, 14d22’46.40"S)
Lower Right (-168.0027937, -14.6337019) (168d 0’10.06"W, 14d38’ 1.33"S)
Center (-168.1344451, -14.5066281) (168d 8’ 4.00"W, 14d30’23.86"S)
…

$ gdalinfo rose_atoll_40m_image.tif
…
Corner Coordinates:
Upper Left ( -14638.981,-1598935.292) (Invalid angle,Invalid angle)
Lower Left ( -14638.981,-1627195.244) (Invalid angle,Invalid angle)
Upper Right ( 14638.981,-1598935.292) (Invalid angle,Invalid angle)
Lower Right ( 14638.981,-1627195.244) (Invalid angle,Invalid angle)
Center ( -0.000,-1613065.268) ( 0d 0’ 0.00"W,Invalid angle)
…

4

If you are specifying a projection different from the input grid, then I think you also have to specify the output region.

5

This is a great suggestion - unfortunately I am still losing longitude information with grdimage even after adding -R and -J values.

In the mean time I’ve found a workaround by generating a .png (no spatial information) and then converting it to a geotiff with gdal_translate:

gdal_translate -of GTiff -a_srs EPSG:4326 -a_ullr -168.2660964 -14.3795543 -168.0027937 -14.6337019 rose_atoll_40m_image.png rose_atoll_40m_image.tif

6

I use the GeoTIFF output from both GMT and PyGMT (recently in the developer version) all the time. In GMT, I don’t do it directly from the grdimage program with -A but do a regular EPS output and then run gmt psconvert -A -P -Tt -W+g $mapOut after completing the plot. This works well, even plots with different map projections.