Interpolation between curves

The audience approached me with the following msg: 

- Hey, you! Give us something else but points, even just for an instance!

Well, here comes. One project had an interesting problem: Given two ships starting approximately from the same point, what could be a way to interpolate their routes? A simple answer is  to use positions after a time shift so that the interpolated movement starts from a nice initial point $p_{ab}= (p_{a0}+p_{b0})/2$, where $p_{a0}=p_a(t_a)$ is a chosen reference point of a ship $a$ at a moment $t=t_a$, $p_{b0}= p_b(t_b)$ is the closest point a ship $b$ gets to $p_{a0}$ on its own route at $t=t_b$, and $t_{ab}=(t_a+t_b)/2$, see Fig.0. Then:\[p(t-t_{ab})+p_{ab}=w_a p_a(t+t_a) + w_b p_b(t+t_b),\] where $w_a,w_b, w_a+w_b= 1$ are the interpolation weights and $p(0)= 0$. But this leads to some awkward outcomes when ships move at different speed. It is possible to get loops and back and worth movements, for example. 

Fig. 0: lineasr interpolation between ship positions.

Another way is to consider routes without speeds and consider the steering process (rudder movements and such, this all being summarized as orientation change per cruising length) as the signal to be interpolated. The speeds can be processed separatedly, then. Let us define an operation $K= C(P)$ which extracts an oriented curvature $K= \{\kappa(s)\,|s\in S\}$ along a curve $P= \{p(s)\,| s\in S\}$ where $S$ is the curve length parameter $S= [0,L]\subset\mathbb{R}$. This is similar (but not equal) to differentiation twice. Then one needs an inverse operation $P= C^{-1}(K, p_0,\beta_0)$ which is similar (but not equal) to integration twice with two boundary conditions: $p_0= p(t=0)$ and $\beta_0$ is the initial tangential direction of $P$ at $t=0$. With this tool set, one can express the curvature interpolation of curves:\[P_{ab}=C^{-1}\left( w_a C(P_a) + w_b C(P_b), p_{ab},\beta_{ab}\right)\]

Fig. 1: Indexing of curve segments. All is quite normal, eh?

Just have to define $C$ and its inverse, anymore. The curvature\[\kappa(s_i)= 2\frac{\beta_i - \beta_{i-1}}{l_i+ l_{i-1}},\] where the indexing of edge lengths $l_i$ and tangential orientations $\beta_i$ are shown in Fig. 1. The inverse process requires an initial point $p_i=p_{ab}$ with $i= 1$ and initial tangential orientation $\beta_i= \beta_{ab}$. Then an Euler -like segment-by-segment construction rule is: \[ \beta_i= \beta_{i-1} + l_{i-1}\kappa_{i-1}\\ p_i= p_{i-1} + l_i \left(\cos(\beta_i),\sin(\beta_i) \right)\].

I mentioned Euler? It means Euler integration, which is notorious in producing cumulative $\mathcal{O}(\epsilon)$ error due to discretization size $\epsilon$. But this can be alleviated with better integration procedures, e.g. 2-3 degree Runge-Kutta, see Hairer & Wanner 1993. But we demonstrate here the results with much too sparse segmentation just to show both that the method works, and what are the perilous effects of noise and errors in the path description. 

Fig. 2 shows how the curve P3 (the rightmost one) has a small detail in the top which has a prominent effeft to interpolations, since also P2 has a similar, yet larger part. 

Fig. 2. Interpolation between curves P2 and P3

Fig. 3. shows a large shape error when replicating the shape of P1 with the weight option $w_1= 1,\,w_3=0$. This means that P1 has too few points to properly discretize the curvature in the beginning, whereas errors do not seem to cumulate in the end part.  

         Fig. 3: Interpolation between curves P1 and P3

Fig. 4 shows a drastic shoot-off in the length of the interpolated curve, which happens when $w_1= 0.05 ...0.15$. There are several possible control mechanisms for cases, where there are very short and very long curves to be interpolated. One is to use absolute lengths (as is shown here), or relative lengths. Absolute lengths scheme requires referencing only to those curves which have still length for each $p(s)$ target point. So, the longest route dominates the end shape, but the orientation has been defined by the common interplay of all curves.  

 

Fig. 4: Interpolation between curves P1 and P2. 

Finally, the Fig 5. has a case where $w_1= w_2=w_3= 1/3$. The initial contours counter each other but the kink to the left in the end is common to all and is shown in the final result. Note how the straightened length reaches further than the interpolant curves themselves. 

Fig. 5: The mean of all 3 curves. 

If this feels somewhat odd and awkward, maybe it is because it it. It is a novel way to interpolate over several curves. There are some similar approaches (see e.g. the Abode curve drawing method (2017)), but they do not focus in the interpolation problem. And then there are moving window approaches, which do not utilize the curvature representation. 

How common is an isosceles triangle?

 Isosceles triangle $t$ is one which has two sides with equal lengths, But it should not be equilateral. (Logically speaking every equilateral triangle is also isosceles in three different ways, but let us not be that logical). And talking about equality, we need a relative tolerance $\epsilon_{12}$ for any two sides and another relative tolerance $\epsilon_3$ for any third side. And it will be hard to dictate the mostly psychological limit when a triangle stops being equilateral and starts being an isosceles. 

Anyways, we are going towards a predicate $P(t,\epsilon_{12},\epsilon_3)$, which is true when the triangle is isosceles with given two limits.\[P(t,\epsilon_{12},\epsilon_3) \equiv |l_1 - l_2| < \epsilon_{12} l_{12} \text{ and } \epsilon_3 l_{12} < |l_{12} - l_3|\], where $l_{12}= (l_1+l_2)/2$ is the mean of any two edge lengths and $\epsilon_{12} \ll \epsilon_3$ and $\{l_1,l_2,l_3\}$ are the edge lenghts of the triangle $t$.  In plain English, two sides of $t$ are of an approx. equal length with a difference of $\epsilon_{12}$ at most, and the third edge length is a significantly different from those two. So, how often one encounters such a triangle? 

Img 1. An approximately isosceles triangle with shape tolerances $\epsilon_{12}$ and $\epsilon_3$.
The tolerance condition visualized as a red segment. 

To fix the psychological ratio value $\epsilon_3/\epsilon_{12}$ I conducted an experiment with a set of researchers (5 persons) in a distinguished University. It seems that $\epsilon_3\approx 5\epsilon_{12}$ is a rather good choice. Have to note that equilateral triangles are of second order rarity (prob(isosceles):prob(equilater) $\approx \mathcal{O}(\epsilon_{12}):\mathcal{O}(\epsilon_{12}^2)$.

Another formulation would define a predicate $P_3(.)$ for an equilateral triangle $t$:\[P_3(t,\epsilon_{12},\epsilon_3)\equiv \max_i(l_i) - \min_i(l_i)  < \epsilon_3 \text{mean}_i \, l_i\].Then, an isosceles would have a predicate $P_{12}(.)$:\[P_{12}(t,\epsilon_{12}) \equiv |l_1 - l_2| < \epsilon_{12} l_{12} \text{ and not } P_3(t,\epsilon_3)\] .  

One almost misses probabilistic geometry here... (Just as there is geometric algebra and algebraic geometry, there is also probabilistic geometry and geometric probability. The latter is -uhm- somewhat constructive and experimental, and the former theoretically solid, but somewhat obscure). So, let's go along the geometric probability road and make a choice of the topology of the probabilistic measure. Delaunay triangulation it is, and with it, the evident boundary effect. And the boundary effect depends on the number of points, when the PC shape is fixed (in this case, fixed to a unit square with a boundary of length 4). We have already outlined how to measure the 'boundariousness' of a PC, and how to define a boundary point. But now we do some census over the inhabitants of the Delaunay triangulation, approximating an infinite 2D plane with random points, at certain densities. 

But about the results: Below left is the ratio $R$ of observed isosceles from a set $P$ of points in a unit square. The small size of $P$ is revealed by fluctuations in different computations (top right corner). An indication that the $R$ really is of first order of the shape accuracy limit $\epsilon_{12}$, is at the right detail, where the ratio $R/\epsilon_{12}$ is shown to stay within a limited range while $\epsilon_{12}$ has a thousandfold range. 

Img. 2. Ratio $R$ of isosceles triangles in a Delaunay triangulation (left).  Visual proof of  $R\in\mathcal{O}(\epsilon_{12})$ behavior (right). The range of values stays within 1...45 while $\epsilon_{12}$ range is 1...1000.

A simplistic explanation for the category of $R$ is that the isosceles top vertex inhabits a narrow vertical band, when the best candidate for the base edge is scaled to unit length, see a sketch below. The band increases its width asymptotically towards infinity, and has some slight complexities at the spots where equilater triangles are. And long needle like triangles (those with the unit scaled presentation having the top vertex close to infinity) are more rare than others, but this does not get analyzed this time. The previous image is a better verification of the $T\in\mathcal{O}(\epsilon)$ behavior. 

Img. 3. Isosceles triangles (striped zone) with the base edge as a unit segment. The part below is symmetrical and omitted. A missing spot is the domain of equilateral triangles. 

There is something to scrutinize in Fig. 2. The Delauny process (should be familiar for the reader at this point) favors close to equilateral triangles because of the so called Delaunay property (every triangle $t$ defines an enclosing circle which has no other points within it than those of $t$), and maybe produces more isoclines, too? This is true, but requires a bit more finesse, and is left out of this text. And extension to 3D is obvious, but rather fruitless.