physics/src/G4RDFluoData.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 //
0028 // Author: Elena Guardincerri (Elena.Guardincerri@ge.infn.it)
0029 //
0030 // History:
0031 // -----------
0032 // 16 Sept 2001  First committed to cvs
0033 //
0034 // -------------------------------------------------------------------
0035
0036 #include <fstream>
0037 #include <sstream>
0038
0039 #include "G4RDFluoData.hh"
0040 #include "G4SystemOfUnits.hh"
0041 #include "G4DataVector.hh"
0042 #include "G4RDFluoTransition.hh"
0043
0044 G4RDFluoData::G4RDFluoData()
0045 {
0046   numberOfVacancies=0;
0047 }
0048
0049 G4RDFluoData::~G4RDFluoData()
0050 {
0051  std::map<G4int,G4DataVector*,std::less<G4int> >::iterator pos;
0052
0053   for (pos = idMap.begin(); pos != idMap.end(); ++pos)
0054     {
0055       G4DataVector* dataSet = (*pos).second;
0056       delete dataSet;
0057     }
0058   for (pos = energyMap.begin(); pos != energyMap.end(); ++pos)
0059     {
0060       G4DataVector* dataSet = (*pos).second;
0061       delete dataSet;
0062     }
0063  for (pos = probabilityMap.begin(); pos != probabilityMap.end(); ++pos)
0064     {
0065       G4DataVector* dataSet = (*pos).second;
0066       delete dataSet;
0067     }
0068 }
0069
0070 size_t G4RDFluoData::NumberOfVacancies() const
0071 {
0072   return numberOfVacancies;
0073 }
0074
0075 G4int G4RDFluoData::VacancyId(G4int vacancyIndex) const
0076 {
0077   G4int n = -1;
0078   if (vacancyIndex<0 || vacancyIndex>=numberOfVacancies)
0079     {G4Exception("G4RDFluoData::VacancyId()", "OutOfRange",
0080                  FatalException, "VacancyIndex outside boundaries!");}
0081   else
0082     {
0083       std::map<G4int,G4DataVector*,std::less<G4int> >::const_iterator pos;
0084       pos = idMap.find(vacancyIndex);
0085       if (pos!= idMap.end())
0086     { G4DataVector dataSet = (*(*pos).second);
0087     n = (G4int) dataSet[0];
0088
0089     }
0090     }
0091   return n;
0092 }
0093
0094 size_t G4RDFluoData::NumberOfTransitions(G4int vacancyIndex) const
0095 {
0096   G4int n = 0;
0097   if (vacancyIndex<0 || vacancyIndex>=numberOfVacancies)
0098     {G4Exception("G4RDFluoData::NumberOfTransitions()", "OutOfRange",
0099                  FatalException, "VacancyIndex outside boundaries!");}
0100   else
0101     {
0102       n = nInitShells[vacancyIndex]-1;
0103       //-1 is necessary because the elements of the vector nInitShells
0104       //include also the vacancy shell:
0105       // -1 subtracts this last one
0106   }
0107  return n;
0108 }
0109 G4int G4RDFluoData::StartShellId(G4int initIndex, G4int vacancyIndex) const
0110 {
0111  G4int n = -1;
0112
0113  if (vacancyIndex<0 || vacancyIndex>=numberOfVacancies)
0114     {G4Exception("G4RDFluoData::StartShellId()", "OutOfRange",
0115                  FatalException, "VacancyIndex outside boundaries!");}
0116  else
0117    {
0118      std::map<G4int,G4DataVector*,std::less<G4int> >::const_iterator pos;
0119
0120      pos = idMap.find(vacancyIndex);
0121
0122      G4DataVector dataSet = *((*pos).second);
0123
0124      G4int nData = dataSet.size();
0125      //The first Element of idMap's dataSets is the original shell of the vacancy,
0126      //so we must start from the first element of dataSet
0127  if (initIndex >= 0 && initIndex < nData)
0128         {
0129           n =  (G4int) dataSet[initIndex+1];
0130
0131         }
0132    }
0133  return n;
0134 }
0135
0136 G4double G4RDFluoData::StartShellEnergy(G4int initIndex, G4int vacancyIndex) const
0137 {
0138   G4double n = -1;
0139
0140   if (vacancyIndex<0 || vacancyIndex>=numberOfVacancies)
0141     {G4Exception("G4RDFluoData::StartShellEnergy()", "OutOfRange",
0142                  FatalException, "VacancyIndex outside boundaries!");}
0143  else
0144    {
0145      std::map<G4int,G4DataVector*,std::less<G4int> >::const_iterator pos;
0146
0147      pos = energyMap.find(vacancyIndex);
0148
0149      G4DataVector dataSet = *((*pos).second);
0150
0151      G4int nData = dataSet.size();
0152      if (initIndex >= 0 && initIndex < nData)
0153        {
0154      n =  dataSet[initIndex];
0155
0156        }
0157    }
0158   return n;
0159 }
0160
0161 G4double G4RDFluoData::StartShellProb(G4int initIndex, G4int vacancyIndex) const
0162 {
0163   G4double n = -1;
0164
0165   if (vacancyIndex<0 || vacancyIndex>=numberOfVacancies)
0166     {G4Exception("G4RDFluoData::StartShellProb()", "OutOfRange",
0167                  FatalException, "VacancyIndex outside boundaries!");}
0168   else
0169     {
0170      std::map<G4int,G4DataVector*,std::less<G4int> >::const_iterator pos;
0171
0172      pos = probabilityMap.find(vacancyIndex);
0173
0174      G4DataVector dataSet = *((*pos).second);
0175
0176      G4int nData = dataSet.size();
0177      if (initIndex >= 0 && initIndex < nData)
0178        {
0179      n =  dataSet[initIndex];
0180
0181        }
0182     }
0183   return n;
0184 }
0185
0186 void G4RDFluoData::LoadData(G4int Z)
0187 {
0188   // Build the complete string identifying the file with the data set
0189
0190   std::ostringstream ost;
0191   if(Z != 0){
0192     ost << "fl-tr-pr-"<< Z << ".dat";
0193   }
0194   else{
0195     ost << "fl-tr-pr-"<<".dat";
0196   }
0197   G4String name(ost.str());
0198
0199   const char* path = G4FindDataDir("G4LEDATA");
0200   if (!path)
0201     {
0202       G4String excep("G4LEDATA environment variable not set!");
0203       G4Exception("G4RDEMDataSet::LoadData()", "InvalidSetup",
0204                   FatalException, excep);
0205     }
0206
0207   G4String pathString(path);
0208   G4String fluor("/fluor/");
0209   G4String dirFile = pathString + fluor + name;
0210   std::ifstream file(dirFile);
0211   std::filebuf* lsdp = file.rdbuf();
0212
0213   if (! (lsdp->is_open()) )
0214     {
0215       G4String excep = "Data file: " + dirFile + " not found";
0216       G4Exception("G4RDEMDataSet::LoadData()", "DataNotFound",
0217                   FatalException, excep);
0218     }
0219
0220   G4double a = 0;
0221   G4int k = 1;
0222   G4int s = 0;
0223
0224   G4int vacIndex = 0;
0225   G4DataVector* initIds = new G4DataVector;
0226   G4DataVector* transEnergies = new G4DataVector;
0227   G4DataVector* transProbabilities = new G4DataVector;
0228
0229   do {
0230     file >> a;
0231     G4int nColumns = 3;
0232     if (a == -1)
0233       {
0234     if (s == 0)
0235       {
0236         // End of a shell data set
0237         idMap[vacIndex] = initIds;
0238             energyMap[vacIndex] = transEnergies;
0239         probabilityMap[vacIndex] = transProbabilities;
0240         //      G4double size=transProbabilities->size();
0241             G4int n = initIds->size();
0242
0243         nInitShells.push_back(n);
0244         numberOfVacancies++;
0245         // Start of new shell data set
0246         initIds = new G4DataVector;
0247             transEnergies = new G4DataVector;
0248         transProbabilities = new G4DataVector;
0249             vacIndex++;
0250       }
0251     s++;
0252     if (s == nColumns)
0253       {
0254         s = 0;
0255       }
0256       }
0257     else if (a == -2)
0258       {
0259     // End of file; delete the empty vectors created
0260     //when encountering the last -1 -1 row
0261     delete initIds;
0262     delete transEnergies;
0263     delete transProbabilities;
0264       }
0265     else
0266       {
0267
0268     if(k%nColumns == 2)
0269       {
0270         // 2nd column is transition  probabilities
0271
0272        if (a != -1) transProbabilities->push_back(a);
0273
0274         k++;
0275       }
0276     else if (k%nColumns == 1)
0277       {
0278         // 1st column is shell id
0279         // if this is the first data of the shell, all the colums are equal
0280         // to the shell Id; so we skip the next colums ang go to the next row
0281         if(initIds->size() == 0) {
0282           if (a != -1) initIds->push_back((G4int)a);
0283           file >> a;
0284           file >> a;
0285           k=k+2;
0286         }
0287         else{
0288           if (a != -1) initIds->push_back(a);
0289         }
0290         k++;
0291       }
0292     else if (k%nColumns == 0)
0293
0294       {//third column is transition energies
0295
0296         if (a != -1)
0297           {G4double e = a * MeV;
0298           transEnergies->push_back(e);}
0299
0300         k=1;
0301       }
0302       }
0303   }
0304   while (a != -2); // end of file
0305   file.close();
0306 }
0307
0308 void G4RDFluoData::PrintData()
0309 {
0310
0311   for (G4int i = 0; i <numberOfVacancies; i++)
0312     {
0313       G4cout << "---- TransitionData for the vacancy nb "
0314          <<i
0315          <<" ----- "
0316          <<G4endl;
0317
0318       for (size_t k = 0; k<NumberOfTransitions(i); k++)
0319     {
0320       G4int id = StartShellId(k,i);
0321     // let's start from 1 because the first (index = 0) element of the vector
0322     // is the id of the intial vacancy
0323       G4double e = StartShellEnergy(k,i) /MeV;
0324       G4double p = StartShellProb(k,i);
0325       G4cout << k <<") Shell id: " << id <<G4endl;
0326       G4cout << " - Transition energy = " << e << " MeV "<<G4endl;
0327       G4cout   << " - Transition probability = " << p <<G4endl;
0328
0329     }
0330       G4cout << "-------------------------------------------------"
0331          << G4endl;
0332     }
0333 }