New parameterization functions for Cone and ConeSurface (#1020)

* Test change to radius from cos to tan

* Drop unneeded type spec

* Add new param functions for Cone and ConeSurface

* Bugfix - missed name change

* Bugfix - updated changed name

* Bugfix in point-vector-svector operations

* Adjust formatting to match prior style

* Add tests for ConeSurface parameterization

* Bugfix in test - wrong name

* Add DomainError test for Cone param

* Bugfix - wrong name in test

* Bugfix - wrong name in test

* Bugfix - promote tuple to Point for comparison

* Add tests for Cone param function

* Apply suggestions from code review

Co-authored-by: Júlio Hoffimann <[email protected]>

* Update test/primitives.jl

* Update test/primitives.jl

* Rename arguments and internal names, use some abstract API

* Update variable name

* Rename arguments and internal names, use some abstract APIs

* Bugfix - update name

* Conversion to abstract methods

* Revert test changes to simplify merge

* Re-add tests after merge

* Final adjustments

---------

Co-authored-by: Júlio Hoffimann <[email protected]>

Author: mikeingold

src/boundingboxes.jl

@@ -55,7 +55,7 @@ end
 
 function boundingbox(c::ConeSurface)
   us = (0, 1 / 4, 1 / 2, 3 / 4)
-  vs = (1,)
+  vs = (0,)
   ps = [c(u, v) for (u, v) in Iterators.product(us, vs)]
   boundingbox([ps; apex(c)])
 end

src/geometries/primitives/cone.jl

@@ -34,3 +34,12 @@ halfangle(c::Cone) = atan(radius(base(c)), height(c))
 
 Base.isapprox(c₁::Cone, c₂::Cone; atol=atol(lentype(c₁)), kwargs...) =
   isapprox(boundary(c₁), boundary(c₂); atol, kwargs...)
+
+function (c::Cone)(r, φ, h)
+  if (r < 0 || r > 1) || (φ < 0 || φ > 1) || (h < 0 || h > 1)
+    throw(DomainError((r, φ, h), "c(r, φ, h) is not defined for r, φ, h outside [0, 1]³."))
+  end
+  a = c.apex
+  b = c.base(r, φ)
+  Segment(b, a)(h)
+end

src/geometries/primitives/conesurface.jl

@@ -36,12 +36,7 @@ function (c::ConeSurface)(φ, h)
   if (φ < 0 || φ > 1) || (h < 0 || h > 1)
     throw(DomainError((φ, h), "c(φ, h) is not defined for φ, h outside [0, 1]²."))
   end
-  n = -normal(c.base)
-  v = c.base(T(0), T(0)) - c.apex
-  l = norm(v)
-  θ = ∠(n, v)
-  o = c.apex + T(h) * v
-  r = T(h) * l * cos(θ)
-  s = Circle(Plane(o, n), r)
-  s(T(φ))
+  a = c.apex
+  b = c.base(one(T), φ)
+  Segment(b, a)(h)
 end

src/sampling/regular.jl

@@ -75,9 +75,9 @@ firstoffset(::CylinderSurface) = (n -> zero(n), n -> zero(n))
 lastoffset(::CylinderSurface) = (n -> inv(n), n -> zero(n))
 extrapoints(c::CylinderSurface, sz) = (bottom(c)(0, 0), top(c)(0, 0))
 
-firstoffset(::ConeSurface) = (n -> zero(n), n -> inv(n))
-lastoffset(::ConeSurface) = (n -> inv(n), n -> zero(n))
-extrapoints(c::ConeSurface, sz) = (apex(c), base(c)(0, 0))
+firstoffset(::ConeSurface) = (n -> zero(n), n -> zero(n))
+lastoffset(::ConeSurface) = (n -> inv(n), n -> inv(n))
+extrapoints(c::ConeSurface, sz) = (base(c)(0, 0), apex(c))
 
 firstoffset(::FrustumSurface) = (n -> zero(n), n -> zero(n))
 lastoffset(::FrustumSurface) = (n -> inv(n), n -> zero(n))

test/primitives.jl

@@ -1046,15 +1046,19 @@ end
   @test crs(c) <: Cartesian{NoDatum}
   @test Meshes.lentype(c) == ℳ
   @test boundary(c) == ConeSurface(d, a)
+  @test_throws DomainError c(T(0), T(0), nextfloat(T(1)))
 
   p = Plane(cart(0, 0, 0), vector(0, 0, 1))
   d = Disk(p, T(2))
-  a = T.((0, 0, 1))
+  a = cart(0, 0, 1)
   c = Cone(d, a)
   @test embeddim(c) == 3
   @test paramdim(c) == 3
   @test crs(c) <: Cartesian{NoDatum}
   @test Meshes.lentype(c) == ℳ
+  @test c(T(0), T(0), T(0)) ≈ centroid(p)
+  @test c(T(0), T(0), T(1)) ≈ a
+  @test c(T(1.0), T(0.25), T(0.0)) ≈ cart(0, 2, 0)
 
   p = Plane(cart(0, 0, 0), vector(0, 0, 1))
   d = Disk(p, T(2))
@@ -1120,15 +1124,19 @@ end
   @test crs(s) <: Cartesian{NoDatum}
   @test Meshes.lentype(s) == ℳ
   @test isnothing(boundary(s))
+  @test_throws DomainError s(T(0), nextfloat(T(1)))
 
   p = Plane(cart(0, 0, 0), vector(0, 0, 1))
   d = Disk(p, T(2))
-  a = T.((0, 0, 1))
+  a = cart(0, 0, 1)
   c = ConeSurface(d, a)
   @test embeddim(c) == 3
   @test paramdim(c) == 2
   @test crs(c) <: Cartesian{NoDatum}
   @test Meshes.lentype(c) == ℳ
+  @test c(T(0), T(0)) ≈ cart(2, 0, 0)
+  @test c(T(0), T(1)) ≈ a
+  @test c(T(0.25), T(0)) ≈ cart(0, 2, 0)
 
   p = Plane(cart(0, 0, 0), vector(0, 0, 1))
   d = Disk(p, T(2))