Wrong point location in inset

[eyalshimony]

Description of the problem

When plotting a point in an inset, it is drawn in the wrong geographic location.

Full code that generated the error

import pygmt

region = [-117.628, -117.478, 35.605, 35.801]
fig = pygmt.Figure()
pygmt.makecpt(
    cmap='etopo1',
    series='0/2800/200',
    continuous=True
)
topo_data = "@earth_relief_01s"
fig.grdimage(
    grid=topo_data,
    region=region,
    projection='M4i',
    shading=True,
    frame=True
)

with fig.inset(position="jTR+w3c", box="+pblack"):
    # Use a plotting function to create a figure inside the inset
    fig.coast(
        region=[-125, -110, 30, 45],
        projection="M3c",
        land="gray",
        borders=[1, 2],
        shorelines="1/thin",
        water="white"
    )
    fig.plot(x=-117.55, y=35.7, size=[0.2], style="c", color="black")
fig.show()

System information

Please paste the output of python -c "import pygmt; pygmt.show_versions()":

PyGMT information:
version: v0.6.0
System information:
python: 3.10.0 (tags/v3.10.0:b494f59, Oct 4 2021, 19:00:18) [MSC v.1929 64 bit (AMD64)]
executable: C:\Users\eyals\AppData\Local\Programs\Python\Python310\python.exe
machine: Windows-10-10.0.19044-SP0
Dependency information:
numpy: 1.21.2
pandas: 1.3.3
xarray: 0.19.0
netCDF4: 1.5.7
packaging: 21.0
ghostscript: 9.54.0
gmt: 6.3.0
GMT library information:
binary dir: C:/Users/eyals/AppData/Local/Programs/Python/Python310
cores: 4
grid layout: rows
library path: c:/programs/gmt6/bin/gmt_w64.dll
padding: 2
plugin dir: c:/programs/gmt6/bin/gmt_plugins
share dir: c:/programs/gmt6/share
version: 6.3.0

[weiji14]

Hi @eyalshimony, thanks for trying PyGMT! So it's a bit unintuitive, but it's best not to set the projection width in the inset plot. Try this instead:

with fig.inset(position="jTR+w3c", box="+pblack"):
    # Use a plotting function to create a figure inside the inset
    fig.coast(
        region=[-125, -110, 30, 45],
        projection="M?c",
        land="gray",
        borders=[1, 2],
        shorelines="1/thin",
        water="white"
    )
    fig.plot(x=-117.55, y=35.7, size=[0.2], style="c", frame="afg", color="black")

Specifically, use projection="M?c" instead of projection="M3c", which means the size is scaled automatically to the inset box. This should produce:

You may need to adjust the width in fig.inset(position="jTR+w3c", ...) to make the box be the correct size.

[eyalshimony]

Thanks for the quick response!

It may be helpful for other users if you'll correct the example in https://www.pygmt.org/v0.3.0/tutorials/inset.html accordingly.

[weiji14]

Thanks for the quick response!

It may be helpful for other users if you'll correct the example in https://www.pygmt.org/v0.3.0/tutorials/inset.html accordingly.

Yes that would be a good idea. Would you like to try making a quick change at https://github.com/GenericMappingTools/pygmt/edit/main/examples/tutorials/advanced/insets.py?

[eyalshimony]

I submitted the change

[eyalshimony]

I am sorry, I checked it and it actually ruins some stuff in the example. Most importantly, the highlighting stops working. I do not know enough about this in order to fix it, so I think it is best be left to the experts.

[michaelgrund]

I submitted the change

Where did you submit it? Just take a look at our contributing guide on how to create a PR for fixing issues 😉.

[weiji14]

I see you've made the change at https://github.com/eyalshimony/pygmt/tree/patch-1. If you click on the 'Contribute' button and then Open Pull Request, we'll be able to see it.

[eyalshimony]

Yes, I understood it afterwards but as I said above changing the projection along your recommendations ruins the highlighting. I have no idea why, so I think it would be better if you'd fix it.

[weiji14]

Hmm yeah, I just tried using M?c for that example but don't understand why the red highlight is incorrect either 😅 You're still welcome to open the Pull Request if you'd like, that way it will be faster for someone to check it and give suggestions on what to do. But I understand if you're short on time.

[seisman]

I'm reopening the issue to track the potential issue in https://github.com/GenericMappingTools/pygmt/blob/main/examples/tutorials/advanced/insets.py

[seisman]

Hmm yeah, I just tried using M?c for that example but don't understand why the red highlight is incorrect either sweat_smile You're still welcome to open the Pull Request if you'd like, that way it will be faster for someone to check it and give suggestions on what to do. But I understand if you're short on time.

We have to use M?, not M?c, because ? will be repalced by an automatic width like 1.5i, so M?c becomes M1.5ic, which is invalid.

[eyalshimony]

With M? the highlighting exists, but it is in the wrong location

[seisman]

With M? the highlighting exists, but it is in the wrong location

Yes, I can confirm the issue. After some tries, I think the only one that works is projection="M?+du".

[seisman]

Here is a GMT bash script to reproduce the issue:

gmt begin map
gmt coast -R-74/-69.5/41/43 -JM15c -Baf -Wthin -Clightblue -Glightyellow
gmt inset begin -DjBL+w3c+o0.5c/0.2c -F+pblack+gwhite
    gmt coast -R-80/-65/35/50 -JM?+du -Ggray -Cwhite -EUS.MA+gred -W1/thin -N1 -N2
gmt inset end
gmt end show

For the coast call in the inset, users have to use -JM?+du so that the highlight is placed at the right position. Other forms like -JM? and -JM?+dw don't work.

@PaulWessel The use of -JM?+du is not intuitive. Is there anything GMT can do to improve it?

[seisman]

Ping @PaulWessel on #1930 (comment)

[PaulWessel]

Yep, saw this issue but up to my neck in issues (not all GMT), so I will get to it when I can.

[PaulWessel]

The key thing that happens is the unfortunate selection of -R and inset size (square). That map is not square, so GMT finds that the hight with -JM3c is 4.08581069804c which exceeds the square so it applies the scale 3/4.08581069804 so that the max dimension is 3c. That means the width is actually 2.20274522369c which is what is used when -JM? is converted to -JM2.20274522369c. However, coast passes the buck to the function gmt_DCW_operation that plots the -E stuff and clearly there is some issues with what projection parameters are used. I will have to debug to see why, and once I understand what happens it should be doable to fix it.

I am not sure what to do for the unfortunate selection of -R and a square inset. It is not easy for the user to know what is the best -R that gives a square map for some projection. Perhaps we need to come up with a better implementation here. Maybe the inset command should be +w3c/? to indicate we want GMT to compute the height of the inset since it is unknown. Or similarly +w?/3c to fix the height and let the width be determined. To me, that seems reasonable: The user selects the region and we compute the missing dimension.

[PaulWessel]

And gmt_DCW_operation ends up making a virtual dataset and passing it via a call to psxy where the command string contains -J. So I think it picks up the -Jm3c since the change to -Jm2.2c happens internally and the history is not updated. I think if I update the history it may be work if it happens before the cal.

[seisman]

@eyalshimony Here is a better solution for you case. So instead of specifying +w3c in position, you can set region and projection parameter so that GMT can automatically determine both width and height of the inset for you.

with fig.inset(position="jTR", box="+pblack", region=[-125, -110, 30, 45], projection="M3c"):
    fig.coast(
        land="gray",
        borders=[1, 2],
        shorelines="1/thin",
        water="white"
    )
    fig.plot(x=-117.55, y=35.7, size=[0.2], style="c", color="black")

Unfortunately, inset's region and projection parameters are not aliased in PyGMT v0.6.1. As a workaround, you can use the following code instead:

with fig.inset(position="jTR", box="+pblack", R="-125/-110/30/45", J="M3c"):