composite_calorimeter/src/CCalMaterial.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 ///////////////////////////////////////////////////////////////////////////////
0027 // File: CCalMaterial.cc
0028 // Description: CCalMaterial transient class to store information from
0029 //              material table (database) which is used to make a G4Material
0030 ///////////////////////////////////////////////////////////////////////////////
0031 #include "CCalMaterial.hh"
0032
0033 //#define debug
0034
0035 CCalMaterial::CCalMaterial(G4String mat, G4double dens, G4int nconst,
0036                            CCalMaterial** constituents, G4double* weights,
0037                            FractionType ft): name(mat), density(dens) {
0038   nElem = 0;
0039
0040   G4int i=0;
0041   for (i=0; i<nconst; i++)
0042     nElem += constituents[i]->NElements();
0043
0044   theElements = new G4String[nElem];
0045   theWeights  = new G4double[nElem];
0046
0047   G4double factor;
0048   G4int nelem=0;
0049   for (i=0; i<nconst; i++) {
0050     if (ft==FTWeight)
0051       factor=1.0;
0052     else
0053       factor=constituents[i]->Density();
0054     for (G4int j=0; j<constituents[i]->NElements(); j++) {
0055       theElements[nelem] = constituents[i]->Element(j);
0056       theWeights[nelem]  = constituents[i]->Weight(j)* weights[i] * factor;
0057       nelem++;
0058     }
0059   }
0060
0061   if (density<0) { //Let's compute density
0062     computeDensity(nconst, constituents, weights, ft);
0063   }
0064   closeMaterial();
0065 }
0066
0067 CCalMaterial::CCalMaterial(const CCalMaterial& mat):
0068   name(mat.name), density(mat.density), nElem(mat.nElem) {
0069   theElements = new G4String[nElem];
0070   theWeights  = new G4double[nElem];
0071   for (G4int i=0; i<nElem; i++){
0072     theElements[i]=mat.theElements[i];
0073     theWeights[i]=mat.theWeights[i];
0074   }
0075 }
0076
0077 CCalMaterial::~CCalMaterial() {
0078   if (theElements)
0079     delete[] theElements;
0080   if (theWeights)
0081     delete[] theWeights;
0082 }
0083
0084 void CCalMaterial::computeDensity(G4int nconst,
0085                                   CCalMaterial** constituents,
0086                                   G4double* weights, FractionType ft) {
0087   G4double mass=0;
0088   G4double volume=0;
0089   for (G4int i=0; i<nconst; i++) {
0090     if (ft==FTWeight) {
0091       mass+=weights[i];
0092       volume+=(weights[i]/constituents[i]->Density());
0093     }
0094     else { //by volume
0095       mass+=(weights[i]*constituents[i]->Density());
0096       volume+=weights[i];
0097     }
0098   }
0099   density=mass/volume;
0100 }
0101
0102 CCalMaterial& CCalMaterial::operator=(const CCalMaterial& mat){
0103   if(theElements)
0104     delete[] theElements;
0105   if(theWeights)
0106     delete[] theWeights;
0107
0108   name=mat.name;
0109   density=mat.density;
0110   nElem=mat.nElem;
0111
0112   theElements = new G4String[nElem];
0113   theWeights  = new G4double[nElem];
0114   for (G4int i=0; i<nElem; i++){
0115     theElements[i]=mat.theElements[i];
0116     theWeights[i]=mat.theWeights[i];
0117   }
0118   return *this;
0119 }
0120
0121 G4bool CCalMaterial::operator==(const CCalMaterial& mat) const{
0122   return (name==mat.name);
0123 }
0124
0125 G4bool CCalMaterial::operator!=(const CCalMaterial& mat) const{
0126   return (name!=mat.name);
0127 }
0128
0129 void CCalMaterial::closeMaterial() {
0130   G4int trueConst=0;
0131
0132   G4double norm=0;
0133
0134   for (G4int i=0; i<nElem; i++) {
0135     norm+=theWeights[i];
0136     if (theElements[i]!="") {
0137       trueConst++;
0138       for (G4int j=i+1; j<nElem; j++) {
0139         if(theElements[i]==theElements[j]){
0140           theWeights[i]+=theWeights[j];
0141           theElements[j]="";
0142         }
0143       }//for j
0144     } //if
0145   }//for i
0146
0147   if (trueConst != nElem) {
0148     G4String* newConst = new G4String[trueConst];
0149     G4double* newWeight = new G4double[trueConst];
0150
0151     G4int newi=0;
0152     for(G4int i=0; i<nElem; i++){
0153       if (theElements[i]!="") {
0154         newConst[newi]  = theElements[i];
0155         newWeight[newi] = theWeights[i]/norm;
0156         newi++;
0157       }
0158     }
0159
0160 #ifdef debug
0161     G4cout << "\tGoing from " << nElem <<" constituents to " << trueConst << G4endl;
0162 #endif
0163     nElem=trueConst;
0164
0165     delete[] theElements;
0166     delete[] theWeights;
0167
0168     theElements=newConst;
0169     theWeights=newWeight;
0170   }
0171   else { //Let's normalize the weights
0172     for (G4int i=0; i<nElem; i++)
0173       theWeights[i] = theWeights[i]/norm;
0174   }
0175 }
0176
0177 std::ostream& operator<<(std::ostream& os, const CCalMaterial& mat) {
0178   os << mat.name << G4endl;
0179   os << "Density= " << mat.density << " g/cm3. Number of Elements: "
0180      << mat.nElem << G4endl;
0181   for (G4int i=0; i<mat.nElem; i++)
0182     os << '\t' << mat.theElements[i] << '\t' << mat.theWeights[i] << G4endl;
0183   return os;
0184 }