CLASS
1.1
|
00001 #include "EvolutiveProduct.hxx" 00002 #include "IsotopicVector.hxx" 00003 #include "ZAI.hxx" 00004 #include "LogFile.hxx" 00005 #include "Defines.hxx" 00006 #include "StringLine.hxx" 00007 00008 #include <TGraphErrors.h> 00009 #include <TString.h> 00010 #include <string> 00011 #include <iostream> 00012 #include <fstream> 00013 #include <algorithm> 00014 #include <map> 00015 #include <TGraphErrors.h> 00016 00017 using namespace std; 00018 //________________________________________________________________________ 00019 // 00020 // EvolutiveProduct 00021 // 00022 // 00023 // 00024 // 00025 //________________________________________________________________________ 00026 00027 EvolutiveProduct::EvolutiveProduct(LogFile* Log, int Z, int A, int I, string DBindexfile) 00028 { 00029 DBGL; 00030 fLog = Log; 00031 ifstream DB_index(DBindexfile.c_str()); 00032 if( !DB_index) 00033 { 00034 cout << "!!!EVOLUTIVE DB!!!! Can't open \"" << DBindexfile << "\"" << endl; 00035 fLog->fLog << "!!!EVOLUTIVE DB!!!! Can't open \"" << DBindexfile << "\"" << endl; 00036 exit (1); 00037 } 00038 bool zaifind = false; 00039 string tmp; 00040 getline(DB_index,tmp); // Read first line 00041 while (!DB_index.eof()) 00042 { 00043 string line; 00044 int start=0; 00045 getline(DB_index,line); // Read first line 00046 string first=StringLine::NextWord(line,start); // Read first word 00047 00048 if(first.size()==0) break; // If First word is null.... quit 00049 00050 int rZ=atoi(first.c_str()); // Get Z 00051 int rA=atoi(StringLine::NextWord(line,start).c_str()); // Get A 00052 int rI=atoi(StringLine::NextWord(line,start).c_str()); // Get Isomeric State 00053 00054 if(rZ == Z && rA == A && rI == I) 00055 { 00056 string file_name = StringLine::NextWord(line,start); 00057 ReadDB(file_name); // Read Decay produc DB file name 00058 zaifind = true; 00059 } 00060 } 00061 if(zaifind == false) 00062 { 00063 #pragma omp critical(LOGupdate) 00064 { 00065 fLog->fLog << "!!Warning!! !!!EVOLUTIVE DB!!! Oups... Can't Find the ZAI : " 00066 << Z << " " << A << " " << I << "!!! It will be considered as stable !!" << endl; 00067 AddAsStable(Z, A, I); 00068 } 00069 } 00070 DBGL; 00071 00072 } 00073 00074 //________________________________________________________________________ 00075 00076 EvolutiveProduct::EvolutiveProduct(string DB_reactor_file) 00077 { 00078 DBGL; 00079 string file_name = DB_reactor_file; 00080 ReadDB(file_name); // Read Evolution Produc DB file name 00081 DBGL; 00082 00083 } 00084 00085 //________________________________________________________________________ 00086 EvolutiveProduct::~EvolutiveProduct() 00087 { 00088 DBGL; 00089 DBGL; 00090 } 00091 00092 //________________________________________________________________________ 00093 //________________________________________________________________________ 00094 00095 void EvolutiveProduct::AddAsStable(int Z,int A,int I) 00096 { 00097 DBGL; 00098 double time[2] = {0, (int)(500*365.4*3600*24)}; 00099 double Err[2] = {0, 0}; 00100 double quantity[2] = {1., 1.}; 00101 ZAI zaitmp(Z, A, I); 00102 00103 fEvolutiveProduct.insert(pair<ZAI ,TGraphErrors* >(zaitmp, new TGraphErrors(2, time, quantity, Err, Err) ) ); 00104 DBGL; 00105 } 00106 //________________________________________________________________________ 00107 00108 //________________________________________________________________________ 00109 void EvolutiveProduct::ReadDB(string DBfile) 00110 { 00111 DBGL; 00112 ifstream DecayDB(DBfile.c_str()); // Open the File 00113 if(!DecayDB) 00114 { 00115 cout << "!!Warning!! !!!EvolutiveProduct!!! \n Can't open \"" << DBfile << "\"\n" << endl; 00116 fLog->fLog << "!!Warning!! !!!EvolutiveProduct!!! \n Can't open \"" << DBfile << "\"\n" << endl; 00117 } 00118 vector<double> vTime; 00119 vector<double> vTimeErr; 00120 00121 string line; 00122 int start = 0; 00123 00124 getline(DecayDB, line); // Nuclei is given with "A Z" 00125 if( StringLine::NextWord(line, start, ' ') != "time") 00126 { 00127 cout << "!!Bad Trouble!! !!!EvolutiveProduct!!! Bad Database file : " << DBfile << endl; 00128 fLog->fLog << "!!Bad Trouble!! !!!EvolutiveProduct!!! Bad Database file : " << DBfile << endl; 00129 exit (1); 00130 } 00131 00132 while(start < (int)line.size()) 00133 { 00134 vTime.push_back(atof(StringLine::NextWord(line, start, ' ').c_str())); 00135 vTimeErr.push_back(0); // Set to Zero because not yet error in DB 00136 00137 } 00138 00139 double Time[vTime.size()]; 00140 double TimeErr[vTime.size()]; 00141 for(int i=0; i < (int)vTime.size();i++) 00142 {Time[i] = vTime[i]; TimeErr[i] = vTimeErr[i];} 00143 00144 00145 while (!DecayDB.eof()) 00146 { 00147 double DPQuantity[vTime.size()]; 00148 double DPQuantityErr[vTime.size()]; 00149 00150 getline(DecayDB, line); // Nuclei is given with "A Z" 00151 start = 0; 00152 int Z = atoi(StringLine::NextWord(line, start, ' ').c_str()); 00153 int A = atoi(StringLine::NextWord(line, start, ' ').c_str()); 00154 int I = atoi(StringLine::NextWord(line, start, ' ').c_str()); 00155 00156 if(A!=0 && Z!=0) 00157 { 00158 00159 ZAI zaitmp(Z, A, I); 00160 int i=0; 00161 while(start < (int)line.size()) 00162 { 00163 long double DPQuantityTmp = atof(StringLine::NextWord(line, start, ' ').c_str()); 00164 DPQuantity[i] = (double)DPQuantityTmp; 00165 DPQuantityErr[i]=0; // Set to Zero because not yet error in DB 00166 i++; 00167 00168 } 00169 fEvolutiveProduct.insert(pair<ZAI ,TGraphErrors* >(zaitmp, new TGraphErrors(vTime.size(), Time, DPQuantity, TimeErr, DPQuantityErr) ) ); 00170 } 00171 00172 } 00173 DBGL; 00174 } 00175 00176 //________________________________________________________________________ 00177 Double_t EvolutiveProduct::Interpolate(double t, TGraphErrors& EvolutionGraph) 00178 { 00179 DBGL; 00180 TString fOption; 00181 return (double)EvolutionGraph.Eval(t,0x0,fOption); 00182 DBGL; 00183 } 00184 00185 //________________________________________________________________________ 00186 TGraphErrors* EvolutiveProduct::GetEvolutionTGraphErrors(const ZAI& zai) 00187 { 00188 DBGL; 00189 map<ZAI ,TGraphErrors *>::iterator it = GetEvolutiveProduct().find(zai) ; 00190 00191 if ( it != GetEvolutiveProduct().end() ) 00192 return it->second; 00193 else 00194 return new TGraphErrors(); 00195 00196 DBGL; 00197 } 00198 00199 //________________________________________________________________________ 00200 IsotopicVector EvolutiveProduct::GetIsotopicVectorAt(double t) 00201 { 00202 DBGL; 00203 IsotopicVector IsotopicVectorTmp; 00204 map<ZAI ,TGraphErrors* >::iterator it; 00205 for( it = fEvolutiveProduct.begin(); it != fEvolutiveProduct.end(); it++ ) 00206 { 00207 IsotopicVectorTmp.Add( (*it).first, Interpolate(t, *((*it).second)) ); 00208 } 00209 00210 DBGL; 00211 return IsotopicVectorTmp; 00212 } 00213 00214 //________________________________________________________________________ 00215 00216 00217 00218 00219 00220