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

Fix plotting issues by mapping from lon/lat manually #14

Merged
merged 1 commit into from
Feb 6, 2016

Conversation

xylar
Copy link
Collaborator

@xylar xylar commented Feb 6, 2016

A number of plotting problems are fixed by first mapping coordinates
from lon/lat to projection coordinates and then plotting, rather than
using the latlon argument to the basemap plot routine.

Added some output so users are aware of what is being plotted and
get a sense of how long it will take.

A number of plotting problems are fixed by first mapping coordinates
from lon/lat to projection coordinates and then plotting, rather than
using the latlon argument to the basemap plot routine.

Added some output so users are aware of what is being plotted and
get a sense of how long it will take.
@xylar
Copy link
Collaborator Author

xylar commented Feb 6, 2016

@pwolfram, let me know if these fixes also work for you. I only tested the default 4 and the new 5 projections (so there are still 8 I didn't test).

@pwolfram
Copy link
Contributor

pwolfram commented Feb 6, 2016

@xylar, thanks for figuring out the solution. I was surprised this was a simple solution and looked into it further.

This appears to be an unanticipated side-effect of https://github.com/matplotlib/basemap/blob/master/lib/mpl_toolkits/basemap/__init__.py#L510. The coordinates are being shifted for the option latlon=True which I think is probably causing the problem because they also are automatically converting to map coordinates in https://github.com/matplotlib/basemap/blob/master/lib/mpl_toolkits/basemap/__init__.py#L543, as your fix does singularly. I don't think this necessarily is a bug on basemap's side but if you think it is I will submit a bug report. Regardless, thanks for figuring this out.

Using your target regions listed above, here is the test for ['atlantic','pacific', 'europe', 'northamerica']

screenshot 2016-02-06 15 34 54

The plots are not quite working correctly if the points are outside the map, but do work if they are contained within the map (e.g., https://github.com/MPAS-Dev/geometric_features/blob/master/landice/region/southGreenland1/region.geojson).

The test for ['ortho','robin', 'robin2', 'hammer'] yields

screenshot 2016-02-06 15 40 31

The results plotted appear to be quite reasonable, with the exception of the Rhone_1 region on the robin2 and hammer plots. We wouldn't expect regions to plot correctly for all the projections because some points will be outside the plotted region. Because of this I think it is fair to merge this PR because your code adds value and fixes an outstanding issue. If these other map types need fixed some time in the future I would say this should be a new PR.

pwolfram added a commit that referenced this pull request Feb 6, 2016
@pwolfram pwolfram merged commit 35af93d into MPAS-Dev:master Feb 6, 2016
@pwolfram
Copy link
Contributor

pwolfram commented Feb 6, 2016

@xylar, as you may have surmised I cannot delete your branch either... thanks for the really clean solution here.

@xylar xylar deleted the fix_plotting_at_poles branch February 7, 2016 09:40
@xylar
Copy link
Collaborator Author

xylar commented Feb 7, 2016

@pwolfram, thanks for reviewing this. It seems there's still some work to do.

By the way, do you know why are the 'europe' and 'northamerica' plots named in that way? It seems like the 'europe' plot shows the Americas and the 'northamerica' plot shows Asia and the Indian Ocean.

@xylar
Copy link
Collaborator Author

xylar commented Feb 7, 2016

@pwolfram, by the way, I do think there's a bug in shiftdata. It tries to find the second biggest jump between longitudes and use that as a threshold for determining where the jump in longitude happens:
https://github.com/matplotlib/basemap/blob/master/lib/mpl_toolkits/basemap/__init__.py#L4837
The problem is that our polygons sometimes have several points that pass the inequality on this line, resulting in a list of indices rather than a single index.

When I get a chance, I'll put together a simpler example that reproduces the problem and file a bug report.

@xylar
Copy link
Collaborator Author

xylar commented Feb 7, 2016

@pwolfram, this seems to be a known issue that you have also commented on: matplotlib/basemap#126

It seems like latlon=True is only useful for longitude data that is monotonic (except that it's periodic), which is obviously not the case for our polygons.

@pelson
Copy link

pelson commented Feb 8, 2016

I don't know how important it is to you to be able to plot geometries on arbitrary maps, but you may be interested in cartopy, whose bread-and-butter is geometry transformations of this kind (e.g. http://scitools.org.uk/cartopy/docs/latest/examples/rotated_pole.html). For sure Basemap is more mature as a package, but Cartopy takes the approach you've taken here @xylar, in that it is cartopy's responsibility to transform the coordinates, not the user. This means it has all of the context it needs to understand datelines, poles etc.

Anyway, nice plots! 👍

@xylar
Copy link
Collaborator Author

xylar commented Feb 8, 2016

@pelson, Thank for the tip. I had heard of cartopy but haven't used it before. It sounds like that's what we should have been using from the start. I'll look into it further.

@pwolfram
Copy link
Contributor

pwolfram commented Feb 9, 2016

Our colleagues @milenaveneziani and @akturner have been using it for their global ocean / sea-ice analysis metrics and I think they have found it useful.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants