Curve Functions

functions for SRVF curve manipulations

moduleauthor:: J. Derek Tucker <jdtuck@sandia.gov>

curve_functions.Basis_Normal_A(q)[source]

Find Normal Basis

Parameters:q – numpy ndarray (n,T) defining T points on n dimensional SRVF

:rtype list :return delg: basis

curve_functions.calc_j(basis)[source]

Calculates Jacobian matrix from normal basis

Parameters:basis – list of numpy ndarray of shape (2,M) of M samples basis
Return type:numpy ndarray
Return j:Jacobian
curve_functions.calculate_variance(beta)[source]

This function calculates variance of curve beta

Parameters:beta – numpy ndarray of shape (2,M) of M samples
Return type:numpy ndarray
Return variance:
 variance
curve_functions.calculatecentroid(beta)[source]

This function calculates centroid of a parameterized curve

Parameters:beta – numpy ndarray of shape (2,M) of M samples
Return type:numpy ndarray
Return centroid:
 center coordinates
curve_functions.curve_to_q(beta, mode='O')[source]

This function converts curve beta to srvf q

Parameters:
  • beta – numpy ndarray of shape (2,M) of M samples
  • mode – Open (‘O’) or closed curve (‘C’) (default ‘O’)
Return type:

numpy ndarray

Return q:

srvf of curve

Return lenb:

length of curve

Return lenq:

length of srvf

curve_functions.curve_zero_crossing(Y, beta, bt, y_max, y_min, gmax, gmin)[source]

finds zero-crossing of optimal gamma, gam = s*gmax + (1-s)*gmin from elastic curve regression model

Parameters:
  • Y – response
  • beta – predicitve function
  • bt – basis function
  • y_max – maximum repsonse for warping function gmax
  • y_min – minimum response for warping function gmin
  • gmax – max warping function
  • gmin – min warping fucntion
Return type:

numpy array

Return gamma:

optimal warping function

Return O_hat:

rotation matrix

curve_functions.elastic_distance_curve(beta1, beta2, closed=0, rotation=True, scale=False, method='DP')[source]

Calculates the two elastic distances between two curves :param beta1: numpy ndarray of shape (2,M) of M samples :param beta2: numpy ndarray of shape (2,M) of M samples :param closed: open (0) or closed (1) curve (default=0) :param rotation: compute optimal rotation (default=True) :param scale: include scale (default=False) :param method: method to apply optimization (default=”DP”) options are “DP” or “RBFGS”

Return type:tuple
Return dist:shape distance
Return dx:phase distance
curve_functions.elastic_shooting(q1, v)[source]

Calculates shooting vector from v to q1

Parameters:
  • q1 – vector of srvf
  • v – shooting vector

:rtype numpy ndarray :return q2n: vector of srvf

curve_functions.find_basis_normal(q)[source]

Finds the basis normal to the srvf

Parameters:q1 – numpy ndarray of shape (2,M) of M samples
Return type:list of numpy ndarray
Return basis:list containing basis vectors
curve_functions.find_best_rotation(q1, q2, allow_reflection=False, only_xy=False)[source]

This function calculates the best rotation between two srvfs using procustes rigid alignment

Parameters:
  • q1 – numpy ndarray of shape (2,M) of M samples
  • q2 – numpy ndarray of shape (2,M) of M samples
  • allow_reflection – bool indicating if reflection is allowed (i.e. if the determinant of the optimal rotation can be -1)
  • only_xy – bool indicating if rotation should only be allowed in the first two dimensions of the space
Return type:

numpy ndarray

Return q2new:

optimal rotated q2 to q1

Return R:

rotation matrix

curve_functions.find_rotation_and_seed_coord(beta1, beta2, closed=0, rotation=True, method='DP')[source]

This function returns a candidate list of optimally oriented and registered (seed) shapes w.r.t. beta1

Parameters:
  • beta1 – numpy ndarray of shape (2,M) of M samples
  • beta2 – numpy ndarray of shape (2,M) of M samples
  • closed – Open (0) or Closed (1)
  • rotation – find rotation (default=True)
  • method – method to apply optimization (default=”DP”) options are “DP” or “RBFGS”
Return type:

numpy ndarray

Return beta2new:
 

optimal rotated beta2 to beta1

Return O:

rotation matrix

Return tau:

seed

curve_functions.find_rotation_and_seed_q(q1, q2)[source]

This function returns a candidate list of optimally oriented and registered (seed) shapes w.r.t. beta1

Parameters:
  • q1 – numpy ndarray of shape (2,M) of M samples
  • q2 – numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return beta2new:
 

optimal rotated beta2 to beta1

Return O:

rotation matrix

Return tau:

seed

curve_functions.find_rotation_and_seed_unique(q1, q2, closed=0, rotation=True, method='DP')[source]

This function returns a candidate list of optimally oriented and registered (seed) shapes w.r.t. beta1

Parameters:
  • beta1 – numpy ndarray of shape (2,M) of M samples
  • beta2 – numpy ndarray of shape (2,M) of M samples
  • closed – Open (0) or Closed (1)
  • rotation – find rotation (default=True)
  • method – method to apply optimization (default=”DP”) options are “DP” or “RBFGS”
Return type:

numpy ndarray

Return beta2new:
 

optimal rotated beta2 to beta1

Return O:

rotation matrix

Return tau:

seed

curve_functions.gram_schmidt(basis)[source]

Performs Gram Schmidt Orthogonlization of a basis_o

param basis:list of numpy ndarray of shape (2,M) of M samples
rtype:list of numpy ndarray
return basis_o:orthogonlized basis
curve_functions.group_action_by_gamma(q, gamma)[source]

This function reparamerized srvf q by gamma

Parameters:
  • f – numpy ndarray of shape (2,M) of M samples
  • gamma – numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return qn:

reparatermized srvf

curve_functions.group_action_by_gamma_coord(f, gamma)[source]

This function reparamerized curve f by gamma

Parameters:
  • f – numpy ndarray of shape (2,M) of M samples
  • gamma – numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return fn:

reparatermized curve

curve_functions.innerprod_q2(q1, q2)[source]

This function calculates the inner product in srvf space

Parameters:
  • q1 – numpy ndarray of shape (2,M) of M samples
  • q2 – numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return val:

inner product

curve_functions.inverse_exp(q1, q2, beta2)[source]

Calculate the inverse exponential to obtain a shooting vector from q1 to q2 in shape space of open curves

Parameters:
  • q1 – numpy ndarray of shape (2,M) of M samples
  • q2 – numpy ndarray of shape (2,M) of M samples
  • beta2 – numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return v:

shooting vectors

curve_functions.inverse_exp_coord(beta1, beta2, closed=0, method='DP')[source]

Calculate the inverse exponential to obtain a shooting vector from beta1 to beta2 in shape space of open curves

Parameters:
  • beta1 – numpy ndarray of shape (2,M) of M samples
  • beta2 – numpy ndarray of shape (2,M) of M samples
  • closed – open (0) or closed (1) curve
  • method – method to apply optimization (default=”DP”) options are “DP” or “RBFGS”
Return type:

numpy ndarray

Return v:

shooting vectors

Return dist:

distance

curve_functions.optimum_reparam_curve(q1, q2, lam=0.0, method='DP')[source]

calculates the warping to align srsf q2 to q1

Parameters:
  • q1 – matrix of size nxN or array of NxM samples of first SRVF
  • time – vector of size N describing the sample points
  • q2 – matrix of size nxN or array of NxM samples samples of second SRVF
  • lam – controls the amount of elasticity (default = 0.0)
  • method – method to apply optimization (default=”DP”) options are “DP” or “RBFGS”
Return type:

vector

Return gam:

describing the warping function used to align q2 with q1

curve_functions.parallel_translate(w, q1, q2, basis, mode=0)[source]

parallel translates q1 and q2 along manifold

Parameters:
  • w – numpy ndarray of shape (2,M) of M samples
  • q1 – numpy ndarray of shape (2,M) of M samples
  • q2 – numpy ndarray of shape (2,M) of M samples
  • basis – list of numpy ndarray of shape (2,M) of M samples
  • mode – open 0 or closed curves 1 (default 0)
Return type:

numpy ndarray

Return wbar:

translated vector

curve_functions.pre_proc_curve(beta, T=100)[source]

This function prepcoessed a curve beta to set of closed curves

Parameters:
  • beta – numpy ndarray of shape (2,M) of M samples
  • T – number of samples (default = 100)
Return type:

numpy ndarray

Return betanew:

projected beta

Return qnew:

projected srvf

Return A:

alignment matrix (not used currently)

curve_functions.project_curve(q)[source]

This function projects srvf q to set of close curves

Parameters:q – numpy ndarray of shape (2,M) of M samples
Return type:numpy ndarray
Return qproj:project srvf
curve_functions.project_tangent(w, q, basis)[source]

projects srvf to tangent space w using basis

Parameters:
  • w – numpy ndarray of shape (2,M) of M samples
  • q – numpy ndarray of shape (2,M) of M samples
  • basis – list of numpy ndarray of shape (2,M) of M samples
Return type:

numpy ndarray

Return wproj:

projected q

curve_functions.psi(x, a, q)[source]

This function formats variance output

Parameters:
  • x – numpy ndarray of shape (2,M) of M samples curve
  • a – numpy ndarray of shape (2,1) mean
  • q – numpy ndarray of shape (2,M) of M samples srvf
Return type:

numpy ndarray

Return psi1:

variance

Return psi2:

cross variance

Return psi3:

curve end

Return psi4:

curve end

curve_functions.q_to_curve(q, scale=1)[source]

This function converts srvf to beta

Parameters:
  • q – numpy ndarray of shape (n,M) of M samples
  • scale – scale of curve
Return type:

numpy ndarray

Return beta:

parameterized curve

curve_functions.resamplecurve(x, N=100, time=None, mode='O')[source]

This function resamples a curve to have N samples

Parameters:
  • x – numpy ndarray of shape (2,M) of M samples
  • N – Number of samples for new curve (default = 100)
  • time – timing vector (Default=None)
  • mode – Open (‘O’) or closed curve (‘C’) (default ‘O’)
Return type:

numpy ndarray

Return xn:

resampled curve

curve_functions.scale_curve(beta)[source]

scales curve to length 1

Parameters:beta – numpy ndarray of shape (2,M) of M samples
Return type:numpy ndarray
Return beta_scaled:
 scaled curve
Return scale:scale factor used
curve_functions.shift_f(f, tau)[source]

shifts a curve f by tau

Parameters:
  • f – numpy ndarray of shape (2,M) of M samples
  • tau – scalar
Return type:

numpy ndarray

Return fn:

shifted curve