TinySpline
Spline Library for a Multitude of Programming Languages
v0.3.0
Classes | Macros | Typedefs | Enumerations | Functions
tinyspline.h File Reference
#include <stddef.h>
+ Include dependency graph for tinyspline.h:
+ This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

struct  tsStatus
 
struct  tsBSpline
 
struct  tsDeBoorNet
 

Macros

#define TS_MAX_NUM_KNOTS   10000
 
#define TS_MIN_KNOT_VALUE   0.0f
 
#define TS_MAX_KNOT_VALUE   1.0f
 
#define TS_KNOT_EPSILON   1e-4
 
#define TS_TRY(label, error, status)
 
#define TS_CALL(label, error, call)
 
#define TS_CATCH(error)   } if ((error)) {
 
#define TS_FINALLY   } {
 
#define TS_END_TRY
 
#define TS_END_TRY_RETURN(error)   TS_END_TRY return (error);
 
#define TS_CALL_ROE(error, call)
 
#define TS_RETURN_SUCCESS(status)
 
#define TS_RETURN_0(status, error, msg)
 
#define TS_RETURN_1(status, error, msg, arg1)
 
#define TS_RETURN_2(status, error, msg, arg1, arg2)
 
#define TS_RETURN_3(status, error, msg, arg1, arg2, arg3)
 
#define TS_THROW_0(label, error, status, val, msg)
 
#define TS_THROW_1(label, error, status, val, msg, arg1)
 
#define TS_THROW_2(label, error, status, val, msg, arg1, arg2)
 
#define TS_THROW_3(label, error, status, val, msg, arg1, arg2, arg3)
 

Typedefs

typedef double tsReal
 

Enumerations

enum  tsError {
  TS_SUCCESS = 0, TS_MALLOC = -1, TS_DIM_ZERO = -2, TS_DEG_GE_NCTRLP = -3,
  TS_U_UNDEFINED = -4, TS_MULTIPLICITY = -5, TS_KNOTS_DECR = -6, TS_NUM_KNOTS = -7,
  TS_UNDERIVABLE = -8, TS_LCTRLP_DIM_MISMATCH = -10, TS_IO_ERROR = -11, TS_PARSE_ERROR = -12,
  TS_INDEX_ERROR = -13, TS_NO_RESULT = -14
}
 
enum  tsBSplineType { TS_OPENED = 0, TS_CLAMPED = 1, TS_BEZIERS = 2 }
 

Functions

size_t ts_bspline_degree (const tsBSpline *spline)
 
tsError ts_bspline_set_degree (tsBSpline *spline, size_t deg, tsStatus *status)
 
size_t ts_bspline_order (const tsBSpline *spline)
 
tsError ts_bspline_set_order (tsBSpline *spline, size_t order, tsStatus *status)
 
size_t ts_bspline_dimension (const tsBSpline *spline)
 
tsError ts_bspline_set_dimension (tsBSpline *spline, size_t dim, tsStatus *status)
 
size_t ts_bspline_len_control_points (const tsBSpline *spline)
 
size_t ts_bspline_num_control_points (const tsBSpline *spline)
 
size_t ts_bspline_sof_control_points (const tsBSpline *spline)
 
tsError ts_bspline_control_points (const tsBSpline *spline, tsReal **ctrlp, tsStatus *status)
 
tsError ts_bspline_control_point_at (const tsBSpline *spline, size_t index, tsReal **ctrlp, tsStatus *status)
 
tsError ts_bspline_set_control_points (tsBSpline *spline, const tsReal *ctrlp, tsStatus *status)
 
tsError ts_bspline_set_control_point_at (tsBSpline *spline, size_t index, const tsReal *ctrlp, tsStatus *status)
 
size_t ts_bspline_num_knots (const tsBSpline *spline)
 
size_t ts_bspline_sof_knots (const tsBSpline *spline)
 
tsError ts_bspline_knots (const tsBSpline *spline, tsReal **knots, tsStatus *status)
 
tsError ts_bspline_set_knots (tsBSpline *spline, const tsReal *knots, tsStatus *status)
 
tsReal ts_deboornet_knot (const tsDeBoorNet *net)
 
size_t ts_deboornet_index (const tsDeBoorNet *net)
 
size_t ts_deboornet_multiplicity (const tsDeBoorNet *net)
 
size_t ts_deboornet_num_insertions (const tsDeBoorNet *net)
 
size_t ts_deboornet_dimension (const tsDeBoorNet *net)
 
size_t ts_deboornet_len_points (const tsDeBoorNet *net)
 
size_t ts_deboornet_num_points (const tsDeBoorNet *net)
 
size_t ts_deboornet_sof_points (const tsDeBoorNet *net)
 
tsError ts_deboornet_points (const tsDeBoorNet *net, tsReal **points, tsStatus *status)
 
size_t ts_deboornet_len_result (const tsDeBoorNet *net)
 
size_t ts_deboornet_num_result (const tsDeBoorNet *net)
 
size_t ts_deboornet_sof_result (const tsDeBoorNet *net)
 
tsError ts_deboornet_result (const tsDeBoorNet *net, tsReal **result, tsStatus *status)
 
tsBSpline ts_bspline_init ()
 
tsError ts_bspline_new (size_t num_control_points, size_t dimension, size_t degree, tsBSplineType type, tsBSpline *spline, tsStatus *status)
 
tsError ts_bspline_copy (const tsBSpline *src, tsBSpline *dest, tsStatus *status)
 
void ts_bspline_move (tsBSpline *src, tsBSpline *dest)
 
void ts_bspline_free (tsBSpline *spline)
 
tsDeBoorNet ts_deboornet_init ()
 
tsError ts_deboornet_copy (const tsDeBoorNet *src, tsDeBoorNet *dest, tsStatus *status)
 
void ts_deboornet_move (tsDeBoorNet *src, tsDeBoorNet *dest)
 
void ts_deboornet_free (tsDeBoorNet *net)
 
tsError ts_bspline_interpolate_cubic (const tsReal *points, size_t n, size_t dim, tsBSpline *_spline_, tsStatus *status)
 
size_t ts_bspline_num_distinct_knots (const tsBSpline *spline)
 
tsError ts_bspline_eval (const tsBSpline *spline, tsReal u, tsDeBoorNet *_deBoorNet_, tsStatus *status)
 
tsError ts_bspline_eval_all (const tsBSpline *spline, const tsReal *us, size_t num, tsReal **points, tsStatus *status)
 
tsError ts_bspline_sample (const tsBSpline *spline, size_t num, tsReal **points, size_t *actual_num, tsStatus *status)
 
tsError ts_bspline_bisect (const tsBSpline *spline, tsReal value, tsReal epsilon, int persnickety, size_t index, int ascending, size_t max_iter, tsDeBoorNet *net, tsStatus *status)
 
void ts_bspline_domain (const tsBSpline *spline, tsReal *min, tsReal *max)
 
tsError ts_bspline_is_closed (const tsBSpline *spline, tsReal epsilon, int *closed, tsStatus *status)
 
tsError ts_bspline_derive (const tsBSpline *spline, size_t n, tsBSpline *_derivative_, tsStatus *status)
 
tsError ts_bspline_insert_knot (const tsBSpline *spline, tsReal u, size_t n, tsBSpline *_result_, size_t *_k_, tsStatus *status)
 
tsError ts_bspline_split (const tsBSpline *spline, tsReal u, tsBSpline *_split_, size_t *_k_, tsStatus *status)
 
tsError ts_bspline_tension (const tsBSpline *spline, tsReal tension, tsBSpline *out, tsStatus *status)
 
tsError ts_bspline_to_beziers (const tsBSpline *spline, tsBSpline *_beziers_, tsStatus *status)
 
tsError ts_bspline_to_json (const tsBSpline *spline, char **_json_, tsStatus *status)
 
tsError ts_bspline_from_json (const char *json, tsBSpline *_spline_, tsStatus *status)
 
tsError ts_bspline_save_json (const tsBSpline *spline, const char *path, tsStatus *status)
 
tsError ts_bspline_load_json (const char *path, tsBSpline *_spline_, tsStatus *status)
 
int ts_knots_equal (tsReal x, tsReal y)
 
void ts_arr_fill (tsReal *arr, size_t num, tsReal val)
 
tsReal ts_distance (const tsReal *x, const tsReal *y, size_t dim)
 

Macro Definition Documentation

◆ TS_CALL

#define TS_CALL (   label,
  error,
  call 
)
Value:
(error) = (call); \
if ((error)) goto __ ## label ## __;

◆ TS_CALL_ROE

#define TS_CALL_ROE (   error,
  call 
)
Value:
{ \
(error) = (call); \
if ((error)) return error; \
}

◆ TS_END_TRY

#define TS_END_TRY
Value:
} \
}

◆ TS_KNOT_EPSILON

#define TS_KNOT_EPSILON   1e-4

Used to check whether two knots are equal.

◆ TS_MAX_KNOT_VALUE

#define TS_MAX_KNOT_VALUE   1.0f

The default maximum knot value of a spline.

◆ TS_MAX_NUM_KNOTS

#define TS_MAX_NUM_KNOTS   10000

The maximum number of knots of a spline.

◆ TS_MIN_KNOT_VALUE

#define TS_MIN_KNOT_VALUE   0.0f

The default minimum knot value of a spline.

◆ TS_RETURN_0

#define TS_RETURN_0 (   status,
  error,
  msg 
)
Value:
{ \
if ((status) != NULL) { \
(status)->code = error; \
sprintf((status)->message, msg); \
} \
return error; \
}

◆ TS_RETURN_1

#define TS_RETURN_1 (   status,
  error,
  msg,
  arg1 
)
Value:
{ \
if ((status) != NULL) { \
(status)->code = error; \
sprintf((status)->message, msg, arg1); \
} \
return error; \
}

◆ TS_RETURN_2

#define TS_RETURN_2 (   status,
  error,
  msg,
  arg1,
  arg2 
)
Value:
{ \
if ((status) != NULL) { \
(status)->code = error; \
sprintf((status)->message, msg, arg1, arg2); \
} \
return error; \
}

◆ TS_RETURN_3

#define TS_RETURN_3 (   status,
  error,
  msg,
  arg1,
  arg2,
  arg3 
)
Value:
{ \
if ((status) != NULL) { \
(status)->code = error; \
sprintf((status)->message, msg, arg1, arg2, arg3); \
} \
return error; \
}

◆ TS_RETURN_SUCCESS

#define TS_RETURN_SUCCESS (   status)
Value:
{ \
if ((status) != NULL) { \
(status)->code = TS_SUCCESS; \
(status)->message[0] = '\0'; \
} \
return TS_SUCCESS; \
}
Definition: tinyspline.h:99

◆ TS_THROW_0

#define TS_THROW_0 (   label,
  error,
  status,
  val,
  msg 
)
Value:
{ \
(error) = val; \
if ((status) != NULL) { \
(status)->code = val; \
sprintf((status)->message, msg); \
} \
goto __ ## label ## __; \
}

◆ TS_THROW_1

#define TS_THROW_1 (   label,
  error,
  status,
  val,
  msg,
  arg1 
)
Value:
{ \
(error) = val; \
if ((status) != NULL) { \
(status)->code = val; \
sprintf((status)->message, msg, arg1); \
} \
goto __ ## label ## __; \
}

◆ TS_THROW_2

#define TS_THROW_2 (   label,
  error,
  status,
  val,
  msg,
  arg1,
  arg2 
)
Value:
{ \
(error) = val; \
if ((status) != NULL) { \
(status)->code = val; \
sprintf((status)->message, msg, arg1, arg2); \
} \
goto __ ## label ## __; \
}

◆ TS_THROW_3

#define TS_THROW_3 (   label,
  error,
  status,
  val,
  msg,
  arg1,
  arg2,
  arg3 
)
Value:
{ \
(error) = val; \
if ((status) != NULL) { \
(status)->code = val; \
sprintf((status)->message, msg, arg1, arg2, arg3); \
} \
goto __ ## label ## __; \
}

◆ TS_TRY

#define TS_TRY (   label,
  error,
  status 
)
Value:
{ \
(error) = TS_SUCCESS; \
if ((status) != NULL) { \
(status)->code = TS_SUCCESS; \
(status)->message[0] = '\0'; \
} \
__ ## label ## __: \
if (!(error)) {
Definition: tinyspline.h:99

Enumeration Type Documentation

◆ tsBSplineType

Describes the structure of the knot vector of a NURBS/B-Spline. For more details, see:

www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/spline/B-spline/bspline-curve.html

◆ tsError

enum tsError

Defines different error codes.

Enumerator
TS_SUCCESS 

No error.

TS_MALLOC 

Unable to allocate memory.

TS_DIM_ZERO 

The dimension of the control points are 0.

TS_DEG_GE_NCTRLP 

degree >= num_control_points.

TS_U_UNDEFINED 

Undefined knot value.

TS_MULTIPLICITY 

s(u) > order

TS_KNOTS_DECR 

Decreasing knot vector.

TS_NUM_KNOTS 

Unexpected number of knots.

TS_UNDERIVABLE 

Spline is not derivable.

TS_LCTRLP_DIM_MISMATCH 

len_control_points % dim != 0.

TS_IO_ERROR 

Error while reading/writing a file.

TS_PARSE_ERROR 

Error while parsing a serialized entity.

TS_INDEX_ERROR 

Index does not exist.

TS_NO_RESULT 

Function returns without result.

Function Documentation

◆ ts_arr_fill()

void ts_arr_fill ( tsReal *  arr,
size_t  num,
tsReal  val 
)

Fills the given array arr with val from arr+0 to arr+ num (exclusive).

◆ ts_bspline_bisect()

tsError ts_bspline_bisect ( const tsBSpline spline,
tsReal  value,
tsReal  epsilon,
int  persnickety,
size_t  index,
int  ascending,
size_t  max_iter,
tsDeBoorNet net,
tsStatus status 
)

Tries to find a point P on spline such that:

ts_distance(P[index], value, 1) <= fabs(epsilon)

This function is using the bisection method to determine P. Accordingly, it is expected that the control points of spline are sorted at component index either in ascending order (if ascending != 0) or in descending order (if ascending == 0). If the control points of spline are not sorted at component index, the behaviour of this function is undefined. For the sake of fail-safeness, the distance of P[index] and value is compared with the absolute value of epsilon (using fabs).

The bisection method is an iterative approach, minimizing the error (epsilon) with each iteration step until an "optimum" was found. However, there may be no point P satisfying the distance condition. Thus, the number of iterations must be limited (max_iter). Depending on the domain of the control points of spline at component index and epsilon, max_iter ranges from 7 to 50. In most cases max_iter == 30 should be fine though. The parameter persnickety allows to define the behaviour of this function is case no point was found after max_iter iterations. If enabled (!= 0), TS_NO_RESULT is returned. If disabled (== 0), the best fitting point is returned.

Parameters
[in]splineThe spline to evaluate
[in]valueThe value (point at component index) to find.
[in]epsilonThe maximum distance (inclusive).
[in]persnicketyIndicates whether TS_NO_RESULT should be returned if there is no point P satisfying the distance condition (!= 0 to enable, == 0 to disable). If disabled, the best fitting point is returned.
[in]indexThe point's component.
[in]ascendingIndicates whether the control points of spline are sorted in ascending (!= 0) or in descending (== 0) order at component index.
[in]max_iterThe maximum number of iterations (16 is a sane default value).
[out]netThe output parameter.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_INDEX_ERROR If the dimension of the control points of spline <= index.
TS_NO_RESULT If persnickety is enabled (!= 0) and there is no point P satisfying the distance condition.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_control_point_at()

tsError ts_bspline_control_point_at ( const tsBSpline spline,
size_t  index,
tsReal **  ctrlp,
tsStatus status 
)

Returns a deep copy of the control point of spline at index (index 0 is the first control points, index 1 is the second control point, and so on).

Parameters
[in]splineThe spline whose control point is read.
[in]indexThe zero based index of the requested control point.
[out]ctrlpThe output array.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_INDEX_ERROR If there is no control point at index.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_control_points()

tsError ts_bspline_control_points ( const tsBSpline spline,
tsReal **  ctrlp,
tsStatus status 
)

Returns a deep copy of the control points of spline.

Parameters
[in]splineThe spline whose control points are read.
[out]ctrlpThe output array.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_copy()

tsError ts_bspline_copy ( const tsBSpline src,
tsBSpline dest,
tsStatus status 
)

Creates a deep copy of src and stores the copied values in dest. Does nothing, if src == dest.

Parameters
[in]srcThe spline to deep copy.
[out]destThe output spline.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_degree()

size_t ts_bspline_degree ( const tsBSpline spline)

Returns the degree of spline.

Parameters
[in]splineThe spline whose degree is read.
Returns
The degree of spline.

◆ ts_bspline_derive()

tsError ts_bspline_derive ( const tsBSpline spline,
size_t  n,
tsBSpline _derivative_,
tsStatus status 
)

Returns the n'th derivative of spline and stores the result in _derivative_. Creates a deep copy of spline, if spline != _result_. For more details, see:

http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/spline/B-spline/bspline-derv.html

The derivative of a spline s of degree d (d > 0) with m control points and n knots is another spline s' of degree d-1 with m-1 control points and n-2 knots, defined over s as:

\begin{eqnarray*} s'(u) &=& \sum_{i=0}^{n-1} N_{i+1,p-1}(u) * (P_{i+1} - P_{i}) * p / (u_{i+p+1}-u_{i+1}) \\ &=& \sum_{i=1}^{n} N_{i,p-1}(u) * (P_{i} - P_{i-1}) * p / (u_{i+p}-u_{i}) \end{eqnarray*}

If s has a clamped knot vector, it can be shown that:

\begin{eqnarray*} s'(u) &=& \sum_{i=0}^{n-1} N_{i,p-1}(u) * (P_{i+1} - P_{i}) * p / (u_{i+p+1}-u_{i+1}) \end{eqnarray*}

where the multiplicity of the first and the last knot value u is p rather than p+1. The derivative of a point (degree == 0) is another point with coordinates 0.

Parameters
splineThe spline to derive.
_derivative_The output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_UNDERIVABLE If the multiplicity of an internal knot of spline is greater than the degree of spline. NOTE: This will be fixed in the future.
TS_MALLOC If allocating memory failed.

< Stores intermediate results.

< Pointer to the control points of worker.

< Pointer to the knots of worker.

< Used in for loops.

< Used to swap worker and derivative.

◆ ts_bspline_dimension()

size_t ts_bspline_dimension ( const tsBSpline spline)

Returns the dimension of spline. The dimension of a spline describes the number of components for each point in ts_bspline_get_control_points(spline). One-dimensional splines are possible, albeit their benefit might be questionable.

Parameters
[in]splineThe spline whose dimension is read.
Returns
The dimension of spline.

◆ ts_bspline_domain()

void ts_bspline_domain ( const tsBSpline spline,
tsReal *  min,
tsReal *  max 
)

Returns the domain of spline.

Parameters
[in]splineThe spline to query.
[out]minThe lower bound of the domain of spline.
[out]maxThe upper bound of the domain of spline.

◆ ts_bspline_eval()

tsError ts_bspline_eval ( const tsBSpline spline,
tsReal  u,
tsDeBoorNet _deBoorNet_,
tsStatus status 
)

Evaluates spline at knot value u and stores the result in _deBoorNet_.

Parameters
splineThe spline to evaluate.
uThe knot value to evaluate.
_deBoorNet_The output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_U_UNDEFINED If spline is not defined at knot value u.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_eval_all()

tsError ts_bspline_eval_all ( const tsBSpline spline,
const tsReal *  us,
size_t  num,
tsReal **  points,
tsStatus status 
)

Evaluates spline at the given knot values us and stores the result points in points. If us contains one or more knot values where spline is discontinuous at, only the first point of the corresponding evaluation result is taken. After calling this function points contains exactly num * ts_bspline_dimension(spline) values.

This function is in particular useful in cases where a multitude of knots need to be evaluated, because only a single instance of tsDeBoorNet is created and reused for all evaluation tasks (therefore the memory footprint is reduced to a minimum).

Parameters
[in]splineThe spline to evaluate.
[in]usThe knot values to evaluate.
[in]numThe number of knots in us.
[out]pointsThe output parameter.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_U_UNDEFINED If spline is not defined at one of the knot values in us.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_free()

void ts_bspline_free ( tsBSpline spline)

Frees the data of spline. After calling this function, the data of spline points to NULL.

Parameters
[out]splineThe spline to free.

◆ ts_bspline_from_json()

tsError ts_bspline_from_json ( const char *  json,
tsBSpline _spline_,
tsStatus status 
)

Parses json and stores the result in _spline_.

Parameters
jsonThe JSON string to parse.
<em>spline</em>The output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_PARSE_ERROR If an error occurred while parsing json.
TS_DIM_ZERO If the dimension is 0.
TS_LCTRLP_DIM_MISMATCH If the length of the control point vector modulo dimension is not 0.
TS_DEG_GE_NCTRLP If the degree is greater or equals to the number of control points.
TS_NUM_KNOTS If the number of knots stored in json does not match to the number of control points and the degree of the spline.
TS_KNOTS_DECR If the knot vector is decreasing.
TS_MULTIPLICITY If there is a knot with multiplicity greater than order.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_init()

tsBSpline ts_bspline_init ( )

Creates a new spline whose data points to NULL.

Returns
A new spline whose data points to NULL.

◆ ts_bspline_insert_knot()

tsError ts_bspline_insert_knot ( const tsBSpline spline,
tsReal  u,
size_t  n,
tsBSpline _result_,
size_t *  _k_,
tsStatus status 
)

Inserts the knot value u up to n times into spline and stores the result in _result_. Creates a deep copy of spline, if spline != _result_.

Parameters
splineThe spline whose knot vector is extended.
uThe knot value to insert.
nHow many times u should be inserted.
_result_The output parameter.
_k_Stores the last index of u in _result_.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_interpolate_cubic()

tsError ts_bspline_interpolate_cubic ( const tsReal *  points,
size_t  n,
size_t  dim,
tsBSpline _spline_,
tsStatus status 
)

Performs a cubic spline interpolation using the thomas algorithm, see:

https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm
http://www.math.ucla.edu/~baker/149.1.02w/handouts/dd_splines.pdf
http://www.bakoma-tex.com/doc/generic/pst-bspline/pst-bspline-doc.pdf

The resulting spline is a sequence of bezier curves connecting each point in points. Each bezier curve b is of degree 3 with dim being the dimension of the each control point in b. The total number of control points is (n - 1) * 4.

On error, all values of _spline_ are set to 0/NULL.

Note: n is the number of points in points and not the length of points. For instance, the following point vector yields to n = 4 and dim = 2:

[x0, y0, x1, y1, x2, y2, x3, y3]
Parameters
pointsThe points to interpolate.
nThe number of points in points.
dimThe dimension of each control point in _spline_.
_spline_The output parameter storing the result of this function.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_DIM_ZERO If dim == 0.
TS_DEG_GE_NCTRLP If n < 2.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_is_closed()

tsError ts_bspline_is_closed ( const tsBSpline spline,
tsReal  epsilon,
int *  closed,
tsStatus status 
)

Returns whether the distance of the first and the last point (not necessarily the control points) of spline is less than or equal to epsilon. The distance of possible inner gaps is not considered.

Parameters
splineThe spline to query.
epsilonThe maximum distance.
closedThe output parameter. Is set to 1 if the distance of the first and the last point (not necessarily the control points) of spline is less than or equal to epsilon. Is set to 0 otherwise.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_knots()

tsError ts_bspline_knots ( const tsBSpline spline,
tsReal **  knots,
tsStatus status 
)

Returns a deep copy of the knots of spline.

Parameters
[in]splineThe spline whose knots are read.
[out]knotsThe output array.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_len_control_points()

size_t ts_bspline_len_control_points ( const tsBSpline spline)

Returns the length of the control point array of spline.

Parameters
[in]splineThe spline with its control point array whose length is read.
Returns
The length of the control point array of spline.

◆ ts_bspline_load_json()

tsError ts_bspline_load_json ( const char *  path,
tsBSpline _spline_,
tsStatus status 
)

Loads the contents of path and stores the result in _spline_.

Parameters
pathPath of the JSON file.
<em>spline</em>The output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_IO_ERROR If path does not exist.
TS_PARSE_ERROR If an error occurred while parsing the contents of path.
TS_DIM_ZERO If the dimension is 0.
TS_LCTRLP_DIM_MISMATCH If the length of the control point vector modulo dimension is not 0.
TS_DEG_GE_NCTRLP If the degree is greater or equals to the number of control points.
TS_NUM_KNOTS If the number of knots stored in json does not match to the number of control points and the degree of the spline.
TS_KNOTS_DECR If the knot vector is decreasing.
TS_MULTIPLICITY If there is a knot with multiplicity greater than order.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_move()

void ts_bspline_move ( tsBSpline src,
tsBSpline dest 
)

Moves the ownership of the data of src to dest. After calling this function, the data of src points to NULL. Does not free the data of dest. Does nothing, if src == dest.

Parameters
[out]srcThe spline whose values are moved to dest.
[out]destThe spline that receives the values of src.

◆ ts_bspline_new()

tsError ts_bspline_new ( size_t  num_control_points,
size_t  dimension,
size_t  degree,
tsBSplineType  type,
tsBSpline spline,
tsStatus status 
)

Creates a new spline and stores the result in spline.

Parameters
[in]num_control_pointsThe number of control points of spline.
[in]dimensionThe dimension of each control point of spline.
[in]degreeThe degree of spline.
[in]typeHow to setup the knot vector of spline.
[out]splineThe output spline.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_DIM_ZERO If dimension == 0.
TS_DEG_GE_NCTRLP If degree >= num_control_points.
TS_NUM_KNOTS If type == ::TS_BEZIERS and (num_control_points % degree + 1) != 0.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_num_control_points()

size_t ts_bspline_num_control_points ( const tsBSpline spline)

Returns the number of control points of spline.

Parameters
[in]splineThe spline whose number of control points is read.
Returns
The number of control points of spline.

◆ ts_bspline_num_distinct_knots()

size_t ts_bspline_num_distinct_knots ( const tsBSpline spline)

Returns the number of distinct knots.

Parameters
splineThe spline to query.
Returns
The number of distinct knots.

◆ ts_bspline_num_knots()

size_t ts_bspline_num_knots ( const tsBSpline spline)

Returns the number of knots of spline.

Parameters
[in]splineThe spline whose number of knots is read.
Returns
The number of knots of spline.

◆ ts_bspline_order()

size_t ts_bspline_order ( const tsBSpline spline)

Returns the order (degree + 1) of spline.

Parameters
[in]splineThe spline whose order is read.
Returns
The order of spline.

◆ ts_bspline_sample()

tsError ts_bspline_sample ( const tsBSpline spline,
size_t  num,
tsReal **  points,
size_t *  actual_num,
tsStatus status 
)

Generates a sequence of num different knots (The knots are equally distributed between the minimum and the maximum of the domain of spline), passes this sequence to ts_bspline_eval_all, and stores the result points in points. If num is 0, the following default value is taken as fallback:

num = (ts_bspline_num_control_points(spline) -
    ts_bspline_degree(spline)) * 30;

That is, the fallback generates 30 knots per Bezier segment. For the sake of stability regarding future changes, the actual number of generated knots (which only differs from num if num is 0) is stored in actual_num. If num is 1, the point located at the minimum of the domain of spline is evaluated.

Parameters
splineThe spline to evaluate.
numThe number of knots to generate.
pointsThe output parameter.
actual_numThe actual number of generated knots. Differs from num only if num is 0. Must not be NULL.
statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_save_json()

tsError ts_bspline_save_json ( const tsBSpline spline,
const char *  path,
tsStatus status 
)

Saves spline as JSON ASCII file.

Parameters
splineThe spline to save.
pathPath of the JSON file.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_IO_ERROR If an error occurred while saving spline.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_set_control_point_at()

tsError ts_bspline_set_control_point_at ( tsBSpline spline,
size_t  index,
const tsReal *  ctrlp,
tsStatus status 
)

Sets the control point of spline at index. Creates a deep copy of ctrlp.

Parameters
[out]splineThe spline whose control point is set.
[in]indexThe zero based index of the control point to set.
[in]ctrlpThe values to deep copy.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_INDEX_ERROR If there is no control point at index.

◆ ts_bspline_set_control_points()

tsError ts_bspline_set_control_points ( tsBSpline spline,
const tsReal *  ctrlp,
tsStatus status 
)

Sets the control points of spline. Creates a deep copy of ctrlp.

Parameters
[out]splineThe spline whose control points are set.
[in]ctrlpThe values to deep copy.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.

◆ ts_bspline_set_degree()

tsError ts_bspline_set_degree ( tsBSpline spline,
size_t  deg,
tsStatus status 
)

Sets the degree of spline.

Parameters
[out]splineThe spline whose degree is set.
[in]degThe degree to be set.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_DEG_GE_NCTRLP If degree >= ts_bspline_get_control_points(spline).

◆ ts_bspline_set_dimension()

tsError ts_bspline_set_dimension ( tsBSpline spline,
size_t  dim,
tsStatus status 
)

Sets the dimension of spline. The following conditions must be satisfied:

(1) dim >= 1
(2) len_control_points(spline) % dim == 0

with len_control_points being the length of the control point array of spline. The dimension of a spline describes the number of components for each point in ts_bspline_get_control_points(spline). One-dimensional splines are possible, albeit their benefit might be questionable.

Parameters
[out]splineThe spline whose dimension is set.
[in]dimThe dimension to be set.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_DIM_ZERO If dimension == 0.
TS_LCTRLP_DIM_MISMATCH If len_control_points(spline) % dim != 0

◆ ts_bspline_set_knots()

tsError ts_bspline_set_knots ( tsBSpline spline,
const tsReal *  knots,
tsStatus status 
)

Sets the knots of spline. Creates a deep copy of knots.

Parameters
[out]splineThe spline whose knots are set.
[in]knotsThe values to deep copy and scale.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_KNOTS_DECR If the knot vector is decreasing.
TS_MULTIPLICITY If there is a knot with multiplicity > order

◆ ts_bspline_set_order()

tsError ts_bspline_set_order ( tsBSpline spline,
size_t  order,
tsStatus status 
)

Sets the order (degree + 1) of spline.

Parameters
[out]splineThe spline whose order is set.
[in]orderThe order to be set.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_DEG_GE_NCTRLP If order > ts_bspline_get_control_points(spline) or if order == 0 ( due to the underflow resulting from: order - 1 => 0 - 1 => INT_MAX which always is >= ts_bspline_get_control_points(spline) ).

◆ ts_bspline_sof_control_points()

size_t ts_bspline_sof_control_points ( const tsBSpline spline)

Returns the size of the control point array of spline. This function may be useful when copying control points using memcpy or memmove.

Parameters
[in]splineThe spline with its control point array whose size is read.
Returns
The size of the control point array of spline.

◆ ts_bspline_sof_knots()

size_t ts_bspline_sof_knots ( const tsBSpline spline)

Returns the size of the knot array of spline. This function may be useful when copying knots using memcpy or memmove.

Parameters
[in]splineThe spline with its knot array whose size is read.
Returns
TS_SUCCESS The size of the knot array of spline.

◆ ts_bspline_split()

tsError ts_bspline_split ( const tsBSpline spline,
tsReal  u,
tsBSpline _split_,
size_t *  _k_,
tsStatus status 
)

Splits spline at knot value u and stores the result in _split_. That is, u is inserted n times such that the multiplicity of u in spline is equal to the spline's order. Creates a deep copy of spline, if spline != _split_.

Parameters
splineThe spline to split.
uThe split point.
_split_The output parameter.
_k_Output parameter. Stores the last index of u in _split_.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_U_UNDEFINED If spline is not defined at knot value u.
TS_MALLOC If allocating memory failed.

◆ ts_bspline_tension()

tsError ts_bspline_tension ( const tsBSpline spline,
tsReal  tension,
tsBSpline out,
tsStatus status 
)

Sets the control points of spline so that their tension corresponds the given tension factor (0 => yields to a line connecting the first and the last control point; 1 => keeps the original shape). If tension < 0 or if tension > 1, the behaviour of this function is undefined, though, it will not result in an error. Creates a deep copy of spline if spline != out.

This function is based on:

 Holten, Danny. "Hierarchical edge bundles: Visualization of adjacency
 relations in hierarchical data." Visualization and Computer Graphics,
 IEEE Transactions on 12.5 (2006): 741-748.

Holten calls it "straightening" (page 744, equation 1).

Parameters
[in]splineThe input spline.
[in]tensionThe tension factor (0 <= tension <= 1).
[out]outThe output spline.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

< The straightening factor.

< Pointer to the control points of out.

< Used in for loops.

◆ ts_bspline_to_beziers()

tsError ts_bspline_to_beziers ( const tsBSpline spline,
tsBSpline _beziers_,
tsStatus status 
)

Subdivides spline into a sequence of Bezier curves by splitting it at each internal knot value. Creates a deep copy of spline, if spline != _beziers_.

Parameters
splineThe spline to subdivide.
_beziers_The output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

< Number of control points to add/remove.

< Index of the split knot value.

< Minimum of the knot values.

< Maximum of the knot values.

< Temporarily stores the result.

< Pointer to the knots of tmp.

< Number of knots in tmp.

◆ ts_bspline_to_json()

tsError ts_bspline_to_json ( const tsBSpline spline,
char **  _json_,
tsStatus status 
)

Serializes spline to a null-terminated JSON string and stores the result in _json_.

Parameters
splineThe spline to serialize.
jsonThe output parameter.
statusOutput parameter. Store the returned error code and a descriptive error message. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_deboornet_copy()

tsError ts_deboornet_copy ( const tsDeBoorNet src,
tsDeBoorNet dest,
tsStatus status 
)

Creates a deep copy of src and stores the copied values in dest. Does nothing, if src == dest.

Parameters
[in]srcThe net to deep copy.
[out]destThe output net.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_deboornet_dimension()

size_t ts_deboornet_dimension ( const tsDeBoorNet net)

Returns the dimension of net. The dimension of a net describes the number of components for each point in ts_bspline_get_points(spline). One-dimensional nets are possible, albeit their benefit might be questionable.

Parameters
[in]netThe net whose dimension is read.
Returns
The dimension of net.

◆ ts_deboornet_free()

void ts_deboornet_free ( tsDeBoorNet net)

Frees the data of net. After calling this function, the data of net points to NULL.

Parameters
[out]netThe net to free.

◆ ts_deboornet_index()

size_t ts_deboornet_index ( const tsDeBoorNet net)

Returns the index [u_k, u_k+1) with u being the knot of net.

Parameters
[in]netThe net whose index is read.
Returns
The index [u_k, u_k+1) with u being the knot of net.

◆ ts_deboornet_init()

tsDeBoorNet ts_deboornet_init ( )

Creates a new net whose data points to NULL.

Returns
A new net whose data points to NULL.

◆ ts_deboornet_knot()

tsReal ts_deboornet_knot ( const tsDeBoorNet net)

Returns the knot (sometimes also called 'u' or 't') of net.

Parameters
[in]netThe net whose knot is read.
Returns
The knot of net.

◆ ts_deboornet_len_points()

size_t ts_deboornet_len_points ( const tsDeBoorNet net)

Returns the length of the point array of net.

Parameters
[in]netThe net with its point array whose length is read.
Returns
The length of the point array of net.

◆ ts_deboornet_len_result()

size_t ts_deboornet_len_result ( const tsDeBoorNet net)

Returns the length of the result array of net.

Parameters
[in]netThe net with its result array whose length is read.
Returns
The length of the result array of net.

◆ ts_deboornet_move()

void ts_deboornet_move ( tsDeBoorNet src,
tsDeBoorNet dest 
)

Moves the ownership of the data of src to dest. After calling this function, the data of src points to NULL. Does not free the data of dest. Does nothing, if src == dest.

Parameters
[out]srcThe net whose values are moved to dest.
[out]destThe net that receives the values of src.

◆ ts_deboornet_multiplicity()

size_t ts_deboornet_multiplicity ( const tsDeBoorNet net)

Returns the multiplicity of the knot of net.

Parameters
[in]netThe net whose multiplicity is read.
Returns
The multiplicity of the knot of net.

◆ ts_deboornet_num_insertions()

size_t ts_deboornet_num_insertions ( const tsDeBoorNet net)

Returns the number of insertion that were necessary to evaluate the knot of net.

Parameters
[in]netThe net with its knot whose number of insertions is read.
Returns
The number of insertions that were necessary to evaluate the knot of net.

◆ ts_deboornet_num_points()

size_t ts_deboornet_num_points ( const tsDeBoorNet net)

Returns the number of points of net.

Parameters
[in]netThe net whose number of points is read.
Returns
The number of points of net.

◆ ts_deboornet_num_result()

size_t ts_deboornet_num_result ( const tsDeBoorNet net)

Returns the number of points in the result array of net (1 <= num_result <= 2).

Parameters
[in]netThe net with its result array whose number of points is read.
Returns
The number of points in the result array of net.

◆ ts_deboornet_points()

tsError ts_deboornet_points ( const tsDeBoorNet net,
tsReal **  points,
tsStatus status 
)

Returns a deep copy of the points of net.

Parameters
[in]netThe net whose points is read.
[out]pointsThe output array.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_deboornet_result()

tsError ts_deboornet_result ( const tsDeBoorNet net,
tsReal **  result,
tsStatus status 
)

Returns a deep copy of the result of net.

Parameters
[in]netThe net whose result is read.
[out]resultThe output array.
[out]statusThe status of this function. May be NULL.
Returns
TS_SUCCESS On success.
TS_MALLOC If allocating memory failed.

◆ ts_deboornet_sof_points()

size_t ts_deboornet_sof_points ( const tsDeBoorNet net)

Returns the size of the point array of net. This function may be useful when copying points using memcpy or memmove.

Parameters
[in]netThe net with its point array whose size is read.
Returns
The size of the point array of net.

◆ ts_deboornet_sof_result()

size_t ts_deboornet_sof_result ( const tsDeBoorNet net)

Returns the size of the result array of net. This function may be useful when copying results using memcpy or memmove.

Parameters
[in]netThe net with its result array whose size is read.
Returns
TS_SUCCESS The size of the result array of net.

◆ ts_distance()

tsReal ts_distance ( const tsReal *  x,
const tsReal *  y,
size_t  dim 
)

Returns the euclidean distance of the points x and y.

Parameters
xThe x value.
yThe y value.
dimThe dimension of x and y.
Returns
The euclidean distance of x and y.

◆ ts_knots_equal()

int ts_knots_equal ( tsReal  x,
tsReal  y 
)

Compares the knots x and y using the epsilon environment TS_KNOT_EPSILON.

Parameters
xThe x value to compare.
yThe y value to compare.
Returns
1 If x is equal to y with respect to the epsilon environment TS_KNOT_EPSILON.
0 Otherwise.