Test Coverage for CartesianPoint: Part 2

src/ConstructiveSolidGeometry/PointsAndVectors/Vectors.jl

@@ -22,6 +22,21 @@ struct CartesianVector{T} <: AbstractCoordinateVector{T, Cartesian}
     x::T
     y::T
     z::T
+
+    CartesianVector{T}(x::T, y::T, z::T) where {T<:AbstractFloat} = new(x, y, z)
+    CartesianVector{T}(x::Real, y::Real, z::Real) where {T} = new(T(x), T(y), T(z))
+end
+
+function CartesianVector(x, y, z)
+    for (name, val) in zip((:x, :y, :z), (x, y, z))
+        (isa(val, Real) || isa(val, Unitful.Length)) ||
+            throw(ArgumentError("Expected $(name) to be a Real or Unitful.Quantity, got $(typeof(val))"))
+    end
+    x_val = to_internal_units(x)
+    y_val = to_internal_units(y)
+    z_val = to_internal_units(z)
+
+    return CartesianVector(x_val, y_val, z_val)
 end
 
 #Type promotion happens here
@@ -41,14 +56,16 @@ CartesianVector{T}(; x = 0, y = 0, z = 0) where {T} = CartesianVector{T}(T(x),T(
 
 # Need to specialize multiplication and division with units for CartesianVector, otherwise result would just be an SVector:
 Base.:(*)(v::CartesianVector{<:Real}, u::Unitful.Units{<:Any,Unitful.𝐋}) = CartesianVector(v.x * u, v.y * u, v.z * u)
-Base.:(/)(v::CartesianVector{<:Quantity{<:Real, Unitful.𝐋}}, u::Unitful.Units) = CartesianVector(v.x / u, v.y / u, v.z / u)
+#Base.:(/)(v::CartesianVector{<:Quantity{<:Real, Unitful.𝐋}}, u::Unitful.Units) = CartesianVector(v.x / u, v.y / u, v.z / u)
 
 # For user convenience, to enable constructs like `ustrip(u"mm", v)` and `NoUnits(v / u"mm")`, we'll support uconvert/ustrip
 # for CartesianVector. Unitful doesn't encourage defining those for collections, but we'll view cartesian vectors as single
 # mathematical objects in regard to units:
-Unitful.uconvert(u::Unitful.Units{<:Any,Unitful.𝐋}, v::CartesianVector{<:Quantity{<:Real, Unitful.𝐋}}) = CartesianVector(uconvert(u, v.x), uconvert(u, v.y), uconvert(u, v.z))
-Unitful.uconvert(u::Unitful.Units{<:Any,Unitful.NoDims}, v::CartesianVector{<:Quantity{<:Real, Unitful.NoDims}}) = CartesianVector(uconvert(u, v.x), uconvert(u, v.y), uconvert(u, v.z))
-Unitful.ustrip(v::CartesianVector) = CartesianVector(ustrip(v.x), ustrip(v.y), ustrip(v.z))
+
+#ToDo: Uncomment when units are supported
+#Unitful.uconvert(u::Unitful.Units{<:Any,Unitful.𝐋}, v::CartesianVector{<:Quantity{<:Real, Unitful.𝐋}}) = CartesianVector(uconvert(u, v.x), uconvert(u, v.y), uconvert(u, v.z))
+#Unitful.uconvert(u::Unitful.Units{<:Any,Unitful.NoDims}, v::CartesianVector{<:Quantity{<:Real, Unitful.NoDims}}) = CartesianVector(uconvert(u, v.x), uconvert(u, v.y), uconvert(u, v.z))
+#Unitful.ustrip(v::CartesianVector) = CartesianVector(ustrip(v.x), ustrip(v.y), ustrip(v.z))
 
 
 # @inline rotate(pt::CartesianPoint{T}, r::RotMatrix{3,T,TT}) where {T, TT} = r.mat * pt

test/ConstructiveSolidGeometry/CSG_coordinates.jl

@@ -7,7 +7,7 @@ using SolidStateDetectors.ConstructiveSolidGeometry: CartesianPoint, CartesianVe
 using StaticArrays: Size, SVector, SMatrix
 using InverseFunctions: inverse
 
-import Unitful: @u_str
+using Unitful
 
 @testset "points_and_vectors" begin
     @testset "cartesian" begin
@@ -127,3 +127,50 @@ import Unitful: @u_str
         @test_throws ArgumentError CylindricalPoint(1u"m", 2u"rad", 3u"rad")
     end
 end
+
+@testset "CartesianVector" begin
+
+    @test CartesianVector(1, 2, 3) isa CartesianVector{Float64}
+    @test CartesianVector(1, 2, 3f0) isa CartesianVector{Float32}
+    @test CartesianVector(1, 2.0f0, Float16(3)) isa CartesianVector{Float32}
+    @test CartesianVector(1.0, 2.0f0, Float16(3)) isa CartesianVector{Float64}
+    @test CartesianVector(1, 2, Float16(3)) isa CartesianVector{Float16}
+
+    @test CartesianVector(1.0u"m", 2.0u"m", 3.0f0u"m") isa CartesianVector{Float64}
+    @test CartesianVector(1.0u"m", 2.0f0u"m", 3.0f0u"m") isa CartesianVector{Float64}
+    @test CartesianVector(1.0f0u"m", 2.0f0u"m", 3.0f0u"m") isa CartesianVector{Float32}
+    @test CartesianVector(1u"m", 2u"m", 3u"m") isa CartesianVector{Float64}
+
+    @test CartesianVector(1.0u"mm", 2.0u"cm", 3.0f0u"m") isa CartesianVector{Float64}
+    @test CartesianVector(1.0f0u"mm", 2.0f0u"mm", 3.0u"m") isa CartesianVector{Float64}
+    @test CartesianVector(1.0u"mm", 2.0f0u"cm", 3.0f0u"m") isa CartesianVector{Float64}
+    @test CartesianVector(1.0f0u"mm", 2.0f0u"cm", 3.0f0u"m") isa CartesianVector{Float32}
+    @test CartesianVector(1u"mm", 2u"cm", 3u"m") isa CartesianVector{Float64}
+
+    @test CartesianVector(1u"mm", 0, 0) isa CartesianVector{Float64}
+    @test CartesianVector(1u"m", 2u"m", 3f0u"m") isa CartesianVector{Float32}
+
+    @test_throws ArgumentError CartesianVector(1u"m", 2u"rad", 3u"m")
+    @test_throws ArgumentError CartesianVector(1u"s", 2u"m", 3u"m")
+    @test_throws ArgumentError CartesianVector(1u"m", 2u"m", 3u"kg")
+
+    v1 = CartesianVector(1, 2, 3)
+    @test v1 == CartesianVector(1.0, 2.0, 3.0)
+
+    v2 = CartesianVector{Float32}(x=1, y=2, z=3)
+    @test v2 == CartesianVector(1f0, 2f0, 3f0)
+
+    z = zero(CartesianVector{Float64})
+    @test z == CartesianVector(0.0, 0.0, 0.0)
+
+    v = CartesianVector(1.0, 2.0, 3.0)
+    v_mul = v * u"m"
+    @test v_mul == CartesianVector(1.0u"m", 2.0u"m", 3.0u"m")
+
+    v = CartesianVector(1f0, 2f0, 3f0)
+    @test -v == CartesianVector(-1f0, -2f0, -3f0)
+
+    cz1 = CartesianZero{Float32}()
+    cz2 = CartesianZero{Float64}()
+    @test -(cz1, cz2) == CartesianVector(0.0, 0.0, 0.0)
+end

test/ConstructiveSolidGeometry/CSG_primitives.jl

@@ -3,6 +3,7 @@ using SolidStateDetectors
 using StaticArrays
 using LinearAlgebra
 using Rotations
+using Unitful
 
 import SolidStateDetectors.ConstructiveSolidGeometry as CSG
 import SolidStateDetectors.ConstructiveSolidGeometry: Geometry, Dictionary
@@ -91,6 +92,99 @@ no_translations = (rotation = one(SMatrix{3, 3, T, 9}), translation = zero(Carte
         tuple_cone = @inferred CSG.Cone{Float64}(φ=(π/4,3*π/4), rotation=SMatrix{3}(0,0,-1,0,1,0,1,0,0))
         rot_cone = @inferred CSG.Cone{Float64}(φ=π/2, rotation=SMatrix{3}(0,0,-1,0,1,0,1,0,0) * SMatrix{3}(cos(π/4),sin(π/4),0,-sin(π/4),cos(π/4),0,0,0,1))
         @test tuple_cone ==rot_cone
+
+        # PartialCylinder (OpenPrimitive)
+        pc_open = CSG.Cone(CSG.OpenPrimitive, r = convert(T, 1.0), φ = convert(T, π/2), hZ = convert(T, 1.0))
+        @test CSG._in(CartesianPoint{T}(0.5, 0.2, 0.0), pc_open)
+        @test !CSG._in(CartesianPoint{T}(1.1, 0.0, 0.0), pc_open)
+        @test !CSG._in(CartesianPoint{T}(-0.5, -0.2, 0.0), pc_open)
+        @test !CSG._in(CartesianPoint{T}(0.5, 0.0, 1.0), pc_open)
+
+        # VaryingCylinder (ClosedPrimitive)
+        vc_closed = CSG.Cone{T}(r = ((convert(T,1.0),), (convert(T,3.0),)), φ = nothing, hZ = convert(T,2.0))
+        @test CSG._in(CartesianPoint{T}(1.9, 0.0, 0.0), vc_closed)
+        @test !CSG._in(CartesianPoint{T}(2.1, 0.0, 0.0), vc_closed)
+        @test CSG._in(CartesianPoint{T}(1.0, 0.0, -2.0), vc_closed)
+        @test CSG._in(CartesianPoint{T}(3.0, 0.0, 2.0), vc_closed)
+
+        # VaryingCylinder (OpenPrimitive)
+        vc_open = CSG.Cone(CSG.OpenPrimitive, r = ((convert(T,1.0),), (convert(T,3.0),)), φ = nothing, hZ = convert(T,2.0))
+        @test !CSG._in(CartesianPoint{T}(2.0, 0.0, 0.0), vc_open)
+        @test CSG._in(CartesianPoint{T}(1.9, 0.0, 0.0), vc_open)
+        @test !CSG._in(CartesianPoint{T}(3.0, 0.0, 2.0), vc_open)
+
+        # PartialVaryingCylinder (ClosedPrimitive)
+        pvc_closed = CSG.Cone{T}(r = ((convert(T,1.0),), (convert(T,2.0),)), φ = convert(T, π), hZ = convert(T, 1.0))
+        @test CSG._in(CartesianPoint{T}(1.4, 0.1, 0.0), pvc_closed)
+        @test !CSG._in(CartesianPoint{T}(1.6, 0.0, 0.0), pvc_closed)
+        @test !CSG._in(CartesianPoint{T}(1.0, -0.1, 0.0), pvc_closed)
+
+        # PartialVaryingCylinder (OpenPrimitive)
+        pvc_open = CSG.Cone(CSG.OpenPrimitive, r = ((convert(T,1.0),), (convert(T,2.0),)), φ = convert(T, π), hZ = convert(T, 1.0))
+        @test !CSG._in(CartesianPoint{T}(1.5, 0.0, 0.0), pvc_open)
+        @test CSG._in(CartesianPoint{T}(1.4, 0.1, 0.0), pvc_open)
+
+        # PartialVaryingTube (OpenPrimitive)
+        pvt_open = CSG.Cone(CSG.OpenPrimitive, r = ((convert(T,1.0), convert(T,2.0)), (convert(T,1.5), convert(T,2.5))), φ = convert(T, π/2), hZ = convert(T, 1.0))
+        @test CSG._in(CartesianPoint{T}(1.6, 0.2, 0.0), pvt_open)
+        @test !CSG._in(CartesianPoint{T}(1.1, 0.0, 0.0), pvt_open)
+        @test !CSG._in(CartesianPoint{T}(2.6, 0.0, 0.0), pvt_open)
+        @test !CSG._in(CartesianPoint{T}(1.6, -0.2, 0.0), pvt_open)
+
+        # Upward cones
+        r_up_closed = ((nothing, convert(T,0.0)), nothing)
+        r_up_open   = ((nothing, convert(T,0.1)), nothing)
+        φ_up = nothing
+        h_up = convert(T, 1.0)
+        origin = CartesianPoint{T}(0.0, 0.0, 0.0)
+        rotation = SMatrix{3,3,T}(I)
+
+        # Closed UpwardCone
+        upward_closed = CSG.UpwardCone{T, CSG.ClosedPrimitive}(r_up_closed, φ_up, h_up, origin, rotation)
+        @test CSG._in(CartesianPoint{T}(0.0, 0.0, 0.5), upward_closed)
+        @test !CSG._in(CartesianPoint{T}(0.0, 0.0, 1.1), upward_closed)
+        @test !CSG._in(CartesianPoint{T}(2.0, 0.0, 0.5), upward_closed)
+        @test length(CSG.surfaces(upward_closed)) > 0
+
+        # Open UpwardCone
+        upward_open = CSG.UpwardCone{T, CSG.OpenPrimitive}(r_up_open, φ_up, h_up, origin, rotation)
+        @test CSG._in(CartesianPoint{T}(0.0, 0.0, 0.5), upward_open)
+        @test !CSG._in(CartesianPoint{T}(0.0, 0.0, 1.0), upward_open)
+        @test !CSG._in(CartesianPoint{T}(2.0, 0.0, 0.5), upward_open)
+        @test length(CSG.surfaces(upward_open)) > 0
+
+        # PartialUpwardCone
+        φ_puc = convert(T, π)
+        h_puc = convert(T, 1.0)
+        r_puc_closed = ((nothing, convert(T,0.1)), nothing)
+        r_puc_open   = ((nothing, convert(T,0.5)), nothing)
+
+        puc_closed = CSG.PartialUpwardCone{T, CSG.ClosedPrimitive}(r_puc_closed, φ_puc, h_puc, origin, rotation)
+        puc_open   = CSG.PartialUpwardCone{T, CSG.OpenPrimitive}(r_puc_open, φ_puc, h_puc, origin, rotation)
+
+        z_test = convert(T, 0.1)
+        r_max_closed = CSG.radius_at_z(puc_closed.hZ, puc_closed.r[1][2], zero(T), z_test)
+        r_max_open   = CSG.radius_at_z(puc_open.hZ,   puc_open.r[1][2], zero(T), z_test)
+
+        # Points inside the cone
+        r_inside_closed = r_max_closed * 0.5
+        r_inside_open   = r_max_open   * 0.5
+
+        # ---- Closed cone tests ----
+        @test CSG._in(CartesianPoint{T}(r_inside_closed, 0.0, z_test), puc_closed)
+        @test !CSG._in(CartesianPoint{T}(r_max_closed * 1.1, 0.0, z_test), puc_closed)
+        @test !CSG._in(CartesianPoint{T}(r_inside_closed, -0.1, z_test), puc_closed)
+
+        # ---- Open cone tests ----
+        x_inside_open = r_inside_open * cos(T(0.1))  # small angle inside φ range
+        y_inside_open = r_inside_open * sin(T(0.1))
+        @test CSG._in(CartesianPoint{T}(x_inside_open, y_inside_open, z_test), puc_open)
+        @test !CSG._in(CartesianPoint{T}(r_max_open * 1.1, 0.0, z_test), puc_open)
+        @test !CSG._in(CartesianPoint{T}(r_inside_open/2, -0.1, z_test), puc_open)
+
+        surfs = CSG.surfaces(puc_open)
+        @test length(surfs) > 0
+        @test all(x -> x !== nothing, surfs)
     end
     @testset "Torus" begin
         for r_tube in (2.0, Dict("from" => 1.0, "to" => 2.0)), 
@@ -276,6 +370,15 @@ no_translations = (rotation = one(SMatrix{3, 3, T, 9}), translation = zero(Carte
         @test !in(CartesianPoint{Float64}(1,1,3),box_open_trafo)
         @test !in(CartesianPoint{Float64}(1,1,3+tol),box_closed_trafo)
         @test in(CartesianPoint{Float64}(1,1,3-tol),box_open_trafo)
+
+        # Test Volume Primitives
+        pt_outside = CartesianPoint{Float32}(3, 0, 0)
+        pt_inside  = CartesianPoint{Float32}(0.5, 0, 0)
+        # Expected distances
+        expected_distance_outside = 2.0f0  # nearest face along x
+        expected_distance_inside  = 0.5f0  # nearest face along x
+        @test CSG.distance(pt_outside, box2) ≈ expected_distance_outside
+        @test CSG.distance(pt_inside, box2) ≈ expected_distance_inside
     end
     @testset "RegularPrism" begin
         prism1 = @inferred CSG.RegularPrism{3}(CSG.ClosedPrimitive,r=1f0, hZ = 2f0)
@@ -401,9 +504,9 @@ no_translations = (rotation = one(SMatrix{3, 3, T, 9}), translation = zero(Carte
         cm_out = CSG.ConeMantle{Float32}( :outwards; r=(1f0, 2f0), φ=nothing, hZ=1f0 )
         cm_in  = CSG.flip(cm_out)
 
-        pt_top    = CartesianPoint{Float32}(0, 1.5f0, 1f0)      # near top radius
-        pt_bottom = CartesianPoint{Float32}(0, 0.5f0, -1f0)     # near bottom radius
-        pt_mid    = CartesianPoint{Float32}(1f0, 0, 0f0)        # mid cone height
+        pt_top    = CartesianPoint{Float32}(0, 1.5f0, 1.0f0)      # near top radius
+        pt_bottom = CartesianPoint{Float32}(0, 0.5f0, -1.0f0)     # near bottom radius
+        pt_mid    = CartesianPoint{Float32}(1.0f0, 0, 0.0f0)        # mid cone height
 
         normal_out = CSG.normal(cm_out, pt_mid)
         normal_in  = CSG.normal(cm_in, pt_mid)
@@ -439,9 +542,9 @@ no_translations = (rotation = one(SMatrix{3, 3, T, 9}), translation = zero(Carte
         end
     end
     @testset "Ellipse" begin
-        ell1 = @inferred CSG.Ellipse{Float32}(r = 1f0, φ=10f0)
+        ell1 = @inferred CSG.Ellipse{Float32}(r = 1.0f0, φ=10.0f0)
         @inferred CSG.Ellipse{T}(r = (1,2))
-        ell2 = @inferred CSG.Ellipse(r = 1f0, φ=10f0) 
+        ell2 = @inferred CSG.Ellipse(r = 1.0f0, φ=10.0f0) 
         @test ell1 == ell2
     end
     @testset "Line" begin
@@ -614,3 +717,71 @@ end
     c = sim.detector.contacts[2].geometry
     @test Geometry(T, Dictionary(c), default_units, no_translations) == c
 end
+
+@testset "Test geometry rounding" begin
+    pt = CartesianPoint{Float64}(1.23456789, 0.0, 1e-13)
+    rounded_pt = CSG.geom_round(pt)
+
+    @test rounded_pt.x == 1.23456789
+    @test rounded_pt.y == 0.0
+    @test rounded_pt.z == 0.0
+
+    pt32 = CartesianPoint{Float32}(1.2345678f0, Float32(1e-7), Float32(-1e-7))
+    rounded_pt32 = CSG.geom_round(pt32)
+
+    @test rounded_pt32.x ≈ 1.234568f0
+    @test rounded_pt32.y == 0.0f0
+    @test rounded_pt32.z == 0.0f0
+
+    pts = [CartesianPoint{Float64}(1e-13, 2.3456789, 0.0),
+           CartesianPoint{Float64}(0.123456789, 1e-14, -0.987654321)]
+    rounded_pts = CSG.geom_round(pts)
+
+    @test rounded_pts[1].x == 0.0
+    @test rounded_pts[1].y ≈ 2.3456789
+    @test rounded_pts[1].z == 0.0
+    @test rounded_pts[2].x ≈ 0.123456789
+    @test rounded_pts[2].y == 0.0
+    @test rounded_pts[2].z ≈ -0.987654321
+end
+
+@testset "Geometry Test" begin
+    outer_transformations = (
+        rotation    = SMatrix{3,3,T}(I),
+        translation = CartesianVector{T}(1.0, 2.0, 3.0)
+    )
+    dict = Dict(
+        "Z" => π/2,
+        "box" => Dict(
+            "hX"        => 1.0,
+            "hY"        => 2.0,
+            "hZ"        => 3.0,
+            "translation" => CartesianVector{T}(4.0, 5.0, 6.0)
+        )
+    )
+    input_units = (angle  = u"rad", length = u"m")
+    geom = CSG.Geometry(T, CSG.Rotation, dict, input_units, outer_transformations)
+
+    primitive_translation_rotated = geom.rotation * dict["box"]["translation"]
+
+    expected_origin = CartesianPoint(
+        outer_transformations.translation[1] + primitive_translation_rotated[1],
+        outer_transformations.translation[2] + primitive_translation_rotated[2],
+        outer_transformations.translation[3] + primitive_translation_rotated[3]
+    )
+    @test geom.origin == expected_origin
+
+    expected_rot = SMatrix{3,3,T,9}(RotZ(π/2))
+    @test geom.rotation ≈ expected_rot
+
+    # --- Test + operator for AbstractConstructiveGeometry ---
+    v = CartesianVector{T}(1.0, -2.0, 0.5)
+    geom_translated = geom + v
+    @test geom_translated.origin == geom.origin + v
+    @test geom_translated.rotation ≈ geom.rotation
+    # --- Test + operator for AbstractPrimitive ---
+    prim = CSG.Box(CSG.ClosedPrimitive, hX=1.0, hY=2.0, hZ=3.0)
+    prim_translated = prim + v
+    @test prim_translated.origin == prim.origin + v
+    @test prim_translated.rotation ≈ prim.rotation
+end

test/test_charge_drift_models.jl

@@ -635,6 +635,9 @@ Base.in(::CartesianPoint{T}, ::OutsideTestVolume{T}) where {T} = false
     dead_vol = SolidStateDetectors.DeadVolume{T, typeof(geom)}("dead", geom)
     result = SolidStateDetectors.modulate_driftvector(sv, pt, dead_vol)
     @test result == CartesianVector{T}(0,0,0)
+
+    result = SolidStateDetectors.modulate_driftvector(sv, pt, missing)
+    @test result == sv
 end
 
 @timed_testset "NBodyChargeCloud Units" begin

test/test_real_detectors.jl

@@ -119,6 +119,23 @@ end
     @test idx_x == SolidStateDetectors.searchsortednearest(ticks, pt.x)
     @test idx_y == SolidStateDetectors.searchsortednearest(ticks, pt.y)
     @test idx_z == SolidStateDetectors.searchsortednearest(ticks, pt.z)
+
+    # Additional tests: find_closest_gridpoint
+    # CartesianPoint version
+    idxs_cart = SolidStateDetectors.find_closest_gridpoint(pt, grid)
+    @test idxs_cart[1] == SolidStateDetectors.searchsortednearest(ticks, pt.x)
+    @test idxs_cart[2] == SolidStateDetectors.searchsortednearest(ticks, pt.y)
+    @test idxs_cart[3] == SolidStateDetectors.searchsortednearest(ticks, pt.z)
+
+    # CylindricalPoint version
+    pt_cyl = CylindricalPoint{T}(0.33f0, π/3, 0.12f0)
+    idxs_cyl = SolidStateDetectors.find_closest_gridpoint(pt_cyl, grid)
+
+    # Convert Cylindrical to Cartesian for comparison
+    pt_cyl_cart = CartesianPoint(pt_cyl)
+    @test idxs_cyl[1] == SolidStateDetectors.searchsortednearest(ticks, pt_cyl_cart.x)
+    @test idxs_cyl[2] == SolidStateDetectors.searchsortednearest(ticks, pt_cyl_cart.y)
+    @test idxs_cyl[3] == SolidStateDetectors.searchsortednearest(ticks, pt_cyl_cart.z)
 end
 @timed_testset "HexagonalPrism" begin
     sim = Simulation{T}(SSD_examples[:Hexagon])