Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4Voxelizer.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id:$
27 //
28 // --------------------------------------------------------------------
29 // GEANT 4 class header file
30 //
31 // G4Voxelizer implementation
32 //
33 // History:
34 // 19.10.12 Marek Gayer, created
35 // --------------------------------------------------------------------
36 
37 #include <iostream>
38 #include <iomanip>
39 #include <sstream>
40 #include <algorithm>
41 #include <set>
42 
43 #include "G4VSolid.hh"
44 
45 #include "G4Orb.hh"
46 #include "G4Voxelizer.hh"
47 #include "G4SolidStore.hh"
48 #include "Randomize.hh"
49 #include "G4PhysicalConstants.hh"
50 #include "G4GeometryTolerance.hh"
51 #include "G4CSGSolid.hh"
52 #include "G4Orb.hh"
53 #include "G4Types.hh"
54 #include "geomdefs.hh"
55 
56 using namespace std;
57 
59 
60 //______________________________________________________________________________
62  : fBoundingBox("VoxBBox", 1, 1, 1)
63 {
65 
67 
69 
71 }
72 
73 //______________________________________________________________________________
75 {
76 }
77 
78 //______________________________________________________________________________
80 {
81  // by reserving the size of candidates, we would avoid reallocation of
82  // the vector which could cause fragmentation
83  //
84  std::vector<G4int> xyz(3), max(3), candidates(fTotalCandidates);
85  const std::vector<G4int> empty(0);
86 
87  for (G4int i = 0; i <= 2; i++) max[i] = fBoundaries[i].size();
88  unsigned int size = max[0] * max[1] * max[2];
89 
90  fEmpty.Clear();
91  fEmpty.ResetBitNumber(size-1);
92  fEmpty.ResetAllBits(true);
93 
94  for (xyz[2] = 0; xyz[2] < max[2]; ++xyz[2])
95  {
96  for (xyz[1] = 0; xyz[1] < max[1]; ++xyz[1])
97  {
98  for (xyz[0] = 0; xyz[0] < max[0]; ++xyz[0])
99  {
100  if (GetCandidatesVoxelArray(xyz, candidates))
101  {
102  G4int index = GetVoxelsIndex(xyz);
103  fEmpty.SetBitNumber(index, false);
104 
105  // rather than assigning directly with:
106  // "fCandidates[index] = candidates;", in an effort to ensure that
107  // capacity would be just exact, we rather use following 3 lines
108  //
109  std::vector<G4int> &c = (fCandidates[index] = empty);
110  c.reserve(candidates.size());
111  c.assign(candidates.begin(), candidates.end());
112  }
113  }
114  }
115  }
116 #ifdef G4SPECSDEBUG
117  G4cout << "Non-empty voxels count: " << fCandidates.size() << G4endl;
118 #endif
119 }
120 
121 //______________________________________________________________________________
122 void G4Voxelizer::BuildVoxelLimits(std::vector<G4VSolid*>& solids,
123  std::vector<G4Transform3D>& transforms)
124 {
125  G4Rotate3D rot;
126  G4Translate3D transl ;
128 
129  // "BuildVoxelLimits"'s aim is to store the coordinates of the origin as
130  // well as the half lengths related to the bounding box of each node.
131  // These quantities are stored in the array "fBoxes" (6 different values per
132  // node
133  //
134  if (G4int numNodes = solids.size()) // Number of nodes in "multiUnion"
135  {
136  fBoxes.resize(numNodes); // Array which will store the half lengths
137  fNPerSlice = 1 + (fBoxes.size() - 1) / (8 * sizeof(unsigned int));
138 
139  // related to a particular node, but also
140  // the coordinates of its origin
141 
143 
144  for (G4int i = 0; i < numNodes; ++i)
145  {
146  G4VSolid& solid = *solids[i];
147  G4Transform3D transform = transforms[i];
149 
150  solid.BoundingLimits(min, max);
151  if (solid.GetEntityType() == "G4Orb")
152  {
153  G4Orb& orb = *(G4Orb*) &solid;
154  G4ThreeVector orbToleranceVector;
155  G4double tolerance = orb.GetRadialTolerance() / 2.0;
156  orbToleranceVector.set(tolerance,tolerance,tolerance);
157  min -= orbToleranceVector;
158  max += orbToleranceVector;
159  }
160  else
161  {
162  min -= toleranceVector;
163  max += toleranceVector;
164  }
165  TransformLimits(min, max, transform);
166  fBoxes[i].hlen = (max - min) / 2;
167  transform.getDecomposition(scale,rot,transl);
168  fBoxes[i].pos = transl.getTranslation();
169  }
170  fTotalCandidates = fBoxes.size();
171  }
172 }
173 
174 //______________________________________________________________________________
175 void G4Voxelizer::BuildVoxelLimits(std::vector<G4VFacet *> &facets)
176 {
177  // "BuildVoxelLimits"'s aim is to store the coordinates of the origin as well
178  // as the half lengths related to the bounding box of each node.
179  // These quantities are stored in the array "fBoxes" (6 different values per
180  // node.
181 
182  if (G4int numNodes = facets.size()) // Number of nodes
183  {
184  fBoxes.resize(numNodes); // Array which will store the half lengths
185  fNPerSlice = 1+(fBoxes.size()-1)/(8*sizeof(unsigned int));
186 
187  G4ThreeVector toleranceVector(10*fTolerance, 10*fTolerance, 10*fTolerance);
188 
189  for (G4int i = 0; i < numNodes; ++i)
190  {
191  G4VFacet &facet = *facets[i];
193  G4ThreeVector x(1,0,0), y(0,1,0), z(0,0,1);
194  G4ThreeVector extent;
195  max.set (facet.Extent(x), facet.Extent(y), facet.Extent(z));
196  min.set (-facet.Extent(-x), -facet.Extent(-y), -facet.Extent(-z));
197  min -= toleranceVector;
198  max += toleranceVector;
199  G4ThreeVector hlen = (max - min) / 2;
200  fBoxes[i].hlen = hlen;
201  fBoxes[i].pos = min + hlen;
202  }
203  fTotalCandidates = fBoxes.size();
204  }
205 }
206 
207 //______________________________________________________________________________
209 {
210  // "DisplayVoxelLimits" displays the dX, dY, dZ, pX, pY and pZ for each node
211 
212  G4int numNodes = fBoxes.size();
213  G4int oldprec = G4cout.precision(16);
214  for(G4int i = 0; i < numNodes; ++i)
215  {
216  G4cout << setw(10) << setiosflags(ios::fixed) <<
217  " -> Node " << i+1 << ":\n" <<
218  "\t * [x,y,z] = " << fBoxes[i].hlen <<
219  "\t * [x,y,z] = " << fBoxes[i].pos << "\n";
220  }
221  G4cout.precision(oldprec);
222 }
223 
224 //______________________________________________________________________________
225 void G4Voxelizer::CreateSortedBoundary(std::vector<G4double> &boundary,
226  G4int axis)
227 {
228  // "CreateBoundaries"'s aim is to determine the slices induced by the
229  // bounding fBoxes, along each axis. The created boundaries are stored
230  // in the array "boundariesRaw"
231 
232  G4int numNodes = fBoxes.size(); // Number of nodes in structure
233 
234  // Determination of the boundaries along x, y and z axis
235  //
236  for(G4int i = 0 ; i < numNodes; ++i)
237  {
238  // For each node, the boundaries are created by using the array "fBoxes"
239  // built in method "BuildVoxelLimits"
240  //
241  G4double p = fBoxes[i].pos[axis], d = fBoxes[i].hlen[axis];
242 
243  // x boundaries
244  //
245 #ifdef G4SPECSDEBUG
246  G4cout << "Boundary " << p - d << " - " << p + d << G4endl;
247 #endif
248  boundary[2*i] = p - d;
249  boundary[2*i+1] = p + d;
250  }
251  std::sort(boundary.begin(), boundary.end());
252 }
253 
254 //______________________________________________________________________________
256 {
257  // "SortBoundaries" orders the boundaries along each axis (increasing order)
258  // and also does not take into account redundant boundaries, i.e. if two
259  // boundaries are separated by a distance strictly inferior to "tolerance".
260  // The sorted boundaries are respectively stored in:
261  // * boundaries[0..2]
262 
263  // In addition, the number of elements contained in the three latter arrays
264  // are precise thanks to variables: boundariesCountX, boundariesCountY and
265  // boundariesCountZ.
266 
267  if (G4int numNodes = fBoxes.size())
268  {
269  const G4double tolerance = fTolerance / 100.0;
270  // Minimal distance to discriminate two boundaries.
271 
272  std::vector<G4double> sortedBoundary(2*numNodes);
273 
274  G4int considered;
275 
276  for (G4int j = 0; j <= 2; ++j)
277  {
278  CreateSortedBoundary(sortedBoundary, j);
279  std::vector<G4double> &boundary = fBoundaries[j];
280  boundary.clear();
281 
282  considered = 0;
283 
284  for(G4int i = 0 ; i < 2*numNodes; ++i)
285  {
286  G4double newBoundary = sortedBoundary[i];
287 #ifdef G4SPECSDEBUG
288  if (j == 0) G4cout << "Examining " << newBoundary << "..." << G4endl;
289 #endif
290  G4int size = boundary.size();
291  if(!size || std::abs(boundary[size-1] - newBoundary) > tolerance)
292  {
293  considered++;
294  {
295 #ifdef G4SPECSDEBUG
296  if (j == 0) G4cout << "Adding boundary " << newBoundary << "..."
297  << G4endl;
298 #endif
299  boundary.push_back(newBoundary);
300  continue;
301  }
302  }
303  // If two successive boundaries are too close from each other,
304  // only the first one is considered
305  }
306 
307  G4int n = boundary.size();
308  G4int max = 100000;
309  if (n > max/2)
310  {
311  G4int skip = n / (max /2); // n has to be 2x bigger then 50.000.
312  // therefore only from 100.000 reduced
313  std::vector<G4double> reduced;
314  for (G4int i = 0; i < n; ++i)
315  {
316  // 50 ok for 2k, 1000, 2000
317  G4int size = boundary.size();
318  if (i % skip == 0 || i == 0 || i == size - 1)
319  {
320  // this condition of merging boundaries was wrong,
321  // it did not count with right part, which can be
322  // completely ommited and not included in final consideration.
323  // Now should be OK
324  //
325  reduced.push_back(boundary[i]);
326  }
327  }
328  boundary = reduced;
329  }
330  }
331  }
332 }
333 
334 //______________________________________________________________________________
336 {
337  char axis[3] = {'X', 'Y', 'Z'};
338  for (G4int i = 0; i <= 2; ++i)
339  {
340  G4cout << " * " << axis[i] << " axis:" << G4endl << " | ";
342  }
343 }
344 
345 //______________________________________________________________________________
346 void G4Voxelizer::DisplayBoundaries(std::vector<G4double> &boundaries)
347 {
348  // Prints the positions of the boundaries of the slices on the three axes
349 
350  G4int count = boundaries.size();
351  G4int oldprec = G4cout.precision(16);
352  for(G4int i = 0; i < count; ++i)
353  {
354  G4cout << setw(10) << setiosflags(ios::fixed) << boundaries[i];
355  if(i != count-1) G4cout << "-> ";
356  }
357  G4cout << "|" << G4endl << "Number of boundaries: " << count << G4endl;
358  G4cout.precision(oldprec);
359 }
360 
361 //______________________________________________________________________________
362 void G4Voxelizer::BuildBitmasks(std::vector<G4double> boundaries[],
363  G4SurfBits bitmasks[], G4bool countsOnly)
364 {
365  // "BuildListNodes" stores in the bitmasks solids present in each slice
366  // along an axis.
367 
368  G4int numNodes = fBoxes.size();
369  G4int bitsPerSlice = GetBitsPerSlice();
370 
371  for (G4int k = 0; k < 3; ++k)
372  {
373  G4int total = 0;
374  std::vector<G4double> &boundary = boundaries[k];
375  G4int voxelsCount = boundary.size() - 1;
376  G4SurfBits &bitmask = bitmasks[k];
377 
378  if (!countsOnly)
379  {
380  bitmask.Clear();
381 #ifdef G4SPECSDEBUG
382  G4cout << "Allocating bitmask..." << G4endl;
383 #endif
384  bitmask.SetBitNumber(voxelsCount*bitsPerSlice-1, false);
385  // it is here so we can set the maximum number of bits. this line
386  // will rellocate the memory and set all to zero
387  }
388  std::vector<G4int> &candidatesCount = fCandidatesCounts[k];
389  candidatesCount.resize(voxelsCount);
390 
391  for(G4int i = 0 ; i < voxelsCount; ++i) { candidatesCount[i] = 0; }
392 
393  // Loop on the nodes, number of slices per axis
394  //
395  for(G4int j = 0 ; j < numNodes; ++j)
396  {
397  // Determination of the minimum and maximum position along x
398  // of the bounding boxe of each node
399  //
400  G4double p = fBoxes[j].pos[k], d = fBoxes[j].hlen[k];
401 
402  G4double min = p - d; // - localTolerance;
403  G4double max = p + d; // + localTolerance;
404 
405  G4int i = BinarySearch(boundary, min);
406  if (i < 0) { i = 0; }
407 
408  do // Loop checking, 13.08.2015, G.Cosmo
409  {
410  if (!countsOnly)
411  {
412  bitmask.SetBitNumber(i*bitsPerSlice+j);
413  }
414  candidatesCount[i]++;
415  total++;
416  i++;
417  }
418  while (max > boundary[i] && i < voxelsCount);
419  }
420  }
421 #ifdef G4SPECSDEBUG
422  G4cout << "Build list nodes completed." << G4endl;
423 #endif
424 }
425 
426 //______________________________________________________________________________
428 {
429  // Decodes the candidates in mask as G4String.
430 
431  stringstream ss;
432  G4int numNodes = fBoxes.size();
433 
434  for(G4int i=0; i<numNodes; ++i)
435  {
436  if (bits.TestBitNumber(i)) { ss << i+1 << " "; }
437  }
438  return ss.str();
439 }
440 
441 //______________________________________________________________________________
443 {
444  // Prints which solids are present in the slices previously elaborated.
445 
446  char axis[3] = {'X', 'Y', 'Z'};
447  G4int size=8*sizeof(G4int)*fNPerSlice;
448  G4SurfBits bits(size);
449 
450  for (G4int j = 0; j <= 2; ++j)
451  {
452  G4cout << " * " << axis[j] << " axis:" << G4endl;
453  G4int count = fBoundaries[j].size();
454  for(G4int i=0; i < count-1; ++i)
455  {
456  G4cout << " Slice #" << i+1 << ": [" << fBoundaries[j][i]
457  << " ; " << fBoundaries[j][i+1] << "] -> ";
458  bits.set(size,(const char *)fBitmasks[j].fAllBits+i
459  *fNPerSlice*sizeof(G4int));
461  G4cout << "[ " << result.c_str() << "] " << G4endl;
462  }
463  }
464 }
465 
466 //______________________________________________________________________________
468 {
469  G4ThreeVector min(fBoundaries[0].front(),
470  fBoundaries[1].front(),
471  fBoundaries[2].front());
472  G4ThreeVector max(fBoundaries[0].back(),
473  fBoundaries[1].back(),
474  fBoundaries[2].back());
475  BuildBoundingBox(min, max);
476 }
477 
478 //______________________________________________________________________________
480  G4ThreeVector& amax,
481  G4double tolerance)
482 {
483  for (G4int i = 0; i <= 2; ++i)
484  {
485  G4double min = amin[i];
486  G4double max = amax[i];
487  fBoundingBoxSize[i] = (max - min) / 2 + tolerance * 0.5;
488  fBoundingBoxCenter[i] = min + fBoundingBoxSize[i];
489  }
493 }
494 
495 // algorithm -
496 
497 // in order to get balanced voxels, merge should always unite those regions,
498 // where the number of voxels is least the number.
499 // We will keep sorted list (std::set) with all voxels. there will be
500 // comparator function between two voxels, which will tell if voxel is less
501 // by looking at his right neighbor.
502 // First, we will add all the voxels into the tree.
503 // We will be pick the first item in the tree, merging it, adding the right
504 // merged voxel into the a list for future reduction (fBitmasks will be
505 // rebuilded later, therefore they need not to be updated).
506 // The merged voxel need to be added to the tree again, so it's position
507 // would be updated.
508 
509 //______________________________________________________________________________
511  G4ThreeVector &reductionRatio)
512 {
513  G4double maxTotal = (G4double) fCandidatesCounts[0].size()
514  * fCandidatesCounts[1].size() * fCandidatesCounts[2].size();
515 
516  if (maxVoxels > 0 && maxVoxels < maxTotal)
517  {
518  G4double ratio = (G4double) maxVoxels / maxTotal;
519  ratio = std::pow(ratio, 1./3.);
520  if (ratio > 1) { ratio = 1; }
521  reductionRatio.set(ratio,ratio,ratio);
522  }
523 }
524 
525 //______________________________________________________________________________
526 void G4Voxelizer::BuildReduceVoxels(std::vector<G4double> boundaries[],
527  G4ThreeVector reductionRatio)
528 {
529  for (G4int k = 0; k <= 2; ++k)
530  {
531  std::vector<G4int> &candidatesCount = fCandidatesCounts[k];
532  G4int max = candidatesCount.size();
533  std::vector<G4VoxelInfo> voxels(max);
534  G4VoxelComparator comp(voxels);
535  std::set<G4int, G4VoxelComparator> voxelSet(comp);
536  std::vector<G4int> mergings;
537 
538  for (G4int j = 0; j < max; ++j)
539  {
540  G4VoxelInfo &voxel = voxels[j];
541  voxel.count = candidatesCount[j];
542  voxel.previous = j - 1;
543  voxel.next = j + 1;
544  voxels[j] = voxel;
545  }
546 
547  for (G4int j = 0; j < max - 1; ++j) { voxelSet.insert(j); }
548  // we go to size-1 to make sure we will not merge the last element
549 
550  G4double reduction = reductionRatio[k];
551  if (reduction != 0)
552  {
553  G4int count = 0, currentCount;
554  while ((currentCount = voxelSet.size()) > 2)
555  {
556  G4double currentRatio = 1 - (G4double) count / max;
557  if ((currentRatio <= reduction) && (currentCount <= 1000))
558  break;
559  const G4int pos = *voxelSet.begin();
560  mergings.push_back(pos + 1);
561 
562  G4VoxelInfo& voxel = voxels[pos];
563  G4VoxelInfo& nextVoxel = voxels[voxel.next];
564 
565  if (voxelSet.erase(pos) != 1)
566  {
567  ;// k = k;
568  }
569  if (voxel.next != max - 1)
570  if (voxelSet.erase(voxel.next) != 1)
571  {
572  ;// k = k;
573  }
574  if (voxel.previous != -1)
575  if (voxelSet.erase(voxel.previous) != 1)
576  {
577  ;// k = k;
578  }
579  nextVoxel.count += voxel.count;
580  voxel.count = 0;
581  nextVoxel.previous = voxel.previous;
582 
583  if (voxel.next != max - 1)
584  voxelSet.insert(voxel.next);
585 
586  if (voxel.previous != -1)
587  {
588  voxels[voxel.previous].next = voxel.next;
589  voxelSet.insert(voxel.previous);
590  }
591  count++;
592  } // Loop checking, 13.08.2015, G.Cosmo
593  }
594 
595  if (mergings.size())
596  {
597  std::sort(mergings.begin(), mergings.end());
598 
599  const std::vector<G4double>& boundary = boundaries[k];
600  int mergingsSize = mergings.size();
601  vector<G4double> reducedBoundary;
602  G4int skip = mergings[0], i = 0;
603  max = boundary.size();
604  for (G4int j = 0; j < max; ++j)
605  {
606  if (j != skip)
607  {
608  reducedBoundary.push_back(boundary[j]);
609  }
610  else if (++i < mergingsSize)
611  {
612  skip = mergings[i];
613  }
614  }
615  boundaries[k] = reducedBoundary;
616  }
617 /*
618  G4int count = 0;
619  while (true) // Loop checking, 13.08.2015, G.Cosmo
620  {
621  G4double reduction = reductionRatio[k];
622  if (reduction == 0)
623  break;
624  G4int currentCount = voxelSet.size();
625  if (currentCount <= 2)
626  break;
627  G4double currentRatio = 1 - (G4double) count / max;
628  if (currentRatio <= reduction && currentCount <= 1000)
629  break;
630  const G4int pos = *voxelSet.begin();
631  mergings.push_back(pos);
632 
633  G4VoxelInfo &voxel = voxels[pos];
634  G4VoxelInfo &nextVoxel = voxels[voxel.next];
635 
636  voxelSet.erase(pos);
637  if (voxel.next != max - 1) { voxelSet.erase(voxel.next); }
638  if (voxel.previous != -1) { voxelSet.erase(voxel.previous); }
639 
640  nextVoxel.count += voxel.count;
641  voxel.count = 0;
642  nextVoxel.previous = voxel.previous;
643 
644  if (voxel.next != max - 1)
645  voxelSet.insert(voxel.next);
646 
647  if (voxel.previous != -1)
648  {
649  voxels[voxel.previous].next = voxel.next;
650  voxelSet.insert(voxel.previous);
651  }
652  count++;
653  }
654 
655  if (mergings.size())
656  {
657  std::sort(mergings.begin(), mergings.end());
658 
659  std::vector<G4double> &boundary = boundaries[k];
660  std::vector<G4double> reducedBoundary(boundary.size() - mergings.size());
661  G4int skip = mergings[0] + 1, cur = 0, i = 0;
662  max = boundary.size();
663  for (G4int j = 0; j < max; ++j)
664  {
665  if (j != skip)
666  {
667  reducedBoundary[cur++] = boundary[j];
668  }
669  else
670  {
671  if (++i < (G4int)mergings.size()) { skip = mergings[i] + 1; }
672  }
673  }
674  boundaries[k] = reducedBoundary;
675  }
676 */
677  }
678 }
679 
680 //______________________________________________________________________________
681 void G4Voxelizer::BuildReduceVoxels2(std::vector<G4double> boundaries[],
682  G4ThreeVector reductionRatio)
683 {
684  for (G4int k = 0; k <= 2; ++k)
685  {
686  std::vector<G4int> &candidatesCount = fCandidatesCounts[k];
687  G4int max = candidatesCount.size();
688  G4int total = 0;
689  for (G4int i = 0; i < max; ++i) total += candidatesCount[i];
690 
691  G4double reduction = reductionRatio[k];
692  if (reduction == 0)
693  break;
694 
695  G4int destination = (G4int) (reduction * max) + 1;
696  if (destination > 1000) destination = 1000;
697  if (destination < 2) destination = 2;
698  G4double average = ((G4double)total / max) / reduction;
699 
700  std::vector<G4int> mergings;
701 
702  std::vector<G4double> &boundary = boundaries[k];
703  std::vector<G4double> reducedBoundary(destination);
704 
705  G4int sum = 0, cur = 0;
706  for (G4int i = 0; i < max; ++i)
707  {
708  sum += candidatesCount[i];
709  if (sum > average * (cur + 1) || i == 0)
710  {
711  G4double val = boundary[i];
712  reducedBoundary[cur] = val;
713  cur++;
714  if (cur == destination)
715  break;
716  }
717  }
718  reducedBoundary[destination-1] = boundary[max];
719  boundaries[k] = reducedBoundary;
720  }
721 }
722 
723 //______________________________________________________________________________
724 void G4Voxelizer::Voxelize(std::vector<G4VSolid*>& solids,
725  std::vector<G4Transform3D>& transforms)
726 {
727  BuildVoxelLimits(solids, transforms);
728  BuildBoundaries();
731  BuildEmpty(); // this does not work well for multi-union,
732  // actually only makes performance slower,
733  // these are only pre-calculated but not used by multi-union
734 
735  for (G4int i = 0; i < 3; ++i)
736  {
737  fCandidatesCounts[i].resize(0);
738  }
739 }
740 
741 //______________________________________________________________________________
742 void G4Voxelizer::CreateMiniVoxels(std::vector<G4double> boundaries[],
743  G4SurfBits bitmasks[])
744 {
745  std::vector<G4int> voxel(3), maxVoxels(3);
746  for (G4int i = 0; i <= 2; ++i) maxVoxels[i] = boundaries[i].size();
747 
748  G4ThreeVector point;
749  for (voxel[2] = 0; voxel[2] < maxVoxels[2] - 1; ++voxel[2])
750  {
751  for (voxel[1] = 0; voxel[1] < maxVoxels[1] - 1; ++voxel[1])
752  {
753  for (voxel[0] = 0; voxel[0] < maxVoxels[0] - 1; ++voxel[0])
754  {
755  std::vector<G4int> candidates;
756  if (GetCandidatesVoxelArray(voxel, bitmasks, candidates, 0))
757  {
758  // find a box for corresponding non-empty voxel
759  G4VoxelBox box;
760  for (G4int i = 0; i <= 2; ++i)
761  {
762  G4int index = voxel[i];
763  const std::vector<G4double> &boundary = boundaries[i];
764  G4double hlen = 0.5 * (boundary[index+1] - boundary[index]);
765  box.hlen[i] = hlen;
766  box.pos[i] = boundary[index] + hlen;
767  }
768  fVoxelBoxes.push_back(box);
769  std::vector<G4int>(candidates).swap(candidates);
770  fVoxelBoxesCandidates.push_back(candidates);
771  }
772  }
773  }
774  }
775 }
776 
777 //______________________________________________________________________________
778 void G4Voxelizer::Voxelize(std::vector<G4VFacet*>& facets)
779 {
780  G4int maxVoxels = fMaxVoxels;
781  G4ThreeVector reductionRatio = fReductionRatio;
782 
783  G4int size = facets.size();
784  if (size < 10)
785  {
786  for (G4int i = 0; i < (G4int) facets.size(); ++i)
787  {
788  if (facets[i]->GetNumberOfVertices() > 3) size++;
789  }
790  }
791 
792  if ((size >= 10 || maxVoxels > 0) && maxVoxels != 0 && maxVoxels != 1)
793  {
794 #ifdef G4SPECSDEBUG
795  G4cout << "Building voxel limits..." << G4endl;
796 #endif
797 
798  BuildVoxelLimits(facets);
799 
800 #ifdef G4SPECSDEBUG
801  G4cout << "Building boundaries..." << G4endl;
802 #endif
803 
804  BuildBoundaries();
805 
806 #ifdef G4SPECSDEBUG
807  G4cout << "Building bitmasks..." << G4endl;
808 #endif
809 
810  BuildBitmasks(fBoundaries, 0, true);
811 
812  if (maxVoxels < 0 && reductionRatio == G4ThreeVector())
813  {
814  maxVoxels = fTotalCandidates;
815  if (fTotalCandidates > 1000000) maxVoxels = 1000000;
816  }
817 
818  SetReductionRatio(maxVoxels, reductionRatio);
819 
821 
822 #ifdef G4SPECSDEBUG
823  G4cout << "Total number of voxels: " << fCountOfVoxels << G4endl;
824 #endif
825 
826  BuildReduceVoxels2(fBoundaries, reductionRatio);
827 
829 
830 #ifdef G4SPECSDEBUG
831  G4cout << "Total number of voxels after reduction: "
832  << fCountOfVoxels << G4endl;
833 #endif
834 
835 #ifdef G4SPECSDEBUG
836  G4cout << "Building bitmasks..." << G4endl;
837 #endif
838 
840 
841  G4ThreeVector reductionRatioMini;
842 
843  G4SurfBits bitmasksMini[3];
844 
845  // section for building mini voxels
846 
847  std::vector<G4double> miniBoundaries[3];
848 
849  for (G4int i = 0; i <= 2; ++i) { miniBoundaries[i] = fBoundaries[i]; }
850 
851  G4int voxelsCountMini = (fCountOfVoxels >= 1000)
852  ? 100 : fCountOfVoxels / 10;
853 
854  SetReductionRatio(voxelsCountMini, reductionRatioMini);
855 
856 #ifdef G4SPECSDEBUG
857  G4cout << "Building reduced voxels..." << G4endl;
858 #endif
859 
860  BuildReduceVoxels(miniBoundaries, reductionRatioMini);
861 
862 #ifdef G4SPECSDEBUG
863  G4int total = CountVoxels(miniBoundaries);
864  G4cout << "Total number of mini voxels: " << total << G4endl;
865 #endif
866 
867 #ifdef G4SPECSDEBUG
868  G4cout << "Building mini bitmasks..." << G4endl;
869 #endif
870 
871  BuildBitmasks(miniBoundaries, bitmasksMini);
872 
873 #ifdef G4SPECSDEBUG
874  G4cout << "Creating Mini Voxels..." << G4endl;
875 #endif
876 
877  CreateMiniVoxels(miniBoundaries, bitmasksMini);
878 
879 #ifdef G4SPECSDEBUG
880  G4cout << "Building bounding box..." << G4endl;
881 #endif
882 
884 
885 #ifdef G4SPECSDEBUG
886  G4cout << "Building empty..." << G4endl;
887 #endif
888 
889  BuildEmpty();
890 
891 #ifdef G4SPECSDEBUG
892  G4cout << "Deallocating unnecessary fields during runtime..." << G4endl;
893 #endif
894  // deallocate fields unnecessary during runtime
895  //
896  fBoxes.resize(0);
897  for (G4int i = 0; i < 3; ++i)
898  {
899  fCandidatesCounts[i].resize(0);
900  fBitmasks[i].Clear();
901  }
902  }
903 }
904 
905 
906 //______________________________________________________________________________
907 void G4Voxelizer::GetCandidatesVoxel(std::vector<G4int> &voxels)
908 {
909  // "GetCandidates" should compute which solids are possibly contained in
910  // the voxel defined by the three slices characterized by the passed indexes.
911 
912  G4cout << " Candidates in voxel [" << voxels[0] << " ; " << voxels[1]
913  << " ; " << voxels[2] << "]: ";
914  std::vector<G4int> candidates;
915  G4int count = GetCandidatesVoxelArray(voxels, candidates);
916  G4cout << "[ ";
917  for (G4int i = 0; i < count; ++i) G4cout << candidates[i];
918  G4cout << "] " << G4endl;
919 }
920 
921 //______________________________________________________________________________
923  std::vector<G4int> &list, G4int i)
924 {
925  for (G4int byte = 0; byte < (G4int) (sizeof(unsigned int)); byte++)
926  {
927  if (G4int maskByte = mask & 0xFF)
928  {
929  for (G4int bit = 0; bit < 8; bit++)
930  {
931  if (maskByte & 1)
932  { list.push_back(8*(sizeof(unsigned int)*i+ byte) + bit); }
933  if (!(maskByte >>= 1)) break;
934  }
935  }
936  mask >>= 8;
937  }
938 }
939 
940 //______________________________________________________________________________
942  const G4Transform3D& transformation) const
943 {
944  // The goal of this method is to convert the quantities min and max
945  // (representing the bounding box of a given solid in its local frame)
946  // to the main frame, using "transformation"
947 
948  G4ThreeVector vertices[8] = // Detemination of the vertices thanks to
949  { // the extension of each solid:
950  G4ThreeVector(min.x(), min.y(), min.z()), // 1st vertice:
951  G4ThreeVector(min.x(), max.y(), min.z()), // 2nd vertice:
952  G4ThreeVector(max.x(), max.y(), min.z()),
953  G4ThreeVector(max.x(), min.y(), min.z()),
954  G4ThreeVector(min.x(), min.y(), max.z()),
955  G4ThreeVector(min.x(), max.y(), max.z()),
956  G4ThreeVector(max.x(), max.y(), max.z()),
957  G4ThreeVector(max.x(), min.y(), max.z())
958  };
959 
962 
963  // Loop on th vertices
964  G4int limit = sizeof(vertices) / sizeof(G4ThreeVector);
965  for (G4int i = 0 ; i < limit; ++i)
966  {
967  // From local frame to the global one:
968  // Current positions on the three axis:
969  G4ThreeVector current = GetGlobalPoint(transformation, vertices[i]);
970 
971  // If need be, replacement of the min & max values:
972  if (current.x() > max.x()) max.setX(current.x());
973  if (current.x() < min.x()) min.setX(current.x());
974 
975  if (current.y() > max.y()) max.setY(current.y());
976  if (current.y() < min.y()) min.setY(current.y());
977 
978  if (current.z() > max.z()) max.setZ(current.z());
979  if (current.z() < min.z()) min.setZ(current.z());
980  }
981 }
982 
983 //______________________________________________________________________________
985  std::vector<G4int> &list, G4SurfBits *crossed) const
986 {
987  // Method returning the candidates corresponding to the passed point
988 
989  list.clear();
990 
991  for (G4int i = 0; i <= 2; ++i)
992  {
993  if(point[i] < fBoundaries[i].front() || point[i] >= fBoundaries[i].back())
994  return 0;
995  }
996 
997  if (fTotalCandidates == 1)
998  {
999  list.push_back(0);
1000  return 1;
1001  }
1002  else
1003  {
1004  if (fNPerSlice == 1)
1005  {
1006  unsigned int mask = 0xFFffFFff;
1007  G4int slice;
1008  if (fBoundaries[0].size() > 2)
1009  {
1010  slice = BinarySearch(fBoundaries[0], point.x());
1011  if (!(mask = ((unsigned int*) fBitmasks[0].fAllBits)[slice]))
1012  return 0;
1013  }
1014  if (fBoundaries[1].size() > 2)
1015  {
1016  slice = BinarySearch(fBoundaries[1], point.y());
1017  if (!(mask &= ((unsigned int*) fBitmasks[1].fAllBits)[slice]))
1018  return 0;
1019  }
1020  if (fBoundaries[2].size() > 2)
1021  {
1022  slice = BinarySearch(fBoundaries[2], point.z());
1023  if (!(mask &= ((unsigned int*) fBitmasks[2].fAllBits)[slice]))
1024  return 0;
1025  }
1026  if (crossed && (!(mask &= ~((unsigned int*)crossed->fAllBits)[0])))
1027  return 0;
1028 
1029  FindComponentsFastest(mask, list, 0);
1030  }
1031  else
1032  {
1033  unsigned int* masks[3], mask; // masks for X,Y,Z axis
1034  for (G4int i = 0; i <= 2; ++i)
1035  {
1036  G4int slice = BinarySearch(fBoundaries[i], point[i]);
1037  masks[i] = ((unsigned int*) fBitmasks[i].fAllBits)
1038  + slice * fNPerSlice;
1039  }
1040  unsigned int* maskCrossed = crossed
1041  ? (unsigned int*)crossed->fAllBits : 0;
1042 
1043  for (G4int i = 0 ; i < fNPerSlice; ++i)
1044  {
1045  // Logic "and" of the masks along the 3 axes x, y, z:
1046  // removing "if (!" and ") continue" => slightly slower
1047  //
1048  if (!(mask = masks[0][i])) continue;
1049  if (!(mask &= masks[1][i])) continue;
1050  if (!(mask &= masks[2][i])) continue;
1051  if (maskCrossed && !(mask &= ~maskCrossed[i])) continue;
1052 
1053  FindComponentsFastest(mask, list, i);
1054  }
1055  }
1056 /*
1057  if (fNPerSlice == 1)
1058  {
1059  unsigned int mask;
1060  G4int slice = BinarySearch(fBoundaries[0], point.x());
1061  if (!(mask = ((unsigned int *) fBitmasks[0].fAllBits)[slice]
1062  )) return 0;
1063  slice = BinarySearch(fBoundaries[1], point.y());
1064  if (!(mask &= ((unsigned int *) fBitmasks[1].fAllBits)[slice]
1065  )) return 0;
1066  slice = BinarySearch(fBoundaries[2], point.z());
1067  if (!(mask &= ((unsigned int *) fBitmasks[2].fAllBits)[slice]
1068  )) return 0;
1069  if (crossed && (!(mask &= ~((unsigned int *)crossed->fAllBits)[0])))
1070  return 0;
1071 
1072  FindComponentsFastest(mask, list, 0);
1073  }
1074  else
1075  {
1076  unsigned int *masks[3], mask; // masks for X,Y,Z axis
1077  for (G4int i = 0; i <= 2; ++i)
1078  {
1079  G4int slice = BinarySearch(fBoundaries[i], point[i]);
1080  masks[i] = ((unsigned int *) fBitmasks[i].fAllBits) + slice*fNPerSlice;
1081  }
1082  unsigned int *maskCrossed =
1083  crossed ? (unsigned int *)crossed->fAllBits : 0;
1084 
1085  for (G4int i = 0 ; i < fNPerSlice; ++i)
1086  {
1087  // Logic "and" of the masks along the 3 axes x, y, z:
1088  // removing "if (!" and ") continue" => slightly slower
1089  //
1090  if (!(mask = masks[0][i])) continue;
1091  if (!(mask &= masks[1][i])) continue;
1092  if (!(mask &= masks[2][i])) continue;
1093  if (maskCrossed && !(mask &= ~maskCrossed[i])) continue;
1094 
1095  FindComponentsFastest(mask, list, i);
1096  }
1097  }
1098 */
1099  }
1100  return list.size();
1101 }
1102 
1103 //______________________________________________________________________________
1104 G4int
1105 G4Voxelizer::GetCandidatesVoxelArray(const std::vector<G4int> &voxels,
1106  const G4SurfBits bitmasks[],
1107  std::vector<G4int> &list,
1108  G4SurfBits *crossed) const
1109 {
1110  list.clear();
1111 
1112  if (fTotalCandidates == 1)
1113  {
1114  list.push_back(0);
1115  return 1;
1116  }
1117  else
1118  {
1119  if (fNPerSlice == 1)
1120  {
1121  unsigned int mask;
1122  if (!(mask = ((unsigned int *) bitmasks[0].fAllBits)[voxels[0]]))
1123  return 0;
1124  if (!(mask &= ((unsigned int *) bitmasks[1].fAllBits)[voxels[1]]))
1125  return 0;
1126  if (!(mask &= ((unsigned int *) bitmasks[2].fAllBits)[voxels[2]]))
1127  return 0;
1128  if (crossed && (!(mask &= ~((unsigned int *)crossed->fAllBits)[0])))
1129  return 0;
1130 
1131  FindComponentsFastest(mask, list, 0);
1132  }
1133  else
1134  {
1135  unsigned int *masks[3], mask; // masks for X,Y,Z axis
1136  for (G4int i = 0; i <= 2; ++i)
1137  {
1138  masks[i] = ((unsigned int *) bitmasks[i].fAllBits)
1139  + voxels[i]*fNPerSlice;
1140  }
1141  unsigned int *maskCrossed = crossed
1142  ? (unsigned int *)crossed->fAllBits : 0;
1143 
1144  for (G4int i = 0 ; i < fNPerSlice; ++i)
1145  {
1146  // Logic "and" of the masks along the 3 axes x, y, z:
1147  // removing "if (!" and ") continue" => slightly slower
1148  //
1149  if (!(mask = masks[0][i])) continue;
1150  if (!(mask &= masks[1][i])) continue;
1151  if (!(mask &= masks[2][i])) continue;
1152  if (maskCrossed && !(mask &= ~maskCrossed[i])) continue;
1153 
1154  FindComponentsFastest(mask, list, i);
1155  }
1156  }
1157  }
1158  return list.size();
1159 }
1160 
1161 //______________________________________________________________________________
1162 G4int
1163 G4Voxelizer::GetCandidatesVoxelArray(const std::vector<G4int> &voxels,
1164  std::vector<G4int> &list, G4SurfBits *crossed) const
1165 {
1166  // Method returning the candidates corresponding to the passed point
1167 
1168  return GetCandidatesVoxelArray(voxels, fBitmasks, list, crossed);
1169 }
1170 
1171 //______________________________________________________________________________
1173 {
1174  for (G4int i = 0; i < 3; ++i)
1175  {
1176  if (point[i] < fBoundaries[i].front() || point[i] > fBoundaries[i].back())
1177  return false;
1178  }
1179  return true;
1180 }
1181 
1182 //______________________________________________________________________________
1183 G4double
1185  const G4ThreeVector &direction) const
1186 {
1187  G4ThreeVector pointShifted = point - fBoundingBoxCenter;
1188  G4double shift = fBoundingBox.DistanceToIn(pointShifted, direction);
1189  return shift;
1190 }
1191 
1192 //______________________________________________________________________________
1193 G4double
1195 {
1196  G4ThreeVector pointShifted = point - fBoundingBoxCenter;
1197  G4double shift = MinDistanceToBox(pointShifted, fBoundingBoxSize);
1198  return shift;
1199 }
1200 
1201 //______________________________________________________________________________
1202 G4double
1204  const G4ThreeVector &f)
1205 {
1206  // Estimates the isotropic safety from a point outside the current solid to
1207  // any of its surfaces. The algorithm may be accurate or should provide a
1208  // fast underestimate.
1209 
1210  G4double safe, safx, safy, safz;
1211  safe = safx = -f.x() + std::abs(aPoint.x());
1212  safy = -f.y() + std::abs(aPoint.y());
1213  if ( safy > safe ) safe = safy;
1214  safz = -f.z() + std::abs(aPoint.z());
1215  if ( safz > safe ) safe = safz;
1216  if (safe < 0.0) return 0.0; // point is inside
1217 
1218  G4double safsq = 0.0;
1219  G4int count = 0;
1220  if ( safx > 0 ) { safsq += safx*safx; count++; }
1221  if ( safy > 0 ) { safsq += safy*safy; count++; }
1222  if ( safz > 0 ) { safsq += safz*safz; count++; }
1223  if (count == 1) return safe;
1224  return std::sqrt(safsq);
1225 }
1226 
1227 //______________________________________________________________________________
1228 G4double
1230  const G4ThreeVector &direction,
1231  std::vector<G4int> &curVoxel) const
1232 {
1233  G4double shift = kInfinity;
1234 
1235  G4int cur = 0; // the smallest index, which would be than increased
1236  for (int i = 0; i <= 2; ++i)
1237  {
1238  // Looking for the next voxels on the considered direction X,Y,Z axis
1239  //
1240  const std::vector<G4double>& boundary = fBoundaries[i];
1241  G4int index = curVoxel[i];
1242  if (direction[i] >= 1e-10)
1243  {
1244  ++index;
1245  }
1246  else
1247  {
1248  if (direction[i] > -1e-10)
1249  continue;
1250  }
1251  G4double dif = boundary[index] - point[i];
1252  G4double distance = dif / direction[i];
1253 
1254  if (shift > distance)
1255  {
1256  shift = distance;
1257  cur = i;
1258  }
1259  }
1260 
1261  if (shift != kInfinity)
1262  {
1263  // updating current voxel using the index corresponding
1264  // to the closest voxel boundary on the ray
1265 
1266  if (direction[cur] > 0)
1267  {
1268  if (++curVoxel[cur] >= (G4int) fBoundaries[cur].size() - 1)
1269  shift = kInfinity;
1270  }
1271  else
1272  {
1273  if (--curVoxel[cur] < 0)
1274  shift = kInfinity;
1275  }
1276  }
1277 
1278 /*
1279  for (G4int i = 0; i <= 2; ++i)
1280  {
1281  // Looking for the next voxels on the considered direction X,Y,Z axis
1282  //
1283  const std::vector<G4double> &boundary = fBoundaries[i];
1284  G4int cur = curVoxel[i];
1285  if(direction[i] >= 1e-10)
1286  {
1287  if (boundary[++cur] - point[i] < fTolerance) // make sure shift would
1288  if (++cur >= (G4int) boundary.size()) // be non-zero
1289  continue;
1290  }
1291  else
1292  {
1293  if(direction[i] <= -1e-10)
1294  {
1295  if (point[i] - boundary[cur] < fTolerance) // make sure shift would
1296  if (--cur < 0) // be non-zero
1297  continue;
1298  }
1299  else
1300  continue;
1301  }
1302  G4double dif = boundary[cur] - point[i];
1303  G4double distance = dif / direction[i];
1304 
1305  if (shift > distance)
1306  shift = distance;
1307  }
1308 */
1309  return shift;
1310 }
1311 
1312 //______________________________________________________________________________
1313 G4bool
1315  const G4ThreeVector &direction,
1316  std::vector<G4int> &curVoxel) const
1317 {
1318  for (G4int i = 0; i <= 2; ++i)
1319  {
1320  G4int index = curVoxel[i];
1321  const std::vector<G4double> &boundary = fBoundaries[i];
1322 
1323  if (direction[i] > 0)
1324  {
1325  if (point[i] >= boundary[++index])
1326  if (++curVoxel[i] >= (G4int) boundary.size() - 1)
1327  return false;
1328  }
1329  else
1330  {
1331  if (point[i] < boundary[index])
1332  if (--curVoxel[i] < 0)
1333  return false;
1334  }
1335 #ifdef G4SPECSDEBUG
1336  G4int indexOK = BinarySearch(boundary, point[i]);
1337  if (curVoxel[i] != indexOK)
1338  curVoxel[i] = indexOK; // put breakpoint here
1339 #endif
1340  }
1341  return true;
1342 }
1343 
1344 //______________________________________________________________________________
1346 {
1347  fMaxVoxels = max;
1348  fReductionRatio.set(0,0,0);
1349 }
1350 
1351 //______________________________________________________________________________
1352 void G4Voxelizer::SetMaxVoxels(const G4ThreeVector &ratioOfReduction)
1353 {
1354  fMaxVoxels = -1;
1355  fReductionRatio = ratioOfReduction;
1356 }
1357 
1358 //______________________________________________________________________________
1360 {
1361  fDefaultVoxelsCount = count;
1362 }
1363 
1364 //______________________________________________________________________________
1366 {
1367  return fDefaultVoxelsCount;
1368 }
1369 
1370 //______________________________________________________________________________
1372 {
1373  G4int size = fEmpty.GetNbytes();
1374  size += fBoxes.capacity() * sizeof(G4VoxelBox);
1375  size += sizeof(G4double) * (fBoundaries[0].capacity()
1376  + fBoundaries[1].capacity() + fBoundaries[2].capacity());
1377  size += sizeof(G4int) * (fCandidatesCounts[0].capacity()
1378  + fCandidatesCounts[1].capacity() + fCandidatesCounts[2].capacity());
1379  size += fBitmasks[0].GetNbytes() + fBitmasks[1].GetNbytes()
1380  + fBitmasks[2].GetNbytes();
1381 
1382  G4int csize = fCandidates.size();
1383  for (G4int i = 0; i < csize; ++i)
1384  {
1385  size += sizeof(vector<G4int>) + fCandidates[i].capacity() * sizeof(G4int);
1386  }
1387 
1388  return size;
1389 }
Float_t x
Definition: compare.C:6
void SetXHalfLength(G4double dx)
Definition: G4Box.cc:136
void set(double x, double y, double z)
G4double total(Particle const *const p1, Particle const *const p2)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
long long fCountOfVoxels
Definition: G4Voxelizer.hh:237
CLHEP::Hep3Vector getTranslation() const
CLHEP::Hep3Vector G4ThreeVector
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4Box.cc:336
G4ThreeVector fBoundingBoxSize
Definition: G4Voxelizer.hh:257
std::vector< G4double > fBoundaries[3]
Definition: G4Voxelizer.hh:244
void ResetBitNumber(unsigned int bitnumber)
Definition: G4SurfBits.hh:161
void Voxelize(std::vector< G4VSolid * > &solids, std::vector< G4Transform3D > &transforms)
Definition: G4Voxelizer.cc:724
G4bool Contains(const G4ThreeVector &point) const
static const G4double pos
G4SurfBits fEmpty
Definition: G4Voxelizer.hh:265
static const G4double kInfinity
Definition: geomdefs.hh:42
void BuildReduceVoxels(std::vector< G4double > fBoundaries[], G4ThreeVector reductionRatio)
Definition: G4Voxelizer.cc:526
G4ThreeVector GetGlobalPoint(const G4Transform3D &trans, const G4ThreeVector &lpoint) const
G4ThreeVector fReductionRatio
Definition: G4Voxelizer.hh:259
G4int GetCandidatesVoxelArray(const G4ThreeVector &point, std::vector< G4int > &list, G4SurfBits *crossed=0) const
Definition: G4Voxelizer.cc:984
static ush mask[]
Definition: csz_inflate.cc:317
virtual void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:625
#define G4endl
Definition: G4ios.hh:61
Float_t y
Definition: compare.C:6
const char * p
Definition: xmltok.h:285
std::map< G4int, std::vector< G4int > > fCandidates
Definition: G4Voxelizer.hh:233
unsigned int GetNbytes() const
Definition: G4SurfBits.hh:102
void BuildVoxelLimits(std::vector< G4VSolid * > &solids, std::vector< G4Transform3D > &transforms)
Definition: G4Voxelizer.cc:122
void SetZHalfLength(G4double dz)
Definition: G4Box.cc:183
void DisplayVoxelLimits() const
Definition: G4Voxelizer.cc:208
Double_t z
double z() const
G4double fTolerance
Definition: G4Voxelizer.hh:263
for(int i=0;i< 401;i++)
void DisplayBoundaries()
Definition: G4Voxelizer.cc:335
unsigned char * fAllBits
Definition: G4SurfBits.hh:113
#define G4ThreadLocal
Definition: tls.hh:69
G4int previous
Definition: G4Voxelizer.hh:65
G4ThreeVector hlen
Definition: G4Voxelizer.hh:58
void setX(double)
void Clear()
Definition: G4SurfBits.cc:92
void BuildBitmasks(std::vector< G4double > fBoundaries[], G4SurfBits bitmasks[], G4bool countsOnly=false)
Definition: G4Voxelizer.cc:362
std::vector< G4int > fCandidatesCounts[3]
Definition: G4Voxelizer.hh:247
void SetMaxVoxels(G4int max)
void DisplayListNodes() const
Definition: G4Voxelizer.cc:442
void SetReductionRatio(G4int maxVoxels, G4ThreeVector &reductionRatio)
Definition: G4Voxelizer.cc:510
G4double DistanceToFirst(const G4ThreeVector &point, const G4ThreeVector &direction) const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
Double_t scale
std::vector< std::vector< G4int > > fVoxelBoxesCandidates
Definition: G4Voxelizer.hh:232
void setZ(double)
G4String GetCandidatesAsString(const G4SurfBits &bits) const
Definition: G4Voxelizer.cc:427
std::vector< G4VoxelBox > fVoxelBoxes
Definition: G4Voxelizer.hh:231
G4int fNPerSlice
Definition: G4Voxelizer.hh:239
G4SurfBits fBitmasks[3]
Definition: G4Voxelizer.hh:251
void GetCandidatesVoxel(std::vector< G4int > &voxels)
Definition: G4Voxelizer.cc:907
Float_t d
void ResetAllBits(G4bool value=false)
Definition: G4SurfBits.cc:154
Definition: G4Orb.hh:62
static G4double MinDistanceToBox(const G4ThreeVector &aPoint, const G4ThreeVector &f)
static void FindComponentsFastest(unsigned int mask, std::vector< G4int > &list, G4int i)
Definition: G4Voxelizer.cc:922
G4ThreeVector pos
Definition: G4Voxelizer.hh:59
void set(unsigned int nbits, const char *array)
Definition: G4SurfBits.cc:174
void CreateSortedBoundary(std::vector< G4double > &boundaryRaw, G4int axis)
Definition: G4Voxelizer.cc:225
G4bool TestBitNumber(unsigned int bitnumber) const
Definition: G4SurfBits.hh:148
G4double G4ParticleHPJENDLHEData::G4double result
G4int fTotalCandidates
Definition: G4Voxelizer.hh:249
G4int AllocatedMemory()
long long CountVoxels(std::vector< G4double > boundaries[]) const
void SetYHalfLength(G4double dy)
Definition: G4Box.cc:160
static G4SolidStore * GetInstance()
G4double DistanceToNext(const G4ThreeVector &point, const G4ThreeVector &direction, std::vector< G4int > &curVoxel) const
void BuildEmpty()
Definition: G4Voxelizer.cc:79
G4Box fBoundingBox
Definition: G4Voxelizer.hh:255
int G4int
Definition: G4Types.hh:78
std::vector< G4VoxelBox > fBoxes
Definition: G4Voxelizer.hh:241
G4int fMaxVoxels
Definition: G4Voxelizer.hh:261
G4GLOB_DLL std::ostream G4cout
double x() const
G4int GetBitsPerSlice() const
G4ThreeVector fBoundingBoxCenter
Definition: G4Voxelizer.hh:253
G4bool UpdateCurrentVoxel(const G4ThreeVector &point, const G4ThreeVector &direction, std::vector< G4int > &curVoxel) const
static G4ThreadLocal G4int fDefaultVoxelsCount
Definition: G4Voxelizer.hh:229
static G4GeometryTolerance * GetInstance()
Char_t n[5]
void BuildBoundingBox()
Definition: G4Voxelizer.cc:467
void SetBitNumber(unsigned int bitnumber, G4bool value=true)
Definition: G4SurfBits.hh:123
G4double DistanceToBoundingBox(const G4ThreeVector &point) const
void BuildReduceVoxels2(std::vector< G4double > fBoundaries[], G4ThreeVector reductionRatio)
Definition: G4Voxelizer.cc:681
void getDecomposition(Scale3D &scale, Rotate3D &rotation, Translate3D &translation) const
Definition: Transform3D.cc:174
static void DeRegister(G4VSolid *pSolid)
double y() const
void TransformLimits(G4ThreeVector &min, G4ThreeVector &max, const G4Transform3D &transformation) const
Definition: G4Voxelizer.cc:941
virtual G4GeometryType GetEntityType() const =0
G4int GetVoxelsIndex(G4int x, G4int y, G4int z) const
void CreateMiniVoxels(std::vector< G4double > fBoundaries[], G4SurfBits bitmasks[])
Definition: G4Voxelizer.cc:742
static G4int BinarySearch(const std::vector< T > &vec, T value)
void BuildBoundaries()
Definition: G4Voxelizer.cc:255
virtual G4double Extent(const G4ThreeVector)=0
G4double GetSurfaceTolerance() const
void setY(double)
static G4int GetDefaultVoxelsCount()
static void SetDefaultVoxelsCount(G4int count)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double GetRadialTolerance() const