X-Git-Url: http://git.salome-platform.org/gitweb/?a=blobdiff_plain;f=src%2FINTERP_KERNEL%2FTransformedTriangleIntersect.cxx;h=36727ef2624ebfce5ed18c859c807574e8966583;hb=1b746b38f3cdeae6654a9501f37fde5e56e59288;hp=fa4c8a53581e8b5498cd5766f531a38ea0964b78;hpb=56fddf07c0b7170f79791d38e2b909a8a5b0b872;p=tools%2Fmedcoupling.git diff --git a/src/INTERP_KERNEL/TransformedTriangleIntersect.cxx b/src/INTERP_KERNEL/TransformedTriangleIntersect.cxx index fa4c8a535..36727ef26 100644 --- a/src/INTERP_KERNEL/TransformedTriangleIntersect.cxx +++ b/src/INTERP_KERNEL/TransformedTriangleIntersect.cxx @@ -1,4 +1,4 @@ -// Copyright (C) 2007-2015 CEA/DEN, EDF R&D +// Copyright (C) 2007-2024 CEA, EDF // // This library is free software; you can redistribute it and/or // modify it under the terms of the GNU Lesser General Public @@ -28,10 +28,10 @@ namespace INTERP_KERNEL { // ---------------------------------------------------------------------------------- - // Correspondance tables describing all the variations of formulas. + // Correspondence tables describing all the variations of formulas. // ---------------------------------------------------------------------------------- - /// \brief Correspondance between facets and double products. + /// \brief Correspondence between facets and double products. /// /// This table encodes Grandy, table IV. Use 3*facet + {0,1,2} as index const TransformedTriangle::DoubleProduct TransformedTriangle::DP_FOR_SEG_FACET_INTERSECTION[12] = @@ -78,7 +78,7 @@ namespace INTERP_KERNEL 9, 10, 11 // XYZ }; - /// \brief Correspondance edge - corners. + /// \brief Correspondence edge - corners. /// /// Gives the two corners associated with each edge /// Use 2*edge + {0, 1} as index @@ -92,7 +92,7 @@ namespace INTERP_KERNEL Z, X // ZX }; - /// \brief Correspondance edge - facets. + /// \brief Correspondence edge - facets. /// /// Gives the two facets shared by and edge. Use 2*facet + {0, 1} as index const TransformedTriangle::TetraFacet TransformedTriangle::FACET_FOR_EDGE[12] = @@ -105,7 +105,7 @@ namespace INTERP_KERNEL OZX, XYZ // ZX }; - /// \brief Correspondance corners - edges. + /// \brief Correspondence corners - edges. /// /// Gives edges meeting at a given corner. Use 3*corner + {0,1,2} as index const TransformedTriangle::TetraEdge TransformedTriangle::EDGES_FOR_CORNER[12] = @@ -168,13 +168,12 @@ namespace INTERP_KERNEL const double alpha = tA / (tA - tB); // calculate point - LOG(4, "corner A = " << corners[0] << " corner B = " << corners[1] ); + LOG(4, "corner A = " << strTC(corners[0]) << " corner B = " << strTC(corners[1]) ); LOG(4, "tA = " << tA << " tB = " << tB << " alpha= " << alpha ); for(int i = 0; i < 3; ++i) { - pt[i] = (1 - alpha) * COORDS_TET_CORNER[3*corners[0] + i] + - alpha * COORDS_TET_CORNER[3*corners[1] + i]; + pt[i] = (1 - alpha)*COORDS_TET_CORNER[3*corners[0] + i] + alpha*COORDS_TET_CORNER[3*corners[1] + i]; #if 0 pt[i] = (1 - alpha) * getCoordinateForTetCorner() + alpha * getCoordinateForTetCorner(); @@ -234,7 +233,7 @@ namespace INTERP_KERNEL /** * Tests if the given segment of the triangle intersects the given edge of the tetrahedron (Grandy, eq. [20] - * If the OPTIMIZE is defined, it does not do the test the double product that should be zero. + * * @param seg segment of the triangle * @param edge edge of tetrahedron * @return true if the segment intersects the edge @@ -249,13 +248,13 @@ namespace INTERP_KERNEL for(int i = 0 ; i < 2 ; ++i) { facet[i] = FACET_FOR_EDGE[2*edge + i]; - + // find the two c-values -> the two for the other edges of the facet int idx1 = 0 ; int idx2 = 1; DoubleProduct dp1 = DP_FOR_SEG_FACET_INTERSECTION[3*facet[i] + idx1]; DoubleProduct dp2 = DP_FOR_SEG_FACET_INTERSECTION[3*facet[i] + idx2]; - + if(dp1 == DoubleProduct( edge )) { idx1 = 2; @@ -266,7 +265,7 @@ namespace INTERP_KERNEL idx2 = 2; dp2 = DP_FOR_SEG_FACET_INTERSECTION[3*facet[i] + idx2]; } - + const double c1 = SIGN_FOR_SEG_FACET_INTERSECTION[3*facet[i] + idx1]*calcStableC(seg, dp1); const double c2 = SIGN_FOR_SEG_FACET_INTERSECTION[3*facet[i] + idx2]*calcStableC(seg, dp2); @@ -343,28 +342,24 @@ namespace INTERP_KERNEL for(int j = 0 ; j < 2; ++j) { const int dpIdx = DP_INDEX[3*facets[j] + i]; - const DoubleProduct dp = DP_FOR_SEG_FACET_INTERSECTION[dpIdx]; - const double sign = SIGN_FOR_SEG_FACET_INTERSECTION[dpIdx]; - c[j] = dpIdx < 0.0 ? 0.0 : sign * calcStableC(seg, dp); + if (dpIdx != -1) + { + const DoubleProduct dp = DP_FOR_SEG_FACET_INTERSECTION[dpIdx]; + const double sign = SIGN_FOR_SEG_FACET_INTERSECTION[dpIdx]; + c[j] = sign * calcStableC(seg, dp); + } + else + c[j] = 0.0; } - - // pt[i] = (c1*s1 + c2*s2) / (s1^2 + s2^2) pt[i] = (c[0] * s[0] + c[1] * s[1]) / denominator; - - // strange bug with -O2 enabled : assertion fails when we don't have the following - // trace - line - //std::cout << "pt[i] = " << pt[i] << std::endl; - //assert(pt[i] >= 0.0); // check we are in tetraeder - //assert(pt[i] <= 1.0); - } } /** * Tests if the given segment of the triangle intersects the given corner of the tetrahedron. - * (Grandy, eq. [21]). If OPTIMIZE is defined, the double products that should be zero are not verified. + * (Grandy, eq. [21]). * * @param seg segment of the triangle * @param corner corner of the tetrahedron @@ -372,8 +367,6 @@ namespace INTERP_KERNEL */ bool TransformedTriangle::testSegmentCornerIntersection(const TriSegment seg, const TetraCorner corner) const { - - // facets meeting at a given corner static const TetraFacet FACETS_FOR_CORNER[12] = { @@ -388,9 +381,7 @@ namespace INTERP_KERNEL { const TetraFacet facet = FACETS_FOR_CORNER[3*corner + i]; if(testSegmentIntersectsFacet(seg, facet)) - { - return true; - } + return true; } return false; @@ -432,11 +423,10 @@ namespace INTERP_KERNEL }; const TetraFacet facet = FACET_FOR_HALFSTRIP_INTERSECTION[edgeIndex]; - // special case : facet H = 0 const bool cond2 = (facet == NO_TET_FACET) ? testSegmentIntersectsHPlane(seg) : testSegmentIntersectsFacet(seg, facet); - LOG(4, "Halfstrip tests (" << seg << ", " << edge << ") : " << (cVals[0]*cVals[1] < 0.0) << ", " << cond2 << ", " << (cVals[2]*cVals[3] > 0.0) ); + LOG(4, "Halfstrip tests (" << strTriS(seg) << ", " << strTE(edge) << ") : " << (cVals[0]*cVals[1] < 0.0) << ", " << cond2 << ", " << (cVals[2]*cVals[3] > 0.0) ); LOG(4, "c2 = " << cVals[2] << ", c3 = " << cVals[3] ); return (cVals[0]*cVals[1] < 0.0) && cond2 && (cVals[2]*cVals[3] > 0.0); @@ -497,7 +487,6 @@ namespace INTERP_KERNEL /** * Tests if the given segment of triangle PQR intersects the ray pointing * in the upwards z - direction from the given corner of the tetrahedron. (Grandy eq. [29]) - * If OPTIMIZE is defined, the double product that should be zero is not verified. * * @param seg segment of the triangle PQR * @param corner corner of the tetrahedron on the h = 0 facet (X, Y, or Z) @@ -573,6 +562,7 @@ namespace INTERP_KERNEL // if two or more c-values are zero we disallow x-edge intersection // Grandy, p.446 + // test for == 0.0 is OK here, if you look at the way double products were pre-comuted: const int numZeros = (cPQ == 0.0 ? 1 : 0) + (cQR == 0.0 ? 1 : 0) + (cRP == 0.0 ? 1 : 0); if(numZeros >= 2 )