Skip to content

Commit 3e2c3dd

Browse files
committed
Move intersectparameters utility to utils/intersect.jl
1 parent 7b015b2 commit 3e2c3dd

File tree

2 files changed

+52
-52
lines changed

2 files changed

+52
-52
lines changed

src/utils/intersect.jl

Lines changed: 52 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,58 @@
22
# Licensed under the MIT License. See LICENSE in the project root.
33
# ------------------------------------------------------------------
44

5+
"""
6+
intersectparameters(a, b, c, d)
7+
8+
Compute the parameters `λ₁` and `λ₂` of the lines
9+
`a + λ₁ ⋅ v⃗₁`, with `v⃗₁ = b - a` and
10+
`c + λ₂ ⋅ v⃗₂`, with `v⃗₂ = d - c` spanned by the input
11+
points `a`, `b` resp. `c`, `d` such that to yield line
12+
points with minimal distance or the intersection point
13+
(if lines intersect).
14+
15+
Furthermore, the ranks `r` of the matrix of the linear
16+
system `A ⋅ λ⃗ = y⃗`, with `A = [v⃗₁ -v⃗₂], y⃗ = c - a`
17+
and the rank `rₐ` of the augmented matrix `[A y⃗]` are
18+
calculated in order to identify the intersection type:
19+
20+
- Intersection: r == rₐ == 2
21+
- Collinear: r == rₐ == 1
22+
- No intersection: r != rₐ
23+
- No intersection and parallel: r == 1, rₐ == 2
24+
- No intersection, skew lines: r == 2, rₐ == 3
25+
"""
26+
function intersectparameters(a::Point, b::Point, c::Point, d::Point)
27+
# augmented linear system
28+
A = ustrip.([(b - a) (c - d) (c - a)])
29+
30+
# normalize by maximum absolute coordinate
31+
A = A / maximum(abs, A)
32+
33+
# numerical tolerance
34+
T = eltype(A)
35+
τ = atol(T)
36+
37+
# check if a vector is non zero
38+
isnonzero(v) = !isapprox(v, zero(v), atol=τ)
39+
40+
# calculate ranks by checking the zero rows of
41+
# the factor R in the QR matrix factorization
42+
_, R = qr(A)
43+
r = sum(isnonzero, eachrow(R[:, SVector(1, 2)]))
44+
rₐ = sum(isnonzero, eachrow(R))
45+
46+
# calculate parameters of intersection
47+
if r 2
48+
λ = A[:, SVector(1, 2)] \ A[:, 3]
49+
λ₁, λ₂ = λ[1], λ[2]
50+
else # parallel or collinear
51+
λ₁, λ₂ = zero(T), zero(T)
52+
end
53+
54+
λ₁, λ₂, r, rₐ
55+
end
56+
557
"""
658
pairwiseintersect(segments; [digits])
759

src/utils/misc.jl

Lines changed: 0 additions & 52 deletions
Original file line numberDiff line numberDiff line change
@@ -55,55 +55,3 @@ function svdbasis(p::AbstractVector{<:Point})
5555
n = Vec(zero(ℒ), zero(ℒ), oneunit(ℒ))
5656
isnegative((u × v) n) ? (v, u) : (u, v)
5757
end
58-
59-
"""
60-
intersectparameters(a, b, c, d)
61-
62-
Compute the parameters `λ₁` and `λ₂` of the lines
63-
`a + λ₁ ⋅ v⃗₁`, with `v⃗₁ = b - a` and
64-
`c + λ₂ ⋅ v⃗₂`, with `v⃗₂ = d - c` spanned by the input
65-
points `a`, `b` resp. `c`, `d` such that to yield line
66-
points with minimal distance or the intersection point
67-
(if lines intersect).
68-
69-
Furthermore, the ranks `r` of the matrix of the linear
70-
system `A ⋅ λ⃗ = y⃗`, with `A = [v⃗₁ -v⃗₂], y⃗ = c - a`
71-
and the rank `rₐ` of the augmented matrix `[A y⃗]` are
72-
calculated in order to identify the intersection type:
73-
74-
- Intersection: r == rₐ == 2
75-
- Collinear: r == rₐ == 1
76-
- No intersection: r != rₐ
77-
- No intersection and parallel: r == 1, rₐ == 2
78-
- No intersection, skew lines: r == 2, rₐ == 3
79-
"""
80-
function intersectparameters(a::Point, b::Point, c::Point, d::Point)
81-
# augmented linear system
82-
A = ustrip.([(b - a) (c - d) (c - a)])
83-
84-
# normalize by maximum absolute coordinate
85-
A = A / maximum(abs, A)
86-
87-
# numerical tolerance
88-
T = eltype(A)
89-
τ = atol(T)
90-
91-
# check if a vector is non zero
92-
isnonzero(v) = !isapprox(v, zero(v), atol=τ)
93-
94-
# calculate ranks by checking the zero rows of
95-
# the factor R in the QR matrix factorization
96-
_, R = qr(A)
97-
r = sum(isnonzero, eachrow(R[:, SVector(1, 2)]))
98-
rₐ = sum(isnonzero, eachrow(R))
99-
100-
# calculate parameters of intersection
101-
if r 2
102-
λ = A[:, SVector(1, 2)] \ A[:, 3]
103-
λ₁, λ₂ = λ[1], λ[2]
104-
else # parallel or collinear
105-
λ₁, λ₂ = zero(T), zero(T)
106-
end
107-
108-
λ₁, λ₂, r, rₐ
109-
end

0 commit comments

Comments
 (0)