How to show large number of lines using Julia, GLMakie

I try to show the road network, which has over 400,000 segments. I first tried the package OpenStreetMapXPlot.jl. However, it takes very long time and cannot plot. Then, I tried to use GLMakie.jl

using Makie, GeometryTypes

scene = Scene()
for r in roads   # over 400,000 roads
     xs = ...
     ys = ...
     lines!(scene, xs, ys, color = :blue, linewidth = 1)
end

It takes a long time to plot 10,000 road segments. It seems also impossible to show the large network.
I use the software QGIS, which is developed using C++, to show the road network and it works very smooth.

Is GLMakie.jl the proper package to show the large road network? Or some other technology is better?
Thanks!

2 Likes

I tried similar code using Plots.jl. The running time for plotting 1000 roads is,

Makie: 6.929608583 s
Plots: 1.138019358 s

The package Plots.jl seems several times quicker. For 400k roads, Plots still takes very long time. Is there other method to show large amount of lines using Julia? Thanks!

Maybe you can avoid the calls to the plotting functions inside the loop? Have you tried calling lines! only once on a large vector of xs and ys, potentially separated by NaNs?

For example, I get a 200x speed up doing that with Plots.jl on this example:

xs = [rand(100) for _ in 1:1000]
ys = [rand(100) for _ in 1:1000]

function plot1(xs, ys)
    plt = plot()
    for (x, y) in zip(xs, ys)
        plot!(plt, x, y)
    end
    plt
end

function plot2(xs, ys)
    plot(reduce(vcat, xs), reduce(vcat, ys))
end

@btime plot1($xs, $ys) ; # 431.335 ms (2231129 allocations: 121.87 MiB)
@btime plot2($xs, $ys) ; # 2.844 ms (2530 allocations: 3.20 MiB)
# BTW in this case a simple call to plot(xs, ys) works, but it slow for me:
@btime plot($xs, $ys) ; # 205.678 ms (1374088 allocations: 72.29 MiB)
2 Likes

Thanks for your reply!
plot2() is different from plot1(). plot1() draws 1000 separated lines, while plot2() draw one line, which links 1000 lines as one. This would be the reason while plot2() is much quicker.
The following are figures of three lines by these two function.

plot1()

plot2()

I think if you do reduce(hcat, ...) instead, they will be plotted as separate series. Not sure how that affects timing, but might be worth trying

Thank you so much! It plots the correct figure. However, the running time is only slight quicker, which is as follows.

  plot1() running time:  380.361 ms (2491155 allocations: 129.92 MiB)
  plot2() running time:  209.021 ms (1610075 allocations: 79.49 MiB)

QGIS, written in c++, can render this network smoothly. Is there other package to draw large network? There should be some method to draw using Julia.

you can try using GR

with this you can plot millions of points [lines].

shadelines

is what you are looking for.

1 Like

I think GLMakie works well with this, as long as you plot everything as one line (use NaN points as gaps).

Here are random “roads” with 900,000 segments

using GLMakie

roads = [i % 10 == 0 ? Point2f0(NaN) : 0.1f0 * randn(Point2f0) .+ Point2f0(x, y) for x in 1:300 for y in 1:300 for i in 1:10]

lines(roads)
julia> @time lines(roads)
  0.022839 seconds (44.10 k allocations: 2.294 MiB)
6 Likes

So amazing! It works! The large urban traffic network, with over 400,000 road segments, can be shown within 1 second!

415.336 ms (534413 allocations: 77.38 MiB)

Thank you so much!

Awesome! The general rule is to have as few high level plot objects as possible, so that the gpu receives only a couple of large arrays. Then it will be fast. 100000 line plot objects have too much overhead

2 Likes

No, it’s the single call to plot that makes it quicker. As others have also said, you just need to concatenate your roads in a single array and separate the segments that need separation by NaNs and it should be fast with Plots. (I did not do the part adding the NaNs in my example.)


[EDIT] For example, to sort-of match the GLMakie one,

xs = [0.1f0 * randn() .+ x for x in 1:300 for y in 1:300 for i in 1:10]
ys = [0.1f0 * randn() .+ y for x in 1:300 for y in 1:300 for i in 1:10]
roads = [i % 10 == 0 for x in 1:300 for y in 1:300 for i in 1:10]
isnotroad = [i % 10 == 0 for x in 1:300 for y in 1:300 for i in 1:10]
xs[isnotroad] .= NaN
ys[isnotroad] .= NaN
@btime plot($xs, $ys) # 23.933 ms for me

and these show what you want as well:

plot(plot(xs, ys), plot(xs, ys, xlim=(100, 110), ylim=(120, 130)))

gives

4 Likes

Thank you so much! I will have a try. It would have similar speed as the Jules’s suggestion.

1 Like

I wonder if high-level plot objects can be automatically “rewritten” to this more efficient representation.

3 Likes

For the record, the inspectdr() backend seems to plot at same speed as gr() in this case, while allowing interactive zooming in and out. However it is 30% slower here than Makie. Makie zoom is nicer as scales are updated, while inspectdr() provides the coordinates of the cursor.

using Plots; inspectdr(legend=false)
xs = [i % 10 == 0 ? NaN : 0.1f0 * randn() .+ x for x in 1:300 for y in 1:300 for i in 1:10]
ys = [i % 10 == 0 ? NaN : 0.1f0 * randn() .+ y for x in 1:300 for y in 1:300 for i in 1:10]
@btime Plots.plot($xs, $ys)  # 18 ms (2436 allocations: 13.88 MiB)
# draw rectangle with RMB to zoom-in, CTRL+f to unzoom
2 Likes

This works very well. Thank you.