home *** CD-ROM | disk | FTP | other *** search
- /******************************************************************************
- * SBsp-Aux.c - Bspline surface auxilary routines. *
- *******************************************************************************
- * Written by Gershon Elber, July. 90. *
- ******************************************************************************/
-
- #ifdef __MSDOS__
- #include <stdlib.h>
- #endif /* __MSDOS__ */
-
- #include <ctype.h>
- #include <stdio.h>
- #include <string.h>
- #include "cagd_loc.h"
-
- /* Define some marcos to make some of the routines below look better. They */
- /* calculate the index of the U, V point of the control mesh in Points. */
- #define DERIVED_SRF(U, V) CAGD_MESH_UV(DerivedSrf, U, V)
- #define RAISED_SRF(U, V) CAGD_MESH_UV(RaisedSrf, U, V)
- #define SRF(U, V) CAGD_MESH_UV(Srf, U, V)
-
- /******************************************************************************
- * Given a bspline surface - subdivide it into two at given parametric value. *
- * Returns pointer to first surface in a list of two srfs (subdivided ones). *
- * The subdivision is exact result of evaluating the surface int a curve at t *
- * using the recursive algorithm - the left resulting points is left surface, *
- * and the right resulting points is right surface (left is below t). *
- ******************************************************************************/
- CagdSrfStruct *BspSrfSubdivAtParam(CagdSrfStruct *Srf, CagdRType t,
- CagdSrfDirType Dir)
- {
- int Row, Col,
- ULength1 = 0,
- ULength2 = 0,
- VLength1 = 0,
- VLength2 = 0,
- ULength = Srf -> ULength,
- VLength = Srf -> VLength,
- UOrder = Srf -> UOrder,
- VOrder = Srf -> VOrder;
- CagdCrvStruct *Crv, *LCrv, *RCrv;
- CagdSrfStruct
- *RSrf = NULL,
- *LSrf = NULL;
-
- switch (Dir) {
- case CAGD_CONST_U_DIR:
- for (Row = 0; Row < VLength; Row++) {
- Crv = BspSrfCrvFromMesh(Srf, Row, CAGD_CONST_V_DIR);
- LCrv = BspCrvSubdivAtParam(Crv, t);
- RCrv = LCrv -> Pnext;
-
- if (Row == 0) { /* Figure out surfaces size and allocate. */
- ULength1 = LCrv -> Length;
- ULength2 = RCrv -> Length;
- LSrf = BspSrfNew(ULength1, VLength, UOrder, VOrder,
- Srf -> PType);
- RSrf = BspSrfNew(ULength2, VLength, UOrder, VOrder,
- Srf -> PType);
- CagdFree((VoidPtr) RSrf -> UKnotVector);
- CagdFree((VoidPtr) RSrf -> VKnotVector);
- CagdFree((VoidPtr) LSrf -> UKnotVector);
- CagdFree((VoidPtr) LSrf -> VKnotVector);
- RSrf -> UKnotVector = BspKnotCopy(RCrv -> KnotVector,
- RCrv -> Length + RCrv -> Order);
- RSrf -> VKnotVector = BspKnotCopy(Srf -> VKnotVector,
- Srf -> VLength + Srf -> VOrder);
- LSrf -> UKnotVector = BspKnotCopy(LCrv -> KnotVector,
- LCrv -> Length + LCrv -> Order);
- LSrf -> VKnotVector = BspKnotCopy(Srf -> VKnotVector,
- Srf -> VLength + Srf -> VOrder);
- }
-
- CagdCrvToMesh(LCrv, Row, CAGD_CONST_V_DIR, LSrf);
- CagdCrvToMesh(RCrv, Row, CAGD_CONST_V_DIR, RSrf);
-
- CagdCrvFree(Crv);
- CagdCrvFree(LCrv);
- CagdCrvFree(RCrv);
- }
- break;
- case CAGD_CONST_V_DIR:
- for (Col = 0; Col < ULength; Col++) {
- Crv = BspSrfCrvFromMesh(Srf, Col, CAGD_CONST_U_DIR);
- LCrv = BspCrvSubdivAtParam(Crv, t);
- RCrv = LCrv -> Pnext;
-
- if (Col == 0) { /* Figure out surfaces size and allocate. */
- VLength1 = LCrv -> Length;
- VLength2 = RCrv -> Length;
- LSrf = BspSrfNew(ULength, VLength1, UOrder, VOrder,
- Srf -> PType);
- RSrf = BspSrfNew(ULength, VLength2, UOrder, VOrder,
- Srf -> PType);
- CagdFree((VoidPtr) RSrf -> UKnotVector);
- CagdFree((VoidPtr) RSrf -> VKnotVector);
- CagdFree((VoidPtr) LSrf -> UKnotVector);
- CagdFree((VoidPtr) LSrf -> VKnotVector);
- RSrf -> UKnotVector = BspKnotCopy(Srf -> UKnotVector,
- Srf -> ULength + Srf -> UOrder);
- RSrf -> VKnotVector = BspKnotCopy(RCrv -> KnotVector,
- RCrv -> Length + RCrv -> Order);
- LSrf -> UKnotVector = BspKnotCopy(Srf -> UKnotVector,
- Srf -> ULength + Srf -> UOrder);
- LSrf -> VKnotVector = BspKnotCopy(LCrv -> KnotVector,
- LCrv -> Length + LCrv -> Order);
- }
-
- CagdCrvToMesh(LCrv, Col, CAGD_CONST_U_DIR, LSrf);
- CagdCrvToMesh(RCrv, Col, CAGD_CONST_U_DIR, RSrf);
-
- CagdCrvFree(Crv);
- CagdCrvFree(LCrv);
- CagdCrvFree(RCrv);
- }
- break;
- default:
- FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV);
- break;
- }
-
- LSrf -> Pnext = RSrf;
- return LSrf;
- }
-
- /******************************************************************************
- * Return a new surface, identical to the original but with one degree higher *
- * in the given direction. *
- ******************************************************************************/
- CagdSrfStruct *BspSrfDegreeRaise(CagdSrfStruct *Srf, CagdSrfDirType Dir)
- {
- FATAL_ERROR(CAGD_ERR_NOT_IMPLEMENTED);
-
- return NULL;
- }
-
- /******************************************************************************
- * Return a new surface equal to the derived surface in the direction Dir. *
- * A test is made to make sure the original surface is not rational. *
- * Let old control polygon be P(i), i = 0 to k-1, and Q(i) be new one then: *
- * Q(i) = (k - 1) * P(i+1) - P(i), i = 0 to k-2. *
- * This is applied to all rows/cols of the surface. *
- ******************************************************************************/
- CagdSrfStruct *BspSrfDerive(CagdSrfStruct *Srf, CagdSrfDirType Dir)
- {
- CagdBType
- IsNotRational = !CAGD_IS_RATIONAL_SRF(Srf);
- int i, j, Row, Col,
- ULength = Srf -> ULength,
- VLength = Srf -> VLength,
- UOrder = Srf -> UOrder,
- VOrder = Srf -> VOrder,
- MaxCoord = CAGD_NUM_OF_PT_COORD(Srf -> PType);
- CagdSrfStruct
- *DerivedSrf = NULL;
-
-
- switch (Dir) {
- case CAGD_CONST_V_DIR:
- if (!IsNotRational || UOrder < 3)
- FATAL_ERROR(CAGD_ERR_RAT_LIN_NO_SUPPORT);
-
- DerivedSrf = BspSrfNew(ULength - 1, VLength,
- UOrder - 1, VOrder, Srf -> PType);
- GEN_COPY(DerivedSrf -> UKnotVector,
- &Srf -> UKnotVector[1],
- sizeof(CagdRType) * (ULength + UOrder - 2));
- GEN_COPY(DerivedSrf -> VKnotVector,
- Srf -> VKnotVector,
- sizeof(CagdRType) * (VLength + VOrder));
-
- for (Row = 0; Row < VLength; Row++)
- for (i = 0; i < ULength - 1; i++)
- for (j = IsNotRational; j <= MaxCoord; j++)
- DerivedSrf -> Points[j][DERIVED_SRF(i, Row)] =
- (ULength - 1) *
- (Srf -> Points[j][SRF(i + 1, Row)] -
- Srf -> Points[j][SRF(i, Row)]);
- break;
- case CAGD_CONST_U_DIR:
- if (!IsNotRational || VOrder < 3)
- FATAL_ERROR(CAGD_ERR_RAT_LIN_NO_SUPPORT);
-
- DerivedSrf = BspSrfNew(ULength, VLength - 1,
- UOrder, VOrder - 1, Srf -> PType);
- GEN_COPY(DerivedSrf -> UKnotVector,
- Srf -> UKnotVector,
- sizeof(CagdRType) * (ULength + UOrder));
- GEN_COPY(DerivedSrf -> VKnotVector,
- &Srf -> VKnotVector[1],
- sizeof(CagdRType) * (VLength + VOrder - 2));
-
- for (Col = 0; Col < ULength; Col++)
- for (i = 0; i < VLength - 1; i++)
- for (j = IsNotRational; j <= MaxCoord; j++)
- DerivedSrf -> Points[j][DERIVED_SRF(Col, i)] =
- (VLength - 1) *
- (Srf -> Points[j][SRF(Col, i + 1)] -
- Srf -> Points[j][SRF(Col, i)]);
- break;
- default:
- FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV);
- break;
- }
-
- return DerivedSrf;
- }
-
- /******************************************************************************
- * Evaluate the tangent to a surface at a given point and given direction. *
- ******************************************************************************/
- CagdVecStruct *BspSrfTangent(CagdSrfStruct *Srf, CagdRType u, CagdRType v,
- CagdSrfDirType Dir)
- {
- CagdVecStruct
- *Tangent = NULL;
- CagdCrvStruct *Crv;
-
- switch (Dir) {
- case CAGD_CONST_V_DIR:
- Crv = BspSrfCrvFromSrf(Srf, v, Dir);
- Tangent = BspCrvTangent(Crv, u);
- CagdCrvFree(Crv);
- break;
- case CAGD_CONST_U_DIR:
- Crv = BspSrfCrvFromSrf(Srf, u, Dir);
- Tangent = BspCrvTangent(Crv, v);
- CagdCrvFree(Crv);
- break;
- default:
- FATAL_ERROR(CAGD_ERR_DIR_NOT_CONST_UV);
- break;
- }
-
- return Tangent;
- }
-
- /******************************************************************************
- * Evaluate the normal of a surface at a given point. *
- ******************************************************************************/
- CagdVecStruct *BspSrfNormal(CagdSrfStruct *Srf, CagdRType u, CagdRType v)
- {
- static CagdVecStruct Normal;
- CagdVecStruct *V, T1, T2;
-
- V = BspSrfTangent(Srf, u, v, CAGD_CONST_U_DIR);
- CAGD_COPY_VECTOR(T1, *V);
-
- V = BspSrfTangent(Srf, u, v, CAGD_CONST_V_DIR);
- CAGD_COPY_VECTOR(T2, *V);
-
- /* The normal is the cross product of T1 and T2: */
- Normal.Vec[0] = T1.Vec[1] * T2.Vec[2] - T1.Vec[2] * T2.Vec[1];
- Normal.Vec[1] = T1.Vec[2] * T2.Vec[0] - T1.Vec[0] * T2.Vec[2];
- Normal.Vec[2] = T1.Vec[0] * T2.Vec[1] - T1.Vec[1] * T2.Vec[0];
-
- CAGD_NORMALIZE_VECTOR(Normal); /* Normalize the vector. */
-
- return &Normal;
- }
-
- /******************************************************************************
- * Evaluate the normals of a surface at a mesh defined by subdividing the *
- * parametric space int a grid of size UFineNess by VFineNess. *
- * The normals are saved in a linear CagdVecStruct vector which is allocated *
- * dynamically. Data is saved u inc. first. *
- * This routine much faster than evaluating normal per each point of such mesh.*
- ******************************************************************************/
- CagdVecStruct *BspSrfMeshNormals(CagdSrfStruct *Srf, int UFineNess,
- int VFineNess)
- {
- int i, j;
- CagdRType u, v, UMin, UMax, VMin, VMax;
- CagdVecStruct *Normals, *NPtr, *T, T1, T2;
- CagdCrvStruct **UCurves, **VCurves, *UCrv, *VCrv;
-
- BspSrfDomain(Srf, &UMin, &UMax, &VMin, &VMax);
- Normals = (CagdVecStruct *) CagdMalloc(sizeof(CagdVecStruct) * UFineNess *
- VFineNess);
-
- UCurves = (CagdCrvStruct **) CagdMalloc(sizeof(CagdCrvStruct *) *
- UFineNess);
- VCurves = (CagdCrvStruct **) CagdMalloc(sizeof(CagdCrvStruct *) *
- VFineNess);
-
- UFineNess--;
- VFineNess--;
- for (i = 0; i <= UFineNess; i++) /* Prepare Iso U curves. */
- UCurves[i] = BspSrfCrvFromSrf(Srf,
- UMin + (UMax - UMin) * ((CagdRType) i) / UFineNess,
- CAGD_CONST_U_DIR);
- for (j = 0; j <= VFineNess; j++) /* Prepare Iso V curves. */
- VCurves[j] = BspSrfCrvFromSrf(Srf,
- VMin + (VMax - VMin) * ((CagdRType) j) / VFineNess,
- CAGD_CONST_V_DIR);
-
- NPtr = Normals;
- for (i = 0; i <= UFineNess; i++) {
- UCrv = UCurves[i];
- u = UMin + (UMax - UMin) * ((CagdRType) i) / UFineNess;
-
- for (j = 0; j <= VFineNess; j++) {
- VCrv = VCurves[j];
- v = VMin + (VMax - VMin) * ((CagdRType) j) / VFineNess;
-
- /* We need to copy the tangents as BspCrvTangent save it in as */
- /* static so the second call will overwrite first call value. */
- T = BspCrvTangent(UCrv, v);
- CAGD_COPY_VECTOR(T1, *T);
- T = BspCrvTangent(VCrv, u);
- CAGD_COPY_VECTOR(T2, *T);
-
- /* The normal is the cross product of T1 and T2: */
- NPtr -> Vec[0] = T1.Vec[1] * T2.Vec[2] - T1.Vec[2] * T2.Vec[1];
- NPtr -> Vec[1] = T1.Vec[2] * T2.Vec[0] - T1.Vec[0] * T2.Vec[2];
- NPtr -> Vec[2] = T1.Vec[0] * T2.Vec[1] - T1.Vec[1] * T2.Vec[0];
-
- CAGD_NORMALIZE_VECTOR(*NPtr); /* Normalize the vector. */
- NPtr++;
- }
- }
-
- for (i = 0; i <= UFineNess; i++) CagdCrvFree(UCurves[i]);
- CagdFree(UCurves);
- for (j = 0; j <= VFineNess; j++) CagdCrvFree(VCurves[j]);
- CagdFree(VCurves);
-
- return Normals;
- }
-