| 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 // G4AdjointCrossSurfChecker class implementat << 26 // $Id: G4AdjointCrossSurfChecker.cc 66872 2013-01-15 01:25:57Z japost $
27 // 27 //
28 // Author: L. Desorgher, SpaceIT GmbH << 28 /////////////////////////////////////////////////////////////////////////////
29 // Contract: ESA contract 21435/08/NL/AT << 29 // Class Name: G4AdjointCrossSurfChecker
30 // Customer: ESA/ESTEC << 30 // Author: L. Desorgher
31 // ------------------------------------------- << 31 // Organisation: SpaceIT GmbH
>> 32 // Contract: ESA contract 21435/08/NL/AT
>> 33 // Customer: ESA/ESTEC
>> 34 /////////////////////////////////////////////////////////////////////////////
32 35
33 #include "G4AdjointCrossSurfChecker.hh" 36 #include "G4AdjointCrossSurfChecker.hh"
34 <<
35 #include "G4AffineTransform.hh" <<
36 #include "G4PhysicalConstants.hh" 37 #include "G4PhysicalConstants.hh"
37 #include "G4PhysicalVolumeStore.hh" << 38 #include "G4SystemOfUnits.hh"
38 #include "G4Step.hh" 39 #include "G4Step.hh"
39 #include "G4StepPoint.hh" 40 #include "G4StepPoint.hh"
40 #include "G4SystemOfUnits.hh" << 41 #include "G4PhysicalVolumeStore.hh"
41 #include "G4VSolid.hh" 42 #include "G4VSolid.hh"
>> 43 #include "G4AffineTransform.hh"
42 44
43 ////////////////////////////////////////////// 45 //////////////////////////////////////////////////////////////////////////////
44 // << 46 //
45 G4ThreadLocal G4AdjointCrossSurfChecker* G4Adj << 47 G4ThreadLocal G4AdjointCrossSurfChecker* G4AdjointCrossSurfChecker::instance = 0;
46 48
47 ////////////////////////////////////////////// 49 //////////////////////////////////////////////////////////////////////////////
48 // 50 //
49 G4AdjointCrossSurfChecker::~G4AdjointCrossSurf << 51 G4AdjointCrossSurfChecker::G4AdjointCrossSurfChecker()
50 << 52 {;
51 ////////////////////////////////////////////// << 53 }
>> 54 ///////////////////////////////////////////////////////////////////////////////
>> 55 //
>> 56 G4AdjointCrossSurfChecker::~G4AdjointCrossSurfChecker()
>> 57 {
>> 58 delete instance;
>> 59 }
>> 60 ////////////////////////////////////////////////////////////////////////////////
52 // 61 //
53 G4AdjointCrossSurfChecker* G4AdjointCrossSurfC 62 G4AdjointCrossSurfChecker* G4AdjointCrossSurfChecker::GetInstance()
54 { 63 {
55 if (instance == nullptr) instance = new G4Ad << 64 if (!instance) instance = new G4AdjointCrossSurfChecker();
56 return instance; 65 return instance;
57 } << 66 }
58 << 67 /////////////////////////////////////////////////////////////////////////////////
59 ////////////////////////////////////////////// <<
60 // 68 //
61 G4bool G4AdjointCrossSurfChecker::CrossingASph << 69 G4bool G4AdjointCrossSurfChecker::CrossingASphere(const G4Step* aStep,G4double sphere_radius, G4ThreeVector sphere_center,G4ThreeVector& crossing_pos, G4double& cos_th , G4bool& GoingIn)
62 G4ThreeVector sphere_center, G4ThreeVector& <<
63 { 70 {
64 G4ThreeVector pos1 = aStep->GetPreStepPoint( << 71 G4ThreeVector pos1= aStep->GetPreStepPoint()->GetPosition() - sphere_center;
65 G4ThreeVector pos2 = aStep->GetPostStepPoint << 72 G4ThreeVector pos2= aStep->GetPostStepPoint()->GetPosition() - sphere_center;
66 G4double r1 = pos1.mag(); << 73 G4double r1= pos1.mag();
67 G4double r2 = pos2.mag(); << 74 G4double r2= pos2.mag();
68 G4bool did_cross = false; << 75 G4bool did_cross =false;
69 << 76
70 if (r1 <= sphere_radius && r2 > sphere_radiu << 77 if (r1<=sphere_radius && r2>sphere_radius){
71 did_cross = true; << 78 did_cross=true;
72 GoingIn = false; << 79 GoingIn=false;
73 } << 80 }
74 else if (r2 <= sphere_radius && r1 > sphere_ << 81 else if (r2<=sphere_radius && r1>sphere_radius){
75 did_cross = true; << 82 did_cross=true;
76 GoingIn = true; << 83 GoingIn=true;
77 } << 84 }
78 << 85
79 if (did_cross) { << 86 if (did_cross) {
80 G4ThreeVector dr = pos2 - pos1; << 87
81 G4double r12 = r1 * r1; << 88 G4ThreeVector dr=pos2-pos1;
82 G4double rdr = dr.mag(); << 89 G4double r12 = r1*r1;
83 G4double a, b, c, d; << 90 G4double rdr = dr.mag();
84 a = rdr * rdr; << 91 G4double a,b,c,d;
85 b = 2. * pos1.dot(dr); << 92 a = rdr*rdr;
86 c = r12 - sphere_radius * sphere_radius; << 93 b = 2.*pos1.dot(dr);
87 d = std::sqrt(b * b - 4. * a * c); << 94 c = r12-sphere_radius*sphere_radius;
88 G4double l = (-b + d) / 2. / a; << 95 d=std::sqrt(b*b-4.*a*c);
89 if (l > 1.) l = (-b - d) / 2. / a; << 96 G4double l= (-b+d)/2./a;
90 crossing_pos = pos1 + l * dr; << 97 if (l > 1.) l=(-b-d)/2./a;
91 cos_th = std::abs(dr.cosTheta(crossing_pos << 98 crossing_pos=pos1+l*dr;
>> 99 cos_th = std::abs(dr.cosTheta(crossing_pos));
>> 100
>> 101
92 } 102 }
93 return did_cross; 103 return did_cross;
94 } 104 }
95 << 105 /////////////////////////////////////////////////////////////////////////////////
96 ////////////////////////////////////////////// <<
97 // 106 //
98 G4bool G4AdjointCrossSurfChecker::GoingInOrOut << 107 G4bool G4AdjointCrossSurfChecker::GoingInOrOutOfaVolume(const G4Step* aStep,const G4String& volume_name, G4double& , G4bool& GoingIn) //from external surface
99 const G4String& volume_name, G4double&, G4bo <<
100 { 108 {
101 G4bool step_at_boundary = (aStep->GetPostSte 109 G4bool step_at_boundary = (aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary);
102 G4bool did_cross = false; << 110 G4bool did_cross =false;
103 if (step_at_boundary) { << 111 if (step_at_boundary){
104 const G4VTouchable* postStepTouchable = aS << 112 const G4VTouchable* postStepTouchable = aStep->GetPostStepPoint()->GetTouchable();
105 const G4VTouchable* preStepTouchable = aSt << 113 const G4VTouchable* preStepTouchable = aStep->GetPreStepPoint()->GetTouchable();
106 if ((preStepTouchable != nullptr) && (post << 114 if (preStepTouchable && postStepTouchable && postStepTouchable->GetVolume() && preStepTouchable->GetVolume()){
107 (postStepTouchable->GetVolume() != nul << 115 G4String post_vol_name = postStepTouchable->GetVolume()->GetName();
108 { << 116 G4String pre_vol_name = preStepTouchable->GetVolume()->GetName();
109 G4String post_vol_name = postStepTouchab << 117
110 G4String pre_vol_name = preStepTouchable << 118 if (post_vol_name == volume_name ){
111 << 119 GoingIn=true;
112 if (post_vol_name == volume_name) { << 120 did_cross=true;
113 GoingIn = true; << 121 }
114 did_cross = true; << 122 else if (pre_vol_name == volume_name){
115 } << 123 GoingIn=false;
116 else if (pre_vol_name == volume_name) { << 124 did_cross=true;
117 GoingIn = false; << 125
118 did_cross = true; << 126 }
119 } << 127
120 } << 128 }
121 } 129 }
122 return did_cross; // still need to compute << 130 return did_cross;
>> 131 //still need to compute the cosine of the direction
123 } 132 }
124 << 133 /////////////////////////////////////////////////////////////////////////////////
125 ////////////////////////////////////////////// <<
126 // 134 //
127 G4bool G4AdjointCrossSurfChecker::GoingInOrOut << 135 G4bool G4AdjointCrossSurfChecker::GoingInOrOutOfaVolumeByExtSurface(const G4Step* aStep,const G4String& volume_name, const G4String& mother_logical_vol_name, G4double& , G4bool& GoingIn) //from external surface
128 const G4String& volume_name, const G4String& <<
129 G4bool& GoingIn) // from external surf. <<
130 { 136 {
131 G4bool step_at_boundary = (aStep->GetPostSte 137 G4bool step_at_boundary = (aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary);
132 G4bool did_cross = false; << 138 G4bool did_cross =false;
133 if (step_at_boundary) { << 139 if (step_at_boundary){
134 const G4VTouchable* postStepTouchable = aS << 140 const G4VTouchable* postStepTouchable = aStep->GetPostStepPoint()->GetTouchable();
135 const G4VTouchable* preStepTouchable = aSt << 141 const G4VTouchable* preStepTouchable = aStep->GetPreStepPoint()->GetTouchable();
136 const G4VPhysicalVolume* postVol = << 142 if (preStepTouchable && postStepTouchable && postStepTouchable->GetVolume() && preStepTouchable->GetVolume()){
137 (postStepTouchable != nullptr) ? postSte << 143 G4String post_vol_name = postStepTouchable->GetVolume()->GetName();
138 const G4VPhysicalVolume* preVol = << 144 G4String post_log_vol_name = postStepTouchable->GetVolume()->GetLogicalVolume()->GetName();
139 (preStepTouchable != nullptr) ? preStepT << 145 G4String pre_vol_name = preStepTouchable->GetVolume()->GetName();
140 if (preStepTouchable != nullptr && postSte << 146 G4String pre_log_vol_name = preStepTouchable->GetVolume()->GetLogicalVolume()->GetName();
141 preVol != nullptr) << 147 if (post_vol_name == volume_name && pre_log_vol_name == mother_logical_vol_name){
142 { << 148 GoingIn=true;
143 G4String post_vol_name = postVol->GetNam << 149 did_cross=true;
144 G4String post_log_vol_name = postVol->Ge << 150 }
145 G4String pre_vol_name = preVol->GetName( << 151 else if (pre_vol_name == volume_name && post_log_vol_name == mother_logical_vol_name ){
146 G4String pre_log_vol_name = preVol->GetL << 152 GoingIn=false;
147 if (post_vol_name == volume_name && pre_ << 153 did_cross=true;
148 GoingIn = true; << 154
149 did_cross = true; << 155 }
150 } << 156
151 else if (pre_vol_name == volume_name && << 157 }
152 GoingIn = false; <<
153 did_cross = true; <<
154 } <<
155 } <<
156 } <<
157 return did_cross; // still need to compute <<
158 } <<
159 <<
160 ////////////////////////////////////////////// <<
161 // <<
162 G4bool G4AdjointCrossSurfChecker::CrossingAGiv <<
163 const G4String& surface_name, G4ThreeVector& <<
164 G4bool& GoingIn) <<
165 { <<
166 G4int ind = FindRegisteredSurface(surface_na <<
167 G4bool did_cross = false; <<
168 if (ind >= 0) { <<
169 did_cross = CrossingAGivenRegisteredSurfac <<
170 } 158 }
171 return did_cross; 159 return did_cross;
>> 160 //still need to compute the cosine of the direction
172 } 161 }
173 << 162 /////////////////////////////////////////////////////////////////////////////////
174 ////////////////////////////////////////////// <<
175 // 163 //
176 G4bool G4AdjointCrossSurfChecker::CrossingAGiv << 164 G4bool G4AdjointCrossSurfChecker::CrossingAGivenRegisteredSurface(const G4Step* aStep,const G4String& surface_name,G4ThreeVector& crossing_pos, G4double& cos_to_surface, G4bool& GoingIn)
177 G4ThreeVector& crossing_pos, G4double& cos_t <<
178 { 165 {
179 G4String surf_type = ListOfSurfaceType[ind]; << 166 G4int ind = FindRegisteredSurface(surface_name);
180 G4double radius = ListOfSphereRadius[ind]; <<
181 G4ThreeVector center = ListOfSphereCenter[in <<
182 G4String vol1 = ListOfVol1Name[ind]; <<
183 G4String vol2 = ListOfVol2Name[ind]; <<
184 <<
185 G4bool did_cross = false; 167 G4bool did_cross = false;
186 if (surf_type == "Sphere") { << 168 if (ind >=0){
187 did_cross = CrossingASphere(aStep, radius, << 169 did_cross = CrossingAGivenRegisteredSurface(aStep, ind, crossing_pos,cos_to_surface, GoingIn);
188 } <<
189 else if (surf_type == "ExternalSurfaceOfAVol <<
190 did_cross = GoingInOrOutOfaVolumeByExtSurf <<
191 crossing_pos = aStep->GetPostStepPoint()-> <<
192 } <<
193 else if (surf_type == "BoundaryBetweenTwoVol <<
194 did_cross = CrossingAnInterfaceBetweenTwoV <<
195 aStep, vol1, vol2, crossing_pos, cos_to_ <<
196 } 170 }
197 return did_cross; 171 return did_cross;
198 } 172 }
199 << 173 /////////////////////////////////////////////////////////////////////////////////
200 ////////////////////////////////////////////// <<
201 // 174 //
202 G4bool G4AdjointCrossSurfChecker::CrossingOneO << 175 G4bool G4AdjointCrossSurfChecker::CrossingAGivenRegisteredSurface(const G4Step* aStep, int ind,G4ThreeVector& crossing_pos, G4double& cos_to_surface, G4bool& GoingIn)
203 G4String& surface_name, G4ThreeVector& cross <<
204 { 176 {
205 for (std::size_t i = 0; i < ListOfSurfaceNam << 177 G4String surf_type = ListOfSurfaceType[ind];
206 if (CrossingAGivenRegisteredSurface(aStep, << 178 G4double radius = ListOfSphereRadius[ind];
207 surface_name = ListOfSurfaceName[i]; << 179 G4ThreeVector center = ListOfSphereCenter[ind];
208 return true; << 180 G4String vol1 = ListOfVol1Name[ind];
209 } << 181 G4String vol2 = ListOfVol2Name[ind];
210 } << 182
211 return false; << 183 G4bool did_cross = false;
>> 184 if (surf_type == "Sphere"){
>> 185 did_cross = CrossingASphere(aStep, radius, center,crossing_pos, cos_to_surface, GoingIn);
>> 186 }
>> 187 else if (surf_type == "ExternalSurfaceOfAVolume"){
>> 188
>> 189 did_cross = GoingInOrOutOfaVolumeByExtSurface(aStep, vol1, vol2, cos_to_surface, GoingIn);
>> 190 crossing_pos= aStep->GetPostStepPoint()->GetPosition();
>> 191
>> 192 }
>> 193 else if (surf_type == "BoundaryBetweenTwoVolumes"){
>> 194 did_cross = CrossingAnInterfaceBetweenTwoVolumes(aStep, vol1, vol2,crossing_pos, cos_to_surface, GoingIn);
>> 195 }
>> 196 return did_cross;
>> 197
>> 198
>> 199 }
>> 200 /////////////////////////////////////////////////////////////////////////////////
>> 201 //
>> 202 //
>> 203 G4bool G4AdjointCrossSurfChecker::CrossingOneOfTheRegisteredSurface(const G4Step* aStep,G4String& surface_name,G4ThreeVector& crossing_pos, G4double& cos_to_surface, G4bool& GoingIn)
>> 204 {
>> 205 for (size_t i=0;i <ListOfSurfaceName.size();i++){
>> 206 if (CrossingAGivenRegisteredSurface(aStep, int(i),crossing_pos, cos_to_surface, GoingIn)){
>> 207 surface_name = ListOfSurfaceName[i];
>> 208 return true;
>> 209 }
>> 210 }
>> 211 return false;
212 } 212 }
213 << 213 /////////////////////////////////////////////////////////////////////////////////
214 ////////////////////////////////////////////// <<
215 // 214 //
216 G4bool G4AdjointCrossSurfChecker::CrossingAnIn << 215 G4bool G4AdjointCrossSurfChecker::CrossingAnInterfaceBetweenTwoVolumes(const G4Step* aStep,const G4String& vol1_name,const G4String& vol2_name,G4ThreeVector& , G4double& , G4bool& GoingIn)
217 const G4String& vol1_name, const G4String& v <<
218 { 216 {
219 G4bool step_at_boundary = (aStep->GetPostSte 217 G4bool step_at_boundary = (aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary);
220 G4bool did_cross = false; << 218 G4bool did_cross =false;
221 if (step_at_boundary) { << 219 if (step_at_boundary){
222 const G4VTouchable* postStepTouchable = aS << 220 const G4VTouchable* postStepTouchable = aStep->GetPostStepPoint()->GetTouchable();
223 const G4VTouchable* preStepTouchable = aSt << 221 const G4VTouchable* preStepTouchable = aStep->GetPreStepPoint()->GetTouchable();
224 if ((preStepTouchable != nullptr) && (post << 222 if (preStepTouchable && postStepTouchable){
225 G4String post_vol_name = postStepTouchab << 223
226 if (post_vol_name.empty()) { << 224 G4String post_vol_name = postStepTouchable->GetVolume()->GetName();
227 post_vol_name = postStepTouchable->Get << 225 if (post_vol_name =="") post_vol_name = postStepTouchable->GetVolume()->GetLogicalVolume()->GetName();
228 } << 226 G4String pre_vol_name = preStepTouchable->GetVolume()->GetName();
229 G4String pre_vol_name = preStepTouchable << 227 if (pre_vol_name =="") pre_vol_name = preStepTouchable->GetVolume()->GetLogicalVolume()->GetName();
230 if (pre_vol_name.empty()) { << 228
231 pre_vol_name = preStepTouchable->GetVo << 229
232 } << 230 if ( pre_vol_name == vol1_name && post_vol_name == vol2_name){
233 if (pre_vol_name == vol1_name && post_vo << 231 GoingIn=true;
234 GoingIn = true; << 232 did_cross=true;
235 did_cross = true; << 233 }
236 } << 234 else if (pre_vol_name == vol2_name && post_vol_name == vol1_name){
237 else if (pre_vol_name == vol2_name && po << 235 GoingIn=false;
238 GoingIn = false; << 236 did_cross=true;
239 did_cross = true; << 237 }
240 } << 238
241 } << 239 }
242 } 240 }
243 return did_cross; // still need to compute << 241 return did_cross;
>> 242 //still need to compute the cosine of the direction
244 } 243 }
245 244
246 ////////////////////////////////////////////// << 245 /////////////////////////////////////////////////////////////////////////////////
247 // << 246 //
248 G4bool G4AdjointCrossSurfChecker::AddaSpherica << 247 G4bool G4AdjointCrossSurfChecker::AddaSphericalSurface(const G4String& SurfaceName, G4double radius, G4ThreeVector pos, G4double& Area)
249 const G4String& SurfaceName, G4double radius <<
250 { 248 {
251 G4int ind = FindRegisteredSurface(SurfaceNam 249 G4int ind = FindRegisteredSurface(SurfaceName);
252 Area = 4. * pi * radius * radius; << 250 Area= 4.*pi*radius*radius;
253 if (ind >= 0) { << 251 if (ind>=0) {
254 ListOfSurfaceType[ind] = "Sphere"; << 252 ListOfSurfaceType[ind]="Sphere";
255 ListOfSphereRadius[ind] = radius; << 253 ListOfSphereRadius[ind]=radius;
256 ListOfSphereCenter[ind] = pos; << 254 ListOfSphereCenter[ind]=pos;
257 ListOfVol1Name[ind] = ""; << 255 ListOfVol1Name[ind]="";
258 ListOfVol2Name[ind] = ""; << 256 ListOfVol2Name[ind]="";
259 AreaOfSurface[ind] = Area; << 257 AreaOfSurface[ind]=Area;
260 } 258 }
261 else { 259 else {
262 ListOfSurfaceName.push_back(SurfaceName); << 260 ListOfSurfaceName.push_back(SurfaceName);
263 ListOfSurfaceType.emplace_back("Sphere"); << 261 ListOfSurfaceType.push_back("Sphere");
264 ListOfSphereRadius.push_back(radius); << 262 ListOfSphereRadius.push_back(radius);
265 ListOfSphereCenter.push_back(pos); << 263 ListOfSphereCenter.push_back(pos);
266 ListOfVol1Name.emplace_back(""); << 264 ListOfVol1Name.push_back("");
267 ListOfVol2Name.emplace_back(""); << 265 ListOfVol2Name.push_back("");
268 AreaOfSurface.push_back(Area); << 266 AreaOfSurface.push_back(Area);
269 } << 267 }
270 return true; 268 return true;
271 } 269 }
272 << 270 /////////////////////////////////////////////////////////////////////////////////
273 ////////////////////////////////////////////// << 271 //
274 // << 272 G4bool G4AdjointCrossSurfChecker::AddaSphericalSurfaceWithCenterAtTheCenterOfAVolume(const G4String& SurfaceName, G4double radius, const G4String& volume_name, G4ThreeVector & center, G4double& area)
275 G4bool G4AdjointCrossSurfChecker::AddaSpherica << 273 {
276 const G4String& SurfaceName, G4double radius << 274
277 G4double& area) << 275 G4VPhysicalVolume* thePhysicalVolume = 0;
278 { << 276 G4PhysicalVolumeStore* thePhysVolStore =G4PhysicalVolumeStore::GetInstance();
279 G4VPhysicalVolume* thePhysicalVolume = nullp << 277 for ( unsigned int i=0; i< thePhysVolStore->size();i++){
280 G4PhysicalVolumeStore* thePhysVolStore = G4P << 278 if ((*thePhysVolStore)[i]->GetName() == volume_name){
281 thePhysicalVolume = thePhysVolStore->GetVolu << 279 thePhysicalVolume = (*thePhysVolStore)[i];
282 if (thePhysicalVolume != nullptr) { << 280 };
283 G4VPhysicalVolume* daughter = thePhysicalV << 281
284 G4LogicalVolume* mother = thePhysicalVolum << 282 }
285 G4AffineTransform theTransformationFromPhy << 283 if (thePhysicalVolume){
286 while (mother != nullptr) { << 284 G4VPhysicalVolume* daughter =thePhysicalVolume;
287 theTransformationFromPhysVolToWorld *= << 285 G4LogicalVolume* mother = thePhysicalVolume->GetMotherLogical();
288 G4AffineTransform(daughter->GetFrameRo << 286 G4AffineTransform theTransformationFromPhysVolToWorld = G4AffineTransform();
289 for (std::size_t i = 0; i < thePhysVolSt << 287 while (mother){
290 if ((*thePhysVolStore)[i]->GetLogicalV << 288 theTransformationFromPhysVolToWorld *=
291 daughter = (*thePhysVolStore)[i]; << 289 G4AffineTransform(daughter->GetFrameRotation(),daughter->GetObjectTranslation());
292 mother = daughter->GetMotherLogical( << 290 /*G4cout<<"Mother "<<mother->GetName()<<std::endl;
293 break; << 291 G4cout<<"Daughter "<<daughter->GetName()<<std::endl;
294 } << 292 G4cout<<daughter->GetObjectTranslation()<<std::endl;
295 } << 293 G4cout<<theTransformationFromPhysVolToWorld.NetTranslation()<<std::endl;*/
296 } << 294 for ( unsigned int i=0; i< thePhysVolStore->size();i++){
297 center = theTransformationFromPhysVolToWor << 295 if ((*thePhysVolStore)[i]->GetLogicalVolume() == mother){
298 G4cout << "Center of the spherical surface << 296 daughter = (*thePhysVolStore)[i];
299 << G4endl; << 297 mother =daughter->GetMotherLogical();
>> 298 break;
>> 299 };
>> 300 }
>> 301
>> 302 }
>> 303 center = theTransformationFromPhysVolToWorld.NetTranslation();
>> 304 G4cout<<"Center of the spherical surface is at the position: "<<center/cm<<" cm"<<std::endl;
>> 305
300 } 306 }
301 else { 307 else {
302 return false; << 308 G4cout<<"The physical volume with name "<<volume_name<<" does not exist!!"<<std::endl;
>> 309 return false;
303 } 310 }
304 return AddaSphericalSurface(SurfaceName, rad 311 return AddaSphericalSurface(SurfaceName, radius, center, area);
305 } <<
306 312
307 ////////////////////////////////////////////// << 313 }
>> 314 /////////////////////////////////////////////////////////////////////////////////
308 // 315 //
309 G4bool G4AdjointCrossSurfChecker::AddanExtSurf << 316 G4bool G4AdjointCrossSurfChecker::AddanExtSurfaceOfAvolume(const G4String& SurfaceName, const G4String& volume_name, G4double& Area)
310 const G4String& SurfaceName, const G4String& <<
311 { 317 {
312 G4int ind = FindRegisteredSurface(SurfaceNam 318 G4int ind = FindRegisteredSurface(SurfaceName);
313 319
314 G4VPhysicalVolume* thePhysicalVolume = nullp << 320 G4VPhysicalVolume* thePhysicalVolume = 0;
315 G4PhysicalVolumeStore* thePhysVolStore = G4P << 321 G4PhysicalVolumeStore* thePhysVolStore =G4PhysicalVolumeStore::GetInstance();
316 thePhysicalVolume = thePhysVolStore->GetVolu << 322 for ( unsigned int i=0; i< thePhysVolStore->size();i++){
317 if (thePhysicalVolume == nullptr) { << 323 if ((*thePhysVolStore)[i]->GetName() == volume_name){
318 return false; << 324 thePhysicalVolume = (*thePhysVolStore)[i];
319 } << 325 };
>> 326
>> 327 }
>> 328 if (!thePhysicalVolume){
>> 329 G4cout<<"The physical volume with name "<<volume_name<<" does not exist!!"<<std::endl;
>> 330 return false;
>> 331 }
320 Area = thePhysicalVolume->GetLogicalVolume() 332 Area = thePhysicalVolume->GetLogicalVolume()->GetSolid()->GetSurfaceArea();
321 G4String mother_vol_name = ""; << 333 G4String mother_vol_name = "";
322 G4LogicalVolume* theMother = thePhysicalVolu 334 G4LogicalVolume* theMother = thePhysicalVolume->GetMotherLogical();
323 << 335
324 if (theMother != nullptr) mother_vol_name = << 336 if (theMother) mother_vol_name= theMother->GetName();
325 if (ind >= 0) { << 337 if (ind>=0) {
326 ListOfSurfaceType[ind] = "ExternalSurfaceO << 338 ListOfSurfaceType[ind]="ExternalSurfaceOfAVolume";
327 ListOfSphereRadius[ind] = 0.; << 339 ListOfSphereRadius[ind]=0.;
328 ListOfSphereCenter[ind] = G4ThreeVector(0. << 340 ListOfSphereCenter[ind]=G4ThreeVector(0.,0.,0.);
329 ListOfVol1Name[ind] = volume_name; << 341 ListOfVol1Name[ind]=volume_name;
330 ListOfVol2Name[ind] = std::move(mother_vol << 342 ListOfVol2Name[ind]=mother_vol_name;
331 AreaOfSurface[ind] = Area; << 343 AreaOfSurface[ind]=Area;
332 } 344 }
333 else { 345 else {
334 ListOfSurfaceName.push_back(SurfaceName); << 346 ListOfSurfaceName.push_back(SurfaceName);
335 ListOfSurfaceType.emplace_back("ExternalSu << 347 ListOfSurfaceType.push_back("ExternalSurfaceOfAVolume");
336 ListOfSphereRadius.push_back(0.); << 348 ListOfSphereRadius.push_back(0.);
337 ListOfSphereCenter.emplace_back(0., 0., 0. << 349 ListOfSphereCenter.push_back(G4ThreeVector(0.,0.,0.));
338 ListOfVol1Name.push_back(volume_name); << 350 ListOfVol1Name.push_back(volume_name);
339 ListOfVol2Name.push_back(std::move(mother_ << 351 ListOfVol2Name.push_back(mother_vol_name);
340 AreaOfSurface.push_back(Area); << 352 AreaOfSurface.push_back(Area);
341 } 353 }
342 return true; 354 return true;
343 } 355 }
344 << 356 /////////////////////////////////////////////////////////////////////////////////
345 ////////////////////////////////////////////// <<
346 // 357 //
347 G4bool G4AdjointCrossSurfChecker::AddanInterfa << 358 G4bool G4AdjointCrossSurfChecker::AddanInterfaceBetweenTwoVolumes(const G4String& SurfaceName, const G4String& volume_name1, const G4String& volume_name2,G4double& Area)
348 const G4String& volume_name1, const G4String <<
349 { 359 {
350 G4int ind = FindRegisteredSurface(SurfaceNam 360 G4int ind = FindRegisteredSurface(SurfaceName);
351 Area = -1.; // the way to compute the surfa << 361 Area=-1.; //the way to compute the surface is not known yet
352 if (ind >= 0) { << 362 if (ind>=0) {
353 ListOfSurfaceType[ind] = "BoundaryBetweenT << 363 ListOfSurfaceType[ind]="BoundaryBetweenTwoVolumes";
354 ListOfSphereRadius[ind] = 0.; << 364 ListOfSphereRadius[ind]=0.;
355 ListOfSphereCenter[ind] = G4ThreeVector(0. << 365 ListOfSphereCenter[ind]=G4ThreeVector(0.,0.,0.);
356 ListOfVol1Name[ind] = volume_name1; << 366 ListOfVol1Name[ind]=volume_name1;
357 ListOfVol2Name[ind] = volume_name2; << 367 ListOfVol2Name[ind]=volume_name2;
358 AreaOfSurface[ind] = Area; << 368 AreaOfSurface[ind]=Area;
>> 369
359 } 370 }
360 else { 371 else {
361 ListOfSurfaceName.push_back(SurfaceName); << 372 ListOfSurfaceName.push_back(SurfaceName);
362 ListOfSurfaceType.emplace_back("BoundaryBe << 373 ListOfSurfaceType.push_back("BoundaryBetweenTwoVolumes");
363 ListOfSphereRadius.push_back(0.); << 374 ListOfSphereRadius.push_back(0.);
364 ListOfSphereCenter.emplace_back(0., 0., 0. << 375 ListOfSphereCenter.push_back(G4ThreeVector(0.,0.,0.));
365 ListOfVol1Name.push_back(volume_name1); << 376 ListOfVol1Name.push_back(volume_name1);
366 ListOfVol2Name.push_back(volume_name2); << 377 ListOfVol2Name.push_back(volume_name2);
367 AreaOfSurface.push_back(Area); << 378 AreaOfSurface.push_back(Area);
368 } 379 }
369 return true; 380 return true;
370 } 381 }
371 << 382 /////////////////////////////////////////////////////////////////////////////////
372 ////////////////////////////////////////////// <<
373 // 383 //
374 void G4AdjointCrossSurfChecker::ClearListOfSel 384 void G4AdjointCrossSurfChecker::ClearListOfSelectedSurface()
375 { 385 {
376 ListOfSurfaceName.clear(); 386 ListOfSurfaceName.clear();
377 ListOfSurfaceType.clear(); 387 ListOfSurfaceType.clear();
378 ListOfSphereRadius.clear(); 388 ListOfSphereRadius.clear();
379 ListOfSphereCenter.clear(); 389 ListOfSphereCenter.clear();
380 ListOfVol1Name.clear(); 390 ListOfVol1Name.clear();
381 ListOfVol2Name.clear(); 391 ListOfVol2Name.clear();
382 } 392 }
383 << 393 /////////////////////////////////////////////////////////////////////////////////
384 ////////////////////////////////////////////// <<
385 // 394 //
386 G4int G4AdjointCrossSurfChecker::FindRegistere 395 G4int G4AdjointCrossSurfChecker::FindRegisteredSurface(const G4String& name)
387 { 396 {
388 for (std::size_t i = 0; i < ListOfSurfaceNam << 397 G4int ind=-1;
389 if (name == ListOfSurfaceName[i]) return G << 398 for (size_t i = 0; i<ListOfSurfaceName.size();i++){
390 } << 399 if (name == ListOfSurfaceName[i]) {
391 return -1; << 400 ind = int (i);
>> 401 return ind;
>> 402 }
>> 403 }
>> 404 return ind;
392 } 405 }
>> 406
393 407