Add voronoi grid

[virgile-baudrot]

The following lines provide a simple way to make grid landscape using voronoi tesselation. It's working but need to be reframe to fit package requirement.

Computing the voronoi is done using VoronoiCells.jl and then the rasterization is by using ArchGDAL.jl.

Digging into both dependencies, I should be able to make things clearer.

using VoronoiCells
using GeometryBasics
using ArchGDAL
using Distributions
using Plots

const VC=VoronoiCells
const GB=GeometryBasics
const AG=ArchGDAL
const DB=Distributions

I created 3 functions. The 2 first create the set of polygones

function voronoitess(min_x, max_x, min_y, max_y, nCell::Int)
    rect = VC.Rectangle(GB.Point2(min_x, min_y), GB.Point2(max_x, max_y))
    points = [GB.Point2(rand(DB.Uniform(min_x, max_x)), rand(DB.Uniform(min_y, max_y))) for _ in 1:nCell]
    return (tess=VC.voronoicells(points, rect),rect,points)
end

vortess = voronoitess(0, 1, 0, 1, 100)

function polygonset(tess)
    tc = tess.tess.Cells
    poly = [ArchGDAL.createpolygon([(tc[i][j][1], tc[i][j][2]) for j in vcat(1:length(tc[i]),1)])
    for i in 1:length(tc)]
    return (poly, tess.rect, tess.points)
end

polyset = polygonset(vortess)
plot(polyset[1])

And the last is making a rasterization providing some features:

# Create a grid over the map
function rasterize(polygonset, nrow, ncol,feature)
    w = rect.UpperRight[1] - rect.LowerLeft[1]
    h = rect.UpperRight[2] - rect.LowerLeft[2]
    gridcell_x=w/ncol
    gridcell_y=h/nrow
    gridX=(rect.LowerLeft[1]:gridcell_x:(rect.UpperRight[1]-gridcell_x)) .+ (gridcell_x/2)
    gridY=(rect.LowerLeft[2]:gridcell_y:(rect.UpperRight[2]-gridcell_y)) .+ (gridcell_y/2)

    gridX_seq=repeat(gridX,inner=nrow)
    gridY_seq=repeat(gridY,outer=ncol)

    gridCell = fill(0,nrow,ncol)
    for j in 1:length(gridCell)
        for i in 1:length(polygonset[1])
            if AG.contains(polygonset[1][i],
                AG.createpoint((gridX_seq[j], gridY_seq[j]))
            ) 
            gridCell[j] = feature[i]
            end        
        end
    end
    return gridCell
end

feature = rand(0:3,length(polyset[1]))
nrow, ncol = 200, 10

gridCell = rasterize(polyset, nrow, ncol, feature)
Plots.heatmap(gridCell)

[mkborregaard]

So is this a suggestion to add the functionality or to put this in the docs?
I wouldn't be keen to add an ArchGDAL dep to this repo. VoroniCells might be OK (it adds Colors, GeometryBasics, VoronoDelaunay, GeometricalPredicates and RecipesBase).

[virgile-baudrot]

Well I don't know what you can do with that 😄 . Just that voronoi tesselation would be nice to have so I propose to work on it if it's something you think is relevant.
Basically, Voronoi is only made from GeometryBasics and VoronoDelaunay, so I just need to remove dependency with ArchGDAL, finding a better algorithm to rasterize.

[rafaqz]

@virgile-baudrot you see my dependency concerns are widespread in the Julia world ;)

But this would be great, but without ArchGDAL.

[mkborregaard]

Just for clarity - isn't rasterized voronoi just nearest-neighbor from a grid to the set of points? We have nearestneighbors.jl as a dep already.

[mkborregaard]

E.g.

using NearestNeighbors
pts = rand(1.:1000, 2, 30)
tr = KDTree(pts)
mat = [knn(tr, [i, j], 1)[1][1] for i in 1:1000, j in 1:1000]

The heatmap:

[tpoisot]

I think it's a good improvement if we don't need ArchGDAL - I think @mkborregaard suggestion makes sense!

[tpoisot]

Wait, so is it NearestNeighborElement? I've added an option to set various values of k, but for k=1, this is a Voronoi tessellation?

[mkborregaard]

I was wondering if it wasn't just that yes. But couldn't tell from the code last night - I was heading home late on the train and was beer coding 😅

[mkborregaard]

But voronoi is just a vector format nearest neighbor yes

[virgile-baudrot]

Oh great. I didn't know that. So fantastic if it's already there!

[tpoisot]

Let's be explicit and add a VoronoiGrid type or something, as a thin wrapper on NearestNeighborElement? @virgile-baudrot do you want to do this? Happy to give feedback.

[virgile-baudrot]

ok, I'll try to implement it and PR

[tpoisot]

This was closed by #46