projects / occt.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
0023995: GeomAPI_ExtremaCurveSurface : wrong result between a curve and a plane
[occt.git] / src / Voxel / Voxel_FastConverter.cxx
... / ...
CommitLineData
1// Created on: 2008-05-30
2// Created by: Vladislav ROMASHKO
3// Copyright (c) 2008-2012 OPEN CASCADE SAS
4//
5// The content of this file is subject to the Open CASCADE Technology Public
6// License Version 6.5 (the "License"). You may not use the content of this file
7// except in compliance with the License. Please obtain a copy of the License
8// at http://www.opencascade.org and read it completely before using this file.
9//
10// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
11// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
12//
13// The Original Code and all software distributed under the License is
14// distributed on an "AS IS" basis, without warranty of any kind, and the
15// Initial Developer hereby disclaims all such warranties, including without
16// limitation, any warranties of merchantability, fitness for a particular
17// purpose or non-infringement. Please see the License for the specific terms
18// and conditions governing the rights and limitations under the License.
19
20
21#include <Voxel_FastConverter.ixx>
22
23#include <Bnd_Box.hxx>
24#include <BRep_Tool.hxx>
25#include <BRepBndLib.hxx>
26#include <BRepMesh.hxx>
27
28#include <TopoDS.hxx>
29#include <TopoDS_Face.hxx>
30#include <TopExp_Explorer.hxx>
31
32#include <gp_Lin2d.hxx>
33#include <gce_MakePln.hxx>
34
35#include <ElSLib.hxx>
36#include <Poly_Triangulation.hxx>
37#include <IntAna2d_AnaIntersection.hxx>
38
39// Printing the progress in stdout.
40//#define CONV_DUMP
41
42Voxel_FastConverter::Voxel_FastConverter(const TopoDS_Shape& shape,
43 Voxel_ROctBoolDS& voxels,
44 const Standard_Real deflection,
45 const Standard_Integer nbx,
46 const Standard_Integer nby,
47 const Standard_Integer nbz,
48 const Standard_Integer nbthreads,
49 const Standard_Boolean useExistingTriangulation)
50:myShape(shape),myVoxels(&voxels),
51 myDeflection(deflection),
52 myNbX(nbx),myNbY(nby),myNbZ(nbz),
53 myNbThreads(nbthreads),myIsBool(2),
54 myNbTriangles(0),
55 myUseExistingTriangulation(useExistingTriangulation)
56{
57 Init();
58}
59
60Voxel_FastConverter::Voxel_FastConverter(const TopoDS_Shape& shape,
61 Voxel_BoolDS& voxels,
62 const Standard_Real deflection,
63 const Standard_Integer nbx,
64 const Standard_Integer nby,
65 const Standard_Integer nbz,
66 const Standard_Integer nbthreads,
67 const Standard_Boolean useExistingTriangulation)
68:myShape(shape),myVoxels(&voxels),
69 myDeflection(deflection),
70 myNbX(nbx),myNbY(nby),myNbZ(nbz),
71 myNbThreads(nbthreads),myIsBool(1),
72 myNbTriangles(0),
73 myUseExistingTriangulation(useExistingTriangulation)
74{
75 Init();
76}
77
78Voxel_FastConverter::Voxel_FastConverter(const TopoDS_Shape& shape,
79 Voxel_ColorDS& voxels,
80 const Standard_Real deflection,
81 const Standard_Integer nbx,
82 const Standard_Integer nby,
83 const Standard_Integer nbz,
84 const Standard_Integer nbthreads,
85 const Standard_Boolean useExistingTriangulation)
86:myShape(shape),myVoxels(&voxels),
87 myDeflection(deflection),
88 myNbX(nbx),myNbY(nby),myNbZ(nbz),
89 myNbThreads(nbthreads),myIsBool(0),
90 myNbTriangles(0),
91 myUseExistingTriangulation(useExistingTriangulation)
92{
93 Init();
94}
95
96void Voxel_FastConverter::Init()
97{
98 if (myShape.IsNull())
99 return;
100 if (myNbThreads < 1)
101 return;
102
103 // Check number of splits.
104 Voxel_DS* voxels = (Voxel_DS*) myVoxels;
105 if (voxels->GetNbX() != myNbX || voxels->GetNbY() != myNbY || voxels->GetNbZ() != myNbZ)
106 {
107 // Compute boundary box of the shape
108 Bnd_Box box;
109 BRepBndLib::Add(myShape, box);
110
111 // Define the voxel model by means of the boundary box of shape
112 Standard_Real xmin, ymin, zmin, xmax, ymax, zmax;
113 box.Get(xmin, ymin, zmin, xmax, ymax, zmax);
114
115 // Initialize the voxels.
116 if (myIsBool == 2)
117 ((Voxel_ROctBoolDS*) voxels)->Init(xmin, ymin, zmin, xmax - xmin, ymax - ymin, zmax - zmin, myNbX, myNbY, myNbZ);
118 else if (myIsBool == 1)
119 ((Voxel_BoolDS*) voxels)->Init(xmin, ymin, zmin, xmax - xmin, ymax - ymin, zmax - zmin, myNbX, myNbY, myNbZ);
120 else if (myIsBool == 0)
121 ((Voxel_ColorDS*) voxels)->Init(xmin, ymin, zmin, xmax - xmin, ymax - ymin, zmax - zmin, myNbX, myNbY, myNbZ);
122 }
123
124 // Check presence of triangulation.
125 TopLoc_Location L;
126 Standard_Boolean triangulate = Standard_False;
127 TopExp_Explorer expl(myShape, TopAbs_FACE);
128 if(myUseExistingTriangulation == Standard_False)
129 {
130 for (; expl.More(); expl.Next())
131 {
132 const TopoDS_Face & F = TopoDS::Face(expl.Current());
133 Handle(Poly_Triangulation) T = BRep_Tool::Triangulation(F, L);
134 if (T.IsNull() || (T->Deflection() > myDeflection))
135 {
136 triangulate = Standard_True;
137 break;
138 }
139 }
140 }
141
142 // Re-create the triangulation.
143 if (triangulate)
144 {
145 BRepMesh::Mesh(myShape, myDeflection);
146 }
147
148 // Compute the number of triangles.
149 myNbTriangles = 0;
150 expl.Init(myShape, TopAbs_FACE);
151 for (; expl.More(); expl.Next())
152 {
153 const TopoDS_Face & F = TopoDS::Face(expl.Current());
154 Handle(Poly_Triangulation) T = BRep_Tool::Triangulation(F, L);
155 if (T.IsNull() == Standard_False)
156 myNbTriangles += T->NbTriangles();
157 }
158}
159
160// Destructor
161void Voxel_FastConverter::Destroy()
162{
163
164}
165
166Standard_Boolean Voxel_FastConverter::Convert(Standard_Integer& progress,
167 const Standard_Integer ithread)
168{
169 if (ithread == 1)
170 progress = 0;
171#ifdef CONV_DUMP
172 if (ithread == 1)
173 printf("Progress = %d \r", progress);
174#endif
175
176 if (myNbX <= 0 || myNbY <= 0 || myNbZ <= 0)
177 return Standard_False;
178
179 if(myNbTriangles == 0)
180 return Standard_False;
181
182 // Half of diagonal of a voxel
183 Voxel_DS* ds = (Voxel_DS*) myVoxels;
184 Standard_Real dx = ds->GetXLen() / (Standard_Real) ds->GetNbX(),
185 dy = ds->GetYLen() / (Standard_Real) ds->GetNbY(),
186 dz = ds->GetZLen() / (Standard_Real) ds->GetNbZ();
187 Standard_Real hdiagonal = sqrt(dx * dx + dy * dy + dz * dz);
188 hdiagonal /= 2.0;
189
190 // Compute the scope of triangles for current thread
191 Standard_Integer start_thread_triangle = 1, end_thread_triangle = myNbTriangles, ithread_triangle = 0;
192 start_thread_triangle = (ithread - 1) * (myNbTriangles / myNbThreads) + 1;
193 end_thread_triangle = (ithread - 0) * (myNbTriangles / myNbThreads);
194
195 // Convert
196 TopLoc_Location L;
197 Standard_Integer iprogress = 0, prev_progress = 0;
198 Standard_Integer n1, n2, n3;
199 Standard_Integer ixmin, iymin, izmin, ixmax, iymax, izmax;
200 Standard_Real xmin, ymin, zmin, xmax, ymax, zmax;
201 TopExp_Explorer expl(myShape, TopAbs_FACE);
202 for (; expl.More(); expl.Next())
203 {
204 const TopoDS_Face & F = TopoDS::Face(expl.Current());
205 Handle(Poly_Triangulation) T = BRep_Tool::Triangulation(F, L);
206 if (T.IsNull())
207 continue;
208
209 gp_Trsf trsf;
210 Standard_Boolean transform = !L.IsIdentity();
211 if (transform)
212 trsf = L.Transformation();
213
214 const TColgp_Array1OfPnt& nodes = T->Nodes();
215 const Poly_Array1OfTriangle& triangles = T->Triangles();
216 Standard_Integer itriangle = triangles.Lower(), nb_triangles = triangles.Upper();
217 for (; itriangle <= nb_triangles; itriangle++)
218 {
219 ithread_triangle++;
220 if (ithread_triangle < start_thread_triangle )
221 continue;
222 if (ithread_triangle > end_thread_triangle)
223 {
224 if (ithread == 1)
225 progress = 100;
226#ifdef CONV_DUMP
227 if (ithread == 1)
228 printf("Progress = %d \r", progress);
229#endif
230 return Standard_True;
231 }
232
233 const Poly_Triangle& t = triangles.Value(itriangle);
234 t.Get(n1, n2, n3);
235 gp_Pnt p1 = nodes.Value(n1);
236 gp_Pnt p2 = nodes.Value(n2);
237 gp_Pnt p3 = nodes.Value(n3);
238 if (transform)
239 {
240 p1.Transform(trsf);
241 p2.Transform(trsf);
242 p3.Transform(trsf);
243 }
244
245 // Get boundary box of the triangle
246 GetBndBox(p1, p2, p3, xmin, ymin, zmin, xmax, ymax, zmax);
247
248 // Find the range of voxels inside the boudary box of the triangle.
249 if (!ds->GetVoxel(xmin, ymin, zmin, ixmin, iymin, izmin))
250 continue;
251 if (!ds->GetVoxel(xmax, ymax, zmax, ixmax, iymax, izmax))
252 continue;
253
254 // Refuse voxels for whom distance from their center to plane of triangle is greater than half of diagonal.
255 // Make a line from center of each voxel to the center of triangle and
256 // compute intersection of the line with sides of triangle.
257 // Refuse the voxel in case of intersection.
258 gce_MakePln mkPlane(p1, p2, p3);
259 if (!mkPlane.IsDone())
260 continue;
261 gp_Pln plane = mkPlane.Value();
262 ComputeVoxelsNearTriangle(plane, p1, p2, p3, hdiagonal, ixmin, iymin, izmin, ixmax, iymax, izmax);
263
264 // Progress
265 if (ithread == 1)
266 {
267 iprogress++;
268 progress = (Standard_Integer) ( (Standard_Real) iprogress / (Standard_Real) myNbTriangles * 100.0 );
269 }
270#ifdef CONV_DUMP
271 if (ithread == 1 && prev_progress != progress)
272 {
273 printf("Progress = %d \r", progress);
274 prev_progress = progress;
275 }
276#endif
277
278 } // iteration of triangles
279 } // iteration of faces
280
281 if (ithread == 1)
282 progress = 100;
283#ifdef CONV_DUMP
284 if (ithread == 1)
285 printf("Progress = %d \r", progress);
286#endif
287 return Standard_True;
288}
289
290Standard_Boolean Voxel_FastConverter::FillInVolume(const Standard_Byte inner,
291 const Standard_Integer ithread)
292{
293 Voxel_DS* ds = (Voxel_DS*) myVoxels;
294 Standard_Integer ix, iy, iz, nbx = ds->GetNbX(), nby = ds->GetNbY(), nbz = ds->GetNbZ();
295 Standard_Boolean prev_surface, surface, volume;
296
297 if (inner)
298 {
299 // Fill-in internal voxels by the value "inner"
300 for (ix = 0; ix < nbx; ix++)
301 {
302 for (iy = 0; iy < nby; iy++)
303 {
304 // Check existence of volume.
305 volume = Standard_False;
306 surface = Standard_False;
307 prev_surface = Standard_False;
308 for (iz = 0; iz < nbz; iz++)
309 {
310 surface = (myIsBool == 1) ?
311 ((Voxel_BoolDS*)myVoxels)->Get(ix, iy, iz) == Standard_True :
312 ((Voxel_ColorDS*)myVoxels)->Get(ix, iy, iz) > 0;
313 if (prev_surface && !surface)
314 {
315 volume = !volume;
316 }
317 prev_surface = surface;
318 }
319 if (volume)
320 continue;
321
322 // Fill-in the volume.
323 volume = Standard_False;
324 surface = Standard_False;
325 prev_surface = Standard_False;
326 for (iz = 0; iz < nbz; iz++)
327 {
328 surface = (myIsBool == 1) ?
329 ((Voxel_BoolDS*)myVoxels)->Get(ix, iy, iz) == Standard_True :
330 ((Voxel_ColorDS*)myVoxels)->Get(ix, iy, iz) > 0;
331 if (prev_surface && !surface)
332 {
333 volume = !volume;
334 }
335 if (volume && !surface)
336 {
337 (myIsBool == 1) ? ((Voxel_BoolDS*)myVoxels)->Set(ix, iy, iz, inner) :
338 ((Voxel_ColorDS*)myVoxels)->Set(ix, iy, iz, inner);
339 }
340 prev_surface = surface;
341 }
342 }
343 }
344 }
345 else
346 {
347 // Set value of interbal voxels to 0 ("inner" = 0)
348 Standard_Boolean next_surface;
349 for (ix = 0; ix < nbx; ix++)
350 {
351 for (iy = 0; iy < nby; iy++)
352 {
353 volume = Standard_False;
354 surface = Standard_False;
355 prev_surface = Standard_False;
356 next_surface = Standard_False;
357 for (iz = 0; iz < nbz; iz++)
358 {
359 surface = (myIsBool == 1) ?
360 ((Voxel_BoolDS*)myVoxels)->Get(ix, iy, iz) == Standard_True :
361 ((Voxel_ColorDS*)myVoxels)->Get(ix, iy, iz) > 0;
362 if (prev_surface != surface)
363 {
364 volume = !volume;
365 }
366 if (volume && iz + 1 < nbz)
367 {
368 next_surface = (myIsBool == 1) ?
369 ((Voxel_BoolDS*)myVoxels)->Get(ix, iy, iz + 1) == Standard_True :
370 ((Voxel_ColorDS*)myVoxels)->Get(ix, iy, iz + 1) > 0;
371 }
372 if (volume && prev_surface == surface && next_surface)
373 {
374 (myIsBool == 1) ? ((Voxel_BoolDS*)myVoxels)->Set(ix, iy, iz, inner) :
375 ((Voxel_ColorDS*)myVoxels)->Set(ix, iy, iz, inner);
376 }
377 prev_surface = surface;
378 }
379 }
380 }
381 }
382
383 return Standard_True;
384}
385
386void Voxel_FastConverter::GetBndBox(const gp_Pnt& p1,
387 const gp_Pnt& p2,
388 const gp_Pnt& p3,
389 Standard_Real& xmin,
390 Standard_Real& ymin,
391 Standard_Real& zmin,
392 Standard_Real& xmax,
393 Standard_Real& ymax,
394 Standard_Real& zmax) const
395{
396 // P1:
397 xmin = p1.X();
398 ymin = p1.Y();
399 zmin = p1.Z();
400 xmax = p1.X();
401 ymax = p1.Y();
402 zmax = p1.Z();
403 // P2:
404 if (xmin > p2.X())
405 xmin = p2.X();
406 if (ymin > p2.Y())
407 ymin = p2.Y();
408 if (zmin > p2.Z())
409 zmin = p2.Z();
410 if (xmax < p2.X())
411 xmax = p2.X();
412 if (ymax < p2.Y())
413 ymax = p2.Y();
414 if (zmax < p2.Z())
415 zmax = p2.Z();
416 // P3:
417 if (xmin > p3.X())
418 xmin = p3.X();
419 if (ymin > p3.Y())
420 ymin = p3.Y();
421 if (zmin > p3.Z())
422 zmin = p3.Z();
423 if (xmax < p3.X())
424 xmax = p3.X();
425 if (ymax < p3.Y())
426 ymax = p3.Y();
427 if (zmax < p3.Z())
428 zmax = p3.Z();
429}
430
431// This method is copied from Voxel_ShapeIntersector.cxx
432static Standard_Boolean mayIntersect(const gp_Pnt2d& p11, const gp_Pnt2d& p12,
433 const gp_Pnt2d& p21, const gp_Pnt2d& p22)
434{
435 if (p11.X() > p21.X() && p11.X() > p22.X() && p12.X() > p21.X() && p12.X() > p22.X())
436 return Standard_False;
437 if (p11.X() < p21.X() && p11.X() < p22.X() && p12.X() < p21.X() && p12.X() < p22.X())
438 return Standard_False;
439 if (p11.Y() > p21.Y() && p11.Y() > p22.Y() && p12.Y() > p21.Y() && p12.Y() > p22.Y())
440 return Standard_False;
441 if (p11.Y() < p21.Y() && p11.Y() < p22.Y() && p12.Y() < p21.Y() && p12.Y() < p22.Y())
442 return Standard_False;
443 return Standard_True;
444}
445
446void Voxel_FastConverter::ComputeVoxelsNearTriangle(const gp_Pln& plane,
447 const gp_Pnt& p1,
448 const gp_Pnt& p2,
449 const gp_Pnt& p3,
450 const Standard_Real hdiagonal,
451 const Standard_Integer ixmin,
452 const Standard_Integer iymin,
453 const Standard_Integer izmin,
454 const Standard_Integer ixmax,
455 const Standard_Integer iymax,
456 const Standard_Integer izmax) const
457{
458 gp_Pnt pc;
459 Standard_Real xc, yc, zc, uc, vc, u1, v1, u2, v2, u3, v3;
460 Standard_Integer ix, iy, iz;
461 IntAna2d_AnaIntersection intersector2d;
462
463 // Project points of triangle onto the plane
464 ElSLib::Parameters(plane, p1, u1, v1);
465 ElSLib::Parameters(plane, p2, u2, v2);
466 ElSLib::Parameters(plane, p3, u3, v3);
467
468 // Make lines of triangle
469 gp_Pnt2d p2d1(u1, v1), p2d2(u2, v2), p2d3(u3, v3), p2dt((u1+u2+u3)/3.0,(v1+v2+v3)/3.0), p2dc;
470 gp_Vec2d v2d12(p2d1, p2d2), v2d23(p2d2, p2d3), v2d31(p2d3, p2d1);
471 gp_Lin2d L1(p2d1, v2d12), L2(p2d2, v2d23), L3(p2d3, v2d31), Lv;
472 Standard_Real d1 = p2d1.Distance(p2d2) - Precision::Confusion(),
473 d2 = p2d2.Distance(p2d3) - Precision::Confusion(),
474 d3 = p2d3.Distance(p2d1) - Precision::Confusion(), dv;
475
476 Voxel_DS* ds = (Voxel_DS*) myVoxels;
477 for (ix = ixmin; ix <= ixmax; ix++)
478 {
479 for (iy = iymin; iy <= iymax; iy++)
480 {
481 for (iz = izmin; iz <= izmax; iz++)
482 {
483 ds->GetCenter(ix, iy, iz, xc, yc, zc);
484 pc.SetCoord(xc, yc, zc);
485 if (plane.Distance(pc) < hdiagonal)
486 {
487 ElSLib::Parameters(plane, pc, uc, vc);
488 p2dc.SetCoord(uc, vc);
489
490 gp_Vec2d v2dct(p2dc, p2dt);
491 dv = v2dct.Magnitude() - Precision::Confusion();
492 Lv.SetLocation(p2dc);
493 Lv.SetDirection(v2dct);
494
495 // Side 1:
496 if (mayIntersect(p2d1, p2d2, p2dc, p2dt))
497 {
498 intersector2d.Perform(Lv, L1);
499 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
500 {
501 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
502 Standard_Real param1 = i2d.ParamOnFirst();
503 Standard_Real param2 = i2d.ParamOnSecond();
504 if (param1 > Precision::Confusion() && param1 < dv &&
505 param2 > Precision::Confusion() && param2 < d1)
506 {
507 continue;
508 }
509 }
510 }
511
512 // Side 2:
513 if (mayIntersect(p2d2, p2d3, p2dc, p2dt))
514 {
515 intersector2d.Perform(Lv, L2);
516 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
517 {
518 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
519 Standard_Real param1 = i2d.ParamOnFirst();
520 Standard_Real param2 = i2d.ParamOnSecond();
521 if (param1 > Precision::Confusion() && param1 < dv &&
522 param2 > Precision::Confusion() && param2 < d2)
523 {
524 continue;
525 }
526 }
527 }
528
529 // Side 3:
530 if (mayIntersect(p2d3, p2d1, p2dc, p2dt))
531 {
532 intersector2d.Perform(Lv, L3);
533 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
534 {
535 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
536 Standard_Real param1 = i2d.ParamOnFirst();
537 Standard_Real param2 = i2d.ParamOnSecond();
538 if (param1 > Precision::Confusion() && param1 < dv &&
539 param2 > Precision::Confusion() && param2 < d3)
540 {
541 continue;
542 }
543 }
544 }
545
546 // Set positive value to this voxel:
547 switch (myIsBool)
548 {
549 case 0:
550 ((Voxel_ColorDS*) myVoxels)->Set(ix, iy, iz, 15);
551 break;
552 case 1:
553 ((Voxel_BoolDS*) myVoxels)->Set(ix, iy, iz, Standard_True);
554 break;
555 case 2:
556 {
557 //((Voxel_ROctBoolDS*) myVoxels)->Set(ix, iy, iz, Standard_True);
558
559 // Check intersection between the triangle & sub-voxels of the voxel.
560 Standard_Real hdiagonal2 = hdiagonal / 2.0, hdiagonal4 = hdiagonal / 4.0;
561 for (Standard_Integer i = 0; i < 8; i++)
562 {
563 ((Voxel_ROctBoolDS*) myVoxels)->GetCenter(ix, iy, iz, i, xc, yc, zc);
564 pc.SetCoord(xc, yc, zc);
565 if (plane.Distance(pc) < hdiagonal2)
566 {
567 ElSLib::Parameters(plane, pc, uc, vc);
568 p2dc.SetCoord(uc, vc);
569
570 gp_Vec2d v2dct(p2dc, p2dt);
571 dv = v2dct.Magnitude() - Precision::Confusion();
572 Lv.SetLocation(p2dc);
573 Lv.SetDirection(v2dct);
574
575 // Side 1:
576 if (mayIntersect(p2d1, p2d2, p2dc, p2dt))
577 {
578 intersector2d.Perform(Lv, L1);
579 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
580 {
581 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
582 Standard_Real param1 = i2d.ParamOnFirst();
583 Standard_Real param2 = i2d.ParamOnSecond();
584 if (param1 > Precision::Confusion() && param1 < dv &&
585 param2 > Precision::Confusion() && param2 < d1)
586 {
587 continue;
588 }
589 }
590 }
591
592 // Side 2:
593 if (mayIntersect(p2d2, p2d3, p2dc, p2dt))
594 {
595 intersector2d.Perform(Lv, L2);
596 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
597 {
598 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
599 Standard_Real param1 = i2d.ParamOnFirst();
600 Standard_Real param2 = i2d.ParamOnSecond();
601 if (param1 > Precision::Confusion() && param1 < dv &&
602 param2 > Precision::Confusion() && param2 < d2)
603 {
604 continue;
605 }
606 }
607 }
608
609 // Side 3:
610 if (mayIntersect(p2d3, p2d1, p2dc, p2dt))
611 {
612 intersector2d.Perform(Lv, L3);
613 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
614 {
615 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
616 Standard_Real param1 = i2d.ParamOnFirst();
617 Standard_Real param2 = i2d.ParamOnSecond();
618 if (param1 > Precision::Confusion() && param1 < dv &&
619 param2 > Precision::Confusion() && param2 < d3)
620 {
621 continue;
622 }
623 }
624 }
625
626 //((Voxel_ROctBoolDS*) myVoxels)->Set(ix, iy, iz, i, Standard_True);
627
628 // Check intersection between the triangle & sub-voxels of the sub-voxel.
629 for (Standard_Integer j = 0; j < 8; j++)
630 {
631 ((Voxel_ROctBoolDS*) myVoxels)->GetCenter(ix, iy, iz, i, j, xc, yc, zc);
632 pc.SetCoord(xc, yc, zc);
633 if (plane.Distance(pc) < hdiagonal4)
634 {
635 ElSLib::Parameters(plane, pc, uc, vc);
636 p2dc.SetCoord(uc, vc);
637
638 gp_Vec2d v2dct(p2dc, p2dt);
639 dv = v2dct.Magnitude() - Precision::Confusion();
640 Lv.SetLocation(p2dc);
641 Lv.SetDirection(v2dct);
642
643 // Side 1:
644 if (mayIntersect(p2d1, p2d2, p2dc, p2dt))
645 {
646 intersector2d.Perform(Lv, L1);
647 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
648 {
649 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
650 Standard_Real param1 = i2d.ParamOnFirst();
651 Standard_Real param2 = i2d.ParamOnSecond();
652 if (param1 > Precision::Confusion() && param1 < dv &&
653 param2 > Precision::Confusion() && param2 < d1)
654 {
655 continue;
656 }
657 }
658 }
659
660 // Side 2:
661 if (mayIntersect(p2d2, p2d3, p2dc, p2dt))
662 {
663 intersector2d.Perform(Lv, L2);
664 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
665 {
666 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
667 Standard_Real param1 = i2d.ParamOnFirst();
668 Standard_Real param2 = i2d.ParamOnSecond();
669 if (param1 > Precision::Confusion() && param1 < dv &&
670 param2 > Precision::Confusion() && param2 < d2)
671 {
672 continue;
673 }
674 }
675 }
676
677 // Side 3:
678 if (mayIntersect(p2d3, p2d1, p2dc, p2dt))
679 {
680 intersector2d.Perform(Lv, L3);
681 if (intersector2d.IsDone() && !intersector2d.ParallelElements() && intersector2d.NbPoints())
682 {
683 const IntAna2d_IntPoint& i2d = intersector2d.Point(1);
684 Standard_Real param1 = i2d.ParamOnFirst();
685 Standard_Real param2 = i2d.ParamOnSecond();
686 if (param1 > Precision::Confusion() && param1 < dv &&
687 param2 > Precision::Confusion() && param2 < d3)
688 {
689 continue;
690 }
691 }
692 }
693
694 ((Voxel_ROctBoolDS*) myVoxels)->Set(ix, iy, iz, i, j, Standard_True);
695 }
696 } // End of "Check level 2".
697
698 }
699 } // End of "Check level 1".
700
701 break;
702 }
703 }
704 }
705 }
706 }
707 }
708}
Open CASCADE Technology repository