diff --git a/Project.toml b/Project.toml index b32b40ff..78605096 100644 --- a/Project.toml +++ b/Project.toml @@ -4,6 +4,7 @@ authors = ["SimonDanisch "] version = "0.2.15" [deps] +EarCut_jll = "5ae413db-bbd1-5e63-b57d-d24a61df00f5" IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" diff --git a/src/GeometryBasics.jl b/src/GeometryBasics.jl index b8fc88ee..183a4cf1 100644 --- a/src/GeometryBasics.jl +++ b/src/GeometryBasics.jl @@ -1,6 +1,7 @@ module GeometryBasics using StaticArrays, Tables, StructArrays, IterTools, LinearAlgebra + using EarCut_jll using Base: @propagate_inbounds @@ -12,15 +13,15 @@ module GeometryBasics include("viewtypes.jl") include("geometry_primitives.jl") include("rectangles.jl") - include("triangulation.jl") include("meshes.jl") + include("triangulation.jl") include("lines.jl") include("boundingboxes.jl") export AbstractGeometry, GeometryPrimitive export Mat, Point, Vec export LineFace, Polytope, Line, NgonFace, convert_simplex - export LineString, Polygon, MultiPoint, MultiLineString, MultiPolygon + export LineString, AbstractPolygon, Polygon, MultiPoint, MultiLineString, MultiPolygon export Simplex, connect, Triangle, NSimplex, Tetrahedron export QuadFace, metafree, coordinates, TetrahedronFace export TupleView, SimplexFace, Mesh, meta @@ -36,7 +37,7 @@ module GeometryBasics export AbstractMesh, Mesh, TriangleMesh export GLNormalMesh2D, PlainTriangleMesh export MetaT, meta_table - + # all the different predefined mesh types # Note: meshes can contain arbitrary meta information, export AbstractMesh, TriangleMesh, PlainMesh, GLPlainMesh, GLPlainMesh2D, GLPlainMesh3D diff --git a/src/basic_types.jl b/src/basic_types.jl index cce77f80..ef3055aa 100644 --- a/src/basic_types.jl +++ b/src/basic_types.jl @@ -92,7 +92,6 @@ function Polytope(::Type{<: NNgon{N}}, P::Type{<: AbstractPoint{NDim, T}}) where return Ngon{NDim, T, N, P} end - const LineP{Dim, T, P <: AbstractPoint{Dim, T}} = Ngon{Dim, T, 2, P} const Line{Dim, T} = LineP{Dim, T, Point{Dim, T}} @@ -111,6 +110,17 @@ const Quadrilateral{Dim, T} = Ngon{Dim, T, 4, P} where P <: AbstractPoint{Dim, T Base.show(io::IO, x::Quadrilateral) = print(io, "Quad(", join(x, ", "), ")") Base.summary(io::IO, x::Type{<: Quadrilateral}) = print(io, "Quad") +function coordinates(lines::AbstractArray{LineP{Dim, T, PointType}}) where {Dim, T, PointType} + if lines isa Base.ReinterpretArray + return coordinates(lines.parent) + else + result = PointType[] + for line in lines + append!(result, coordinates(line)) + end + return result + end +end """ A `Simplex` is a generalization of an N-dimensional tetrahedra and can be thought @@ -181,11 +191,11 @@ struct LineString{ points::V end -coordinates(x::LineString) = x.points +coordinates(x::LineString) = coordinates(x.points) Base.copy(x::LineString) = LineString(copy(x.points)) -Base.size(x::LineString) = size(coordinates(x)) -Base.getindex(x::LineString, i) = getindex(coordinates(x), i) +Base.size(x::LineString) = size(getfield(x, :points)) +Base.getindex(x::LineString, i) = getindex(getfield(x, :points), i) function LineString(points::AbstractVector{LineP{Dim, T, P}}) where {Dim, T, P} return LineString{Dim, T, P, typeof(points)}(points) @@ -276,6 +286,27 @@ function Polygon(exterior::AbstractVector{P}, faces::AbstractVector{<: LineFace} return Polygon(LineString(exterior, faces)) end +function Polygon(exterior::AbstractVector{P}, interior::AbstractVector{<:AbstractVector{P}}) where P <: AbstractPoint{Dim, T} where {Dim, T} + ext = LineString(exterior) + # We need to take extra care for empty interiors, since + # if we just map over it, it won't infer the element type correctly! + int = typeof(ext)[] + foreach(x-> push!(int, LineString(x)), interior) + return Polygon(ext, int) +end + +function coordinates(polygon::Polygon{N, T, PointType}) where {N, T, PointType} + exterior = coordinates(polygon.exterior) + if isempty(polygon.interiors) + return exterior + else + result = PointType[] + append!(result, exterior) + foreach(x-> append!(result, coordinates(x)), polygon.interiors) + return result + end +end + """ MultiPolygon(polygons::AbstractPolygon) """ diff --git a/src/interfaces.jl b/src/interfaces.jl index 90992f5e..16a5e662 100644 --- a/src/interfaces.jl +++ b/src/interfaces.jl @@ -75,13 +75,12 @@ faces(tesselation::Tesselation) = faces(tesselation.primitive, nvertices(tessela normals(tesselation::Tesselation) = normals(tesselation.primitive, nvertices(tesselation)) texturecoordinates(tesselation::Tesselation) = texturecoordinates(tesselation.primitive, nvertices(tesselation)) - ## Decompose methods # Dispatch type to make `decompose(UV{Vec2f0}, primitive)` work # and to pass through tesselation information # Types that can be converted to a mesh via the functions below -const Meshable{Dim, T} = Union{Tesselation{Dim, T}, Mesh{Dim, T}, +const Meshable{Dim, T} = Union{Tesselation{Dim, T}, Mesh{Dim, T}, AbstractPolygon{Dim, T}, GeometryPrimitive{Dim, T}, AbstractVector{<: AbstractPoint{Dim, T}}} struct UV{T} end @@ -108,23 +107,14 @@ function decompose(::Type{Point}, primitive::Meshable{Dim, T}) where {Dim, T} return collect_with_eltype(Point{Dim, T}, metafree(coordinates(primitive))) end -function decompose(::Type{T}, primitive) where {T} - return collect_with_eltype(T, primitive) +function decompose(::Type{Point}, primitive::LineString{Dim, T}) where {Dim, T} + return collect_with_eltype(Point{Dim, T}, metafree(coordinates(primitive))) end -function decompose(::Type{P}, pol::Polygon) where {P<:AbstractPoint} - if isempty(pol.interiors) - return decompose(P, pol.exterior) - else - arr = copy(decompose(P, pol.exterior)) - for i in pol.interiors - append!(arr, decompose(P, i)) - end - return arr - end +function decompose(::Type{T}, primitive) where {T} + return collect_with_eltype(T, primitive) end -decompose(::Type{P}, ls::LineString) where {P<:AbstractPoint} = ls.points.parent.data decompose_uv(primitive) = decompose(UV(), primitive) decompose_uvw(primitive) = decompose(UVW(), primitive) decompose_normals(primitive) = decompose(Normal(), primitive) diff --git a/src/meshes.jl b/src/meshes.jl index 5c1821fd..9631b53c 100644 --- a/src/meshes.jl +++ b/src/meshes.jl @@ -14,7 +14,6 @@ function normals(mesh::AbstractMesh) return nothing end - const GLTriangleElement = Triangle{3, Float32} const GLTriangleFace = TriangleFace{GLIndex} const PointWithUV{Dim, T} = PointMeta{Dim, T, Point{Dim, T}, (:uv,), Tuple{Vec{2, T}}} @@ -144,6 +143,12 @@ Polygon triangluation! function mesh(polygon::AbstractVector{P}; pointtype=P, facetype=GLTriangleFace, normaltype=nothing) where {P<:AbstractPoint{2}} + return mesh(Polygon(polygon); pointtype, facetype, normaltype) +end + +function mesh(polygon::AbstractPolygon{Dim, T}; pointtype=Point{Dim, T}, facetype=GLTriangleFace, + normaltype=nothing) where {Dim, T} + faces = decompose(facetype, polygon) positions = decompose(pointtype, polygon) @@ -151,7 +156,6 @@ function mesh(polygon::AbstractVector{P}; pointtype=P, facetype=GLTriangleFace, n = normals(positions, faces; normaltype=normaltype) positions = meta(positions; normals=n) end - return Mesh(positions, faces) end diff --git a/src/triangulation.jl b/src/triangulation.jl index 9b1e14a1..51efa8b2 100644 --- a/src/triangulation.jl +++ b/src/triangulation.jl @@ -169,3 +169,66 @@ function decompose(::Type{FaceType}, points::AbstractArray{P}) where {P<:Abstrac return result end + + +function earcut_triangulate(polygon::Vector{Vector{Point{2, Float64}}}) + lengths = map(x-> UInt32(length(x)), polygon) + len = UInt32(length(lengths)) + array = ccall( + (:u32_triangulate_f64, libearcut), + Tuple{Ptr{GLTriangleFace}, Cint}, + (Ptr{Ptr{Float64}}, Ptr{UInt32}, UInt32), + polygon, lengths, len + ) + return unsafe_wrap(Vector{GLTriangleFace}, array[1], array[2]) +end + +function earcut_triangulate(polygon::Vector{Vector{Point{2, Float32}}}) + lengths = map(x-> UInt32(length(x)), polygon) + len = UInt32(length(lengths)) + array = ccall( + (:u32_triangulate_f32, libearcut), + Tuple{Ptr{GLTriangleFace}, Cint}, + (Ptr{Ptr{Float32}}, Ptr{UInt32}, UInt32), + polygon, lengths, len + ) + return unsafe_wrap(Vector{GLTriangleFace}, array[1], array[2]) +end + +function earcut_triangulate(polygon::Vector{Vector{Point{2, Int64}}}) + lengths = map(x-> UInt32(length(x)), polygon) + len = UInt32(length(lengths)) + array = ccall( + (:u32_triangulate_i64, libearcut), + Tuple{Ptr{GLTriangleFace}, Cint}, + (Ptr{Ptr{Int64}}, Ptr{UInt32}, UInt32), + polygon, lengths, len + ) + return unsafe_wrap(Vector{GLTriangleFace}, array[1], array[2]) +end + +function earcut_triangulate(polygon::Vector{Vector{Point{2, Int32}}}) + lengths = map(x-> UInt32(length(x)), polygon) + len = UInt32(length(lengths)) + array = ccall( + (:u32_triangulate_i32, libearcut), + Tuple{Ptr{GLTriangleFace}, Cint}, + (Ptr{Ptr{Int32}}, Ptr{UInt32}, UInt32), + polygon, lengths, len + ) + return unsafe_wrap(Vector{GLTriangleFace}, array[1], array[2]) +end + +best_earcut_eltype(x) = Float64 +best_earcut_eltype(::Type{Float64}) = Float64 +best_earcut_eltype(::Type{<:AbstractFloat}) = Float32 +best_earcut_eltype(::Type{Int64}) = Int64 +best_earcut_eltype(::Type{Int32}) = Int32 +best_earcut_eltype(::Type{<:Integer}) = Int64 + +function faces(polygon::Polygon{Dim, T}) where {Dim, T} + PT = Point{Dim, best_earcut_eltype(T)} + points = [decompose(PT, polygon.exterior)] + foreach(x-> push!(points, decompose(PT, x)), polygon.interiors) + return earcut_triangulate(points) +end diff --git a/src/viewtypes.jl b/src/viewtypes.jl index e1e4d2eb..7b4fdd46 100644 --- a/src/viewtypes.jl +++ b/src/viewtypes.jl @@ -36,6 +36,8 @@ struct TupleView{T, N, Skip, A} <: AbstractVector{T} connect::Bool end +coordinates(tw::TupleView) = coordinates(tw.data) + function Base.size(x::TupleView{T, N, M}) where {T, N, M} nitems = length(x.data) รท (N - (N - M)) nitems = nitems - max(N - M, 0) diff --git a/test/runtests.jl b/test/runtests.jl index 080e3df2..f20e73e3 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -56,7 +56,7 @@ using GeometryBasics: attributes end end - + @testset "polygon with metadata" begin polys = [Polygon(rand(Point{2, Float32}, 20)) for i in 1:10] pnames = [randstring(4) for i in 1:10] @@ -68,7 +68,7 @@ using GeometryBasics: attributes multipoly = MultiPolygonMeta(polys, name = pnames, value = numbers, category = bin) @test multipoly isa AbstractVector @test poly isa GeometryBasics.AbstractPolygon - + @test GeometryBasics.getcolumn(poly, :name) == pnames[1] @test GeometryBasics.MetaFree(PolygonMeta) == Polygon @@ -92,7 +92,7 @@ using GeometryBasics: attributes @test metafree(pm) === p @test propertynames(pm) == (:position, :a, :b) end - + @testset "MultiPoint with metadata" begin p = collect(Point{2, Float64}(x, x+1) for x in 1:5) @test p isa AbstractVector @@ -144,10 +144,10 @@ end multipoly = MetaT(multipol, name = pnames, value = numbers, category = bin) @test multipoly isa MetaT @test poly isa MetaT - + @test GeometryBasics.getcolumn(poly, :name) == pnames[1] @test GeometryBasics.getcolumn(multipoly, :name) == pnames - + meta_p = MetaT(polys[1], boundingbox=Rect(0, 0, 2, 2)) @test meta_p.boundingbox === Rect(0, 0, 2, 2) @test GeometryBasics.metafree(meta_p) == polys[1] @@ -155,7 +155,7 @@ end @test GeometryBasics.metafree(multipoly) == multipol @test GeometryBasics.meta(meta_p) == (boundingbox = GeometryBasics.HyperRectangle{2,Int64}([0, 0], [2, 2]),) @test GeometryBasics.meta(poly) == (name = pnames[1], value = 0.0, category = bin[1]) - @test GeometryBasics.meta(multipoly) == (name = pnames, value = numbers, category = bin) + @test GeometryBasics.meta(multipoly) == (name = pnames, value = numbers, category = bin) end @testset "MetaT{Point}" begin @@ -171,7 +171,7 @@ end @test GeometryBasics.metafree(pm) == p @test GeometryBasics.meta(pm) == (a = 1, b = 2) end - + @testset "MetaT{MultiPoint}" begin p = collect(Point{2, Float64}(x, x+1) for x in 1:5) @test p isa AbstractVector @@ -623,7 +623,7 @@ end @testset "MetaT and heterogeneous data" begin ls = [LineString([Point(i, (i+1)^2/6), Point(i*0.86,i+5), Point(i/3, i/7)]) for i in 1:10] - mls = MultiLineString([LineString([Point(i+1, (i)^2/6), Point(i*0.75,i+8), Point(i/2.5, i/6.79)]) for i in 5:10]) + mls = MultiLineString([LineString([Point(i+1, (i)^2/6), Point(i*0.75,i+8), Point(i/2.5, i/6.79)]) for i in 5:10]) poly = Polygon(Point{2, Int}[(40, 40), (20, 45), (45, 30), (40, 40)]) geom = [ls..., mls, poly] prop = Any[(country_states = "India$(i)", rainfall = (i*9)/2) for i in 1:11] @@ -637,17 +637,14 @@ end @test propertynames(sa) === (:main, :country_states, :rainfall) @test getproperty(sa, :country_states) isa Array{Any} @test getproperty(sa, :main) == geom - - @test GeometryBasics.getnamestypes(typeof(feat[1])) == + + @test GeometryBasics.getnamestypes(typeof(feat[1])) == (LineString{2,Float64,Point{2,Float64},Base.ReinterpretArray{GeometryBasics.Ngon{2,Float64,2,Point{2,Float64}},1,Tuple{Point{2,Float64},Point{2,Float64}},TupleView{Tuple{Point{2,Float64},Point{2,Float64}},2,1,Array{Point{2,Float64},1}}}}, (:country_states, :rainfall), Tuple{String,Float64}) - - @test StructArrays.staticschema(typeof(feat[1])) == - NamedTuple{(:main, :country_states, :rainfall),Tuple{LineString{2,Float64,Point{2,Float64},Base.ReinterpretArray{GeometryBasics.Ngon{2,Float64,2,Point{2,Float64}},1,Tuple{Point{2,Float64},Point{2,Float64}},TupleView{Tuple{Point{2,Float64},Point{2,Float64}},2,1,Array{Point{2,Float64},1}}}}, - String,Float64}} + @test StructArrays.createinstance(typeof(feat[1]), LineString([Point(1, (2)^2/6), Point(1*0.86,6), Point(1/3, 1/7)]), "Mumbai", 100) isa typeof(feat[1]) - + @test Base.getindex(feat[1], 1) isa Line @test Base.size(feat[1]) == (2,) end