Cartographic projections using RecipesBase

Specific Domains Visualization
1

Hi to everybody,

In the past years, I have developed a library that implements Healpix projection functions: GitHub - ziotom78/Healpix.jl: Healpix library written in Julia. Unlike Michael Eastwood’s Libhealpix.jl, this library is 100% pure Julia, because I needed to make it run under Windows as well (Libhealpix is a wrapper to the C++ Healpix library, which does not support Windows yet).

So far, I have managed to produce cartographic projections using Cairo to create RGB bitmaps, which look nice in Jupyter but are not shown in the REPL. As the dependency on Cairo, Images, and ImageMagick is quite heavy and does not work well for those applications that need to run on clusters, I would like to move to RecipesBase.

Using Cairo has allowed me to render images similar to the one here:

800×600 255 KB

This is a Mollweide projection of the celestial sphere, where pixels along the Galactic plane have been masked and are shown in gray.

My question is: how am I supposed to produce this kind of images using RecipesBase? My understanding is that I should write a @recipe returning a matrix containing the 2D projection of the map, using missing for the pixels outside the ellipse. This matrix should contain floating-point values, which would be mapped to colors by a user-specified color scale, like the ones defined in ColorSchemes.jl. But how could I specify that «masked» pixels must be painted using a dedicated color?

The Python bindings to Healpix use Matplotlib’s features to provide custom axis that take care of the projection (healpy/projaxes.py#L719); from my understanding, in this case pixels get a gray color if they map to some special values (like numpy.NaN), and pixels outside the ellipse are never painted at all.

The best solution would allow to add layers of vector elements: the most common case is a graticule showing parallels and meridians, and points and lines indicating sky sources and regions.
(P.S. I have made the case of Mollweide projections, but the library should provide other projections as well, of corse. So I am looking for some general way to implement this, which does not assume that the pixels will be always plotted within an ellipse.)

Thanks a lot!
Maurizio.

1 Like
2

It’d be a lot easier to answer with some specific code examples. But yes, briefly, you’l make a heatmap of a 2-d Matrix{Float64} with NaNs as transparent fields. Masking and adding graticules etc is easiest done by adding additional @series where primary = false (so they don’t show up in your color legend).

2 Likes
3

Many thanks @mkborregaard, the reason why I did not submit any code (as I usually do) is that I’m actually in the planning stage and have yet to write anything. Your idea to use @series is perfect, many thanks!

1 Like
4

Great - you’re welcome to post code for comments once you get to that stage.

5

After a few trials, I have been able to produce something really close to what I want. But I am struggling to understand how to properly implement masks.

This short script reproduces the problem:

using RecipesBase

# This type represents the map
struct T end

# Very silly, "cross-like" projection. The function
# returns the bitmap to display and a 2D bitmap
# marking the position of masked pixels (to be
# shown using a uniform gray shade)
function projection()
    # NaN marks the pixels outside the map area
    img  = [[NaN  1. NaN];
            [2.   2.  2.];
            [NaN  3. NaN]]

    # NaN means "unmasked" (default), 0 means "masked"
    # Only the pixel at the center of the cross should be masked
    mask = [[NaN NaN NaN];
            [NaN  0. NaN];
            [NaN NaN NaN]]

    (img, mask)
end

@recipe function plot(::T)
    seriestype := :heatmap
    aspect_ratio := 1
    colorbar := :bottom
    framestyle := :none
    
    img, mask = projection()

    @series begin
        seriestype := :heatmap
        primary := false
        color_palette := [:gray] # Masked pixels should be gray
        mask
    end
    
    img
end

using Plots
gr()

plot(T())

The result is the following:

cross

The shape of the cross and the shades look as expected, but the center of the cross has not been painted with gray. I don’t understand if the mask is being ignored, or if it has been drawn below the map itself.

May you provide me some tip to debug this issue, please?

6

It’s plotting below, yes, but unfortunately it seems Plots has issues with overplotting different heatmaps with different color gradients :frowning:

You could implement your own seriestype that draws the mask as a polygon…

7

unfortunately it seems Plots has issues with overplotting different heatmaps

That’s a pity to hear. I have looked for an issue of this kind in Issues · JuliaPlots/Plots.jl · GitHub, but I found nothing. Should I report this somewhere?

In the meantime, I will try to follow the polygon route…

8

Fantastic, I have been able to put together a solution that works pretty well. Have a look at this:

masked_map
masked_map936×457 12.7 KB

If you are curious, the implementation of the recipe is in file src/projections.jl, and a quite sketchy documentation is here.

Using polygons was easier than I thought. Many thanks, @mkborregaard, I would have never been able to do this without your suggestions!

7 Likes
9

Very nice!

1 Like