Useful for visualizing on a sphere? #72

dlfivefifty · 2018-03-02T20:08:06Z

I'm looking for a "good" method to do visualizations on the sphere, including heat maps and vector fields.

Is this possible in Makie?

SimonDanisch · 2018-03-02T20:45:22Z

Can you give an example? how i imagine it, mapping something to the surface of a sphere, it should be possible ;)

dlfivefifty · 2018-03-02T20:59:44Z

For the density plot, @MikaelSlevinsky has some nice examples he did using PlotlyJS (https://arxiv.org/pdf/1801.04902.pdf)

For vector fields on a sphere, the simplest is just to plot a lot of arrows, a bit nicer would be a stream plot but where the stream lines are on the sphere, with colour used to indicate magnitude.

(In my case the vector fields are restricted to the tangent space.)

SimonDanisch · 2018-03-05T11:42:48Z

Is the code somewhere, so I can get an idea of the interface?

dlfivefifty · 2018-03-05T11:45:16Z

Yes: https://github.com/JuliaApproximation/SpectralTimeStepping.jl/blob/master/NonlocalSphere/plot_Allen_Cahn.jl
The key line is

  s = surface(x=x, y=y, z=z, colorscale = "Viridis", surfacecolor = F, cmin = -1.0, cmax = 1.0, showscale = false)

where F is a Matrix of height values of the function.

dlfivefifty · 2018-03-05T11:49:21Z

For the arrow plot case, the interface should match quiver: something like

quiver(pts, quiver=df)

where pts is a Vector{<:Point} or Vector{<:Tuple} and df is also a Vector{<:Point} or Vector{<:Tuple}.

SimonDanisch · 2018-03-05T11:52:05Z

Do you expect the arrows do bend according to the surface of the sphere?

SimonDanisch · 2018-03-05T11:53:12Z

and will the points / vectors be 2d, and you will indicate with a keyword argument, that they should be on a sphere?

dlfivefifty · 2018-03-05T11:53:30Z

Do you expect the arrows do bend according to the surface of the sphere?

I "expect" them to not because that's hard to implement, but I wish they would.

dlfivefifty · 2018-03-05T11:53:45Z

will the points / vectors be 2d

No, 3D. (In principle I'd expect the interface to work for other surfaces too.)

SimonDanisch · 2018-03-05T12:34:00Z

Surface works!

I will make a pr to make it a bit more user friendly, since right now, you need to map F to an array of colors yourself... But that actually fits well with implementing image(x::Matrix{Float64}, colormap = ..., color_norm = (-1, 1)) to also work as intended ;)

I just recently added a prototype for quiver, I'll take a look how far it is ;)

SimonDanisch · 2018-03-05T12:40:00Z

Do you also have a piece of code that gives me pts and df?

dlfivefifty · 2018-03-05T14:12:36Z

For quiver, probably something like:

f   = (x,y,z) -> x*exp(cos(y)*z)
∇f  = (x,y,z) -> SVector(exp(cos(y)*z), -sin(y)*z*x*exp(cos(y)*z), x*cos(y)*exp(cos(y)*z))
∇ˢf = (x,y,z) -> ∇f(x,y,z) - SVector(x,y,z)*dot(SVector(x,y,z), ∇f(x,y,z))

θ = [0;(0.5:n-0.5)/n;1]
φ = [(0:2n-2)*2/(2n-1);2]
x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
z = [cospi(θ) for θ in θ, φ in φ]

pts = vec(SVector.(x, y, z))
∇ˢF = vec(∇ˢf.(x, y, z))

quiver(pts, df=∇ˢF)

SimonDanisch · 2018-03-05T15:48:17Z

Something like this?

dlfivefifty · 2018-03-05T15:54:37Z

Perfect! (I assume all the arrows are tangential to the sphere... hard to tell without being able to rotate and zoom.)

If the arrows could wrap around the sphere it would be even nicer... but I guess that will take some work (and may be too slow in practice)...

SimonDanisch · 2018-03-05T16:08:39Z

Slowness shouldn't be a problem, but someone will need to put some time into emiting the correct line paths... I'm guessing one would need to project the 3D vectors on the surface of the sphere?

They're not really tangential, since pts .+ ∇ˢF doesn't result on the points being on the sphere. Maybe I'm interpreting ∇ˢF wrong? I'm using pts[i] => pts[i] + ∇ˢF[i] to generate line segments

dlfivefifty · 2018-03-05T16:32:50Z

It should be true that dot(∇ˢF[i], pts[i]) ≈ 0 for all i by construction (since norm(SVector(x,y,z)) == 1). In other words, ∇ˢF[i] is tangential to the normal at pts[i].

This means the arrow connecting pts[i] and ∇ˢF[i] should line in the tangent plane (that is if you draw a plane tangential to the sphere at pts[i])

Note I'm assuming this is a unit sphere, so that (x,y,z) is both the point in space and the normal vector to the surface at that point.

dlfivefifty · 2018-03-05T16:34:47Z

They're not really tangential, since pts .+ ∇ˢF doesn't result on the points being on the sphere

To clarify further: pts .+ ∇ˢF should not be on the sphere.

SimonDanisch · 2018-03-05T16:45:17Z

ok!
Well the code is easy enough to just paste it here, if you want to play around with it:

using GeometryTypes, Colors, Makie

function arrows(
        parent, points::AbstractVector{Pair{Point3f0, Point3f0}};
        arrowhead = Pyramid(Point3f0(0, 0, -0.5), 1f0, 1f0), arrowtail = nothing, arrowsize = 0.3,
        linecolor = :black, arrowcolor = linecolor, linewidth = 1,
        linestyle = nothing, scale = Vec3f0(1)
    )
    linesegment(parent, points, color = linecolor, linewidth = linewidth, linestyle = linestyle, scale = scale)
    rotations = map(points) do p
        p1, p2 = p
        dir = p2 .- p1
        GLVisualize.rotation_between(Vec3f0(0, 0, 1), Vec3f0(dir))
    end
    meshscatter(
        last.(points), marker = arrowhead, markersize = arrowsize, color = arrowcolor,
        rotations = rotations, scale = scale
    )
end
n = 20
f   = (x,y,z) -> x*exp(cos(y)*z)
∇f  = (x,y,z) -> Point3f0(exp(cos(y)*z), -sin(y)*z*x*exp(cos(y)*z), x*cos(y)*exp(cos(y)*z))
∇ˢf = (x,y,z) -> ∇f(x,y,z) - Point3f0(x,y,z)*dot(Point3f0(x,y,z), ∇f(x,y,z))

θ = [0;(0.5:n-0.5)/n;1]
φ = [(0:2n-2)*2/(2n-1);2]
x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
z = [cospi(θ) for θ in θ, φ in φ]

pts = vec(Point3f0.(x, y, z))
∇ˢF = vec(∇ˢf.(x, y, z))

arr = map(pts, ∇ˢF) do p, dir
    p => p .+ (dir .* 0.2f0)
end
scene = Scene();
surface(scene , x, y, z)
arrows(scene, arr, arrowsize = 0.03, linecolor = :gray, linewidth = 3)

dlfivefifty · 2018-03-05T16:47:13Z

Cool! That's good enough for me.

dlfivefifty · 2018-03-05T16:50:06Z

PS Is it easy to create animated GIFs like on the README?

I'll send you some nice movies of evolving functions / vector fields on the sphere once we get it working.

MikaelSlevinsky · 2018-03-05T17:01:00Z

This is really cool! Care to give ballpark timings? (PlotlyJS worked great up to about 512 x 1024 grids)

More of a transforms-related issue, but the grids probably shouldn't cluster at the poles for a quiver plot.

dlfivefifty · 2018-03-05T17:02:26Z

Is there Julia code for generating uniformly distributed points on the sphere?

MikaelSlevinsky · 2018-03-05T17:14:01Z

There's Matlab code for minimum energy or maximum determinant nodes https://github.com/gradywright/spherepts/blob/master/docs/NodesGridsMeshesOnSphere.pdf

but perhaps HEALPix would be better transforms-wise.

MikaelSlevinsky · 2018-03-05T17:17:38Z

Oh wait, the min energy and max det use huge look-up tables. https://github.com/gradywright/spherepts/tree/master/nodes

dlfivefifty · 2018-03-05T17:20:54Z

It’s MIT license so we can just copy the tables

SebastianM-C · 2018-03-05T18:54:34Z

You can generate uniformly distributed points oh the unit sphere with something like this:

function SphericalToCartesian(r::T,θ::T,ϕ::T) where T<:AbstractArray
    x = @.r*sin(θ)*cos(ϕ)
    y = @.r*sin(θ)*sin(ϕ)
    z = @.r*cos(θ)
    
    x,y,z
end

n = 1000 # number of points to generate

r = ones(n);
θ = acos.(1 .- 2 .* rand(n));
ϕ = 2π * rand(n);

x,y,z = SphericalToCartesian(r,θ,ϕ);

dlfivefifty · 2018-03-05T19:25:40Z

I think random points are not a good idea since the arrows will be hard to read. Also you’ll have clustering of points.

PS here’s another random points on the sphere which is coordinate independent:

qr(randn(3,3))[1][:,1]

SimonDanisch · 2018-03-05T22:34:36Z

So this works now on master:

using GeometryTypes, Colors, Makie
# the surface
scene = Scene()
s = Makie.surface(x, y, z, image = F, colormap = :viridis, colornorm = (-1.0, 1.0))
s[:image] = F .+ 0.1 # update image

# the quiver

n = 20
f   = (x,y,z) -> x*exp(cos(y)*z)
∇f  = (x,y,z) -> Point3f0(exp(cos(y)*z), -sin(y)*z*x*exp(cos(y)*z), x*cos(y)*exp(cos(y)*z))
∇ˢf = (x,y,z) -> ∇f(x,y,z) - Point3f0(x,y,z)*dot(Point3f0(x,y,z), ∇f(x,y,z))

θ = [0;(0.5:n-0.5)/n;1]
φ = [(0:2n-2)*2/(2n-1);2]
x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
z = [cospi(θ) for θ in θ, φ in φ]

pts = vec(Point3f0.(x, y, z))
∇ˢF = vec(∇ˢf.(x, y, z))
scene = Scene();

surface(x, y, z)
arr = Makie.arrows(
    scene,
    pts, ∇ˢF,
    arrowsize = 0.03, linecolor = (:gray, 0.7), linewidth = 3
)
arr[:lengthscale] = 0.3

center!(scene)
# record a video
io = VideoStream(scene, joinpath(homedir(), "Desktop"), "test")
arr[:directions] = ∇ˢF
for i = 1:100
    arr[:directions] = (to_value(arr, :directions) .+ (rand(Float32, length(∇ˢF)) * 0.05))
    recordframe!(io)
    sleep(1/30)
end
finish(io, "mp4") # could also be gif, webm or mkv

SimonDanisch · 2018-03-06T09:29:39Z

Please let me know about any problems in a new issue ;)

dlfivefifty · 2018-03-06T10:11:29Z

The code work, and I checked and all the arrows are in the tangent plane, yay! It's a bit clearer when we normalize and scale:

scene = Scene(); surface(scene , x, y, z); arr = Makie.arrows(
                  scene,
                  pts, normalize.(∇ˢF)./15,
                  arrowsize = 0.03, linecolor = (:gray, 0.7), linewidth = 3
              )

Some questions:

You never gave the code for the density map plot
For the arrows on the sphere, does the fact that they are (circular) arcs help? I.e., is there an arc analogue to linesegment (used to make the arrows)?

SimonDanisch · 2018-03-06T10:44:38Z

The code for the density map is above:

using GeometryTypes, Colors, Makie
# the surface
scene = Scene()
s = Makie.surface(x, y, z, image = F, colormap = :viridis, colornorm = (-1.0, 1.0))
s[:image] = F .+ 0.1 # update image

dlfivefifty · 2018-03-06T10:45:06Z

Forgot one question:

\3. How do you clear the scene?

SimonDanisch · 2018-03-06T10:45:55Z

I usually just make a new scene! In the future, I plan to offer a more Plots.jl like api, where plot (without !) will just create a new scene!

SimonDanisch · 2018-03-06T10:47:18Z

Btw, there are arguments for normalization and scaling:
https://github.com/JuliaPlots/Makie.jl/blob/master/src/plotsbase/basic_recipes.jl#L14
Not sure if lengthscale is a particularly nice name^^

SimonDanisch · 2018-03-06T10:56:18Z

2.) For the arrows on the sphere, does the fact that they are (circular) arcs help? I.e., is there an arc analogue to linesegment (used to make the arrows)?

We can definitely just draw lines with a high enough granularity, to look like they're bend along the sphere ;) Could you supply me with an example code, using lines(arc::Vector{Point3f0}), to show how to translate ∇ˢF into a vector field of arcs?

dlfivefifty · 2018-03-06T11:20:07Z

At that point we might as well do streamlines. Here's a first attempt (it works but is dreadfully slow, and needs a better point distribution to look good):

function sphere_streamline(∇ˢf, pt, h=0.01f0, n=5)    
    ret = [pt]
    df = normalize(∇ˢf(pt...))
    push!(ret, normalize(pt .+ h*df))
    for k=2:n
        cur_pt = last(ret)
        push!(ret, cur_pt)
        df = normalize(∇ˢf(cur_pt...))
        push!(ret, normalize(cur_pt .+ h*df))
    end
    ret
end

function plotstreamline(scene, ret)
    linecolor = :gray
    linewidth = 3
    scale = Vec3f0(1)
    linesegment(scene, ret, color = linecolor, linewidth = linewidth, linestyle = nothing, scale = scale)
    Makie.arrows(scene, [ret[end-1]], [ret[end]-ret[end-1]], arrowsize = 0.01, linecolor = (:gray, 0.7), linewidth = 3)
end

lns = sphere_streamline.(∇ˢf, pts, 0.03f0)

scene = Scene();
    for ln in lns
        plotstreamline(scene, ln)
    end

SimonDanisch · 2018-03-06T11:48:28Z

ok, let's make it fast :)

SimonDanisch · 2018-03-06T11:49:53Z

Wow, it's slower than I imagined ;) Good news is, it should be easy to fix!

SimonDanisch · 2018-03-06T11:51:28Z

I didn't expect it to be fast, but there must also be a performance bug, since it seems to become at least quadratically slower

SimonDanisch · 2018-03-06T12:28:35Z

This seems to be reasonably fast:

using Makie, GeometryTypes
function sphere_streamline(linebuffer, ∇ˢf, pt, h=0.01f0, n=5)
    push!(linebuffer, pt)
    df = normalize(∇ˢf(pt[1], pt[2], pt[3]))
    push!(linebuffer, normalize(pt .+ h*df))
    for k=2:n
        cur_pt = last(linebuffer)
        push!(linebuffer, cur_pt)
        df = normalize(∇ˢf(cur_pt...))
        push!(linebuffer, normalize(cur_pt .+ h*df))
    end
    return
end

function streamlines(
        scene, ∇ˢf, pts::AbstractVector{T};
        h=0.01f0, n=5, color = :black, linewidth = 1
    ) where T
    linebuffer = T[]
    sub = Scene(
        scene,
        h = h, n = 5, color = :black, linewidth = 1
    )
    lines = lift_node(to_node(∇ˢf), to_node(pts), sub[:h], sub[:n]) do ∇ˢf, pts, h, n
        empty!(linebuffer)
        for point in pts
            sphere_streamline(linebuffer, ∇ˢf, point, h, n)
        end
        linebuffer
    end
    linesegment(sub, lines, color = sub[:color], linewidth = sub[:linewidth])
    sub
end

# needs to be in a function for ∇ˢf to be fast and inferable
function test()
    n = 20
    f   = (x,y,z) -> x*exp(cos(y)*z)
    ∇f  = (x,y,z) -> Point3f0(exp(cos(y)*z), -sin(y)*z*x*exp(cos(y)*z), x*cos(y)*exp(cos(y)*z))
    ∇ˢf = (x,y,z) -> ∇f(x,y,z) - Point3f0(x,y,z)*dot(Point3f0(x,y,z), ∇f(x,y,z))
    θ = [0;(0.5:n-0.5)/n;1]
    φ = [(0:2n-2)*2/(2n-1);2]
    x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
    y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
    z = [cospi(θ) for θ in θ, φ in φ]

    pts = vec(Point3f0.(x, y, z))
    scene = Scene()
    lns = streamlines(scene, ∇ˢf, pts)
    # those can be changed interactively:
    lns[:color] = :black
    lns[:h] = 0.06
    lns[:linewidth] = 1.0
    for i = linspace(0.01, 0.1, 100)
        lns[:h] = i
        yield()
    end
end
test()

loicspace · 2018-07-09T23:40:31Z

Hi Hi,

I'm trying to play some with the arrows, and I'm using the code Simon put in there,

using GeometryTypes, Colors, Makie 
# the quiver
n = 20
f   = (x,y,z) -> x*exp(cos(y)*z)
∇f  = (x,y,z) -> Point3f0(exp(cos(y)*z), -sin(y)*z*x*exp(cos(y)*z), x*cos(y)*exp(cos(y)*z))
∇ˢf = (x,y,z) -> ∇f(x,y,z) - Point3f0(x,y,z)*dot(Point3f0(x,y,z), ∇f(x,y,z))

θ = [0;(0.5:n-0.5)/n;1]
φ = [(0:2n-2)*2/(2n-1);2]
x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
z = [cospi(θ) for θ in θ, φ in φ]

pts = vec(Point3f0.(x, y, z))
∇ˢF = vec(∇ˢf.(x, y, z))
scene = Scene();

surface(x, y, z)
arr = Makie.arrows(scene,pts, ∇ˢF)

But I am getting the following error:

DimensionMismatch("tuples could not be broadcast to a common size")
Arrows{...}(::AbstractPlotting.Scene, ::AbstractPlotting.Attributes, ::Tuple{AbstractPlotting.Scene,Array{GeometryTypes.Point{3,Float32},1},Array{GeometryTypes.Point{3,Float32},1}}) at interfaces.jl:288
plot!(::AbstractPlotting.Scene, ::Type{Arrows{...}}, ::AbstractPlotting.Attributes, ::AbstractPlotting.Scene, ::Vararg{Any,N} where N) at interfaces.jl:382
#plot!#175(::Array{Any,1}, ::Function, ::AbstractPlotting.Scene, ::Type{T} where T, ::AbstractPlotting.Scene, ::Vararg{Any,N} where N) at interfaces.jl:364
plot!(::AbstractPlotting.Scene, ::Type{T} where T, ::AbstractPlotting.Scene, ::Array{GeometryTypes.Point{3,Float32},1}, ::Vararg{Array{GeometryTypes.Point{3,Float32},1},N} where N) at interfaces.jl:364
#plot#169(::Array{Any,1}, ::Function, ::Type{T} where T, ::AbstractPlotting.Scene, ::Vararg{Any,N} where N) at interfaces.jl:351
plot(::Type{T} where T, ::AbstractPlotting.Scene, ::Array{GeometryTypes.Point{3,Float32},1}, ::Vararg{Array{GeometryTypes.Point{3,Float32},1},N} where N) at interfaces.jl:351
#arrows#258(::Array{Any,1}, ::Function, ::AbstractPlotting.Scene, ::Vararg{Any,N} where N) at recipes.jl:77
arrows(::AbstractPlotting.Scene, ::Array{GeometryTypes.Point{3,Float32},1}, ::Vararg{Array{GeometryTypes.Point{3,Float32},1},N} where N) at recipes.jl:77
include_string(::String, ::String) at loading.jl:522
include_string(::String, ::String, ::Int64) at eval.jl:30
include_string(::Module, ::String, ::String, ::Int64, ::Vararg{Int64,N} where N) at eval.jl:34
(::Atom.##102#107{String,Int64,String})() at eval.jl:82
withpath(::Atom.##102#107{String,Int64,String}, ::String) at utils.jl:30
withpath(::Function, ::String) at eval.jl:38
hideprompt(::Atom.##101#106{String,Int64,String}) at repl.jl:62
macro expansion at eval.jl:80 [inlined]
(::Atom.##100#105{Dict{String,Any}})() at task.jl:80

I have both Makie and AbstractPlotting checked out. What am I missing?
Cheers,

Loic

mateuszbaran · 2021-07-17T21:02:15Z

In case someone needs an updates version:

using GLMakie, Makie, LinearAlgebra

n = 10

θ = [0;(0.5:n-0.5)/n;1]
φ = [(0:2n-2)*2/(2n-1);2]
x = [cospi(φ)*sinpi(θ) for θ in θ, φ in φ]
y = [sinpi(φ)*sinpi(θ) for θ in θ, φ in φ]
z = [cospi(θ) for θ in θ, φ in φ]

f2(x, y, z) = cross([x, y, z], [1.0, 0.0, 0.0])
tans = f2.(vec(x), vec(y), vec(z))
u = [a[1] for a in tans]
v = [a[2] for a in tans]
w = [a[3] for a in tans]

scene = Scene();

surface(x, y, z)
arr = Makie.arrows!(
    vec(x), vec(y), vec(z), u, v, w;
    arrowsize = 0.1, linecolor = (:gray, 0.7), linewidth = 0.02, lengthscale = 0.1
)

SimonDanisch closed this as completed Mar 6, 2018

ffreyer mentioned this issue May 23, 2019

Streamplot in Makie #355

Closed

hacklint mentioned this issue Jun 5, 2021

Keyboard controlled 3d camera #1024

Closed

19 tasks

jakebolewski mentioned this issue Oct 4, 2021

Add ability to plot unstructred 2D mesh on the sphere CliMA/ClimaCore.jl#224

Closed

Useful for visualizing on a sphere? #72

Useful for visualizing on a sphere? #72

Comments

dlfivefifty commented Mar 2, 2018

SimonDanisch commented Mar 2, 2018

dlfivefifty commented Mar 2, 2018 • edited Loading

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018 • edited Loading

dlfivefifty commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018 • edited Loading

SimonDanisch commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018 • edited Loading

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018 • edited Loading

dlfivefifty commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

MikaelSlevinsky commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

MikaelSlevinsky commented Mar 5, 2018

MikaelSlevinsky commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

SebastianM-C commented Mar 5, 2018

dlfivefifty commented Mar 5, 2018

SimonDanisch commented Mar 5, 2018

SimonDanisch commented Mar 6, 2018

dlfivefifty commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

dlfivefifty commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

dlfivefifty commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

SimonDanisch commented Mar 6, 2018

loicspace commented Jul 9, 2018

mateuszbaran commented Jul 17, 2021

dlfivefifty commented Mar 2, 2018 •

edited

Loading

dlfivefifty commented Mar 5, 2018 •

edited

Loading

dlfivefifty commented Mar 5, 2018 •

edited

Loading

dlfivefifty commented Mar 5, 2018 •

edited

Loading

dlfivefifty commented Mar 5, 2018 •

edited

Loading