## Synopsis

```#include <cpml-1/cpml.h>

void                cpml_curve_put_pair_at_time         (`const CpmlPrimitive *curve`,
`double t`,
`CpmlPair *pair`);
void                cpml_curve_put_vector_at_time       (`const CpmlPrimitive *curve`,
`double t`,
`CpmlVector *vector`);
```

## Description

The following functions manipulate `CPML_CURVE` CpmlPrimitive. No validation is made on the input so use the following methods only when you are sure the `primitive` argument is effectively a cubic Bézier curve.

### TODO

• the `get_length()` method must be implemented;
• actually the `put_extents()` method is implemented by computing the bounding box of the control polygon and this will likely include some empty space: there is room for improvements;
• the `put_pair_at()` method must be implemented;
• the `put_vector_at()` method must be implemented;
• the `get_closest_pos()` method must be implemented;
• the `put_intersections()` method must be implemented;
• by default, the offset curve is calculated by using the point at t=0.5 as reference: use a better candidate;
• in the `offset()` implementation, when the equations are inconsistent, the alternative approach performs very bad if `v0` and `v3` are opposite or staggered.

### Offseting algorithm

Given a cubic Bézier primitive, it must be found the approximated Bézier curve parallel to the original one at a specific distance (the so called "offset curve"). The four points needed to build the new curve must be returned.

To solve the offset problem, a custom algorithm is used. First, the resulting curve MUST have the same slope at the start and end point. These constraints are not sufficient to resolve the system, so I let the curve pass thought a given point (pm, known and got from the original curve) at a given time (m, now hardcoded to 0.5).

Firstly, some useful variables are defined:

 ```1 2 3 4``` ```v0 = unitvector(p[1] − p[0]) × offset; v3 = unitvector(p[3] − p[2]) × offset; p0 = p[0] + normal v0; p3 = p[3] + normal v3.```

The resulting curve must have the same slopes than the original one at the start and end points. Forcing the same slopes means:

 `1` `p1 = p0 + k0 v0.`

where `k0` is an arbitrary factor. Decomposing for `x` and `y`:

 ```1 2``` ```p1.x = p0.x + k0 v0.x; p1.y = p0.y + k0 v0.y.```

and doing the same for the end point:

 ```1 2``` ```p2.x = p3.x + k3 v3.x; p2.y = p3.y + k3 v3.y.```

This does not give a resolvable system, though. The curve will be interpolated by forcing its path to pass throught `pm` when `time` is `m`, where `0 ≤ m ≤ 1`. Knowing the function of the cubic Bézier:

 `1` `C(t) = (1 − t)³p0 + 3t(1 − t)²p1 + 3t²(1 − t)p2 + t³p3.`

and forcing `t = m` and `C(t) = pm`:

 ```1 2 3``` ```pm = (1 − m)³p0 + 3m(1 − m)²p1 + 3m²(1 − m)p2 + m³p3. (1 − m) p1 + m p2 = (pm − (1 − m)³p0 − m³p3) / (3m (1 − m)).```

gives this final system:

 ```1 2 3 4 5 6``` ```p1.x = p0.x + k0 v0.x; p1.y = p0.y + k0 v0.y; p2.x = p3.x + k3 v3.x; p2.y = p3.y + k3 v3.y; (1 − m) p1.x + m p2.x = (pm.x − (1 − m)³p0.x − m³p3.x) / (3m (1 − m)); (1 − m) p1.y + m p2.y = (pm.y − (1 − m)³p0.y − m³p3.y) / (3m (1 − m)).```

Substituting and resolving for `k0` and `k3`:

 ```1 2 3 4 5``` ```(1 − m) k0 v0.x + m k3 v3.x = (pm.x − (1 − m)³p0.x − m³p3.x) / (3m (1 − m)) − (1 − m) p0.x − m p3.x; (1 − m) k0 v0.y + m k3 v3.y = (pm.y − (1 − m)³p0.y − m³p3.y) / (3m (1 − m)) − (1 − m) p0.y − m p3.y. (1 − m) k0 v0.x + m k3 v3.x = (pm.x − (1 − m)²(1+2m) p0.x − m²(3 − 2m) p3.x) / (3m (1 − m)); (1 − m) k0 v0.y + m k3 v3.y = (pm.y − (1 − m)²(1+2m) p0.y − m²(3 − 2m) p3.y) / (3m (1 − m)).```

Letting:

 `1` `pk = (pm − (1 − m)²(1+2m) p0 − m²(3 − 2m) p3) / (3m (1 − m)).`

reduces the above to this final equations:

 ```1 2``` ```(1 − m) k0 v0.x + m k3 v3.x = pk.x; (1 − m) k0 v0.y + m k3 v3.y = pk.y.```

If `v0.x ≠ 0`, the system can be resolved for `k0` and `k3` calculated accordingly:

 ```1 2 3 4 5 6 7 8``` ```k0 = (pk.x − m k3 v3.x) / ((1 − m) v0.x); (pk.x − m k3 v3.x) v0.y / v0.x + m k3 v3.y = pk.y. k0 = (pk.x − m k3 v3.x) / ((1 − m) v0.x); k3 m (v3.y − v3.x v0.y / v0.x) = pk.y − pk.x v0.y / v0.x. k3 = (pk.y − pk.x v0.y / v0.x) / (m (v3.y − v3.x v0.y / v0.x)); k0 = (pk.x − m k3 v3.x) / ((1 − m) v0.x).```

Otherwise, if `v3.x ≠ 0`, the system can be solved for `k3` and `k0` calculated accordingly:

 ```1 2 3 4 5 6 7 8``` ```k3 = (pk.x − (1 − m) k0 v0.x) / (m v3.x); (1 − m) k0 v0.y + (pk.x − (1 − m) k0 v0.x) v3.y / v3.x = pk.y. k3 = (pk.x − (1 − m) k0 v0.x) / (m v3.x); k0 (1 − m) (v0.y − k0 v0.x v3.y / v3.x) = pk.y − pk.x v3.y / v3.x. k0 = (pk.y − pk.x v3.y / v3.x) / ((1 − m) (v0.y − v0.x v3.y / v3.x)); k3 = (pk.x − (1 − m) k0 v0.x) / (m v3.x).```

The whole process must be guarded against division by 0 exceptions. If either `v0.x` and `v3.x` are `0`, the first equation will be inconsistent. More in general, the `v0.x × v3.y = v3.x × v3.y` condition must be avoided. This is the first situation to avoid, in which case an alternative approach should be used.

## Details

### cpml_curve_put_pair_at_time ()

```void                cpml_curve_put_pair_at_time         (`const CpmlPrimitive *curve`,
`double t`,
`CpmlPair *pair`);```

Given the `curve` Bézier cubic, finds the coordinates at time `t` (where 0 is the start and 1 is the end) and stores the result in `pair`. Keep in mind `t` is not homogeneous, so 0.5 does not necessarily means the mid point.

 `curve` : the CpmlPrimitive curve data `t` : the "time" value `pair` : the destination pair

Since 1.0

### cpml_curve_put_vector_at_time ()

```void                cpml_curve_put_vector_at_time       (`const CpmlPrimitive *curve`,
`double t`,
`CpmlVector *vector`);```

Given the `curve` Bézier cubic, finds the slope at time `t` (where 0 is the start and 1 is the end) and stores the result in `vector`. Keep in mind `t` is not homogeneous, so 0.5 does not necessarily means the mid point.

 `curve` : the CpmlPrimitive curve data `t` : the "time" value `vector` : the destination vector

Since 1.0