Source code for pyspline.pyCurve

# External modules
import numpy as np
from scipy.sparse import linalg

# Local modules
from . import libspline
from .utils import Error, _assembleMatrix, checkInput, closeTecplot, openTecplot, writeTecplot1D


[docs] class Curve(object): """ Create an instance of a b-spline curve. There are two ways to initialize the class * **Creation**: Create an instance of the spline class directly by supplying the required information. kwargs MUST contain the following information: ``k, t, coef``. * **LMS/Interpolation**: Create an instance of the spline class by using an interpolating spline to given data points or a LMS spline. The following keyword argument information is required: 1. ``k`` Spline Order 2. ``X`` real array size (N, nDim) of data to fit. **OR** 1. ``x`` (1D) and ``s`` for 1D 2. ``x`` (1D) and ``y`` (1D) for 2D spatial curve 3. ``x`` (1D) and ``y``` (1D) and ``z`` 1D for 3D spatial curve 3. ``s`` real array of size (N). Optional for nDim >= 2 Parameters ---------- k : {2, 3, 4} Order for spline. A spline with order :math:`n` has at most :math:`C^{n-2}` continuity. The :math:`n-1` derivative will be piecewise constant. nCtl : int Number of control points. If this is specified then LMS will be used instead of interpolation. t : array, list Knot vector. Used only for creation. Must be size ``nCtl + k`` coef : array, size (nCtl, nDim) Coefficients to use. Only used for creation. The second dimension determine the spatial order of the spline X : arary, size (N, ndim) Full array of data to interpolate/fit x : array, size (N) Just x coordinates of data to fit y : array, size (N) Just y coordinates of data to fit z : array, size (N) Just z coordinates of data to fit s : array, size(N) Optional parameter to use. Not required for nDim >=2 nIter : int The number of Hoscheks parameter correction iterations to run. Only used for LMS fitting. weight : array, size(N) A array of weighting factors for each fitting point. A value of -1 can be used to exactly constrain a point. By default, all weights are 1.0 deriv : array, size (len(derivPtr), nDim) Array containing derivative information the user wants to use at a particular point. **EXPERIMENTAL** derivPtr : int, array, size (nDerivPtr) Array of indices pointing to the index of points in X (or x,y,z), for which the user has supplied a derivative for in ``deriv``. **EXPERIMENTAL** derivWeights : array size(nDerivPtr) Optional array of the weighting to use for derivatives. A value of -1 can be used to exactly constrain a derivative. **EXPERIMENTAL** Examples -------- >>> x = [0, 0.5, 1.0] >>> y = [0, 0.25, 1.0] >>> s = [0., 0.5, 1.0] >>> # Spatial interpolated seg (k=2 makes this a straight line) >>> line_seg = Curve(x=x, y=y, k=2) >>> # With explicit parameter values >>> line_seg = Curve(x=x, y=y, k=2, s=s) >>> #Interpolate parabolic curve >>> parabola = Curve(x=x, y=y, k=3) >>> #Interpolate parabolic curve with parameter values >>> parabola = Curve(x=x, y=y, k=3, s=s) """ def __init__(self, **kwargs): self.length = None self.gpts = None self.data = None self.sdata = None self.localInterp = False # We have provided information to create curve directly if "k" in kwargs and "t" in kwargs and "coef" in kwargs: self.s = None self.X = None self.N = None self.k = checkInput(kwargs["k"], "k", int, 0) self.coef = np.atleast_2d(kwargs["coef"]) self.nCtl = self.coef.shape[0] self.nDim = self.coef.shape[1] self.t = checkInput(kwargs["t"], "t", float, 1, (self.nCtl + self.k,)) self.origData = False self.calcGrevillePoints() elif "localInterp" in kwargs: # Local, non-global interpolation. We could use this for # second-order (linear interpolation) but there is not # much point, since it would be the same as 'interp' # below. self.localInterp = True self.k = 4 self.origData = True if "X" in kwargs: self.X = np.atleast_2d(kwargs["X"]) if np.ndim(kwargs["X"]) == 1: self.X = np.transpose(self.X) elif "x" in kwargs and "y" in kwargs and "z" in kwargs: self.X = np.vstack([kwargs["x"], kwargs["y"], kwargs["z"]]).T elif "x" in kwargs and "y" in kwargs: self.X = np.vstack([kwargs["x"], kwargs["y"]]).T elif "x" in kwargs: self.X = np.transpose(np.atleast_2d(kwargs["x"])) # enf if self.nDim = self.X.shape[1] self.N = self.X.shape[0] if "s" in kwargs: self.s = checkInput(kwargs["s"], "s", float, 1, (self.N,)) self.length = 1.0 # Assume it is 0->1 else: if self.nDim <= 1: raise Error("For 1D splines, the basis, 's' must be given.") self._getParameterization() # Now we have the data we need to generate the local # interpolation T = np.zeros((self.N, self.nDim)) # Compute tangents qq = np.zeros_like(self.X) T = np.zeros_like(self.X) deltaS = np.zeros(self.N) for i in range(1, self.N): deltaS[i] = self.s[i] - self.s[i - 1] qq[i, :] = self.X[i] - self.X[i - 1] for i in range(1, self.N - 1): a = deltaS[i] / (deltaS[i] + deltaS[i + 1]) T[i] = (1 - a) * qq[i] + a * qq[i + 1] # Do the start and end points: (eqn: 9.32, The NURBS book) T[0] = 2 * qq[1] / deltaS[1] - T[1] T[-1] = 2 * qq[-1] / deltaS[-1] - T[-2] # Normalize for i in range(self.N): T[i] = T[i] / np.linalg.norm(T[i]) # Final coefficients and t self.coef = np.zeros((2 * (self.N - 1) + 2, self.nDim)) self.t = np.zeros(len(self.coef) + self.k) self.nCtl = len(self.coef) # End Pts self.coef[0] = self.X[0].copy() self.coef[-1] = self.X[-1].copy() # Interior coefficients for i in range(self.N - 1): a = self.length * (self.s[i + 1] - self.s[i]) self.coef[2 * i + 1] = self.X[i] + a / 3.0 * T[i] self.coef[2 * i + 2] = self.X[i + 1] - a / 3.0 * T[i + 1] # Knots self.t[-4:] = 1.0 u = np.zeros(self.N) for i in range(0, self.N - 1): u[i + 1] = u[i] + np.linalg.norm(self.coef[2 * i + 2] - self.coef[2 * i + 1]) for i in range(1, self.N - 1): self.t[2 * i + 2] = u[i] / u[self.N - 1] self.t[2 * i + 3] = u[i] / u[self.N - 1] else: # lms or interpolate function if not ("k" in kwargs and ("X" in kwargs or "x" in kwargs)): raise ValueError( "At least spline order, k and X (or x=, y=) MUST be defined for (interpolation) spline creation. " + "nCtl=<number of control points> must be specified for a LMS fit" ) self.origData = True if "X" in kwargs: self.X = np.atleast_2d(kwargs["X"]) if np.ndim(kwargs["X"]) == 1: self.X = np.transpose(self.X) elif "x" in kwargs and "y" in kwargs and "z" in kwargs: self.X = np.vstack([kwargs["x"], kwargs["y"], kwargs["z"]]).T elif "x" in kwargs and "y" in kwargs: self.X = np.vstack([kwargs["x"], kwargs["y"]]).T elif "x" in kwargs: self.X = np.transpose(np.atleast_2d(kwargs["x"])) # enf if self.nDim = self.X.shape[1] self.N = self.X.shape[0] self.k = checkInput(kwargs["k"], "k", int, 0) self.t = None if "nIter" in kwargs: self.nIter = checkInput(kwargs["nIter"], "nIter", int, 0) else: self.nIter = 1 if "s" in kwargs: self.s = checkInput(kwargs["s"], "s", float, 1, (self.N,)) else: if self.nDim <= 1: raise Error("For 1D splines, the basis, 's' must be given.") self._getParameterization() if "weights" in kwargs: self.weights = checkInput(kwargs["weights"], "weights", float, 1, (self.N,)) else: self.weights = np.ones(self.N) if "deriv" in kwargs and "derivPtr" in kwargs: self.deriv = checkInput(kwargs["deriv"], "deriv", float, 2) self.derivPtr = checkInput(kwargs["derivPtr"], "derivPtr", int, 1, (len(self.deriv),)) else: self.deriv = None self.derivPtr = np.array([]) if "derivWeights" in kwargs and self.deriv is not None: self.derivWeights = checkInput(kwargs["derivWeights"], "derivWeights", float, 1, (len(self.derivPtr),)) else: if self.deriv is not None: self.derivWeights = np.ones(len(self.deriv)) else: self.derivWeights = None # We are doing an interpolation...set all weights to 1 if "nCtl" not in kwargs: self.interp = True self.weights[:] = 1 if self.derivWeights is not None: self.derivWeights[:] = 1 else: self.interp = False self.nCtl = checkInput(kwargs["nCtl"], "nCtl", int, 0) self.recompute(self.nIter, computeKnots=True)
[docs] def recompute(self, nIter, computeKnots=True): """ Run iterations of Hoscheks Parameter Correction on the current curve Parameters ---------- nIter : int The number of parameter correction iterations to run computeKnots : bool Flag whether or not the knots should be recomputed """ # Return if we don't have original data to fit if not self.origData or self.localInterp: return # Do the separation between the constrained and unconstrained: # u -> unconstrained # s -> constrained suSelect = np.where(self.weights > 0.0) scSelect = np.where(self.weights <= 0.0) S = self.X[suSelect] su = self.s[suSelect] T = self.X[scSelect] sc = self.s[scSelect] weights = self.weights[np.where(self.weights > 0.0)] nu = len(S) nc = len(T) # And the derivative info if self.deriv is not None: sduSelect = np.where(self.derivWeights > 0.0) sdcSelect = np.where(self.derivWeights <= 0.0) S = np.vstack((S, self.deriv[sduSelect])) sdu = self.s[self.derivPtr][sduSelect] T = np.vstack((T, self.deriv[sdcSelect])) sdc = self.s[self.derivPtr][sdcSelect] weights = np.append(weights, self.derivWeights[np.where(self.derivWeights > 0.0)]) ndu = len(sdu) ndc = len(sdc) else: sdu = np.array([], "d") sdc = np.array([], "d") ndu = 0 ndc = 0 W = _assembleMatrix(weights, np.arange(len(weights)), np.arange(len(weights) + 1), (len(weights), len(weights))) if self.interp: self.nCtl = nu + nc + ndu + ndc self.nIter = 1 # Sanity check to make sure k is ok if nu + nc + ndu + ndc < self.k: self.k = nu + nc + ndu + ndc if computeKnots: # Generate the knot vector, if necessary greville points and # empty coefficients if self.interp: self.t = libspline.knots_interp(self.s, self.derivPtr, self.k) else: self.t = libspline.knots_lms(self.s, self.nCtl, self.k) self.calcGrevillePoints() self.coef = np.zeros((self.nCtl, self.nDim), "d") # Get the 'N' jacobian nVals = np.zeros((nu + ndu) * self.k) # | nRowPtr = np.zeros(nu + ndu + 1, "intc") # | -> CSR formulation nColInd = np.zeros((nu + ndu) * self.k, "intc") # | libspline.curve_jacobian_wrap(su, sdu, self.t, self.k, self.nCtl, nVals, nRowPtr, nColInd) N = _assembleMatrix(nVals, nColInd, nRowPtr, (nu + ndu, self.nCtl)).tocsc() if self.interp: # Factorize once for efficiency solve = linalg.factorized(N) for idim in range(self.nDim): self.coef[:, idim] = solve(S[:, idim]) return # If we do NOT have an interpolation: length = libspline.poly_length(self.X.T) for _i in range(nIter): su = self.s[suSelect] sc = self.s[scSelect] if self.deriv is not None: sdu = self.s[self.derivPtr][sduSelect] sdc = self.s[self.derivPtr][sdcSelect] libspline.curve_jacobian_wrap(su, sdu, self.t, self.k, self.nCtl, nVals, nRowPtr, nColInd) NTWN = (N.transpose() * W * N).tocsc() # We need this either way if nc + ndc == 0: # We are doing LMS but no # constraints...just a straight weighted # LMS # Factorize once for efficiency solve = linalg.factorized(NTWN) for idim in range(self.nDim): self.coef[:, idim] = solve(N.transpose() * W * S[:, idim]) else: # Now its more complicated since we have # constraints --only works with scipy Sparse # matrices mVals = np.zeros((nc + ndc) * self.k) # | mRowPtr = np.zeros(nc + ndc + 1, "intc") # | -> CSR mColInd = np.zeros((nc + ndc) * self.k, "intc") # | libspline.curve_jacobian_wrap(sc, sdc, self.t, self.k, self.nCtl, mVals, mRowPtr, mColInd) M = _assembleMatrix(mVals, mColInd, mRowPtr, (nc + ndc, self.nCtl)) # Now we must assemble the constrained jacobian # [ N^T*W*T M^T][P] = [ N^T*W*S] # [ M 0 ][R] [ T ] MT = M.transpose().tocsr() jVal, jColInd, jRowPtr = libspline.constr_jac( NTWN.data, NTWN.indptr, NTWN.indices, MT.data, MT.indptr, MT.indices, M.data, M.indptr, M.indices, self.nCtl, ) # Create sparse csr matrix and factorize J = _assembleMatrix(jVal, jColInd, jRowPtr, (self.nCtl + nc + ndc, self.nCtl + nc + ndc)) # Factorize once for efficiency solve = linalg.factorized(J) for idim in range(self.nDim): rhs = np.hstack((N.transpose() * W * S[:, idim], T[:, idim])) self.coef[:, idim] = solve(rhs)[0 : self.nCtl] # end if (constr - not constrained # Run para correction libspline.curve_para_corr(self.t, self.k, self.s, self.coef.T, length, self.X.T) # end for (iter loop) # Check the RMS rms = 0.0 for idim in range(self.nDim): rms += np.linalg.norm(N * self.coef[:, idim] - S[:, idim]) ** 2 rms = np.sqrt(rms / self.N)
def _getParameterization(self): """Compute a parametrization for the curve based on an arc-length formulation """ self.s = np.zeros(self.N, "d") for i in range(self.N - 1): dist = 0 for idim in range(self.nDim): dist += (self.X[i + 1, idim] - self.X[i, idim]) ** 2 self.s[i + 1] = self.s[i] + np.sqrt(dist) self.length = self.s[-1] self.s = self.s / self.s[-1]
[docs] def reverse(self): """ Reverse the direction of this curve """ self.coef = self.coef[::-1, :] self.t = 1 - self.t[::-1]
[docs] def insertKnot(self, u, r): """ Insert at knot in the curve at parametric position u Parameters ---------- u : float Parametric position to split at r : int Number of time to insert. Should be > 0. Returns ------- actualR : int The number of times the knot was **actually** inserted breakPt : int Index in the knot vector of the new knot(s) """ u = checkInput(u, "u", float, 0) r = checkInput(r, "r", int, 0) if u <= 0: return if u >= 1.0: return actualR, tNew, coefNew, breakPt = libspline.insertknot(u, r, self.t, self.k, self.coef.T) self.t = tNew[0 : self.nCtl + self.k + actualR] self.coef = coefNew[:, 0 : self.nCtl + actualR].T self.nCtl = self.nCtl + actualR # break_pt is converted to zero based ordering here!!! return actualR, breakPt - 1
[docs] def splitCurve(self, u): """ Split the curve at parametric position u. This uses the :func:`insertKnot` routine Parameters ---------- u : float Parametric position to insert knot. Returns ------- curve1 : pySpline curve object Curve from s=[0, u] curve2 : pySpline curve object Curve from s=[u, 1] Notes ----- curve1 and curve2 may be None if the parameter u is outside the range of (0, 1) """ u = checkInput(u, "u", float, 0) if u <= 0.0: return None, Curve(t=self.t.copy(), k=self.k, coef=self.coef.copy()) if u >= 1.0: return Curve(t=self.t.copy(), k=self.k, coef=self.coef.copy()), None r, breakPt = self.insertKnot(u, self.k - 1) # Break point is now at the right so we need to adjust the # counter to the left breakPt = breakPt - self.k + 2 # Process knot vectors: uu = self.t[breakPt] t1 = np.hstack((self.t[0 : breakPt + self.k - 1].copy(), uu)) / uu t2 = (np.hstack((uu, self.t[breakPt:].copy())) - uu) / (1.0 - uu) coef1 = self.coef[0:breakPt, :].copy() coef2 = self.coef[breakPt - 1 :, :].copy() return Curve(t=t1, k=self.k, coef=coef1), Curve(t=t2, k=self.k, coef=coef2)
[docs] def windowCurve(self, uLow, uHigh): """ Compute the segment of the curve between the two parameter values Parameters ---------- uLow : float Lower bound for the clip uHigh : float Upper bound for the clip """ __, c = self.splitCurve(uLow) c, __ = c.splitCurve((uHigh - uLow) / (1.0 - uLow)) return c
[docs] def getLength(self): """Compute the length of the curve using the Euclidean Norm Returns ------- length : float Physical length of curve """ points = self.getValue(self.gpts) length = 0 for i in range(len(points) - 1): length += np.linalg.norm(points[i] - points[i + 1]) return length
[docs] def calcGrevillePoints(self): """Calculate the Greville points""" self.gpts = np.zeros(self.nCtl) for i in range(self.nCtl): for n in range(self.k - 1): # degree loop self.gpts[i] += self.t[i + n + 1] self.gpts[i] = self.gpts[i] / (self.k - 1)
[docs] def calcInterpolatedGrevillePoints(self): """ Compute greville points, but with additional interpolated knots """ self.calcGrevillePoints() s = [self.gpts[0]] N = 2 for i in range(len(self.gpts) - 1): for j in range(N): s.append((self.gpts[i + 1] - self.gpts[i]) * (j + 1) / (N + 1) + self.gpts[i]) s.append(self.gpts[i + 1]) self.sdata = np.array(s)
def __call__(self, s): """ Equivalent to getValue() """ return self.getValue(s)
[docs] def getValue(self, s): """ Evaluate the spline at parametric position, s Parameters ---------- s : float or array Parametric position(s) at which to evaluate the curve. Returns ------- values : array of size nDim or an array of size (N, 3) The evaluated points. If a single scalar s was given, the result with be an array of length nDim (or a scalar if nDim=1). If a vector of s values were given it will be an array of size (N, 3) (or size (N) if ndim=1) """ s = np.array(s).T if self.coef.dtype == np.dtype("d"): vals = libspline.eval_curve(np.atleast_1d(s), self.t, self.k, self.coef.T) else: vals = libspline.eval_curve_c(np.atleast_1d(s).astype("D"), self.t, self.k, self.coef.T) return vals.squeeze().T
[docs] def getDerivative(self, s): """ Evaluate the derivative of the spline at parametric position, s Parameters ---------- s : float Parametric location for derivative Returns ------- ds : array The first derivative. This is an array of size nDim """ if self.coef.dtype == np.dtype("d"): ds = libspline.eval_curve_deriv(s, self.t, self.k, self.coef.T).squeeze() else: ds = libspline.eval_curve_deriv_c(s, self.t, self.k, self.coef.T).squeeze() return ds
[docs] def getSecondDerivative(self, s): """ Evaluate the second derivative of the spline at parametric position, s Parameters ---------- s : float Parametric location for derivative Returns ------- d2s : array The second derivative. This is an array of size nDim """ if self.coef.dtype == np.dtype("d"): d2s = libspline.eval_curve_deriv2(s, self.t, self.k, self.coef.T).squeeze() else: d2s = libspline.eval_curve_deriv2_c(s, self.t, self.k, self.coef.T).squeeze() return d2s
[docs] def projectPoint(self, x0, nIter=25, eps=1e-10, **kwargs): """ Perform a point inversion algorithm. Attempt to find the closest parameter values to the given points x0. Parameters ---------- x0 : array A point or list of points in nDim space for which the minimum distance to the curve is sought. nIter : int Maximum number of Newton iterations to perform. eps : float Desired parameter tolerance. Returns ------- s : float or array Solution to the point inversion. s are the parametric locations yielding the minimum distance to points x0 D : float or array Physical distances between the points and the curve. This is simply ||curve(s) - X0||_2. """ x0 = np.atleast_2d(x0) if "s" in kwargs: s = np.atleast_1d(kwargs["s"]) else: s = -1 * np.ones(len(x0)) if len(x0) != len(s): raise Error("projectPoint: The length of x0 and s must be the same") # If necessary get brute-force starting point if np.any(s < 0) or np.any(s > 1): self.computeData() s = libspline.point_curve_start(x0.T, self.sdata, self.data.T) D = np.zeros_like(x0) for i in range(len(x0)): s[i], D[i] = libspline.point_curve(x0[i], self.t, self.k, self.coef.T, nIter, eps, s[i]) return s.squeeze(), D.squeeze()
[docs] def projectCurve(self, inCurve, nIter=25, eps=1e-10, **kwargs): """ Find the minimum distance between this curve (self) and a second curve passed in (inCurve) Parameters ---------- inCurve : pySpline.curve object Other curve to use nIter : int Maximum number of Newton iterations to perform. eps : float Desired parameter tolerance. s : float Initial guess for curve1 (this curve class) t : float Initial guess for inCurve (curve passed in ) Returns ------- s : float Parametric position on curve1 (this class) t : float Parametric position on curve2 (inCurve) D : float Minimum distance between this curve and inCurve. It is equivalent to ||self(s) - inCurve(t)||_2. """ s = -1 t = -1 if "s" in kwargs: s = checkInput(kwargs["s"], "s", float, 0) if "t" in kwargs: t = checkInput(kwargs["t"], "t", float, 0) eps = checkInput(eps, "eps", float, 0) if s < 0 or s > 1 or t < 0 or t > 1: self.computeData() inCurve.computeData() s, t = libspline.curve_curve_start(self.data.T, self.sdata, inCurve.data.T, inCurve.sdata) s, t, Diff = libspline.curve_curve( self.t, self.k, self.coef.T, inCurve.t, inCurve.k, inCurve.coef.T, nIter, eps, s, t ) return s, t, Diff
[docs] def projectCurveMultiSol(self, inCurve, nIter=25, eps=1e-10, **kwargs): """ Find the minimum distance between this curve (self) and a second curve passed in (inCurve). Try to find as many local solutions as possible Parameters ---------- inCurve : pySpline.curve object Other curve to use nIter : int Maximum number of Newton iterations to perform. eps : float Desired parameter tolerance. s : float Initial guess for curve1 (this curve class) t : float Initial guess for inCurve (curve passed in ) Returns ------- s : array Parametric position(s) on curve1 (this class) t : float Parametric position)s_ on curve2 (inCurve) D : float Minimum distance(s) between this curve and inCurve. It is equivalent to ||self(s) - inCurve(t)||_2. """ s = -1 t = -1 if "s" in kwargs: s = checkInput(kwargs["s"], "s", float, 0) if "t" in kwargs: t = checkInput(kwargs["t"], "t", float, 0) eps = checkInput(eps, "eps", float, 0) self.computeData() inCurve.computeData() uSol = [] tSol = [] diff = [] for i in range(len(self.sdata)): for j in range(len(inCurve.sdata)): s, t, Diff = libspline.curve_curve( self.t, self.k, self.coef.T, inCurve.t, inCurve.k, inCurve.coef.T, nIter, eps, self.sdata[i], inCurve.sdata[j], ) if np.linalg.norm(Diff) < eps: # Its a solution. Check it it is already in list: if len(uSol) == 0: uSol.append(s) tSol.append(t) diff.append(Diff) for ii in range(len(uSol)): if abs(uSol[ii] - s) < eps and abs(tSol[ii] - t) < eps: pass else: uSol.append(s) tSol.append(t) diff.append(Diff) return np.array(uSol), np.array(tSol), np.array(diff)
[docs] def computeData(self, recompute=False): """ Compute discrete data that is used for the Tecplot Visualization as well as the data for doing the brute-force checks Parameters ---------- recompute : bool If True, recompute the data even if it has already been computed. """ # We will base the data on interpolated greville points if self.data is None or recompute: self.calcInterpolatedGrevillePoints() self.data = self.getValue(self.sdata)
[docs] def writeTecplot(self, fileName, curve=True, coef=True, orig=True): """ Write the curve to a tecplot .dat file Parameters ---------- fileName : str File name for tecplot file. Should have .dat extension curve : bool Flag to write discrete approximation of the actual curve coef : bool Flag to write b-spline coefficients orig : bool Flag to write original data (used for fitting) if it exists """ f = openTecplot(fileName, self.nDim) if curve: self.computeData() writeTecplot1D(f, "interpolated", self.data) if coef: writeTecplot1D(f, "control_pts", self.coef) if orig and self.origData: writeTecplot1D(f, "orig_data", self.X) closeTecplot(f)
[docs] def writeIGES_directory(self, handle, Dcount, Pcount, twoD=False): """ Write the IGES file header information (Directory Entry Section) for this curve. DO NOT MODIFY ANYTHING HERE UNLESS YOU KNOW **EXACTLY** WHAT YOU ARE DOING! """ if self.nDim != 3: raise Error("Must have 3 dimensions to write to IGES file") paraEntries = 6 + len(self.t) + self.nCtl + 3 * self.nCtl + 5 paraLines = (paraEntries - 11) // 3 + 2 if np.mod(paraEntries - 11, 3) != 0: paraLines += 1 if twoD: handle.write(" 126%8d 0 0 1 0 0 001010501D%7d\n" % (Pcount, Dcount)) handle.write( " 126 0 2%8d 0 0D%7d\n" % (paraLines, Dcount + 1) ) else: handle.write(" 126%8d 0 0 1 0 0 000000001D%7d\n" % (Pcount, Dcount)) handle.write( " 126 0 2%8d 0 0D%7d\n" % (paraLines, Dcount + 1) ) Dcount += 2 Pcount += paraLines return Pcount, Dcount
[docs] def writeIGES_parameters(self, handle, Pcount, counter): """Write the IGES parameter information for this curve. DO NOT MODIFY ANYTHING HERE UNLESS YOU KNOW **EXACTLY** WHAT YOU ARE DOING! """ handle.write( "%10d,%10d,%10d,0,0,0,0, %7dP%7d\n" % (126, self.nCtl - 1, self.k - 1, Pcount, counter) ) counter += 1 pos_counter = 0 for i in range(len(self.t)): pos_counter += 1 handle.write("%20.12g," % (np.real(self.t[i]))) if np.mod(pos_counter, 3) == 0: handle.write(" %7dP%7d\n" % (Pcount, counter)) counter += 1 pos_counter = 0 for _i in range(self.nCtl): pos_counter += 1 handle.write("%20.12g," % (1.0)) if np.mod(pos_counter, 3) == 0: handle.write(" %7dP%7d\n" % (Pcount, counter)) counter += 1 pos_counter = 0 for i in range(self.nCtl): for idim in range(3): pos_counter += 1 handle.write("%20.12g," % (np.real(self.coef[i, idim]))) if np.mod(pos_counter, 3) == 0: handle.write(" %7dP%7d\n" % (Pcount, counter)) counter += 1 pos_counter = 0 if pos_counter == 1: handle.write("%s %7dP%7d\n" % (" " * 40, Pcount, counter)) counter += 1 elif pos_counter == 2: handle.write("%s %7dP%7d\n" % (" " * 20, Pcount, counter)) counter += 1 # Ouput the ranges handle.write("%20.12g,%20.12g,0.0,0.0,0.0; " % (np.min(self.t), np.max(self.t))) handle.write(" %7dP%7d\n" % (Pcount, counter)) counter += 1 Pcount += 2 return Pcount, counter