Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Sg/in on updates #37

Closed
wants to merge 39 commits into from
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
39 commits
Select commit Hold shift + click to select a range
aba090f
Fix up intersection point base calculation
skygering Oct 2, 2023
ab16785
Update intersects and add line tests
skygering Oct 4, 2023
7e75e86
Add more tests and debug intersects
skygering Oct 4, 2023
a7a7367
Add comments to point_in_poly
skygering Oct 4, 2023
b99e37d
Remove CairoMakie
skygering Oct 4, 2023
319bd88
Update equals and overlaps
skygering Oct 18, 2023
0b9799d
Remove use of findfirst for 1.6 compat
skygering Oct 18, 2023
90fff1c
Updated geom, multi-geom equality
skygering Oct 20, 2023
ab211f3
Merge branch 'asinghvi17:main' into sg/update_intersects
skygering Oct 20, 2023
6184971
Merge pull request #1 from Caltech-OCTO/sg/update_intersects
skygering Oct 20, 2023
37e0fc7
Add files for new functions
skygering Oct 26, 2023
92778b0
Start cleaning up geom on and in geom
skygering Oct 27, 2023
bce6381
Uncomment polygon_in_polygon
skygering Oct 27, 2023
98af38e
Add more geom in geom
skygering Oct 28, 2023
5034be2
Update line in curve and add comments
skygering Nov 8, 2023
e7e5df7
Generalized line in curve
skygering Nov 9, 2023
0aa1529
Working line line processor
skygering Nov 28, 2023
b557520
Finish within tests
skygering Nov 30, 2023
c8d71ad
Add multipolygon tests
skygering Dec 1, 2023
57a31d4
Add within comments
skygering Dec 1, 2023
3ddeb78
Add disjoint comments and multipolygon code
skygering Dec 1, 2023
da0cc2b
Add processor comments and clean
skygering Dec 13, 2023
be5c24f
Generalize processor
skygering Dec 19, 2023
c8af05a
Started poly-poly processor generalization
skygering Dec 20, 2023
d491545
Add tests and reorganize
skygering Dec 21, 2023
4443707
Merge branch 'main' of github.com:Caltech-OCTO/GeometryOps.jl
skygering Dec 21, 2023
4be9961
Continue to debug and test
skygering Dec 22, 2023
d9d636a
Update crosses algorithm
skygering Dec 26, 2023
ff4b9e5
New crosses code
skygering Dec 27, 2023
d21698f
Merge branch 'main' of github.com:Caltech-OCTO/GeometryOps.jl
skygering Dec 27, 2023
e456796
Add crosses tests
skygering Dec 28, 2023
8bd648a
Update line line processors for crosses and overlaps
skygering Dec 28, 2023
920a24b
Clean up and comment 5 of 10 functions
skygering Dec 29, 2023
380a3a9
Done with source code other than overlaps
skygering Dec 29, 2023
e88b4ca
Fixed all source code
skygering Jan 2, 2024
053e3bc
Cleanup and comment more
skygering Jan 2, 2024
c4a178f
Merge branch 'main' of github.com:Caltech-OCTO/GeometryOps.jl
skygering Jan 2, 2024
ce4f9d2
Merge branch 'main' into sg/in_on_updates
skygering Jan 2, 2024
d001c5f
Fix bugs
skygering Jan 2, 2024
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .gitignore
Original file line number Diff line number Diff line change
@@ -1,3 +1,3 @@
/Manifest.toml
/docs/build/
.vscode/
.vscode/
1 change: 0 additions & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,6 @@ Proj = "c94c279d-25a6-4763-9509-64d165bea63e"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[compat]
GeometryBasics = "0.4.7"
Proj = "1"
julia = "1.3"

Expand Down
18 changes: 11 additions & 7 deletions src/GeometryOps.jl
Original file line number Diff line number Diff line change
Expand Up @@ -23,16 +23,20 @@ include("methods/bools.jl")
include("methods/distance.jl")
include("methods/area.jl")
include("methods/centroid.jl")
include("methods/intersects.jl")
include("methods/contains.jl")
include("methods/crosses.jl")
include("methods/disjoint.jl")
include("methods/overlaps.jl")
include("methods/within.jl")
include("methods/geom_relations/contains.jl")
include("methods/geom_relations/coveredby.jl")
include("methods/geom_relations/covers.jl")
include("methods/geom_relations/crosses.jl")
include("methods/geom_relations/disjoint.jl")
include("methods/geom_relations/geom_geom_processors.jl")
include("methods/geom_relations/intersects.jl")
include("methods/geom_relations/overlaps.jl")
include("methods/geom_relations/within.jl")
include("methods/polygonize.jl")
include("methods/barycentric.jl")
include("methods/equals.jl")

include("methods/orientation.jl")
include("methods/geom_relations/touches.jl")
include("transformations/extent.jl")
include("transformations/flip.jl")
include("transformations/simplify.jl")
Expand Down
296 changes: 0 additions & 296 deletions src/methods/bools.jl
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
# # Boolean conditions

export isclockwise, isconcave
export point_on_line, point_in_polygon, point_in_ring
export line_on_line, line_in_polygon, polygon_in_polygon

# These are all adapted from Turf.jl.

Expand Down Expand Up @@ -89,297 +87,3 @@ function isconcave(poly)::Bool

return false
end

# """
# isparallel(line1::LineString, line2::LineString)::Bool

# Return `true` if each segment of `line1` is parallel to the correspondent segment of `line2`

# ## Examples
# ```julia
# import GeoInterface as GI, GeometryOps as GO
# julia> line1 = GI.LineString([(9.170356, 45.477985), (9.164434, 45.482551), (9.166644, 45.484003)])
# GeoInterface.Wrappers.LineString{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}([(9.170356, 45.477985), (9.164434, 45.482551), (9.166644, 45.484003)], nothing, nothing)

# julia> line2 = GI.LineString([(9.169356, 45.477985), (9.163434, 45.482551), (9.165644, 45.484003)])
# GeoInterface.Wrappers.LineString{false, false, Vector{Tuple{Float64, Float64}}, Nothing, Nothing}([(9.169356, 45.477985), (9.163434, 45.482551), (9.165644, 45.484003)], nothing, nothing)

# julia>
# GO.isparallel(line1, line2)
# true
# ```
# """
# function isparallel(line1, line2)::Bool
# seg1 = linesegment(line1)
# seg2 = linesegment(line2)

# for i in eachindex(seg1)
# coors2 = nothing
# coors1 = seg1[i]
# coors2 = seg2[i]
# _isparallel(coors1, coors2) == false && return false
# end
# return true
# end

# @inline function _isparallel(p1, p2)
# slope1 = bearing_to_azimuth(rhumb_bearing(GI.x(p1), GI.x(p2)))
# slope2 = bearing_to_azimuth(rhumb_bearing(GI.y(p1), GI.y(p2)))

# return slope1 === slope2
# end


"""
point_on_line(point::Point, line::LineString; ignore_end_vertices::Bool=false)::Bool

Return true if a point is on a line. Accept a optional parameter to ignore the
start and end vertices of the linestring.

## Examples

```jldoctest
import GeoInterface as GI, GeometryOps as GO

point = (1, 1)
line = GI.LineString([(0, 0), (3, 3), (4, 4)])
GO.point_on_line(point, line)

# output
true
```
"""
function point_on_line(point, line; ignore_end_vertices::Bool=false)::Bool
line_points = tuple_points(line)
n = length(line_points)

exclude_boundary = :none
for i in 1:n - 1
if ignore_end_vertices
if i === 1
exclude_boundary = :start
elseif i === n - 2
exclude_boundary = :end
elseif (i === 1 && i + 1 === n - 1)
exclude_boundary = :both
end
end
if point_on_segment(point, (line_points[i], line_points[i + 1]); exclude_boundary)
return true
end
end
return false
end

function point_on_seg(point, start, stop)
# Parse out points
x, y = GI.x(point), GI.y(point)
x1, y1 = GI.x(start), GI.y(start)
x2, y2 = GI.x(stop), GI.y(stop)
Δxl = x2 - x1
Δyl = y2 - y1
# Determine if point is on segment
cross = (x - x1) * Δyl - (y - y1) * Δxl
if cross == 0 # point is on line extending to infinity
# is line between endpoints
if abs(Δxl) >= abs(Δyl) # is line between endpoints
return Δxl > 0 ? x1 <= x <= x2 : x2 <= x <= x1
else
return Δyl > 0 ? y1 <= y <= y2 : y2 <= y <= y1
end
end
return false
end

function point_on_segment(point, (start, stop); exclude_boundary::Symbol=:none)::Bool
x, y = GI.x(point), GI.y(point)
x1, y1 = GI.x(start), GI.y(start)
x2, y2 = GI.x(stop), GI.y(stop)

dxc = x - x1
dyc = y - y1
dx1 = x2 - x1
dy1 = y2 - y1

# TODO use better predicate for crossing here
cross = dxc * dy1 - dyc * dx1
cross != 0 && return false

# Will constprop optimise these away?
if exclude_boundary === :none
if abs(dx1) >= abs(dy1)
return dx1 > 0 ? x1 <= x && x <= x2 : x2 <= x && x <= x1
end
return dy1 > 0 ? y1 <= y && y <= y2 : y2 <= y && y <= y1
elseif exclude_boundary === :start
if abs(dx1) >= abs(dy1)
return dx1 > 0 ? x1 < x && x <= x2 : x2 <= x && x < x1
end
return dy1 > 0 ? y1 < y && y <= y2 : y2 <= y && y < y1
elseif exclude_boundary === :end
if abs(dx1) >= abs(dy1)
return dx1 > 0 ? x1 <= x && x < x2 : x2 < x && x <= x1
end
return dy1 > 0 ? y1 <= y && y < y2 : y2 < y && y <= y1
elseif exclude_boundary === :both
if abs(dx1) >= abs(dy1)
return dx1 > 0 ? x1 < x && x < x2 : x2 < x && x < x1
end
return dy1 > 0 ? y1 < y && y < y2 : y2 < y && y < y1
end
return false
end

"""
point_in_polygon(point::Point, polygon::Union{Polygon, MultiPolygon}, ignore_boundary::Bool=false)::Bool

Take a Point and a Polygon and determine if the point
resides inside the polygon. The polygon can be convex or concave. The function accounts for holes.

## Examples

```jldoctest
import GeoInterface as GI, GeometryOps as GO

point = (-77.0, 44.0)
poly = GI.Polygon([[(-81, 41), (-81, 47), (-72, 47), (-72, 41), (-81, 41)]])
GO.point_in_polygon(point, poly)

# output
true
```
"""
point_in_polygon(point, polygon; kw...)::Bool =
point_in_polygon(GI.trait(point), point, GI.trait(polygon), polygon; kw...)
function point_in_polygon(
::PointTrait, point,
::PolygonTrait, poly;
ignore_boundary::Bool=false,
check_extent::Bool=false,
)::Bool
# Cheaply check that the point is inside the polygon extent
if check_extent
point_in_extent(point, GI.extent(poly)) || return false
end

# Then check the point is inside the exterior ring
point_in_polygon(
point,GI.getexterior(poly);
ignore_boundary, check_extent=false,
) || return false

# Finally make sure the point is not in any of the holes,
# flipping the boundary condition
for ring in GI.gethole(poly)
point_in_polygon(
point, ring;
ignore_boundary=!ignore_boundary,
) && return false
end
return true
end

function point_in_polygon(
::PointTrait, pt,
::Union{LineStringTrait,LinearRingTrait}, ring;
ignore_boundary::Bool=false,
check_extent::Bool=false,
)::Bool
x, y = GI.x(pt), GI.y(pt)
# Cheaply check that the point is inside the ring extent
if check_extent
point_in_extent(point, GI.extent(ring)) || return false
end
# Then check the point is inside the ring
inside = false
n = GI.npoint(ring)
p_start = GI.getpoint(ring, 1)
p_end = GI.getpoint(ring, n)
# Handle closed vs opne rings
if GI.x(p_start) == GI.x(p_end) && GI.y(p_start) == GI.y(p_end)
n -= 1
end
# Loop over all points in the ring
for i in 1:(n - 1)
# First point on edge
p_i = GI.getpoint(ring, i)
xi, yi = GI.x(p_i), GI.y(p_i)
# Second point on edge (j = i + 1)
p_j = GI.getpoint(ring, i + 1)
xj, yj = GI.x(p_j), GI.y(p_j)
# Check if point is on the ring boundary
on_boundary = ( # vertex to point has same slope as edge
yi * (xj - x) + yj * (x - xi) == y * (xj - xi) &&
(xi - x) * (xj - x) <= 0 && # x is between xi and xj
(yi - y) * (yj - y) <= 0 # y is between yi and yj
)
on_boundary && return !ignore_boundary
# Check if ray from point passes through edge
intersects = (
(yi > y) !== (yj > y) &&
(x < (xj - xi) * (y - yi) / (yj - yi) + xi)
)
if intersects
inside = !inside
end
end
return inside
end

function point_in_extent(p, extent::Extents.Extent)
(x1, x2), (y1, y1) = extent.X, extent.Y
return x1 <= GI.x(p) && y1 <= GI.y(p) && x2 >= GI.x(p) && y2 >= GI.y(p)
end

line_on_line(line1, line2) = line_on_line(trait(line1), line1, trait(line2), line2)
function line_on_line(t1::GI.AbstractCurveTrait, line1, t2::AbstractCurveTrait, line2)
for p in GI.getpoint(line1)
# FIXME: all points being on the line doesn't
# actually mean the whole line is on the line...
point_on_line(p, line2) || return false
end
return true
end

line_in_polygon(line, poly) = line_in_polygon(trait(line), line, trait(poly), poly)

function line_in_polygon(
::AbstractCurveTrait, line,
::Union{AbstractPolygonTrait,LinearRingTrait}, poly
)
Extents.intersects(GI.extent(poly), GI.extent(line)) || return false

inside = false
for i in 1:GI.npoint(line) - 1
p = GI.getpoint(line, i)
p2 = GI.getpoint(line, i + 1)
point_in_polygon(p, poly) || return false
if !inside
inside = point_in_polygon(p, poly; ignore_boundary=true)
end
# FIXME This seems like a hack, we should check for intersections rather than midpoint??
if !inside
mid = ((GI.x(p) + GI.x(p2)) / 2, (GI.y(p) + GI.y(p2)) / 2)
inside = point_in_polygon(mid, poly; ignore_boundary=true)
end
end
return inside
end

function polygon_in_polygon(poly1, poly2)
# edges1, edges2 = to_edges(poly1), to_edges(poly2)
# extent1, extent2 = to_extent(edges1), to_extent(edges2)
# Check the extents intersect
Extents.intersects(GI.extent(poly1), GI.extent(poly2)) || return false

# Check all points in poly1 are in poly2
for point in GI.getpoint(poly1)
point_in_polygon(point, poly2) || return false
end

# Check the line of poly1 does not intersect the line of poly2
#intersects(poly1, poly2) && return false

# poly1 must be in poly2
return true
end
25 changes: 0 additions & 25 deletions src/methods/contains.jl

This file was deleted.

Loading
Loading