|
#include "EquivalenceModel.hxx"
|
|
#include "EQM_CANDU2.hxx"
|
|
#include "CLASSLogger.hxx"
|
|
#include "StringLine.hxx"
|
|
|
|
#include <string>
|
|
#include <iostream>
|
|
#include <fstream>
|
|
#include <algorithm>
|
|
#include <cmath>
|
|
#include <cassert>
|
|
|
|
#include "TSystem.h"
|
|
#include "TMVA/Reader.h"
|
|
#include "TMVA/Tools.h"
|
|
#include "TMVA/MethodCuts.h"
|
|
|
|
|
|
//________________________________________________________________________
|
|
//
|
|
// EQM_CANDU2
|
|
//
|
|
// Equivalenve Model based on multi layer perceptron from TMVA (root cern)
|
|
// For REP MOX use
|
|
//
|
|
//________________________________________________________________________
|
|
|
|
EQM_CANDU2::EQM_CANDU2(string TMVAWeightPath):EquivalenceModel(new CLASSLogger("EQM_CANDU2.log"))
|
|
{
|
|
|
|
WARNING << "yes you are a CANDU2 first constructor"<<endl;
|
|
fTMVAWeightPath = TMVAWeightPath;
|
|
|
|
ZAI U8(92,238,0);
|
|
ZAI U5(92,235,0);
|
|
ZAI U6(92,236,0);
|
|
double U5_enrich = 0.0025;
|
|
|
|
ZAI Pu8(94,238,0);
|
|
ZAI Pu9(94,239,0);
|
|
ZAI Pu0(94,240,0);
|
|
ZAI Pu1(94,241,0);
|
|
ZAI Pu2(94,242,0);
|
|
|
|
fFissileList = U5*1+U6*1+U8*1;
|
|
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
|
|
|
|
SetBuildFuelFirstGuess(0.04);
|
|
|
|
INFO << "__An equivalence model of PWR MOX has been define__" << endl;
|
|
INFO << "\tThis model is based on a multi layer perceptron" << endl;
|
|
INFO << "\t\tThe TMVA weight file is :" << endl;
|
|
INFO << "\t\t\t" << fTMVAWeightPath << endl;
|
|
EquivalenceModel::PrintInfo();
|
|
|
|
}
|
|
|
|
//________________________________________________________________________
|
|
EQM_CANDU2::EQM_CANDU2(CLASSLogger* log, string TMVAWeightPath):EquivalenceModel(log)
|
|
{
|
|
WARNING << "yes you are a CANDU2 second constructor"<<endl;
|
|
fTMVAWeightPath = TMVAWeightPath;
|
|
|
|
ZAI U8(92,238,0);
|
|
ZAI U2(92,232,0);
|
|
ZAI U3(92,233,0);
|
|
ZAI U4(92,234,0);
|
|
ZAI U5(92,235,0);
|
|
ZAI U6(92,236,0);
|
|
double U5_enrich = 0.0025;
|
|
|
|
ZAI Pu8(94,238,0);
|
|
ZAI Pu9(94,239,0);
|
|
ZAI Pu0(94,240,0);
|
|
ZAI Pu1(94,241,0);
|
|
ZAI Pu2(94,242,0);
|
|
|
|
// fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
|
|
|
|
fFissileList = U2*1+U3*1+U4*1+U5*1+U6*1+U8*1;
|
|
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
|
|
|
|
SetBuildFuelFirstGuess(0.04);
|
|
|
|
INFO << "__An equivalence model of PWR MOX has been defined__" << endl;
|
|
INFO << "\tThis model is based on a multi layer perceptron" << endl;
|
|
INFO << "\t\tThe TMVA weight file is :" << endl;
|
|
INFO << "\t\t\t" << fTMVAWeightPath << endl;
|
|
EquivalenceModel::PrintInfo();
|
|
|
|
}
|
|
|
|
//________________________________________________________________________
|
|
TTree* EQM_CANDU2::CreateTMVAInputTree(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp)
|
|
{
|
|
TTree* InputTree = new TTree("EQTMP", "EQTMP");
|
|
float Pu8 = 0;
|
|
float Pu9 = 0;
|
|
float Pu10 = 0;
|
|
float Pu11 = 0;
|
|
float Pu12 = 0;
|
|
float Am1 = 0;
|
|
float U5_enrichment = 0;
|
|
float BU = 0;
|
|
|
|
InputTree->Branch( "Pu8" ,&Pu8 ,"Pu8/F" );
|
|
InputTree->Branch( "Pu9" ,&Pu9 ,"Pu9/F" );
|
|
InputTree->Branch( "Pu10" ,&Pu10 ,"Pu10/F" );
|
|
InputTree->Branch( "Pu11" ,&Pu11 ,"Pu11/F" );
|
|
InputTree->Branch( "Pu12" ,&Pu12 ,"Pu12/F" );
|
|
InputTree->Branch( "Am1" ,&Am1 ,"Am1/F" );
|
|
InputTree->Branch( "U5_enrichment" ,&U5_enrichment ,"U5_enrichment/F" );
|
|
InputTree->Branch( "BU" ,&BU ,"BU/F" );
|
|
|
|
|
|
float U8 = Fertil.GetZAIIsotopicQuantity(92,238,0);
|
|
float U5 = Fertil.GetZAIIsotopicQuantity(92,235,0);
|
|
float U4 = Fertil.GetZAIIsotopicQuantity(92,234,0);
|
|
|
|
float UTOT = U8 + U5 + U4;
|
|
|
|
Pu8 = Fissil.GetZAIIsotopicQuantity(94,238,0);
|
|
Pu9 = Fissil.GetZAIIsotopicQuantity(94,239,0);
|
|
Pu10 = Fissil.GetZAIIsotopicQuantity(94,240,0);
|
|
Pu11 = Fissil.GetZAIIsotopicQuantity(94,241,0);
|
|
Pu12 = Fissil.GetZAIIsotopicQuantity(94,242,0);
|
|
Am1 = Fissil.GetZAIIsotopicQuantity(95,241,0);
|
|
|
|
double TOTPU = (Pu8+Pu9+Pu10+Pu11+Pu12+Am1);
|
|
|
|
Pu8 = Pu8 / TOTPU;
|
|
Pu9 = Pu9 / TOTPU;
|
|
Pu10 = Pu10 / TOTPU;
|
|
Pu11 = Pu11 / TOTPU;
|
|
Pu12 = Pu12 / TOTPU;
|
|
Am1 = Am1 / TOTPU;
|
|
|
|
U5_enrichment = U5 / UTOT;
|
|
|
|
BU = BurnUp;
|
|
if(Pu8 + Pu9 + Pu10 + Pu11 + Pu12 + Am1 > 1.00001 )//?????1.00001??? I don't know it! goes in condition if = 1 !! may be float/double issue ...
|
|
{
|
|
ERROR << Pu8 << " " << Pu9 << " " << Pu10 << " " << Pu11 << " " << Pu12 << " " << Am1 << endl;
|
|
exit(0);
|
|
}
|
|
// All value are molar (!weight)
|
|
|
|
InputTree->Fill();
|
|
return InputTree;
|
|
}
|
|
//________________________________________________________________________
|
|
double EQM_CANDU2::ExecuteTMVA(TTree* theTree)
|
|
{
|
|
// --- 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
|
|
Float_t Pu8,Pu9,Pu10,Pu11,Pu12,Am1,BU,U5_enrichment;
|
|
|
|
reader->AddVariable( "BU" ,&BU );
|
|
reader->AddVariable( "U5_enrichment",&U5_enrichment );
|
|
reader->AddVariable( "Pu8" ,&Pu8 );
|
|
reader->AddVariable( "Pu9" ,&Pu9 );
|
|
reader->AddVariable( "Pu10" ,&Pu10);
|
|
reader->AddVariable( "Pu11" ,&Pu11);
|
|
reader->AddVariable( "Pu12" ,&Pu12);
|
|
reader->AddVariable( "Am1" ,&Am1 );
|
|
|
|
// --- Book the MVA methods
|
|
|
|
// Book method MLP
|
|
TString methodName = "MLP method";
|
|
reader->BookMVA( methodName, fTMVAWeightPath );
|
|
theTree->SetBranchAddress( "BU" ,&BU );
|
|
theTree->SetBranchAddress( "U5_enrichment" ,&U5_enrichment ) ;
|
|
theTree->SetBranchAddress( "Pu8" ,&Pu8 );
|
|
theTree->SetBranchAddress( "Pu9" ,&Pu9 );
|
|
theTree->SetBranchAddress( "Pu10" ,&Pu10 );
|
|
theTree->SetBranchAddress( "Pu11" ,&Pu11 );
|
|
theTree->SetBranchAddress( "Pu12" ,&Pu12 );
|
|
theTree->SetBranchAddress( "Am1" ,&Am1 );
|
|
theTree->GetEntry(0);
|
|
|
|
Float_t val = (reader->EvaluateRegression( methodName ))[0];
|
|
|
|
delete reader;
|
|
delete theTree;
|
|
|
|
return (double)val; //retourne teneur
|
|
}
|
|
//________________________________________________________________________
|
|
double EQM_CANDU2::GetFissileMolarFraction(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp)
|
|
{DBGL
|
|
//double retourTMVA=ExecuteTMVA(CreateTMVAInputTree(Fissil,Fertil,BurnUp));
|
|
//INFO<<"FissileMolarFraction ="<<retourTMVA<<endl;
|
|
return 1.0;
|
|
}
|