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 67009 2013-01-29 16:00:21Z gcosmo $ 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 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