Cartopy projections?

[nvaytet]

@Adnan-Ali-Ahmad
In addition, we could also add cartopy projections of data. I realize this sounds a bit like scope creep, and maintaining the code base is daunting task, but I'm always happy to help

@nvaytet
Do you have an example of an application for this? I'm trying to understand what you need it for.

@Adnan-Ali-Ahmad
Fig.6 of your 2018 paper is a great example of this. Although at the time you may have used basemap, it is now deprecated in favor of cartopy.

[Adnan-Ali-Ahmad]

Hello, I've started working on implementing this in osyris and I have a question.

I first select the cells that are in contact with a sphere of radius "r". After this, we simply have to project the cell values unto a cartopy map, but I'm currently torn between two methods; the first consists of interpolating the projected data unto a 2D grid of latitude-longitude and plot the data with pcolormesh (or countourf), which allows the user to specify his desired resolution. However, there also exists tricontourf, which doesn't need any interpolation and allows directly plotting the projection. The latter benefits from simplicity (and looks good), but it might be computationally expensive at higher resolution simulations (but then again, so will the interpolation). Here's what this looks like when plotting a density sphere 200 AU away from the origin (of max density) with the sample data:

The first three plots are done using pcolormesh rendering of the interpolated data, whereas tricontourf is a direct render of the cell values. For reference, here's what the 3D sphere looks like (the scatter points are cell centers):

[nvaytet]

I suspect tricontourf is using something like https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.griddata.html
I.e. it first computes an unstructured mesh of triangles (probably Delaunay) and then does its interpolations.
I was using this in earlier versions of osiris for making the maps, but it got very slow at high res.

So I would be in favour of doing our own interpolations, using numba as I did for the maps.

Finally, we should probably think about how this task is related to #52 , because maybe it could involve converting to spherical coords, making a normal map in theta and phi, and then just applying the projection?

[Adnan-Ali-Ahmad]

I suppose so too. Yes, the latitude-longitudes are computed using a cartesian to spherical transformation. I've also coded a simple 3D cartesian rotation.

Also, how would we make these plots? osyris.cartopy() ? osyris.plot_surface() ?

[nvaytet]

I was just thinking:

1. convert to spherical coordinates using some utility (which we yet have to define)
2. make a (theta, phi) map at r = 200 AU using map(..., projection='hammer') or something like that.

If you don't specify the projection, the map will just be rectangular, like a normal map, but still represent theta and phi on the axes.
You could also envisage making slices in (r, theta) at a constant phi, and use a projection that would draw a wedge, that would look something like:

In fact, we might even be able to use the thick maps to make azimuthal averages around the phi dimension.

We have to think of the (r, theta, phi) dimensions as logical dimensions, just like (x,y,z), and just let matplotlib handle the projections.

I think it is very much worth thinking for a while about the best way to do this, rather than jump straight in and start coding because this initial effort will most likely be for nothing, and wasting developer time is not a good thing.