4 // Author: I. Hrivnacova
6 // Class TG4GeometryServices
7 // -------------------------
8 // See the class description in the header file.
10 #include "TG4GeometryServices.h"
11 #include "TG4Limits.h"
12 #include "TG4IntMap.h"
13 #include "TG4NameMap.h"
14 #include "TG4G3Units.h"
15 #include "TG4G3ControlVector.h"
16 #include "TG4Globals.h"
18 #include <G4LogicalVolumeStore.hh>
19 #include <G4LogicalVolume.hh>
20 #include <G4Material.hh>
21 #include <G4Element.hh>
22 #include <G4UserLimits.hh>
24 #include <G3EleTable.hh>
25 #include <g4std/vector>
29 TG4GeometryServices* TG4GeometryServices::fgInstance = 0;
30 const G4double TG4GeometryServices::fgkAZTolerance = 0.001;
31 const G4double TG4GeometryServices::fgkDensityTolerance = 0.005;
33 //_____________________________________________________________________________
34 TG4GeometryServices::TG4GeometryServices(TG4IntMap* mediumMap,
36 : fMediumMap(mediumMap),
41 TG4Globals::Exception(
42 "TG4GeometryServices: attempt to create two instances of singleton.");
48 //_____________________________________________________________________________
49 TG4GeometryServices::TG4GeometryServices() {
51 TG4Globals::Exception(
52 "TG4GeometryServices default constructor is protected.");
56 //_____________________________________________________________________________
57 TG4GeometryServices::TG4GeometryServices(const TG4GeometryServices& right) {
59 TG4Globals::Exception(
60 "Attempt to copy TG4GeometryServices singleton.");
64 //_____________________________________________________________________________
65 TG4GeometryServices::~TG4GeometryServices() {
69 //=============================================================================
73 //=============================================================================
75 //_____________________________________________________________________________
77 TG4GeometryServices::operator=(const TG4GeometryServices& right)
79 // check assignement to self
80 if (this == &right) return *this;
82 TG4Globals::Exception(
83 "Attempt to assign TG4GeometryServices singleton.");
88 //=============================================================================
92 //=============================================================================
94 //_____________________________________________________________________________
95 G4bool TG4GeometryServices::IsG3Volume(const G4String& lvName) const
97 // Returns true if the logical volume of given volumeName
98 // was not created by Gsposp method with a generic name
102 if (lvName.contains(gSeparator))
108 //_____________________________________________________________________________
109 G4bool TG4GeometryServices::CompareElement(G4double a, G4double z,
110 const G4Element* element) const
112 // Compares given parameters with those of a given element,
113 // returns true if they are equal, false otherwise.
116 G4double ae = element->GetA()/TG4G3Units::AtomicWeight();
117 G4double ze = element->GetZ();
119 // g3tog4 can redefine A
120 G4double ax = G3Ele.GetEle(z)->GetA()/TG4G3Units::AtomicWeight();
122 if (abs(ax - ae) < fgkAZTolerance && abs(z - ze) < fgkAZTolerance)
128 //_____________________________________________________________________________
129 G4bool TG4GeometryServices::CompareMaterial(G4int nofElements, G4double density,
130 const G4Material* material) const
132 // Compares given density with those of a given material,
133 // returns true if they are equal, false otherwise.
136 G4double dm = material->GetDensity()/TG4G3Units::MassDensity();
137 G4int ne = material->GetNumberOfElements();
139 // density percentual difference
140 G4double diff = abs(density - dm)/(density + dm)*2.;
142 if (nofElements == ne && diff < fgkDensityTolerance)
148 //_____________________________________________________________________________
149 G4double* TG4GeometryServices::ConvertAtomWeight(G4int nmat,
150 G4double* a, G4double* wmat) const
152 // In case of proportions given in atom counts (nmat<0),
153 // the wmat[i] are converted to weight fractions.
154 // (From g3tog4 G4gsmixt.)
155 // The new array has to be delete by client.
158 G4double* weight = new G4double[abs(nmat)];
163 for (i=0; i<abs(nmat); i++) {
164 // total molecular weight
165 aMol += wmat[i]*a[i];
168 G4String text = "TG4GeometryServices::ConvertAtomWeight:\n";
169 text = text + " Total molecular weight = 0";
170 TG4Globals::Warning(text);
172 for (G4int i=0; i<abs(nmat); i++) {
174 weight[i] = wmat[i]*a[i]/aMol;
178 for (G4int j=0; j<nmat; j++) weight[j] = wmat[j];
183 //=============================================================================
187 //=============================================================================
189 //_____________________________________________________________________________
190 G4double* TG4GeometryServices::CreateG4doubleArray(Float_t* array,
193 // Converts Float_t* array to G4double*,
194 // !! The new array has to be deleted by user.
197 G4double* doubleArray;
199 doubleArray = new G4double[size];
200 for (G4int i=0; i<size; i++) doubleArray[i] = array[i];
208 //_____________________________________________________________________________
209 G4String TG4GeometryServices::CutName(const char* name) const
211 // Removes spaces after the name if present.
214 G4String cutName = name;
215 G4int i = cutName.length();
216 while (cutName(--i) == ' ') cutName = cutName(0,i);
221 //_____________________________________________________________________________
222 G4String TG4GeometryServices::CutMaterialName(const char* name) const
224 // Removes the $ with precedent spaces at the name if present.
227 G4String cutName = name;
228 cutName = cutName.substr(0,cutName.find('$'));
230 return CutName(cutName);
233 //_____________________________________________________________________________
234 G4String TG4GeometryServices::G4ToG3VolumeName(const G4String& name) const
236 // Cuts _copyNo extension added to logical volume name in case
237 // the logical volume was created by Gsposp method.
240 G4String cutName = name;
241 if (cutName.contains(gSeparator))
242 cutName = cutName(0,cutName.first(gSeparator));
247 //_____________________________________________________________________________
248 G4String TG4GeometryServices::GenerateLimitsName(G4int id,
249 const G4String& medName,
250 const G4String& matName) const
252 // Generate unique name for user limits composed from the tracking medium id,
253 // name and its material name.
256 TG4Globals::AppendNumberToString(name, id);
257 name = name + "__med_" + medName + "__mat_" + matName;
262 //_____________________________________________________________________________
263 G4Material* TG4GeometryServices::MixMaterials(G4String name, G4double density,
264 const TG4StringVector& matNames,
265 const TG4doubleVector& matWeights)
267 // Creates a mixture of selected materials
270 // number of materials to be mixed
271 G4int nofMaterials = matNames.size();
272 if (nofMaterials != matWeights.size()) {
273 G4String text = "TG4GeometryServices::MixMaterials: ";
274 text = text + "different number of material names and weigths.";
275 TG4Globals::Exception(text);
278 // G4cout << "Nof of materials to be mixed: " << nofMaterials << G4endl;
280 // fill vector of materials
281 G4std::vector<G4Material*> matVector;
283 for (im=0; im< nofMaterials; im++) {
285 G4Material* material = G4Material::GetMaterial(matNames[im]);
286 matVector.push_back(material);
289 // create the mixed material
290 G4Material* mixture = new G4Material(name, density, nofMaterials);
291 for (im=0; im< nofMaterials; im++) {
292 G4Material* material = matVector[im];
293 G4double fraction = matWeights[im];
294 mixture->AddMaterial(material, fraction);
300 //_____________________________________________________________________________
301 void TG4GeometryServices::PrintNameMap() const
303 // Prints the map of volumes names to second names.
306 fNameMap->PrintAll();
309 //_____________________________________________________________________________
310 void TG4GeometryServices::PrintLimits(const G4String& name) const
312 // Finds the limits with the specified name and prints them.
315 TG4Limits* limits = FindLimits(name, true);
317 if (limits) limits->Print();
320 //_____________________________________________________________________________
321 void TG4GeometryServices::PrintVolumeLimits(const G4String& volumeName) const
323 // Finds a logical volume with the specified name and prints
327 G4LogicalVolume* lv = FindLogicalVolume(volumeName, false);
330 TG4Limits* limits = GetLimits(lv->GetUserLimits());
331 G4cout << volumeName << " ";
335 G4cout << "has not the limits set." << G4endl;
339 //_____________________________________________________________________________
340 void TG4GeometryServices::PrintStatistics(G4bool open, G4bool close) const
342 // Print G4 geometry statistics
346 if (open) TG4Globals::PrintStars(true);
348 G4cout << " GEANT4 Geometry statistics: " << G4endl
349 << " " << NofG4LogicalVolumes()
350 << " logical volumes" << G4endl
352 << NofG4PhysicalVolumes()
353 << " physical volumes" << G4endl
355 << G4Material::GetNumberOfMaterials()
356 << " materials" << G4endl
358 << TG4Limits::GetNofLimits()
359 << " user limits" << G4endl;
361 if (close) TG4Globals::PrintStars(false);
364 //_____________________________________________________________________________
365 Int_t TG4GeometryServices::NofG3Volumes() const
367 // Returns the total number of logical volumes corresponding
369 // The logical volume that were created by Gsposp method
370 // with a generic name (name_copyNo) are NOT included.
373 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
376 for (G4int i=0; i<lvStore->size(); i++) {
377 G4LogicalVolume* lv = (*lvStore)[i];
378 if (IsG3Volume(lv->GetName())) counter++;
384 //_____________________________________________________________________________
385 Int_t TG4GeometryServices::NofG4LogicalVolumes() const
387 // Returns the total number of logical volumes in the geometry.
390 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
391 return lvStore->size();
394 //_____________________________________________________________________________
395 Int_t TG4GeometryServices::NofG4PhysicalVolumes() const
397 // Returns the total number of physical volumes in the geometry.
400 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
403 for (G4int i=0; i<lvStore->size(); i++) {
404 counter += ((*lvStore)[i])->GetNoDaughters();
411 //______________________________________________________________________________
412 G4bool TG4GeometryServices::IsSpecialControls() const
414 // Returns true if a process control in some limits instance is set.
417 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
419 for (G4int i=0; i<lvStore->size(); i++) {
420 TG4Limits* limits = GetLimits((*lvStore)[i]->GetUserLimits());
421 if (limits && limits->IsControl()) return true;
427 //_____________________________________________________________________________
428 TG4Limits* TG4GeometryServices::GetLimits(G4UserLimits* limits) const
430 // Checks and converts type of the given limits.
433 if (!limits) return 0;
435 TG4Limits* tg4Limits = dynamic_cast<TG4Limits*> (limits);
438 G4Exception("TG4GeometryServices::GetLimits: Wrong limits type.");
445 //_____________________________________________________________________________
446 const G4String& TG4GeometryServices::GetMapSecond(const G4String& name)
448 // Returns the second string associated with the name in
452 return fNameMap->GetSecond(name);
455 //_____________________________________________________________________________
457 TG4GeometryServices::FindLogicalVolume(const G4String& name, G4bool silent) const
459 // Finds a logical volume with the specified name in G4LogicalVolumeStore.
462 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
464 for (G4int i=0; i<lvStore->size(); i++) {
465 G4LogicalVolume* lv = (*lvStore)[i];
466 if (lv->GetName() == name) return lv;
470 G4String text = "TG4GeometryServices:: FindLogicalVolume:\n";
471 text = text + " Logical volume " + name + " does not exist.";
472 TG4Globals::Warning(text);
478 //_____________________________________________________________________________
480 TG4GeometryServices::FindLimits(const G4String& name, G4bool silent) const
482 // Finds limits with the specified name.
485 G4LogicalVolumeStore* lvStore = G4LogicalVolumeStore::GetInstance();
487 for (G4int i=0; i<lvStore->size(); i++) {
488 G4LogicalVolume* lv = (*lvStore)[i];
489 TG4Limits* limits = GetLimits(lv->GetUserLimits());
490 if (limits && limits->GetName() == name) return limits;
494 G4String text = "TG4GeometryServices:: FindLimits:\n";
495 text = text + " Logical volume " + name + " does not exist.";
496 TG4Globals::Warning(text);
501 //_____________________________________________________________________________
502 G4int TG4GeometryServices::GetMediumId(G4LogicalVolume* lv) const
504 // Returns the second index for materials (having its origin in
505 // G4 tracking media concept)
508 return fMediumMap->GetSecond(lv->GetName());
511 //_____________________________________________________________________________
512 G4double TG4GeometryServices::GetEffA(G4Material* material) const
514 // Returns A or effective A=sum(pi*Ai) (if compound/mixture)
515 // of given material.
519 G4int nofElements = material->GetNumberOfElements();
520 if (nofElements > 1) {
521 G4String text = "Effective A for material mixture (";
522 text = text + material->GetName();
523 text = text + ") is used.";
524 //TG4Globals::Warning(text);
526 for (G4int i=0; i<nofElements; i++) {
527 G4double aOfElement = material->GetElement(i)->GetA();
528 G4double massFraction = material->GetFractionVector()[i];
529 a += aOfElement*massFraction /(TG4G3Units::AtomicWeight());
533 a = material->GetA();
534 a /= TG4G3Units::AtomicWeight();
539 //_____________________________________________________________________________
540 G4double TG4GeometryServices::GetEffZ(G4Material* material) const
542 // Returns Z or effective Z=sum(pi*Zi) (if compound/mixture)
543 // of given material.
547 G4int nofElements = material->GetNumberOfElements();
548 if (nofElements > 1) {
549 G4String text = "Effective Z for material mixture (";
550 text = text + material->GetName();
551 text = text + ") is used.";
552 //TG4Globals::Warning(text);
554 for (G4int i=0; i<nofElements; i++) {
555 G4double zOfElement = material->GetElement(i)->GetZ();
556 G4double massFraction = material->GetFractionVector()[i];
557 z += zOfElement*massFraction;
561 z = material->GetZ();
566 //_____________________________________________________________________________
567 G4Material* TG4GeometryServices::FindMaterial(G4double a, G4double z,
568 G4double density) const
570 // Finds material in G4MaterialTable with specified parameters,
571 // returns 0 if not found.
574 const G4MaterialTable* kpMatTable = G4Material::GetMaterialTable();
576 for (G4int i=0; i<G4Material::GetNumberOfMaterials(); i++) {
578 G4Material* material = (*kpMatTable)[i];
580 if (CompareElement(a, z, material->GetElement(0)) &&
581 CompareMaterial(1, density, material))
589 //_____________________________________________________________________________
590 G4Material* TG4GeometryServices::FindMaterial(G4double* a, G4double* z,
592 G4int nmat, G4double* wmat) const
594 // Finds material in G4MaterialTable with specified parameters,
595 // returns 0 if not found.
598 G4double* weight = ConvertAtomWeight(nmat, a, wmat);
600 // loop over materials
601 G4Material* found = 0;
602 for (G4int i=0; i<G4Material::GetNumberOfMaterials(); i++) {
604 G4Material* material = (*G4Material::GetMaterialTable())[i];
605 G4int nm = material->GetNumberOfElements();
607 if (CompareMaterial(nm, density, material)) {
609 // loop over elements
611 for (G4int ie=0; ie<nm; ie++) {
613 G4double we = (material->GetFractionVector())[ie];
615 if (!CompareElement(a[ie], z[ie], material->GetElement(ie)) ||
616 abs(weight[ie] - we) > fgkAZTolerance) {