diff --git a/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx b/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx
index aae0883a9a5320b8cc3af827a6158ec5fd76db50..e5ec0fd535b0f0443ed48e00b8bc9293416f1a7f 100644
--- a/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx
+++ b/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx
@@ -2,6 +2,7 @@
#include "XSModel.hxx"
#include "XSM_MLP.hxx"
#include "CLASSLogger.hxx"
+#include "CLASSMethod.hxx"
#include "StringLine.hxx"
#include "TMVA/Reader.h"
@@ -34,51 +35,74 @@
//________________________________________________________________________
XSM_MLP::XSM_MLP(string TMVA_Weight_Directory,string InformationFile, bool IsTimeStep):XSModel(new CLASSLogger("XSM_MLP.log"))
{
-
+
fIsStepTime=IsTimeStep;
fTMVAWeightFolder = TMVA_Weight_Directory;
if(InformationFile=="")
- fMLPInformationFile = TMVA_Weight_Directory+"/Data_Base_Info.nfo";
+ fInformationFile = TMVA_Weight_Directory+"/Data_Base_Info.nfo";
else
- fMLPInformationFile=fTMVAWeightFolder+InformationFile;
-
+ fInformationFile = fTMVAWeightFolder+InformationFile;
+
GetMLPWeightFiles();
- GetDataBaseInformation();
-
- INFO<<"__A cross section interpolator using" <<endl;
- INFO<<"Multi Layer Perceptron has been define__"<<endl;
+
+ INFO << "__A cross section interpolator using" << endl;
+ INFO << "Multi Layer Perceptron has been define__" << endl;
INFO << " \t His TMVA folder is : \" " << fTMVAWeightFolder << "\"" << endl;
-
+
+ LoadKeyword();
+ ReadNFO();
+
}
+
//________________________________________________________________________
XSM_MLP::XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string InformationFile, bool IsTimeStep):XSModel(Log)
{
-
- fIsStepTime=IsTimeStep;
+
+ fIsStepTime = IsTimeStep;
fTMVAWeightFolder = TMVA_Weight_Directory;
- if(InformationFile=="")
- fMLPInformationFile = TMVA_Weight_Directory+"/Data_Base_Info.nfo";
+ if( InformationFile == "" )
+ fInformationFile = TMVA_Weight_Directory + "/Data_Base_Info.nfo";
else
- fMLPInformationFile=fTMVAWeightFolder+InformationFile;
-
+ fInformationFile = fTMVAWeightFolder+InformationFile;
+
GetMLPWeightFiles();
- GetDataBaseInformation();
-
- INFO<<"__A cross section interpolator using" <<endl;
- INFO<<"Multi Layer Perceptron has been define__"<<endl;
+
+ INFO << "__A cross section interpolator using" << endl;
+ INFO << "Multi Layer Perceptron has been define__" << endl;
INFO << " \t His TMVA folder is : \" " << fTMVAWeightFolder << "\"" << endl;
-
+
+ LoadKeyword();
+ ReadNFO();
+
}
+
//________________________________________________________________________
XSM_MLP::~XSM_MLP()
{
fMapOfTMVAVariableNames.clear();
}
+
+
+
//________________________________________________________________________
-void XSM_MLP::GetDataBaseInformation()
+void XSM_MLP::ReadLine(string line)
{
- ifstream FILE(fMLPInformationFile.c_str());
+ DBGL
+
+ int start = 0;
+ string keyword = tlc(StringLine::NextWord(line, start, ' '));
+ (this->*fKeyword[ keyword ])(line);
+
+ DBGL
+}
+
+
+//________________________________________________________________________
+void XSM_MLP::GetDataBaseInformation()
+{
+ ifstream FILE(fInformationFile.c_str());
+
if(FILE.good())
{
while(!FILE.eof())
@@ -92,13 +116,13 @@ void XSM_MLP::GetDataBaseInformation()
size_t foundTime = line.find("Time (s) :");
size_t foundZAI = line.find("Z A I Name (input MLP) :");
size_t foundDomain = line.find("Fuel range (Z A I min max) :");
-
+
int pos=0;
if(foundRType != std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDBRType = atof( (StringLine::NextWord(line,pos,':')).c_str() );
}
- pos=0;
+ pos=0;
if(foundFType != std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDBFType = atof( (StringLine::NextWord(line,pos,':')).c_str() );
@@ -113,14 +137,14 @@ void XSM_MLP::GetDataBaseInformation()
{ StringLine::NextWord(line,pos,':');
fDBPower = atof( (StringLine::NextWord(line,pos,':') ).c_str() );
}
- pos=0;
+ pos=0;
if(foundTime!=std::string::npos)
{
StringLine::NextWord(line,pos,':');
while( pos< (int)line.size() )
fMLP_Time.push_back( atof( (StringLine::NextWord(line,pos,' ')).c_str() ));
}
- pos=0;
+ pos=0;
if(foundZAI != std::string::npos)
{ string Z;
string A;
@@ -134,14 +158,14 @@ void XSM_MLP::GetDataBaseInformation()
ssline<<line;
ssline>>Z>>A>>I>>Name;
if(StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) )
- {
+ {
fMapOfTMVAVariableNames.insert( pair<ZAI,string>(ZAI(atoi(Z.c_str()),atoi(A.c_str()),atoi(I.c_str())),Name) );
}
-
+
}while((StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I)) && !FILE.eof());
-
+
FILE.seekg(posoflinebeforbadline); //return one line before
-
+
}
if(foundDomain != std::string::npos)
{ string Z;
@@ -157,24 +181,24 @@ void XSM_MLP::GetDataBaseInformation()
ssline<<line;
ssline>>Z>>A>>I>>min>>max;
if(StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) && StringLine::IsDouble(min) && StringLine::IsDouble(max) )
- {
+ {
fZAILimits.insert( pair<ZAI,pair<double,double> >(ZAI(atoi(Z.c_str()),atoi(A.c_str()),atoi(I.c_str())),make_pair(atof(min.c_str()),atof(max.c_str()))) );
}
-
+
}
while((StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) )&& !FILE.eof());
FILE.seekg(posoflinebeforbadline); //return one line before
-
+
}
-
+
}
}
else
{
- ERROR << "Can't find/open file " << fMLPInformationFile << endl;
+ ERROR << "Can't find/open file " << fInformationFile << endl;
exit(0);
}
-
+
/********DEBUG*************************************/
INFO<<"\tMLP XS Data Base Information : "<<endl;
INFO<<"\t\tHeavy Metal (t) :"<<fDBHMMass<<endl;
@@ -183,21 +207,22 @@ void XSM_MLP::GetDataBaseInformation()
for (int i = 0; i < (int)fMLP_Time.size(); ++i)
INFO<<"\t\t\t"<<fMLP_Time[i]<<endl;
INFO<<"\t\tZ A I Name (input MLP) :"<<endl;
-
+
map<ZAI ,string >::iterator it;
-
+
for (it= fMapOfTMVAVariableNames.begin();it!=fMapOfTMVAVariableNames.end();it++)
INFO<<"\t\t\t"<< it->first.Z()<<" "<<it->first.A()<<" "<<it->second<<endl;
-
+
INFO<<"\t\tFuel range"<<endl;
for (map<ZAI,pair<double,double> >::iterator it_dom = fZAILimits.begin();it_dom!=fZAILimits.end();it_dom++)
INFO<<"\t\t\t"<< it_dom->second.first<<" <= "<<it_dom->first.Z()<<" "<<it_dom->first.A()<<" "<<it_dom->first.I()<<" <= "<<it_dom->second.second<<endl;;
-
-
+
+
}
//________________________________________________________________________
void XSM_MLP::GetMLPWeightFiles()
-{DBGL
+{
+ DBGL
/**********Get All TMVA weight files*******************/
//check if the folder containing weights exists
DIR* rep = NULL;
@@ -205,10 +230,10 @@ void XSM_MLP::GetMLPWeightFiles()
rep = opendir(fTMVAWeightFolder.c_str());
if (rep == NULL)
{
- ERROR<<" Reading error for TMVA weight folder "<<fTMVAWeightFolder.c_str()<<" : "<<strerror(errno)<<endl;
+ ERROR<<" Reading error for TMVA weight folder "<<fTMVAWeightFolder.c_str()<<" : "<<strerror(errno)<<endl;
exit(1);
}
-
+
/**Save file names of TMVA weights*/
fWeightFiles.resize(0);
while ((fichierLu = readdir(rep)) != NULL)
@@ -218,7 +243,7 @@ void XSM_MLP::GetMLPWeightFiles()
{
if(FileName[FileName.size()-3]=='x' && FileName[FileName.size()-2]=='m' && FileName[FileName.size()-1]=='l' && FileName[0]!='.' )
fWeightFiles.push_back(FileName);
-
+
}
}
DBGL
@@ -235,71 +260,71 @@ void XSM_MLP::ReadWeightFile(string Filename, int &Z, int &A, int &I, int &React
A=-1;
I=-1;
Reaction=-1;
-
+
size_t found = Filename.find("XS");
-
+
string NameJOB;
NameJOB=Filename.substr(found);
int pos=0;
-
+
StringLine::NextWord(NameJOB, pos, '_');
-
+
Z = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
-
+
A = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
-
+
I = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
-
- string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
- size_t foundext = sReaction.find(".weights.xml");
- sReaction = sReaction.substr(0,foundext);
-
+
+ string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
+ size_t foundext = sReaction.find(".weights.xml");
+ sReaction = sReaction.substr(0,foundext);
+
if(sReaction=="fis")
Reaction=0;
if(sReaction=="cap")
Reaction=1;
if(sReaction=="n2n")
Reaction=2;
-
+
if(Z<=0 || A<=0 || I<0 || Reaction==-1)
{
ERROR << " wrong TMVA weight format " << endl;
exit(0);
}
-
+
}
//________________________________________________________________________
TTree* XSM_MLP::CreateTMVAInputTree(IsotopicVector isotopicvector,int TimeStep)
{
/******Create Input data tree to be interpreted by TMVA::Reader***/
TTree* InputTree = new TTree("InTMP", "InTMP");
-
+
vector<float> InputTMVA;
for(int i = 0 ; i< (int)fMapOfTMVAVariableNames.size() ; i++)
InputTMVA.push_back(0);
-
+
float Time=0;
-
+
IsotopicVector IVInputTMVA;
map<ZAI ,string >::iterator it;
int j=0;
-
+
for( it = fMapOfTMVAVariableNames.begin() ; it!=fMapOfTMVAVariableNames.end() ; it++)
{
- InputTree->Branch( ((*it).second).c_str() ,&InputTMVA[j], ((*it).second + "/F").c_str());
+ InputTree->Branch( ((*it).second).c_str() ,&InputTMVA[j], ((*it).second + "/F").c_str());
IVInputTMVA+= ((*it).first)*1;
j++;
}
-
+
if( !fIsStepTime)
InputTree->Branch( "Time" ,&Time ,"Time/F" );
-
+
IsotopicVector IVAccordingToUserInfoFile = isotopicvector.GetThisComposition(IVInputTMVA);
-
+
double Ntot = IVAccordingToUserInfoFile.GetSumOfAll();
-
+
IVAccordingToUserInfoFile = IVAccordingToUserInfoFile/Ntot;
-
+
j=0;
map<ZAI ,string >::iterator it2;
DBGV("INPUT TMVA");
@@ -309,11 +334,11 @@ TTree* XSM_MLP::CreateTMVAInputTree(IsotopicVector isotopicvector,int TimeStep)
DBGV((*it2).first.Z()<<" "<<(*it2).first.A()<<" "<<InputTMVA[j]);
j++;
}
-
+
Time=fMLP_Time[TimeStep];
-
+
InputTree->Fill();
-
+
return InputTree;
}
//________________________________________________________________________
@@ -322,14 +347,14 @@ double XSM_MLP::ExecuteTMVA(string WeightFile,TTree* InputTree)
DBGV( "File :" << WeightFile);
// --- Create the Reader object
TMVA::Reader *reader = new TMVA::Reader( "Silent" );
-
+
// Create a set of variables and declare them to the reader
// - the variable names MUST corresponds in name and type to those given in the weight file(s) used
vector<float> InputTMVA;
for(int i = 0 ; i< (int)fMapOfTMVAVariableNames.size() ; i++)
InputTMVA.push_back(0);
Float_t Time;
-
+
map<ZAI ,string >::iterator it;
int j=0;
for( it = fMapOfTMVAVariableNames.begin() ; it!=fMapOfTMVAVariableNames.end() ; it++)
@@ -338,18 +363,18 @@ double XSM_MLP::ExecuteTMVA(string WeightFile,TTree* InputTree)
}
if(!fIsStepTime)
reader->AddVariable( "Time" ,&Time);
-
+
// --- Book the MVA methods
-
+
string dir = fTMVAWeightFolder;
if(dir[dir.size()-1]!='/')
- dir+="/";
-
+ dir+="/";
+
// Book method MLP
TString methodName = "MLP method";
TString weightpath = dir + WeightFile ;
reader->BookMVA( methodName, weightpath );
-
+
map<ZAI ,string >::iterator it2;
j=0;
for( it2 = fMapOfTMVAVariableNames.begin() ; it2!=fMapOfTMVAVariableNames.end() ; it2++)
@@ -357,24 +382,25 @@ double XSM_MLP::ExecuteTMVA(string WeightFile,TTree* InputTree)
InputTree->SetBranchAddress(( (*it2).second ).c_str(),&InputTMVA[j]);
j++;
}
-
+
if(!fIsStepTime)
InputTree->SetBranchAddress( "Time" ,&Time );
-
+
InputTree->GetEntry(0);
Float_t val = (reader->EvaluateRegression( methodName ))[0];
-
+
delete reader;
DBGL
-
+
return (double)val;
}
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSectionsTime(IsotopicVector IV)
-{DBGL
-
+{
+ DBGL
+
EvolutionData EvolutionDataFromMLP = EvolutionData();
-
+
map<ZAI,TGraph*> ExtrapolatedXS[3];
/*************DATA BASE INFO****************/
EvolutionDataFromMLP.SetReactorType(fDBRType);
@@ -390,44 +416,44 @@ EvolutionData XSM_MLP::GetCrossSectionsTime(IsotopicVector IV)
int Reaction=-2;
ReadWeightFile( fWeightFiles[i], Z, A, I, Reaction);
if( Z >= GetZAIThreshold() )
- {
+ {
for(int TimeStep=0;TimeStep<int(fMLP_Time.size());TimeStep++)
{
TTree* InputTree = CreateTMVAInputTree(IV,TimeStep);
-
+
pair< map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = ExtrapolatedXS[Reaction].insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph() ) );
-
+
double XSValue = ExecuteTMVA(fWeightFiles[i],InputTree );
if(IResult.second )
{
(IResult.first)->second->SetPoint(0, (double)fMLP_Time[TimeStep], XSValue );
-
+
}
else
{
(IResult.first)->second->SetPoint( (IResult.first)->second->GetN(), (double)fMLP_Time[TimeStep], XSValue );
}
-
+
delete InputTree;
- }
- }
+ }
+ }
}
-
+
/**********Sorting TGraph*********/
for(int x=0;x<3;x++)
{ map<ZAI,TGraph*>::iterator it;
for(it = ExtrapolatedXS[x].begin(); it != ExtrapolatedXS[x].end(); it++)
it->second->Sort();
-
+
}
/**********Filling Matrices*/
EvolutionDataFromMLP.SetFissionXS(ExtrapolatedXS[0]);
EvolutionDataFromMLP.SetCaptureXS(ExtrapolatedXS[1]);
EvolutionDataFromMLP.Setn2nXS(ExtrapolatedXS[2]);
-
-
+
+
DBGL
return EvolutionDataFromMLP;
}
@@ -443,29 +469,29 @@ void XSM_MLP::ReadWeightFileStep(string Filename, int &Z, int &A, int &I, int &R
A=-1;
I=-1;
Reaction=-1;
-
+
size_t found = Filename.find("XS");
-
+
string NameJOB;
NameJOB=Filename.substr(found);
int pos=0;
-
+
StringLine::NextWord(NameJOB, pos, '_');
-
+
Z = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
A = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
I = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
-
- string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
+
+ string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
if(sReaction=="fis")
Reaction=0;
if(sReaction=="cap")
Reaction=1;
if(sReaction=="n2n")
Reaction=2;
-
+
TimeStep = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
-
+
if(Z==-1 || A==-1 || I==-1 || Reaction==-1 || TimeStep==-1)
{
ERROR << " wrong TMVA weight format " << endl;
@@ -475,20 +501,21 @@ void XSM_MLP::ReadWeightFileStep(string Filename, int &Z, int &A, int &I, int &R
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSectionsStep(IsotopicVector IV)
-{DBGL
+{
+ DBGL
TTree* InputTree=CreateTMVAInputTree(IV);
-
+
EvolutionData EvolutionDataFromMLP = EvolutionData();
-
+
map<ZAI,TGraph*> ExtrapolatedXS[3];
/*************DATA BASE INFO****************/
EvolutionDataFromMLP.SetReactorType("PWR");
EvolutionDataFromMLP.SetFuelType("MOX");
EvolutionDataFromMLP.SetPower(fDBPower);
EvolutionDataFromMLP.SetHeavyMetalMass(fDBHMMass);
-
+
/************* The Cross sections***********/
-
+
for(int i=0;i<int(fWeightFiles.size());i++)
{
int Z=-2;
@@ -498,13 +525,13 @@ EvolutionData XSM_MLP::GetCrossSectionsStep(IsotopicVector IV)
int TimeStep=-2;
ReadWeightFileStep( fWeightFiles[i], Z, A, I, Reaction, TimeStep);
if( Z >= GetZAIThreshold() )
- {
+ {
ZAI zaitmp = ZAI(Z,A,I);
-
+
pair< map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = ExtrapolatedXS[Reaction].insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph() ) );
-
+
if( IResult.second )
{
(IResult.first)->second->SetPoint(0, (double)fMLP_Time[TimeStep], ExecuteTMVA(fWeightFiles[i],InputTree) );
@@ -515,7 +542,7 @@ EvolutionData XSM_MLP::GetCrossSectionsStep(IsotopicVector IV)
}
}
}
-
+
/**********Sorting TGraph*********/
for(int x=0;x<3;x++)
{ map<ZAI,TGraph*>::iterator it;
@@ -526,21 +553,22 @@ EvolutionData XSM_MLP::GetCrossSectionsStep(IsotopicVector IV)
EvolutionDataFromMLP.SetFissionXS(ExtrapolatedXS[0]);
EvolutionDataFromMLP.SetCaptureXS(ExtrapolatedXS[1]);
EvolutionDataFromMLP.Setn2nXS(ExtrapolatedXS[2]);
-
+
delete InputTree;
DBGL
return EvolutionDataFromMLP;
}
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSections(IsotopicVector IV ,double t)
-{DBGL
+{
+ DBGL
if(t!=0)
WARNING << " Argument t has non effect here " << endl;
-
+
EvolutionData EV;
if(fIsStepTime)
EV=GetCrossSectionsStep(IV);
-
+
else
EV=GetCrossSectionsTime(IV);
diff --git a/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx b/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx
index eda1d1faa8a7c35569f88d335b87fb41be37de5b..4142b3cad9ec9bc6dc705d92cd44971c796cc657 100644
--- a/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx
+++ b/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx
@@ -6,8 +6,8 @@
/*!
\file
\brief Header file for XSM_MLP class.
-
-
+
+
@authors BLG
@version 1.0
*/
@@ -21,14 +21,19 @@
typedef long long int cSecond;
using namespace std;
+
+class XSM_MLP;
+#ifndef __CINT__
+typedef void (XSM_MLP::*DMthPtr)( const string & ) ;
+#endif
//-----------------------------------------------------------------------------//
//! Defines a XSModel getting mean cross sections from neural network execution
/*!
Define a XSM_MLP.
- This is the class to predict cross sections with a
+ This is the class to predict cross sections with a
set of Multi Layer Perceptrons (MLP)
-
+
@authors BLG
@version 1.0
*/
@@ -38,19 +43,19 @@ using namespace std;
class XSM_MLP : public XSModel
{
public :
-
+
/*!
\name Constructor/Desctructor
*/
//@{
-
+
//{
/// Normal Constructor
/*!
\param TMVA_Weight_Directory : The directory where all the TMVA weight are located
\param InformationFile : Name of the information file located in TMVA_Weight_Directory (default : Data_Base_Info.nfo)
\param IsTimeStep : if true , one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
-
+
*/
XSM_MLP(string TMVA_Weight_Directory,string InformationFile="",bool IsTimeStep=false);
//}
@@ -69,40 +74,43 @@ class XSM_MLP : public XSModel
~XSM_MLP();
//@}
-
-
- virtual EvolutionData GetCrossSections(IsotopicVector IV,double t=0); //!< Return calculated cross section by the MLP regression
-
-
+
+ void LoadKeyword() {}
+
+ EvolutionData GetCrossSections(IsotopicVector IV,double t=0); //!< Return calculated cross section by the MLP regression
+
+ void ReadLine(string line);
+
private :
void GetDataBaseInformation(); //!< Read information file and fill Reactor Type, Fuel type, HM mass, Power, time vector, and TMVA input variables names
-
- void GetMLPWeightFiles(); //!< Find all .xml file in TMVA_Weight_Directory
+
+ void GetMLPWeightFiles(); //!< Find all .xml file in TMVA_Weight_Directory
EvolutionData GetCrossSectionsStep(IsotopicVector IV); //!< Return calculated cross section by the MLP regression when fIsTimeStep==true
EvolutionData GetCrossSectionsTime(IsotopicVector IV); //!< Return calculated cross section by the MLP regression when fIsTimeStep==false
void ReadWeightFile(string Filename, int &Z, int &A, int &I, int &Reaction) ; //!< Select the reaction according to the weight file name
void ReadWeightFileStep(string Filename, int &Z, int &A, int &I, int &Reaction, int &TimeStep);; //!< Select the reaction according to the weight file name
-
-
-
+
+
+
double ExecuteTMVA(string WeightFile, TTree* InputTree); //!<Execute the MLP according to the input tree created
TTree* CreateTMVAInputTree(IsotopicVector isotopicvector,int TimeStep=0); //!<Create input tmva tree to be read by ExecuteTMVA
-
-
- vector<double> fMLP_Time; //!< Time vector of the data base
- vector<string> fWeightFiles; //!< All the weight file contains in fTMVAWeightFolder
+
+
+ vector<double> fMLP_Time; //!< Time vector of the data base
+ vector<string> fWeightFiles; //!< All the weight file contains in fTMVAWeightFolder
string fTMVAWeightFolder; //!< folder containing all the weight file
- string fMLPInformationFile; //!< file containing Reactor Type, Fuel type, HM mass, Power, time vector, and TMVA input variables names (looks the manual for format details)
- bool fIsStepTime; //!< true if one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
-
- map<ZAI,string> fMapOfTMVAVariableNames;//!< List of TMVA input variable names (read from fMLPInformationFile ) , name depends on the training step
+ bool fIsStepTime; //!< true if one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
+ map<ZAI,string> fMapOfTMVAVariableNames;//!< List of TMVA input variable names (read from fMLPInformationFile ) , name depends on the training step
+#ifndef __CINT__
+ map<string, DMthPtr> fDKeyword;
+#endif
};
#endif
diff --git a/source/branches/BaM_Dev/include/XSModel.hxx b/source/branches/BaM_Dev/include/XSModel.hxx
index bd522cce4a8782c2dde356352ba73983eff9ba13..719fe067239d904c2d628624eb2a6c834c0dea6a 100644
--- a/source/branches/BaM_Dev/include/XSModel.hxx
+++ b/source/branches/BaM_Dev/include/XSModel.hxx
@@ -30,14 +30,14 @@ typedef void (XSModel::*MthPtr)( const string & ) ;
//! Defines a mean cross section predictor
/*!
-This is the mother class for methods related to XS prediction
-
-\warning
+ This is the mother class for methods related to XS prediction
+
+ \warning
Never instantiate XSModel in your CLASS input but it's derivated class
-
+
@see XSM_CLOSEST
@see XSM_MLP
-
+
@author BaM
@author BLG
@version 1.0
@@ -47,19 +47,19 @@ This is the mother class for methods related to XS prediction
class XSModel : public CLASSObject
{
-
- public :
-
+
+ public :
+
/*!
\name Constructor/Desctructor
*/
//@{
-
+
XSModel(); //!<Default constructor
XSModel(CLASSLogger* log); //!<Logger constructor
-
+
//@}
@@ -90,23 +90,29 @@ class XSModel : public CLASSObject
virtual bool isIVInDomain(IsotopicVector IV) ;
//@}
+ void ReadNFO();
+ virtual void ReadLine(string line);
+
+
void ReadZAIlimits(const string &line);
void ReadType(const string &line);
void ReadRParam(const string &line);
- void LoadKeyword();
+ virtual void LoadKeyword();
void SetZAIThreshold(int Z_Threshold){fZAIThreshold = Z_Threshold;}//!< Set the Z threshold : ZAI with Z < fZAIThreshold are not manage by CLASS
int GetZAIThreshold(){return fZAIThreshold;}//!< Get the Z threshold
-
+
protected :
+ bool freaded;
+ string fInformationFile; //!< file containing Reactor Type, Fuel type, HM mass, Power, time vector, and TMVA input variables names (looks the manual for format details)
double fDBPower; //!< Power of the data base (read from fMLPInformationFile )
double fDBHMMass; //!< Heavy metal mass of the data base (read from fMLPInformationFile )
string fDBFType; //!< Fuel Type (e.g MOX, UOX, ThU, ThPu ...)
string fDBRType; //!< Reactor Type (e.g PWR, FBR-Na, ADS..)
-
+
map< ZAI, pair<double,double> > fZAILimits; //!< Fresh fuel range : map<ZAI<min edge ,max edge >>
#ifndef __CINT__
diff --git a/source/branches/BaM_Dev/src/EvolutionData.cxx b/source/branches/BaM_Dev/src/EvolutionData.cxx
index 6be228a73bd75dc7130dd5719c19c368bcb245b0..e06985a6b6775fd1f2de5ad9621f2cf9b8752c70 100755
--- a/source/branches/BaM_Dev/src/EvolutionData.cxx
+++ b/source/branches/BaM_Dev/src/EvolutionData.cxx
@@ -13,64 +13,64 @@
#include <fstream>
#include <algorithm>
- //________________________________________________________________________
- //
- // EvolutionData
- //
- //
- //
- //
- //________________________________________________________________________
+//________________________________________________________________________
+//
+// EvolutionData
+//
+//
+//
+//
+//________________________________________________________________________
double Distance(IsotopicVector IV1, EvolutionData Evd1 )
{
-
+
double d2=0;
IsotopicVector IV2 = Evd1.GetIsotopicVectorAt(0.);
-
+
IsotopicVector IVtmp = IV1;
map<ZAI ,double> IVtmpIsotopicQuantity = IVtmp.GetIsotopicQuantity();
map<ZAI ,double >::iterator it;
-
+
double SumOfXs = 0;
for( it = IVtmpIsotopicQuantity.begin(); it != IVtmpIsotopicQuantity.end(); it++)
{
-
+
SumOfXs += Evd1.GetXSForAt(0., (*it).first, 1);
SumOfXs += Evd1.GetXSForAt(0., (*it).first, 2);
SumOfXs += Evd1.GetXSForAt(0., (*it).first, 3);
-
+
}
-
-
+
+
for( it = IVtmpIsotopicQuantity.begin(); it != IVtmpIsotopicQuantity.end(); it++)
{
double Z1 = 0.0;
double Z2 = 0.0;
double XS = 0.0;
-
+
Z1 = IV1.GetZAIIsotopicQuantity( (*it).first );
Z2 = IV2.GetZAIIsotopicQuantity( (*it).first );
XS = Evd1.GetXSForAt(0., (*it).first, 1)
- + Evd1.GetXSForAt(0., (*it).first, 2)
- + Evd1.GetXSForAt(0., (*it).first, 3);
-
+ + Evd1.GetXSForAt(0., (*it).first, 2)
+ + Evd1.GetXSForAt(0., (*it).first, 3);
+
d2 += pow(Z1-Z2 , 2 ) * pow(XS,2);
-
+
}
return sqrt(d2)/SumOfXs;
-
+
}
double Distance(EvolutionData Evd1, IsotopicVector IV1 )
{
return Distance(IV1,Evd1);
}
- //________________________________________________________________________
- //________________________________________________________________________
- //________________________________________________________________________
+//________________________________________________________________________
+//________________________________________________________________________
+//________________________________________________________________________
EvolutionData operator*(EvolutionData const& evol, double F)
{
@@ -96,29 +96,29 @@ EvolutionData operator*(EvolutionData const& evol, double F)
evoltmp.SetFissionXS(evol.GetFissionXS());
evoltmp.SetCaptureXS(evol.GetCaptureXS());
evoltmp.Setn2nXS(evol.Getn2nXS());
-
-
+
+
return evoltmp;
}
//________________________________________________________________________
EvolutionData operator*(double F, EvolutionData const& evol)
{
-
+
return evol*F;
-
+
}
//________________________________________________________________________
EvolutionData operator/(EvolutionData const& evol, double F)
{
-
+
return evol*(1./F);
-
+
}
EvolutionData Multiply(EvolutionData const& evol, double F)
{
-
+
EvolutionData evoltmp;
map<ZAI ,TGraph* > EvolutionData = evol.GetInventoryEvolution();
map<ZAI ,TGraph* >::iterator it;
@@ -126,7 +126,7 @@ EvolutionData Multiply(EvolutionData const& evol, double F)
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
-
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -134,16 +134,16 @@ EvolutionData Multiply(EvolutionData const& evol, double F)
Y[i] = y*F;
}
evoltmp.NucleiInsert( pair<ZAI, TGraph*> ( (*it).first,new TGraph((*it).second->GetN(), X, Y) ) );
-
+
}
-
-
+
+
EvolutionData = evol.GetFissionXS();
for(it = EvolutionData.begin(); it != EvolutionData.end(); it++)
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
-
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -151,15 +151,15 @@ EvolutionData Multiply(EvolutionData const& evol, double F)
Y[i] = y*F;
}
evoltmp.FissionXSInsert( pair<ZAI, TGraph*> ( (*it).first,new TGraph((*it).second->GetN(), X, Y) ) );
-
+
}
-
+
EvolutionData = evol.GetCaptureXS();
for(it = EvolutionData.begin(); it != EvolutionData.end(); it++)
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
-
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -174,7 +174,7 @@ EvolutionData Multiply(EvolutionData const& evol, double F)
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
-
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -182,23 +182,23 @@ EvolutionData Multiply(EvolutionData const& evol, double F)
Y[i] = y*F;
}
evoltmp.n2nXSInsert( pair<ZAI, TGraph*> ( (*it).first,new TGraph((*it).second->GetN(), X, Y) ) );
-
+
}
-
+
evoltmp.SetPower(evol.GetPower()*F);
-
-
-
+
+
+
return evoltmp;
-
+
}
//________________________________________________________________________
EvolutionData Multiply(double F, EvolutionData const& evol)
{
-
+
return Multiply(evol,F);
-
+
}
EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
@@ -207,19 +207,19 @@ EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
map<ZAI ,TGraph* > EvolutionData1 = EvolSum.GetInventoryEvolution();
map<ZAI ,TGraph* > EvolutionData2 = evol2.GetInventoryEvolution();
map<ZAI ,TGraph* >::iterator it;
-
+
for(it = EvolutionData2.begin(); it != EvolutionData2.end(); it++)
{
pair<map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = EvolutionData1.insert( pair<ZAI, TGraph*> ( *it ) );
if(!(IResult.second) )
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
map<ZAI ,TGraph* >::iterator it2 = EvolutionData1.find( (*it).first );
-
-
+
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -229,27 +229,27 @@ EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
IResult.first->second = new TGraph((*it).second->GetN(), X, Y);
}
-
+
}
EvolSum.SetInventoryEvolution(EvolutionData1);
-
-
+
+
EvolutionData1 = evol1.GetFissionXS();
EvolutionData2 = evol2.GetFissionXS();
-
+
for(it = EvolutionData2.begin(); it != EvolutionData2.end(); it++)
{
pair<map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = EvolutionData1.insert( pair<ZAI, TGraph*> ( *it ) );
-
+
if(!(IResult.second) )
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
map<ZAI ,TGraph* >::iterator it2 = EvolutionData1.find( (*it).first );
-
-
+
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -260,24 +260,24 @@ EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
}
}
EvolSum.SetFissionXS(EvolutionData1);
-
-
+
+
EvolutionData1 = EvolSum.GetCaptureXS();
EvolutionData2 = evol2.GetCaptureXS();
-
+
for(it = EvolutionData2.begin(); it != EvolutionData2.end(); it++)
{
pair<map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = EvolutionData1.insert( pair<ZAI, TGraph*> ( *it ) );
-
+
if(!(IResult.second) )
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
map<ZAI ,TGraph* >::iterator it2 = EvolutionData1.find( (*it).first );
-
-
+
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -288,24 +288,24 @@ EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
}
}
EvolSum.SetCaptureXS(EvolutionData1);
-
-
+
+
EvolutionData1 = EvolSum.Getn2nXS();
EvolutionData2 = evol2.Getn2nXS();
-
+
for(it = EvolutionData2.begin(); it != EvolutionData2.end(); it++)
{
pair<map<ZAI, TGraph*>::iterator, bool> IResult;
-
+
IResult = EvolutionData1.insert( pair<ZAI, TGraph*> ( *it ) );
-
+
if(!(IResult.second) )
{
double X[(*it).second->GetN()];
double Y[(*it).second->GetN()];
map<ZAI ,TGraph* >::iterator it2 = EvolutionData1.find( (*it).first );
-
-
+
+
for(int i = 0; i < (*it).second->GetN(); i++)
{
double y;
@@ -316,7 +316,7 @@ EvolutionData Sum(EvolutionData const& evol1, EvolutionData const& evol2)
}
}
EvolSum.Setn2nXS(EvolutionData1);
-
+
return EvolSum;
}
@@ -352,7 +352,7 @@ EvolutionData::EvolutionData():CLASSObject()
fFlux = 0;
}
- //________________________________________________________________________
+//________________________________________________________________________
EvolutionData::EvolutionData(CLASSLogger* log):CLASSObject(log)
{
fIsCrossSection = false;
@@ -362,7 +362,7 @@ EvolutionData::EvolutionData(CLASSLogger* log):CLASSObject(log)
fFlux = 0;
}
- //________________________________________________________________________
+//________________________________________________________________________
EvolutionData::EvolutionData(CLASSLogger* log, string DB_file, bool oldread, ZAI zai):CLASSObject(log)
{
@@ -383,17 +383,17 @@ EvolutionData::EvolutionData(CLASSLogger* log, string DB_file, bool oldread, ZAI
}
- //________________________________________________________________________
+//________________________________________________________________________
EvolutionData::~EvolutionData()
{
-
+
}
//________________________________________________________________________
void EvolutionData::DeleteEvolutionData()
{
-
+
map<ZAI ,TGraph* >::iterator it_del;
-
+
for( it_del = fInventoryEvolution.begin(); it_del != fInventoryEvolution.end(); it_del++)
{
delete (*it_del).second;
@@ -414,19 +414,19 @@ void EvolutionData::DeleteEvolutionData()
delete (*it_del).second;
(*it_del).second = 0;
}
-
-
+
+
delete fKeff;
delete fFlux;
-
+
fInventoryEvolution.clear();
fFissionXS.clear();
fCaptureXS.clear();
fn2nXS.clear();
-
+
fFlux = 0;
fKeff = 0;
-
+
}
@@ -466,7 +466,7 @@ bool EvolutionData::n2nXSInsert(pair<ZAI, TGraph*> zaitoinsert)
}
- //________________________________________________________________________
+//________________________________________________________________________
void EvolutionData::AddAsStable(ZAI zai)
{
@@ -477,7 +477,7 @@ void EvolutionData::AddAsStable(ZAI zai)
}
- //________________________________________________________________________
+//________________________________________________________________________
Double_t EvolutionData::Interpolate(double t, TGraph& EvolutionGraph)
{
@@ -486,7 +486,7 @@ Double_t EvolutionData::Interpolate(double t, TGraph& EvolutionGraph)
}
- //________________________________________________________________________
+//________________________________________________________________________
TGraph* EvolutionData::GetEvolutionTGraph(const ZAI& zai)
{
@@ -500,7 +500,7 @@ TGraph* EvolutionData::GetEvolutionTGraph(const ZAI& zai)
}
- //________________________________________________________________________
+//________________________________________________________________________
IsotopicVector EvolutionData::GetIsotopicVectorAt(double t)
{
@@ -515,7 +515,7 @@ IsotopicVector EvolutionData::GetIsotopicVectorAt(double t)
return IsotopicVectorTmp;
}
- //________________________________________________________________________
+//________________________________________________________________________
double EvolutionData::GetXSForAt(double t, ZAI zai, int ReactionId)
{
@@ -552,9 +552,9 @@ double EvolutionData::GetXSForAt(double t, ZAI zai, int ReactionId)
- //________________________________________________________________________
+//________________________________________________________________________
- //________________________________________________________________________
+//________________________________________________________________________
@@ -564,18 +564,18 @@ void EvolutionData::ReadDB(string DBfile, bool oldread)
if(oldread)
{
-
+
OldReadDB(DBfile);
return;
}
-
+
ReadInfo(); // Read the .info associeted
-
+
ifstream DecayDB(DBfile.c_str()); // Open the File
if(!DecayDB) //check if file is correctly open
{
ERROR << " \n Can't open \"" << DBfile << "\"\n" << endl;
- ERROR <<"\t-> Hint : If loading .dat files using a .idx file (like for a decay data base)\nmake sure the paths in it are correct";
+ ERROR << "\t-> Hint : If loading .dat files using a .idx file (like for a decay data base)\nmake sure the paths in it are correct" << endl;;
exit(1);
}
vector<double> vTime;
@@ -599,7 +599,7 @@ void EvolutionData::ReadDB(string DBfile, bool oldread)
for(int i=0; i < (int)vTime.size();i++)
Time[i] = vTime[i];
const int NTimeStep = sizeof(Time)/sizeof(double);
-
+
enum { Keff=1, Flux, Inv, XSFis, XSCap, XSn2n };
@@ -610,7 +610,7 @@ void EvolutionData::ReadDB(string DBfile, bool oldread)
keyword_map["xscap"] = XSCap;
keyword_map["xsn2n"] = XSn2n;
keyword_map["inv"] = Inv;
-
+
getline(DecayDB, line);
do
{
@@ -624,7 +624,7 @@ void EvolutionData::ReadDB(string DBfile, bool oldread)
case Flux:
ReadFlux(line, Time, NTimeStep);
break;
-
+
case Inv:
ReadInv(line, Time, NTimeStep);
break;
@@ -647,22 +647,22 @@ void EvolutionData::ReadDB(string DBfile, bool oldread)
}
getline(DecayDB, line);
-
+
}while ( !DecayDB.eof() );
-
+
fHeavyMetalMass = cZAIMass.GetMass( GetIsotopicVectorAt(0.).GetActinidesComposition() );
-
+
DecayDB.close();
-
-
+
+
}
void EvolutionData::ReadKeff(string line, double* time, int NTimeStep)
{
if(fKeff != 0)
delete fKeff;
-
+
int start = 0;
if( tlc(StringLine::NextWord(line, start, ' ')) != "keff" ) // Check the keyword
{
@@ -672,15 +672,15 @@ void EvolutionData::ReadKeff(string line, double* time, int NTimeStep)
double Keff[NTimeStep];
-
+
int i = 0;
while(start < (int)line.size() && i < NTimeStep ) // Read the Data
{
Keff[i] = atof(StringLine::NextWord(line, start, ' ').c_str()) ;
i++;
}
-
-
+
+
fKeff = new TGraph(NTimeStep, time, Keff); // Add the TGraph
@@ -693,7 +693,7 @@ void EvolutionData::ReadFlux(string line, double* time, int NTimeStep)
if(fFlux != 0)
delete fFlux;
-
+
int start = 0;
if( tlc(StringLine::NextWord(line, start, ' ')) != "flux" ) // Check the keyword
@@ -714,7 +714,7 @@ void EvolutionData::ReadFlux(string line, double* time, int NTimeStep)
fFlux = new TGraph(NTimeStep, time, Flux); // Add the TGraph
-
+
}
@@ -729,7 +729,7 @@ void EvolutionData::ReadInv(string line, double* time, int NTimeStep)
ERROR << " Bad keyword : \"inv\" not found !" << endl;
exit(1);
}
- // Read the Z A I
+ // Read the Z A I
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
@@ -744,7 +744,7 @@ void EvolutionData::ReadInv(string line, double* time, int NTimeStep)
Inv[i] = atof(StringLine::NextWord(line, start, ' ').c_str()) ;
i++;
}
- // Add the TGraph
+ // Add the TGraph
fInventoryEvolution.insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph(NTimeStep, time, Inv) ) );
}
@@ -755,14 +755,14 @@ void EvolutionData::ReadInv(string line, double* time, int NTimeStep)
void EvolutionData::ReadXSFis(string line, double* time, int NTimeStep)
{
-
+
int start = 0;
if( tlc(StringLine::NextWord(line, start, ' ')) != "xsfis" ) // Check the keyword
{
ERROR << " Bad keyword : \"xsfis\" not found !" << endl;
exit(1);
}
- // Read the Z A I
+ // Read the Z A I
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
@@ -777,8 +777,8 @@ void EvolutionData::ReadXSFis(string line, double* time, int NTimeStep)
XSFis[i] = atof(StringLine::NextWord(line, start, ' ').c_str()) ;
i++;
}
-
- // Add the TGraph
+
+ // Add the TGraph
fFissionXS.insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph(NTimeStep, time, XSFis) ) );
}
@@ -788,14 +788,14 @@ void EvolutionData::ReadXSFis(string line, double* time, int NTimeStep)
void EvolutionData::ReadXSCap(string line, double* time, int NTimeStep)
{
-
+
int start = 0;
if( tlc(StringLine::NextWord(line, start, ' ')) != "xscap" ) // Check the keyword
{
ERROR << " Bad keyword : \"xscap\" not found !" << endl;
exit(1);
}
- // Read the Z A I
+ // Read the Z A I
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
@@ -810,8 +810,8 @@ void EvolutionData::ReadXSCap(string line, double* time, int NTimeStep)
XSCap[i] = atof(StringLine::NextWord(line, start, ' ').c_str()) ;
i++;
}
-
- // Add the TGraph
+
+ // Add the TGraph
fCaptureXS.insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph(NTimeStep, time, XSCap) ) );
}
@@ -821,14 +821,14 @@ void EvolutionData::ReadXSCap(string line, double* time, int NTimeStep)
void EvolutionData::ReadXSn2n(string line, double* time, int NTimeStep)
{
-
+
int start = 0;
if( tlc(StringLine::NextWord(line, start, ' ')) != "xsn2n" ) // Check the keyword
{
ERROR << " Bad keyword : \"xsn2n\" not found !" << endl;
exit(1);
}
- // Read the Z A I
+ // Read the Z A I
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
@@ -843,7 +843,7 @@ void EvolutionData::ReadXSn2n(string line, double* time, int NTimeStep)
XSn2n[i] = atof(StringLine::NextWord(line, start, ' ').c_str()) ;
i++;
}
- // Add the TGraph
+ // Add the TGraph
fn2nXS.insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph(NTimeStep, time, XSn2n) ) );
}
@@ -857,27 +857,27 @@ void EvolutionData::ReadInfo()
string InfoDBFile = fDB_file.erase(fDB_file.size()-3,fDB_file.size());
InfoDBFile += "Info";
ifstream InfoDB_tmp(InfoDBFile.c_str()); // Open the File
-
+
if(!InfoDB_tmp)
{
InfoDBFile = InfoDBFile.erase(InfoDBFile.size()-4,InfoDBFile.size());
InfoDBFile += "info";
-
+
}
InfoDB_tmp.close();
-
+
ifstream InfoDB(InfoDBFile.c_str()); // Open the File
if(!InfoDB)
{
WARNING << "!!ERROR!! !!!EvolutionData!!! \n Can't open \"" << InfoDBFile << "\"\n" << endl;
}
-
+
int start = 0;
string line;
getline(InfoDB, line);
if ( tlc(StringLine::NextWord(line, start, ' ')) == "reactor")
fReactorType = StringLine::NextWord(line, start, ' ');
-
+
start = 0;
getline(InfoDB, line);
if (tlc(StringLine::NextWord(line, start, ' ')) == "fueltype")
@@ -896,32 +896,32 @@ void EvolutionData::ReadInfo()
void EvolutionData::OldReadDB(string DBfile)
{
-
-
+
+
ifstream DecayDB(DBfile.c_str()); // Open the File
if(!DecayDB)
{
ERROR << " Can't open \"" << DBfile << "\"\n" << endl;
ERROR <<"\t-> Hint : If loading .dat files using a .idx file (like for a decay data base)\nmake sure the paths in it are correct";
-
+
exit(1);
}
vector<double> vTime;
vector<double> vTimeErr;
-
+
string line;
int start = 0;
-
+
getline(DecayDB, line);
if( StringLine::NextWord(line, start, ' ') != "time")
{
ERROR << " Bad Database file : " << DBfile << endl;
exit (1);
}
-
+
while(start < (int)line.size())
vTime.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
-
+
fFinalTime = vTime.back();
double Time[vTime.size()];
for(int i=0; i < (int)vTime.size();i++)
@@ -935,48 +935,48 @@ void EvolutionData::OldReadDB(string DBfile)
vector<double> vKeff;
while(start < (int)line.size())
vKeff.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
-
+
double Keff[vKeff.size()];
for(int i=0; i < (int)vKeff.size();i++)
Keff[i] = vKeff[i];
-
+
fKeff = new TGraph(vTime.size(), Time, Keff);
-
+
start = 0;
getline(DecayDB, line);
if (StringLine::NextWord(line, start, ' ') == "flux")
{
-
-
+
+
while(start < (int)line.size())
vFlux.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
-
+
double Flux[vFlux.size()];
for(int i=0; i < (int)vFlux.size();i++)
Flux[i] = vFlux[i];
-
+
fFlux = new TGraph(vTime.size(), Time, Flux);
-
+
}
}
-
-
+
+
do
{
-
-
+
+
start = 0;
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
-
+
if(A!=0 && Z!=0)
{
double DPQuantity[vTime.size()];
for(int k = 0; k < (int)vTime.size(); k++ )
DPQuantity[k] = 0;
-
-
+
+
ZAI zaitmp(Z, A, I);
int i=0;
while(start < (int)line.size())
@@ -984,39 +984,39 @@ void EvolutionData::OldReadDB(string DBfile)
double DPQuantityTmp = atof(StringLine::NextWord(line, start, ' ').c_str());
DPQuantity[i] = (double)DPQuantityTmp;
i++;
-
+
}
TGraph* tgraphtmp = new TGraph((int)vTime.size()-1, Time, DPQuantity);
fInventoryEvolution.insert(pair<ZAI ,TGraph* >(zaitmp, tgraphtmp) );
}
-
+
getline(DecayDB, line);
if(line == "" || line == "CrossSection" ) break;
}while (!DecayDB.eof() );
-
+
if(line == "CrossSection")
{
fIsCrossSection = true;
getline(DecayDB, line);
-
+
if (line == "Fission")
{
getline(DecayDB, line);
-
+
do
{
double DPQuantity[vTime.size()];
for(int k = 0; k < (int)vTime.size(); k++ )
DPQuantity[k] = 0;
-
+
start = 0;
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
if(A!=0 && Z!=0)
{
-
-
+
+
ZAI zaitmp(Z, A, I);
int i=0;
while(start < (int)line.size())
@@ -1024,7 +1024,7 @@ void EvolutionData::OldReadDB(string DBfile)
long double DPQuantityTmp = atof(StringLine::NextWord(line, start, ' ').c_str());
DPQuantity[i] = (double)DPQuantityTmp;
i++;
-
+
}
fFissionXS.insert(pair<ZAI ,TGraph* >(zaitmp, new TGraph(vTime.size()-1, Time, DPQuantity) ) );
}
@@ -1032,27 +1032,27 @@ void EvolutionData::OldReadDB(string DBfile)
if(line == "" || line == "Capture" ) break;
}while ( !DecayDB.eof() );
}
-
+
if (line == "Capture")
{
getline(DecayDB, line); // Nuclei is given with "A Z"
-
+
do
{
double DPQuantity[vTime.size()];
for(int k = 0; k < (int)vTime.size(); k++ )
DPQuantity[k] = 0;
-
-
+
+
start = 0;
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
-
+
if(A!=0 && Z!=0)
{
-
-
+
+
ZAI zaitmp(Z, A, I);
int i=0;
while(start < (int)line.size())
@@ -1060,36 +1060,36 @@ void EvolutionData::OldReadDB(string DBfile)
long double DPQuantityTmp = atof(StringLine::NextWord(line, start, ' ').c_str());
DPQuantity[i] = (double)DPQuantityTmp;
i++;
-
+
}
fCaptureXS.insert(pair<ZAI ,TGraph* >(zaitmp, new TGraph(vTime.size()-1, Time, DPQuantity) ) );
}
getline(DecayDB, line); // Nuclei is given with "A Z"
if(line == "" || line == "n2n" ) break;
}while ( !DecayDB.eof() );
-
+
}
-
+
if (line == "n2n")
{
-
+
getline(DecayDB, line); // Nuclei is given with "A Z"
-
+
do
{
double DPQuantity[vTime.size()];
for(int k = 0; k < (int)vTime.size(); k++ )
DPQuantity[k] = 0;
-
+
start = 0;
int Z = atoi(StringLine::NextWord(line, start, ' ').c_str());
int A = atoi(StringLine::NextWord(line, start, ' ').c_str());
int I = atoi(StringLine::NextWord(line, start, ' ').c_str());
-
+
if(A!=0 && Z!=0)
{
-
-
+
+
ZAI zaitmp(Z, A, I);
int i=0;
while(start < (int)line.size())
@@ -1097,20 +1097,20 @@ void EvolutionData::OldReadDB(string DBfile)
long double DPQuantityTmp = atof(StringLine::NextWord(line, start, ' ').c_str());
DPQuantity[i] = (double)DPQuantityTmp;
i++;
-
+
}
fn2nXS.insert(pair<ZAI ,TGraph* >(zaitmp, new TGraph(vTime.size()-1, Time, DPQuantity) ) );
}
getline(DecayDB, line); // Nuclei is given with "A Z"
if(line == "" ) break;
-
+
}while ( !DecayDB.eof() );
}
-
+
}
DecayDB.close();
start = 0;
-
+
string InfoDBFile = DBfile.erase(DBfile.size()-3,DBfile.size());
InfoDBFile += "info";
ifstream InfoDB(InfoDBFile.c_str()); // Open the File
@@ -1119,7 +1119,7 @@ void EvolutionData::OldReadDB(string DBfile)
INFO << " Can't open \"" << InfoDBFile << "\"\n" << endl;
return;
}
-
+
start = 0;
getline(InfoDB, line);
if (StringLine::NextWord(line, start, ' ') == "Reactor")
@@ -1148,7 +1148,7 @@ void EvolutionData::OldReadDB(string DBfile)
double Flux[vFlux.size()];
for(int i=0; i < (int)vFlux.size();i++)
Flux[i] = vFlux[i];
-
+
fFlux = new TGraph(vTime.size()-1, Time, Flux);
}
InfoDB.close();
diff --git a/source/branches/BaM_Dev/src/XSModel.cxx b/source/branches/BaM_Dev/src/XSModel.cxx
index 1522534f530316e928731576658ce14480ee2661..d83064ca38c775368571797bc1e05c2e080ba02a 100644
--- a/source/branches/BaM_Dev/src/XSModel.cxx
+++ b/source/branches/BaM_Dev/src/XSModel.cxx
@@ -17,23 +17,71 @@ using namespace std;
XSModel::XSModel():CLASSObject()
{
- LoadKeyword();
+ XSModel::LoadKeyword();
}
XSModel::XSModel(CLASSLogger* log):CLASSObject(log)
{
- LoadKeyword();
-}
+
+ XSModel::LoadKeyword();
+}
+
+
+void XSModel::ReadNFO()
+{
+ DBGL
+ ifstream NFO(fInformationFile.c_str());
+
+ if(!NFO)
+ {
+ ERROR << "Can't find/open file " << fInformationFile << endl;
+ exit(0);
+ }
+
+ do
+ {
+ string line;
+ getline(NFO,line);
+
+ XSModel::ReadLine(line);
+
+ } while(NFO.eof());
+
+ DBGL
+}
+
+//________________________________________________________________________
+void XSModel::ReadLine(string line)
+{
+ DBGL
+
+ if (!freaded)
+ {
+ int start = 0;
+ string keyword = tlc(StringLine::NextWord(line, start, ' '));
+
+ (this->*fKeyword[ keyword ])(line);
+ freaded = true;
+
+ ReadLine(line);
+
+ }
+
+ freaded = false;
+
+ DBGL
+}
+
void XSModel::LoadKeyword()
{
DBGL
fKeyword.insert( pair<string, MthPtr>( "k_zail", & XSModel::ReadZAIlimits));
- fKeyword.insert( pair<string, MthPtr>( "k_reactor", & XSModel::ReadType) );
- fKeyword.insert( pair<string, MthPtr>( "k_fuel", & XSModel::ReadType) );
- fKeyword.insert( pair<string, MthPtr>( "k_mass", & XSModel::ReadRParam) );
- fKeyword.insert( pair<string, MthPtr>( "k_power", & XSModel::ReadRParam) );
+ fKeyword.insert( pair<string, MthPtr>( "k_reactor", & XSModel::ReadType) );
+ fKeyword.insert( pair<string, MthPtr>( "k_fuel", & XSModel::ReadType) );
+ fKeyword.insert( pair<string, MthPtr>( "k_mass", & XSModel::ReadRParam) );
+ fKeyword.insert( pair<string, MthPtr>( "k_power", & XSModel::ReadRParam) );
DBGL
}
@@ -41,15 +89,15 @@ void XSModel::ReadRParam(const string &line)
{
DBGL
int start = 0;
- string type = tlc(StringLine::NextWord(line, start, ' '));
- if( type != "k_power" || type != "k_mass" ) // Check the keyword
+ string keyword = tlc(StringLine::NextWord(line, start, ' '));
+ if( keyword != "k_power" || keyword != "k_mass" ) // Check the keyword
{
- ERROR << " Bad keyword : " << type << " Not found !" << endl;
+ ERROR << " Bad keyword : " << keyword << " Not found !" << endl;
exit(1);
}
- if( type == "k_mass" )
+ if( keyword == "k_mass" )
fDBHMMass = atof(StringLine::NextWord(line, start, ' ').c_str());
- else if( type == "k_mass" )
+ else if( keyword == "k_mass" )
fDBPower = atof(StringLine::NextWord(line, start, ' ').c_str());
DBGL
@@ -60,17 +108,17 @@ void XSModel::ReadType(const string &line)
{
DBGL
int start = 0;
- string type = tlc(StringLine::NextWord(line, start, ' '));
- if( type != "k_fuel" || type != "k_reactor" ) // Check the keyword
+ string keyword = tlc(StringLine::NextWord(line, start, ' '));
+ if( keyword != "k_fuel" || keyword != "k_reactor" ) // Check the keyword
{
- ERROR << " Bad keyword : " << type << " Not found !" << endl;
+ ERROR << " Bad keyword : " << keyword << " Not found !" << endl;
exit(1);
}
- if( type == "k_fuel" )
+ if( keyword == "k_fuel" )
fDBFType = StringLine::NextWord(line, start, ' ');
- else if( type == "k_reactor" )
+ else if( keyword == "k_reactor" )
fDBRType = StringLine::NextWord(line, start, ' ');
-
+
DBGL
}
@@ -79,7 +127,8 @@ void XSModel::ReadZAIlimits(const string &line)
{
DBGL
int start = 0;
- if( tlc(StringLine::NextWord(line, start, ' ')) != "k_zail" ) // Check the keyword
+ string keyword = tlc(StringLine::NextWord(line, start, ' '));
+ if( keyword != "k_zail" ) // Check the keyword
{
ERROR << " Bad keyword : \"k_zail\" not found !" << endl;
exit(1);
@@ -98,9 +147,9 @@ void XSModel::ReadZAIlimits(const string &line)
bool XSModel::isIVInDomain(IsotopicVector IV)
{
-DBGL
+ DBGL
bool IsInDomain=true;
-
+
if(fZAILimits.empty())
{
WARNING << "Fresh Fuel variation domain is not set" << endl;
@@ -108,19 +157,19 @@ DBGL
WARNING << "Proceed finger crossed !!" << endl;
return true;
}
-
+
else
- {
+ {
IsotopicVector IVNorm = IV /IV.GetSumOfAll();
for (map< ZAI,pair<double,double> >::iterator Domain_it=fZAILimits.begin(); Domain_it!=fZAILimits.end(); Domain_it++)
- {
+ {
double ThatZAIProp = IVNorm.GetIsotopicQuantity()[Domain_it->first] ;
double ThatZAIMin = Domain_it->second.first;
double ThatZAIMax = Domain_it->second.second;
if( (ThatZAIProp > ThatZAIMax) || (ThatZAIProp < ThatZAIMin) )
- {
- IsInDomain = false;
-
+ {
+ IsInDomain = false;
+
WARNING<<"Fresh fuel out of model range"<<endl;
WARNING<<"\t AT LEAST this ZAI is accused to be outrange :"<<endl;
WARNING<<"\t\t"<<Domain_it->first.Z()<<" "<<Domain_it->first.A()<<" "<<Domain_it->first.I()<<endl;
@@ -129,9 +178,9 @@ DBGL
WARNING<<IVNorm.sPrint()<<endl;
break;
}
- }
+ }
}
-DBGL
-return IsInDomain;
-
+ DBGL
+ return IsInDomain;
+
}
\ No newline at end of file