| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // class G4PVPlacement Implementation
27 //
28 // 24.07.95 P.Kent, First non-stub version.
29 // -------------------------------------------
30
31 #include "G4PVPlacement.hh"
32 #include "G4AffineTransform.hh"
33 #include "G4UnitsTable.hh"
34 #include "G4LogicalVolume.hh"
35 #include "G4VSolid.hh"
36
37 // -------------------------------------------
38 // Constructor
39 //
40 G4PVPlacement::G4PVPlacement( G4RotationMatrix
41 const G4ThreeVector& t
42 const G4String& pName,
43 G4LogicalVolume*
44 G4VPhysicalVolum
45 G4bool pMany,
46 G4int pCopyNo,
47 G4bool pSurfChk
48 : G4VPhysicalVolume(pRot, tlate, pName, pLog
49 fmany(pMany), fcopyNo(pCopyNo)
50 {
51 if (pMother != nullptr)
52 {
53 G4LogicalVolume* motherLogical = pMother->
54 if (pLogical == motherLogical)
55 {
56 G4Exception("G4PVPlacement::G4PVPlacemen
57 FatalException, "Cannot plac
58 }
59 SetMotherLogical(motherLogical);
60 motherLogical->AddDaughter(this);
61 if (pSurfChk) { CheckOverlaps(); }
62 }
63 }
64
65 // -------------------------------------------
66 // Constructor
67 //
68 G4PVPlacement::G4PVPlacement( const G4Transfor
69 const G4String&
70 G4LogicalV
71 G4VPhysica
72 G4bool pMa
73 G4int pCop
74 G4bool pSu
75 : G4VPhysicalVolume(NewPtrRotMatrix(Transfor
76 Transform3D.getTranslati
77 fmany(pMany), fcopyNo(pCopyNo)
78 {
79 fallocatedRotM = (GetRotation() != nullptr);
80 if (pMother != nullptr)
81 {
82 G4LogicalVolume* motherLogical = pMother->
83 if (pLogical == motherLogical)
84 G4Exception("G4PVPlacement::G4PVPlacemen
85 FatalException, "Cannot plac
86 SetMotherLogical(motherLogical);
87 motherLogical->AddDaughter(this);
88 if (pSurfChk) { CheckOverlaps(); }
89 }
90 }
91
92 // -------------------------------------------
93 // Constructor
94 //
95 // The logical volume of the mother is utilise
96 //
97 G4PVPlacement::G4PVPlacement( G4RotationMatrix
98 const G4ThreeVector& t
99 G4LogicalVolume*
100 const G4String& pName,
101 G4LogicalVolume*
102 G4bool pMany,
103 G4int pCopyNo,
104 G4bool pSurfChk
105 : G4VPhysicalVolume(pRot, tlate, pName, pCur
106 fmany(pMany), fcopyNo(pCopyNo)
107 {
108 if (pCurrentLogical == pMotherLogical)
109 {
110 G4Exception("G4PVPlacement::G4PVPlacement(
111 FatalException, "Cannot place
112 }
113 SetMotherLogical(pMotherLogical);
114 if (pMotherLogical != nullptr) { pMotherLogi
115 if ((pSurfChk) && ((pMotherLogical) != nullp
116 }
117
118 // -------------------------------------------
119 // Constructor
120 //
121 G4PVPlacement::G4PVPlacement( const G4Transfor
122 G4LogicalV
123 const G4String&
124 G4LogicalV
125 G4bool pMa
126 G4int pCop
127 G4bool pSu
128 : G4VPhysicalVolume(nullptr, Transform3D.get
129 pName, pCurrentLogical,
130 fmany(pMany), fcopyNo(pCopyNo)
131 {
132 if (pCurrentLogical == pMotherLogical)
133 {
134 G4Exception("G4PVPlacement::G4PVPlacement(
135 FatalException, "Cannot place
136 }
137 SetRotation( NewPtrRotMatrix(Transform3D.get
138 fallocatedRotM = (GetRotation() != nullptr);
139 SetMotherLogical(pMotherLogical);
140 if (pMotherLogical != nullptr) { pMotherLogi
141 if ((pSurfChk) && ((pMotherLogical) != nullp
142 }
143
144 // -------------------------------------------
145 // Fake default constructor - sets only member
146 // for usage restri
147 //
148 G4PVPlacement::G4PVPlacement( __void__& a )
149 : G4VPhysicalVolume(a)
150 {
151 }
152
153 // -------------------------------------------
154 // Destructor
155 //
156 G4PVPlacement::~G4PVPlacement()
157 {
158 if( fallocatedRotM ){ delete this->GetRotati
159 }
160
161 // -------------------------------------------
162 // IsMany
163 //
164 G4bool G4PVPlacement::IsMany() const
165 {
166 return fmany;
167 }
168
169 // -------------------------------------------
170 // SetCopyNo
171 //
172 void G4PVPlacement::SetCopyNo(G4int newCopyNo)
173 {
174 fcopyNo = newCopyNo;
175 }
176
177 // -------------------------------------------
178 // IsReplicated
179 //
180 G4bool G4PVPlacement::IsReplicated() const
181 {
182 return false;
183 }
184
185 // -------------------------------------------
186 // IsParameterised
187 //
188 G4bool G4PVPlacement::IsParameterised() const
189 {
190 return false;
191 }
192
193 // -------------------------------------------
194 // GetParameterisation
195 //
196 G4VPVParameterisation* G4PVPlacement::GetParam
197 {
198 return nullptr;
199 }
200
201 // -------------------------------------------
202 // GetReplicationData
203 //
204 void G4PVPlacement::
205 GetReplicationData( EAxis&, G4int&, G4double&,
206 {
207 // No-operations
208 }
209
210 // -------------------------------------------
211 // IsRegularRepeatedStructure
212 //
213 // This is for specialised repeated volumes (r
214 //
215 G4bool G4PVPlacement::IsRegularStructure() con
216 {
217 return false;
218 }
219
220 // -------------------------------------------
221 // IsRegularRepeatedStructure
222 //
223 // This is for specialised repeated volumes (r
224 //
225 G4int G4PVPlacement::GetRegularStructureId() c
226 {
227 return 0;
228 }
229
230 // -------------------------------------------
231 // VolumeType
232 //
233 // Information to help identify sub-navigator
234 //
235 EVolume G4PVPlacement::VolumeType() const
236 {
237 return kNormal;
238 }
239
240 // -------------------------------------------
241 // CheckOverlaps
242 //
243 G4bool G4PVPlacement::CheckOverlaps(G4int res,
244 G4bool ver
245 {
246 if (res <= 0) { return false; }
247
248 G4VSolid* solid = GetLogicalVolume()->GetSol
249 G4LogicalVolume* motherLog = GetMotherLogica
250 if (motherLog == nullptr) { return false; }
251
252 G4int trials = 0;
253 G4bool retval = false;
254
255 if (verbose)
256 {
257 G4cout << "Checking overlaps for volume "
258 << GetName() << ':' << GetCopyNo()
259 << " (" << solid->GetEntityType() <
260 }
261
262 // Check that random points are gererated co
263 //
264 G4ThreeVector ptmp = solid->GetPointOnSurfac
265 if (solid->Inside(ptmp) != kSurface)
266 {
267 G4String position[3] = { "outside", "surfa
268 std::ostringstream message;
269 message << "Sample point is not on the sur
270 << " The issue is detecte
271 << GetName() << ':' << GetCopyNo()
272 << " (" << solid->GetEntityType()
273 << " generated point " <<
274 << " is " << position[solid->Insid
275 G4Exception("G4PVPlacement::CheckOverlaps(
276 "GeomVol1002", JustWarning, me
277 return false;
278 }
279
280 // Generate random points on the surface of
281 // transform them into the mother volume coo
282 // and find the bonding box
283 //
284 std::vector<G4ThreeVector> points(res);
285 G4double xmin = kInfinity, ymin = kInfinit
286 G4double xmax = -kInfinity, ymax = -kInfinit
287 G4AffineTransform Tm(GetRotation(), GetTrans
288 for (G4int i = 0; i < res; ++i)
289 {
290 points[i] = Tm.TransformPoint(solid->GetPo
291 xmin = std::min(xmin, points[i].x());
292 ymin = std::min(ymin, points[i].y());
293 zmin = std::min(zmin, points[i].z());
294 xmax = std::max(xmax, points[i].x());
295 ymax = std::max(ymax, points[i].y());
296 zmax = std::max(zmax, points[i].z());
297 }
298 G4ThreeVector scenter(0.5*(xmax+xmin), 0.5*(
299 G4double sradius = 0.5*G4ThreeVector(xmax-xm
300
301 // Check overlap with the mother volume
302 //
303 G4int overlapCount = 0;
304 G4double overlapSize = -kInfinity;
305 G4ThreeVector overlapPoint;
306 G4VSolid* motherSolid = motherLog->GetSolid(
307 for (G4int i = 0; i < res; ++i)
308 {
309 G4ThreeVector mp = points[i];
310 if (motherSolid->Inside(mp) != kOutside) c
311 G4double distin = motherSolid->DistanceToI
312 if (distin < tol) continue; // too small o
313 ++overlapCount;
314 if (distin <= overlapSize) continue;
315 overlapSize = distin;
316 overlapPoint = mp;
317 }
318
319 // Print information on overlap
320 //
321 if (overlapCount > 0)
322 {
323 ++trials;
324 retval = true;
325 std::ostringstream message;
326 message << "Overlap with mother volume !"
327 << " Overlap is detected
328 << GetName() << ':' << GetCopyNo()
329 << " (" << solid->GetEntityType()
330 << " with its mother volume " << m
331 << " (" << motherSolid->GetEntityT
332 << " protrusion at mother
333 << " by " << G4BestUnit(overlapSiz
334 << " (max of " << overlapCount <<
335 if (trials >= maxErr)
336 {
337 message << G4endl
338 << "NOTE: Reached maximum fixed
339 << "- of overlaps reports for th
340 }
341 G4Exception("G4PVPlacement::CheckOverlaps(
342 "GeomVol1002", JustWarning, me
343 if (trials >= maxErr) { return true; }
344 }
345
346 // Checking overlaps with each 'sister' volu
347 //
348 G4VSolid* previous = nullptr;
349 G4ThreeVector pmin_local(0.,0.,0.), pmax_loc
350
351 for (std::size_t k = 0; k < motherLog->GetNo
352 {
353 G4VPhysicalVolume* daughter = motherLog->G
354 if (daughter == this) continue;
355 G4bool check_encapsulation = true;
356
357 G4AffineTransform Td(daughter->GetRotation
358 G4VSolid* daughterSolid = daughter->GetLog
359 if (previous != daughterSolid)
360 {
361 daughterSolid->BoundingLimits(pmin_local
362 previous = daughterSolid;
363 }
364 overlapCount = 0;
365 overlapSize = -kInfinity;
366 if (!Td.IsRotated()) { // no rotation, onl
367 G4ThreeVector offset = Td.NetTranslation
368 G4ThreeVector pmin(pmin_local + offset);
369 G4ThreeVector pmax(pmax_local + offset);
370 if (pmin.x() >= xmax) continue;
371 if (pmin.y() >= ymax) continue;
372 if (pmin.z() >= zmax) continue;
373 if (pmax.x() <= xmin) continue;
374 if (pmax.y() <= ymin) continue;
375 if (pmax.z() <= zmin) continue;
376 for (G4int i = 0; i < res; ++i)
377 {
378 G4ThreeVector p = points[i];
379 if (p.x() <= pmin.x()) continue;
380 if (p.x() >= pmax.x()) continue;
381 if (p.y() <= pmin.y()) continue;
382 if (p.y() >= pmax.y()) continue;
383 if (p.z() <= pmin.z()) continue;
384 if (p.z() >= pmax.z()) continue;
385 G4ThreeVector md = p - offset;
386 if (daughterSolid->Inside(md) == kInsi
387 {
388 check_encapsulation = false;
389 G4double distout = daughterSolid->Di
390 if (distout < tol) continue; // too
391 ++overlapCount;
392 if (distout <= overlapSize) continue
393 overlapSize = distout;
394 overlapPoint = md;
395 }
396 }
397 }
398 else // transformation with rotation
399 {
400 G4ThreeVector pmin(pmin_local), pmax(pma
401 G4ThreeVector dcenter = Td.TransformPoin
402 G4double dradius = 0.5*((pmax - pmin).ma
403 if ((scenter - dcenter).mag2() >= (sradi
404 if (dcenter.x() - dradius >= xmax) conti
405 if (dcenter.y() - dradius >= ymax) conti
406 if (dcenter.z() - dradius >= zmax) conti
407 if (dcenter.x() + dradius <= xmin) conti
408 if (dcenter.y() + dradius <= ymin) conti
409 if (dcenter.z() + dradius <= zmin) conti
410
411 G4ThreeVector pbox[8] = {
412 G4ThreeVector(pmin.x(), pmin.y(), pmin
413 G4ThreeVector(pmax.x(), pmin.y(), pmin
414 G4ThreeVector(pmin.x(), pmax.y(), pmin
415 G4ThreeVector(pmax.x(), pmax.y(), pmin
416 G4ThreeVector(pmin.x(), pmin.y(), pmax
417 G4ThreeVector(pmax.x(), pmin.y(), pmax
418 G4ThreeVector(pmin.x(), pmax.y(), pmax
419 G4ThreeVector(pmax.x(), pmax.y(), pmax
420 };
421 G4double dxmin = kInfinity, dymin = kI
422 G4double dxmax = -kInfinity, dymax = -kI
423 for (const auto & i : pbox)
424 {
425 G4ThreeVector p = Td.TransformPoint(i)
426 dxmin = std::min(dxmin, p.x());
427 dymin = std::min(dymin, p.y());
428 dzmin = std::min(dzmin, p.z());
429 dxmax = std::max(dxmax, p.x());
430 dymax = std::max(dymax, p.y());
431 dzmax = std::max(dzmax, p.z());
432 }
433 if (dxmin >= xmax) continue;
434 if (dymin >= ymax) continue;
435 if (dzmin >= zmax) continue;
436 if (dxmax <= xmin) continue;
437 if (dymax <= ymin) continue;
438 if (dzmax <= zmin) continue;
439 for (G4int i = 0; i < res; ++i)
440 {
441 G4ThreeVector p = points[i];
442 if (p.x() >= dxmax) continue;
443 if (p.x() <= dxmin) continue;
444 if (p.y() >= dymax) continue;
445 if (p.y() <= dymin) continue;
446 if (p.z() >= dzmax) continue;
447 if (p.z() <= dzmin) continue;
448 G4ThreeVector md = Td.InverseTransform
449 if (daughterSolid->Inside(md) == kInsi
450 {
451 check_encapsulation = false;
452 G4double distout = daughterSolid->Di
453 if (distout < tol) continue; // too
454 ++overlapCount;
455 if (distout <= overlapSize) continue
456 overlapSize = distout;
457 overlapPoint = md;
458 }
459 }
460 }
461
462 // Print information on overlap
463 //
464 if (overlapCount > 0)
465 {
466 ++trials;
467 retval = true;
468 std::ostringstream message;
469 message << "Overlap with volume already
470 << " Overlap is detecte
471 << GetName() << ':' << GetCopyNo
472 << " (" << solid->GetEntityType(
473 << daughter->GetName() << ':' <<
474 << " (" << daughterSolid->GetEnt
475 << " overlap at local p
476 << " by " << G4BestUnit(overlapS
477 << " (max of " << overlapCount <
478 if (trials >= maxErr)
479 {
480 message << G4endl
481 << "NOTE: Reached maximum fixe
482 << "- of overlaps reports for
483 }
484 G4Exception("G4PVPlacement::CheckOverlap
485 "GeomVol1002", JustWarning,
486 if (trials >= maxErr) { return true; }
487 }
488 else if (check_encapsulation)
489 {
490 // Now checking that 'sister' volume is
491 // and overlapping. Generate a single po
492 // the 'sister' volume and verify that t
493 // the current volume
494 //
495 G4ThreeVector pSurface = daughterSolid->
496 G4ThreeVector normal = daughterSolid->Su
497 G4ThreeVector pInside = pSurface - norma
498 G4ThreeVector dPoint = (daughterSolid->I
499 pInside : pSurface;
500
501 // Transform the generated point to the
502 // and then to current volume's coordina
503 //
504 G4ThreeVector mp2 = Td.TransformPoint(dP
505 G4ThreeVector msi = Tm.InverseTransformP
506
507 if (solid->Inside(msi) == kInside)
508 {
509 ++trials;
510 retval = true;
511 std::ostringstream message;
512 message << "Overlap with volume alread
513 << " Overlap is detec
514 << GetName() << ':' << GetCopy
515 << " (" << solid->GetEntityTyp
516 << " apparently fully
517 << daughter->GetName() << ':'
518 << " (" << daughterSolid->GetE
519 << " at the same level!";
520 if (trials >= maxErr)
521 {
522 message << G4endl
523 << "NOTE: Reached maximum fi
524 << "- of overlaps reports fo
525 }
526 G4Exception("G4PVPlacement::CheckOverl
527 "GeomVol1002", JustWarning
528 if (trials >= maxErr) { return true;
529 }
530 }
531 }
532
533 if (verbose && trials == 0) { G4cout << "OK!
534 return retval;
535 }
536
537 // -------------------------------------------
538 // NewPtrRotMatrix
539 //
540 // Auxiliary function for 2nd & 4th constructo
541 // Creates a new rotation matrix on the heap (
542 // argument into it.
543 //
544 // NOTE: Ownership of the returned pointer is
545 // No entity is currently responsible to
546 //
547 G4RotationMatrix*
548 G4PVPlacement::NewPtrRotMatrix(const G4Rotatio
549 {
550 G4RotationMatrix* pRotMatrix;
551 if ( RotMat.isIdentity() )
552 {
553 pRotMatrix = nullptr;
554 }
555 else
556 {
557 pRotMatrix = new G4RotationMatrix(RotMat)
558 }
559 return pRotMatrix;
560 }
561