source: trunk/examples/advanced/composite_calorimeter/src/CCalMaterialFactory.cc @ 1304

Last change on this file since 1304 was 807, checked in by garnier, 16 years ago

update

File size: 12.3 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer                                           *
4// *                                                                  *
5// * The  Geant4 software  is  copyright of the Copyright Holders  of *
6// * the Geant4 Collaboration.  It is provided  under  the terms  and *
7// * conditions of the Geant4 Software License,  included in the file *
8// * LICENSE and available at  http://cern.ch/geant4/license .  These *
9// * include a list of copyright holders.                             *
10// *                                                                  *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work  make  any representation or  warranty, express or implied, *
14// * regarding  this  software system or assume any liability for its *
15// * use.  Please see the license in the file  LICENSE  and URL above *
16// * for the full disclaimer and the limitation of liability.         *
17// *                                                                  *
18// * This  code  implementation is the result of  the  scientific and *
19// * technical work of the GEANT4 collaboration.                      *
20// * By using,  copying,  modifying or  distributing the software (or *
21// * any work based  on the software)  you  agree  to acknowledge its *
22// * use  in  resulting  scientific  publications,  and indicate your *
23// * acceptance of all terms of the Geant4 Software license.          *
24// ********************************************************************
25//
26///////////////////////////////////////////////////////////////////////////////
27// File: CCalMaterialFactory.cc
28// Description: CCalMaterialFactory is a factory class to vuild G4Material
29//              from CCalMaterial and CCalAmaterial
30///////////////////////////////////////////////////////////////////////////////
31#include "CCalMaterialFactory.hh"
32#include "CCalutils.hh"
33#include <fstream>
34#include <stdlib.h>
35
36#include "G4Material.hh"
37
38//#define ddebug
39//#define debug
40
41typedef CCalMaterial*  ptrCCalMaterial;
42typedef CCalAMaterial* ptrCCalAMaterial;
43
44
45CCalMaterialFactory * CCalMaterialFactory::instance = 0;
46G4String CCalMaterialFactory::elementfile = "";
47G4String CCalMaterialFactory::mixturefile = "";
48
49
50CCalMaterialFactory* CCalMaterialFactory::getInstance(const G4String& matfile,
51                                                      const G4String& mixfile){
52  if ((matfile=="" || matfile==elementfile) &&
53      (mixfile=="" || mixfile==mixturefile))
54    return getInstance();
55  else if ((matfile != "" && elementfile != "" && matfile != elementfile) ||
56           (mixfile != "" && mixturefile != "" && mixfile != mixturefile)) {
57    G4cerr << "ERROR: Trying to get materials from " << matfile << " and " 
58         << mixfile << " while previously were retrieved from " 
59         << elementfile << " and " << mixturefile << "." << G4endl;
60    return 0;
61  } else {
62    if (elementfile == "") 
63      elementfile=matfile;
64    if (mixturefile == "") 
65      mixturefile=mixfile;
66    return getInstance();
67  }
68}
69
70
71CCalMaterialFactory* CCalMaterialFactory::getInstance(const G4String& matfile){
72  return getInstance(matfile,matfile);
73}
74
75
76CCalMaterialFactory* CCalMaterialFactory::getInstance(){
77  if (elementfile=="" || mixturefile=="") {
78    G4cerr << "ERROR: You haven't defined files to be used for materials in "
79         << "CCalMaterialFactory::getInstance(const G4String&,const G4String&)"
80         << G4endl;
81    return 0;
82  }
83
84  if (instance==0) {
85    instance = new CCalMaterialFactory;
86    return instance;
87  }
88  else
89    return instance;
90}
91
92
93CCalMaterialFactory::~CCalMaterialFactory(){
94  CCalMaterialTable::iterator ite;
95  for(ite = theCCalMaterials.begin(); ite != theCCalMaterials.end(); ite++ ){
96    delete *ite;
97  }
98  theCCalMaterials.clear();
99  CCalAMaterialTable::iterator itea;
100  for(itea = theCCalAMaterials.begin(); itea != theCCalAMaterials.end(); 
101      itea++ ){
102    delete *itea;
103  }
104  theCCalAMaterials.clear();
105}
106
107 
108G4Material* CCalMaterialFactory::findMaterial(const G4String & mat) const {
109  G4Material* theMat=findG4Material(mat);
110
111  if (theMat) {
112#ifdef ddebug
113    G4cout << "Material " << mat << " already defined. Returning previous "
114         << "instance." << G4endl;
115#endif
116    return theMat;
117  } else { 
118    CCalMaterial* CCalmat=findCCalMaterial(mat);
119    if (CCalmat){
120      G4Material* G4Mat = new G4Material(CCalmat->Name(),
121                                         CCalmat->Density()*g/cm3, 
122                                         CCalmat->NElements());
123      for(G4int i=0; i<CCalmat->NElements(); i++) {
124        G4Element* elem = findElement(CCalmat->Element(i));
125        if (!elem) {
126          G4cerr << "       Could not build material " << mat << "." << G4endl;
127          exit(-10);
128        }
129        G4Mat->AddElement(elem, CCalmat->Weight(i));
130      }
131#ifdef ddebug
132    G4cout << "Material " << mat << " has been built successfully." << G4endl;
133#endif
134      return G4Mat;
135    } else {
136      G4cerr << "ERROR: Material " << mat << " not found in CCal database!!!" 
137           << G4endl;
138      return 0;
139    }
140  }
141}
142
143
144G4Element* CCalMaterialFactory::findElement(const G4String & mat) const {
145  const G4ElementTable  theElements = *(G4Element::GetElementTable());
146  for (unsigned int i=0; i<theElements.size(); i++)
147    if (theElements[i]->GetName()==mat){
148#ifdef ddebug
149      G4cout << "Element " << mat << " found!" << G4endl;
150#endif
151      return theElements[i];
152    }
153  return 0;
154}
155
156
157G4Element* CCalMaterialFactory::addElement(const G4String & name,
158                                           const G4String & symbol,
159                                           G4double Z, G4double A,
160                                           G4double density) {
161
162  G4Element* theEl = new G4Element(name, symbol, Z, A*g/mole);
163  //Make it also as a material.
164  CCalAMaterial* theMat = new CCalAMaterial(name,A,density);
165  theCCalAMaterials.push_back(theMat);
166
167#ifdef ddebug
168  G4cout << "Element " << name << " created!" << G4endl;
169#endif
170  return theEl;
171}
172
173
174G4Material* CCalMaterialFactory::addMaterial(const G4String& name,
175                                             G4double density,
176                                             G4int nconst,
177                                             G4String mats[],
178                                             G4double prop[],
179                                             MatDescription md){
180  addCCalMaterial(name, density, nconst, mats, prop, md);
181  return findMaterial(name);
182}
183
184
185void CCalMaterialFactory::readElements(const G4String& matfile) {
186
187  G4String path = getenv("CCAL_GLOBALPATH");
188  G4cout << " ==> Opening file " << matfile << " to read elements..." << G4endl;
189  std::ifstream is;
190  bool ok = openGeomFile(is, path, matfile);
191  if (!ok) {
192    G4cerr << "ERROR: Could not open file " << matfile << G4endl;
193    return;
194  }
195
196  // Find *DO GMAT
197  findDO(is, G4String("GMAT"));
198 
199  readElements(is);
200 
201  is.close();
202}
203
204
205void CCalMaterialFactory::readMaterials(const G4String& matfile) {
206
207  G4String path = getenv("CCAL_GLOBALPATH");
208  G4cout << " ==> Opening file " << matfile << " to read materials..." << G4endl;
209  std::ifstream is;
210  bool ok = openGeomFile(is, path, matfile);
211  if (!ok) {
212    G4cerr << "ERROR: Could not open file " << matfile << G4endl;
213    return;
214  }
215
216  // Find *DO GMIX
217  findDO(is, G4String("GMIX"));
218
219  readMaterials(is);
220
221  is.close();
222}
223
224
225//===========================================================================
226// Protected & private methods ==============================================
227
228
229G4Material* CCalMaterialFactory::findG4Material(const G4String & mat) const {
230  const G4MaterialTable theG4Materials = *(G4Material::GetMaterialTable());
231  for (unsigned int i=0; i<theG4Materials.size(); i++) {
232    if (theG4Materials[i]->GetName()==mat){
233      return theG4Materials[i];
234    }
235  }
236  return 0;
237}
238
239
240CCalMaterial* CCalMaterialFactory::findCCalMaterial(const G4String & mat) 
241  const {
242  for (unsigned int i=0; i<theCCalMaterials.size(); i++)
243    if (theCCalMaterials[i]->Name()==mat){
244#ifdef ddebug
245      G4cout << "CCalMaterial " << mat << " found!" << G4endl;
246#endif
247      return theCCalMaterials[i];
248    }
249  return (CCalMaterial*) findCCalAMaterial(mat);
250}
251
252
253CCalAMaterial* CCalMaterialFactory::findCCalAMaterial(const G4String & mat) 
254  const {
255  for (unsigned int i=0; i<theCCalAMaterials.size(); i++)
256    if (theCCalAMaterials[i]->Name()==mat){
257#ifdef ddebug
258      G4cout << "CCalMaterial " << mat << " found!" << G4endl;
259#endif
260      return theCCalAMaterials[i];
261    }
262  return 0;
263}
264
265
266CCalMaterial* CCalMaterialFactory::addCCalMaterial(const G4String& name, 
267                                                   G4double density,
268                                                   G4int nconst,
269                                                   G4String mats[], 
270                                                   G4double prop[],
271                                                   MatDescription md){
272  ptrCCalMaterial* matcol=0;
273  ptrCCalAMaterial* amatcol=0;
274
275  if (md==byAtomic)
276    amatcol = new ptrCCalAMaterial[nconst];
277  else
278    matcol = new ptrCCalMaterial[nconst];
279
280  for (G4int i=0; i<nconst; i++){
281    if (md==byAtomic) {
282      CCalAMaterial* amat = findCCalAMaterial(mats[i]);
283      if (amat)
284        amatcol[i]=amat;
285      else {
286        G4cerr << "ERROR: Trying to build" << name << " out of unknown " 
287             << mats[i] << "." << G4endl
288             << "Skiping this material!" << G4endl;
289        delete[] amatcol;
290        return 0;
291      }
292    } //by Atomic fractions
293    else {
294      CCalMaterial* mat = findCCalMaterial(mats[i]);
295      if (mat)
296        matcol[i]=mat;
297      else {
298        G4cerr << "ERROR: Trying to build" <<name << " out of unknown " 
299             << mats[i] << "." << G4endl
300             << "Skiping this material!" << G4endl;
301        delete[] matcol;
302        return 0;
303      }
304    }
305  } //for
306
307  //Let's do the CCalMaterial!
308  if (md==byAtomic) {
309    CCalAMaterial* amaterial = new CCalAMaterial(name, density, nconst, 
310                                                 amatcol, prop);
311    delete[] amatcol;
312    theCCalAMaterials.push_back(amaterial);
313#ifdef ddebug
314    G4cout << *amaterial << G4endl;
315#endif
316    return amaterial;
317  } else {
318    CCalMaterial::FractionType ft;
319    if (md == byWeight)
320      ft=CCalMaterial::FTWeight;
321    else
322      ft=CCalMaterial::FTVolume;
323    CCalMaterial* material = new CCalMaterial(name, density, nconst, 
324                                            matcol, prop, ft);
325    delete[] matcol;
326    theCCalMaterials.push_back(material);
327#ifdef ddebug
328    G4cout << *material << G4endl;
329#endif
330    return material;
331  } 
332}
333
334
335void CCalMaterialFactory::readElements(std::ifstream& is){
336  G4String name, symbol;
337 
338  G4cout << "     ==> Reading elements... " << G4endl;
339#ifdef debug
340  G4cout << "       Element    \tsymbol\tA\tZ\tdensity\tX_0          abs_l"<< G4endl;
341#endif
342  //There should come the list of materials. #. Defines a comment
343  //*DO defines the beguining of the Mixes block.
344
345  readName(is,name);
346  while (name != "*ENDDO") {
347    //It should be an element definition
348    G4double A, Z, density;
349    is >> symbol >> A >> Z >> density >> jump;   
350#ifdef debug
351    G4cout << "       " << name << "    \t" << symbol << "\t" 
352         << A << "\t" << Z << "\t" << density << G4endl;
353#endif   
354    addElement(name, symbol, Z, A, density);
355    readName(is,name);
356  };
357  G4cout << "     " << G4Element::GetElementTable()->size() 
358       << " elements read from file" << G4endl << G4endl;
359}
360
361
362void CCalMaterialFactory::readMaterials(std::ifstream& is){
363  G4String name, matname;
364
365  G4cout << "     ==> Reading materials... " << G4endl;
366
367  //Take into account the special case of vacuum...
368#ifdef debug
369  G4cout <<"       \"Vacuum\"" << G4endl;
370#endif
371  G4double density  = universe_mean_density;   //from PhysicalConstants.h
372  G4double pressure = 1.E-19*pascal;
373  G4double temperature = 0.1*kelvin;
374  new G4Material("Vacuum", /*Z=*/ 1., /*A=*/ 1.01*g/mole,
375                 density, kStateGas, temperature, pressure);
376
377  //There should come the list of materials. #. Defines a comment
378  //*ENDDO defines the block.
379  readName(is,name);
380  while (name != "*ENDDO") {
381    //It should be a material definition
382    matname=name;
383    G4int nElem;
384    G4double density;
385    is >> nElem >> density >> jump;
386
387#ifdef debug
388    G4cout <<"       " << matname
389           << " made of " << nElem
390           << " elements. Density=" << density
391           << G4endl;
392#endif
393
394    G4int absnelem = std::abs(nElem);
395
396    G4String* mats    = new G4String[absnelem];
397    G4double* weights = new G4double[absnelem];
398   
399    G4double prop;
400    for(int i=0; i<absnelem; i++) {
401      readName(is, name);
402      is >> prop >> jump;
403      mats[i]=name;
404      weights[i]=std::abs(prop);
405    } //for...
406    MatDescription md;
407    if (nElem>0 && prop<0)
408      md = byAtomic;
409    else if (nElem>0)
410      md = byWeight;
411    else
412      md = byVolume;
413
414    addCCalMaterial(matname, density, absnelem, mats, weights, md);
415    delete[] mats;
416    delete[] weights;
417   
418    readName(is,name);
419  };  //while
420
421  G4cout << "     " << theCCalMaterials.size() << " materials read from " 
422       << mixturefile << G4endl << G4endl;
423}
424
425
426CCalMaterialFactory::CCalMaterialFactory() {
427  readElements (elementfile);
428  readMaterials(mixturefile);
429}
Note: See TracBrowser for help on using the repository browser.