/****************************************************************************** * Cagd_aux.c - auxiliary routine to interface to different free from types. * ******************************************************************************* * (C) Gershon Elber, Technion, Israel Institute of Technology * ******************************************************************************* * Written by Gershon Elber, July. 90. * ******************************************************************************/ #include "cagd_loc.h" #define VEC_FIELD_TRIES 10 #define VEC_FIELD_START_STEP 1e-6 /***************************************************************************** * DESCRIPTION: M * Returns the parametric domain of a curve. M * * * PARAMETERS: M * Crv: To get its parametric domain. M * TMin: Where to put the minimal domain's boundary. M * TMax: Where to put the maximal domain's boundary. M * * * RETURN VALUE: M * void M * * * SEE ALSO: M * BspCrvDomain M * * * KEYWORDS: M * CagdCrvDomain, domain, parametric domain M *****************************************************************************/ void CagdCrvDomain(CagdCrvStruct *Crv, CagdRType *TMin, CagdRType *TMax) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: case CAGD_CPOWER_TYPE: *TMin = 0.0; *TMax = 1.0; break; case CAGD_CBSPLINE_TYPE: BspCrvDomain(Crv, TMin, TMax); break; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); break; } } /***************************************************************************** * DESCRIPTION: M * Affinely reset the parametric domain of a curve, in place. M * * * PARAMETERS: M * Crv: To reset its parametric domain. M * TMin: Minimal domain's new boundary. M * TMax: Maximal domain's new boundary. M * * * RETURN VALUE: M * CagdCrvStruct *: Modified curve, in place. M * * * SEE ALSO: M * BspCrvDomain, BspKnotAffineTrans2, CagdSrfSetDomain M * * * KEYWORDS: M * CagdCrvSetDomain, domain, parametric domain M *****************************************************************************/ CagdCrvStruct *CagdCrvSetDomain(CagdCrvStruct *Crv, CagdRType TMin, CagdRType TMax) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: if (APX_EQ(TMin, 0.0) && APX_EQ(TMax, 1.0)) { return Crv; } else { /* Convert to a B-spline curve. */ Crv -> Order = Crv -> Length; Crv -> KnotVector = BspKnotUniformOpen(Crv -> Length, Crv -> Order, NULL); Crv -> GType = CAGD_CBSPLINE_TYPE; } case CAGD_CBSPLINE_TYPE: BspKnotAffineTrans2(Crv -> KnotVector, Crv -> Order + Crv -> Length, TMin, TMax); break; case CAGD_CPOWER_TYPE: default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); break; } return Crv; } /***************************************************************************** * DESCRIPTION: M * Given a curve and parameter value t, evaluate the curve at t. M * * * PARAMETERS: M * Crv: To evaluate at the given parametric location t. M * t: The parameter value at which the curve Crv is to be evaluated. M * * * RETURN VALUE: M * CagdRType *: A vector holding all the coefficients of all components M * of curve Crv's point type. If for example the curve's M * point type is P2, the W, X, and Y will be saved in the M * first three locations of the returned vector. The first M * location (index 0) of the returned vector is reserved for M * the rational coefficient W and XYZ always starts at second M * location of the returned vector (index 1). M * * * SEE ALSO: M * BzrCrvEvalAtParam, BspCrvEvalAtParam, BzrCrvEvalVecAtParam, M * BspCrvEvalVecAtParam, BspCrvEvalCoxDeBoor, CagdCrvEvalToPolyline M * * * KEYWORDS: M * CagdCrvEval, evaluation M *****************************************************************************/ CagdRType *CagdCrvEval(CagdCrvStruct *Crv, CagdRType t) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvEvalAtParam(Crv, t); case CAGD_CBSPLINE_TYPE: return BspCrvEvalAtParam(Crv, t); case CAGD_CPOWER_TYPE: return PwrCrvEvalAtParam(Crv, t); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns the parametric domain of a surface. M * * * PARAMETERS: M * Srf: To get its parametric domain. M * UMin: Where to put the minimal U domain's boundary. M * UMax: Where to put the maximal U domain's boundary. M * VMin: Where to put the minimal V domain's boundary. M * VMax: Where to put the maximal V domain's boundary. M * * * RETURN VALUE: M * void M * * * SEE ALSO: M * BspSrfDomain M * * * KEYWORDS: M * CagdSrfDomain, domain, parametric domain M *****************************************************************************/ void CagdSrfDomain(CagdSrfStruct *Srf, CagdRType *UMin, CagdRType *UMax, CagdRType *VMin, CagdRType *VMax) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: case CAGD_SPOWER_TYPE: *UMin = 0.0; *UMax = 1.0; *VMin = 0.0; *VMax = 1.0; break; case CAGD_SBSPLINE_TYPE: BspSrfDomain(Srf, UMin, UMax, VMin, VMax); break; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); break; } } /***************************************************************************** * DESCRIPTION: M * Affinely reset the parametric domain of a surface, in place. M * * * PARAMETERS: M * Srf: To reset its parametric domain. M * UMin: Minimal domain's new U boundary. M * UMax: Maximal domain's new U boundary. M * VMin: Minimal domain's new V boundary. M * VMax: Maximal domain's new V boundary. M * * * RETURN VALUE: M * CagdSrfStruct *: Modified surface, in place. M * * * SEE ALSO: M * BspSrfDomain, BspKnotAffineTrans2, CagdCrvSetDomain M * * * KEYWORDS: M * CagdSrfSetDomain, domain, parametric domain M *****************************************************************************/ CagdSrfStruct *CagdSrfSetDomain(CagdSrfStruct *Srf, CagdRType UMin, CagdRType UMax, CagdRType VMin, CagdRType VMax) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: if (APX_EQ(UMin, 0.0) && APX_EQ(UMax, 1.0) && APX_EQ(VMin, 0.0) && APX_EQ(VMax, 1.0)) { return Srf; } else { /* Convert to a B-spline surface. */ Srf -> UOrder = Srf -> ULength; Srf -> VOrder = Srf -> VLength; Srf -> UKnotVector = BspKnotUniformOpen(Srf -> ULength, Srf -> UOrder, NULL); Srf -> VKnotVector = BspKnotUniformOpen(Srf -> VLength, Srf -> VOrder, NULL); Srf -> GType = CAGD_SBSPLINE_TYPE; } case CAGD_SBSPLINE_TYPE: BspKnotAffineTrans2(Srf -> UKnotVector, Srf -> UOrder + Srf -> ULength, UMin, UMax); BspKnotAffineTrans2(Srf -> VKnotVector, Srf -> VOrder + Srf -> VLength, VMin, VMax); break; case CAGD_SPOWER_TYPE: default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); break; } return Srf; } /***************************************************************************** * DESCRIPTION: M * Given a surface and parameter values u, v, evaluate the surface at (u, v). M * * * PARAMETERS: M * Srf: To evaluate at the given parametric location (u, v). M * u, v: The parameter values at which the curve Crv is to be evaluated. M * * * RETURN VALUE: M * CagdRType *: A vector holding all the coefficients of all components M * of surface Srf's point type. If for example the surface's M * point type is P2, the W, X, and Y will be saved in the M * first three locations of the returned vector. The first M * location (index 0) of the returned vector is reserved for M * the rational coefficient W and XYZ always starts at second M * location of the returned vector (index 1). M * * * SEE ALSO: M * CagdCrvEval, BspSrfEvalAtParam, BzrSrfEvalAtParam, M * BspSrfEvalAtParam2, TrimSrfEval M * * * KEYWORDS: M * CagdSrfEval, evaluation M *****************************************************************************/ CagdRType *CagdSrfEval(CagdSrfStruct *Srf, CagdRType u, CagdRType v) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfEvalAtParam(Srf, u, v); case CAGD_SBSPLINE_TYPE: return BspSrfEvalAtParam(Srf, u, v); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Routine to convert a single freeform surface to set of triangles M * approximating it. FineNess is a fineness control on result and the larger M t is more triangles may result. A value of 10 is a good start value. M * NULL is returned in case of an error, otherwise list of CagdPolygonStruct. M * * * PARAMETERS: M * Srf: To approximate into triangles. M * FineNess: Control on accuracy, the higher the finer. M * ComputeNormals: If TRUE, normal information is also computed. M * FourPerFlat: If TRUE, four triangles are created per flat surface. M * If FALSE, only 2 triangles are created. M * ComputeUV: If TRUE, UV values are stored and returned as well. M * * * RETURN VALUE: M * CagdPolygonStruct *: A list of polygons with optional normal and/or M * UV parametric information. M * NULL is returned in case of an error. M * * * SEE ALSO: M * BzrSrf2Polygons, BspSrf2Polygons, CagdCrv2Polyline, CagdSrf2Polylines, M * CagdSrf2PolygonStrip, CagdSrf2PolygonsN M * * * KEYWORDS: M * CagdSrf2Polygons, evaluation, polygonal approximation M *****************************************************************************/ CagdPolygonStruct *CagdSrf2Polygons(CagdSrfStruct *Srf, int FineNess, CagdBType ComputeNormals, CagdBType FourPerFlat, CagdBType ComputeUV) { CagdPolygonStruct *Pls; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrf2Polygons(Srf, FineNess, ComputeNormals, FourPerFlat, ComputeUV); case CAGD_SBSPLINE_TYPE: return BspSrf2Polygons(Srf, FineNess, ComputeNormals, FourPerFlat, ComputeUV); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); Pls = BzrSrf2Polygons(Srf, FineNess, ComputeNormals, FourPerFlat, ComputeUV); CagdSrfFree(Srf); return Pls; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Routine to convert a single freeform surface to set of triangles M * approximating it using a uniform fixed resolution of Nu x Nv. M * NULL is returned in case of an error, otherwise list of CagdPolygonStruct. M * * * PARAMETERS: M * Srf: To approximate into triangles. M * Nu, Nv: The number of uniform samples in U and V of surface. M * ComputeNormals: If TRUE, normal information is also computed. M * FourPerFlat: If TRUE, four triangles are created per flat surface. M * If FALSE, only 2 triangles are created. M * ComputeUV: If TRUE, UV values are stored and returned as well. M * * * RETURN VALUE: M * CagdPolygonStruct *: A list of polygons with optional normal and/or M * UV parametric information. M * NULL is returned in case of an error. M * * * SEE ALSO: M * BzrSrf2Polygons, BspSrf2Polygons, CagdCrv2Polyline, CagdSrf2Polylines, M * CagdSrf2PolygonStrip, CagdSrf2Polygons M * * * KEYWORDS: M * CagdSrf2PolygonsN, evaluation, polygonal approximation M *****************************************************************************/ CagdPolygonStruct *CagdSrf2PolygonsN(CagdSrfStruct *Srf, int Nu, int Nv, CagdBType ComputeNormals, CagdBType FourPerFlat, CagdBType ComputeUV) { CagdPolygonStruct *Pls; if (Nu < 2 || Nv < 2) { CAGD_FATAL_ERROR(CAGD_ERR_WRONG_SIZE); return FALSE; } switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrf2PolygonsN(Srf, Nu, Nv, ComputeNormals, FourPerFlat, ComputeUV); case CAGD_SBSPLINE_TYPE: return BspSrf2PolygonsN(Srf, Nu, Nv, ComputeNormals, FourPerFlat, ComputeUV); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); Pls = BzrSrf2PolygonsN(Srf, Nu, Nv, ComputeNormals, FourPerFlat, ComputeUV); CagdSrfFree(Srf); return Pls; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Routine to convert a single Bspline surface to NumOfIsolines polylines M * in each parametric direction with SamplesPerCurve in each isoparametric M * curve. M * Polyline are always E3 of CagdPolylineStruct type. M * Iso parametric curves are sampled equally spaced in parametric space. M * NULL is returned in case of an error, otherwise list of M * CagdPolylineStruct. Attempt is made to extract isolines along C1 M * discontinuities first. M * * * PARAMETERS: M * Srf: Srf to extract isoparametric curves from. M * NumOfIsocurves: To extarct from Srf in each (U or V) direction. M * SamplesPerCurve: Fineness control on piecewise linear curve M * approximation. M * * * RETURN VALUE: M * CagdPolylineStruct *: List of polygons representing a piecewise linear M * approximation of the extracted isoparamteric M * curves or NULL is case of an error. M * * * SEE ALSO: M * BspCrv2Polyline, BzrSrf2Polylines, IritSurface2Polylines, M * IritTrimSrf2Polylines, SymbSrf2Polylines, TrimSrf2Polylines M * * * KEYWORDS: M * CagdSrf2Polylines, polylines, isoparametric curves M *****************************************************************************/ CagdPolylineStruct *CagdSrf2Polylines(CagdSrfStruct *Srf, int NumOfIsocurves[2], int SamplesPerCurve) { CagdPolylineStruct *Pls; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrf2Polylines(Srf, NumOfIsocurves, SamplesPerCurve); case CAGD_SBSPLINE_TYPE: return BspSrf2Polylines(Srf, NumOfIsocurves, SamplesPerCurve); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); Pls = BzrSrf2Polylines(Srf, NumOfIsocurves, SamplesPerCurve); CagdSrfFree(Srf); return Pls; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Routine to extract from a surface NumOfIsoline isocurve list M * in each param. direction. M * Iso parametric curves are sampled equally spaced in parametric space. M * NULL is returned in case of an error, otherwise list of CagdCrvStruct. M * * * PARAMETERS: M * Srf: To extract isoparametric curves from. M * NumOfIsocurves: In each (U or V) direction. M * * * RETURN VALUE: M * CagdCrvStruct *: List of extracted isoparametric curves. These curves M * inherit the order and continuity of the original Srf. M * NULL is returned in case of an error. M * * * SEE ALSO: M * BspSrf2PCurves, BzrSrf2Curves, SymbSrf2Curves M * * * KEYWORDS: M * CagdSrf2Curves, curves, isoparametric curves M *****************************************************************************/ CagdCrvStruct *CagdSrf2Curves(CagdSrfStruct *Srf, int NumOfIsocurves[2]) { CagdCrvStruct *Crvs; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrf2Curves(Srf, NumOfIsocurves); case CAGD_SBSPLINE_TYPE: return BspSrf2Curves(Srf, NumOfIsocurves); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); Crvs = BzrSrf2Curves(Srf, NumOfIsocurves); CagdSrfFree(Srf); return Crvs; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Evaluates a vector field surface to a unit size vector. If fails, moves a M * tad until success. Useful for normal field evaluations. M * * * PARAMETERS: M * Vec: Where resulting unit length vector is to be saved. M * VecFieldSrf: A surface representing a vector field. M * U, V: Parameter locations. M * * * RETURN VALUE: M * void M * * * KEYWORDS: M * CagdEvaluateSurfaceVecField, normal, vector field M *****************************************************************************/ void CagdEvaluateSurfaceVecField(CagdVType Vec, CagdSrfStruct *VecFieldSrf, CagdRType U, CagdRType V) { CagdRType *R = CagdSrfEval(VecFieldSrf, U, V); CagdCoerceToE3(Vec, &R, -1, VecFieldSrf -> PType); if (PT_SQR_LENGTH(Vec) < SQR(IRIT_UEPS)) { int i = 0; CagdRType UMin, UMax, VMin, VMax, UMid, VMid, Step = VEC_FIELD_START_STEP; CagdSrfDomain(VecFieldSrf, &UMin, &UMax, &VMin, &VMax); UMid = (UMin + UMax) * 0.5; VMid = (VMin + VMax) * 0.5; while (PT_SQR_LENGTH(Vec) < SQR(IRIT_UEPS) && i++ < VEC_FIELD_TRIES) { U += U < UMid ? Step : -Step; V += V < VMid ? Step : -Step; Step *= 2.0; R = CagdSrfEval(VecFieldSrf, U, V); CagdCoerceToE3(Vec, &R, -1, VecFieldSrf -> PType); } } PT_NORMALIZE(Vec); } /***************************************************************************** * DESCRIPTION: M * Given a curve, computes its derivative curve (Hodograph). M * * * PARAMETERS: M * Crv: To compute its Hodograph curve. M * * * RETURN VALUE: M * CagdCrvStruct *: Resulting hodograph. M * * * SEE ALSO: M * BzrCrvDerive, BspCrvDerive, BzrCrvDeriveRational, BspCrvDeriveRational M * CagdCrvDeriveScalar M * * * KEYWORDS: M * CagdCrvDerive, derivatives, Hodograph M *****************************************************************************/ CagdCrvStruct *CagdCrvDerive(CagdCrvStruct *Crv) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvDerive(Crv); case CAGD_CBSPLINE_TYPE: return BspCrvDerive(Crv); case CAGD_CPOWER_TYPE: return PwrCrvDerive(Crv); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve, computes the derivative of all is scalar components. For a M * Euclidean curve this is the same as CagdCrvDerive but for a rational curve M * the returned curves is not the vector field but simply the derivatives of M * all the curve's coefficients, including the weights. M * * * PARAMETERS: M * Crv: To compute derivatives of all its components. M * * * RETURN VALUE: M * CagdCrvStruct *: Resulting derivative. M * * * SEE ALSO: M * BzrCrvDerive, BspCrvDerive, BzrCrvDeriveRational, BspCrvDeriveRational M * CagdCrvDerive, BzrCrvDeriveScalar, BspCrvDeriveScalar M * * * KEYWORDS: M * CagdCrvDeriveScalar, derivatives, Hodograph M *****************************************************************************/ CagdCrvStruct *CagdCrvDeriveScalar(CagdCrvStruct *Crv) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvDeriveScalar(Crv); case CAGD_CBSPLINE_TYPE: return BspCrvDeriveScalar(Crv); case CAGD_CPOWER_TYPE: return PwrCrvDeriveScalar(Crv); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve, compute its integral curve. M * * * PARAMETERS: M * Crv: To compute its integral curve. M * * * RETURN VALUE: M * CagdCrvStruct *: Resulting integral curve. M * * * SEE ALSO: M * BzrCrvIntegrate, BspCrvIntegrate M * * * KEYWORDS: M * CagdCrvIntegrate, integrals M *****************************************************************************/ CagdCrvStruct *CagdCrvIntegrate(CagdCrvStruct *Crv) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvIntegrate(Crv); case CAGD_CBSPLINE_TYPE: return BspCrvIntegrate(Crv); case CAGD_CPOWER_TYPE: return PwrCrvIntegrate(Crv); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve, compute its moebius transformation. M * * * PARAMETERS: M * Crv: To compute its moebius transformation. M * c: The scaling coefficient - c^n is the ratio between the first M * and last weight of the curve. M * If c == 0, the first and last weights are made equal. M * * * RETURN VALUE: M * CagdCrvStruct *: Resulting curve after the moebius transformation. M * * * SEE ALSO: M * BzrCrvMoebiusTransform, BspCrvMoebiusTransform M * * * KEYWORDS: M * CagdCrvMoebiusTransform, moebius transformation M *****************************************************************************/ CagdCrvStruct *CagdCrvMoebiusTransform(CagdCrvStruct *Crv, CagdRType c) { CagdCrvStruct *TCrv; switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvMoebiusTransform(Crv, c); case CAGD_CBSPLINE_TYPE: return BspCrvMoebiusTransform(Crv, c); case CAGD_CPOWER_TYPE: Crv = CnvrtPower2BezierCrv(Crv); TCrv = BzrCrvMoebiusTransform(Crv, c); CagdCrvFree(Crv); Crv = CnvrtBezier2PowerCrv(TCrv); CagdCrvFree(TCrv); return Crv; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface, computes its partial derivative in the prescibed M * direction Dir. M * * * PARAMETERS: M * Srf: To compute its derivative surface in direction Dir. M * Dir: Direction of differentiation. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: Resulting partial derivative surface. M * * * SEE ALSO: M * BzrSrfDerive, BspSrfDerive, BzrSrfDeriveRational, BspSrfDeriveRational M * CagdSrfDeriveScalar M * * * KEYWORDS: M * CagdSrfDerive, derivatives, partial derivatives M *****************************************************************************/ CagdSrfStruct *CagdSrfDerive(CagdSrfStruct *Srf, CagdSrfDirType Dir) { CagdSrfStruct *TSrf; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfDerive(Srf, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfDerive(Srf, Dir); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); TSrf = BzrSrfDerive(Srf, Dir); CagdSrfFree(Srf); Srf = CnvrtBezier2PowerSrf(TSrf); CagdSrfFree(TSrf); return Srf; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface, computes its partial derivative in the prescibed M * direction Dir of all is scalar components. M * For a Euclidean surface this is the same as CagdSrfDerive but for a M * rational surface the returned surfaces is not the vector field but simply M * the derivatives of all the surface's coefficients, including the weights. M * * * PARAMETERS: M * Srf: To compute derivatives of all its components. M * Dir: Direction of differentiation. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: Resulting derivative. M * * * SEE ALSO: M * BzrSrfDerive, BspSrfDerive, BzrSrfDeriveRational, BspSrfDeriveRational M * CagdSrfDerive, BzrSrfDeriveScalar, BspSrfDeriveScalar M * * * KEYWORDS: M * CagdSrfDeriveScalar, derivatives, Hodograph M *****************************************************************************/ CagdSrfStruct *CagdSrfDeriveScalar(CagdSrfStruct *Srf, CagdSrfDirType Dir) { CagdSrfStruct *TSrf; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfDeriveScalar(Srf, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfDeriveScalar(Srf, Dir); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); TSrf = BzrSrfDeriveScalar(Srf, Dir); CagdSrfFree(Srf); Srf = CnvrtBezier2PowerSrf(TSrf); CagdSrfFree(TSrf); return Srf; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface, compute its integral surface. M * * * PARAMETERS: M * Srf: To compute its integral surface. M * Dir: Direction of integration. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: Resulting integral surface. M * * * SEE ALSO: M * BzrSrfIntegrate, BspSrfIntegrate M * * * KEYWORDS: M * CagdSrfIntegrate, integrals M *****************************************************************************/ CagdSrfStruct *CagdSrfIntegrate(CagdSrfStruct *Srf, CagdSrfDirType Dir) { CagdSrfStruct *TSrf; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfIntegrate(Srf, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfIntegrate(Srf, Dir); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); TSrf = BzrSrfIntegrate(Srf, Dir); CagdSrfFree(Srf); Srf = CnvrtBezier2PowerSrf(TSrf); CagdSrfFree(TSrf); return Srf; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface, compute its moebius transformation. M * * * PARAMETERS: M * Srf: Surface to apply the Moebius transformation to. M * c: The scaling coefficient - c^n is the ratio between the first M * and last weight of the surface, along each row or column. M * If c == 0, the first and last weights are made equal, in the M * first row/column. M * Dir: Direction to apply the Moebius transformation, row or col. M * If Dir == CAGD_BOTH_DIR, the transformation is applied to M * both the row and column directions, in this order. M * * * RETURN VALUE: M * CagdSrfStruct *: Resulting surface after the moebius transformation. M * * * SEE ALSO: M * BzrSrfMoebiusTransform, BspSrfMoebiusTransform M * * * KEYWORDS: M * CagdSrfMoebiusTransform, moebius transformation M *****************************************************************************/ CagdSrfStruct *CagdSrfMoebiusTransform(CagdSrfStruct *Srf, CagdRType c, CagdSrfDirType Dir) { CagdSrfStruct *TSrf; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfMoebiusTransform(Srf, c, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfMoebiusTransform(Srf, c, Dir); case CAGD_SPOWER_TYPE: Srf = CnvrtPower2BezierSrf(Srf); TSrf = BzrSrfMoebiusTransform(Srf, c, Dir); CagdSrfFree(Srf); Srf = CnvrtBezier2PowerSrf(TSrf); CagdSrfFree(TSrf); return Srf; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve - subdivides it into two curves at the given parameter M * value t. M * Returns pointer to first curve in a list of two subdivided curves. M * * * PARAMETERS: M * Crv: To subdivide at the prescibed parameter value t. M * t: The parameter to subdivide the curve Crv at. M * * * RETURN VALUE: M * CagdCrvStruct *: A list of the two curves resulting from the process M * of subdivision. M * * * SEE ALSO: M * CagdCrvSubdivAtParams M * * * KEYWORDS: M * CagdCrvSubdivAtParam, subdivision M *****************************************************************************/ CagdCrvStruct *CagdCrvSubdivAtParam(CagdCrvStruct *Crv, CagdRType t) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvSubdivAtParam(Crv, t); case CAGD_CBSPLINE_TYPE: return BspCrvSubdivAtParam(Crv, t); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve - subdivides it into curves at all the given parameter M * values Pts. Pts is assumed to hold the parameters in order in the first M * point coordinate. M * Returns pointer to first curve in a list of subdivided curves. M * M * * * PARAMETERS: M * Crv: Curve to split at all parameter values as prescribed by Pts. M * Bezier curves are promoted to Bspline curves in this function. M * Pts: Ordered list of parameter values (first coordinate of point) to M * split curve Crv at. M * Eps: parameter closer than Eps and/or closer to boundary than Eps are M 8 ignored. M * * * RETURN VALUE: M * CagdCrvStruct *: List of splitted curves, in order. M * * * SEE ALSO: M * CagdCrvSubdivAtParam M * * * KEYWORDS: M * CagdCrvSubdivAtParams M *****************************************************************************/ CagdCrvStruct *CagdCrvSubdivAtParams(CagdCrvStruct *Crv, CagdPtStruct *Pts, CagdRType Eps) { CagdRType t, TMax; CagdCrvStruct *Crvs, *OrigCrv = Crv, *RetList = NULL; CagdCrvDomain(Crv, &t, &TMax); switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: Crv = CnvrtBezier2BsplineCrv(Crv); break; case CAGD_CBSPLINE_TYPE: Crv = CagdCrvCopy(Crv); break; case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } while (Pts != NULL) { if (Pts -> Pt[0] >= TMax || (Eps != 0.0 && APX_EQ_EPS(Pts -> Pt[0], TMax, Eps))) break; /* Verify monotonicity of parameter values. */ if (Pts -> Pt[0] > t && (Eps == 0.0 || !APX_EQ_EPS(Pts -> Pt[0], t, Eps))) { Crvs = CagdCrvSubdivAtParam(Crv, Pts -> Pt[0]); if (Crv != OrigCrv) CagdCrvFree(Crv); Crv = Crvs -> Pnext; LIST_PUSH(Crvs, RetList); t = Pts -> Pt[0]; } Pts = Pts -> Pnext; } LIST_PUSH(Crv, RetList); return CagdListReverse(RetList); } /***************************************************************************** * DESCRIPTION: M * Given a curve - extracts a sub-region within the domain specified by t1 M * and t2. M * * * PARAMETERS: M * Crv: To extract a sub-region from. M * t1, t2: Parametric domain boundaries of sub-region. M * * * RETURN VALUE: M * CagdCrvStruct *: Sub-region extracted from Crv from t1 to t2. M * * * KEYWORDS: M * CagdCrvRegionFromCrv, regions, subdivision M *****************************************************************************/ CagdCrvStruct *CagdCrvRegionFromCrv(CagdCrvStruct *Crv, CagdRType t1, CagdRType t2) { CagdRType TMin, TMax; CagdCrvStruct *Crvs; CagdBType BezCrv = FALSE, OpenEnd = TRUE, NewCrv = FALSE; switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: BezCrv = TRUE; break; case CAGD_CBSPLINE_TYPE: OpenEnd = BspCrvHasOpenEC(Crv); break; case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } CagdCrvDomain(Crv, &TMin, &TMax); if (!CAGD_IS_BEZIER_CRV(Crv)) { /* No limits on Bezier curves! */ CAGD_DOMAIN_T_VERIFY(t1, TMin, TMax); CAGD_DOMAIN_T_VERIFY(t2, TMin, TMax); } if (APX_EQ(t1, TMin) && APX_EQ(t2, TMax) && OpenEnd) { return CagdCrvCopy(Crv); } if (t1 > t2) SWAP(CagdRType, t1, t2); if (!OpenEnd && CAGD_IS_PERIODIC_CRV(Crv)) { Crv = CnvrtPeriodic2FloatCrv(Crv); NewCrv = TRUE; } if (!APX_EQ(t1, TMin) || !OpenEnd) { Crvs = CagdCrvSubdivAtParam(Crv, t1); if (NewCrv) CagdCrvFree(Crv); Crv = Crvs -> Pnext; Crvs -> Pnext = NULL; CagdCrvFree(Crvs); /* Free the first region. */ NewCrv = TRUE; } if (APX_EQ(t2, TMax) && OpenEnd) return NewCrv ? Crv : CagdCrvCopy(Crv); else { if (BezCrv) t2 = (t2 - t1) / (TMax - t1); Crvs = CagdCrvSubdivAtParam(Crv, t2); if (NewCrv) CagdCrvFree(Crv); CagdCrvFree(Crvs -> Pnext); /* Free the second region. */ Crvs -> Pnext = NULL; return Crvs; /* Returns the first region. */ } } /***************************************************************************** * DESCRIPTION: M * Given a curve - refines it at the given n knots as defined by vector t. M * If Replace is TRUE, the values in t replaces current knot vector. M * Returns pointer to refined surface (Note a Bezier curve will be converted M * into a Bspline curve). M * * * PARAMETERS: M * Crv: To refine. M * Replace: If TRUE, t holds knots in exactly the same length as the M * length of the knot vector of Crv and t simply replaces the M * knot vector. M * t: Vector of knots with length of n. M * n: Length of vector t. M * * * RETURN VALUE: M * CagdCrvStruct *: A refined curve of Crv after insertion of all the M * knots as specified by vector t of length n. M * * * KEYWORDS: M * CagdCrvRefineAtParams, refinement, subdivision M *****************************************************************************/ CagdCrvStruct *CagdCrvRefineAtParams(CagdCrvStruct *Crv, CagdBType Replace, CagdRType *t, int n) { CagdCrvStruct *BspCrv, *TCrv; switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: BspCrv = CnvrtBezier2BsplineCrv(Crv); TCrv = BspCrvKnotInsertNDiff(BspCrv, Replace, t, n); CagdCrvFree(BspCrv); return TCrv; case CAGD_CBSPLINE_TYPE: return BspCrvKnotInsertNDiff(Crv, Replace, t, n); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns a new curve that is the reversed curve of Crv by reversing the M * control polygon and the knot vector of Crv is a Bspline curve. M * See also BspKnotReverse. M * * * PARAMETERS: M * Crv: To be reversed. M * * * RETURN VALUE: M * CagdCrvStruct *: Reversed curve of Crv. M * * * KEYWORDS: M * CagdCrvReverse, reverse M *****************************************************************************/ CagdCrvStruct *CagdCrvReverse(CagdCrvStruct *Crv) { CagdBType IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv); int i, Len, Col, Length = Crv -> Length, MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType); CagdCrvStruct *ReversedCrv = CagdCrvCopy(Crv); CagdRType *KV, **Points = ReversedCrv -> Points; switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: case CAGD_CBSPLINE_TYPE: break; case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } /* Reverse the Ctl Polygon: */ Len = (Length >> 1); for (Col = 0; Col < Len; Col++) for (i = IsNotRational; i <= MaxCoord; i++) SWAP(CagdRType, Points[i][Col], Points[i][Length - Col - 1]); /* Reverse the knot vector if it exists: */ if (Crv -> GType == CAGD_CBSPLINE_TYPE && Crv -> KnotVector != NULL) { KV = BspKnotReverse(Crv -> KnotVector, Crv -> Order + CAGD_CRV_PT_LST_LEN(Crv)); IritFree(ReversedCrv -> KnotVector); ReversedCrv -> KnotVector = KV; } return ReversedCrv; } /***************************************************************************** * DESCRIPTION: M * Returns a new curve representing the same curve as Crv but with its degree M * raised by one. M * * * PARAMETERS: M * Crv: To raise its degree. M * * * RETURN VALUE: M * CagdCrvStruct *: A curve with same geometry as Crv but with one degree M * higher. M * * * KEYWORDS: M * CagdCrvDegreeRaise, degree raising M *****************************************************************************/ CagdCrvStruct *CagdCrvDegreeRaise(CagdCrvStruct *Crv) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvDegreeRaise(Crv); case CAGD_CBSPLINE_TYPE: return CagdCrvBlossomDegreeRaise(Crv); case CAGD_CPOWER_TYPE: return PwrCrvDegreeRaise(Crv); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns a new curve representing the same curve as Crv but with its degree M * raised by one. M * * * PARAMETERS: M * Crv: To raise its degree. M * * * RETURN VALUE: M * CagdCrvStruct *: A curve with same geometry as Crv but with one degree M * higher. M * * * KEYWORDS: M * CagdCrvDegreeReduce, degree raising M *****************************************************************************/ CagdCrvStruct *CagdCrvDegreeReduce(CagdCrvStruct *Crv) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvDegreeReduce(Crv); case CAGD_CBSPLINE_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_BSPLINE_NO_SUPPORT); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns a new curve representing the same curve as Crv but with its degree M * raised to NewOrder M * * * PARAMETERS: M * Crv: To raise its degree. M * NewOrder: Expected new order of the raised curve. M * * * RETURN VALUE: M * CagdCrvStruct *: A curve with same geometry as Crv but with order that M * is equal to NewOrder. M * * * KEYWORDS: M * CagdCrvDegreeRaiseN, degree raising M *****************************************************************************/ CagdCrvStruct *CagdCrvDegreeRaiseN(CagdCrvStruct *Crv, int NewOrder) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvDegreeRaiseN(Crv, NewOrder); case CAGD_CBSPLINE_TYPE: return BspCrvDegreeRaiseN(Crv, NewOrder); case CAGD_CPOWER_TYPE: return PwrCrvDegreeRaiseN(Crv, NewOrder); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns a new surface representing the same surface as Srf but with its M * degree raised by one. M * * * PARAMETERS: M * Srf: To raise its degree. M * Dir: Direction of degree raising. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: A surface with same geometry as Srf but with one M * degree higher. M * * * SEE ALSO: M * BzrSrfDegreeRaise, BspSrfDegreeRaise, TrimSrfDegreeRaise M * * * KEYWORDS: M * CagdSrfDegreeRaise, degree raising M *****************************************************************************/ CagdSrfStruct *CagdSrfDegreeRaise(CagdSrfStruct *Srf, CagdSrfDirType Dir) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfDegreeRaise(Srf, Dir); case CAGD_SBSPLINE_TYPE: return CagdSrfBlossomDegreeRaise(Srf, Dir); case CAGD_SPOWER_TYPE: return PwrSrfDegreeRaise(Srf, Dir); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Returns a new surface, identical to the original but with higher degrees, M * as prescribed by NewUOrder, NewVOrder. M * * * PARAMETERS: M * Srf: To raise its degree. M * NewUOrder: New U order of Srf. M * NewVOrder: New V order of Srf. M * * * RETURN VALUE: M * CagdSrfStruct *: A surface with higher degrees as prescribed by M * NewUOrder/NewVOrder. M * * * SEE ALSO: M * CagdSrfDegreeRaise, BzrSrfDegreeRaise, TrimSrfDegreeRaise M * BspSrfDegreeRaise, BzrSrfDegreeRaiseN, BspSrfDegreeRaiseN M ** * * KEYWORDS: M * CagdSrfDegreeRaiseN, degree raising M *****************************************************************************/ CagdSrfStruct *CagdSrfDegreeRaiseN(CagdSrfStruct *Srf, int NewUOrder, int NewVOrder) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfDegreeRaiseN(Srf, NewUOrder, NewVOrder); case CAGD_SBSPLINE_TYPE: return BspSrfDegreeRaiseN(Srf, NewUOrder, NewVOrder); case CAGD_SPOWER_TYPE: return PwrSrfDegreeRaiseN(Srf, NewUOrder, NewVOrder); default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Extracts an isoparametric curve from the surface Srf in direction Dir at M * the parameter value of t. M * * * PARAMETERS: M * Srf: To extract an isoparametric curve from. M * t: Parameter value of extracted isoparametric curve. M * Dir: Direction of extracted isoparametric curve. Either U or V. M * * * RETURN VALUE: M * CagdCrvStruct *: An isoparametric curve of Srf. This curve inherit the M * order and continuity of surface Srf in direction Dir. M * * * SEE ALSO: M * BzrSrfCrvFromSrf, BspSrfCrvFromSrf, CagdCrvFromMesh, BzrSrfCrvFromMesh, M * BspSrfCrvFromMesh, TrngCrvFromTriSrf M * * * KEYWORDS: M * CagdCrvFromSrf, isoparametric curves, curve from surface M *****************************************************************************/ CagdCrvStruct *CagdCrvFromSrf(CagdSrfStruct *Srf, CagdRType t, CagdSrfDirType Dir) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfCrvFromSrf(Srf, t, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfCrvFromSrf(Srf, t, Dir); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Extracts the four boundary curves of the given surface. M * * * PARAMETERS: M * Srf: To extract the boundary from. M * * * RETURN VALUE: M * CagdCrvStruct **: A pointer to a static vector of four curve pointers, M * representing the four boundaries of surface Srf M * in order of UMin, UMax, VMin, VMax. M * * * KEYWORDS: M * CagdBndryCrvsFromSrf, isoparametric curves, curve from surface M *****************************************************************************/ CagdCrvStruct **CagdBndryCrvsFromSrf(CagdSrfStruct *Srf) { static CagdCrvStruct *Crvs[4]; CagdRType UMin, UMax, VMin, VMax; CagdSrfDomain(Srf, &UMin, &UMax, &VMin, &VMax); Crvs[0] = CagdCrvFromSrf(Srf, UMin, CAGD_CONST_U_DIR); Crvs[1] = CagdCrvFromSrf(Srf, UMax, CAGD_CONST_U_DIR); Crvs[2] = CagdCrvFromSrf(Srf, VMin, CAGD_CONST_V_DIR); Crvs[3] = CagdCrvFromSrf(Srf, VMax, CAGD_CONST_V_DIR); return Crvs; } /***************************************************************************** * DESCRIPTION: M * Extracts a curve from the mesh of surface Srf in direction Dir at index M * Index. M * * * PARAMETERS: M * Srf: To extract a curve from. M * Index: Index along the mesh of Srf to extract the curve from. M * Dir: Direction of extracted curve. Either U or V. M * * * RETURN VALUE: M * CagdCrvStruct *: A curve from Srf. This curve inherit the order and M * continuity of surface Srf in direction Dir. However, M * thiscurve is not on surface Srf, in general. M * * * SEE ALSO: M * CagdCrvFromSrf, BzrSrfCrvFromSrf, BspSrfCrvFromSrf, BzrSrfCrvFromMesh, M * BspSrfCrvFromMesh, CagdCrvToMesh M * * * KEYWORDS: M * CagdCrvFromMesh, isoparametric curves, curve from mesh M *****************************************************************************/ CagdCrvStruct *CagdCrvFromMesh(CagdSrfStruct *Srf, int Index, CagdSrfDirType Dir) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfCrvFromMesh(Srf, Index, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfCrvFromMesh(Srf, Index, Dir); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Substitutes a row/column of surface Srf from the given curve Crv at M * surface direction Dir and mesh index Index. Curve must have the same M * PtType/Length as the surface in the selected direction. M * * * PARAMETERS: M * Crv: To substitute into the surface Srf. M * Index: Of mesh where the curve Crv should be substituted in. M * Dir: Either U or V. M * Srf: That a row or a column of should be replaced by Crv. M * * * RETURN VALUE: M * void M * * * SEE ALSO: M * CagdCrvFromSrf, CagdCrvFromMesh M * * * KEYWORDS: M * CagdCrvToMesh, curve from mesh M *****************************************************************************/ void CagdCrvToMesh(CagdCrvStruct *Crv, int Index, CagdSrfDirType Dir, CagdSrfStruct *Srf) { CagdBType IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf); int i, j, Length = Crv -> Length, ULength = Srf -> ULength, VLength = Srf -> VLength, MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType); CagdRType *CrvP, *SrfP; if (Crv -> PType != Srf -> PType || Length != (Dir == CAGD_CONST_U_DIR ? VLength : ULength)) CAGD_FATAL_ERROR(CAGD_ERR_PT_OR_LEN_MISMATCH); switch (Dir) { case CAGD_CONST_U_DIR: if (Index + 1 > ULength) CAGD_FATAL_ERROR(CAGD_ERR_INDEX_NOT_IN_MESH); for (i = IsNotRational; i <= MaxCoord; i++) { CrvP = Crv -> Points[i]; SrfP = Srf -> Points[i] + Index * CAGD_NEXT_U(Srf); for (j = 0; j < Length; j++) { *SrfP = *CrvP++; SrfP += CAGD_NEXT_V(Srf); } } break; case CAGD_CONST_V_DIR: if (Index + 1 > VLength) CAGD_FATAL_ERROR(CAGD_ERR_INDEX_NOT_IN_MESH); for (i = IsNotRational; i <= MaxCoord; i++) { CrvP = Crv -> Points[i]; SrfP = Srf -> Points[i] + Index * CAGD_NEXT_V(Srf); for (j = 0; j < Length; j++) { *SrfP = *CrvP++; SrfP += CAGD_NEXT_U(Srf); } } break; default: CAGD_FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV); break; } } /***************************************************************************** * DESCRIPTION: M * Given a surface - subdivides it into two sub-surfaces at given parametric M * value t in the given direction Dir. M * Returns pointer to first surface in a list of two subdivided surfaces. M * * * PARAMETERS: M * Srf: To subdivide at the prescibed parameter value t. M * t: The parameter to subdivide the curve Crv at. M * Dir: Direction of subdivision. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: A list of the two surfaces resulting from the process M * of subdivision. M * * * KEYWORDS: M * CagdSrfSubdivAtParam, subdivision M *****************************************************************************/ CagdSrfStruct *CagdSrfSubdivAtParam(CagdSrfStruct *Srf, CagdRType t, CagdSrfDirType Dir) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfSubdivAtParam(Srf, t, Dir); case CAGD_SBSPLINE_TYPE: return BspSrfSubdivAtParam(Srf, t, Dir); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface - extracts a sub-region within the domain specified by t1 M * and t2, in the direction Dir. M * * * PARAMETERS: M * Srf: To extract a sub-region from. M * t1, t2: Parametric domain boundaries of sub-region. M * Dir: Direction of region extraction. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: Sub-region extracted from Srf from t1 to t2. M * * * KEYWORDS: M * CagdSrfRegionFromSrf, regions, subdivision M *****************************************************************************/ CagdSrfStruct *CagdSrfRegionFromSrf(CagdSrfStruct *Srf, CagdRType t1, CagdRType t2, CagdSrfDirType Dir) { CagdRType TMin, TMax, R1, R2; CagdSrfStruct *Srfs; CagdBType OpenEnd = FALSE, NewSrf = FALSE; if (t1 > t2) SWAP(CagdRType, t1, t2); if (Dir == CAGD_CONST_U_DIR) CagdSrfDomain(Srf, &TMin, &TMax, &R1, &R2); else CagdSrfDomain(Srf, &R1, &R2, &TMin, &TMax); if (!CAGD_IS_BEZIER_SRF(Srf)) { /* No limits on Bezier surfaces! */ CAGD_DOMAIN_T_VERIFY(t1, TMin, TMax); CAGD_DOMAIN_T_VERIFY(t2, TMin, TMax); } switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: /* Update t2 to be between t1 and TMax as it will come back */ /* after the first subdivision to be between zero and one. */ t2 = 1.0 - (1.0 - t2) / (1.0 - t1); break; case CAGD_SBSPLINE_TYPE: OpenEnd = BspSrfHasOpenECDir(Srf, Dir); break; case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } if (!OpenEnd && (Dir == CAGD_CONST_U_DIR ? CAGD_IS_UPERIODIC_SRF(Srf) : CAGD_IS_VPERIODIC_SRF(Srf))) { Srf = CnvrtPeriodic2FloatSrf(Srf); NewSrf = TRUE; } if (!APX_EQ(t1, TMin) || !OpenEnd) { Srfs = CagdSrfSubdivAtParam(Srf, t1, Dir); if (NewSrf) CagdSrfFree(Srf); Srf = Srfs -> Pnext; Srfs -> Pnext = NULL; CagdSrfFree(Srfs); /* Free the first region. */ NewSrf = TRUE; } if (APX_EQ(t2, TMax) && OpenEnd) return NewSrf ? Srf : CagdSrfCopy(Srf); else { Srfs = CagdSrfSubdivAtParam(Srf, t2, Dir); if (NewSrf) CagdSrfFree(Srf); CagdSrfFree(Srfs -> Pnext); /* Free the second region. */ Srfs -> Pnext = NULL; return Srfs; /* Returns the first region. */ } } /***************************************************************************** * DESCRIPTION: M * Given a surface - refines it at the given n knots as defined by vector t. M * If Replace is TRUE, the values in t replaces current knot vector. M * Returns pointer to refined surface (Note a Bezier surface will be M * converted into a Bspline surface). M * * * PARAMETERS: M * Srf: To refine. M * Dir: Direction of refinement. Either U or V. M * Replace: If TRUE, t holds knots in exactly the same length as the M * length of the knot vector of Srf and t simply replaces the M * knot vector. M * t: Vector of knots with length of n. M * n: Length of vector t. M * * * RETURN VALUE: M * CagdSrfStruct *: A refined curve of Srf after insertion of all the M * knots as specified by vector t of length n. M * * * KEYWORDS: M * CagdSrfRefineAtParams, refinement, subdivision M *****************************************************************************/ CagdSrfStruct *CagdSrfRefineAtParams(CagdSrfStruct *Srf, CagdSrfDirType Dir, CagdBType Replace, CagdRType *t, int n) { CagdSrfStruct *BspSrf, *TSrf; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: BspSrf = CnvrtBezier2BsplineSrf(Srf); TSrf = BspSrfKnotInsertNDiff(BspSrf, Dir, Replace, t, n); CagdSrfFree(BspSrf); return TSrf; case CAGD_SBSPLINE_TYPE: return BspSrfKnotInsertNDiff(Srf, Dir, Replace, t, n); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve Crv and a parameter value t, returns the (unit) tangent M * direction of Crv at t. M * The unnormalized normal does not equal dC/dt in its magnitude, only in M * its direction. M * * * PARAMETERS: M * Crv: To compute (unit) tangent vector for. M * t: Location where to evaluate the tangent of Crv. M * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, returned is an unnormalized vector in the right M * direction of the tangent. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit tanegnt M * information. M * * * KEYWORDS: M * CagdCrvTangent, tangent M *****************************************************************************/ CagdVecStruct *CagdCrvTangent(CagdCrvStruct *Crv, CagdRType t, CagdBType Normalize) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvTangent(Crv, t, Normalize); case CAGD_CBSPLINE_TYPE: return BspCrvTangent(Crv, t, Normalize); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve Crv and a parameter value t, returns the (unit) binormal M * direction of Crv at t. M * * * PARAMETERS: M * Crv: To compute (unit) binormal vector for. M * t: Location where to evaluate the binormal of Crv. M * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, length is a function of given parametrization. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit binormal M * information. M * * * KEYWORDS: M * CagdCrvBiNormal, binormal M *****************************************************************************/ CagdVecStruct *CagdCrvBiNormal(CagdCrvStruct *Crv, CagdRType t, CagdBType Normalize) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvBiNormal(Crv, t, Normalize); case CAGD_CBSPLINE_TYPE: return BspCrvBiNormal(Crv, t, Normalize); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve Crv and a parameter value t, returns the (unit) normal M * direction of Crv at t. M * * * PARAMETERS: M * Crv: To compute (unit) normal vector for. M * t: Location where to evaluate the normal of Crv. M * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, length is a function of given parametrization. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit normal M * information. M * * * KEYWORDS: M * CagdCrvNormal, normal M *****************************************************************************/ CagdVecStruct *CagdCrvNormal(CagdCrvStruct *Crv, CagdRType t, CagdBType Normalize) { switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: return BzrCrvNormal(Crv, t, Normalize); case CAGD_CBSPLINE_TYPE: return BspCrvNormal(Crv, t, Normalize); case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a curve Crv and a parameter value t, returns the (unit) normal M * direction of Crv at t, that is consistent over inflection points. M * That is, this normal is not flipped over inflection points and is always M * 90 rotation from the tangent vector. M * Needless to say, this function is for two dimensional palanr curves. M * * * PARAMETERS: M * Crv: To compute (unit) normal vector for. M * t: Location where to evaluate the normal of Crv. M * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, length is a function of given parametrization. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit normal M * information. M * * * KEYWORDS: M * CagdCrvNormalXY, normal M *****************************************************************************/ CagdVecStruct *CagdCrvNormalXY(CagdCrvStruct *Crv, CagdRType t, CagdBType Normalize) { CagdRType R; CagdVecStruct *Vec = NULL; switch (Crv -> GType) { case CAGD_CBEZIER_TYPE: Vec = BzrCrvTangent(Crv, t, FALSE); break; case CAGD_CBSPLINE_TYPE: Vec = BspCrvTangent(Crv, t, FALSE); break; case CAGD_CPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_CRV); return NULL; } /* Rotate 90 degrees: (x, y) -> (y, -x). */ R = Vec -> Vec[0]; Vec -> Vec[0] = Vec -> Vec[1]; Vec -> Vec[1] = -R; if (Normalize) PT_NORMALIZE(Vec -> Vec); return Vec; } /***************************************************************************** * DESCRIPTION: M * Given a surface Srf and a parameter values u, v, returns the (unit) M * tangent vector of Srf in direction Dir. M * * * PARAMETERS: M * Srf: To compute (unit) tangent vector for. M * u, v: Location where to evaluate the tangent of Srf. M * Dir: Direction of tangent. Either U or V. * * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, length is a function of given parametrization. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit tangent M * information. M * * * SEE ALSO: M * BzrSrfTangent, BspSrfTangent M * * * KEYWORDS: M * CagdSrfTangent, tangent M *****************************************************************************/ CagdVecStruct *CagdSrfTangent(CagdSrfStruct *Srf, CagdRType u, CagdRType v, CagdSrfDirType Dir, CagdBType Normalize) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfTangent(Srf, u, v, Dir, Normalize); case CAGD_SBSPLINE_TYPE: return BspSrfTangent(Srf, u, v, Dir, Normalize); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Given a surface Srf and a parameter values u, v, returns the (unit) normal M * vector of Srf. M * * * PARAMETERS: M * Srf: To compute (unit) normal vector for. M * u, v: Location where to evaluate the normal of Srf. M * Normalize: If TRUE, attempt is made to normalize the returned vector. M * If FALSE, length is a function of given parametrization. M * * * RETURN VALUE: M * CagdVecStruct *: A pointer to a static vector holding the unit normal M * information. M * * * SEE ALSO: M * BzrSrfNormal, BspSrfNormal, SymbSrfNormalSrf, TrngTriSrfNrml, M * BspSrfMeshNormals M * * * KEYWORDS: M * CagdSrfNormal, normal M *****************************************************************************/ CagdVecStruct *CagdSrfNormal(CagdSrfStruct *Srf, CagdRType u, CagdRType v, CagdBType Normalize) { switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: return BzrSrfNormal(Srf, u, v, Normalize); case CAGD_SBSPLINE_TYPE: return BspSrfNormal(Srf, u, v, Normalize); case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } } /***************************************************************************** * DESCRIPTION: M * Computes a new parametric direction OrthoUVDir that is orthogonal, in M * Eculidean space, to given parametric UVDir at parametric position UV of M * Srf. M * Clearly both Srf(OrthoUVDir) and Srf(UVDir) are in the tangent space. M * * * PARAMETERS: M * Srf: Surface to compute orthogonal direction to, in its tangent plane. M * UV: Location on Srf where to compute the orthogonal direction. M * UVDir: Direction to compute its orthogonal direction in Euclidean space. M * * * RETURN VALUE: M * CagdUVType *: UV direction that is orthogonal to UVDir, in Euclidean M * space, allocated statically. NULL, if error. M * * * KEYWORDS: M * CagdSrfUVDirOrthoE3 M *****************************************************************************/ CagdUVType *CagdSrfUVDirOrthoE3(CagdSrfStruct *Srf, CagdUVType *UV, CagdUVType *UVDir) { STATIC_DATA CagdUVType OrthoUVDir; int i; CagdRType Desc; CagdVecStruct *V; CagdVType dSdU, dSdV, Nrml, E3UVDir, E3OrthoUVDir; /* Compute the two tangent vectors and normal of Srf at UV. */ V = CagdSrfTangent(Srf, (*UV)[0], (*UV)[1], CAGD_CONST_V_DIR, FALSE); VEC_COPY(dSdU, V -> Vec); V = CagdSrfTangent(Srf, (*UV)[0], (*UV)[1], CAGD_CONST_U_DIR, FALSE); VEC_COPY(dSdV, V -> Vec); CROSS_PROD(Nrml, dSdU, dSdV); /* Find UVDir and its orthogonal in Srf tangent/Euclidean space. */ for (i = 0; i < 3; i++) E3UVDir[i] = (*UVDir)[0] * dSdU[i] + (*UVDir)[1] * dSdV[i]; CROSS_PROD(E3OrthoUVDir, Nrml, E3UVDir); if (PT_APX_EQ_ZERO_EPS(E3OrthoUVDir, IRIT_UEPS)) return NULL; /* We now need to solve for the parameteric direction that would yield */ /* E3OrthUVDir in the Euclidean space, at location UV on the surface. */ /* We have three dependent equations in X, Y, and Z and just two */ /* unknowns, U and V. Select the two equations that are the largest. */ if (FABS(Nrml[2]) > FABS(Nrml[1]) && FABS(Nrml[2]) > FABS(Nrml[0])) { /* Use the XY equations. */ if (FABS(Desc = dSdU[0] * dSdV[1] - dSdU[1] * dSdV[0]) < IRIT_UEPS) return NULL; OrthoUVDir[0] = (E3OrthoUVDir[0] * dSdV[1] - E3OrthoUVDir[1] * dSdV[0]) / Desc; OrthoUVDir[1] = (dSdU[0] * E3OrthoUVDir[1] - dSdU[1] * E3OrthoUVDir[0]) / Desc; } else if (FABS(Nrml[1]) > FABS(Nrml[2]) && FABS(Nrml[1]) > FABS(Nrml[0])) { /* Use the XZ equations. */ if (FABS(Desc = dSdU[0] * dSdV[2] - dSdU[2] * dSdV[0]) < IRIT_UEPS) return NULL; OrthoUVDir[0] = (E3OrthoUVDir[0] * dSdV[2] - E3OrthoUVDir[2] * dSdV[0]) / Desc; OrthoUVDir[1] = (dSdU[0] * E3OrthoUVDir[2] - dSdU[2] * E3OrthoUVDir[0]) / Desc; } else { /* Use the YZ equations. */ if (FABS(Desc = dSdU[1] * dSdV[2] - dSdU[2] * dSdV[1]) < IRIT_UEPS) return NULL; OrthoUVDir[0] = (E3OrthoUVDir[1] * dSdV[2] - E3OrthoUVDir[2] * dSdV[1]) / Desc; OrthoUVDir[1] = (dSdU[1] * E3OrthoUVDir[2] - dSdU[2] * E3OrthoUVDir[1]) / Desc; } #ifdef DEBUG /* Check we got it right... */ for (i = 0; i < 3; i++) E3OrthoUVDir[i] = OrthoUVDir[0] * dSdU[i] + OrthoUVDir[1] * dSdV[i]; VEC_NORMALIZE(E3OrthoUVDir); VEC_NORMALIZE(E3UVDir); if (FABS(DOT_PROD(E3OrthoUVDir, E3UVDir)) > IRIT_EPS) fprintf(stderr, "CagdSrfFindUVDirOrthoE3 orthogonality failed\n"); #endif /* DEBUG */ VEC2D_NORMALIZE(OrthoUVDir); return &OrthoUVDir; } /***************************************************************************** * DESCRIPTION: M * Returns a new surface that is the reversed surface of Srf by reversing the M * control mesh and the knot vector (if Bspline surface) of Srf in the U M * direction. See also BspKnotReverse. M * * * PARAMETERS: M * Srf: To be reversed. M * * * RETURN VALUE: M * CagdSrfStruct *: Reversed surface of Srf. M * * * SEE ALSO: M * CagdSrfReverse2, CagdSrfReverseDir M * * * KEYWORDS: M * CagdSrfReverse, reverse M *****************************************************************************/ CagdSrfStruct *CagdSrfReverse(CagdSrfStruct *Srf) { return CagdSrfReverseDir(Srf, CAGD_CONST_U_DIR); } /***************************************************************************** * DESCRIPTION: M * Returns a new surface that is the reversed surface of Srf by reversing the M * control mesh and the knot vector (if Bspline surface) of Srf in the Dir M * direction. See also BspKnotReverse. M * * * PARAMETERS: M * Srf: To be reversed. M * Dir: Direction to reverse the Mesh along. Either U or V. M * * * RETURN VALUE: M * CagdSrfStruct *: Reversed surface of Srf. M * * * SEE ALSO: M * CagdSrfReverse2 M * * * KEYWORDS: M * CagdSrfReverseDir, reverse M *****************************************************************************/ CagdSrfStruct *CagdSrfReverseDir(CagdSrfStruct *Srf, CagdSrfDirType Dir) { CagdBType IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf); int i, Len, Row, Col, ULength = Srf -> ULength, VLength = Srf -> VLength, MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType); CagdSrfStruct *ReversedSrf = CagdSrfCopy(Srf); CagdRType *KV, **Points = ReversedSrf -> Points; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: case CAGD_SBSPLINE_TYPE: break; case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } /* Reverse the Mesh: */ switch (Dir) { case CAGD_CONST_U_DIR: Len = (ULength >> 1); for (Row = 0; Row < VLength; Row++) for (Col = 0; Col < Len; Col++) for (i = IsNotRational; i <= MaxCoord; i++) SWAP(CagdRType, Points[i][Row * ULength + Col], Points[i][Row * ULength + ULength - Col - 1]); break; case CAGD_CONST_V_DIR: Len = (VLength >> 1); for (Col = 0; Col < ULength; Col++) for (Row = 0; Row < Len; Row++) for (i = IsNotRational; i <= MaxCoord; i++) SWAP(CagdRType, Points[i][Row * ULength + Col], Points[i][(VLength - Row - 1) * ULength + Col]); break; default: CAGD_FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV); } /* Reverse the U knot vector if it exists: */ if (Srf -> GType == CAGD_SBSPLINE_TYPE) { switch (Dir) { case CAGD_CONST_U_DIR: KV = BspKnotReverse(Srf -> UKnotVector, Srf -> UOrder + CAGD_SRF_UPT_LST_LEN(Srf)); IritFree(ReversedSrf -> UKnotVector); ReversedSrf -> UKnotVector = KV; break; case CAGD_CONST_V_DIR: KV = BspKnotReverse(Srf -> VKnotVector, Srf -> VOrder + CAGD_SRF_VPT_LST_LEN(Srf)); IritFree(ReversedSrf -> VKnotVector); ReversedSrf -> VKnotVector = KV; break; } } return ReversedSrf; } /***************************************************************************** * DESCRIPTION: M * Returns a new surface that is the reversed surface of Srf by flipping the M * U and the V directions of the surface. M * See also BspKnotReverse. M * * * PARAMETERS: M * Srf: To be reversed. M * * * RETURN VALUE: M * CagdSrfStruct *: Reversed surface of Srf. M * * * SEE ALSO: M * CagdSrfReverse, CagdSrfReverseDir M * * * KEYWORDS: M * CagdSrfReverse2, reverse M *****************************************************************************/ CagdSrfStruct *CagdSrfReverse2(CagdSrfStruct *Srf) { CagdBType IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf); int i, Row, Col, ULength = Srf -> ULength, VLength = Srf -> VLength, MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType); CagdSrfStruct *ReversedSrf = CagdSrfCopy(Srf); CagdRType **Points = Srf -> Points, **RevPoints = ReversedSrf -> Points; switch (Srf -> GType) { case CAGD_SBEZIER_TYPE: case CAGD_SBSPLINE_TYPE: break; case CAGD_SPOWER_TYPE: CAGD_FATAL_ERROR(CAGD_ERR_POWER_NO_SUPPORT); return NULL; default: CAGD_FATAL_ERROR(CAGD_ERR_UNDEF_SRF); return NULL; } /* Reverse the Mesh: */ for (Row = 0; Row < VLength; Row++) for (Col = 0; Col < ULength; Col++) for (i = IsNotRational; i <= MaxCoord; i++) RevPoints[i][Col * VLength + Row] = Points[i][Row * ULength + Col]; /* Swap the U and the V knot vectors if the exists: */ if (Srf -> GType == CAGD_SBSPLINE_TYPE) { SWAP(CagdRType *, ReversedSrf -> UKnotVector, ReversedSrf -> VKnotVector); } /* And swap the orders and lengths. */ SWAP(int, ReversedSrf -> UOrder, ReversedSrf -> VOrder); SWAP(int, ReversedSrf -> ULength, ReversedSrf -> VLength); return ReversedSrf; }