Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/abla/src/G4AblaDataFile.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/hadronic/models/abla/src/G4AblaDataFile.cc (Version 11.3.0) and /processes/hadronic/models/abla/src/G4AblaDataFile.cc (Version 11.0.p4)


  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 // ABLAXX statistical de-excitation model          26 // ABLAXX statistical de-excitation model
 27 // Jose Luis Rodriguez, UDC (translation from  <<  27 // Jose Luis Rodriguez, GSI (translation from ABLA07 and contact person)
 28 // Pekka Kaitaniemi, HIP (initial translation      28 // Pekka Kaitaniemi, HIP (initial translation of ablav3p)
 29 // Aleksandra Kelic, GSI (ABLA07 code)             29 // Aleksandra Kelic, GSI (ABLA07 code)
 30 // Davide Mancusi, CEA (contact person INCL)       30 // Davide Mancusi, CEA (contact person INCL)
 31 // Aatos Heikkinen, HIP (project coordination)     31 // Aatos Heikkinen, HIP (project coordination)
 32 //                                                 32 //
 33                                                    33 
                                                   >>  34 #define ABLAXX_IN_GEANT4_MODE 1
                                                   >>  35 
                                                   >>  36 #include "globals.hh"
                                                   >>  37 
 34 #include "G4AblaDataFile.hh"                       38 #include "G4AblaDataFile.hh"
 35 #include "G4AblaDataDefs.hh"                   <<  39 
                                                   >>  40 #ifdef ABLAXX_IN_GEANT4_MODE
 36 #include "globals.hh"                              41 #include "globals.hh"
                                                   >>  42 #else
                                                   >>  43 #include "G4INCLGeant4Compat.hh"
                                                   >>  44 #endif
 37                                                    45 
 38 #include <cmath>                               << 
 39 #include <cstdlib>                             << 
 40 #include <fstream>                                 46 #include <fstream>
                                                   >>  47 #include <cmath>
 41 #include <iostream>                                48 #include <iostream>
                                                   >>  49 #include <cstdlib>
                                                   >>  50 
                                                   >>  51 #ifdef ABLAXX_IN_GEANT4_MODE
                                                   >>  52 G4AblaDataFile::G4AblaDataFile() {
                                                   >>  53 #else
                                                   >>  54 G4AblaDataFile::G4AblaDataFile(G4INCL::Config *config)
                                                   >>  55   : G4AblaVirtualData(config) {
                                                   >>  56   theConfig = config;
                                                   >>  57 #endif
                                                   >>  58   verboseLevel = 0;
                                                   >>  59 }
 42                                                    60 
 43 G4AblaDataFile::G4AblaDataFile() { verboseLeve <<  61 G4AblaDataFile::~G4AblaDataFile()
                                                   >>  62 {
                                                   >>  63 }
 44                                                    64 
 45 /**                                                65 /**
 46  * Read all data from files.                       66  * Read all data from files.
 47  */                                                67  */
 48 G4bool G4AblaDataFile::readData()              <<  68 bool G4AblaDataFile::readData()
 49 {                                                  69 {
 50     if (!G4FindDataDir("G4ABLADATA"))          <<  70 #ifdef ABLAXX_IN_GEANT4_MODE
 51     {                                          <<  71   if(!std::getenv("G4ABLADATA")) {
 52         G4ExceptionDescription ed;             <<  72     //    throw G4HadronicException(__FILE__, __LINE__, "ERROR: Data
 53         ed << " Data missing: set environment  <<  73     //    missing. Set environment variable G4ABLA3.0 to point to the
 54            << " to point to the directory cont <<  74     //    directory containing data files needed by INCL and ABLA
 55            << " by the ABLA model" << G4endl;  <<  75     //    models.");
 56         G4Exception("G4AblaDataFile::readData( <<  76     // G4String errorMessage1 = "ERROR: Data missing. Set environment variable G4ABLADATA\n";
 57     }                                          <<  77     // G4String errorMessage2 = "\t to point to the directory containing data files needed\n";
 58     G4String dataPath(G4FindDataDir("G4ABLADAT <<  78     // G4String errorMessage3 = "\t by INCL and ABLA models.\n";
 59                                                <<  79     // G4String errorMessage = errorMessage1 + errorMessage2 + errorMessage3;
 60     G4String flAlphaFile(dataPath + "/flalpha. <<  80     // G4Exception(errorMessage);
 61     G4String frldmFile(dataPath + "/frldm.dat" <<  81     G4ExceptionDescription ed;
 62     G4String vgsldFile(dataPath + "/vgsld.dat" <<  82     ed << " Data missing: set environment variable G4ABLADATA\n"
 63     G4String rmsFile(dataPath + "/rms.dat");   <<  83        << " to point to the directory containing data files needed\n"
 64     G4String defoFile(dataPath + "/defo.dat"); <<  84        << " by the ABLA model" << G4endl;
 65     G4String massFile(dataPath + "/mass2020.da <<  85     G4Exception("G4AblaDataFile::readData()","ABLA_001",
 66                                                <<  86                 FatalException, ed);
 67     if (verboseLevel > 1)                      <<  87   }
 68     {                                          <<  88   G4String dataPath(std::getenv("G4ABLADATA"));
 69         // G4cout <<"Data path   = " << dataPa <<  89 #else
 70         // G4cout <<"FlAlphaFile = " << flAlph <<  90   G4String dataPath(theConfig->getABLAXXDataFilePath().c_str());
 71         // G4cout <<"FrldmFile   = " << frldmF <<  91 #endif
 72         // G4cout <<"VgsldFile   = " << vgsldF <<  92   G4String flAlphaFile(dataPath + "/flalpha.dat");
 73     }                                          <<  93   G4String frldmFile(  dataPath + "/frldm.dat");
 74                                                <<  94   G4String vgsldFile(  dataPath + "/vgsld.dat");
 75     std::ifstream flalphain(flAlphaFile.c_str( <<  95   G4String pace2File(  dataPath + "/pace2.dat");
 76     std::ifstream frldmin(frldmFile.c_str());  <<  96   G4String rmsFile(  dataPath + "/rms.dat");
 77     std::ifstream vgsldin(vgsldFile.c_str());  <<  97   G4String defoFile(  dataPath + "/defo.dat");
 78     std::ifstream rmsin(rmsFile.c_str());      <<  98   G4String massFile(  dataPath + "/mass2003.dat");
 79     std::ifstream defoin(defoFile.c_str());    <<  99 
 80     std::ifstream massin(massFile.c_str());    << 100   if(verboseLevel > 1) {
 81                                                << 101     // G4cout <<"Data path   = " << dataPath    << G4endl;
 82     if (!massin.is_open())                     << 102     // G4cout <<"FlAlphaFile = " << flAlphaFile << G4endl;
 83     {                                          << 103     // G4cout <<"FrldmFile   = " << frldmFile   << G4endl;
 84         massFile = dataPath + "/mass2016.dat"; << 104     // G4cout <<"VgsldFile   = " << vgsldFile   << G4endl;
 85         massin.close();                        << 105     // G4cout <<"Pace2File   = " << pace2File   << G4endl;
 86         massin.open(massFile.c_str());         << 106   }
 87         std::cout << "Mass evaluation file mas << 107   
 88                                                << 108   std::ifstream flalphain(flAlphaFile.c_str());
 89         if (!massin.is_open())                 << 109   std::ifstream frldmin(frldmFile.c_str());  
 90         {                                      << 110   std::ifstream vgsldin(vgsldFile.c_str());  
 91             massFile = dataPath + "/mass2003.d << 111   std::ifstream pace2in(pace2File.c_str());
 92             massin.close();                    << 112   std::ifstream rmsin(rmsFile.c_str());
 93             massin.open(massFile.c_str());     << 113   std::ifstream defoin(defoFile.c_str());
 94             std::cout << "Mass evaluation file << 114   std::ifstream massin(massFile.c_str());
 95         }                                      << 115 
 96     }                                          << 116   std::filebuf *buf1 = flalphain.rdbuf();
 97                                                << 117   std::filebuf *buf2 = frldmin.rdbuf();
 98     std::filebuf* buf1 = flalphain.rdbuf();    << 118   std::filebuf *buf3 = vgsldin.rdbuf();
 99     std::filebuf* buf2 = frldmin.rdbuf();      << 119   std::filebuf *buf4 = pace2in.rdbuf();  
100     std::filebuf* buf3 = vgsldin.rdbuf();      << 120   std::filebuf *buf5 = rmsin.rdbuf(); 
101     std::filebuf* buf4 = rmsin.rdbuf();        << 121   std::filebuf *buf6 = defoin.rdbuf();
102     std::filebuf* buf5 = defoin.rdbuf();       << 122   std::filebuf *buf7 = massin.rdbuf();
103     std::filebuf* buf6 = massin.rdbuf();       << 123   if (!((buf1->is_open()) && (buf2->is_open()) && (buf3->is_open()) && (buf4->is_open()) && (buf5->is_open()) && (buf6->is_open()) && (buf7->is_open()))) {
104     if (!((buf1->is_open()) && (buf2->is_open( << 124 #ifdef ABLAXX_IN_GEANT4_MODE
105           (buf6->is_open())))                  << 125     G4ExceptionDescription ed;
106     {                                          << 126     ed << "Data missing: could not find ABLA data file in " << dataPath
107         G4ExceptionDescription ed;             << 127        << "defined by environment variable G4ABLADATA" << G4endl;
108         ed << "Data missing: could not find AB << 128     G4Exception("G4AblaDataFile::readData()", "ABLA", FatalException, ed);
109            << "defined by environment variable << 129 #else
110         G4Exception("G4AblaDataFile::readData( << 130     std::cerr << "Error opening file." << std::endl;
111     }                                          << 131 #endif
112                                                << 132   }
113     G4double fflalpha, ffrldm, fvgsld, frms;   << 133   
114     G4int fj = 0, fk = 0, a2, a3, a4;          << 134   G4double fflalpha, ffrldm, fvgsld, fpace2, frms;
115     G4double fbeta2, fbeta4;                   << 135   int fj,fk,a2,a3,a4;
116     G4double a7;                               << 136   G4double fbeta2,fbeta4;
117     const G4int rows = 99;                     << 137   G4double a7;
118     const G4int cols = 154;                    << 138   const G4int rows = 99;
119     const G4int rowsbeta = 137;                << 139   const G4int cols = 154;
120     const G4int colsbeta = 251;                << 140   const G4int rowsbeta = 137;
121                                                << 141   const G4int colsbeta = 251;
122     for (G4int i = 0; i < zcols; i++)          << 142   const G4int rowsmass = 13;
123     {                                          << 143   const G4int colsmass = 154;
124         for (G4int j = 0; j < nrows; j++)      << 144   const G4int massnumbers = 263;
125         {                                      << 145   for(int i = 0; i < rows; i++) {
126             setAlpha(j, i, 0.0);               << 146     for(int j = 0; j < cols; j++) {
127             setEcnz(j, i, 0.0);                << 147       setAlpha(j, i, 0.0);
128             setVgsld(j, i, 0.0);               << 148       setEcnz( j, i, 0.0);
129             setRms(j, i, 0.0);                 << 149       setVgsld(j, i, 0.0);
130         }                                      << 150       setRms(j, i, 0.0);
131     }                                          << 151     }
132                                                << 152   }
133     for (G4int i = 0; i < rows; i++)           << 153   
134     {                                          << 154   for(int i = 0; i < rows; i++) {
135         for (G4int j = 0; j < cols; j++)       << 155     for(int j = 0; j < cols; j++) {
136         {                                      << 156       flalphain >> fflalpha;
137             flalphain >> fflalpha;             << 157       frldmin >> ffrldm;
138             frldmin >> ffrldm;                 << 158       vgsldin >> fvgsld;  
139             vgsldin >> fvgsld;                 << 159       rmsin >> frms;    
140             rmsin >> frms;                     << 160       setAlpha(j, i, fflalpha);
141             setAlpha(j, i, fflalpha);          << 161       setEcnz( j, i, ffrldm);
142             setEcnz(j, i, ffrldm);             << 162       setVgsld(j, i, fvgsld);
143             setVgsld(j, i, fvgsld);            << 163       setRms(j, i, frms);
144             setRms(j, i, frms);                << 164     }
145         }                                      << 165   }
146     }                                          << 166 
147                                                << 167 
148     for (G4int i = 0; i < rowsbeta; i++)       << 168   for(int i = 0; i < rowsbeta; i++) {
149     {                                          << 169     for(int j = 0; j < colsbeta; j++) {
150         for (G4int j = 0; j < colsbeta; j++)   << 170       setBeta2(j, i, 0.0);
151         {                                      << 171       setBeta4(j, i, 0.0);
152             setBeta2(j, i, 0.0);               << 172     }
153             setBeta4(j, i, 0.0);               << 173   }
154         }                                      << 174 
155     }                                          << 175   for(int i = 0; i < 8983; i++) {
156                                                << 176   defoin >> fj >> fk >> fbeta2 >> fbeta4;
157     defoin >> fj >> fk >> fbeta2 >> fbeta4;    << 177   setBeta2(fk, fj, fbeta2);
158     while (!defoin.eof())                      << 178   setBeta4(fk, fj, fbeta4);
159     {                                          << 179   }
160         setBeta2(fk, fj, fbeta2);              << 180 
161         setBeta4(fk, fj, fbeta4);              << 181   for(int i = 0; i < rowsmass; i++) {
162         defoin >> fj >> fk >> fbeta2 >> fbeta4 << 182     for(int j = 0; j < colsmass; j++) {
163     }                                          << 183       setMexp(j, i, 0.0);
164                                                << 184       setMexpID(j,i,0);
165     for (G4int i = 0; i < zcols; i++)          << 185     }
166     {                                          << 186   }
167         for (G4int j = 0; j < nrows; j++)      << 187    massin >> a2 >> a3 >> a4 >> a7 ;
168         {                                      << 188    while(!massin.eof()){
169             setMexp(j, i, 0.0);                << 189    //
170             setMexpID(j, i, 0);                << 190      if(a3<13.){
171         }                                      << 191      setMexpID(a2,a3,1);
172     }                                          << 192      setMexp(a2,a3,938.7829835*a3+939.5653301*a2-1.*a4*a7/1000.);
173     massin >> a2 >> a3 >> a4 >> a7;            << 193      }
174     while (!massin.eof())                      << 194      massin >> a2 >> a3 >> a4 >> a7 ;
175     {                                          << 195    }
176         //                                     << 196 
177         if (a3 < lpcols)                       << 197   flalphain.close();
178         {                                      << 198   frldmin.close();  
179             setMexpID(a2, a3, 1);              << 199   vgsldin.close();
180             setMexp(a2, a3, 938.7829835 * a3 + << 200   rmsin.close();
181         }                                      << 201   defoin.close();
182         massin >> a2 >> a3 >> a4 >> a7;        << 202   massin.close();
183     }                                          << 203 
184                                                << 204   G4String str1, str2, str3;
185     flalphain.close();                         << 205   for(int i = 0; i < 500; i++) {
186     frldmin.close();                           << 206     for(int j = 0; j < 500; j++) {
187     vgsldin.close();                           << 207       setPace2(i, j, 0.0);
188     rmsin.close();                             << 208     }
189     defoin.close();                            << 209   }
190     massin.close();                            << 210   
                                                   >> 211   int A = 0, Zbegin = 0, Zend = 0;
                                                   >> 212   for(int i = 0; i < massnumbers; i++) {
                                                   >> 213     pace2in >> str1 >> A >> str2 >> Zbegin >> str3 >> Zend;
                                                   >> 214     if(Zbegin >= 0 && Zbegin < getPaceCols() &&
                                                   >> 215        A >= 0 && A < getPaceRows()) {
                                                   >> 216       for(int j = Zbegin; j <= Zend; j++) {
                                                   >> 217   pace2in >> fpace2;
                                                   >> 218   setPace2(A, j, fpace2);
                                                   >> 219       }
                                                   >> 220     }
                                                   >> 221   }
                                                   >> 222   pace2in.close();
                                                   >> 223   if(std::abs(getPace2(A, Zend) - 114516.10) > 1e-6) {
                                                   >> 224     std::cerr << "ERROR: Problem in parsing datafile " + pace2File << std::endl;
                                                   >> 225     return false;
                                                   >> 226   }
191                                                   227 
192     return true;                               << 228   return true;
193 }                                                 229 }
                                                   >> 230 
194                                                   231