0024826: Wrapping of parallelisation algorithms
[occt.git] / src / OpenGl / OpenGl_BVHTreeSelector.cxx
CommitLineData
b7cd4ba7 1// Created on: 2013-12-25
2// Created by: Varvara POSKONINA
3// Copyright (c) 1999-2014 OPEN CASCADE SAS
4//
5// This file is part of Open CASCADE Technology software library.
6//
7// This library is free software; you can redistribute it and/or modify it under
8// the terms of the GNU Lesser General Public License version 2.1 as published
9// by the Free Software Foundation, with special exception defined in the file
10// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11// distribution for complete text of the license and disclaimer of any warranty.
12//
13// Alternatively, this file may be used under the terms of Open CASCADE
14// commercial license or contractual agreement.
15
16#include <OpenGl_BVHTreeSelector.hxx>
17#include <OpenGl_BVHClipPrimitiveSet.hxx>
18
19#include <vector>
20#include <algorithm>
21
22// =======================================================================
23// function : DotProduct
24// purpose : Calculates a dot product of 4-dimensional vectors in homogeneous coordinates
25// =======================================================================
26static Standard_ShortReal DotProduct (const OpenGl_Vec4& theA,
27 const OpenGl_Vec4& theB)
28{
29 return theA.x() * theB.x() + theA.y() * theB.y() + theA.z() * theB.z();
30}
31
32// =======================================================================
33// function : BinarySign
34// purpose :
35// =======================================================================
36static OpenGl_Vec4 BinarySign (const OpenGl_Vec4& theVec)
37{
38 return OpenGl_Vec4 (theVec.x() > 0.0f ? 1.0f : 0.0f,
39 theVec.y() > 0.0f ? 1.0f : 0.0f,
40 theVec.z() > 0.0f ? 1.0f : 0.0f,
41 theVec.w() > 0.0f ? 1.0f : 0.0f);
42}
43
44// =======================================================================
45// function : InversedBinarySign
46// purpose :
47// =======================================================================
48static OpenGl_Vec4 InversedBinarySign (const OpenGl_Vec4& theVec)
49{
50 return OpenGl_Vec4 (theVec.x() > 0.0f ? 0.0f : 1.0f,
51 theVec.y() > 0.0f ? 0.0f : 1.0f,
52 theVec.z() > 0.0f ? 0.0f : 1.0f,
53 theVec.w() > 0.0f ? 0.0f : 1.0f);
54}
55
56// =======================================================================
57// function : OpenGl_BVHTreeSelector
58// purpose :
59// =======================================================================
60OpenGl_BVHTreeSelector::OpenGl_BVHTreeSelector()
61: myIsProjectionParallel (Standard_True),
62 myProjectionState (0),
63 myModelViewState (0)
64{
65 //
66}
67
68// =======================================================================
69// function : SetClipVolume
70// purpose : Retrieves view volume's planes equations and its vertices from projection and modelview matrices.
71// =======================================================================
72void OpenGl_BVHTreeSelector::SetViewVolume (const Handle(Graphic3d_Camera)& theCamera)
73{
74 myIsProjectionParallel = theCamera->IsOrthographic();
75 const OpenGl_Mat4& aProjMat = theCamera->ProjectionMatrixF();
76 const OpenGl_Mat4& aModelMat = theCamera->OrientationMatrixF();
77
78 Standard_ShortReal nLeft = 0.0f, nRight = 0.0f, nTop = 0.0f, nBottom = 0.0f;
79 Standard_ShortReal fLeft = 0.0f, fRight = 0.0f, fTop = 0.0f, fBottom = 0.0f;
80 Standard_ShortReal aNear = 0.0f, aFar = 0.0f;
81 if (!myIsProjectionParallel)
82 {
83 // handle perspective projection
84 aNear = aProjMat.GetValue (2, 3) / (- 1.0f + aProjMat.GetValue (2, 2));
85 aFar = aProjMat.GetValue (2, 3) / ( 1.0f + aProjMat.GetValue (2, 2));
86 // Near plane
87 nLeft = aNear * (aProjMat.GetValue (0, 2) - 1.0f) / aProjMat.GetValue (0, 0);
88 nRight = aNear * (aProjMat.GetValue (0, 2) + 1.0f) / aProjMat.GetValue (0, 0);
89 nTop = aNear * (aProjMat.GetValue (1, 2) + 1.0f) / aProjMat.GetValue (1, 1);
90 nBottom = aNear * (aProjMat.GetValue (1, 2) - 1.0f) / aProjMat.GetValue (1, 1);
91 // Far plane
92 fLeft = aFar * (aProjMat.GetValue (0, 2) - 1.0f) / aProjMat.GetValue (0, 0);
93 fRight = aFar * (aProjMat.GetValue (0, 2) + 1.0f) / aProjMat.GetValue (0, 0);
94 fTop = aFar * (aProjMat.GetValue (1, 2) + 1.0f) / aProjMat.GetValue (1, 1);
95 fBottom = aFar * (aProjMat.GetValue (1, 2) - 1.0f) / aProjMat.GetValue (1, 1);
96 }
97 else
98 {
99 // handle orthographic projection
100 aNear = (1.0f / aProjMat.GetValue (2, 2)) * (aProjMat.GetValue (2, 3) + 1.0f);
101 aFar = (1.0f / aProjMat.GetValue (2, 2)) * (aProjMat.GetValue (2, 3) - 1.0f);
102 // Near plane
103 nLeft = ( 1.0f + aProjMat.GetValue (0, 3)) / (-aProjMat.GetValue (0, 0));
104 fLeft = nLeft;
105 nRight = ( 1.0f - aProjMat.GetValue (0, 3)) / aProjMat.GetValue (0, 0);
106 fRight = nRight;
107 nTop = ( 1.0f - aProjMat.GetValue (1, 3)) / aProjMat.GetValue (1, 1);
108 fTop = nTop;
109 nBottom = (-1.0f - aProjMat.GetValue (1, 3)) / aProjMat.GetValue (1, 1);
110 fBottom = nBottom;
111 }
112
113 OpenGl_Vec4 aLeftTopNear (nLeft, nTop, -aNear, 1.0f), aRightBottomFar (fRight, fBottom, -aFar, 1.0f);
114 OpenGl_Vec4 aLeftBottomNear (nLeft, nBottom, -aNear, 1.0f), aRightTopFar (fRight, fTop, -aFar, 1.0f);
115 OpenGl_Vec4 aRightBottomNear (nRight, nBottom, -aNear, 1.0f), aLeftTopFar (fLeft, fTop, -aFar, 1.0f);
116 OpenGl_Vec4 aRightTopNear (nRight, nTop, -aNear, 1.0f), aLeftBottomFar (fLeft, fBottom, -aFar, 1.0f);
117
118 const OpenGl_Mat4 aViewProj = aModelMat * aProjMat;
119 OpenGl_Mat4 anInvModelView;
120 aModelMat.Inverted(anInvModelView);
121
122 myClipVerts[ClipVert_LeftTopNear] = anInvModelView * aLeftTopNear;
123 myClipVerts[ClipVert_RightBottomFar] = anInvModelView * aRightBottomFar;
124 myClipVerts[ClipVert_LeftBottomNear] = anInvModelView * aLeftBottomNear;
125 myClipVerts[ClipVert_RightTopFar] = anInvModelView * aRightTopFar;
126 myClipVerts[ClipVert_RightBottomNear] = anInvModelView * aRightBottomNear;
127 myClipVerts[ClipVert_LeftTopFar] = anInvModelView * aLeftTopFar;
128 myClipVerts[ClipVert_RightTopNear] = anInvModelView * aRightTopNear;
129 myClipVerts[ClipVert_LeftBottomFar] = anInvModelView * aLeftBottomFar;
130
131 // UNNORMALIZED!
132 myClipPlanes[Plane_Left] = aViewProj.GetRow (3) + aViewProj.GetRow (0);
133 myClipPlanes[Plane_Right] = aViewProj.GetRow (3) - aViewProj.GetRow (0);
134 myClipPlanes[Plane_Top] = aViewProj.GetRow (3) - aViewProj.GetRow (1);
135 myClipPlanes[Plane_Bottom] = aViewProj.GetRow (3) + aViewProj.GetRow (1);
136 myClipPlanes[Plane_Near] = aViewProj.GetRow (3) + aViewProj.GetRow (2);
137 myClipPlanes[Plane_Far] = aViewProj.GetRow (3) - aViewProj.GetRow (2);
138
139 gp_Pnt aPtCenter = theCamera->Center();
140 OpenGl_Vec4 aCenter (static_cast<Standard_ShortReal> (aPtCenter.X()),
141 static_cast<Standard_ShortReal> (aPtCenter.Y()),
142 static_cast<Standard_ShortReal> (aPtCenter.Z()),
143 1.0f);
144
145 for (Standard_Integer aPlaneIter = 0; aPlaneIter < PlanesNB; ++aPlaneIter)
146 {
147 OpenGl_Vec4 anEq = myClipPlanes[aPlaneIter];
148 if (SignedPlanePointDistance (anEq, aCenter) > 0)
149 {
150 anEq *= -1.0f;
151 myClipPlanes[aPlaneIter] = anEq;
152 }
153 }
154}
155
156// =======================================================================
157// function : SignedPlanePointDistance
158// purpose :
159// =======================================================================
160Standard_ShortReal OpenGl_BVHTreeSelector::SignedPlanePointDistance (const OpenGl_Vec4& theNormal,
161 const OpenGl_Vec4& thePnt)
162{
163 const Standard_ShortReal aNormLength = std::sqrt (theNormal.x() * theNormal.x()
164 + theNormal.y() * theNormal.y()
165 + theNormal.z() * theNormal.z());
3c648527 166
167 if (aNormLength < FLT_EPSILON)
168 return 0.0f;
169
b7cd4ba7 170 const Standard_ShortReal anInvNormLength = 1.0f / aNormLength;
171 const Standard_ShortReal aD = theNormal.w() * anInvNormLength;
172 const Standard_ShortReal anA = theNormal.x() * anInvNormLength;
173 const Standard_ShortReal aB = theNormal.y() * anInvNormLength;
174 const Standard_ShortReal aC = theNormal.z() * anInvNormLength;
175 return aD + (anA * thePnt.x() + aB * thePnt.y() + aC * thePnt.z());
176}
177
178// =======================================================================
179// function : CacheClipPtsProjections
180// purpose : Caches view volume's vertices projections along its normals and AABBs dimensions
181// Must be called at the beginning of each BVH tree traverse loop
182// =======================================================================
183void OpenGl_BVHTreeSelector::CacheClipPtsProjections()
184{
185 Standard_ShortReal aProjectedVerts[ClipVerticesNB];
186 for (Standard_Integer aPlaneIter = 0; aPlaneIter < PlanesNB; ++aPlaneIter)
187 {
188 const OpenGl_Vec4 aPlane = myClipPlanes[aPlaneIter];
189 for (Standard_Integer aCornerIter = 0; aCornerIter < ClipVerticesNB; ++aCornerIter)
190 {
191 Standard_ShortReal aProjection = DotProduct (aPlane, myClipVerts[aCornerIter]);
192 aProjectedVerts[aCornerIter] = aProjection;
193 }
194 myMaxClipProjectionPts[aPlaneIter] = *std::max_element (aProjectedVerts, aProjectedVerts + ClipVerticesNB);
195 myMinClipProjectionPts[aPlaneIter] = *std::min_element (aProjectedVerts, aProjectedVerts + ClipVerticesNB);
196 }
197
198 OpenGl_Vec4 aDimensions[3] =
199 {
200 OpenGl_Vec4 (1.0f, 0.0f, 0.0f, 1.0f),
201 OpenGl_Vec4 (0.0f, 1.0f, 0.0f, 1.0f),
202 OpenGl_Vec4 (0.0f, 0.0f, 1.0f, 1.0f)
203 };
204
205 for (Standard_Integer aDim = 0; aDim < 3; ++aDim)
206 {
207 for (Standard_Integer aCornerIter = 0; aCornerIter < ClipVerticesNB; ++aCornerIter)
208 {
209 Standard_ShortReal aProjection = DotProduct (aDimensions[aDim], myClipVerts[aCornerIter]);
210 aProjectedVerts[aCornerIter] = aProjection;
211 }
212 myMaxOrthoProjectionPts[aDim] = *std::max_element (aProjectedVerts, aProjectedVerts + ClipVerticesNB);
213 myMinOrthoProjectionPts[aDim] = *std::min_element (aProjectedVerts, aProjectedVerts + ClipVerticesNB);
214 }
215}
216
217// =======================================================================
218// function : Intersect
219// purpose : Detects if AABB overlaps view volume using separating axis theorem (SAT)
220// =======================================================================
221Standard_Boolean OpenGl_BVHTreeSelector::Intersect (const OpenGl_Vec4& theMinPt,
222 const OpenGl_Vec4& theMaxPt) const
223{
224 // E1
225 // |_ E0
226 // /
227 // E2
228 const OpenGl_Vec4 aShiftedBoxMax = theMaxPt - theMinPt;
229 Standard_ShortReal aBoxProjMax = 0.0f, aBoxProjMin = 0.0f;
230 Standard_ShortReal aFrustumProjMax = 0.0f, aFrustumProjMin = 0.0f;
231
232 // E0 test
233 aBoxProjMax = aShiftedBoxMax.x();
234 aFrustumProjMax = myMaxOrthoProjectionPts[0] - DotProduct (OpenGl_Vec4 (1.0f, 0.0f, 0.0f, 1.0f), theMinPt);
235 aFrustumProjMin = myMinOrthoProjectionPts[0] - DotProduct (OpenGl_Vec4 (1.0f, 0.0f, 0.0f, 1.0f), theMinPt);
236 if (aBoxProjMin > aFrustumProjMax
237 || aBoxProjMax < aFrustumProjMin)
238 {
239 return Standard_False;
240 }
241
242 // E1 test
243 aBoxProjMax = aShiftedBoxMax.y();
244 aFrustumProjMax = myMaxOrthoProjectionPts[1] - DotProduct (OpenGl_Vec4 (0.0f, 1.0f, 0.0f, 1.0f), theMinPt);
245 aFrustumProjMin = myMinOrthoProjectionPts[1] - DotProduct (OpenGl_Vec4 (0.0f, 1.0f, 0.0f, 1.0f), theMinPt);
246 if (aBoxProjMin > aFrustumProjMax
247 || aBoxProjMax < aFrustumProjMin)
248 {
249 return Standard_False;
250 }
251
252 // E2 test
253 aBoxProjMax = aShiftedBoxMax.z();
254 aFrustumProjMax = myMaxOrthoProjectionPts[2] - DotProduct (OpenGl_Vec4 (0.0f, 0.0f, 1.0f, 1.0f), theMinPt);
255 aFrustumProjMin = myMinOrthoProjectionPts[2] - DotProduct (OpenGl_Vec4 (0.0f, 0.0f, 1.0f, 1.0f), theMinPt);
256 if (aBoxProjMin > aFrustumProjMax
257 || aBoxProjMax < aFrustumProjMin)
258 {
259 return Standard_False;
260 }
261
262 const Standard_Integer anIncFactor = myIsProjectionParallel ? 2 : 1;
263 for (Standard_Integer aPlaneIter = 0; aPlaneIter < 5; aPlaneIter += anIncFactor)
264 {
265 OpenGl_Vec4 aPlane = myClipPlanes[aPlaneIter];
266 OpenGl_Vec4 aPt1 (0.0f), aPt2 (0.0f);
267 aPt1 = BinarySign (aPlane) * aShiftedBoxMax;
268 aBoxProjMax = DotProduct (aPlane, aPt1);
269 aFrustumProjMax = myMaxClipProjectionPts[aPlaneIter] - DotProduct (aPlane, theMinPt);
270 aFrustumProjMin = myMinClipProjectionPts[aPlaneIter] - DotProduct (aPlane, theMinPt);
271 if (aFrustumProjMin < aBoxProjMax
272 && aBoxProjMax < aFrustumProjMax)
273 {
274 continue;
275 }
276
277 aPt2 = InversedBinarySign (aPlane) * aShiftedBoxMax;
278 aBoxProjMin = DotProduct (aPlane, aPt2);
279 if (aBoxProjMin > aFrustumProjMax
280 || aBoxProjMax < aFrustumProjMin)
281 {
282 return Standard_False;
283 }
284 }
285
286 return Standard_True;
287}