Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Basemap fails with Mollweide projection across prime meridian or negative longitude defining both central meridian and range #8403

Closed
thomasAmorrow opened this issue Mar 14, 2024 · 5 comments · Fixed by #8404
Labels
bug Something isn't working

Comments

@thomasAmorrow
Copy link

Description of the problem

Basemap fails when creating a Molleweide projection that crosses 0 longitude and when creating a Molleweide projection with a negative longitude central meridian and negative longitude for E/W range arguments.

Interestingly, it succeeds when defining a negative longitude central meridian and a positive E/W range.

This might be two separate issues, one associated with crossing 0 longitude and the other somehow associated with using negative longitude for range E/W, but I can't disentangle them. 6 examples included, the last 3 fail with the same error message.

The actual need I have is to create a map in a Molleweide projection defined by a central meridian west of 0 longitude and with a range that spans across 0 longitude.

Minimal Complete Verifiable Example

# Load the PyGMT package.
import pygmt

fig1 = pygmt.Figure()
fig2 = pygmt.Figure()
fig3 = pygmt.Figure()
fig4 = pygmt.Figure()
fig5 = pygmt.Figure()
fig6 = pygmt.Figure()

# Basemap succeeds with negative central meridian and positive region EW
fig1.basemap(frame=["afg", "+t1"], projection="w-30/0.1", region="260/350/0/90")
fig1.show()

# Basemap succeeds with positive central meridian and positive region EW
fig2.basemap(frame=["afg", "+t2"], projection="w330/0.1", region="260/350/0/90")
fig2.show()

# Basemap succeeds with positive central meridian and negative region EW
fig3.basemap(frame=["afg", "+t3"], projection="w30/0.1", region="-100/-10/0/90")
fig3.show()

# Basemap fails with negative central meridian and negative region EW
fig4.basemap(frame=["afg", "+t4"], projection="w-30/0.1", region="-100/-10/0/90")
fig4.show()

# Basemap fails with negative central meridian and negative region EW across 0
fig5.basemap(frame=["afg", "+t5"], projection="w-30/0.1", region="-100/20/0/90")
fig5.show()

# Basemap fails with positive central meridian and negative region EW across 0
fig6.basemap(frame=["afg", "+t6"], projection="w330/0.1", region="-100/20/0/90")
fig6.show()

Full error message

basemap [ERROR]: Eastern boundary is > 180 degrees from specified central meridian and thus your region is invalid
basemap [ERROR]: General map projection error

---------------------------------------------------------------------------

GMTCLibError                              Traceback (most recent call last)

<ipython-input-30-e9ccff593116> in <cell line: 27>()
     25 
     26 # Basemap fails with negative central meridian and negative region EW
---> 27 fig4.basemap(frame=["afg", "+t4"], projection="w-30/0.1", region="-100/-10/0/90")
     28 
     29 fig4.show()

3 frames

/usr/local/lib/python3.10/site-packages/pygmt/helpers/decorators.py in new_module(*args, **kwargs)
    601                 warnings.warn(msg, category=SyntaxWarning, stacklevel=2)
    602 
--> 603             return module_func(*args, **kwargs)
    604 
    605         new_module.aliases = aliases

/usr/local/lib/python3.10/site-packages/pygmt/helpers/decorators.py in new_module(*args, **kwargs)
    774 
    775             # Execute the original function and return its output
--> 776             return module_func(*bound.args, **bound.kwargs)
    777 
    778         return new_module

/usr/local/lib/python3.10/site-packages/pygmt/src/basemap.py in basemap(self, **kwargs)
     85     kwargs = self._preprocess(**kwargs)
     86     with Session() as lib:
---> 87         lib.call_module(module="basemap", args=build_arg_string(kwargs))

/usr/local/lib/python3.10/site-packages/pygmt/clib/session.py in call_module(self, module, args)
    622         )
    623         if status != 0:
--> 624             raise GMTCLibError(
    625                 f"Module '{module}' failed with status code {status}:\n{self._error_message}"
    626             )

GMTCLibError: Module 'basemap' failed with status code 74:
basemap [ERROR]: Eastern boundary is > 180 degrees from specified central meridian and thus your region is invalid
basemap [ERROR]: General map projection error

System information

PyGMT information:
  version: v0.11.0
System information:
  python: 3.10.12 (main, Nov 20 2023, 15:14:05) [GCC 11.4.0]
  executable: /usr/bin/python3.real
  machine: Linux-6.1.58+-x86_64-with-glibc2.35
Dependency information:
  numpy: 1.25.2
  pandas: 2.2.1
  xarray: 2024.2.0
  netCDF4: 1.6.5
  packaging: 24.0
  contextily: None
  geopandas: 0.13.2
  ipython: None
  rioxarray: None
  ghostscript: 10.02.1
GMT library information:
  binary version: 6.5.0
  cores: 2
  grid layout: rows
  image layout: 
  library path: /usr/local/lib/libgmt.so
  padding: 2
  plugin dir: /usr/local/lib/gmt/plugins
  share dir: /usr/local/share/gmt
  version: 6.5.0
@thomasAmorrow thomasAmorrow added the bug Something isn't working label Mar 14, 2024
Copy link

welcome bot commented Mar 14, 2024

👋 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
Copy link
Member

seisman commented Mar 15, 2024

The equivalent GMT commands also fail. So I'm transfering the issue to the upstream GMT repository. The equivalent GMT commands are:

gmt basemap -Jw-30/0.1 -R-100/-10/0/90 -Bafg -pdf fig4
gmt basemap -Jw-30/0.1 -R-100/20/0/90 -Bafg -pdf fig5
gmt basemap -Jw330/0.1 -R-100/20/0/90 -Bafg -pdf fig6

@seisman seisman transferred this issue from GenericMappingTools/pygmt Mar 15, 2024
Copy link

welcome bot commented Mar 15, 2024

👋 Thanks for opening your first issue here! Please make sure you filled out the template with as much detail as possible. We appreciate that you took the time to contribute!

Please make sure you read our Contributing Guide and abide by our Code of Conduct.

@yvonnefroehlich
Copy link
Member

yvonnefroehlich commented Mar 15, 2024

The actual need I have is to create a map in a Molleweide projection defined by a central meridian west of 0 longitude and with a range that spans across 0 longitude.

Thanks @thomasAmorrow for your detailed report! I can reproduce this issue with GMT 6.5.0, but it worked with GMT 6.4.0. The analogous issue occurs with the Hammer (H) projection. The projections N, Kf, Ks, I, and R seem to work.

Depending on how fast you need this map and if you do not need any features or bug fixes of GMT 6.5.0, you may want to consider setting up a conda environment with PyGMT and GMT 6.4.0 (, even it is generally recommended to use the latest versions).

For me, the following code works with GMT 6.4.0, but fails with GMT 6.5.0 (the same holds for the equivalent GMT commands by @seisman in comment #8403 (comment)):

import pygmt

dict_J_R = {
    "0/180/0/90": "W30/10c",  # Works
    "-180/0/0/90": "W-30/10c",  # Fails with GMT 6.5.0
    "-20/180/0/90": "W30/10c",  # Works
    "-180/20/0/90": "W-30/10c",  # Fails with GMT 6.5.0
}

fig = pygmt.Figure()

for panel in dict_J_R.keys():
    fig.basemap(
        region=panel,
        projection=dict_J_R[panel],
        frame=["afg", f"+tJ={dict_J_R[panel]}, R={panel}"],
    )
    fig.shift_origin(xshift="+w+1c")

fig.show()
# fig.savefig(fname="mollweide_subregion.png")

Output figure:
mollweide_subregion

@joa-quim
Copy link
Member

Hmm, we are having an issue with central longitudes in [0 360] and regions in [-180 180]. For some reason the central longitude is converted to [0 360] and a test saying

/* For regional (<360) areas we cannot have clon > 180 away from either boundary */

fails. For the moment the solution is to give the -R in [0 360]. e.g., this works

gmt basemap -Jw-30/0.1 -R260/350/0/90 -Bafg -png lixo

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working
Projects
None yet
4 participants