PyGMT fails after conda update

[Sistjerne]

Dear PyGMT support,

I get an error in an old code after a coda update.

pygmt.exceptions.GMTCLibError: Module 'coast' failed with status code 74:
coast [ERROR]: Map region exceeds 360 degrees
coast [ERROR]: General map projection error

The example is heavily inspired by: https://nbviewer.org/github/weiji14/nzasc2021/blob/main/key_figure.ipynb#Figure-of-Antarctic-active-subglacial-lake-map and the code is here:

import pygmt
import numpy as np
 
arctic_dem = path_to_arctic_dem
back_img = arctic_dem
figheight = 115  # in mm
xl, xh, yl, yh = [-799890.9468068393, 940434.5640275762, -3464708.338152134, -542960.1486965496]
ratio = (yh - yl) / (figheight / 1000)

# Make a GMT region string and projection strings in both S and Lon/Lat
projX_gl = "s-45/90/70/1:"
regionX  = [xl, xh, yl, yh]
regionr  = "-58/58/15/80r"
projX    = "x1:" + str(ratio)
projS    = projX_gl + str(ratio)

fig = pygmt.Figure()
fig.coast(region=regionX, projection=projS, resolution="f", water=True)
fig.grdimage(
    grid="@earth_relief_03m",
    projection=projS,
    region=regionX,
    cmap="grayC",
    shading=True,
)
fig.coast(Q=True)  # end water clip path

fig.basemap(projection=projS, region=regionX, frame=["xa45g45", "ya10g10"])

pygmt.makecpt(
    series=[-5000, 5000, 1],
    cmap="grayC",
    continuous=True,
    reverse=True,
)
fig.grdimage(grid=back_img, projection=projX, cmap=True, nan_transparent=True, shading=False, transparency=30)

fig.plot(
        x          = np.array([-45, -32]),
        y          = np.array([62, 77]), 
        region     = regionr,
        projection = projS,
        fill       = "red",
        style      = "c0.5c",
        )
fig.savefig("test.png", show=True)

System information:

>>> import pygmt
>>> pygmt.show_versions()
PyGMT information:
  version: v0.11.1.dev9
System information:
  python: 3.12.1 | packaged by conda-forge | (main, Dec 23 2023, 08:03:24) [GCC 12.3.0]
  executable: /home/s3mpc/anaconda3/envs/s3mpc/bin/python
  machine: Linux-3.10.0-1160.42.2.el7.x86_64-x86_64-with-glibc2.17
Dependency information:
  numpy: 1.26.4
  pandas: 2.2.0
  xarray: 2024.1.1
  netCDF4: 1.6.5
  packaging: 23.2
  contextily: None
  geopandas: 0.14.3
  ipython: None
  rioxarray: 0.15.1
  ghostscript: 10.02.1
GMT library information:
  binary version: 6.5.0
  cores: 32
  grid layout: rows
  image layout: 
  library path: /home/s3mpc/anaconda3/envs/s3mpc/lib/libgmt.so
  padding: 2
  plugin dir: /home/s3mpc/anaconda3/envs/s3mpc/lib/gmt/plugins
  share dir: /home/s3mpc/anaconda3/envs/s3mpc/share/gmt
  version: 6.5.0

In the old installation, I got the following Warning:

coast [WARNING]: For a UTM or TM projection, your region -799890.9468068393/940434.5640275762/-3464708.338152134/-542960.1486965496 is too large to be in degrees and thus assumed to be in meters
pygmt-session [WARNING]: For a UTM or TM projection, your region -799890.9468068393/940434.5640275762/-3464708.338152134/-542960.1486965496 is too large to be in degrees and thus assumed to be in meters
grdimage [WARNING]: For a UTM or TM projection, your region -799890.9468068393/940434.5640275762/-3464708.338152134/-542960.1486965496 is too large to be in degrees and thus assumed to be in meters
coast [WARNING]: For a UTM or TM projection, your region -799890.9468068393/940434.5640275762/-3464708.338152134/-542960.1486965496 is too large to be in degrees and thus assumed to be in meters
basemap [WARNING]: For a UTM or TM projection, your region -799890.9468068393/940434.5640275762/-3464708.338152134/-542960.1486965496 is too large to be in degrees and thus assumed to be in meters

This makes sense, but has there been a change to GMT/PyGMT? Do I need to adjust the code, or is it an installation issue?

[welcome]

👋 Thanks for opening your first issue here! Please make sure you filled out the template with as much detail as possible. You might also want to take a look at our contributing guidelines and code of conduct.

[seisman]

What happens if you add region=regionr to the following line?

fig.grdimage(grid=back_img, projection=projX, cmap=True, nan_transparent=True, shading=False, transparency=30)

[Sistjerne]

The code never reaches this line. It fails in the line fig.coast(region=regionX, projection=projS, resolution="f", water=True)

[seisman]

What about changing regionX to regionr in fig.coast(region=regionX, projection=projS, resolution="f", water=True)?

[Sistjerne]

Now I get to another issue (I have experienced after updating my conda environment). I get Segmentation fault (core dumped)
when running

fig.grdimage(
    grid="@earth_relief_03m",
    projection=projS,
    region=regionX,
    cmap="grayC",
    shading=True,
)

It helps if I change grid="@earth_relief_03m" to grid="@earth_relief_06m" but I do not experience this in my old conda environment on the same server.

The next issue arises when I want to plot the line

fig.grdimage(grid=back_img, projection=projX, cmap=True, nan_transparent=True, shading=False, transparency=30)

I get the same error: grdimage [ERROR]: Map region exceeds 360 degrees. back_img is the Arctic DEM given in NSIDC Sea Ice Polar Stereographic projection. Do I need to do some reprojection before hand?

[weiji14]

Good to see someone using my old code 😄 You'll need to use regionr for that fig.grdimage(grid="@earth_relief_03m", ... call. Also, can you provide the output of gmt grdinfo <path_to_arcticdem.tif>, or gdalinfo <path_to_arcticdem.tif>? Or the link to the Arctic DEM file you downloaded? This will help us to determine the projection parameters you'll need for plotting back_img. I managed to plot @earth_relief_03m without crashing on my computer, but cannot reproduce the ArcticDEM part. This is as far as I got:

By the way, the IGPP Earth Relief grids (https://www.generic-mapping-tools.org/remote-datasets/earth-relief.html) based on SRTM15+V2.5.5 actually uses ArcticDEM R7 above 60°N (see https://doi.org/10.1029/2019EA000658), but I'm guessing you're trying to plot a higher resolution version?

[Sistjerne]

Thanks for helping out!

Just for your information, I can run the code above in this setup:

PyGMT information:
  version: v0.8.0
System information:
  python: 3.9.15 | packaged by conda-forge | (main, Nov 22 2022, 08:45:29)  [GCC 10.4.0]
  executable: /g5/procdata/skr/anaconda3/envs/s3mpc/bin/python
  machine: Linux-3.10.0-1062.7.1.el7.x86_64-x86_64-with-glibc2.17
Dependency information:
  numpy: 1.23.1
  pandas: 1.4.3
  xarray: 2022.6.0
  netCDF4: 1.6.0
  packaging: 21.3
  geopandas: 0.11.1
  ghostscript: 9.54.0
GMT library information:
  cores: 56
  grid layout: rows
  library path: /nfs/g5/procdata/skr/anaconda3/envs/s3mpc/lib/libgmt.so
  padding: 2
  plugin dir: /nfs/g5/procdata/skr/anaconda3/envs/s3mpc/lib/gmt/plugins
  share dir: /g5/procdata/skr/anaconda3/envs/s3mpc/share/gmt
  version: 6.3.0

and get this output:

The Arctic DEM I am using can be downloaded here: https://www.dropbox.com/scl/fi/1l3il1wvczofm0fgdpw10/arcticdem_mosaic_1km_v3.0_GL.tif?rlkey=gahc2sqauyoyqdmmrqhqiftee&dl=0

I may be able to use the Earh relief grid for this particular plot, but I would like to make this type of plot with different projections for other cases.

[seisman]

I'm going to close this issue since it's more likely a question rather than a bug report. If you'd like to ask for more support about this question, please post it to the forum https://forum.generic-mapping-tools.org/.