Warning, /include/opencascade/SelectMgr_Frustum.lxx is written in an unsupported language. File is not indexed.
0001 // Created on: 2015-03-16
0002 // Created by: Varvara POSKONINA
0003 // Copyright (c) 2005-2014 OPEN CASCADE SAS
0004 //
0005 // This file is part of Open CASCADE Technology software library.
0006 //
0007 // This library is free software; you can redistribute it and/or modify it under
0008 // the terms of the GNU Lesser General Public License version 2.1 as published
0009 // by the Free Software Foundation, with special exception defined in the file
0010 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
0011 // distribution for complete text of the license and disclaimer of any warranty.
0012 //
0013 // Alternatively, this file may be used under the terms of Open CASCADE
0014 // commercial license or contractual agreement.
0015
0016 #include <gp_Pln.hxx>
0017 #include <NCollection_Vector.hxx>
0018 #include <Poly_Array1OfTriangle.hxx>
0019 #include <Standard_Assert.hxx>
0020 #include <SelectMgr_FrustumBuilder.hxx>
0021
0022 // =======================================================================
0023 // function : isSeparated
0024 // purpose : Checks if AABB and frustum are separated along the given axis.
0025 // =======================================================================
0026 template <int N>
0027 Standard_Boolean SelectMgr_Frustum<N>::isSeparated(const SelectMgr_Vec3& theBoxMin,
0028 const SelectMgr_Vec3& theBoxMax,
0029 const gp_XYZ& theDirect,
0030 Standard_Boolean* theInside) const
0031 {
0032 const Standard_Real aMinB =
0033 theDirect.X() * (theDirect.X() < 0.0 ? theBoxMax.x() : theBoxMin.x())
0034 + theDirect.Y() * (theDirect.Y() < 0.0 ? theBoxMax.y() : theBoxMin.y())
0035 + theDirect.Z() * (theDirect.Z() < 0.0 ? theBoxMax.z() : theBoxMin.z());
0036
0037 const Standard_Real aMaxB =
0038 theDirect.X() * (theDirect.X() < 0.0 ? theBoxMin.x() : theBoxMax.x())
0039 + theDirect.Y() * (theDirect.Y() < 0.0 ? theBoxMin.y() : theBoxMax.y())
0040 + theDirect.Z() * (theDirect.Z() < 0.0 ? theBoxMin.z() : theBoxMax.z());
0041
0042 Standard_ASSERT_RAISE(aMaxB >= aMinB, "Error! Failed to project box");
0043
0044 // frustum projection
0045 Standard_Real aMinF = DBL_MAX;
0046 Standard_Real aMaxF = -DBL_MAX;
0047
0048 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0049 {
0050 const Standard_Real aProj = myVertices[aVertIdx].XYZ().Dot(theDirect);
0051
0052 aMinF = Min(aMinF, aProj);
0053 aMaxF = Max(aMaxF, aProj);
0054
0055 if (aMinF <= aMaxB && aMaxF >= aMinB)
0056 {
0057 if (theInside == NULL || !(*theInside)) // only overlap test
0058 {
0059 return Standard_False;
0060 }
0061 }
0062 }
0063
0064 if (aMinF > aMaxB || aMaxF < aMinB)
0065 {
0066 return Standard_True; // fully separated
0067 }
0068 else if (theInside != NULL) // to check for inclusion?
0069 {
0070 *theInside &= aMinB >= aMinF && aMaxB <= aMaxF;
0071 }
0072
0073 return Standard_False;
0074 }
0075
0076 // =======================================================================
0077 // function : isSeparated
0078 // purpose : Checks if triangle and frustum are separated along the
0079 // given axis
0080 // =======================================================================
0081 template <int N>
0082 Standard_Boolean SelectMgr_Frustum<N>::isSeparated(const gp_Pnt& thePnt1,
0083 const gp_Pnt& thePnt2,
0084 const gp_Pnt& thePnt3,
0085 const gp_XYZ& theAxis) const
0086 {
0087 // frustum projection
0088 Standard_Real aMinF = RealLast();
0089 Standard_Real aMaxF = RealFirst();
0090
0091 // triangle projection
0092 Standard_Real aMinTr = RealLast();
0093 Standard_Real aMaxTr = RealFirst();
0094
0095 Standard_Real aTriangleProj;
0096
0097 aTriangleProj = theAxis.Dot(thePnt1.XYZ());
0098 aMinTr = Min(aMinTr, aTriangleProj);
0099 aMaxTr = Max(aMaxTr, aTriangleProj);
0100
0101 aTriangleProj = theAxis.Dot(thePnt2.XYZ());
0102 aMinTr = Min(aMinTr, aTriangleProj);
0103 aMaxTr = Max(aMaxTr, aTriangleProj);
0104
0105 aTriangleProj = theAxis.Dot(thePnt3.XYZ());
0106 aMinTr = Min(aMinTr, aTriangleProj);
0107 aMaxTr = Max(aMaxTr, aTriangleProj);
0108
0109 for (Standard_Integer aVertIter = 0; aVertIter < N * 2; ++aVertIter)
0110 {
0111 const Standard_Real aProj = myVertices[aVertIter].XYZ().Dot(theAxis);
0112
0113 aMinF = Min(aMinF, aProj);
0114 aMaxF = Max(aMaxF, aProj);
0115
0116 if (aMinF <= aMaxTr && aMaxF >= aMinTr)
0117 {
0118 return Standard_False;
0119 }
0120 }
0121
0122 return aMinF > aMaxTr || aMaxF < aMinTr;
0123 }
0124
0125 // =======================================================================
0126 // function : hasBoxOverlap
0127 // purpose : Returns true if selecting volume is overlapped by
0128 // axis-aligned bounding box with minimum corner at point
0129 // theMinPnt and maximum at point theMaxPnt
0130 // =======================================================================
0131 template <int N>
0132 Standard_Boolean SelectMgr_Frustum<N>::hasBoxOverlap(const SelectMgr_Vec3& theMinPnt,
0133 const SelectMgr_Vec3& theMaxPnt,
0134 Standard_Boolean* theInside) const
0135 {
0136 for (Standard_Integer anAxis = 0; anAxis < 3; ++anAxis)
0137 {
0138 if (theMinPnt[anAxis] > myMaxOrthoVertsProjections[anAxis]
0139 || theMaxPnt[anAxis] < myMinOrthoVertsProjections[anAxis])
0140 {
0141 return Standard_False; // fully separated
0142 }
0143 else if (theInside != NULL) // to check for inclusion?
0144 {
0145 *theInside &= theMinPnt[anAxis] >= myMinOrthoVertsProjections[anAxis]
0146 && theMaxPnt[anAxis] <= myMaxOrthoVertsProjections[anAxis];
0147 }
0148 }
0149
0150 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0151 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0152 {
0153 const gp_XYZ& aPlane = myPlanes[aPlaneIdx].XYZ();
0154
0155 const Standard_Real aBoxProjMin =
0156 aPlane.X() * (aPlane.X() < 0.f ? theMaxPnt.x() : theMinPnt.x())
0157 + aPlane.Y() * (aPlane.Y() < 0.f ? theMaxPnt.y() : theMinPnt.y())
0158 + aPlane.Z() * (aPlane.Z() < 0.f ? theMaxPnt.z() : theMinPnt.z());
0159
0160 const Standard_Real aBoxProjMax =
0161 aPlane.X() * (aPlane.X() < 0.f ? theMinPnt.x() : theMaxPnt.x())
0162 + aPlane.Y() * (aPlane.Y() < 0.f ? theMinPnt.y() : theMaxPnt.y())
0163 + aPlane.Z() * (aPlane.Z() < 0.f ? theMinPnt.z() : theMaxPnt.z());
0164
0165 Standard_ASSERT_RAISE(aBoxProjMax >= aBoxProjMin, "Error! Failed to project box");
0166
0167 if (aBoxProjMin > myMaxVertsProjections[aPlaneIdx]
0168 || aBoxProjMax < myMinVertsProjections[aPlaneIdx])
0169 {
0170 return Standard_False; // fully separated
0171 }
0172 else if (theInside != NULL) // to check for inclusion?
0173 {
0174 *theInside &= aBoxProjMin >= myMinVertsProjections[aPlaneIdx]
0175 && aBoxProjMax <= myMaxVertsProjections[aPlaneIdx];
0176 }
0177 }
0178
0179 for (Standard_Integer aDim = 0; aDim < 3; ++aDim)
0180 {
0181 // the following code performs a speedup of cross-product
0182 // of vector with 1.0 at the position aDim and myEdgeDirs[aVolDir]
0183 const Standard_Integer aNext = (aDim + 1) % 3;
0184 const Standard_Integer aNextNext = (aDim + 2) % 3;
0185 for (Standard_Integer aVolDir = 0, aDirectionsNb = Camera()->IsOrthographic() ? 4 : 6;
0186 aVolDir < aDirectionsNb;
0187 ++aVolDir)
0188 {
0189 gp_XYZ aDirection(DBL_MAX, DBL_MAX, DBL_MAX);
0190 aDirection.ChangeData()[aDim] = 0;
0191 aDirection.ChangeData()[aNext] = -myEdgeDirs[aVolDir].XYZ().GetData()[aNextNext];
0192 aDirection.ChangeData()[aNextNext] = myEdgeDirs[aVolDir].XYZ().GetData()[aNext];
0193
0194 if (isSeparated(theMinPnt, theMaxPnt, aDirection, theInside))
0195 {
0196 return Standard_False;
0197 }
0198 }
0199 }
0200
0201 return Standard_True;
0202 }
0203
0204 // =======================================================================
0205 // function : hasPointOverlap
0206 // purpose : SAT intersection test between defined volume and given point
0207 // =======================================================================
0208 template <int N>
0209 Standard_Boolean SelectMgr_Frustum<N>::hasPointOverlap(const gp_Pnt& thePnt) const
0210 {
0211 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0212 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0213 {
0214 const Standard_Real aPointProj = myPlanes[aPlaneIdx].XYZ().Dot(thePnt.XYZ());
0215
0216 if (aPointProj > myMaxVertsProjections[aPlaneIdx]
0217 || aPointProj < myMinVertsProjections[aPlaneIdx])
0218 {
0219 return Standard_False;
0220 }
0221 }
0222
0223 return Standard_True;
0224 }
0225
0226 // =======================================================================
0227 // function : hasSegmentOverlap
0228 // purpose : SAT intersection test between defined volume and given segment
0229 // =======================================================================
0230 template <int N>
0231 Standard_Boolean SelectMgr_Frustum<N>::hasSegmentOverlap(const gp_Pnt& theStartPnt,
0232 const gp_Pnt& theEndPnt) const
0233 {
0234 const gp_XYZ& aDir = theEndPnt.XYZ() - theStartPnt.XYZ();
0235 if (aDir.Modulus() < Precision::Confusion())
0236 return Standard_True;
0237
0238 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0239 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0240 {
0241 Standard_Real aMinSegm = RealLast(), aMaxSegm = RealFirst();
0242 Standard_Real aMinF = RealLast(), aMaxF = RealFirst();
0243
0244 Standard_Real aProj1 = myPlanes[aPlaneIdx].XYZ().Dot(theStartPnt.XYZ());
0245 Standard_Real aProj2 = myPlanes[aPlaneIdx].XYZ().Dot(theEndPnt.XYZ());
0246 aMinSegm = Min(aProj1, aProj2);
0247 aMaxSegm = Max(aProj1, aProj2);
0248
0249 aMaxF = myMaxVertsProjections[aPlaneIdx];
0250 aMinF = myMinVertsProjections[aPlaneIdx];
0251
0252 if (aMinSegm > aMaxF || aMaxSegm < aMinF)
0253 {
0254 return Standard_False;
0255 }
0256 }
0257
0258 Standard_Real aMin1 = DBL_MAX, aMax1 = -DBL_MAX;
0259 Standard_Real aMin2 = DBL_MAX, aMax2 = -DBL_MAX;
0260 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0261 {
0262 Standard_Real aProjection = aDir.Dot(myVertices[aVertIdx].XYZ());
0263 aMax2 = Max(aMax2, aProjection);
0264 aMin2 = Min(aMin2, aProjection);
0265 }
0266 Standard_Real aProj1 = aDir.Dot(theStartPnt.XYZ());
0267 Standard_Real aProj2 = aDir.Dot(theEndPnt.XYZ());
0268 aMin1 = Min(aProj1, aProj2);
0269 aMax1 = Max(aProj1, aProj2);
0270 if (aMin1 > aMax2 || aMax1 < aMin2)
0271 {
0272 return Standard_False;
0273 }
0274
0275 Standard_Integer aDirectionsNb = Camera()->IsOrthographic() ? 4 : 6;
0276 for (Standard_Integer aEdgeDirIdx = 0; aEdgeDirIdx < aDirectionsNb; ++aEdgeDirIdx)
0277 {
0278 Standard_Real aMinSegm = DBL_MAX, aMaxSegm = -DBL_MAX;
0279 Standard_Real aMinF = DBL_MAX, aMaxF = -DBL_MAX;
0280
0281 const gp_XYZ aTestDir = aDir.Crossed(myEdgeDirs[aEdgeDirIdx].XYZ());
0282
0283 Standard_Real Proj1 = aTestDir.Dot(theStartPnt.XYZ());
0284 Standard_Real Proj2 = aTestDir.Dot(theEndPnt.XYZ());
0285 aMinSegm = Min(Proj1, Proj2);
0286 aMaxSegm = Max(Proj1, Proj2);
0287
0288 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0289 {
0290 Standard_Real aProjection = aTestDir.Dot(myVertices[aVertIdx].XYZ());
0291 aMaxF = Max(aMaxF, aProjection);
0292 aMinF = Min(aMinF, aProjection);
0293 }
0294
0295 if (aMinSegm > aMaxF || aMaxSegm < aMinF)
0296 {
0297 return Standard_False;
0298 }
0299 }
0300
0301 return Standard_True;
0302 }
0303
0304 // =======================================================================
0305 // function : hasPolygonOverlap
0306 // purpose : SAT intersection test between frustum given and planar convex
0307 // polygon represented as ordered point set
0308 // =======================================================================
0309 template <int N>
0310 Standard_Boolean SelectMgr_Frustum<N>::hasPolygonOverlap(const TColgp_Array1OfPnt& theArrayOfPnts,
0311 gp_Vec& theNormal) const
0312 {
0313 Standard_Integer aStartIdx = theArrayOfPnts.Lower();
0314 Standard_Integer anEndIdx = theArrayOfPnts.Upper();
0315
0316 const gp_XYZ& aPnt1 = theArrayOfPnts.Value(aStartIdx).XYZ();
0317 const gp_XYZ& aPnt2 = theArrayOfPnts.Value(aStartIdx + 1).XYZ();
0318 const gp_XYZ& aPnt3 = theArrayOfPnts.Value(aStartIdx + 2).XYZ();
0319 const gp_XYZ aVec1 = aPnt1 - aPnt2;
0320 const gp_XYZ aVec2 = aPnt3 - aPnt2;
0321 theNormal = aVec2.Crossed(aVec1);
0322 const gp_XYZ& aNormal = theNormal.XYZ();
0323 Standard_Real aPolygProjection = aNormal.Dot(aPnt1);
0324
0325 Standard_Real aMax = RealFirst();
0326 Standard_Real aMin = RealLast();
0327 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0328 {
0329 Standard_Real aProjection = aNormal.Dot(myVertices[aVertIdx].XYZ());
0330 aMax = Max(aMax, aProjection);
0331 aMin = Min(aMin, aProjection);
0332 }
0333 if (aPolygProjection > aMax || aPolygProjection < aMin)
0334 {
0335 return Standard_False;
0336 }
0337
0338 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0339 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0340 {
0341 Standard_Real aMaxF = RealFirst();
0342 Standard_Real aMinF = RealLast();
0343 Standard_Real aMaxPolyg = RealFirst();
0344 Standard_Real aMinPolyg = RealLast();
0345 const gp_XYZ& aPlane = myPlanes[aPlaneIdx].XYZ();
0346 for (Standard_Integer aPntIter = aStartIdx; aPntIter <= anEndIdx; ++aPntIter)
0347 {
0348 Standard_Real aProjection = aPlane.Dot(theArrayOfPnts.Value(aPntIter).XYZ());
0349 aMaxPolyg = Max(aMaxPolyg, aProjection);
0350 aMinPolyg = Min(aMinPolyg, aProjection);
0351 }
0352 aMaxF = myMaxVertsProjections[aPlaneIdx];
0353 aMinF = myMinVertsProjections[aPlaneIdx];
0354 if (aMinPolyg > aMaxF || aMaxPolyg < aMinF)
0355 {
0356 return Standard_False;
0357 }
0358 }
0359
0360 Standard_Integer aDirectionsNb = Camera()->IsOrthographic() ? 4 : 6;
0361 for (Standard_Integer aPntsIter = 0,
0362 aLastIdx = anEndIdx - aStartIdx,
0363 aLen = theArrayOfPnts.Length();
0364 aPntsIter <= aLastIdx;
0365 ++aPntsIter)
0366 {
0367 const gp_XYZ aSegmDir = theArrayOfPnts.Value((aPntsIter + 1) % aLen + aStartIdx).XYZ()
0368 - theArrayOfPnts.Value(aPntsIter + aStartIdx).XYZ();
0369 for (Standard_Integer aVolDir = 0; aVolDir < aDirectionsNb; ++aVolDir)
0370 {
0371 Standard_Real aMaxPolyg = RealFirst();
0372 Standard_Real aMinPolyg = RealLast();
0373 Standard_Real aMaxF = RealFirst();
0374 Standard_Real aMinF = RealLast();
0375 const gp_XYZ aTestDir = aSegmDir.Crossed(myEdgeDirs[aVolDir].XYZ());
0376
0377 for (Standard_Integer aPntIter = aStartIdx; aPntIter <= anEndIdx; ++aPntIter)
0378 {
0379 Standard_Real aProjection = aTestDir.Dot(theArrayOfPnts.Value(aPntIter).XYZ());
0380 aMaxPolyg = Max(aMaxPolyg, aProjection);
0381 aMinPolyg = Min(aMinPolyg, aProjection);
0382 }
0383
0384 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0385 {
0386 Standard_Real aProjection = aTestDir.Dot(myVertices[aVertIdx].XYZ());
0387 aMaxF = Max(aMaxF, aProjection);
0388 aMinF = Min(aMinF, aProjection);
0389 }
0390
0391 if (aMinPolyg > aMaxF || aMaxPolyg < aMinF)
0392 {
0393 return Standard_False;
0394 }
0395 }
0396 }
0397
0398 return Standard_True;
0399 }
0400
0401 // =======================================================================
0402 // function : hasTriangleOverlap
0403 // purpose : SAT intersection test between defined volume and given triangle
0404 // =======================================================================
0405 template <int N>
0406 Standard_Boolean SelectMgr_Frustum<N>::hasTriangleOverlap(const gp_Pnt& thePnt1,
0407 const gp_Pnt& thePnt2,
0408 const gp_Pnt& thePnt3,
0409 gp_Vec& theNormal) const
0410 {
0411 const gp_XYZ aTrEdges[3] = {thePnt2.XYZ() - thePnt1.XYZ(),
0412 thePnt3.XYZ() - thePnt2.XYZ(),
0413 thePnt1.XYZ() - thePnt3.XYZ()};
0414
0415 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0416 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0417 {
0418 const gp_XYZ& aPlane = myPlanes[aPlaneIdx].XYZ();
0419 Standard_Real aTriangleProj;
0420
0421 aTriangleProj = aPlane.Dot(thePnt1.XYZ());
0422 Standard_Real aTriangleProjMin = aTriangleProj;
0423 Standard_Real aTriangleProjMax = aTriangleProj;
0424
0425 aTriangleProj = aPlane.Dot(thePnt2.XYZ());
0426 aTriangleProjMin = Min(aTriangleProjMin, aTriangleProj);
0427 aTriangleProjMax = Max(aTriangleProjMax, aTriangleProj);
0428
0429 aTriangleProj = aPlane.Dot(thePnt3.XYZ());
0430 aTriangleProjMin = Min(aTriangleProjMin, aTriangleProj);
0431 aTriangleProjMax = Max(aTriangleProjMax, aTriangleProj);
0432
0433 Standard_Real aFrustumProjMax = myMaxVertsProjections[aPlaneIdx];
0434 Standard_Real aFrustumProjMin = myMinVertsProjections[aPlaneIdx];
0435 if (aTriangleProjMin > aFrustumProjMax || aTriangleProjMax < aFrustumProjMin)
0436 {
0437 return Standard_False;
0438 }
0439 }
0440
0441 theNormal = aTrEdges[2].Crossed(aTrEdges[0]);
0442 if (isSeparated(thePnt1, thePnt2, thePnt3, theNormal.XYZ()))
0443 {
0444 return Standard_False;
0445 }
0446
0447 Standard_Integer aDirectionsNb = myCamera->IsOrthographic() ? 4 : 6;
0448 for (Standard_Integer aTriangleEdgeIdx = 0; aTriangleEdgeIdx < 3; ++aTriangleEdgeIdx)
0449 {
0450 for (Standard_Integer aVolDir = 0; aVolDir < aDirectionsNb; ++aVolDir)
0451 {
0452 const gp_XYZ& aTestDirection = myEdgeDirs[aVolDir].XYZ().Crossed(aTrEdges[aTriangleEdgeIdx]);
0453
0454 if (isSeparated(thePnt1, thePnt2, thePnt3, aTestDirection))
0455 {
0456 return Standard_False;
0457 }
0458 }
0459 }
0460 return Standard_True;
0461 }
0462
0463 // =======================================================================
0464 // function : hasSphereOverlap
0465 // purpose :
0466 // =======================================================================
0467 template <int N>
0468 Standard_Boolean SelectMgr_Frustum<N>::hasSphereOverlap(const gp_Pnt& thePnt,
0469 const Standard_Real theRadius,
0470 Standard_Boolean* theInside) const
0471 {
0472 Standard_Boolean isOverlapFull = Standard_True;
0473 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0474 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N; aPlaneIdx += anIncFactor)
0475 {
0476 const gp_XYZ& aPlane = myPlanes[aPlaneIdx].XYZ();
0477 const Standard_Real aNormVecLen = Sqrt(aPlane.Dot(aPlane));
0478 const Standard_Real aCenterProj = aPlane.Dot(thePnt.XYZ()) / aNormVecLen;
0479 const Standard_Real aMaxDist = myMaxVertsProjections[aPlaneIdx] / aNormVecLen;
0480 const Standard_Real aMinDist = myMinVertsProjections[aPlaneIdx] / aNormVecLen;
0481 if (aCenterProj > (aMaxDist + theRadius) || aCenterProj < (aMinDist - theRadius))
0482 {
0483 return Standard_False; // fully separated
0484 }
0485 else if (theInside)
0486 {
0487 *theInside &= aCenterProj >= (aMinDist + theRadius) && aCenterProj <= (aMaxDist - theRadius);
0488 }
0489 isOverlapFull &= aCenterProj >= (aMinDist + theRadius) && aCenterProj <= (aMaxDist - theRadius);
0490 }
0491 if (theInside || isOverlapFull)
0492 {
0493 return Standard_True;
0494 }
0495 const gp_Vec aVecPlane1(myVertices[0], myVertices[2]);
0496 const gp_Vec aVecPlane2(myVertices[0], myVertices[2 * N - 2]);
0497 if (aVecPlane1.IsParallel(aVecPlane2, Precision::Angular()))
0498 {
0499 return Standard_False;
0500 }
0501 const gp_Dir aNorm(aVecPlane1.Crossed(aVecPlane2));
0502 gp_Pnt aBoundariesCArr[5];
0503 NCollection_Array1<gp_Pnt> aBoundaries(aBoundariesCArr[0], 0, N);
0504 for (Standard_Integer anIdx = 0; anIdx < N * 2; anIdx += 2)
0505 {
0506 aBoundaries.SetValue(anIdx / 2, myVertices[anIdx]);
0507 }
0508 // distance from point(x,y,z) to plane(A,B,C,D) d = | Ax + By + Cz + D | / sqrt (A^2 + B^2 + C^2)
0509 // = aPnt.Dot (Norm) / 1
0510 const gp_Pnt aCenterProj = thePnt.XYZ() - aNorm.XYZ() * thePnt.XYZ().Dot(aNorm.XYZ());
0511 Standard_Boolean isBoundaryInside = Standard_False;
0512 return IsBoundaryIntersectSphere(aCenterProj, theRadius, aNorm, aBoundaries, isBoundaryInside);
0513 }
0514
0515 // =======================================================================
0516 // function : IsDotInside
0517 // purpose :
0518 // =======================================================================
0519 template <int N>
0520 Standard_Boolean SelectMgr_Frustum<N>::isDotInside(const gp_Pnt& thePnt,
0521 const TColgp_Array1OfPnt& theVertices) const
0522 {
0523 Standard_Real anAngle = 0.0;
0524 for (Standard_Integer aVertIdx = 0; aVertIdx < theVertices.Size(); aVertIdx++)
0525 {
0526 const gp_Pnt aVert1 = theVertices[aVertIdx];
0527 const gp_Pnt aVert2 =
0528 (aVertIdx == (theVertices.Size() - 1) ? theVertices[0] : theVertices[aVertIdx + 1]);
0529 const gp_Vec aVec1(thePnt, aVert1);
0530 const gp_Vec aVec2(thePnt, aVert2);
0531 anAngle += aVec1.Angle(aVec2);
0532 }
0533 if (Abs(anAngle - 2.0 * M_PI) < Precision::Angular())
0534 {
0535 return true;
0536 }
0537 return false;
0538 }
0539
0540 // =======================================================================
0541 // function : isSegmentsIntersect
0542 // purpose :
0543 // =======================================================================
0544 template <int N>
0545 Standard_Boolean SelectMgr_Frustum<N>::isSegmentsIntersect(const gp_Pnt& thePnt1Seg1,
0546 const gp_Pnt& thePnt2Seg1,
0547 const gp_Pnt& thePnt1Seg2,
0548 const gp_Pnt& thePnt2Seg2) const
0549 {
0550 const gp_Mat aMatPln(thePnt2Seg1.X() - thePnt1Seg1.X(),
0551 thePnt2Seg1.Y() - thePnt1Seg1.Y(),
0552 thePnt2Seg1.Z() - thePnt1Seg1.Z(),
0553 thePnt1Seg2.X() - thePnt1Seg1.X(),
0554 thePnt1Seg2.Y() - thePnt1Seg1.Y(),
0555 thePnt1Seg2.Z() - thePnt1Seg1.Z(),
0556 thePnt2Seg2.X() - thePnt1Seg1.X(),
0557 thePnt2Seg2.Y() - thePnt1Seg1.Y(),
0558 thePnt2Seg2.Z() - thePnt1Seg1.Z());
0559 if (Abs(aMatPln.Determinant()) > Precision::Confusion())
0560 {
0561 return false;
0562 }
0563
0564 Standard_Real aFst[4] = {thePnt1Seg1.X(), thePnt2Seg1.X(), thePnt1Seg2.X(), thePnt2Seg2.X()};
0565 Standard_Real aSnd[4] = {thePnt1Seg1.Y(), thePnt2Seg1.Y(), thePnt1Seg2.Y(), thePnt2Seg2.Y()};
0566 if (aFst[0] == aFst[2] && aFst[1] == aFst[3])
0567 {
0568 aFst[0] = thePnt1Seg1.Z();
0569 aFst[1] = thePnt2Seg1.Z();
0570 aFst[2] = thePnt1Seg2.Z();
0571 aFst[3] = thePnt2Seg2.Z();
0572 }
0573 if (aSnd[0] == aSnd[2] && aSnd[1] == aSnd[3])
0574 {
0575 aSnd[0] = thePnt1Seg1.Z();
0576 aSnd[1] = thePnt2Seg1.Z();
0577 aSnd[2] = thePnt1Seg2.Z();
0578 aSnd[3] = thePnt2Seg2.Z();
0579 }
0580 const gp_Mat2d aMat(gp_XY(aFst[0] - aFst[1], aSnd[0] - aSnd[1]),
0581 gp_XY(aFst[3] - aFst[2], aSnd[3] - aSnd[2]));
0582
0583 const gp_Mat2d aMatU(gp_XY(aFst[0] - aFst[2], aSnd[0] - aSnd[2]),
0584 gp_XY(aFst[3] - aFst[2], aSnd[3] - aSnd[2]));
0585
0586 const gp_Mat2d aMatV(gp_XY(aFst[0] - aFst[1], aSnd[0] - aSnd[1]),
0587 gp_XY(aFst[0] - aFst[2], aSnd[0] - aSnd[2]));
0588 if (aMat.Determinant() == 0.0)
0589 {
0590 return false;
0591 }
0592
0593 const Standard_Real anU = aMatU.Determinant() / aMat.Determinant();
0594 const Standard_Real aV = aMatV.Determinant() / aMat.Determinant();
0595 if (anU >= 0.0 && anU <= 1.0 && aV >= 0.0 && aV <= 1.0)
0596 {
0597 return true;
0598 }
0599 return false;
0600 }
0601
0602 // =======================================================================
0603 // function : isIntersectCircle
0604 // purpose :
0605 // =======================================================================
0606 template <int N>
0607 Standard_Boolean SelectMgr_Frustum<N>::isIntersectCircle(
0608 const Standard_Real theRadius,
0609 const gp_Pnt& theCenter,
0610 const gp_Trsf& theTrsf,
0611 const TColgp_Array1OfPnt& theVertices) const
0612 {
0613 const gp_Trsf aTrsfInv = theTrsf.Inverted();
0614 const gp_Dir aRayDir = gp_Dir(myEdgeDirs[N == 4 ? 4 : 0]).Transformed(aTrsfInv);
0615 if (aRayDir.Z() == 0.0)
0616 {
0617 return false;
0618 }
0619
0620 for (Standard_Integer anIdx = theVertices.Lower(); anIdx <= theVertices.Upper(); anIdx++)
0621 {
0622 const gp_Pnt aPntStart = theVertices.Value(anIdx).Transformed(aTrsfInv);
0623 const gp_Pnt aPntFinish = anIdx == theVertices.Upper()
0624 ? theVertices.Value(theVertices.Lower()).Transformed(aTrsfInv)
0625 : theVertices.Value(anIdx + 1).Transformed(aTrsfInv);
0626
0627 // Project points on the end face plane
0628 const Standard_Real aParam1 = (theCenter.Z() - aPntStart.Z()) / aRayDir.Z();
0629 const Standard_Real aX1 = aPntStart.X() + aRayDir.X() * aParam1;
0630 const Standard_Real anY1 = aPntStart.Y() + aRayDir.Y() * aParam1;
0631
0632 const Standard_Real aParam2 = (theCenter.Z() - aPntFinish.Z()) / aRayDir.Z();
0633 const Standard_Real aX2 = aPntFinish.X() + aRayDir.X() * aParam2;
0634 const Standard_Real anY2 = aPntFinish.Y() + aRayDir.Y() * aParam2;
0635
0636 // Solving quadratic equation anA * T^2 + 2 * aK * T + aC = 0
0637 const Standard_Real anA = (aX1 - aX2) * (aX1 - aX2) + (anY1 - anY2) * (anY1 - anY2);
0638 const Standard_Real aK = aX1 * (aX2 - aX1) + anY1 * (anY2 - anY1);
0639 const Standard_Real aC = aX1 * aX1 + anY1 * anY1 - theRadius * theRadius;
0640
0641 const Standard_Real aDiscr = aK * aK - anA * aC;
0642 if (aDiscr >= 0.0)
0643 {
0644 const Standard_Real aT1 = (-aK + Sqrt(aDiscr)) / anA;
0645 const Standard_Real aT2 = (-aK - Sqrt(aDiscr)) / anA;
0646 if ((aT1 >= 0 && aT1 <= 1) || (aT2 >= 0 && aT2 <= 1))
0647 {
0648 return true;
0649 }
0650 }
0651 }
0652 return false;
0653 }
0654
0655 // =======================================================================
0656 // function : isInsideCylinderEndFace
0657 // purpose :
0658 // =======================================================================
0659 template <int N>
0660 Standard_Boolean SelectMgr_Frustum<N>::isInsideCylinderEndFace(
0661 const Standard_Real theBottomRad,
0662 const Standard_Real theTopRad,
0663 const Standard_Real theHeight,
0664 const gp_Trsf& theTrsf,
0665 const TColgp_Array1OfPnt& theVertices) const
0666 {
0667 const gp_Trsf aTrsfInv = theTrsf.Inverted();
0668 const gp_Dir aRayDir = gp_Dir(myEdgeDirs[N == 4 ? 4 : 0]).Transformed(aTrsfInv);
0669 if (aRayDir.Z() == 0.0)
0670 {
0671 return false;
0672 }
0673
0674 for (Standard_Integer anIdx = theVertices.Lower(); anIdx <= theVertices.Upper(); anIdx++)
0675 {
0676 const gp_Pnt aLoc = theVertices.Value(anIdx).Transformed(aTrsfInv);
0677
0678 const Standard_Real aTime1 = (0 - aLoc.Z()) / aRayDir.Z();
0679 const Standard_Real aX1 = aLoc.X() + aRayDir.X() * aTime1;
0680 const Standard_Real anY1 = aLoc.Y() + aRayDir.Y() * aTime1;
0681
0682 const Standard_Real aTime2 = (theHeight - aLoc.Z()) / aRayDir.Z();
0683 const Standard_Real aX2 = aLoc.X() + aRayDir.X() * aTime2;
0684 const Standard_Real anY2 = aLoc.Y() + aRayDir.Y() * aTime2;
0685
0686 if (aX1 * aX1 + anY1 * anY1 <= theBottomRad * theBottomRad
0687 && aX2 * aX2 + anY2 * anY2 <= theTopRad * theTopRad)
0688 {
0689 continue;
0690 }
0691
0692 return false;
0693 }
0694 return true;
0695 }
0696
0697 // =======================================================================
0698 // function : hasCylinderOverlap
0699 // purpose :
0700 // =======================================================================
0701 template <int N>
0702 Standard_Boolean SelectMgr_Frustum<N>::hasCylinderOverlap(const Standard_Real theBottomRad,
0703 const Standard_Real theTopRad,
0704 const Standard_Real theHeight,
0705 const gp_Trsf& theTrsf,
0706 const Standard_Boolean theIsHollow,
0707 Standard_Boolean* theInside) const
0708 {
0709 gp_Pnt aVerticesBuf[N];
0710 TColgp_Array1OfPnt aVertices(aVerticesBuf[0], 0, N - 1);
0711 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0712 if (anIncFactor == 2)
0713 {
0714 const Standard_Integer anIndices[] = {0, 2, 6, 4};
0715 for (Standard_Integer anIdx = 0; anIdx < N; anIdx++)
0716 {
0717 aVertices.SetValue(anIdx, myVertices[anIndices[anIdx]]);
0718 }
0719 }
0720 else
0721 {
0722 for (Standard_Integer anIdx = 0; anIdx < N; anIdx++)
0723 {
0724 aVertices.SetValue(anIdx, myVertices[anIdx]);
0725 }
0726 }
0727
0728 if (theIsHollow
0729 && isInsideCylinderEndFace(theBottomRad, theTopRad, theHeight, theTrsf, aVertices))
0730 {
0731 if (theInside != NULL)
0732 {
0733 *theInside = false;
0734 }
0735 return false;
0736 }
0737
0738 const gp_Dir aCylNorm(gp::DZ().Transformed(theTrsf));
0739 const gp_Pnt aBottomCenter(gp::Origin().Transformed(theTrsf));
0740 const gp_Pnt aTopCenter = aBottomCenter.XYZ() + aCylNorm.XYZ() * theHeight;
0741
0742 const gp_Dir aViewRayDir = gp_Dir(myEdgeDirs[N == 4 ? 4 : 0]);
0743 const gp_Pln aPln(myVertices[0], aViewRayDir);
0744 Standard_Real aCoefA, aCoefB, aCoefC, aCoefD;
0745 aPln.Coefficients(aCoefA, aCoefB, aCoefC, aCoefD);
0746
0747 const Standard_Real aTBottom = -(aBottomCenter.XYZ().Dot(aViewRayDir.XYZ()) + aCoefD);
0748 const gp_Pnt aBottomCenterProject(aCoefA * aTBottom + aBottomCenter.X(),
0749 aCoefB * aTBottom + aBottomCenter.Y(),
0750 aCoefC * aTBottom + aBottomCenter.Z());
0751 const Standard_Real aTTop = -(aTopCenter.XYZ().Dot(aViewRayDir.XYZ()) + aCoefD);
0752 const gp_Pnt aTopCenterProject(aCoefA * aTTop + aTopCenter.X(),
0753 aCoefB * aTTop + aTopCenter.Y(),
0754 aCoefC * aTTop + aTopCenter.Z());
0755 gp_Vec aCylNormProject(0, 0, 0);
0756 if (aTopCenterProject.Distance(aBottomCenterProject) > 0.0)
0757 {
0758 aCylNormProject = gp_Vec((aTopCenterProject.XYZ() - aBottomCenterProject.XYZ())
0759 / aTopCenterProject.Distance(aBottomCenterProject));
0760 }
0761
0762 gp_Pnt aPoints[6];
0763 const gp_Dir aDirEndFaces = (aCylNorm.IsParallel(aViewRayDir, Precision::Angular()))
0764 ? gp::DY().Transformed(theTrsf)
0765 : aCylNorm.Crossed(aViewRayDir);
0766
0767 const Standard_Real anAngle = aCylNorm.Angle(aViewRayDir);
0768 aPoints[0] =
0769 aBottomCenterProject.XYZ() - aCylNormProject.XYZ() * theBottomRad * Abs(Cos(anAngle));
0770 aPoints[1] = aBottomCenterProject.XYZ() + aDirEndFaces.XYZ() * theBottomRad;
0771 aPoints[2] = aTopCenterProject.XYZ() + aDirEndFaces.XYZ() * theTopRad;
0772 aPoints[3] = aTopCenterProject.XYZ() + aCylNormProject.XYZ() * theTopRad * Abs(Cos(anAngle));
0773 aPoints[4] = aTopCenterProject.XYZ() - aDirEndFaces.XYZ() * theTopRad;
0774 aPoints[5] = aBottomCenterProject.XYZ() - aDirEndFaces.XYZ() * theBottomRad;
0775 const TColgp_Array1OfPnt aPointsArr(aPoints[0], 0, 5);
0776
0777 for (Standard_Integer anIdx = 0; anIdx < N; anIdx++)
0778 {
0779 if ((aCylNormProject.Dot(aCylNormProject) == 0.0
0780 && aVertices.Value(anIdx).Distance(aPoints[0]) <= Max(theTopRad, theBottomRad))
0781 || isDotInside(aVertices.Value(anIdx), aPointsArr))
0782 {
0783 if (theInside != NULL)
0784 {
0785 *theInside = false;
0786 }
0787 return true;
0788 }
0789 }
0790
0791 for (Standard_Integer anIdx = aVertices.Lower(); anIdx <= aVertices.Upper(); anIdx++)
0792 {
0793 const gp_Pnt aPnt1Seg = aVertices[anIdx];
0794 const gp_Pnt aPnt2Seg =
0795 (anIdx == aVertices.Upper()) ? aVertices[aVertices.Lower()] : aVertices[anIdx + 1];
0796 if (isSegmentsIntersect(aPoints[1], aPoints[2], aPnt1Seg, aPnt2Seg)
0797 || isSegmentsIntersect(aPoints[4], aPoints[5], aPnt1Seg, aPnt2Seg)
0798 || isSegmentsIntersect(aPoints[4], aPoints[2], aPnt1Seg, aPnt2Seg)
0799 || isSegmentsIntersect(aPoints[1], aPoints[5], aPnt1Seg, aPnt2Seg))
0800 {
0801 if (theInside != NULL)
0802 {
0803 *theInside = false;
0804 }
0805 return true;
0806 }
0807 }
0808
0809 if (!theIsHollow
0810 && (isIntersectCircle(theBottomRad, gp_Pnt(0, 0, 0), theTrsf, aVertices)
0811 || isIntersectCircle(theTopRad, gp_Pnt(0, 0, theHeight), theTrsf, aVertices)))
0812 {
0813 if (theInside != NULL)
0814 {
0815 *theInside = false;
0816 }
0817 return true;
0818 }
0819 bool isCylInsideRec = true;
0820 for (int i = 0; i < 6; ++i)
0821 {
0822 isCylInsideRec &= isDotInside(aPoints[i], aVertices);
0823 }
0824 if (theInside != NULL)
0825 {
0826 *theInside &= isCylInsideRec;
0827 }
0828 return isCylInsideRec;
0829 }
0830
0831 // =======================================================================
0832 // function : hasCircleOverlap
0833 // purpose :
0834 // =======================================================================
0835 template <int N>
0836 Standard_Boolean SelectMgr_Frustum<N>::hasCircleOverlap(const Standard_Real theRadius,
0837 const gp_Trsf& theTrsf,
0838 const Standard_Boolean theIsFilled,
0839 Standard_Boolean* theInside) const
0840 {
0841 gp_Pnt aVerticesBuf[N];
0842 TColgp_Array1OfPnt aVertices(aVerticesBuf[0], 0, N - 1);
0843 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0844 if (anIncFactor == 2)
0845 {
0846 const Standard_Integer anIndices[] = {0, 2, 6, 4};
0847 for (Standard_Integer anIdx = 0; anIdx < N; anIdx++)
0848 {
0849 aVertices.SetValue(anIdx, myVertices[anIndices[anIdx]]);
0850 }
0851 }
0852 else
0853 {
0854 for (Standard_Integer anIdx = 0; anIdx < N; anIdx++)
0855 {
0856 aVertices.SetValue(anIdx, myVertices[anIdx]);
0857 }
0858 }
0859
0860 if (isIntersectCircle(theRadius, gp_Pnt(0, 0, 0), theTrsf, aVertices))
0861 {
0862 if (theInside != NULL)
0863 {
0864 *theInside = false;
0865 }
0866 return true;
0867 }
0868
0869 gp_Pnt aCircCenter = gp::Origin(); //.Transformed (theTrsf);
0870 const gp_Dir aViewRayDir = gp_Dir(myEdgeDirs[N == 4 ? 4 : 0]);
0871 const gp_Pln aPln(myVertices[0], aViewRayDir);
0872 Standard_Real aCoefA, aCoefB, aCoefC, aCoefD;
0873 aPln.Coefficients(aCoefA, aCoefB, aCoefC, aCoefD);
0874
0875 const Standard_Real aTCenter = -(aCircCenter.XYZ().Dot(aViewRayDir.XYZ()) + aCoefD);
0876 const gp_Pnt aCenterProject(aCoefA * aTCenter, aCoefB * aTCenter, aCoefC * aTCenter);
0877
0878 const Standard_Boolean isCenterInside = isDotInside(aCenterProject, aVertices);
0879
0880 Standard_Boolean isInside = false;
0881 for (Standard_Integer anIdx = aVertices.Lower(); anIdx <= aVertices.Upper(); anIdx++)
0882 {
0883 if (aVertices.Value(anIdx).Distance(aCenterProject) > theRadius)
0884 {
0885 isInside = true;
0886 break;
0887 }
0888 }
0889
0890 if (theInside != NULL)
0891 {
0892 *theInside = isInside && isCenterInside;
0893 }
0894
0895 return theIsFilled ? !isInside || (isCenterInside && isInside) : isInside && isCenterInside;
0896 }
0897
0898 //=======================================================================
0899 // function : DumpJson
0900 // purpose :
0901 //=======================================================================
0902 template <int N>
0903 void SelectMgr_Frustum<N>::DumpJson(Standard_OStream& theOStream, Standard_Integer theDepth) const
0904 {
0905 OCCT_DUMP_TRANSIENT_CLASS_BEGIN(theOStream)
0906
0907 const Standard_Integer anIncFactor = (Camera()->IsOrthographic() && N == 4) ? 2 : 1;
0908 for (Standard_Integer aPlaneIdx = 0; aPlaneIdx < N + 1; aPlaneIdx += anIncFactor)
0909 {
0910 const gp_Vec& aPlane = myPlanes[aPlaneIdx];
0911 OCCT_DUMP_FIELD_VALUES_DUMPED(theOStream, theDepth, &aPlane)
0912
0913 OCCT_DUMP_FIELD_VALUE_NUMERICAL(theOStream, myMaxVertsProjections[aPlaneIdx])
0914 OCCT_DUMP_FIELD_VALUE_NUMERICAL(theOStream, myMinVertsProjections[aPlaneIdx])
0915 }
0916
0917 for (Standard_Integer aVertIdx = 0; aVertIdx < N * 2; ++aVertIdx)
0918 {
0919 const gp_Pnt& aVertex = myVertices[aVertIdx];
0920 OCCT_DUMP_FIELD_VALUES_DUMPED(theOStream, theDepth, &aVertex)
0921 }
0922
0923 OCCT_DUMP_FIELD_VALUE_NUMERICAL(theOStream, myPixelTolerance)
0924 OCCT_DUMP_FIELD_VALUE_POINTER(theOStream, myBuilder)
0925 OCCT_DUMP_FIELD_VALUE_POINTER(theOStream, myCamera)
0926
0927 for (Standard_Integer anIndex = 0; anIndex < 3; anIndex++)
0928 {
0929 Standard_Real aMaxOrthoVertsProjections = myMaxOrthoVertsProjections[anIndex];
0930 Standard_Real aMinOrthoVertsProjections = myMinOrthoVertsProjections[anIndex];
0931
0932 OCCT_DUMP_FIELD_VALUE_NUMERICAL(theOStream, aMaxOrthoVertsProjections)
0933 OCCT_DUMP_FIELD_VALUE_NUMERICAL(theOStream, aMinOrthoVertsProjections)
0934 }
0935
0936 for (Standard_Integer anIndex = 0; anIndex < 6; anIndex++)
0937 {
0938 const gp_Vec& anEdgeDir = myEdgeDirs[anIndex];
0939 OCCT_DUMP_FIELD_VALUES_DUMPED(theOStream, theDepth, &anEdgeDir)
0940 }
0941 }
