/****************************************************************************** * Bsp2Poly.c - Bezier to polygon/polylines conversion routines. * ******************************************************************************* * Written by Gershon Elber, Mar. 90. * ******************************************************************************/ #include "cagd_loc.h" static CagdPolygonStruct *BspC1Srf2Polygons(CagdSrfStruct *Srf, int FineNess, CagdBType ComputeNormals, CagdBType FourPerFlat, CagdBType ComputeUV); /***************************************************************************** * Routine to convert a single bspline surface to set of triangles * * approximating it. FineNess is a finess control on result and the bigger it * * is more triangles may result. a value of 10 is a good start value. * * NULL is returned in case of an error, otherwise list of CagdPolygonStruct. * * This routine looks for C1 discontinuities in the surface and split it * * into C1 continuous patches to invoke BspC1Srf2Polygons to gen. polygons. * *****************************************************************************/ CagdPolygonStruct *BspSrf2Polygons(CagdSrfStruct *Srf, int FineNess, CagdBType ComputeNormals, CagdBType FourPerFlat, CagdBType ComputeUV) { CagdRType u, v; int UOrder = Srf -> UOrder, VOrder = Srf -> VOrder, ULength = Srf -> ULength, VLength = Srf -> VLength; CagdBType HasUDiscont = BspKnotC1Discont(Srf -> UKnotVector, UOrder, ULength, &u), HasVDiscont = BspKnotC1Discont(Srf -> VKnotVector, VOrder, VLength, &v); CagdPolygonStruct *Poly; if (HasUDiscont || HasVDiscont) { CagdSrfStruct *Srf1 = HasUDiscont ? BspSrfSubdivAtParam(Srf, u, CAGD_CONST_U_DIR) : BspSrfSubdivAtParam(Srf, v, CAGD_CONST_V_DIR), *Srf2 = Srf1 -> Pnext; CagdPolygonStruct *Poly1 = BspSrf2Polygons(Srf1, FineNess, ComputeNormals, FourPerFlat, ComputeUV), *Poly2 = BspSrf2Polygons(Srf2, FineNess, ComputeNormals, FourPerFlat, ComputeUV); CagdSrfFreeList(Srf1); /* Chain the two lists together: */ for (Poly = Poly1; Poly -> Pnext != NULL; Poly = Poly -> Pnext); Poly -> Pnext = Poly2; Poly = Poly1; } else Poly = BspC1Srf2Polygons(Srf, FineNess, ComputeNormals, FourPerFlat, ComputeUV); return Poly; } /***************************************************************************** * Routine to convert a single C1 continuouse bspline srf to set of triangles * * approximating it. FineNess is a finess control on result and the bigger it * * is more triangles may result. a value of 10 is a good start value. * * NULL is returned in case of an error, otherwise list of CagdPolygonStruct. * *****************************************************************************/ static CagdPolygonStruct *BspC1Srf2Polygons(CagdSrfStruct *Srf, int FineNess, CagdBType ComputeNormals, CagdBType FourPerFlat, CagdBType ComputeUV) { int i, j, FineNessU1, FineNessV1, FineNessU, FineNessV, BaseIndex; CagdRType u, v, UMin, UMax, VMin, VMax, *Pt; CagdPointType PType = Srf -> PType; CagdPtStruct PtCenter, *Pt1, *Pt2, *Pt3, *Pt4, *PtMesh, *PtMeshPtr; CagdUVStruct UVCenter, *UVMeshPtr, *UV1 = NULL, *UV2 = NULL, *UV3 = NULL, *UV4 = NULL, *UVMesh = NULL; CagdVecStruct NlCenter, *Nl1 = NULL, *Nl2 = NULL, *Nl3 = NULL, *Nl4 = NULL, *PtNrml = NULL; CagdCrvStruct *Crv; CagdPolygonStruct *Poly, *PolyHead = NULL; if (!CAGD_IS_BSPLINE_SRF(Srf)) return NULL; /* Simple heuristic to estimate how many samples to compute. */ FineNessU = Srf -> ULength * FineNess / 10; FineNessV = Srf -> VLength * FineNess / 10; if (FineNessU < 2) FineNessU = 2; if (FineNessV < 2) FineNessV = 2; switch (_CagdLin2Poly) { case CAGD_REG_POLY_PER_LIN: break; case CAGD_ONE_POLY_PER_LIN: if (Srf -> UOrder == 2) FineNessU = 2; if (Srf -> VOrder == 2) FineNessV = 2; break; case CAGD_ONE_POLY_PER_COLIN: break; } FineNessU1 = FineNessU - 1; FineNessV1 = FineNessV - 1; /* Current to surface property such as curvature is used as subdivison */ /* criterion and the surface is subdivided, equally spaced in parametric */ /* space, using FineNess as number of subdivisions per axis. */ /* Allocate a mesh to hold all vertices so common vertices need not be */ /* Evaluated twice, and evaluate the surface at these mesh points. */ PtMeshPtr = PtMesh = (CagdPtStruct *) IritMalloc(FineNessU * FineNessV * sizeof(CagdPtStruct)); if (ComputeUV) { UVMeshPtr = UVMesh = (CagdUVStruct *) IritMalloc(FineNessU * FineNessV * sizeof(CagdUVStruct)); } BspSrfDomain(Srf, &UMin, &UMax, &VMin, &VMax); for (i = 0; i < FineNessU; i++) { u = UMin + (UMax - UMin) * i / ((CagdRType) FineNessU1); Crv = BspSrfCrvFromSrf(Srf, u, CAGD_CONST_U_DIR); for (j = 0; j < FineNessV; j++, PtMeshPtr++) { v = VMin + (VMax - VMin) * j / ((CagdRType) FineNessV1); Pt = BspCrvEvalAtParam(Crv, v); CagdCoerceToE3(PtMeshPtr -> Pt, &Pt, -1, PType); if (ComputeUV) { UVMeshPtr -> UV[0] = u; UVMeshPtr -> UV[1] = v; UVMeshPtr++; } } CagdCrvFree(Crv); } if (ComputeNormals) PtNrml = BspSrfMeshNormals(Srf, FineNessU, FineNessV); /* Now that we have the mesh, create the polygons. */ for (i = 0; i < FineNessU1; i++) for (j = 0; j < FineNessV1; j++) { BaseIndex = i * FineNessV + j; Pt1 = &PtMesh[BaseIndex]; /* Cache the four flat corners. */ Pt2 = &PtMesh[BaseIndex + 1]; Pt3 = &PtMesh[BaseIndex + FineNessV + 1]; Pt4 = &PtMesh[BaseIndex + FineNessV]; if (ComputeNormals) { Nl1 = &PtNrml[BaseIndex]; Nl2 = &PtNrml[BaseIndex + 1]; Nl3 = &PtNrml[BaseIndex + FineNessV + 1]; Nl4 = &PtNrml[BaseIndex + FineNessV]; } if (ComputeUV) { UV1 = &UVMesh[BaseIndex]; UV2 = &UVMesh[BaseIndex + 1]; UV3 = &UVMesh[BaseIndex + FineNessV + 1]; UV4 = &UVMesh[BaseIndex + FineNessV]; } if (FourPerFlat) { /* Eval middle point and create 4 triangles. */ CAGD_COPY_POINT(PtCenter, *Pt1); CAGD_ADD_POINT(PtCenter, *Pt2); CAGD_ADD_POINT(PtCenter, *Pt3); CAGD_ADD_POINT(PtCenter, *Pt4); CAGD_MULT_POINT(PtCenter, 0.25); if (ComputeNormals) { /* Average the four normals to find the middle one. */ CAGD_COPY_VECTOR(NlCenter, *Nl1); CAGD_ADD_VECTOR(NlCenter, *Nl2); CAGD_ADD_VECTOR(NlCenter, *Nl3); CAGD_ADD_VECTOR(NlCenter, *Nl4); CAGD_NORMALIZE_VECTOR(NlCenter); } if (ComputeUV) { UVCenter.UV[0] = (UV1 -> UV[0] + UV2 -> UV[0] + UV3 -> UV[0] + UV4 -> UV[0]) / 4.0; UVCenter.UV[1] = (UV1 -> UV[1] + UV2 -> UV[1] + UV3 -> UV[1] + UV4 -> UV[1]) / 4.0; } Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt1, Pt2, &PtCenter, Nl1, Nl2, &NlCenter, UV1, UV2, &UVCenter); CAGD_LIST_PUSH(Poly, PolyHead); Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt2, Pt3, &PtCenter, Nl2, Nl3, &NlCenter, UV2, UV3, &UVCenter); CAGD_LIST_PUSH(Poly, PolyHead); Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt3, Pt4, &PtCenter, Nl3, Nl4, &NlCenter, UV3, UV4, &UVCenter); CAGD_LIST_PUSH(Poly, PolyHead); Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt4, Pt1, &PtCenter, Nl4, Nl1, &NlCenter, UV4, UV1, &UVCenter); CAGD_LIST_PUSH(Poly, PolyHead); } else { /* Only two along the diagonal... */ Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt1, Pt2, Pt3, Nl1, Nl2, Nl3, UV1, UV2, UV3); CAGD_LIST_PUSH(Poly, PolyHead); Poly = _CagdMakePolygon(ComputeNormals, ComputeUV, Pt3, Pt4, Pt1, Nl3, Nl4, Nl1, UV3, UV4, UV1); CAGD_LIST_PUSH(Poly, PolyHead); } } IritFree((VoidPtr) PtMesh); if (ComputeNormals) IritFree((VoidPtr) PtNrml); if (ComputeUV) IritFree((VoidPtr) UVMesh); return PolyHead; } /***************************************************************************** * Routine to convert a single bspline surface to NumOfIsolines polylines list* * in each param. direction with SamplesPerCurve in each isoparametric curve. * * Polyline are always E3 of CagdPolylineStruct type. * * Iso parametric curves are sampled equally spaced in parametric space. * * NULL is returned in case of an error, otherwise list of CagdPolylineStruct.* * Attempt is made to extract isolines along C1 discontinuities first. * *****************************************************************************/ CagdPolylineStruct *BspSrf2Polylines(CagdSrfStruct *Srf, int NumOfIsocurves[2], int SamplesPerCurve) { int i, NumC1Disconts, NumOfIsos, IsRegion, n = 1 << SamplesPerCurve, ULength = Srf -> ULength, VLength = Srf -> VLength, UOrder = Srf -> UOrder, VOrder = Srf -> VOrder; CagdRType u, v, UMin, UMax, VMin, VMax, *C1Disconts, *IsoParams, *RefKV, *UKV = Srf -> UKnotVector, *VKV = Srf -> VKnotVector; CagdCrvStruct *Crv; CagdPolylineStruct *Poly, *PolyList = NULL; BspKnotAlphaCoeffType *A; if (!CAGD_IS_BSPLINE_SRF(Srf)) return NULL; /* Make sure the curve is open. We move 2 Epsilons to make sure region */ /* extraction will occur. Otherwise the curve will be copied as is. */ if (!BspKnotHasOpenEC(UKV, ULength, UOrder) || !BspKnotHasOpenEC(VKV, VLength, VOrder)) { CagdSrfStruct *TSrf = CagdSrfRegionFromSrf(Srf, UKV[UOrder - 1] + EPSILON * 2, UKV[ULength] - EPSILON * 2, CAGD_CONST_U_DIR); Srf = CagdSrfRegionFromSrf(TSrf, VKV[VOrder - 1] + EPSILON * 2, VKV[VLength] - EPSILON * 2, CAGD_CONST_V_DIR); CagdSrfFree(TSrf); IsRegion = TRUE; } else IsRegion = FALSE; /* Make sure requested format is something reasonable. */ if (SamplesPerCurve < 1) SamplesPerCurve = 1; if (SamplesPerCurve > CAGD_MAX_BEZIER_CACHE_ORDER) SamplesPerCurve = CAGD_MAX_BEZIER_CACHE_ORDER; if (NumOfIsocurves[0] < 2) NumOfIsocurves[0] = 2; if (NumOfIsocurves[1] <= 0) NumOfIsocurves[1] = NumOfIsocurves[0]; else if (NumOfIsocurves[1] < 2) NumOfIsocurves[1] = 2; BspSrfDomain(Srf, &UMin, &UMax, &VMin, &VMax); /* Compute discontinuities along the u axis and use that to determine */ /* where to extract isolines along u. */ /* Note C1Disconts is freed by BspKnotParamValues. */ /* Add heuristically more samples if surface has interior knots. */ NumOfIsos = NumOfIsocurves[0]; if (UOrder > 2) NumOfIsos += (ULength - UOrder) / 2; C1Disconts = BspKnotAllC1Discont(Srf -> UKnotVector, UOrder, ULength, &NumC1Disconts); IsoParams = BspKnotParamValues(UMin, UMax, NumOfIsos, C1Disconts, NumC1Disconts); RefKV = BspKnotPrepEquallySpaced(MAX(n - VLength, 1), VMin, VMax); A = BspKnotEvalAlphaCoefMerge(VOrder, Srf -> VKnotVector, VLength, RefKV, MAX(n - VLength, 1)); IritFree((VoidPtr) RefKV); for (i = 0; i < NumOfIsos; i++) { u = IsoParams[i]; Crv = BspSrfCrvFromSrf(Srf, u, CAGD_CONST_U_DIR); Poly = BspCrv2Polyline(Crv, SamplesPerCurve, A); Poly -> Pnext = PolyList; PolyList = Poly; CagdCrvFree(Crv); } IritFree((VoidPtr) IsoParams); BspKnotFreeAlphaCoef(A); /* Compute discontinuities along the v axis and use that to determine */ /* where to extract isolines along v. */ /* Note C1Disconts is freed by BspKnotParamValues. */ /* Add heuristically more samples if surface has interior knots. */ NumOfIsos = NumOfIsocurves[1]; if (VOrder > 2) NumOfIsos += (VLength - VOrder) / 2; C1Disconts = BspKnotAllC1Discont(Srf -> VKnotVector, VOrder, VLength, &NumC1Disconts); IsoParams = BspKnotParamValues(VMin, VMax, NumOfIsos, C1Disconts, NumC1Disconts); RefKV = BspKnotPrepEquallySpaced(MAX(n - ULength, 1), UMin, UMax); A = BspKnotEvalAlphaCoefMerge(UOrder, Srf -> UKnotVector, ULength, RefKV, MAX(n - ULength, 1)); IritFree((VoidPtr) RefKV); for (i = 0; i < NumOfIsos; i++) { v = IsoParams[i]; Crv = BspSrfCrvFromSrf(Srf, v, CAGD_CONST_V_DIR); Poly = BspCrv2Polyline(Crv, SamplesPerCurve, A); Poly -> Pnext = PolyList; PolyList = Poly; CagdCrvFree(Crv); } IritFree((VoidPtr) IsoParams); BspKnotFreeAlphaCoef(A); return PolyList; } /***************************************************************************** * Routine to convert a single bspline surface to NumOfIsoline isocurve list * * in each param. direction. * * Iso parametric curves are sampled equally spaced in parametric space. * * NULL is returned in case of an error, otherwise list of CagdCrvStruct. * *****************************************************************************/ CagdCrvStruct *BspSrf2Curves(CagdSrfStruct *Srf, int NumOfIsocurves[2]) { int i, NumC1Disconts, ULength = Srf -> ULength, VLength = Srf -> VLength, UOrder = Srf -> UOrder, VOrder = Srf -> VOrder; CagdRType u, v, UMin, UMax, VMin, VMax, *C1Disconts, *IsoParams; CagdCrvStruct *Crv, *CrvList = NULL; if (!CAGD_IS_BSPLINE_SRF(Srf)) return NULL; /* Make sure requested format is something reasonable. */ if (NumOfIsocurves[0] < 2) NumOfIsocurves[0] = 2; if (NumOfIsocurves[1] <= 0) NumOfIsocurves[1] = NumOfIsocurves[0]; else if (NumOfIsocurves[1] < 2) NumOfIsocurves[1] = 2; BspSrfDomain(Srf, &UMin, &UMax, &VMin, &VMax); /* Compute discontinuities along the u axis and use that to determine */ /* where to extract isolines along u. */ /* Note C1Disconts is freed by BspKnotParamValues. */ C1Disconts = BspKnotAllC1Discont(Srf -> UKnotVector, UOrder, ULength, &NumC1Disconts); IsoParams = BspKnotParamValues(UMin, UMax, NumOfIsocurves[0], C1Disconts, NumC1Disconts); for (i = 0; i < NumOfIsocurves[0]; i++) { u = IsoParams[i]; Crv = BspSrfCrvFromSrf(Srf, u, CAGD_CONST_U_DIR); Crv -> Pnext = CrvList; CrvList = Crv; } IritFree((VoidPtr) IsoParams); /* Compute discontinuities along the v axis and use that to determine */ /* where to extract isolines along v. */ /* Note C1Disconts is freed by BspKnotParamValues. */ C1Disconts = BspKnotAllC1Discont(Srf -> VKnotVector, VOrder, VLength, &NumC1Disconts); IsoParams = BspKnotParamValues(VMin, VMax, NumOfIsocurves[1], C1Disconts, NumC1Disconts); for (i = 0; i < NumOfIsocurves[1]; i++) { v = IsoParams[i]; Crv = BspSrfCrvFromSrf(Srf, v, CAGD_CONST_V_DIR); Crv -> Pnext = CrvList; CrvList = Crv; } IritFree((VoidPtr) IsoParams); return CrvList; } /***************************************************************************** * Routine to approx. a single bspline curve as a polyline with * * 2^SamplesPerCurve samples. Polyline is always E3 CagdPolylineStruct type. * * Curve is refined equally spaced in parametric space, unless the curve is * * linear in which the control polygon is simply being copied. * * If A is specified, it is used to refine the curve. * * NULL is returned in case of an error, otherwise CagdPolylineStruct. * *****************************************************************************/ CagdPolylineStruct *BspCrv2Polyline(CagdCrvStruct *Crv, int SamplesPerCurve, BspKnotAlphaCoeffType *A) { int i, j, n, IsRegion, Order = Crv -> Order, Len = Crv -> Length, IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv), MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType); CagdRType *Polyline[CAGD_MAX_PT_SIZE], Scaler, *KV = Crv -> KnotVector; CagdPtStruct *NewPolyline; CagdPolylineStruct *P; if (!CAGD_IS_BSPLINE_CRV(Crv)) return NULL; /* Make sure the curve is open. We move 2 Epsilons to make sure region */ /* extraction will occur. Otherwise the curve will be copied as is. */ if (!BspKnotHasOpenEC(KV, Len, Order)) { Crv = CagdCrvRegionFromCrv(Crv, KV[Order - 1] + EPSILON * 2, KV[Len] - EPSILON * 2); IsRegion = TRUE; } else IsRegion = FALSE; /* Make sure requested format is something reasonable. */ if (SamplesPerCurve < 1) SamplesPerCurve = 1; if ((1 << SamplesPerCurve) <= Crv -> Length) { /* Make sure 2^SamplesPerCurve can hold the entire control polygon. */ for (i = 1, SamplesPerCurve = 0; i <= Crv -> Length; i <<= 1, SamplesPerCurve++); } n = MAX(A ? A -> RefLength : 0, 1 << SamplesPerCurve); P = CagdPolylineNew(n); NewPolyline = P -> Polyline; /* Allocate temporary memory to hold evaluated curve. */ for (i = 0; i < CAGD_MAX_PT_SIZE; i++) Polyline[i] = (CagdRType *) IritMalloc(sizeof(CagdRType) * n); if (MaxCoord > 3) MaxCoord = 3; n = P -> Length = CagdCrvEvalToPolyline(Crv, A == NULL ? SamplesPerCurve : 0, Polyline, A); for (i = n - 1; i >= 0; i--) { /* Convert to E3 polyline. */ if (IsNotRational) Scaler = 1.0; else Scaler = Polyline[0][i]; for (j = 0; j < MaxCoord; j++) NewPolyline[i].Pt[j] = Polyline[j+1][i] / Scaler; for (j = MaxCoord; j < 3; j++) NewPolyline[i].Pt[j] = 0.0; } for (i = 0; i < CAGD_MAX_PT_SIZE; i++) IritFree((VoidPtr) Polyline[i]); if (IsRegion) CagdCrvFree(Crv); return P; }