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

Last change on this file since 1341 was 1340, checked in by garnier, 15 years ago

update ti head
File size: 45.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.39 2010/09/06 09:39:21 gcosmo Exp $
28// GEANT4 tag $Name: geommng-V09-03-05 $
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 emin = kInfinity, emax = -kInfinity;
369 G4VoxelLimits limits;
370 G4AffineTransform origin;
371 pVolume->GetSolid()->CalculateExtent(axis, limits, origin, emin, emax);
372 if ( (std::fabs((emin-fminExtent)/fminExtent) +
373 std::fabs((emax-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(); )
513 {
514 if (*i==tmpProx)
515 {
516 i = pTestSlices->erase(i);
517 }
518 else
519 {
520 ++i;
521 }
522 }
523 if ( tmpProx ) { delete tmpProx; }
524 }
525 delete pTestSlices;
526 }
527 }
528 }
529 // Check for error case.. when limits already 3d,
530 // so cannot select a new axis
531 //
532 if (!pGoodSlices)
533 {
534 G4cerr << "ERROR - G4SmartVoxelHeader::BuildVoxelsWithinLimits()"
535 << G4endl
536 << " Illegal limits: only 3 dimensions allowed." << G4endl;
537 G4Exception("G4SmartVoxelHeader::BuildVoxelsWithinLimits()",
538 "InvalidSetup", FatalException,
539 "Cannot select more than 3 axis for optimisation.");
540 return;
541 }
542
543 //
544 // We have selected pGoodSlices, with a score testSliceScore
545 //
546
547 // Store chosen axis, slice ptr
548 //
549 fslices=*pGoodSlices; // Set slice information, copy ptrs in collection
550 delete pGoodSlices; // Destroy slices vector, but not contained
551 // proxies or nodes
552 faxis=goodSliceAxis;
553
554#ifdef G4GEOMETRY_VOXELDEBUG
555 G4cout << G4endl << " Volume = " << pVolume->GetName()
556 << G4endl << " Selected axis = " << faxis << G4endl;
557 for (size_t islice=0; islice<fslices.size(); islice++)
558 {
559 G4cout << " Node #" << islice << " = {";
560 for (G4int j=0; j<fslices[islice]->GetNode()->GetNoContained(); j++)
561 {
562 G4cout << " " << fslices[islice]->GetNode()->GetVolume(j);
563 }
564 G4cout << " }" << G4endl;
565 }
566 G4cout << G4endl;
567#endif
568
569 // Calculate and set min and max extents given our axis
570 //
571 G4VSolid* outerSolid = pVolume->GetSolid();
572 const G4AffineTransform origin;
573 if(!outerSolid->CalculateExtent(faxis,pLimits,origin,fminExtent,fmaxExtent))
574 {
575 outerSolid->CalculateExtent(faxis,noLimits,origin,fminExtent,fmaxExtent);
576 }
577
578 // Calculate equivalent nos
579 //
580 BuildEquivalentSliceNos();
581 CollectEquivalentNodes(); // Collect common nodes
582 RefineNodes(pVolume,pLimits); // Refine nodes creating headers
583
584 // No common headers can exist because collapsed by construction
585}
586
587// ***************************************************************************
588// Calculates and stores the minimum and maximum equivalent neighbour
589// values for all slices at our level.
590//
591// Precondition: all slices are nodes.
592// For each potential start of a group of equivalent nodes:
593// o searches forwards in fslices to find group end
594// o loops from start to end setting start and end slices.
595// ***************************************************************************
596//
597void G4SmartVoxelHeader::BuildEquivalentSliceNos()
598{
599 G4int sliceNo, minNo, maxNo, equivNo;
600 G4int maxNode = fslices.size();
601 G4SmartVoxelNode *startNode, *sampleNode;
602 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
603 {
604 minNo = sliceNo;
605
606 // Get first node (see preconditions - will throw exception if a header)
607 //
608 startNode = fslices[minNo]->GetNode();
609
610 // Find max equivalent
611 //
612 for (equivNo=minNo+1; equivNo<maxNode; equivNo++)
613 {
614 sampleNode = fslices[equivNo]->GetNode();
615 if (!((*startNode) == (*sampleNode))) { break; }
616 }
617 maxNo = equivNo-1;
618 if (maxNo != minNo)
619 {
620 // Set min and max nos
621 //
622 for (equivNo=minNo; equivNo<=maxNo; equivNo++)
623 {
624 sampleNode = fslices[equivNo]->GetNode();
625 sampleNode->SetMinEquivalentSliceNo(minNo);
626 sampleNode->SetMaxEquivalentSliceNo(maxNo);
627 }
628 // Advance outer loop to end of equivalent group
629 //
630 sliceNo = maxNo;
631 }
632 }
633}
634
635// ***************************************************************************
636// Collects common nodes at our level, deleting all but one to save
637// memory, and adjusting stored slice pointers appropriately.
638//
639// Preconditions:
640// o the slices have not previously be "collected"
641// o all of the slices are nodes.
642// ***************************************************************************
643//
644void G4SmartVoxelHeader::CollectEquivalentNodes()
645{
646 G4int sliceNo, maxNo, equivNo;
647 G4int maxNode=fslices.size();
648 G4SmartVoxelNode *equivNode;
649 G4SmartVoxelProxy *equivProxy;
650
651 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
652 {
653 equivProxy=fslices[sliceNo];
654
655 // Assumption (see preconditions): all slices are nodes
656 //
657 equivNode = equivProxy->GetNode();
658 maxNo = equivNode->GetMaxEquivalentSliceNo();
659 if (maxNo != sliceNo)
660 {
661#ifdef G4GEOMETRY_VOXELDEBUG
662 G4cout << "**** G4SmartVoxelHeader::CollectEquivalentNodes" << G4endl
663 << " Collecting Nodes = "
664 << sliceNo << " - " << maxNo << G4endl;
665#endif
666 // Do collection between sliceNo and maxNo inclusive
667 //
668 for (equivNo=sliceNo+1; equivNo<=maxNo; equivNo++)
669 {
670 delete fslices[equivNo]->GetNode();
671 delete fslices[equivNo];
672 fslices[equivNo] = equivProxy;
673 }
674 sliceNo = maxNo;
675 }
676 }
677}
678
679// ***************************************************************************
680// Collects common headers at our level, deleting all but one to save
681// memory, and adjusting stored slice pointers appropriately.
682//
683// Preconditions:
684// o if a header forms part of a range of equivalent slices
685// (ie. GetMaxEquivalentSliceNo()>GetMinEquivalentSliceNo()),
686// it is assumed that all slices in the range are headers.
687// o this will be true if a constant Expression is used to evaluate
688// when to refine nodes.
689// ***************************************************************************
690//
691void G4SmartVoxelHeader::CollectEquivalentHeaders()
692{
693 G4int sliceNo, maxNo, equivNo;
694 G4int maxNode = fslices.size();
695 G4SmartVoxelHeader *equivHeader, *sampleHeader;
696 G4SmartVoxelProxy *equivProxy;
697
698 for (sliceNo=0; sliceNo<maxNode; sliceNo++)
699 {
700 equivProxy = fslices[sliceNo];
701 if (equivProxy->IsHeader())
702 {
703 equivHeader = equivProxy->GetHeader();
704 maxNo = equivHeader->GetMaxEquivalentSliceNo();
705 if (maxNo != sliceNo)
706 {
707 // Attempt collection between sliceNo and maxNo inclusive:
708 // look for common headers. All slices between sliceNo and maxNo
709 // are guaranteed to be headers but may not have equal contents
710 //
711#ifdef G4GEOMETRY_VOXELDEBUG
712 G4cout << "**** G4SmartVoxelHeader::CollectEquivalentHeaders" << G4endl
713 << " Collecting Headers =";
714#endif
715 for (equivNo=sliceNo+1; equivNo<=maxNo; equivNo++)
716 {
717 sampleHeader = fslices[equivNo]->GetHeader();
718 if ( (*sampleHeader) == (*equivHeader) )
719 {
720#ifdef G4GEOMETRY_VOXELDEBUG
721 G4cout << " " << equivNo;
722#endif
723 // Delete sampleHeader + proxy and replace with equivHeader/Proxy
724 //
725 delete sampleHeader;
726 delete fslices[equivNo];
727 fslices[equivNo] = equivProxy;
728 }
729 else
730 {
731 // Not equal. Set this header to be
732 // the current header for comparisons
733 //
734 equivProxy = fslices[equivNo];
735 equivHeader = equivProxy->GetHeader();
736 }
737
738 }
739#ifdef G4GEOMETRY_VOXELDEBUG
740 G4cout << G4endl;
741#endif
742 // Skip past examined slices
743 //
744 sliceNo = maxNo;
745 }
746 }
747 }
748}
749
750// ***************************************************************************
751// Builds the nodes corresponding to slices between the specified limits
752// and along the specified axis, using candidate volume no.s in the vector
753// pCandidates. If the `daughters' are replicated volumes (ie. the logical
754// volume has a single replicated/parameterised volume for a daughter)
755// the candidate no.s are interpreted as PARAMETERISED volume no.s &
756// PARAMETERISATIONs are applied to compute transformations & solid
757// dimensions appropriately. The volume must be parameterised - ie. has a
758// parameterisation object & non-consuming) - in this case.
759//
760// Returns pointer to built node "structure" (guaranteed non NULL) consisting
761// of G4SmartVoxelNodeProxies referring to G4SmartVoxelNodes.
762// ***************************************************************************
763//
764G4ProxyVector* G4SmartVoxelHeader::BuildNodes(G4LogicalVolume* pVolume,
765 G4VoxelLimits pLimits,
766 const G4VolumeNosVector* pCandidates,
767 EAxis pAxis)
768{
769 G4double motherMinExtent= kInfinity, motherMaxExtent= -kInfinity,
770 targetMinExtent= kInfinity, targetMaxExtent= -kInfinity;
771 G4VPhysicalVolume *pDaughter=0;
772 G4VPVParameterisation *pParam=0;
773 G4VSolid *targetSolid;
774 G4AffineTransform targetTransform;
775 G4bool replicated;
776 G4int nCandidates = pCandidates->size();
777 G4int nVol, nNode, targetVolNo;
778 G4VoxelLimits noLimits;
779
780#ifdef G4GEOMETRY_VOXELDEBUG
781 G4cout << "**** G4SmartVoxelHeader::BuildNodes" << G4endl
782 << " Limits = " << pLimits << G4endl
783 << " Axis = " << pAxis << G4endl
784 << " Candidates = " << nCandidates << G4endl;
785#endif
786
787 // Compute extent of logical volume's solid along this axis
788 // NOTE: results stored locally and not preserved/reused
789 //
790 G4VSolid* outerSolid = pVolume->GetSolid();
791 const G4AffineTransform origin;
792 if( !outerSolid->CalculateExtent(pAxis, pLimits, origin,
793 motherMinExtent, motherMaxExtent) )
794 {
795 outerSolid->CalculateExtent(pAxis, noLimits, origin,
796 motherMinExtent, motherMaxExtent);
797 }
798 G4VolumeExtentVector minExtents(nCandidates,0.);
799 G4VolumeExtentVector maxExtents(nCandidates,0.);
800
801 if ( (pVolume->GetNoDaughters()==1)
802 && (pVolume->GetDaughter(0)->IsReplicated()==true) )
803 {
804 // Replication data not required: only parameterisation object
805 // and volume no. List used
806 //
807 pDaughter = pVolume->GetDaughter(0);
808 pParam = pDaughter->GetParameterisation();
809 if (!pParam)
810 {
811 G4cerr << "PANIC! - G4SmartVoxelHeader::BuildNodes()" << G4endl
812 << " Replicated volume with no parameterisation object !"
813 << G4endl;
814 G4Exception("G4SmartVoxelHeader::BuildNodes()", "InvalidSetup",
815 FatalException, "Missing parameterisation.");
816 return 0;
817 }
818
819 // Setup daughter's transformations
820 //
821 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
822 pDaughter->GetTranslation());
823 replicated = true;
824 }
825 else
826 {
827 replicated = false;
828 }
829
830 // Compute extents
831 //
832 for (nVol=0; nVol<nCandidates; nVol++)
833 {
834 targetVolNo=(*pCandidates)[nVol];
835 if (replicated == false)
836 {
837 pDaughter=pVolume->GetDaughter(targetVolNo);
838
839 // Setup daughter's transformations
840 //
841 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
842 pDaughter->GetTranslation());
843 // Get underlying (and setup) solid
844 //
845 targetSolid = pDaughter->GetLogicalVolume()->GetSolid();
846 }
847 else
848 {
849 // Find solid
850 //
851 targetSolid = pParam->ComputeSolid(targetVolNo,pDaughter);
852
853 // Setup solid
854 //
855 targetSolid->ComputeDimensions(pParam,targetVolNo,pDaughter);
856
857 // Setup transform
858 //
859 pParam->ComputeTransformation(targetVolNo,pDaughter);
860 targetTransform = G4AffineTransform(pDaughter->GetRotation(),
861 pDaughter->GetTranslation());
862 }
863 // Calculate extents
864 //
865 if(!targetSolid->CalculateExtent(pAxis, pLimits, targetTransform,
866 targetMinExtent, targetMaxExtent))
867 {
868 targetSolid->CalculateExtent(pAxis, noLimits, targetTransform,
869 targetMinExtent,targetMaxExtent);
870 }
871 minExtents[nVol] = targetMinExtent;
872 maxExtents[nVol] = targetMaxExtent;
873
874#ifdef G4GEOMETRY_VOXELDEBUG
875 G4cout << "---------------------------------------------------" << G4endl
876 << " Volume = " << pDaughter->GetName() << G4endl
877 << " Min Extent = " << targetMinExtent << G4endl
878 << " Max Extent = " << targetMaxExtent << G4endl
879 << "---------------------------------------------------" << G4endl;
880#endif
881
882 // Check not entirely outside mother when processing toplevel nodes
883 //
884 if ( (!pLimits.IsLimited()) && ((targetMaxExtent<=motherMinExtent)
885 ||(targetMinExtent>=motherMaxExtent)) )
886 {
887 G4cerr << "PANIC! - G4SmartVoxelHeader::BuildNodes()" << G4endl
888 << " Daughter physical volume "
889 << pDaughter->GetName() << G4endl
890 << " is entirely outside mother logical volume "
891 << pVolume->GetName() << " !!" << G4endl;
892 G4Exception("G4SmartVoxelHeader::BuildNodes()", "InvalidSetup",
893 FatalException, "Overlapping daughter with mother volume.");
894 }
895
896#ifdef G4GEOMETRY_VOXELDEBUG
897 // Check for straddling volumes when debugging.
898 // If a volume is >kStraddlePercent percent over the mother
899 // boundary, print a warning.
900 //
901 if (!pLimits.IsLimited())
902 {
903 G4double width;
904 G4int kStraddlePercent=5;
905 width = maxExtents[nVol]-minExtents[nVol];
906 if ( (((motherMinExtent-minExtents[nVol])*100/width) > kStraddlePercent)
907 ||(((maxExtents[nVol]-motherMaxExtent)*100/width) > kStraddlePercent) )
908 {
909 G4cout << "**** G4SmartVoxelHeader::BuildNodes" << G4endl
910 << " WARNING : Daughter # " << nVol
911 << " name = " << pDaughter->GetName() << G4endl
912 << " Crosses mother boundary of logical volume, name = "
913 << pVolume->GetName() << G4endl
914 << " by more than " << kStraddlePercent
915 << "%" << G4endl;
916 }
917 }
918#endif
919 }
920
921 // Extents of all daughters known
922
923 // Calculate minimum slice width, only including volumes inside the limits
924 //
925 G4double minWidth = kInfinity;
926 G4double currentWidth;
927 for (nVol=0; nVol<nCandidates; nVol++)
928 {
929 // currentWidth should -always- be a positive value. Inaccurate computed extent
930 // from the solid or situations of malformed geometries (overlaps) may lead to
931 // negative values and therefore unpredictable crashes !
932 //
933 currentWidth = std::abs(maxExtents[nVol]-minExtents[nVol]);
934 if ( (currentWidth<minWidth)
935 && (maxExtents[nVol]>=pLimits.GetMinExtent(pAxis))
936 && (minExtents[nVol]<=pLimits.GetMaxExtent(pAxis)) )
937 {
938 minWidth = currentWidth;
939 }
940 }
941
942 // No. of Nodes formula - nearest integer to
943 // mother width/half min daughter width +1
944 //
945 G4double noNodesExactD = ((motherMaxExtent-motherMinExtent)*2.0/minWidth)+1.0;
946
947 // Compare with "smartless quality", i.e. the average number of slices
948 // used per contained volume.
949 //
950 G4double smartlessComputed = noNodesExactD / nCandidates;
951 G4double smartlessUser = pVolume->GetSmartless();
952 G4double smartless = (smartlessComputed <= smartlessUser)
953 ? smartlessComputed : smartlessUser;
954 G4double noNodesSmart = smartless*nCandidates;
955 G4int noNodesExactI = G4int(noNodesSmart);
956 G4int noNodes = ((noNodesSmart-noNodesExactI)>=0.5)
957 ? noNodesExactI+1 : noNodesExactI;
958 if( noNodes == 0 ) { noNodes=1; }
959
960#ifdef G4GEOMETRY_VOXELDEBUG
961 G4cout << " Smartless computed = " << smartlessComputed << G4endl
962 << " Smartless volume = " << smartlessUser
963 << " => # Smartless = " << smartless << G4endl;
964 G4cout << " Min width = " << minWidth
965 << " => # Nodes = " << noNodes << G4endl;
966#endif
967
968 if (noNodes>kMaxVoxelNodes)
969 {
970 noNodes=kMaxVoxelNodes;
971#ifdef G4GEOMETRY_VOXELDEBUG
972 G4cout << " Nodes Clipped to = " << kMaxVoxelNodes << G4endl;
973#endif
974 }
975 G4double nodeWidth = (motherMaxExtent-motherMinExtent)/noNodes;
976
977 // Create G4VoxelNodes. Will Add proxies before setting fslices
978 //
979 G4NodeVector* nodeList = new G4NodeVector();
980 if (!nodeList)
981 {
982 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
983 << " NodeList allocation failed." << G4endl;
984 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
985 FatalException, "NodeList allocation error.");
986 return 0;
987 }
988 nodeList->reserve(noNodes);
989
990 for (nNode=0; nNode<noNodes; nNode++)
991 {
992 G4SmartVoxelNode *pNode;
993 pNode = new G4SmartVoxelNode(nNode);
994 if (!pNode)
995 {
996 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
997 << " Node allocation failed." << G4endl;
998 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
999 FatalException, "Node allocation error.");
1000 return 0;
1001 }
1002 nodeList->push_back(pNode);
1003 }
1004
1005 // All nodes created (empty)
1006
1007 // Fill nodes: Step through extent lists
1008 //
1009 for (nVol=0; nVol<nCandidates; nVol++)
1010 {
1011 G4int nodeNo, minContainingNode, maxContainingNode;
1012 minContainingNode = G4int((minExtents[nVol]-motherMinExtent)/nodeWidth);
1013 maxContainingNode = G4int((maxExtents[nVol]-motherMinExtent)/nodeWidth);
1014
1015 // Only add nodes that are inside the limits of the axis
1016 //
1017 if ( (maxContainingNode>=0) && (minContainingNode<noNodes) )
1018 {
1019 // If max extent is on max boundary => maxContainingNode=noNodes:
1020 // should be one less as nodeList has noNodes entries
1021 //
1022 if (maxContainingNode>=noNodes)
1023 {
1024 maxContainingNode = noNodes-1;
1025 }
1026 //
1027 // Protection against protruding volumes
1028 //
1029 if (minContainingNode<0)
1030 {
1031 minContainingNode=0;
1032 }
1033 for (nodeNo=minContainingNode; nodeNo<=maxContainingNode; nodeNo++)
1034 {
1035 (*nodeList)[nodeNo]->Insert((*pCandidates)[nVol]);
1036 }
1037 }
1038 }
1039
1040 // All nodes filled
1041
1042 // Create proxy List : caller has deletion responsibility
1043 // (but we must delete nodeList *itself* - not the contents)
1044 //
1045 G4ProxyVector* proxyList = new G4ProxyVector();
1046 if (!proxyList)
1047 {
1048 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
1049 << " Proxy list allocation failed." << G4endl;
1050 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
1051 FatalException, "Proxy list allocation error.");
1052 return 0;
1053 }
1054 proxyList->reserve(noNodes);
1055
1056 //
1057 // Fill proxy List
1058 //
1059 for (nNode=0; nNode<noNodes; nNode++)
1060 {
1061 // Get rid of possible excess capacity in the internal node vector
1062 //
1063 ((*nodeList)[nNode])->Shrink();
1064 G4SmartVoxelProxy* pProxyNode = new G4SmartVoxelProxy((*nodeList)[nNode]);
1065 if (!pProxyNode)
1066 {
1067 G4cerr << "ERROR - G4SmartVoxelHeader::BuildNodes()" << G4endl
1068 << " Proxy node allocation failed." << G4endl;
1069 G4Exception("G4SmartVoxelHeader::BuildNodes()", "FatalError",
1070 FatalException, "Proxy node allocation failed.");
1071 return 0;
1072 }
1073 proxyList->push_back(pProxyNode);
1074 }
1075 delete nodeList;
1076 return proxyList;
1077}
1078
1079// ***************************************************************************
1080// Calculate a "quality value" for the specified vector of voxels.
1081// The value returned should be >0 and such that the smaller the number
1082// the higher the quality of the slice.
1083//
1084// Preconditions: pSlice must consist of G4SmartVoxelNodeProxies only
1085// Process:
1086// o Examine each node in turn, summing:
1087// no. of non-empty nodes
1088// no. of volumes in each node
1089// o Calculate Quality=sigma(volumes in nod)/(no. of non-empty nodes)
1090// if all nodes empty, return kInfinity
1091// o Call G4Exception on finding a G4SmartVoxelHeaderProxy
1092// ***************************************************************************
1093//
1094G4double G4SmartVoxelHeader::CalculateQuality(G4ProxyVector *pSlice)
1095{
1096 G4double quality;
1097 G4int nNodes = pSlice->size();
1098 G4int noContained, maxContained=0, sumContained=0, sumNonEmptyNodes=0;
1099 G4SmartVoxelNode *node;
1100
1101 for (G4int i=0; i<nNodes; i++)
1102 {
1103 if ((*pSlice)[i]->IsNode())
1104 {
1105 // Definitely a node. Add info to running totals
1106 //
1107 node = (*pSlice)[i]->GetNode();
1108 noContained = node->GetNoContained();
1109 if (noContained)
1110 {
1111 sumNonEmptyNodes++;
1112 sumContained += noContained;
1113 //
1114 // Calc maxContained for statistics
1115 //
1116 if (noContained>maxContained)
1117 {
1118 maxContained = noContained;
1119 }
1120 }
1121 }
1122 else
1123 {
1124 G4cout << "ERROR - G4SmartVoxelHeader::CalculateQuality()" << G4endl
1125 << " Not defined for sliced volumes." << G4endl;
1126 G4Exception("G4SmartVoxelHeader::CalculateQuality()", "NotApplicable",
1127 FatalException, "Not applicable to replicated volumes.");
1128 }
1129 }
1130
1131 // Calculate quality with protection against no non-empty nodes
1132 //
1133 if (sumNonEmptyNodes)
1134 {
1135 quality = sumContained/sumNonEmptyNodes;
1136 }
1137 else
1138 {
1139 quality = kInfinity;
1140 }
1141
1142#ifdef G4GEOMETRY_VOXELDEBUG
1143 G4cout << "**** G4SmartVoxelHeader::CalculateQuality" << G4endl
1144 << " Quality = " << quality << G4endl
1145 << " Nodes = " << nNodes
1146 << " of which " << sumNonEmptyNodes << " non empty" << G4endl
1147 << " Max Contained = " << maxContained << G4endl;
1148#endif
1149
1150 return quality;
1151}
1152
1153// ***************************************************************************
1154// Examined each contained node, refines (creates a replacement additional
1155// dimension of voxels) when there is more than one voxel in the slice.
1156// Does not refine further if already limited in two dimensions (=> this
1157// is the third level of limits)
1158//
1159// Preconditions: slices (nodes) have been built.
1160// ***************************************************************************
1161//
1162void G4SmartVoxelHeader::RefineNodes(G4LogicalVolume* pVolume,
1163 G4VoxelLimits pLimits)
1164{
1165 G4int refinedDepth=0, minVolumes;
1166 G4int maxNode = fslices.size();
1167
1168 if (pLimits.IsXLimited())
1169 {
1170 refinedDepth++;
1171 }
1172 if (pLimits.IsYLimited())
1173 {
1174 refinedDepth++;
1175 }
1176 if (pLimits.IsZLimited())
1177 {
1178 refinedDepth++;
1179 }
1180
1181 // Calculate minimum number of volumes necessary to refine
1182 //
1183 switch (refinedDepth)
1184 {
1185 case 0:
1186 minVolumes=kMinVoxelVolumesLevel2;
1187 break;
1188 case 1:
1189 minVolumes=kMinVoxelVolumesLevel3;
1190 break;
1191 default:
1192 minVolumes=10000; // catch refinedDepth=3 and errors
1193 break;
1194 }
1195
1196 if (refinedDepth<2)
1197 {
1198 G4int targetNo, noContainedDaughters, minNo, maxNo, replaceNo, i;
1199 G4double sliceWidth = (fmaxExtent-fminExtent)/maxNode;
1200 G4VoxelLimits newLimits;
1201 G4SmartVoxelNode* targetNode;
1202 G4SmartVoxelProxy* targetNodeProxy;
1203 G4SmartVoxelHeader* replaceHeader;
1204 G4SmartVoxelProxy* replaceHeaderProxy;
1205 G4VolumeNosVector* targetList;
1206 G4SmartVoxelProxy* lastProxy;
1207
1208 for (targetNo=0; targetNo<maxNode; targetNo++)
1209 {
1210 // Assume all slices are nodes (see preconditions)
1211 //
1212 targetNodeProxy = fslices[targetNo];
1213 targetNode = targetNodeProxy->GetNode();
1214
1215 if (targetNode->GetNoContained() >= minVolumes)
1216 {
1217 noContainedDaughters = targetNode->GetNoContained();
1218 targetList = new G4VolumeNosVector();
1219 if (!targetList)
1220 {
1221 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1222 << " Target volume node list allocation failed."
1223 << G4endl;
1224 G4Exception("G4SmartVoxelHeader::RefineNodes()",
1225 "FatalError", FatalException,
1226 "Target volume node list allocation error.");
1227 return;
1228 }
1229 targetList->reserve(noContainedDaughters);
1230 for (i=0; i<noContainedDaughters; i++)
1231 {
1232 targetList->push_back(targetNode->GetVolume(i));
1233 }
1234 minNo = targetNode->GetMinEquivalentSliceNo();
1235 maxNo = targetNode->GetMaxEquivalentSliceNo();
1236
1237#ifdef G4GEOMETRY_VOXELDEBUG
1238 G4cout << "**** G4SmartVoxelHeader::RefineNodes" << G4endl
1239 << " Refining nodes " << minNo
1240 << " - " << maxNo << " inclusive" << G4endl;
1241#endif
1242 if (minNo > maxNo) // Delete node and list to be replaced
1243 { // and avoid further action ...
1244 delete targetNode;
1245 delete targetList;
1246 return;
1247 }
1248
1249 // Delete node proxies at start of collected sets of nodes/headers
1250 //
1251 lastProxy=0;
1252 for (replaceNo=minNo; replaceNo<=maxNo; replaceNo++)
1253 {
1254 if (lastProxy != fslices[replaceNo])
1255 {
1256 lastProxy=fslices[replaceNo];
1257 delete lastProxy;
1258 }
1259 }
1260 // Delete node to be replaced
1261 //
1262 delete targetNode;
1263
1264 // Create new headers + proxies and replace in fslices
1265 //
1266 newLimits = pLimits;
1267 newLimits.AddLimit(faxis,fminExtent+sliceWidth*minNo,
1268 fminExtent+sliceWidth*(maxNo+1));
1269 replaceHeader = new G4SmartVoxelHeader(pVolume,newLimits,
1270 targetList,replaceNo);
1271 if (!replaceHeader)
1272 {
1273 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1274 << " Refined VoxelHeader allocation failed." << G4endl;
1275 G4Exception("G4SmartVoxelHeader::RefineNodes()", "FatalError",
1276 FatalException, "Refined VoxelHeader allocation error.");
1277 return;
1278 }
1279 replaceHeader->SetMinEquivalentSliceNo(minNo);
1280 replaceHeader->SetMaxEquivalentSliceNo(maxNo);
1281 replaceHeaderProxy = new G4SmartVoxelProxy(replaceHeader);
1282 if (!replaceHeaderProxy)
1283 {
1284 G4cerr << "ERROR - G4SmartVoxelHeader::RefineNodes()" << G4endl
1285 << " Refined VoxelProxy allocation failed." << G4endl;
1286 G4Exception("G4SmartVoxelHeader::RefineNodes()", "FatalError",
1287 FatalException, "Refined VoxelProxy allocation error.");
1288 return;
1289 }
1290 for (replaceNo=minNo; replaceNo<=maxNo; replaceNo++)
1291 {
1292 fslices[replaceNo] = replaceHeaderProxy;
1293 }
1294 // Finished replacing current `equivalent' group
1295 //
1296 delete targetList;
1297 targetNo=maxNo;
1298 }
1299 }
1300 }
1301}
1302
1303// ***************************************************************************
1304// Returns true if all slices have equal contents.
1305// Preconditions: all equal slices have been collected.
1306// Procedure:
1307// o checks all slice proxy pointers are equal
1308// o returns true if only one slice or all slice proxies pointers equal.
1309// ***************************************************************************
1310//
1311G4bool G4SmartVoxelHeader::AllSlicesEqual() const
1312{
1313 G4int noSlices = fslices.size();
1314 G4SmartVoxelProxy* refProxy;
1315
1316 if (noSlices>1)
1317 {
1318 refProxy=fslices[0];
1319 for (G4int i=1; i<noSlices; i++)
1320 {
1321 if (refProxy!=fslices[i])
1322 {
1323 return false;
1324 }
1325 }
1326 }
1327 return true;
1328}
1329
1330// ***************************************************************************
1331// Streaming operator for debugging.
1332// ***************************************************************************
1333//
1334std::ostream& operator << (std::ostream& s, const G4SmartVoxelHeader& h)
1335{
1336 s << "Axis = " << G4int(h.faxis) << G4endl;
1337 G4SmartVoxelProxy *collectNode=0, *collectHead=0;
1338 G4int collectNodeNo=0;
1339 G4int collectHeadNo=0;
1340 size_t i, j;
1341 G4bool haveHeaders=false;
1342
1343 for (i=0; i<h.fslices.size(); i++)
1344 {
1345 s << "Slice #" << i << " = ";
1346 if (h.fslices[i]->IsNode())
1347 {
1348 if (h.fslices[i]!=collectNode)
1349 {
1350 s << "{";
1351 for (G4int j=0; j<h.fslices[i]->GetNode()->GetNoContained(); j++)
1352 {
1353 s << " " << h.fslices[i]->GetNode()->GetVolume(j);
1354 }
1355 s << " }" << G4endl;
1356 collectNode = h.fslices[i];
1357 collectNodeNo = i;
1358 }
1359 else
1360 {
1361 s << "As slice #" << collectNodeNo << G4endl;
1362 }
1363 }
1364 else
1365 {
1366 haveHeaders=true;
1367 if (h.fslices[i] != collectHead)
1368 {
1369 s << "Header" << G4endl;
1370 collectHead = h.fslices[i];
1371 collectHeadNo = i;
1372 }
1373 else
1374 {
1375 s << "As slice #" << collectHeadNo << G4endl;
1376 }
1377 }
1378 }
1379
1380 if (haveHeaders)
1381 {
1382 collectHead=0;
1383 for (j=0; j<h.fslices.size(); j++)
1384 {
1385 if (h.fslices[j]->IsHeader())
1386 {
1387 s << "Header at Slice #" << j << " = ";
1388 if (h.fslices[j] != collectHead)
1389 {
1390 s << G4endl
1391 << (*(h.fslices[j]->GetHeader()));
1392 collectHead = h.fslices[j];
1393 collectHeadNo = j;
1394 }
1395 else
1396 {
1397 s << "As slice #" << collectHeadNo << G4endl;
1398 }
1399 }
1400 }
1401 }
1402 return s;
1403}
Note: See TracBrowser for help on using the repository browser.