source: trunk/source/geometry/management/src/G4SmartVoxelHeader.cc@ 1192

Last change on this file since 1192 was 1058, checked in by garnier, 17 years ago

file release beta
File size: 44.1 KB
Line 
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//
27// $Id: G4SmartVoxelHeader.cc,v 1.31 2008/01/24 15:46:47 gcosmo Exp $
28// GEANT4 tag $Name: geant4-09-02-ref-02 $
29//
30//
31// class G4SmartVoxelHeader
32//
33// Implementation
34//
35// Define G4GEOMETRY_VOXELDEBUG for debugging information on G4cout
36//
37// History:
38// 29.04.02 Use 3D voxelisation for non consuming replication - G.C.
39// 18.04.01 Migrated to STL vector - G.C.
40// 12.02.99 Introduction of new quality/smartless: max for (slices/candid) S.G.
41// 11.02.99 Voxels at lower levels are now built for collapsed slices S.G.
42// 21.07.95 Full implementation, supporting non divided physical volumes
43// 14.07.95 Initial version - stubb definitions only
44// --------------------------------------------------------------------
45
46#include "G4SmartVoxelHeader.hh"
47
48#include "G4ios.hh"
49
50#include "G4LogicalVolume.hh"
51#include "G4VPhysicalVolume.hh"
52#include "G4VoxelLimits.hh"
53
54#include "voxeldefs.hh"
55#include "G4AffineTransform.hh"
56#include "G4VSolid.hh"
57#include "G4VPVParameterisation.hh"
58
59// ***************************************************************************
60// Constructor for topmost header, to begin voxel construction at a
61// given logical volume.
62// Constructs target List of volumes, calls "Build and refine" constructor.
63// Assumes all daughters represent single volumes (ie. no divisions
64// or parametric)
65// ***************************************************************************
66//
67G4SmartVoxelHeader::G4SmartVoxelHeader(G4LogicalVolume* pVolume,
68 G4int pSlice)
69 : fminEquivalent(pSlice),
70 fmaxEquivalent(pSlice),
71 fparamAxis(kUndefined)
72{
73 G4int nDaughters = pVolume->GetNoDaughters();
74 G4VoxelLimits limits; // Create `unlimited' limits object
75
76 // Determine whether daughter is replicated
77 //
78 if ((nDaughters!=1) || (!pVolume->GetDaughter(0)->IsReplicated()))
79 {
80 // Daughter not replicated => conventional voxel Build
81 // where each daughters extents are computed
82 //
83 BuildVoxels(pVolume);
84 }
85 else
86 {
87 // Single replicated daughter
88 //
89 BuildReplicaVoxels(pVolume);
90 }
91}
92
93// ***************************************************************************
94// Protected constructor:
95// builds and refines voxels between specified limits, considering only
96// the physical volumes numbered `pCandidates'. `pSlice' is used to set max
97// and min equivalent slice nos for the header - they apply to the level
98// of the header, not its nodes.
99// ***************************************************************************
100//
101G4SmartVoxelHeader::G4SmartVoxelHeader(G4LogicalVolume* pVolume,
102 const G4VoxelLimits& pLimits,
103 const G4VolumeNosVector* pCandidates,
104 G4int pSlice)
105 : fminEquivalent(pSlice),
106 fmaxEquivalent(pSlice),
107 fparamAxis(kUndefined)
108{
109#ifdef G4GEOMETRY_VOXELDEBUG
110 G4cout << "**** G4SmartVoxelHeader::G4SmartVoxelHeader" << G4endl
111 << " Limits " << pLimits << G4endl
112 << " Candidate #s = " ;
113 for (size_t i=0;i<pCandidates->size();i++)
114 {
115 G4cout << (*pCandidates)[i] << " ";
116 }
117 G4cout << G4endl;
118#endif
119
120 BuildVoxelsWithinLimits(pVolume,pLimits,pCandidates);
121}
122
123// ***************************************************************************
124// Destructor:
125// deletes all proxies and underlying objects.
126// ***************************************************************************
127//
128G4SmartVoxelHeader::~G4SmartVoxelHeader()
129{
130 // Manually destroy underlying nodes/headers
131 // Delete collected headers and nodes once only
132 //
133 G4int node, proxy, maxNode=fslices.size();
134 G4SmartVoxelProxy *lastProxy=0;
135 G4SmartVoxelNode *dyingNode, *lastNode=0;
136 G4SmartVoxelHeader *dyingHeader, *lastHeader=0;
137
138 for (node=0; node<maxNode; node++)
139 {
140 if (fslices[node]->IsHeader())
141 {
142 dyingHeader = fslices[node]->GetHeader();
143 if (lastHeader!=dyingHeader)
144 {
145 lastHeader = dyingHeader;
146 lastNode = 0;
147 delete dyingHeader;
148 }
149 }
150 else
151 {
152 dyingNode = fslices[node]->GetNode();
153 if (dyingNode!=lastNode)
154 {
155 lastNode=dyingNode;
156 lastHeader=0;
157 delete dyingNode;
158 }
159 }
160 }
161 // Delete proxies
162 //
163 for (proxy=0; proxy<maxNode; proxy++)
164 {
165 if (fslices[proxy]!=lastProxy)
166 {
167 lastProxy = fslices[proxy];
168 delete lastProxy;
169 }
170 }
171 // Don't need to clear slices
172 // fslices.clear();
173}
174
175// ***************************************************************************
176// Equality operator: returns true if contents are equivalent.
177// Implies a deep search through contained nodes/header.
178// Compares headers' axes,sizes,extents. Returns false if different.
179// For each contained proxy, determines whether node/header, compares and
180// returns if different. Compares and returns if proxied nodes/headers
181// are different.
182// ***************************************************************************
183//
184G4bool G4SmartVoxelHeader::operator == (const G4SmartVoxelHeader& pHead) const
185{
186 if ( (GetAxis() == pHead.GetAxis())
187 && (GetNoSlices() == pHead.GetNoSlices())
188 && (GetMinExtent() == pHead.GetMinExtent())
189 && (GetMaxExtent() == pHead.GetMaxExtent()) )
190 {
191 G4int node, maxNode;
192 G4SmartVoxelProxy *leftProxy, *rightProxy;
193 G4SmartVoxelHeader *leftHeader, *rightHeader;
194 G4SmartVoxelNode *leftNode, *rightNode;
195
196 maxNode=GetNoSlices();
197 for (node=0; node<maxNode; node++)
198 {
199 leftProxy = GetSlice(node);
200 rightProxy = pHead.GetSlice(node);
201 if (leftProxy->IsHeader())
202 {
203 if (rightProxy->IsNode())
204 {
205 return false;
206 }
207 else
208 {
209 leftHeader = leftProxy->GetHeader();
210 rightHeader = rightProxy->GetHeader();
211 if (!(*leftHeader==*rightHeader))
212 {
213 return false;
214 }
215 }
216 }
217 else
218 {
219 if (rightProxy->IsHeader())
220 {
221 return false;
222 }
223 else
224 {
225 leftNode = leftProxy->GetNode();
226 rightNode = rightProxy->GetNode();
227 if (!(*leftNode==*rightNode))
228 {
229 return false;
230 }
231 }
232 }
233 }
234 return true;
235 }
236 else
237 {
238 return false;
239 }
240}
241
242// ***************************************************************************
243// Builds voxels for daughters specified volume, in NON-REPLICATED case
244// o Create List of target volume nos (all daughters; 0->noDaughters-1)
245// o BuildWithinLimits does Build & also determines mother dimensions.
246// ***************************************************************************
247//
248void G4SmartVoxelHeader::BuildVoxels(G4LogicalVolume* pVolume)
249{
250 G4VoxelLimits limits; // Create `unlimited' limits object
251 G4int nDaughters = pVolume->GetNoDaughters();
252
253 G4VolumeNosVector targetList;
254 targetList.reserve(nDaughters);
255 for (G4int i=0; i<nDaughters; i++)
256 {
257 targetList.push_back(i);
258 }
259 BuildVoxelsWithinLimits(pVolume, limits, &targetList);
260}
261
262// ***************************************************************************
263// Builds voxels for specified volume containing a single replicated volume.
264// If axis is not specified (i.e. "kUndefined"), 3D voxelisation is applied,
265// and the best axis is determined according to heuristics as for placements.
266// ***************************************************************************
267//
268void G4SmartVoxelHeader::BuildReplicaVoxels(G4LogicalVolume* pVolume)
269{
270 G4VPhysicalVolume *pDaughter=0;
271
272 // Replication data
273 //
274 EAxis axis;
275 G4int nReplicas;
276 G4double width,offset;
277 G4bool consuming;
278
279 // Consistency check: pVolume should contain single replicated volume
280 //
281 if ( (pVolume->GetNoDaughters()==1)
282 && (pVolume->GetDaughter(0)->IsReplicated()==true) )
283 {
284 // Obtain replication data
285 //
286 pDaughter=pVolume->GetDaughter(0);
287 pDaughter->GetReplicationData(axis,nReplicas,width,offset,consuming);
288 fparamAxis = axis;
289 if ( consuming==false )
290 {
291 G4VoxelLimits limits; // Create `unlimited' limits object
292 G4VolumeNosVector targetList;
293 targetList.reserve(nReplicas);
294 for (G4int i=0; i<nReplicas; i++)
295 {
296 targetList.push_back(i);
297 }
298 if (axis != kUndefined)
299 {
300 // Apply voxelisation along the specified axis only
301
302 G4ProxyVector* pSlices=BuildNodes(pVolume,limits,&targetList,axis);
303 faxis = axis;
304 fslices = *pSlices;
305 delete pSlices;
306
307 // Calculate and set min and max extents given our axis
308 //
309 const G4AffineTransform origin;
310 pVolume->GetSolid()->CalculateExtent(faxis, limits, origin,
311 fminExtent, fmaxExtent);
312 // Calculate equivalent nos
313 //
314 BuildEquivalentSliceNos();
315 CollectEquivalentNodes(); // Collect common nodes
316 }
317 else
318 {
319 // Build voxels similarly as for normal placements considering
320 // all three cartesian axes.
321
322 BuildVoxelsWithinLimits(pVolume, limits, &targetList);
323 }
324 }
325 else
326 {
327 // Replication is consuming -> Build voxels directly
328 //
329 // o Cartesian axes - range is -width*nREplicas/2 to +width*nREplicas/2
330 // nReplicas replications result
331 // o Radial axis (rho) = range is 0 to width*nReplicas
332 // nReplicas replications result
333 // o Phi axi - range is offset to offset+width*nReplicas radians
334 //
335 // Equivalent slices no computation & collection not required - all
336 // slices are different
337 //
338 switch (axis)
339 {
340 case kXAxis:
341 case kYAxis:
342 case kZAxis:
343 fminExtent = -width*nReplicas*0.5;
344 fmaxExtent = width*nReplicas*0.5;
345 break;
346 case kRho:
347 fminExtent = offset;
348 fmaxExtent = width*nReplicas+offset;
349 break;
350 case kPhi:
351 fminExtent = offset;
352 fmaxExtent = offset+width*nReplicas;
353 break;
354 default:
355 G4cerr << "ERROR - G4SmartVoxelHeader::BuildReplicaVoxels()"
356 << G4endl
357 << " Illegal axis !" << G4endl;
358 G4Exception("G4SmartVoxelHeader::BuildReplicaVoxels()", "FatalError",
359 FatalException, "Illegal axis.");
360 break;
361 }
362 faxis = axis; // Set axis
363 BuildConsumedNodes(nReplicas);
364 if ( (axis==kXAxis) || (axis==kYAxis) || (axis==kZAxis) )
365 {
366 // Sanity check on extent
367 //
368 G4double min, max;
369 G4VoxelLimits limits;
370 G4AffineTransform origin;
371 pVolume->GetSolid()->CalculateExtent(axis, limits, origin, min, max);
372 if ( (std::fabs((min-fminExtent)/fminExtent) +
373 std::fabs((max-fmaxExtent)/fmaxExtent)) > 0.05)
374 {
375 G4cerr << "ERROR - G4SmartVoxelHeader::BuildReplicaVoxels()"
376 << G4endl
377 << " Replicated geometry, logical volume: "
378 << pVolume->GetName() << G4endl;
379 G4Exception("G4SmartVoxelHeader::BuildReplicaVoxels", "FatalError",
380 FatalException, "Sanity check: wrong solid extent.");
381 }
382 }
383 }
384 }
385 else
386 {
387 G4cerr << "ERROR - G4SmartVoxelHeader::BuildReplicaVoxels()"
388 << G4endl
389 << " There must be a single replicated volume !" << G4endl;
390 G4Exception("G4SmartVoxelHeader::BuildReplicaVoxels", "InvalidSetup",
391 FatalException, "Only one replicated daughter is allowed !");
392 }
393}
394
395// ***************************************************************************
396// Builds `consumed nodes': nReplicas nodes each containing one replication,
397// numbered in sequence 0->nReplicas-1
398// o Modifies fslices `in place'
399// o faxis,fminExtent,fmaxExtent NOT modified.
400// ***************************************************************************
401//
402void G4SmartVoxelHeader::BuildConsumedNodes(G4int nReplicas)
403{
404 G4int nNode, nVol;
405 G4SmartVoxelNode *pNode;
406 G4SmartVoxelProxy *pProxyNode;
407
408 // Create and fill nodes in temporary G4NodeVector (on stack)
409 //
410 G4NodeVector nodeList;
411 nodeList.reserve(nReplicas);
412 for (nNode=0; nNode<nReplicas; nNode++)
413 {
414 pNode=new G4SmartVoxelNode(nNode);
415 if (!pNode)
416 {
417 G4cerr << "ERROR - G4SmartVoxelHeader::BuildConsumedNodes()" << G4endl
418 << " Node allocation failed." << G4endl;
419 G4Exception("G4SmartVoxelHeader::BuildConsumedNodes()", "FatalError",
420 FatalException, "Node allocation error.");
421 }
422 nodeList.push_back(pNode);
423 }
424 for (nVol=0; nVol<nReplicas; nVol++)
425 {
426 nodeList[nVol]->Insert(nVol); // Insert replication of number
427 } // identical to voxel number
428
429 // Create & fill proxy List `in place' by modifying instance data fslices
430 //
431 fslices.clear();
432 for (nNode=0; nNode<nReplicas; nNode++)
433 {
434 pProxyNode = new G4SmartVoxelProxy(nodeList[nNode]);
435 if (!pProxyNode)
436 {
437 G4cerr << "ERROR - G4SmartVoxelHeader::BuildConsumedNodes()" << G4endl
438 << " Proxy node allocation failed." << G4endl;
439 G4Exception("G4SmartVoxelHeader::BuildConsumedNodes()", "FatalError",
440 FatalException, "Proxy node allocation error.");
441 }
442 fslices.push_back(pProxyNode);
443 }
444}
445
446// ***************************************************************************
447// Builds and refines voxels between specified limits, considering only
448// the physical volumes numbered `pCandidates'.
449// o Chooses axis
450// o Determines min and max extents (of mother solid) within limits.
451// ***************************************************************************
452//
453void
454G4SmartVoxelHeader::BuildVoxelsWithinLimits(G4LogicalVolume* pVolume,
455 G4VoxelLimits pLimits,
456 const G4VolumeNosVector* pCandidates)
457{
458 // Choose best axis for slicing by:
459 // 1. Trying all unlimited cartesian axes
460 // 2. Select axis which gives greatest no slices
461
462 G4ProxyVector *pGoodSlices=0, *pTestSlices, *tmpSlices;
463 G4double goodSliceScore=kInfinity, testSliceScore;
464 EAxis goodSliceAxis = kXAxis;
465 EAxis testAxis = kXAxis;
466 G4int node, maxNode, iaxis;
467 G4VoxelLimits noLimits;
468
469 // Try all non-limited cartesian axes
470 //
471 for (iaxis=0; iaxis<3; iaxis++)
472 {
473 switch(iaxis)
474 {
475 case 0:
476 testAxis = kXAxis;
477 break;
478 case 1:
479 testAxis = kYAxis;
480 break;
481 case 2:
482 testAxis = kZAxis;
483 break;
484 }
485 if (!pLimits.IsLimited(testAxis))
486 {
487 pTestSlices = BuildNodes(pVolume,pLimits,pCandidates,testAxis);
488 testSliceScore = CalculateQuality(pTestSlices);
489 if ( (!pGoodSlices) || (testSliceScore<goodSliceScore) )
490 {
491 goodSliceAxis = testAxis;
492 goodSliceScore = testSliceScore;
493 tmpSlices = pGoodSlices;
494 pGoodSlices = pTestSlices;
495 pTestSlices = tmpSlices;
496 }
497 if (pTestSlices)
498 {
499 // Destroy pTestSlices and all its contents
500 //
501 maxNode=pTestSlices->size();
502 for (node=0; node<maxNode; node++)
503 {
504 delete (*pTestSlices)[node]->GetNode();
505 }
506 G4SmartVoxelProxy* tmpProx;
507 while (pTestSlices->size()>0)
508 {
509 tmpProx = pTestSlices->back();
510 pTestSlices->pop_back();
511 for (G4ProxyVector::iterator i=pTestSlices->begin();
512 i!=pTestSlices->end(); i++)
513 {
514 if (*i==tmpProx)
515 {
516 pTestSlices->erase(i); i--;
517 }
518 }
519 if ( tmpProx ) { delete tmpProx; }
520 }
521 delete pTestSlices;
522 }
523 }
524 }
525 // Check for error case.. when limits already 3d,
526 // so cannot select a new axis
527 //
528 if (!pGoodSlices)
529 {
530 G4cerr << "ERROR - G4SmartVoxelHeader::BuildVoxelsWithinLimits()"
531 << G4endl
532 << " Illegal limits: only 3 dimensions allowed." << G4endl;
533 G4Exception("G4SmartVoxelHeader::BuildVoxelsWithinLimits()",
534 "InvalidSetup", FatalException,
535 "Cannot select more than 3 axis for optimisation.");
536 }
537
538 //
539 // We have selected pGoodSlices, with a score testSliceScore
540 //
541
542 // Store chosen axis, slice ptr
543 //
544 fslices=*pGoodSlices; // Set slice information, copy ptrs in collection
545 delete pGoodSlices; // Destroy slices vector, but not contained
546 // proxies or nodes
547 faxis=goodSliceAxis;
548
549#ifdef G4GEOMETRY_VOXELDEBUG
550 G4cout << G4endl << " Selected axis = " << faxis << G4endl;
551 for (size_t islice=0; islice<fslices.size(); islice++)
552 {
553 G4cout << " Node #" << islice << " = {";
554 for (G4int j=0; j<fslices[islice]->GetNode()->GetNoContained(); j++)
555 {
556 G4cout << " " << fslices[islice]->GetNode()->GetVolume(j);
557 }
558 G4cout << " }" << G4endl;
559 }
560 G4cout << G4endl;
561#endif
562
563 // Calculate and set min and max extents given our axis
564 //
565 G4VSolid* outerSolid = pVolume->GetSolid();
566 const G4AffineTransform origin;
567 if(!outerSolid->CalculateExtent(faxis,pLimits,origin,fminExtent,fmaxExtent))
568 {
569 outerSolid->CalculateExtent(faxis,noLimits,origin,fminExtent,fmaxExtent);
570 }
571
572 // Calculate equivalent nos
573 //
574 BuildEquivalentSliceNos();
575 CollectEquivalentNodes(); // Collect common nodes
576 RefineNodes(pVolume,pLimits); // Refine nodes creating headers
577
578 // No common headers can exist because collapsed by construction
579}
580
581// ***************************************************************************
582// Calculates and stores the minimum and maximum equivalent neighbour
583// values for all slices at our level.
584//
585// Precondition: all slices are nodes.
586// For each potential start of a group of equivalent nodes:
587// o searches forwards in fslices to find group end
588// o loops from start to end setting start and end slices.
589// ***************************************************************************
590//
591void G4SmartVoxelHeader::BuildEquivalentSliceNos()
592{
593 G4int sliceNo, minNo, maxNo, equivNo;
594 G4int maxNode = fslices.size();
595 G4SmartVoxelNode *startNode, *sampleNode;
596 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
597 {
598 minNo = sliceNo;
599
600 // Get first node (see preconditions - will throw exception if a header)
601 //
602 startNode = fslices[minNo]->GetNode();
603
604 // Find max equivalent
605 //
606 for (equivNo=minNo+1; equivNo<maxNode; equivNo++)
607 {
608 sampleNode = fslices[equivNo]->GetNode();
609 if (!((*startNode) == (*sampleNode))) { break; }
610 }
611 maxNo = equivNo-1;
612 if (maxNo != minNo)
613 {
614 // Set min and max nos
615 //
616 for (equivNo=minNo; equivNo<=maxNo; equivNo++)
617 {
618 sampleNode = fslices[equivNo]->GetNode();
619 sampleNode->SetMinEquivalentSliceNo(minNo);
620 sampleNode->SetMaxEquivalentSliceNo(maxNo);
621 }
622 // Advance outer loop to end of equivalent group
623 //
624 sliceNo = maxNo;
625 }
626 }
627}
628
629// ***************************************************************************
630// Collects common nodes at our level, deleting all but one to save
631// memory, and adjusting stored slice pointers appropriately.
632//
633// Preconditions:
634// o the slices have not previously be "collected"
635// o all of the slices are nodes.
636// ***************************************************************************
637//
638void G4SmartVoxelHeader::CollectEquivalentNodes()
639{
640 G4int sliceNo, maxNo, equivNo;
641 G4int maxNode=fslices.size();
642 G4SmartVoxelNode *equivNode;
643 G4SmartVoxelProxy *equivProxy;
644
645 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
646 {
647 equivProxy=fslices[sliceNo];
648
649 // Assumption (see preconditions): all slices are nodes
650 //
651 equivNode = equivProxy->GetNode();
652 maxNo = equivNode->GetMaxEquivalentSliceNo();
653 if (maxNo != sliceNo)
654 {
655#ifdef G4GEOMETRY_VOXELDEBUG
656 G4cout << "**** G4SmartVoxelHeader::CollectEquivalentNodes" << G4endl
657 << " Collecting Nodes = "
658 << sliceNo << " - " << maxNo << G4endl;
659#endif
660 // Do collection between sliceNo and maxNo inclusive
661 //
662 for (equivNo=sliceNo+1; equivNo<=maxNo; equivNo++)
663 {
664 delete fslices[equivNo]->GetNode();
665 delete fslices[equivNo];
666 fslices[equivNo] = equivProxy;
667 }
668 sliceNo = maxNo;
669 }
670 }
671}
672
673// ***************************************************************************
674// Collects common headers at our level, deleting all but one to save
675// memory, and adjusting stored slice pointers appropriately.
676//
677// Preconditions:
678// o if a header forms part of a range of equivalent slices
679// (ie. GetMaxEquivalentSliceNo()>GetMinEquivalentSliceNo()),
680// it is assumed that all slices in the range are headers.
681// o this will be true if a constant Expression is used to evaluate
682// when to refine nodes.
683// ***************************************************************************
684//
685void G4SmartVoxelHeader::CollectEquivalentHeaders()
686{
687 G4int sliceNo, maxNo, equivNo;
688 G4int maxNode = fslices.size();
689 G4SmartVoxelHeader *equivHeader, *sampleHeader;
690 G4SmartVoxelProxy *equivProxy;
691
692 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
693 {
694 equivProxy = fslices[sliceNo];
695 if (equivProxy->IsHeader())
696 {
697 equivHeader = equivProxy->GetHeader();
698 maxNo = equivHeader->GetMaxEquivalentSliceNo();
699 if (maxNo != sliceNo)
700 {
701 // Attempt collection between sliceNo and maxNo inclusive:
702 // look for common headers. All slices between sliceNo and maxNo
703 // are guaranteed to be headers but may not have equal contents
704 //
705#ifdef G4GEOMETRY_VOXELDEBUG
706 G4cout << "**** G4SmartVoxelHeader::CollectEquivalentHeaders" << G4endl
707 << " Collecting Headers =";
708#endif
709 for (equivNo=sliceNo+1; equivNo<=maxNo; equivNo++)
710 {
711 sampleHeader = fslices[equivNo]->GetHeader();
712 if ( (*sampleHeader) == (*equivHeader) )
713 {
714#ifdef G4GEOMETRY_VOXELDEBUG
715 G4cout << " " << equivNo;
716#endif
717 // Delete sampleHeader + proxy and replace with equivHeader/Proxy
718 //
719 delete sampleHeader;
720 delete fslices[equivNo];
721 fslices[equivNo] = equivProxy;
722 }
723 else
724 {
725 // Not equal. Set this header to be
726 // the current header for comparisons
727 //
728 equivProxy = fslices[equivNo];
729 equivHeader = equivProxy->GetHeader();
730 }
731
732 }
733#ifdef G4GEOMETRY_VOXELDEBUG
734 G4cout << G4endl;
735#endif
736 // Skip past examined slices
737 //
738 sliceNo = maxNo;
739 }
740 }
741 }
742}
743
744// ***************************************************************************
745// Builds the nodes corresponding to slices between the specified limits
746// and along the specified axis, using candidate volume no.s in the vector
747// pCandidates. If the `daughters' are replicated volumes (ie. the logical
748// volume has a single replicated/parameterised volume for a daughter)
749// the candidate no.s are interpreted as PARAMETERISED volume no.s &
750// PARAMETERISATIONs are applied to compute transformations & solid
751// dimensions appropriately. The volume must be parameterised - ie. has a
752// parameterisation object & non-consuming) - in this case.
753//
754// Returns pointer to built node "structure" (guaranteed non NULL) consisting
755// of G4SmartVoxelNodeProxies referring to G4SmartVoxelNodes.
756// ***************************************************************************
757//
758G4ProxyVector* G4SmartVoxelHeader::BuildNodes(G4LogicalVolume* pVolume,
759 G4VoxelLimits pLimits,
760 const G4VolumeNosVector* pCandidates,
761 EAxis pAxis)
762{
763 G4double motherMinExtent, motherMaxExtent, targetMinExtent, targetMaxExtent;
764 G4VPhysicalVolume *pDaughter=0;
765 G4VPVParameterisation *pParam=0;
766 G4VSolid *targetSolid;
767 G4AffineTransform targetTransform;
768 G4bool replicated;
769 G4int nCandidates = pCandidates->size();
770 G4int nVol, nNode, targetVolNo;
771 G4VoxelLimits noLimits;
772
773#ifdef G4GEOMETRY_VOXELDEBUG
774 G4cout << "**** G4SmartVoxelHeader::BuildNodes" << G4endl
775 << " Limits = " << pLimits << G4endl
776 << " Axis = " << pAxis << G4endl
777 << " Candidates = " << nCandidates << G4endl;
778#endif
779
780 // Compute extent of logical volume's solid along this axis
781 // NOTE: results stored locally and not preserved/reused
782 //
783 G4VSolid* outerSolid = pVolume->GetSolid();
784 const G4AffineTransform origin;
785 if( !outerSolid->CalculateExtent(pAxis, pLimits, origin,
786 motherMinExtent, motherMaxExtent) )
787 {
788 outerSolid->CalculateExtent(pAxis, noLimits, origin,
789 motherMinExtent, motherMaxExtent);
790 }
791 G4VolumeExtentVector minExtents(nCandidates,0.);
792 G4VolumeExtentVector maxExtents(nCandidates,0.);
793
794 if ( (pVolume->GetNoDaughters()==1)
795 && (pVolume->GetDaughter(0)->IsReplicated()==true) )
796 {
797 // Replication data not required: only parameterisation object
798 // and volume no. List used
799 //
800 pDaughter = pVolume->GetDaughter(0);
801 pParam = pDaughter->GetParameterisation();
802 if (!pParam)
803 {
804 G4cerr << "PANIC! - G4SmartVoxelHeader::BuildNodes()" << G4endl
805 << " Replicated volume with no parameterisation object !"
806 << G4endl;
807 G4Exception("G4SmartVoxelHeader::BuildNodes()", "InvalidSetup",
808 FatalException, "Missing parameterisation.");
809 }
810
811 // Setup daughter's transformations
812 //
813 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
814 pDaughter->GetTranslation());
815 replicated = true;
816 }
817 else
818 {
819 replicated = false;
820 }
821
822 // Compute extents
823 //
824 for (nVol=0; nVol<nCandidates; nVol++)
825 {
826 targetVolNo=(*pCandidates)[nVol];
827 if (replicated == false)
828 {
829 pDaughter=pVolume->GetDaughter(targetVolNo);
830
831 // Setup daughter's transformations
832 //
833 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
834 pDaughter->GetTranslation());
835 // Get underlying (and setup) solid
836 //
837 targetSolid = pDaughter->GetLogicalVolume()->GetSolid();
838 }
839 else
840 {
841 // Find solid
842 //
843 targetSolid = pParam->ComputeSolid(targetVolNo,pDaughter);
844
845 // Setup solid
846 //
847 targetSolid->ComputeDimensions(pParam,targetVolNo,pDaughter);
848
849 // Setup transform
850 //
851 pParam->ComputeTransformation(targetVolNo,pDaughter);
852 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
853 pDaughter->GetTranslation());
854 }
855 // Calculate extents
856 //
857 if(!targetSolid->CalculateExtent(pAxis, pLimits, targetTransform,
858 targetMinExtent, targetMaxExtent))
859 {
860 targetSolid->CalculateExtent(pAxis, noLimits, targetTransform,
861 targetMinExtent,targetMaxExtent);
862 }
863 minExtents[nVol] = targetMinExtent;
864 maxExtents[nVol] = targetMaxExtent;
865
866 // Check not entirely outside mother when processing toplevel nodes
867 //
868 if ( (!pLimits.IsLimited()) && ((targetMaxExtent<=motherMinExtent)
869 ||(targetMinExtent>=motherMaxExtent)) )
870 {
871 G4cerr << "PANIC! - G4SmartVoxelHeader::BuildNodes()" << G4endl
872 << " Daughter physical volume "
873 << pDaughter->GetName() << G4endl
874 << " is entirely outside mother logical volume "
875 << pVolume->GetName() << " !!" << G4endl;
876 G4Exception("G4SmartVoxelHeader::BuildNodes()", "InvalidSetup",
877 FatalException, "Overlapping daughter with mother volume.");
878 }
879
880#ifdef G4GEOMETRY_VOXELDEBUG
881 // Check for straddling volumes when debugging.
882 // If a volume is >kStraddlePercent percent over the mother
883 // boundary, print a warning.
884 //
885 if (!pLimits.IsLimited())
886 {
887 G4double width;
888 G4int kStraddlePercent=5;
889 width = maxExtents[nVol]-minExtents[nVol];
890 if ( (((motherMinExtent-minExtents[nVol])*100/width) > kStraddlePercent)
891 ||(((maxExtents[nVol]-motherMaxExtent)*100/width) > kStraddlePercent) )
892 {
893 G4cout << "**** G4SmartVoxelHeader::BuildNodes" << G4endl
894 << " WARNING : Daughter # " << nVol
895 << " name = " << pDaughter->GetName() << G4endl
896 << " Crosses mother boundary of logical volume, name = "
897 << pVolume->GetName() << G4endl
898 << " by more than " << kStraddlePercent
899 << "%" << G4endl;
900 }
901 }
902#endif
903 }
904
905 // Extents of all daughters known
906
907 // Calculate minimum slice width, only including volumes inside the limits
908 //
909 G4double minWidth = kInfinity;
910 G4double currentWidth;
911 for (nVol=0; nVol<nCandidates; nVol++)
912 {
913 // currentWidth should -always- be a positive value. Inaccurate computed extent
914 // from the solid or situations of malformed geometries (overlaps) may lead to
915 // negative values and therefore unpredictable crashes !
916 //
917 currentWidth = std::abs(maxExtents[nVol]-minExtents[nVol]);
918 if ( (currentWidth<minWidth)
919 && (maxExtents[nVol]>=pLimits.GetMinExtent(pAxis))
920 && (minExtents[nVol]<=pLimits.GetMaxExtent(pAxis)) )
921 {
922 minWidth = currentWidth;
923 }
924 }
925
926 // No. of Nodes formula - nearest integer to
927 // mother width/half min daughter width +1
928 //
929 G4double noNodesExactD = ((motherMaxExtent-motherMinExtent)*2.0/minWidth)+1.0;
930
931 // Compare with "smartless quality", i.e. the average number of slices
932 // used per contained volume.
933 //
934 G4double smartlessComputed = noNodesExactD / nCandidates;
935 G4double smartlessUser = pVolume->GetSmartless();
936 G4double smartless = (smartlessComputed <= smartlessUser)
937 ? smartlessComputed : smartlessUser;
938 G4double noNodesSmart = smartless*nCandidates;
939 G4int noNodesExactI = G4int(noNodesSmart);
940 G4int noNodes = ((noNodesSmart-noNodesExactI)>=0.5)
941 ? noNodesExactI+1 : noNodesExactI;
942 if( noNodes == 0 ) { noNodes=1; }
943
944#ifdef G4GEOMETRY_VOXELDEBUG
945 G4cout << " Smartless computed = " << smartlessComputed << G4endl
946 << " Smartless volume = " << smartlessUser
947 << " => # Smartless = " << smartless << G4endl;
948 G4cout << " Min width = " << minWidth
949 << " => # Nodes = " << noNodes << G4endl;
950#endif
951
952 if (noNodes>kMaxVoxelNodes)
953 {
954 noNodes=kMaxVoxelNodes;
955#ifdef G4GEOMETRY_VOXELDEBUG
956 G4cout << " Nodes Clipped to = " << kMaxVoxelNodes << G4endl;
957#endif
958 }
959 G4double nodeWidth = (motherMaxExtent-motherMinExtent)/noNodes;
960
961// Create G4VoxelNodes. Will Add proxies before setting fslices
962//
963 G4NodeVector* nodeList = new G4NodeVector();
964 nodeList->reserve(noNodes);
965 if (!nodeList)
966 {
967 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
968 << " NodeList allocation failed." << G4endl;
969 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
970 FatalException, "NodeList allocation error.");
971 }
972 for (nNode=0; nNode<noNodes; nNode++)
973 {
974 G4SmartVoxelNode *pNode;
975 pNode = new G4SmartVoxelNode(nNode);
976 if (!pNode)
977 {
978 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
979 << " Node allocation failed." << G4endl;
980 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
981 FatalException, "Node allocation error.");
982 }
983 nodeList->push_back(pNode);
984 }
985
986 // All nodes created (empty)
987
988 // Fill nodes: Step through extent lists
989 //
990 for (nVol=0; nVol<nCandidates; nVol++)
991 {
992 G4int nodeNo, minContainingNode, maxContainingNode;
993 minContainingNode = G4int((minExtents[nVol]-motherMinExtent)/nodeWidth);
994 maxContainingNode = G4int((maxExtents[nVol]-motherMinExtent)/nodeWidth);
995
996 // Only add nodes that are inside the limits of the axis
997 //
998 if ( (maxContainingNode>=0) && (minContainingNode<noNodes) )
999 {
1000 // If max extent is on max boundary => maxContainingNode=noNodes:
1001 // should be one less as nodeList has noNodes entries
1002 //
1003 if (maxContainingNode>=noNodes)
1004 {
1005 maxContainingNode = noNodes-1;
1006 }
1007 //
1008 // Protection against protruding volumes
1009 //
1010 if (minContainingNode<0)
1011 {
1012 minContainingNode=0;
1013 }
1014 for (nodeNo=minContainingNode; nodeNo<=maxContainingNode; nodeNo++)
1015 {
1016 (*nodeList)[nodeNo]->Insert((*pCandidates)[nVol]);
1017 }
1018 }
1019 }
1020
1021 // All nodes filled
1022
1023 // Create proxy List : caller has deletion responsibility
1024 // (but we must delete nodeList *itself* - not the contents)
1025 //
1026 G4ProxyVector* proxyList = new G4ProxyVector();
1027 proxyList->reserve(noNodes);
1028 if (!proxyList)
1029 {
1030 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
1031 << " Proxy list allocation failed." << G4endl;
1032 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
1033 FatalException, "Proxy list allocation error.");
1034 }
1035 //
1036 // Fill proxy List
1037 //
1038 for (nNode=0; nNode<noNodes; nNode++)
1039 {
1040 G4SmartVoxelProxy* pProxyNode = new G4SmartVoxelProxy((*nodeList)[nNode]);
1041 if (!pProxyNode)
1042 {
1043 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
1044 << " Proxy node allocation failed." << G4endl;
1045 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
1046 FatalException, "Proxy node allocation failed.");
1047 }
1048 proxyList->push_back(pProxyNode);
1049 }
1050 delete nodeList;
1051 return proxyList;
1052}
1053
1054// ***************************************************************************
1055// Calculate a "quality value" for the specified vector of voxels.
1056// The value returned should be >0 and such that the smaller the number
1057// the higher the quality of the slice.
1058//
1059// Preconditions: pSlice must consist of G4SmartVoxelNodeProxies only
1060// Process:
1061// o Examine each node in turn, summing:
1062// no. of non-empty nodes
1063// no. of volumes in each node
1064// o Calculate Quality=sigma(volumes in nod)/(no. of non-empty nodes)
1065// if all nodes empty, return kInfinity
1066// o Call G4Exception on finding a G4SmartVoxelHeaderProxy
1067// ***************************************************************************
1068//
1069G4double G4SmartVoxelHeader::CalculateQuality(G4ProxyVector *pSlice)
1070{
1071 G4double quality;
1072 G4int nNodes = pSlice->size();
1073 G4int noContained, maxContained=0, sumContained=0, sumNonEmptyNodes=0;
1074 G4SmartVoxelNode *node;
1075
1076 for (G4int i=0; i<nNodes; i++)
1077 {
1078 if ((*pSlice)[i]->IsNode())
1079 {
1080 // Definitely a node. Add info to running totals
1081 //
1082 node = (*pSlice)[i]->GetNode();
1083 noContained = node->GetNoContained();
1084 if (noContained)
1085 {
1086 sumNonEmptyNodes++;
1087 sumContained += noContained;
1088 //
1089 // Calc maxContained for statistics
1090 //
1091 if (noContained>maxContained)
1092 {
1093 maxContained = noContained;
1094 }
1095 }
1096 }
1097 else
1098 {
1099 G4cout << "ERROR - G4SmartVoxelHeader::CalculateQuality()" << G4endl
1100 << " Not defined for sliced volumes." << G4endl;
1101 G4Exception("G4SmartVoxelHeader::CalculateQuality()", "NotApplicable",
1102 FatalException, "Not applicable to replicated volumes.");
1103 }
1104 }
1105
1106 // Calculate quality with protection against no non-empty nodes
1107 //
1108 if (sumNonEmptyNodes)
1109 {
1110 quality = sumContained/sumNonEmptyNodes;
1111 }
1112 else
1113 {
1114 quality = kInfinity;
1115 }
1116
1117#ifdef G4GEOMETRY_VOXELDEBUG
1118 G4cout << "**** G4SmartVoxelHeader::CalculateQuality" << G4endl
1119 << " Quality = " << quality << G4endl
1120 << " Nodes = " << nNodes
1121 << " of which " << sumNonEmptyNodes << " non empty" << G4endl
1122 << " Max Contained = " << maxContained << G4endl;
1123#endif
1124
1125 return quality;
1126}
1127
1128// ***************************************************************************
1129// Examined each contained node, refines (creates a replacement additional
1130// dimension of voxels) when there is more than one voxel in the slice.
1131// Does not refine further if already limited in two dimensions (=> this
1132// is the third level of limits)
1133//
1134// Preconditions: slices (nodes) have been built.
1135// ***************************************************************************
1136//
1137void G4SmartVoxelHeader::RefineNodes(G4LogicalVolume* pVolume,
1138 G4VoxelLimits pLimits)
1139{
1140 G4int refinedDepth=0, minVolumes;
1141 G4int maxNode = fslices.size();
1142
1143 if (pLimits.IsXLimited())
1144 {
1145 refinedDepth++;
1146 }
1147 if (pLimits.IsYLimited())
1148 {
1149 refinedDepth++;
1150 }
1151 if (pLimits.IsZLimited())
1152 {
1153 refinedDepth++;
1154 }
1155
1156 // Calculate minimum number of volumes necessary to refine
1157 //
1158 switch (refinedDepth)
1159 {
1160 case 0:
1161 minVolumes=kMinVoxelVolumesLevel2;
1162 break;
1163 case 1:
1164 minVolumes=kMinVoxelVolumesLevel3;
1165 break;
1166 default:
1167 minVolumes=10000; // catch refinedDepth=3 and errors
1168 break;
1169 }
1170
1171 if (refinedDepth<2)
1172 {
1173 G4int targetNo, noContainedDaughters, minNo, maxNo, replaceNo, i;
1174 G4double sliceWidth = (fmaxExtent-fminExtent)/maxNode;
1175 G4VoxelLimits newLimits;
1176 G4SmartVoxelNode* targetNode;
1177 G4SmartVoxelProxy* targetNodeProxy;
1178 G4SmartVoxelHeader* replaceHeader;
1179 G4SmartVoxelProxy* replaceHeaderProxy;
1180 G4VolumeNosVector* targetList;
1181 G4SmartVoxelProxy* lastProxy;
1182
1183 for (targetNo=0; targetNo<maxNode; targetNo++)
1184 {
1185 // Assume all slices are nodes (see preconditions)
1186 //
1187 targetNodeProxy = fslices[targetNo];
1188 targetNode = targetNodeProxy->GetNode();
1189
1190 if (targetNode->GetNoContained() >= minVolumes)
1191 {
1192 noContainedDaughters = targetNode->GetNoContained();
1193 targetList = new G4VolumeNosVector();
1194 targetList->reserve(noContainedDaughters);
1195 if (!targetList)
1196 {
1197 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1198 << " Target volume node list allocation failed."
1199 << G4endl;
1200 G4Exception("G4SmartVoxelHeader::RefineNodes()",
1201 "FatalError", FatalException,
1202 "Target volume node list allocation error.");
1203 }
1204 for (i=0; i<noContainedDaughters; i++)
1205 {
1206 targetList->push_back(targetNode->GetVolume(i));
1207 }
1208 minNo = targetNode->GetMinEquivalentSliceNo();
1209 maxNo = targetNode->GetMaxEquivalentSliceNo();
1210
1211#ifdef G4GEOMETRY_VOXELDEBUG
1212 G4cout << "**** G4SmartVoxelHeader::RefineNodes" << G4endl
1213 << " Refining nodes " << minNo
1214 << " - " << maxNo << " inclusive" << G4endl;
1215#endif
1216 // Delete node proxies at start of collected sets of nodes/headers
1217 //
1218 lastProxy=0;
1219 for (replaceNo=minNo; replaceNo<=maxNo; replaceNo++)
1220 {
1221 if (lastProxy != fslices[replaceNo])
1222 {
1223 lastProxy=fslices[replaceNo];
1224 delete lastProxy;
1225 }
1226 }
1227 // Delete node to be replaced
1228 //
1229 delete targetNode;
1230
1231 // Create new headers + proxies and replace in fslices
1232 //
1233 newLimits = pLimits;
1234 newLimits.AddLimit(faxis,fminExtent+sliceWidth*minNo,
1235 fminExtent+sliceWidth*(maxNo+1));
1236 replaceHeader = new G4SmartVoxelHeader(pVolume,newLimits,
1237 targetList,replaceNo);
1238 if (!replaceHeader)
1239 {
1240 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1241 << " Refined VoxelHeader allocation failed." << G4endl;
1242 G4Exception("G4SmartVoxelHeader::RefineNodes()", "FatalError",
1243 FatalException, "Refined VoxelHeader allocation error.");
1244 }
1245 replaceHeader->SetMinEquivalentSliceNo(minNo);
1246 replaceHeader->SetMaxEquivalentSliceNo(maxNo);
1247 replaceHeaderProxy = new G4SmartVoxelProxy(replaceHeader);
1248 if (!replaceHeader)
1249 {
1250 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1251 << " Refined VoxelProxy allocation failed." << G4endl;
1252 G4Exception("G4SmartVoxelHeader::RefineNodes()", "FatalError",
1253 FatalException, "Refined VoxelProxy allocation error.");
1254 }
1255 for (replaceNo=minNo; replaceNo<=maxNo; replaceNo++)
1256 {
1257 fslices[replaceNo] = replaceHeaderProxy;
1258 }
1259 // Finished replacing current `equivalent' group
1260 //
1261 delete targetList;
1262 targetNo=maxNo;
1263 }
1264 }
1265 }
1266}
1267
1268// ***************************************************************************
1269// Returns true if all slices have equal contents.
1270// Preconditions: all equal slices have been collected.
1271// Procedure:
1272// o checks all slice proxy pointers are equal
1273// o returns true if only one slice or all slice proxies pointers equal.
1274// ***************************************************************************
1275//
1276G4bool G4SmartVoxelHeader::AllSlicesEqual() const
1277{
1278 G4int noSlices = fslices.size();
1279 G4SmartVoxelProxy* refProxy;
1280
1281 if (noSlices>1)
1282 {
1283 refProxy=fslices[0];
1284 for (G4int i=1; i<noSlices; i++)
1285 {
1286 if (refProxy!=fslices[i])
1287 {
1288 return false;
1289 }
1290 }
1291 }
1292 return true;
1293}
1294
1295// ***************************************************************************
1296// Streaming operator for debugging.
1297// ***************************************************************************
1298//
1299std::ostream& operator << (std::ostream& s, const G4SmartVoxelHeader& h)
1300{
1301 s << "Axis = " << G4int(h.faxis) << G4endl;
1302 G4SmartVoxelProxy *collectNode=0, *collectHead=0;
1303 G4int collectNodeNo=0;
1304 G4int collectHeadNo=0;
1305 size_t i, j;
1306 G4bool haveHeaders=false;
1307
1308 for (i=0; i<h.fslices.size(); i++)
1309 {
1310 s << "Slice #" << i << " = ";
1311 if (h.fslices[i]->IsNode())
1312 {
1313 if (h.fslices[i]!=collectNode)
1314 {
1315 s << "{";
1316 for (G4int j=0; j<h.fslices[i]->GetNode()->GetNoContained(); j++)
1317 {
1318 s << " " << h.fslices[i]->GetNode()->GetVolume(j);
1319 }
1320 s << " }" << G4endl;
1321 collectNode = h.fslices[i];
1322 collectNodeNo = i;
1323 }
1324 else
1325 {
1326 s << "As slice #" << collectNodeNo << G4endl;
1327 }
1328 }
1329 else
1330 {
1331 haveHeaders=true;
1332 if (h.fslices[i] != collectHead)
1333 {
1334 s << "Header" << G4endl;
1335 collectHead = h.fslices[i];
1336 collectHeadNo = i;
1337 }
1338 else
1339 {
1340 s << "As slice #" << collectHeadNo << G4endl;
1341 }
1342 }
1343 }
1344
1345 if (haveHeaders)
1346 {
1347 collectHead=0;
1348 for (j=0; j<h.fslices.size(); j++)
1349 {
1350 if (h.fslices[j]->IsHeader())
1351 {
1352 s << "Header at Slice #" << j << " = ";
1353 if (h.fslices[j] != collectHead)
1354 {
1355 s << G4endl
1356 << (*(h.fslices[j]->GetHeader()));
1357 collectHead = h.fslices[j];
1358 collectHeadNo = j;
1359 }
1360 else
1361 {
1362 s << "As slice #" << collectHeadNo << G4endl;
1363 }
1364 }
1365 }
1366 }
1367 return s;
1368}
Note: See TracBrowser for help on using the repository browser.