Geant4  10.01.p01
UVoxelizer.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * This Software is part of the AIDA Unified Solids Library package *
4 // * See: https://aidasoft.web.cern.ch/USolids *
5 // ********************************************************************
6 //
7 // $Id:$
8 //
9 // --------------------------------------------------------------------
10 //
11 // UVoxelizer
12 //
13 // 19.10.12 Marek Gayer
14 // Created from original implementation in ROOT
15 // --------------------------------------------------------------------
16 
17 #include <iostream>
18 #include <iomanip>
19 #include <sstream>
20 #include <algorithm>
21 #include <set>
22 
23 #include "VUSolid.hh"
24 #include "UUtils.hh"
25 #include "UOrb.hh"
26 #include "UVoxelizer.hh"
27 
28 using namespace std;
29 
30 //______________________________________________________________________________
31 UVoxelizer::UVoxelizer() : fBoundingBox("TessBBox", 1, 1, 1)
32 {
34 
36 
38 }
39 
40 //______________________________________________________________________________
42 {
43 }
44 
45 
47 {
48  // by reserving the size of candidates, we would avoid reallocation of
49  // the vector which could cause fragmentation
50  vector<int> xyz(3), candidates(fTotalCandidates);
51  const vector<int> empty(0);
52  int max[3];
53 
54  for (int i = 0; i <= 2; ++i) max[i] = fBoundaries[i].size();
55  unsigned int size = max[0] * max[1] * max[2];
56 
57  fEmpty.Clear();
58  fEmpty.ResetBitNumber(size - 1);
59  fEmpty.ResetAllBits(true);
60 
61  for (xyz[2] = 0; xyz[2] < max[2]; ++xyz[2])
62  {
63  for (xyz[1] = 0; xyz[1] < max[1]; ++xyz[1])
64  {
65  for (xyz[0] = 0; xyz[0] < max[0]; ++xyz[0])
66  {
67  if (GetCandidatesVoxelArray(xyz, candidates))
68  {
69  int index = GetVoxelsIndex(xyz);
70  fEmpty.SetBitNumber(index, false);
71  // rather than assigning directly with:
72  // "fCandidates[index] = candidates;", in an effort to ensure that
73  // capacity would be just exact, we rather use following 3 lines:
74  vector<int>& c = (fCandidates[index] = empty);
75  c.reserve(candidates.size());
76  c.assign(candidates.begin(), candidates.end());
77  }
78  }
79  }
80  }
81 #ifdef USPECSDEBUG
82  cout << "Non-empty voxels count: " << fCandidates.size() << endl;
83 #endif
84 }
85 
86 #ifdef USOLIDSONLY
87 //______________________________________________________________________________
88 //void UVoxelizer::BuildVoxelLimits(vector<VUSolid*>& solids, vector<UTransform3D*>& transforms)
89 void UVoxelizer::BuildVoxelLimits(vector<VUSolid*>& solids, vector<UTransform3D>& transforms)
90 {
91  // "BuildVoxelLimits"'s aim is to store the coordinates of the origin as well as
92  // the half lengths related to the bounding box of each node.
93  // These quantities are stored in the array "fBoxes" (6 different values per
94  // node.
95  if (int numNodes = solids.size()) // Number of nodes in "multiUnion"
96  {
97  fBoxes.resize(numNodes); // Array which will store the half lengths
98  fNPerSlice = 1 + (fBoxes.size() - 1) / (8 * sizeof(unsigned int));
99 
100  // related to a particular node, but also
101  // the coordinates of its origin
102 
103  UVector3 toleranceVector;
104  toleranceVector.Set(fTolerance);
105 
106  for (int i = 0; i < numNodes; ++i)
107  {
108  VUSolid& solid = *solids[i];
109  //UTransform3D& transform = *transforms[i];
110  UTransform3D transform = transforms[i];
111  UVector3 min, max;
112 
113  solid.Extent(min, max);
114  if (solid.GetEntityType() == "Orb")
115  {
116  UOrb& orb = *(UOrb*) &solid;
117  UVector3 orbToleranceVector;
118  double tolerance = orb.GetRadialTolerance() / 2.0;
119  orbToleranceVector.Set(tolerance);
120  min -= orbToleranceVector;
121  max += orbToleranceVector;
122  }
123  else
124  {
125  min -= toleranceVector;
126  max += toleranceVector;
127  }
128  UUtils::TransformLimits(min, max, transform);
129  fBoxes[i].hlen = (max - min) / 2;
130  fBoxes[i].pos = transform.fTr;
131  }
132  fTotalCandidates = fBoxes.size();
133  }
134 
135 
136 }
137 #endif // USOLIDSONLY
138 
139 void UVoxelizer::BuildVoxelLimits(vector<VUFacet*>& facets)
140 {
141  // "BuildVoxelLimits"'s aim is to store the coordinates of the origin as well as
142  // the half lengths related to the bounding box of each node.
143  // These quantities are stored in the array "fBoxes" (6 different values per
144  // node.
145  if (int numNodes = facets.size()) // Number of nodes in "multiUnion"
146  {
147  fBoxes.resize(numNodes); // Array which will store the half lengths
148  fNPerSlice = 1 + (fBoxes.size() - 1) / (8 * sizeof(unsigned int));
149 
150  UVector3 toleranceVector;
151  toleranceVector.Set(10 * fTolerance);
152 
153  for (int i = 0; i < numNodes; ++i)
154  {
155  VUFacet& facet = *facets[i];
156  UVector3 min, max;
157  UVector3 x(1, 0, 0), y(0, 1, 0), z(0, 0, 1);
158  UVector3 extent;
159  max.Set(facet.Extent(x), facet.Extent(y), facet.Extent(z));
160  min.Set(-facet.Extent(-x), -facet.Extent(-y), -facet.Extent(-z));
161  min -= toleranceVector;
162  max += toleranceVector;
163  UVector3 hlen = (max - min) / 2;
164  fBoxes[i].hlen = hlen;
165  fBoxes[i].pos = min + hlen;
166  }
167  fTotalCandidates = fBoxes.size();
168  }
169 }
170 
171 
172 //______________________________________________________________________________
174 {
175  // "DisplayVoxelLimits" displays the dX, dY, dZ, pX, pY and pZ for each node
176  int numNodes = fBoxes.size();
177  int oldprc = cout.precision();
178  for (int i = 0; i < numNodes; ++i)
179  {
180  cout << setw(10) << setiosflags(ios::fixed) << setprecision(16) <<
181  " -> Node " << i + 1 << ":\n" <<
182  "\t * [x,y,z] = " << fBoxes[i].hlen <<
183  "\t * [x,y,z] = " << fBoxes[i].pos << "\n";
184  }
185  cout.precision(oldprc) ;
186 }
187 
188 //______________________________________________________________________________
189 void UVoxelizer::CreateSortedBoundary(vector<double>& boundary, int axis)
190 {
191  // "CreateBoundaries"'s aim is to determine the slices induced by the bounding fBoxes,
192  // along each axis. The created boundaries are stored in the array "boundariesRaw"
193  int numNodes = fBoxes.size(); // Number of nodes in structure of "UMultiUnion" type
194  // Determination of the boundaries along x, y and z axis:
195  for (int i = 0 ; i < numNodes; ++i)
196  {
197  // For each node, the boundaries are created by using the array "fBoxes"
198  // built in method "BuildVoxelLimits":
199  double p = fBoxes[i].pos[axis], d = fBoxes[i].hlen[axis];
200  // x boundaries
201  boundary[2 * i] = p - d;
202  boundary[2 * i + 1] = p + d;
203  }
204  sort(boundary.begin(), boundary.end());
205 }
206 
208 {
209  // "SortBoundaries" orders the boundaries along each axis (increasing order)
210  // and also does not take into account redundant boundaries, ie if two boundaries
211  // are separated by a distance strictly inferior to "tolerance".
212  // The sorted boundaries are respectively stored in:
213  // * boundaries[0..2]
214 
215  // In addition, the number of elements contained in the three latter arrays are
216  // precise thanks to variables: boundariesCountX, boundariesCountY and boundariesCountZ.
217 
218  if (int numNodes = fBoxes.size())
219  {
220  const double tolerance = fTolerance / 100.0; // Minimal distance to discriminate two boundaries.
221  vector<double> sortedBoundary(2 * numNodes);
222 
223  int considered;
224 
225  for (int j = 0; j <= 2; ++j)
226  {
227  CreateSortedBoundary(sortedBoundary, j);
228  vector<double>& boundary = fBoundaries[j];
229  boundary.clear();
230 
231  considered = 0;
232 
233  for (int i = 0 ; i < 2 * numNodes; ++i)
234  {
235  double newBoundary = sortedBoundary[i];
236  int size = boundary.size();
237  if (!size || abs(boundary[size - 1] - newBoundary) > tolerance)
238  {
239  considered++;
240  {
241  boundary.push_back(newBoundary);
242  continue;
243  }
244  }
245  // If two successive boundaries are too close from each other, only the first one is considered
246  {
247  // cout << "Skipping boundary [" << j << "] : " << i << endl;
248  }
249  }
250 
251  int n = boundary.size();
252  int max = 100000;
253  if (n > max / 2)
254  {
255  int skip = n / (max / 2); // n has to be 2x bigger then 50.000. therefore only from 100.000 re
256  vector<double> reduced;
257  for (int i = 0; i < n; ++i)
258  {
259  // 50 ok for 2k, 1000, 2000
260  int size = boundary.size();
261  if (i % skip == 0 || i == 0 || i == size - 1) // this condition of merging boundaries was wrong, it did not count with right part, which can be completely ommited and not included in final consideration. Now should be OK
262  reduced.push_back(boundary[i]);
263  }
264  boundary = reduced;
265  }
266  }
267  }
268 }
269 
271 {
272  char axis[3] = {'X', 'Y', 'Z'};
273  for (int i = 0; i <= 2; ++i)
274  {
275  cout << " * " << axis[i] << " axis:" << endl << " | ";
277  }
278 }
279 
280 //______________________________________________________________________________
281 void UVoxelizer::DisplayBoundaries(vector<double>& boundaries)
282 {
283  // Prints the positions of the boundaries of the slices on the three axis:
284 
285  int count = boundaries.size();
286  int oldprc = cout.precision();
287  for (int i = 0; i < count; ++i)
288  {
289  cout << setw(10) << setiosflags(ios::fixed) << setprecision(16) << boundaries[i];
290  // printf("%19.15f ", boundaries[0][i]);
291  // cout << boundaries[0][i] << " ";
292  if (i != count - 1) cout << "-> ";
293  }
294  cout << "|" << endl << "Number of boundaries: " << count << endl;
295  cout.precision(oldprc) ;
296 }
297 
298 void UVoxelizer::BuildBitmasks(std::vector<double> boundaries[], UBits bitmasks[])
299 {
300  // "BuildListNodes" stores in the bitmasks solids present in each slice along an axis.
301  // const double localTolerance = 0;
302  int numNodes = fBoxes.size();
303  int bitsPerSlice = GetBitsPerSlice();
304 
305  for (int k = 0; k < 3; k++)
306  {
307  int total = 0;
308  vector<double>& boundary = boundaries[k];
309  int voxelsCount = boundary.size() - 1;
310  UBits& bitmask = bitmasks[k];
311 
312  if (bitmasks != NULL)
313  {
314  bitmask.Clear();
315  bitmask.SetBitNumber(voxelsCount * bitsPerSlice - 1, false); // it is here so we can set the maximum number of bits. this line will rellocate the memory and set all to zero
316  }
317  vector<int>& candidatesCount = fCandidatesCounts[k];
318  candidatesCount.resize(voxelsCount);
319 
320  for (int i = 0 ; i < voxelsCount; ++i) candidatesCount[i] = 0;
321 
322  // Loop on the nodes, number of slices per axis
323  for (int j = 0 ; j < numNodes; j++)
324  {
325  // Determination of the minimum and maximum position along x
326  // of the bounding boxe of each node:
327  double p = fBoxes[j].pos[k], d = fBoxes[j].hlen[k];
328 
329  double min = p - d; // - localTolerance;
330  double max = p + d; // + localTolerance;
331 
332  int i = BinarySearch(boundary, min);
333  if (i < 0) i = 0;
334 
335  do
336  {
337  if (bitmasks != NULL)
338  bitmask.SetBitNumber(i * bitsPerSlice + j);
339 
340  candidatesCount[i]++;
341  total++;
342  i++;
343  }
344  while (max > boundary[i] && i < voxelsCount);
345  }
346 
347 #ifdef DEBUG
348  string fileName = "reduced" + UUtils::ToString(k);
349  UUtils::SaveVectorToExternalFile(candidatesCount, fileName);
350 #endif
351 
352  /*
353  if (voxelsCount < 20000)
354  {
355  int total2 = 0;
356  for(int i = 0 ; i < voxelsCount; ++i)
357  {
358  candidatesCount[i] = 0;
359  // Loop on the nodes, number of slices per axis
360  for(int j = 0 ; j < numNodes; j++)
361  {
362  // Determination of the minimum and maximum position along x
363  // of the bounding boxe of each node:
364  double p = fBoxes[j].pos[k], d = fBoxes[j].hlen[k];
365  double leftBoundary = boundary[i];
366  double rightBoundary = boundary[i+1];
367  double min = p - d; // - localTolerance;
368  double max = p + d; // + localTolerance;
369  if (min < rightBoundary && max > leftBoundary)
370  {
371  // Storage of the number of the node in the array:
372  bitmask.SetBitNumber(i*bitsPerSlice+j);
373  candidatesCount[i]++;
374  total2++;
375  }
376  }
377  }
378  if (total2 != total)
379  total2 = total;
380  }
381  */
382  }
383 #ifdef USPECSDEBUG
384  cout << "Build list nodes completed" << endl;
385 #endif
386 }
387 
388 //______________________________________________________________________________
390 {
391  // Decodes the candidates in mask as string.
392  stringstream ss;
393  int numNodes = fBoxes.size();
394 
395  for (int i = 0; i < numNodes; ++i)
396  if (bits.TestBitNumber(i)) ss << i + 1 << " ";
397  const string& result = ss.str();
398  return result;
399 }
400 
401 
402 //______________________________________________________________________________
404 {
405  char axis[3] = {'X', 'Y', 'Z'};
406  // Prints which solids are present in the slices previously elaborated.
407  int size = 8 * sizeof(int) * fNPerSlice;
408  UBits bits(size);
409 
410  for (int j = 0; j <= 2; ++j)
411  {
412  cout << " * " << axis[j] << " axis:" << endl;
413  int count = fBoundaries[j].size();
414  for (int i = 0; i < count - 1; ++i)
415  {
416  cout << " Slice #" << i + 1 << ": [" << fBoundaries[j][i] << " ; " << fBoundaries[j][i + 1] << "] -> ";
417  bits.Set(size, (const char*)fBitmasks[j].fAllBits + i * fNPerSlice * sizeof(int));
418  string result = GetCandidatesAsString(bits);
419  cout << "[ " << result.c_str() << "] " << endl;
420  }
421  }
422 }
423 
425 {
426  UVector3 min(fBoundaries[0].front(), fBoundaries[1].front(), fBoundaries[2].front());
427  UVector3 max(fBoundaries[0].back(), fBoundaries[1].back(), fBoundaries[2].back());
428  BuildBoundingBox(min, max);
429 }
430 
432 {
433  for (int i = 0; i <= 2; ++i)
434  {
435  double min = amin[i];
436  double max = amax[i];
437  fBoundingBoxSize[i] = (max - min) / 2 + tolerance * 0.5;
438  fBoundingBoxCenter[i] = min + fBoundingBoxSize[i];
439  }
441 }
442 
443 
444 // algorithm -
445 
446 // in order to get balanced voxels, merge should always unite those regions, where the number of voxels is least
447 // the number
448 
449 // we will keep sorted list (std::set) with all voxels. there will be comparator function between two voxels,
450 // which will tell if voxel is less by looking at his right neighbor.
451 // first, we will add all the voxels into the tree.
452 // we will be pick the first item in the tree, merging it, adding the right merged voxel into the a list for
453 // future reduction (fBitmasks will be rebuilded later, therefore they need not to be updated).
454 // the merged voxel need to be added to the tree again, so it's position would be updated
455 
456 
458 {
459  public:
460 
461  vector<UVoxelInfo>& fVoxels;
462 
463  UVoxelComparator(vector<UVoxelInfo>& voxels) : fVoxels(voxels)
464  {
465 
466  }
467 
468  bool operator()(int l, int r)
469  {
470  UVoxelInfo& lv = fVoxels[l], &rv = fVoxels[r];
471  int left = lv.count + fVoxels[lv.next].count;
472  int right = rv.count + fVoxels[rv.next].count;;
473  return (left == right) ? l < r : left < right;
474  }
475 };
476 
477 void UVoxelizer::SetReductionRatio(int maxVoxels, UVector3& reductionRatio)
478 {
479  double maxTotal = (double) fCandidatesCounts[0].size() * fCandidatesCounts[1].size() * fCandidatesCounts[2].size();
480 
481  if (maxVoxels > 0 && maxVoxels < maxTotal)
482  {
483  double ratio = (double) maxVoxels / maxTotal;
484  ratio = pow(ratio, 1. / 3.);
485  if (ratio > 1) ratio = 1;
486  reductionRatio.Set(ratio);
487  }
488 }
489 
490 void UVoxelizer::BuildReduceVoxels(vector<double> boundaries[], UVector3 reductionRatio)
491 {
492  for (int k = 0; k <= 2; ++k)
493  {
494  vector<int>& candidatesCount = fCandidatesCounts[k];
495  int max = candidatesCount.size();
496  vector<UVoxelInfo> voxels(max);
497  UVoxelComparator comp(voxels);
498  set<int, UVoxelComparator> voxelSet(comp);
499  vector<int> mergings;
500 
501  for (int j = 0; j < max; ++j)
502  {
503  UVoxelInfo& voxel = voxels[j];
504  voxel.count = candidatesCount[j];
505  voxel.previous = j - 1;
506  voxel.next = j + 1;
507  voxels[j] = voxel;
508  }
509  // voxels[max - 1].count = 99999999999;
510 
511  for (int j = 0; j < max - 1; ++j) voxelSet.insert(j); // we go to size-1 to make sure we will not merge the last element
512 
513  double reduction = reductionRatio[k];
514  if (reduction != 0)
515  {
516  int count = 0, currentCount;
517  while ((currentCount = voxelSet.size()) > 2)
518  {
519  double currentRatio = 1 - (double) count / max;
520  if (currentRatio <= reduction && currentCount <= 1000)
521  break;
522  const int pos = *voxelSet.begin();
523  mergings.push_back(pos + 1);
524 
525  UVoxelInfo& voxel = voxels[pos];
526  UVoxelInfo& nextVoxel = voxels[voxel.next];
527 
528  if (voxelSet.erase(pos) != 1)
529  {
530  ;// k = k;
531  }
532  if (voxel.next != max - 1)
533  if (voxelSet.erase(voxel.next) != 1)
534  {
535  // k = k;
536  }
537  if (voxel.previous != -1)
538  if (voxelSet.erase(voxel.previous) != 1)
539  {
540  // k = k;
541  }
542  nextVoxel.count += voxel.count;
543  voxel.count = 0;
544  nextVoxel.previous = voxel.previous;
545 
546  if (voxel.next != max - 1)
547  voxelSet.insert(voxel.next);
548 
549  if (voxel.previous != -1)
550  {
551  voxels[voxel.previous].next = voxel.next;
552  voxelSet.insert(voxel.previous);
553  }
554  count++;
555  }
556  }
557 
558  // for (int i = 0; i < max; ++i) cout << voxels[i].count << ", ";
559 
560  if (mergings.size())
561  {
562  sort(mergings.begin(), mergings.end());
563 
564  const vector<double>& boundary = boundaries[k];
565  int mergingsSize = mergings.size();
566  vector<double> reducedBoundary;
567  int skip = mergings[0], i = 0;
568  max = boundary.size();
569  for (int j = 0; j < max; ++j)
570  {
571  if (j != skip)
572  reducedBoundary.push_back(boundary[j]);
573  else if (++i < mergingsSize)
574  {
575  skip = mergings[i];
576  }
577  }
578  boundaries[k] = reducedBoundary;
579  }
580  }
581 }
582 
583 
584 void UVoxelizer::BuildReduceVoxels2(vector<double> boundaries[], UVector3 reductionRatio)
585 {
586  for (int k = 0; k <= 2; ++k)
587  {
588  vector<int>& candidatesCount = fCandidatesCounts[k];
589  int max = candidatesCount.size();
590  int total = 0;
591  for (int i = 0; i < max; ++i) total += candidatesCount[i];
592 
593  double reduction = reductionRatio[k];
594  if (reduction == 0)
595  break;
596 
597  int destination = (int)(reduction * max) + 1;
598  if (destination > 1000) destination = 1000;
599  if (destination < 2) destination = 2;
600  double average = ((double)total / max) / reduction;
601 
602  vector<int> mergings;
603 
604  vector<double>& boundary = boundaries[k];
605  vector<double> reducedBoundary(destination);
606 
607  int sum = 0, cur = 0;
608  for (int i = 0; i < max; ++i)
609 
610  {
611  sum += candidatesCount[i];
612  if (sum > average * (cur + 1) || i == 0)
613  {
614  double val = boundary[i];
615  reducedBoundary[cur] = val;
616  cur++;
617  if (cur == destination)
618  break;
619  }
620  }
621  reducedBoundary[destination - 1] = boundary[max];
622  boundaries[k] = reducedBoundary;
623  }
624 }
625 
626 #ifdef USOLIDSONLY
627 //______________________________________________________________________________
628 void UVoxelizer::Voxelize(vector<VUSolid*>& solids, vector<UTransform3D>& transforms)
629 {
630 
631  BuildVoxelLimits(solids, transforms);
632 
633  BuildBoundaries();
634 
636 
638 
639  BuildEmpty(); // these does not work well for multi-union, actually only makes performance slower, these are only pre-calculated but not used by multi-union
640 
641 // fBoxes.resize(0); // for multiunion, fBoxes are actually needed
642 
643  for (int i = 0; i < 3; ++i) fCandidatesCounts[i].resize(0);
644 }
645 #endif // USOLIDSONLY
646 
647 void UVoxelizer::CreateMiniVoxels(std::vector<double> boundaries[], UBits bitmasks[])
648 {
649  vector<int> voxel(3), maxVoxels(3);
650  for (int i = 0; i <= 2; ++i) maxVoxels[i] = boundaries[i].size();
651 
652  UVector3 point;
653  for (voxel[2] = 0; voxel[2] < maxVoxels[2] - 1; ++voxel[2])
654  {
655  for (voxel[1] = 0; voxel[1] < maxVoxels[1] - 1; ++voxel[1])
656  {
657  for (voxel[0] = 0; voxel[0] < maxVoxels[0] - 1; ++voxel[0])
658  {
659  vector<int> candidates;
660  if (GetCandidatesVoxelArray(voxel, bitmasks, candidates, NULL))
661  {
662  // find a box for corresponding non-empty voxel
663  UVoxelBox box;
664  for (int i = 0; i <= 2; ++i)
665  {
666  int index = voxel[i];
667  const vector<double>& boundary = boundaries[i];
668  double hlen = 0.5 * (boundary[index + 1] - boundary[index]);
669  box.hlen[i] = hlen;
670  box.pos[i] = boundary[index] + hlen;
671  }
672  fVoxelBoxes.push_back(box);
673  vector<int>(candidates).swap(candidates);
674  fVoxelBoxesCandidates.push_back(candidates);
675  }
676  }
677  }
678  }
679 }
680 
681 void UVoxelizer::Voxelize(vector<VUFacet*>& facets)
682 {
683  int maxVoxels = fMaxVoxels;
684  UVector3 reductionRatio = fReductionRatio;
685 
686  int size = facets.size();
687  if (size < 10)
688  for (int i = 0; i < (int) facets.size(); ++i)
689  if (facets[i]->GetNumberOfVertices() > 3) size++;
690 
691  if ((size >= 10 || maxVoxels > 0) && maxVoxels != 0 && maxVoxels != 1)
692  {
693  BuildVoxelLimits(facets);
694 
695  BuildBoundaries();
696 
697  BuildBitmasks(fBoundaries, NULL);
698 
699  if (maxVoxels < 0 && reductionRatio == UVector3())
700  {
701  /*
702  if (fTotalCandidates < 1000)
703  maxVoxels = 10 * fTotalCandidates;
704  else
705  {
706  if (fTotalCandidates < 10000) maxVoxels = 10000;
707  else maxVoxels = fTotalCandidates;
708  }
709  if (fTotalCandidates > 1000000) maxVoxels = 1000000;
710  */
711 
712  maxVoxels = fTotalCandidates;
713  if (fTotalCandidates > 1000000) maxVoxels = 1000000;
714  }
715 
716  SetReductionRatio(maxVoxels, reductionRatio);
717 
719 
720 #ifdef USPECSDEBUG
721  cout << "Total number of voxels: " << fCountOfVoxels << endl;
722 #endif
723 
724  BuildReduceVoxels2(fBoundaries, reductionRatio);
725 
727 
728 #ifdef USPECSDEBUG
729  cout << "Total number of voxels after reduction: " << fCountOfVoxels << endl;
730 #endif
731 
733 
734  UVector3 reductionRatioMini;
735 
736  UBits bitmasksMini[3];
737 
738  // section for building mini voxels
739  std::vector<double> miniBoundaries[3];
740 
741  for (int i = 0; i <= 2; ++i) miniBoundaries[i] = fBoundaries[i];
742 
743  int voxelsCountMini = (fCountOfVoxels >= 1000) ? 100 : fCountOfVoxels / 10;
744  // if (voxelCountMini < 8) voxelCountMini = 8;
745 
746 // voxelsCountMini = 1;
747 
748  SetReductionRatio(voxelsCountMini, reductionRatioMini);
749 
750  BuildReduceVoxels(miniBoundaries, reductionRatioMini);
751 
752 #ifdef USPECSDEBUG
753  int total = CountVoxels(miniBoundaries);
754  cout << "Total number of mini voxels: " << total << endl;
755 #endif
756 
757  BuildBitmasks(miniBoundaries, bitmasksMini);
758 
759  CreateMiniVoxels(miniBoundaries, bitmasksMini);
760 
761  BuildEmpty();
762 
763  // deallocate fields unnecessary during runtime
764  fBoxes.resize(0);
765 
766  for (int i = 0; i < 3; ++i)
767  {
768  fCandidatesCounts[i].resize(0);
769  fBitmasks[i].Clear();
770  }
771  }
772 }
773 
774 
775 //______________________________________________________________________________
776 // "GetCandidates" should compute which solids are possibly contained in
777 // the voxel defined by the three slices characterized by the passed indexes.
778 void UVoxelizer::GetCandidatesVoxel(vector<int>& voxels)
779 {
780  cout << " Candidates in voxel [" << voxels[0] << " ; " << voxels[1] << " ; " << voxels[2] << "]: ";
781  vector<int> candidates;
782  int count = GetCandidatesVoxelArray(voxels, candidates);
783  cout << "[ ";
784  for (int i = 0; i < count; ++i) cout << candidates[i];
785  cout << "] " << endl;
786 }
787 
788 #ifdef USOLIDSONLY
789 
790 inline void findComponents2(unsigned int mask, vector<int>& list, int i)
791 {
792  for (int bit = 0; bit < (int)(8 * sizeof(unsigned int)); bit++)
793  {
794  if (mask & 1)
795  list.push_back(8 * sizeof(unsigned int)*i + bit);
796  if (!(mask >>= 1)) break; // new
797  }
798 }
799 
800 inline void findComponents3(unsigned int mask, vector<int>& list, int i)
801 {
802  // nejrychlejsi asi bude: ve for cyclu traversovat pres vsechny byty:
803  // 1. voendovanou hodnotu pres 0xFF, aby zustal jen posledni byte
804  // 2. pokud 0, continue
805  // 3. najit prvni nastaveny bit pres lookup table
806  // 4. pricist 8*j k nalezenemu bitu, ulozit do seznamu
807  // 5. odecist 1 << bit, dokud neni *-*nula pokracovat na 3
808  // 6. posunout hodnotu do prava >> 8, pokracovat na 1
809 
810  for (int byte = 0; byte < (int)(sizeof(unsigned int)); byte++)
811  {
812  if (int maskByte = mask & 0xFF)
813  {
814  do
815  {
816  static const int firstBits[256] =
817  {
818  8, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
819  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
820  5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
821  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
822  6, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
823  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
824  5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
825  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
826  7, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
827  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
828  5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
829  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
830  6, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
831  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
832  5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
833  4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0
834  };
835 
836  int bit = firstBits[maskByte];
837  list.push_back(8 * (sizeof(unsigned int)*i + byte) + bit);
838  maskByte -= 1 << bit;
839  }
840  while (maskByte);
841  }
842  mask >>= 8;
843  }
844 }
845 
846 #endif // USOLIDSONLY
847 
848 void UVoxelizer::FindComponentsFastest(unsigned int mask, vector<int>& list, int i)
849 {
850  for (int byte = 0; byte < (int)(sizeof(unsigned int)); byte++)
851  {
852  if (int maskByte = mask & 0xFF)
853  {
854  for (int bit = 0; bit < 8; bit++)
855  {
856  if (maskByte & 1)
857  list.push_back(8 * (sizeof(unsigned int)*i + byte) + bit);
858  if (!(maskByte >>= 1)) break;
859  }
860  }
861  mask >>= 8;
862  }
863 }
864 
865 
866 //______________________________________________________________________________
867 // Method returning the candidates corresponding to the passed point
868 int UVoxelizer::GetCandidatesVoxelArray(const UVector3& point, vector<int>& list, UBits* crossed) const
869 {
870  list.clear();
871 
872  for (int i = 0; i <= 2; ++i)
873  if (point[i] < fBoundaries[i].front() || point[i] >= fBoundaries[i].back())
874  return 0;
875 
876  if (fTotalCandidates == 1)
877  {
878  list.push_back(0);
879  return 1;
880  }
881  else
882  {
883  if (fNPerSlice == 1)
884  {
885  unsigned int mask = 0xFFffFFff;
886  int slice;
887  if (fBoundaries[0].size() > 2)
888  {
889  slice = BinarySearch(fBoundaries[0], point.x);
890  if (!(mask = ((unsigned int*) fBitmasks[0].fAllBits)[slice])) return 0;
891  }
892  if (fBoundaries[1].size() > 2)
893  {
894  slice = BinarySearch(fBoundaries[1], point.y);
895  if (!(mask &= ((unsigned int*) fBitmasks[1].fAllBits)[slice]
896  )) return 0;
897  }
898  if (fBoundaries[2].size() > 2)
899  {
900  slice = BinarySearch(fBoundaries[2], point.z);
901  if (!(mask &= ((unsigned int*) fBitmasks[2].fAllBits)[slice]
902  )) return 0;
903  }
904  if (crossed && (!(mask &= ~((unsigned int*)crossed->fAllBits)[0]))) return 0;
905 
906  FindComponentsFastest(mask, list, 0);
907  }
908  else
909  {
910  unsigned int* masks[3], mask; // masks for X,Y,Z axis
911  for (int i = 0; i <= 2; ++i)
912  {
913  int slice = BinarySearch(fBoundaries[i], point[i]);
914  // if (slice < 0 || slice == fBoundaries[i].size()-1) return 0; // not neccesary anymore
915  masks[i] = ((unsigned int*) fBitmasks[i].fAllBits) + slice * fNPerSlice;
916  }
917  unsigned int* maskCrossed = crossed ? (unsigned int*)crossed->fAllBits : NULL;
918 
919  for (int i = 0 ; i < fNPerSlice; ++i)
920  {
921  // Logic "and" of the masks along the 3 axes x, y, z:
922  // removing "if (!" and ") continue" => slightly slower
923  if (!(mask = masks[0][i])) continue;
924  if (!(mask &= masks[1][i])) continue;
925  if (!(mask &= masks[2][i])) continue;
926  if (maskCrossed && !(mask &= ~maskCrossed[i])) continue;
927 
928  FindComponentsFastest(mask, list, i);
929  }
930  }
931  }
932  return list.size();
933 }
934 
935 
936 
937 
938 
939 int UVoxelizer::GetCandidatesVoxelArray(const vector<int>& voxels, const UBits bitmasks[], vector<int>& list, UBits* crossed) const
940 {
941  list.clear();
942 
943  if (fTotalCandidates == 1)
944  {
945  list.push_back(0);
946  return 1;
947  }
948  else
949  {
950  if (fNPerSlice == 1)
951  {
952  unsigned int mask;
953  if (!(mask = ((unsigned int*) bitmasks[0].fAllBits)[voxels[0]]
954  )) return 0;
955  if (!(mask &= ((unsigned int*) bitmasks[1].fAllBits)[voxels[1]]
956  )) return 0;
957  if (!(mask &= ((unsigned int*) bitmasks[2].fAllBits)[voxels[2]]
958  )) return 0;
959  if (crossed && (!(mask &= ~((unsigned int*)crossed->fAllBits)[0]))) return 0;
960 
961  FindComponentsFastest(mask, list, 0);
962  }
963  else
964  {
965  unsigned int* masks[3], mask; // masks for X,Y,Z axis
966  for (int i = 0; i <= 2; ++i)
967  masks[i] = ((unsigned int*) bitmasks[i].fAllBits) + voxels[i] * fNPerSlice;
968 
969  unsigned int* maskCrossed = crossed ? (unsigned int*)crossed->fAllBits : NULL;
970 
971  for (int i = 0 ; i < fNPerSlice; ++i)
972  {
973  // Logic "and" of the masks along the 3 axes x, y, z:
974  // removing "if (!" and ") continue" => slightly slower
975  if (!(mask = masks[0][i])) continue;
976  if (!(mask &= masks[1][i])) continue;
977  if (!(mask &= masks[2][i])) continue;
978  if (maskCrossed && !(mask &= ~maskCrossed[i])) continue;
979 
980  FindComponentsFastest(mask, list, i);
981  }
982  }
983  }
984  return list.size();
985 }
986 
987 
988 
989 //______________________________________________________________________________
990 // Method returning the candidates corresponding to the passed point
991 int UVoxelizer::GetCandidatesVoxelArray(const vector<int>& voxels, vector<int>& list, UBits* crossed) const
992 {
993  return GetCandidatesVoxelArray(voxels, fBitmasks, list, crossed);
994 }
995 
996 
997 
998 bool UVoxelizer::Contains(const UVector3& point) const
999 {
1000  for (int i = 0; i < 3; ++i)
1001  if (point[i] < fBoundaries[i].front() || point[i] > fBoundaries[i].back())
1002  return false;
1003  return true;
1004 }
1005 
1006 
1007 
1008 double UVoxelizer::DistanceToFirst(const UVector3& point, const UVector3& direction) const
1009 {
1010  UVector3 pointShifted = point - fBoundingBoxCenter;
1011  double shift = fBoundingBox.DistanceToIn(pointShifted, direction);
1012  return shift;
1013 }
1014 
1015 double UVoxelizer::SafetyToBoundingBox(const UVector3& point) const
1016 {
1017  UVector3 pointShifted = point - fBoundingBoxCenter;
1018  double shift = MinDistanceToBox(pointShifted, fBoundingBoxSize);
1019  return shift;
1020 }
1021 
1022 double UVoxelizer::MinDistanceToBox(const UVector3& aPoint, const UVector3& f)
1023 {
1024  // Estimates the isotropic safety from a point outside the current solid to any
1025  // of its surfaces. The algorithm may be accurate or should provide a fast
1026  // underestimate.
1027  double safe, safx, safy, safz;
1028  safe = safx = -f.x + std::abs(aPoint.x);
1029  safy = -f.y + std::abs(aPoint.y);
1030  if (safy > safe) safe = safy;
1031  safz = -f.z + std::abs(aPoint.z);
1032  if (safz > safe) safe = safz;
1033  if (safe < 0.0) return 0.0; // point is inside
1034 // if (!aAccurate) return safe;
1035  double safsq = 0.0;
1036  int count = 0;
1037  if (safx > 0)
1038  {
1039  safsq += safx * safx;
1040  count++;
1041  }
1042  if (safy > 0)
1043  {
1044  safsq += safy * safy;
1045  count++;
1046  }
1047  if (safz > 0)
1048  {
1049  safsq += safz * safz;
1050  count++;
1051  }
1052  if (count == 1) return safe;
1053  return std::sqrt(safsq);
1054 }
1055 
1056 
1057 double UVoxelizer::DistanceToNext(const UVector3& point, const UVector3& direction, vector<int>& curVoxel) const
1058 {
1059  double shift = UUtils::kInfinity;
1060  int cur = 0; // the smallest index, which would be than increased
1061 
1062  for (int i = 0; i <= 2; ++i)
1063  {
1064  // Looking for the next voxels on the considered direction X,Y,Z axis
1065  const vector<double>& boundary = fBoundaries[i];
1066  int index = curVoxel[i];
1067  if (direction[i] >= 1e-10)
1068  {
1069  ++index;
1070 // if (boundary[++index] - point[i] < fTolerance)
1071 
1072 // if (++index >= (int) boundary.size() - 1)
1073 // exit = true;
1074 
1075 // make sure shift would be non-zero => not needed anymore
1076 
1077 // if (++index >= (int) boundary.size()) // the shift is not important, only is important to increase the boundary index, even if it would be zero, than we would in next steps always increase another dimension, so it would continue alright
1078 // continue;
1079  }
1080  else
1081  {
1082  if (direction[i] <= -1e-10)
1083  {
1084 // if (point[i] - boundary[index] < fTolerance) // make sure shift would be non-zero
1085 // if (--index < 0)
1086 // continue;
1087 
1088 // if (--index < 0)
1089 // exit = true;
1090  }
1091  else
1092  continue;
1093  }
1094  double dif = boundary[index] - point[i];
1095  double distance = dif / direction[i];
1096 
1097  if (shift > distance)
1098  {
1099  shift = distance;
1100  cur = i;
1101  }
1102  }
1103 
1104  if (shift != UUtils::kInfinity)
1105  {
1106  // updating current voxel using the index corresponding to the closest voxel boundary on the ray
1107  if (direction[cur] > 0)
1108  {
1109  if (++curVoxel[cur] >= (int) fBoundaries[cur].size() - 1)
1110  shift = UUtils::kInfinity;
1111  }
1112  else
1113  {
1114  if (--curVoxel[cur] < 0)
1115  shift = UUtils::kInfinity;
1116  }
1117  }
1118 
1119  return shift;
1120 }
1121 
1122 
1123 
1124 
1126 {
1127  fMaxVoxels = max;
1128  fReductionRatio.Set(0);
1129 }
1130 
1131 void UVoxelizer::SetMaxVoxels(const UVector3& ratioOfReduction)
1132 {
1133  fMaxVoxels = -1;
1134  fReductionRatio = ratioOfReduction;
1135 }
1136 
1138 {
1139  int size = fEmpty.GetNbytes();
1140  size += fBoxes.capacity() * sizeof(UVoxelBox);
1141  size += sizeof(double) * (fBoundaries[0].capacity() + fBoundaries[1].capacity() + fBoundaries[2].capacity());
1142  size += sizeof(int) * (fCandidatesCounts[0].capacity() + fCandidatesCounts[1].capacity() + fCandidatesCounts[2].capacity());
1143  size += fBitmasks[0].GetNbytes() + fBitmasks[1].GetNbytes() + fBitmasks[2].GetNbytes();
1144 
1145  int csize = fCandidates.size();
1146  for (int i = 0; i < csize; ++i)
1147  size += sizeof(vector<int>) + fCandidates[i].capacity() * sizeof(int);
1148 
1149  return size;
1150 }
1151 
1153 {
1154  fDefaultVoxelsCount = count;
1155 }
1156 
1158 {
1159  return fDefaultVoxelsCount;
1160 }
1161 
1163 
1164 #ifdef USOLIDSONLY
1165 
1166 /*
1167 // THIS METHOD IS NOT NECCESSARY. destructors are called automatically. it is enough to call destructor of
1168 // object which contains the voxel finder
1169 //
1170 inline void DeleteObjects()
1171 {
1172 for (int k = 0; k < 3; k++)
1173 {
1174 std::vector<double> &boundary = fBoundaries[k];
1175 boundary.resize(0);
1176 UBits &bitmask = fBitmasks[k];
1177 bitmask.Clear();
1178 bitmask.SetBitNumber(1, false); // it is here so we can set the maximum
1179 }
1180 }
1181 */
1182 
1183 /*
1184 UVoxelCandidatesIterator::UVoxelCandidatesIterator(const UVoxelizer &f, const UVector3 &point) : curInt(-1), curBit((int) (8*sizeof(unsigned int))), nextAvailable(true)
1185 {
1186  carNodes = f.fBoxes.size();
1187  if (carNodes > 1)
1188  {
1189  n = 1+(carNodes-1)/(8*sizeof(unsigned int));
1190  int sliceX = UVoxelizer::BinarySearch(f.fBoundaries[0], point.x);
1191  int sliceY = UVoxelizer::BinarySearch(f.fBoundaries[1], point.y);
1192  int sliceZ = UVoxelizer::BinarySearch(f.fBoundaries[2], point.z);
1193 
1194  unsigned int *maskX = ((unsigned int *) f.fBitmasks[0].fAllBits) + sliceX*n;
1195  maskXLeft = (sliceX && point.x == f.fBoundaries[0][sliceX]) ? maskX - n : NULL;
1196 
1197  unsigned int *maskY = ((unsigned int *) f.fBitmasks[1].fAllBits) + sliceY*n;
1198  maskYLeft = (sliceY && point.y == f.fBoundaries[1][sliceY]) ? maskY - n : NULL;
1199  unsigned int *maskZ = ((unsigned int *) f.fBitmasks[2].fAllBits) + sliceZ*n;
1200  maskZLeft = (sliceZ && point.z == f.fBoundaries[2][sliceZ]) ? maskZ - n : NULL;
1201 
1202  if (sliceX == -1 || sliceX == f.fBoundaries[0].back()) nextAvailable = false;
1203  if (sliceY == -1 || sliceY == f.fBoundaries[1].back()) nextAvailable = false;
1204  if (sliceZ == -1 || sliceZ == f.fBoundaries[2].back()) nextAvailable = false;
1205  }
1206  else
1207  {
1208  if (point.x < f.fBoundaries[0].front() || point.x > f.fBoundaries[0].back()) nextAvailable = false;
1209  if (point.y < f.fBoundaries[1].front() || point.y > f.fBoundaries[1].back()) nextAvailable = false;
1210  if (point.z < f.fBoundaries[2].front() || point.z > f.fBoundaries[2].back()) nextAvailable = false;
1211  }
1212 }
1213 
1214 int UVoxelCandidatesIterator::Next()
1215 {
1216  if (!nextAvailable) return -1;
1217 
1218  if (carNodes == 1)
1219  {
1220  nextAvailable = false;
1221  return 0;
1222  }
1223  else
1224  {
1225  do
1226  {
1227  while (curBit == (int) (8*sizeof(unsigned int)))
1228  {
1229  if (curInt++ >= n)
1230  {
1231  nextAvailable = false;
1232  return -1;
1233  }
1234  // Logic "and" of the masks along the 3 axes x, y, z:
1235  // removing "if (!" and ") continue" => slightly slower
1236  if (!(mask = maskXLeft ? maskX[curInt] | maskXLeft[curInt] : maskX[curInt]) ||
1237  !(mask &= maskYLeft ? maskY[curInt] | maskYLeft[curInt] : maskY[curInt]) ||
1238  !(mask &= maskZLeft ? maskZ[curInt] | maskZLeft[curInt] : maskZ[curInt]))
1239  {
1240  continue;
1241  }
1242  curBit = 0;
1243  }
1244  int shifted = 1 << curBit; // new
1245  if (mask & shifted)
1246  {
1247  int res = 8*sizeof(unsigned int)*curInt+curBit;
1248  if (!(mask -= shifted))
1249  curBit = (int) (8*sizeof(unsigned int));
1250 
1251  return res;
1252  }
1253  curBit++;
1254  }
1255  while (nextAvailable);
1256  return -1;
1257  }
1258 }
1259 */
1260 
1261 #endif // USOLIDSONLY
void BuildVoxelLimits(std::vector< VUSolid * > &solids, std::vector< UTransform3D > &transforms)
static int fDefaultVoxelsCount
Definition: UVoxelizer.hh:241
UVector3 hlen
Definition: UVoxelizer.hh:38
virtual UGeometryType GetEntityType() const =0
void DisplayBoundaries()
Definition: UVoxelizer.cc:270
int GetVoxelsIndex(int x, int y, int z) const
Definition: UVoxelizer.hh:126
double DistanceToNext(const UVector3 &point, const UVector3 &direction, std::vector< int > &curVoxel) const
Definition: UVoxelizer.cc:1057
vector< UVoxelInfo > & fVoxels
Definition: UVoxelizer.cc:461
static void FindComponentsFastest(unsigned int mask, std::vector< int > &list, int i)
Definition: UVoxelizer.cc:848
unsigned int GetNbytes() const
Definition: UBits.hh:124
static ush mask[]
Definition: csz_inflate.cc:332
G4double z
Definition: TRTMaterials.hh:39
void CreateMiniVoxels(std::vector< double > fBoundaries[], UBits bitmasks[])
Definition: UVoxelizer.cc:647
virtual double Extent(const UVector3)=0
void Voxelize(std::vector< VUSolid * > &solids, std::vector< UTransform3D > &transforms)
void SetMaxVoxels(int max)
Definition: UVoxelizer.cc:1125
void BuildEmpty()
Definition: UVoxelizer.cc:46
bool operator()(int l, int r)
Definition: UVoxelizer.cc:468
static const G4double tolerance
int fMaxVoxels
Definition: UVoxelizer.hh:297
void ResetAllBits(bool value=false)
Definition: UBits.cc:386
static double Tolerance()
Definition: VUSolid.hh:127
bool Contains(const UVector3 &point) const
Definition: UVoxelizer.cc:998
double x
Definition: UVector3.hh:136
Definition: UBox.hh:33
std::vector< UVoxelBox > fBoxes
Definition: UVoxelizer.hh:281
int fTotalCandidates
Definition: UVoxelizer.hh:287
double DistanceToFirst(const UVector3 &point, const UVector3 &direction) const
Definition: UVoxelizer.cc:1008
std::vector< double > fBoundaries[3]
Definition: UVoxelizer.hh:283
static int BinarySearch(const std::vector< T > &vec, T value)
Definition: UVoxelizer.hh:57
UVector3 fTr
Definition: UTransform3D.hh:29
static const double kInfinity
Definition: UUtils.hh:53
Definition: UOrb.hh:28
unsigned char * fAllBits
Definition: UBits.hh:42
void BuildReduceVoxels2(std::vector< double > fBoundaries[], UVector3 reductionRatio)
Definition: UVoxelizer.cc:584
void BuildBoundaries()
Definition: UVoxelizer.cc:207
static double MinDistanceToBox(const UVector3 &aPoint, const UVector3 &f)
Definition: UVoxelizer.cc:1022
void Clear()
Definition: UBits.cc:72
void BuildBoundingBox()
Definition: UVoxelizer.cc:424
std::vector< std::vector< int > > fVoxelBoxesCandidates
Definition: UVoxelizer.hh:245
UVector3 pos
Definition: UVoxelizer.hh:39
int GetCandidatesVoxelArray(const UVector3 &point, std::vector< int > &list, UBits *crossed=NULL) const
Definition: UVoxelizer.cc:868
void BuildBitmasks(std::vector< double > fBoundaries[], UBits bitmasks[])
Definition: UVoxelizer.cc:298
const G4int n
void Set(unsigned int nbits, const char *array)
Definition: UBits.cc:407
long long CountVoxels(std::vector< double > boundaries[]) const
Definition: UVoxelizer.hh:189
void SetBitNumber(unsigned int bitnumber, bool value=true)
Definition: UBits.hh:172
double DistanceToIn(const UVector3 &aPoint, const UVector3 &aDirection, double aPstep=UUtils::kInfinity) const
Definition: UBox.cc:117
int previous
Definition: UVoxelizer.hh:45
UVector3 fReductionRatio
Definition: UVoxelizer.hh:295
int fNPerSlice
Definition: UVoxelizer.hh:279
std::string GetCandidatesAsString(const UBits &bits)
Definition: UVoxelizer.cc:389
UBits fBitmasks[3]
Definition: UVoxelizer.hh:289
G4double total(Particle const *const p1, Particle const *const p2)
void Set(double xx, double yy, double zz)
Definition: UVector3.hh:245
double GetRadialTolerance()
Definition: UOrb.hh:74
T max(const T t1, const T t2)
brief Return the largest of the two arguments
UBox fBoundingBox
Definition: UVoxelizer.hh:292
bool TestBitNumber(unsigned int bitnumber) const
Definition: UBits.hh:199
void GetCandidatesVoxel(std::vector< int > &voxels)
Definition: UVoxelizer.cc:778
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static int GetDefaultVoxelsCount()
Definition: UVoxelizer.cc:1157
std::vector< UVoxelBox > fVoxelBoxes
Definition: UVoxelizer.hh:243
double fTolerance
Definition: UVoxelizer.hh:299
long long fCountOfVoxels
Definition: UVoxelizer.hh:251
int GetBitsPerSlice() const
Definition: UVoxelizer.hh:113
void BuildReduceVoxels(std::vector< double > fBoundaries[], UVector3 reductionRatio)
Definition: UVoxelizer.cc:490
std::vector< int > fCandidatesCounts[3]
Definition: UVoxelizer.hh:285
void CreateSortedBoundary(std::vector< double > &boundaryRaw, int axis)
Definition: UVoxelizer.cc:189
UVoxelComparator(vector< UVoxelInfo > &voxels)
Definition: UVoxelizer.cc:463
void DisplayListNodes()
Definition: UVoxelizer.cc:403
void DisplayVoxelLimits()
Definition: UVoxelizer.cc:173
std::map< int, std::vector< int > > fCandidates
Definition: UVoxelizer.hh:247
virtual void Extent(UVector3 &aMin, UVector3 &aMax) const =0
std::string ToString(int number)
Definition: UUtils.cc:77
double z
Definition: UVector3.hh:138
void SetReductionRatio(int maxVoxels, UVector3 &reductionRatio)
Definition: UVoxelizer.cc:477
double SafetyToBoundingBox(const UVector3 &point) const
Definition: UVoxelizer.cc:1015
int AllocatedMemory()
Definition: UVoxelizer.cc:1137
void ResetBitNumber(unsigned int bitnumber)
Definition: UBits.hh:212
static void SetDefaultVoxelsCount(int count)
Definition: UVoxelizer.cc:1152
UVector3 fBoundingBoxSize
Definition: UVoxelizer.hh:293
double y
Definition: UVector3.hh:137
UBits fEmpty
Definition: UVoxelizer.hh:301
static const G4double pos
Definition: UBits.hh:38
void TransformLimits(UVector3 &min, UVector3 &max, const UTransform3D &transformation)
Definition: UUtils.cc:29
UVector3 fBoundingBoxCenter
Definition: UVoxelizer.hh:291