DomainError in intersection(ray, triangle)

Question

DomainError in intersection(ray, triangle)

Closed this issue 6 months ago · 5 comments

I am trying to calculate if triangles intersect with rays and the intersection returns a DomainError when it tries to return with a Touching ray.
MWE:

# Inputs

julia> rt = Ray(Point(-18.51314216447587, -45.79111127739282, 81.42486924309667), Vec(24.856095593894622, 1.563812983514694, -2.174626583982459))
Ray{3,Float64}
├─ p: Point(-18.51314216447587, -45.79111127739282, 81.42486924309667)
└─ v: Vec(24.856095593894622, 1.563812983514694, -2.174626583982459)

julia> rr = Triangle(Point(-18.25, -46.0, 84.282382318787), Point(-18.625, -46.0, 79.9961127052515), Point(-18.66442649342761, -45.373333832178474, 79.9961127052515))
Triangle{3,Float64}
├─ Point(-18.25, -46.0, 84.282382318787)
├─ Point(-18.625, -46.0, 79.9961127052515)
└─ Point(-18.66442649342761, -45.373333832178474, 79.9961127052515)

# the point of the ray is on the center of triangle
julia> mean(coordinates.(vertices(rr)))
Vec(-18.51314216447587, -45.79111127739282, 81.42486924309667)

julia> intersection(rt, rr)
ERROR: DomainError with -2.2956335921835282e-17:
r(t) is not defined for t < 0.
Stacktrace:
 [1] (::Ray{3, Float64})(t::Float64)
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25
 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0
 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101
 [4] top-level scope
   @ REPL[36]:1

I am assuming that this should be Touching intersection, and using the debugger, I can actually confirm it.

I have no idea how to resolve this. I guess this is a numerical error, because the following works fine:

julia> t2 = Triangle(Point(0,0,0), Point(0,1,0), Point(1,1,0))
Triangle{3,Float64}
├─ Point(0.0, 0.0, 0.0)
├─ Point(0.0, 1.0, 0.0)
└─ Point(1.0, 1.0, 0.0)

julia> r2 = Ray(Point(0.4,0.8,0), Vec3(0,0,1))
Ray{3,Float64}
├─ p: Point(0.4, 0.8, 0.0)
└─ v: Vec(0.0, 0.0, 1.0)

julia> intersection(t2, r2)
Intersection{Point3}(Touching, Point(0.4, 0.8, 0.0))

Answer 1 · 2024-03-18T10:31:39.000Z

Cab you confirm that this issue also happens if you scale the coordinates? Try to send the Ray and Triangle to the Scale transform before calling the intersection. Em seg., 18 de mar. de 2024 05:16, cserteGT3 ***@***.***> escreveu:

…

I am trying to calculate if triangles intersect with rays and the intersection returns a DomainError when it tries to return with a Touching ray. MWE: # Inputs julia> rt = Ray(Point(-18.51314216447587, -45.79111127739282, 81.42486924309667), Vec(24.856095593894622, 1.563812983514694, -2.174626583982459)) Ray{3,Float64} ├─ p: Point(-18.51314216447587, -45.79111127739282, 81.42486924309667) └─ v: Vec(24.856095593894622, 1.563812983514694, -2.174626583982459) julia> rr = Triangle(Point(-18.25, -46.0, 84.282382318787), Point(-18.625, -46.0, 79.9961127052515), Point(-18.66442649342761, -45.373333832178474, 79.9961127052515)) Triangle{3,Float64} ├─ Point(-18.25, -46.0, 84.282382318787) ├─ Point(-18.625, -46.0, 79.9961127052515) └─ Point(-18.66442649342761, -45.373333832178474, 79.9961127052515) # the point of the ray is on the center of triangle julia> mean(coordinates.(vertices(rr)))Vec(-18.51314216447587, -45.79111127739282, 81.42486924309667) julia> intersection(rt, rr) ERROR: DomainError with -2.2956335921835282e-17:r(t) is not defined for t < 0. Stacktrace: [1] (::Ray{3, Float64})(t::Float64) @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64}) @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64}) @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101 [4] top-level scope @ REPL[36]:1 I am assuming that this should be Touching intersection, and using the debugger, I can actually confirm it. I have no idea how to resolve this. I guess this is a numerical error, because the following works fine: julia> t2 = Triangle(Point(0,0,0), Point(0,1,0), Point(1,1,0)) Triangle{3,Float64} ├─ Point(0.0, 0.0, 0.0) ├─ Point(0.0, 1.0, 0.0) └─ Point(1.0, 1.0, 0.0) julia> r2 = Ray(Point(0.4,0.8,0), Vec3(0,0,1)) Ray{3,Float64} ├─ p: Point(0.4, 0.8, 0.0) └─ v: Vec(0.0, 0.0, 1.0) julia> intersection(t2, r2)Intersection{Point3}(Touching, Point(0.4, 0.8, 0.0)) — Reply to this email directly, view it on GitHub <#805>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAZQW3PHP36YBW53U6BSQ2DYY2PG7AVCNFSM6AAAAABE3CYWT6VHI2DSMVQWIX3LMV43ASLTON2WKOZSGE4TCNJQHE4TENA> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

Answer 2 · 2024-03-18T11:14:01.000Z

I'm not sure what to expect, but here they are:

julia> rt = Ray(Point(-18.51314216447587, -45.79111127739282, 81.42486924309667), Vec(24.856095593894622, 1.563812983514694, -2.174626583982459))
Ray{3,Float64}
├─ p: Point(-18.51314216447587, -45.79111127739282, 81.42486924309667)
└─ v: Vec(24.856095593894622, 1.563812983514694, -2.174626583982459)

julia> rr = Triangle(Point(-18.25, -46.0, 84.282382318787), Point(-18.625, -46.0, 79.9961127052515), Point(-18.66442649342761, -45.373333832178474, 79.9961127052515))
Triangle{3,Float64}
├─ Point(-18.25, -46.0, 84.282382318787)
├─ Point(-18.625, -46.0, 79.9961127052515)
└─ Point(-18.66442649342761, -45.373333832178474, 79.9961127052515)

julia> intersection(rt, rr)
ERROR: DomainError with -2.2956335921835282e-17:
r(t) is not defined for t < 0.
Stacktrace:
 [1] (::Ray{3, Float64})(t::Float64)
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25
 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0
 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101
 [4] top-level scope
   @ REPL[25]:1

julia> intersection(Scale(1.0)(rt), Scale(1.0)(rr))
ERROR: DomainError with -2.2956335921835282e-17:
r(t) is not defined for t < 0.
Stacktrace:
 [1] (::Ray{3, Float64})(t::Float64)
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25
 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0
 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101
 [4] top-level scope
   @ REPL[26]:1

julia> intersection(Scale(1)(rt), Scale(1)(rr))
ERROR: DomainError with -2.2956335921835282e-17:
r(t) is not defined for t < 0.
Stacktrace:
 [1] (::Ray{3, Float64})(t::Float64)
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25
 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0
 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101
 [4] top-level scope
   @ REPL[27]:1

julia> intersection(Scale(1.0001)(rt), Scale(1.0001)(rr))
Intersection{Point3}(Touching, Point(-18.51499347869232, -45.79569038852056, 81.43301173002098))

julia> intersection(Scale(0.5)(rt), Scale(0.5)(rr))
ERROR: DomainError with -2.2956335921835282e-17:
r(t) is not defined for t < 0.
Stacktrace:
 [1] (::Ray{3, Float64})(t::Float64)
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\primitives\ray.jl:25
 [2] intersection(f::Function, ray::Ray{3, Float64}, tri::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections\rays.jl:0
 [3] intersection(g₁::Ray{3, Float64}, g₂::Triangle{3, Float64})
   @ Meshes C:\Users\csert\.julia\dev\Meshes\src\intersections.jl:101
 [4] top-level scope
   @ REPL[29]:1

Answer 3 · 2024-03-18T11:30:09.000Z

So these algorithms that depend on numerical tolerances usually perform well when coordinates are within [0,1] or [-1,1], etc, i.e., normalized coordinates. If it is an option, you can scale the whole domain by the boundingbox:

b = boundingbox(mesh)
s = maximum(sides(b))
smesh = mesh |> Scale(1 / s)

# intersections should work fine
intersection(..., ...)

We can also consider the normalization inside the intersection method, but it may too computationally expensive to normalize on every call. We can add tips like this to the documentation in case people encounter errors. What do you think?

Answer 4 · 2024-03-18T15:16:42.000Z

I believe that adding this normalization as a preprocessing step is a reasonable approach. There are transforms for that already implemented, we need to document it better though. Feel free to open a PR with improved documentation.

Answer 5 · 2024-03-20T09:49:35.000Z

Thank you for making this clear. I agree that normallization should not be included by default and it's easy to add as a previous step.