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