Add parametric function for Chain (#1155)

* add parametric function for Chain

* add tests

* Apply suggestions from code review

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

* reduce allocations

* Update src/geometries/polytopes.jl

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

* fix

* use measure

* Apply suggestions from code review

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

* Fix type instabilities

* add tests for Cartesian coords

* Apply suggestions from code review

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

* remove approximate latlon tests

* Adjust tests

---------

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>
Co-authored-by: Júlio Hoffimann <[email protected]>

Author: 3 people

src/geometries/polytopes.jl

@@ -8,7 +8,7 @@
 We say that a geometry is a K-polytope when it is a collection of "flat" sides
 that constitute a `K`-dimensional subspace. They are called chain, polygon and
 polyhedron respectively for 1D (`K=1`), 2D (`K=2`) and 3D (`K=3`) subspaces.
-The parameter `K` is also known as the rank or parametric dimension 
+The parameter `K` is also known as the rank or parametric dimension
 of the polytope (<https://en.wikipedia.org/wiki/Abstract_polytope>).
 
 The term polytope expresses a particular combinatorial structure. A polyhedron,
@@ -156,6 +156,23 @@ function angles(c::Chain)
   map(i -> ∠(vs[i - 1], vs[i], vs[i + 1]), i1:i2)
 end
 
+function (c::Chain)(t)
+  if t < 0 || t > 1
+    throw(DomainError(t, "c(t) is not defined for t outside [0, 1]."))
+  end
+  segs = segments(c)
+  sums = cumsum(measure.(segs))
+  sums /= last(sums)
+  # find k such that sums[k] ≤ t < sums[k + 1]
+  k = searchsortedfirst(sums, t) - 1
+  # select segment s at index k
+  s, _ = iterate(segs, k)
+  # reparametrization of t within s
+  ∑ₖ = iszero(k) ? zero(eltype(sums)) : sums[k]
+  ∑ₖ₊₁ = sums[k + 1]
+  s((t - ∑ₖ) / (∑ₖ₊₁ - ∑ₖ))
+end
+
 # implementations of Chain
 include("polytopes/segment.jl")
 include("polytopes/rope.jl")

test/polytopes.jl

@@ -264,6 +264,25 @@ end
   r = Ring(merc.([(0, 0), (1, 0), (1, 1), (0, 1)]))
   @test crs(centroid(r)) === crs(r)
 
+  # parameterization
+  r = boundary(Box(cart(0, 0), cart(1, 1)))
+  @test r(T(0)) == cart(0, 0)
+  @test r(T(0.25)) == cart(1, 0)
+  @test r(T(0.5)) == cart(1, 1)
+  @test r(T(0.75)) == cart(0, 1)
+  r = Rope(cart(0, 0), cart(3, 0), cart(4, 0))
+  @test r(T(0)) == cart(0, 0)
+  @test r(T(0.25)) == cart(1, 0)
+  @test r(T(0.5)) == cart(2, 0)
+  @test r(T(0.75)) == cart(3, 0)
+  @test r(T(1)) == cart(4, 0)
+  r = Ring(latlon(45, 0), latlon(45, 90), latlon(90, 90))
+  @test r(T(0)) == latlon(45, 0)
+  @test r(T(1)) == latlon(45, 0)
+  r = Rope(latlon(45, 0), latlon(45, 90), latlon(90, 90))
+  @test r(T(0)) == latlon(45, 0)
+  @test r(T(1)) == latlon(90, 90)
+
   ri = Ring(cart.([(1, 1), (2, 2), (3, 3)]))
   ro = Rope(cart.([(1, 1), (2, 2), (3, 3)]))
   @test sprint(show, ri) == "Ring((x: 1.0 m, y: 1.0 m), (x: 2.0 m, y: 2.0 m), (x: 3.0 m, y: 3.0 m))"