[FR] Interpolation does not work for the heatmap in Plots.jl

[songhanzhang]

I wish to draw a heatmap for a stress field, and a linear interpolation is preferred (such like the "shading interp" in Matlab) in order to present a smooth distribution. After searching around, I haven't found a helpful answer.

I noticed an earlier issue, where the interpolation is involved in the heatmap automatically.

But unfortunately, I cannot produce the result by using the same code:

z = reshape(1:25,5,5)
heatmap(z)

Instead, the result is as below:

I have tried

heatmap(z, interpolation = true)
heatmap(z, interpolation = :linear)

and so on, but none of them works.

So, can anyone please tell me how to draw a heatmap with interpolation? Thank you!

Version information:

Julia Version 1.4.2
Commit 44fa15b150* (2020-05-23 18:35 UTC)
Platform Info:
OS: Windows (x86_64-w64-mingw32)
CPU: Intel(R) Xeon(R) E-2236 CPU @ 3.40GHz
WORD_SIZE: 64
LIBM: libopenlibm
LLVM: libLLVM-8.0.1 (ORCJIT, skylake)
Environment:
JULIA_NUM_THREADS = 6

[angeris]

Oooh! I was just running into this issue. I am wondering if there is a "canonical" workaround for it or something of the like :)

In particular, one approach is to use the interpolation from contourf by setting levels to be some large number and setting linewidth=0, but this seems like a not-great hack.

Is there an alternative? If so, that would be awesome!

(This might also be a GR backend thing, though I'm not 100% sure.)

Thanks again for all of your great work, Plots.jl team! Love this package :)

[mdmaas]

Hi, I found this issue. Of course, one possible workaround is to interpolate and plot ;-)

using Interpolations

function Interp2D(data, factor)
    
    IC = CubicSplineInterpolation((axes(data,1), axes(data,2)), data)

    finerx = LinRange(firstindex(data,1), lastindex(data,1), size(data,1) * factor)
    finery = LinRange(firstindex(data,2), lastindex(data,2), size(data,2) * factor)
    nx = length(finerx)
    ny = length(finery)

    data_interp = Array{Float64}(undef,nx,ny)
    for i ∈ 1:nx, j ∈ 1:ny
        data_interp[i,j] = IC(finerx[i],finery[j])
    end

    return finerx, finery, data_interp

end

# Test
data = rand(6,6)

# Simple heatmap
heatmap(data)      

# Heatmap of Interpolated data            
heatmap(Interp2D(data, 8))

I'd like to think about how to incorporate this into Plots.jl, but the package is too large, and I can't find a way to start a PR. Any help would be appreciated!

[songhanzhang]

Hi, I found this issue. Of course, one possible workaround is to interpolate and plot ;-)

using Interpolations

function Interp2D(data, factor)
    
    IC = CubicSplineInterpolation((axes(data,1), axes(data,2)), data)

    finerx = LinRange(firstindex(data,1), lastindex(data,1), size(data,1) * factor)
    finery = LinRange(firstindex(data,2), lastindex(data,2), size(data,2) * factor)
    nx = length(finerx)
    ny = length(finery)

    data_interp = Array{Float64}(undef,nx,ny)
    for i ∈ 1:nx, j ∈ 1:ny
        data_interp[i,j] = IC(finerx[i],finery[j])
    end

    return finerx, finery, data_interp

end

# Test
data = rand(6,6)

# Simple heatmap
heatmap(data)      

# Heatmap of Interpolated data            
heatmap(Interp2D(data, 8))

I'd like to think about how to incorporate this into Plots.jl, but the package is too large, and I can't find a way to start a PR. Any help would be appreciated!

Thank you mdmass! This is quite a feasible approach, but I am not sure if a very large factor (which may lead to a huge computational cost) is always needed ensuring that the interpolated heatmap becomes sufficiently smooth. I will try it later on.

[BeastyBlacksmith]

I'd like to think about how to incorporate this into Plots.jl, but the package is too large, and I can't find a way to start a PR. Any help would be appreciated!

Hey @mdmaas if you are still interested in this send me a pm on either discourse or zulip and we'll find a way to get you started.

[mdmaas]

I'd like to think about how to incorporate this into Plots.jl, but the package is too large, and I can't find a way to start a PR. Any help would be appreciated!

Hey @mdmaas if you are still interested in this send me a pm on either discourse or zulip and we'll find a way to get you started.

Oh sure. I was about to create a recipe for a dummy type to incorporate this and a few other settings I'd like to have in extra kwargs. Do you think a PR is a better idea? Maybe I could do the recipe first and discuss later how to incorporate it to the main package.

[BeastyBlacksmith]

Maybe I could do the recipe first and discuss later how to incorporate it to the main package.

Sure, prototyping in a recipe is a good idea.

[Umut-Can-Physics]

Hi, I found this issue. Of course, one possible workaround is to interpolate and plot ;-)

using Interpolations

function Interp2D(data, factor)
    
    IC = CubicSplineInterpolation((axes(data,1), axes(data,2)), data)

    finerx = LinRange(firstindex(data,1), lastindex(data,1), size(data,1) * factor)
    finery = LinRange(firstindex(data,2), lastindex(data,2), size(data,2) * factor)
    nx = length(finerx)
    ny = length(finery)

    data_interp = Array{Float64}(undef,nx,ny)
    for i ∈ 1:nx, j ∈ 1:ny
        data_interp[i,j] = IC(finerx[i],finery[j])
    end

    return finerx, finery, data_interp

end

# Test
data = rand(6,6)

# Simple heatmap
heatmap(data)      

# Heatmap of Interpolated data            
heatmap(Interp2D(data, 8))

I'd like to think about how to incorporate this into Plots.jl, but the package is too large, and I can't find a way to start a PR. Any help would be appreciated!

Thank you mdmass! This is quite a feasible approach, but I am not sure if a very large factor (which may lead to a huge computational cost) is always needed ensuring that the interpolated heatmap becomes sufficiently smooth. I will try it later on.

Hi! Thank you for your response, I had to same issue. But can you generalize to a rectangular random matrix? For example rand(10,15). When I try to interpolate the rectangular random matrix as a heatmap, Plot looks like bullshit :)

[mdmaas]

Hi @Umut-Can-Physics,

The problem was that finerx and finery were in the wrong order.

So this should work:

using Plots
using Interpolations

function Interp2D(data, factor)
    
    IC = CubicSplineInterpolation((axes(data,1), axes(data,2)), data)

    finerx = LinRange(firstindex(data,1), lastindex(data,1), size(data,1) * factor)
    finery = LinRange(firstindex(data,2), lastindex(data,2), size(data,2) * factor)
    nx = length(finerx)
    ny = length(finery)

    data_interp = Array{Float64}(undef,nx,ny)
    for i ∈ 1:nx, j ∈ 1:ny
        data_interp[i,j] = IC(finerx[i],finery[j])
    end

    return finery, finerx, data_interp

end

data = rand(10,16)

heatmap(Interp2D(data,8))

[mdmaas]

Interpolation with a factor of 16 looks pretty good

[mdmaas]

Hm, it looked better on VSCode before saving the PNG... I guess one could save the file with PDF format...

[mdmaas]

Sure, prototyping in a recipe is a good idea.

I've been playing around with this a little bit...

I guess Plots.jl should support both interpolated heatmaps with given vectors x and y, and fall back to using a range of integers if no arguments are given.

The problem I found with using given x and y, is that Interpolations.jl doesn't seem to support 2D Cubic Splines unless the indices are uniform grids. And people might have non-uniform x and y, and still want to interpolate. Additionally, Interpolations.jl is no longer developed, and a PR for gridded cubic splines got stuck in the middle of nowhere

I guess we could rely on DataInterpolations for 1D interp, and build a simple 2d interpolation routine based on it.

Does this seem right to you @BeastyBlacksmith? I guess if this were to be incoporated into Plots.jl you would need an extra dependecy somewhere. Anyway, let me know what you think.

[mdmaas]

I'm prototyping a recipe here: pcolor.

Given that both Matlab and Python have "pcolor" which is an alias for heatmap, I thought it was a good idea to use that name.

So far only linear interpolation is implemented (which looks really bad on random data), and will implement 2d cubic splines shortly.

This might eventually get merged into Plots.jl if the mantainers judge it is convenient, but for the meantime the idea would be to have this available:

using Plots
using pcolor
pcolor(x,y,data,interp=:cubic)

[BeastyBlacksmith]

You might want to have a look at https://github.com/kbarbary/Dierckx.jl

[mdmaas]

You might want to have a look at https://github.com/kbarbary/Dierckx.jl

There, it's done in the pcolor repo

One issue I couldn't figure out is why pcolor(data) doesn't plot anything.

Another nice thing to have would be access to the figure resolution to automatically set the oversampling_factor so we get an interpolated image of the size of the resolution... any clues on how to access fig resolution from within a recipe? Is it possible?

[BeastyBlacksmith]

Isn't that given by the sizes of plotattributes[:z]? Or the x and y component.

[mdmaas]

Oh thanks, I finally figured I could get the figure size with:

plotattributes[:plot_object].attr[:size]

So I think it's best to just interpolate to a number of points equal to the resolution size (600x400 by default) and remove the oversampling_factor option.

[mdmaas]

This now works to full resolution even when you change it with:

pcolor(data, size=(1200,1200), interpolate=:true)

which is nice.

The only thing I haven't figured out how to fix is that pcolor(data) doesn't plot anything, and I indended to forward that to heatmap(data).

[BeastyBlacksmith]

Oh thanks, I finally figured I could get the figure size with:

plotattributes[:plot_object].attr[:size]

So I think it's best to just interpolate to a number of points equal to the resolution size (600x400 by default) and remove the oversampling_factor option.

Note that that is the size of the whole figure, not just the plotarea, including axis, labels, more subplots if this is a layout, etc.

[mdmaas]

Hi! So I went back to this issue.

I got around the first problem I had, with pcolor(data) not doing anything. It turned out that I was defining a @userplot using lower caps, which resulted in the dummy type collisioning with the new plot function.

So, it's working now. I also also had to rename the package to PColorPlot.jl, fortunately it wasn't registered yet. The repo is now here

Note that that is the size of the whole figure, not just the plotarea, including axis, labels, more subplots if this is a layout, etc.

This one I couldn't fix it. I believe that when the plot recipe gets called it doesn't have access to the final layout information. Given that the interpolation library is fast enough, I don't believe this to be much of an issue.

Anyway, do you think this could be merged into Plots.jl somehow (maybe adding the extra flag to heatmap directly?), or prefer this as a separate package? Now we have package extensions, so it could be an option, too.

[BeastyBlacksmith]

We can definetely integrate this into plots. I imagine the interface being a new keyword interpolate/interpolation that accepts a tuple specifying the resolution per dimension, e.g. interpolate = (20, 40) being a 2D interpolation with 20 points in x direction and 40 points in y-direction. This should happen after the recipe pipeline I'd put that here:

Plots.jl/src/pipeline.jl

Line 375 in ec4b7a4

end

[mdmaas]

I'll get a PR started, then, as soon as I have some time for this.