train.C File Reference

#include <cstdlib>
#include <iostream>
#include <map>
#include <string>
#include "TChain.h"
#include "TFile.h"
#include "TTree.h"
#include "TString.h"
#include "TObjString.h"
#include "TSystem.h"
#include "TROOT.h"
#include "TMVA/Tools.h"
#include "TMVA/Factory.h"
#include "TMVA/DataLoader.h"

Go to the source code of this file.

Functions
void	train ()

Function Documentation

void train

(

)

Definition at line 24 of file train.C.

References om::cout, train_bdt::dataloader, allTimeWatchdog::endl, train_bdt::factory, file, lt, pot, string, and art::to_string().

Referenced by cvn::CVNEventDump::analyze(), cvn::CVNCosmicEventDump::analyze(), cvn::CVNEventProngDump::analyze(), and cvn::CVNAddTrainingData::produce().

{
  TFile *file = new TFile("files/hadded_all_files.root");

  // Need to balance the beam MC / cosmic weights using livetimes.
  // Use recent estimates for each
  // Note, overall exposure scale doesn't affect training, but lt/POT does
  TFile *filepot = new TFile("files/hadded_fluxswap_files.root");

  double pot = 9.47964e20;
  double lt  = 438.1584;

  TH1 *hPOT = (TH1F*)filepot->Get("pot");
  double totPOT = hPOT->Integral();
  cout<<"pot from fluxswap is "<<totPOT<<endl;
  delete hPOT;
  TH1 *hLt = (TH1F*)file->Get("livetime");
  double totLt = hLt->Integral();
  delete hLt;

  double nuScale = pot/totPOT;
  double cosScale= lt/totLt;

  cout<<" scale for mc "<<nuScale;
  cout<<" scale for cosmic "<<cosScale;

  TTree *sigTree = (TTree*)file->Get("sigTree");
  TTree *bkgTree = (TTree*)file->Get("bkgTree");

  TFile* outputFile = TFile::Open("files/test_current_bdt.root","recreate");
  TMVA::Factory *factory = new TMVA::Factory("MVAnalysis",outputFile,"V");
  TMVA::DataLoader *dataloader=new TMVA::DataLoader("dataset");

  dataloader->AddVariable("nhit",       'F');
  dataloader->AddVariable("cvnsse",       'F');
  //dataloader->AddVariable("cvne",       'F');
  dataloader->AddVariable("sparseassym",'F');
  //dataloader->AddVariable("pxp",        'F');
  //dataloader->AddVariable("pyp",        'F');
  dataloader->AddVariable("ptp",        'F');
  dataloader->AddVariable("disttop",    'F');
  dataloader->AddVariable("distnottop", 'F');
  //dataloader->AddVariable("vtxx",       'F'); 
  //dataloader->AddVariable("vtxy",       'F');
  //dataloader->AddVariable("vtxz",       'F');
  //dataloader->AddVariable("prongl",       'F');
  //dataloader->AddVariable("inelast",    'F');
  //dataloader->AddVariable("hitsperpl",  'F');
  //dataloader->AddVariable("numofshow",  'I');
  //dataloader->AddVariable("widthofshow",'F');
  //dataloader->AddVariable("costheta",   'F');

  // These are questionable
  //dataloader->AddVariable("calE",       'F');
  //dataloader->AddVariable("costhetaZ",  'F');
  //dataloader->AddVariable("prongl",     'F');
  //dataloader->AddVariable("remid",      'F');

  //dataloader->AddVariable("distfront", 'F');
  //dataloader->AddVariable("distback",  'F');
  //dataloader->AddVariable("disteast",  'F');
  //dataloader->AddVariable("cosang",  'F');
  //dataloader->AddVariable("vtxz",  'F');
  // dataloader->AddVariable("vtxx",  'F');
  // dataloader->AddVariable("vtxy",  'F');
  // dataloader->AddVariable("vtxz",  'F');
  //dataloader->AddVariable("efrac",  'F');
  //dataloader->AddVariable("nplfr",  'F');
  //dataloader->AddVariable("maxy",  'F');
  //dataloader->AddVariable("distwest",  'F');
  //dataloader->AddVariable("distbottom",  'F');
  //dataloader->AddVariable("lid",  'F');
  //dataloader->AddVariable("starty",  'F');
  //dataloader->AddVariable("stopy",  'F');

  // Dumb oscillations, for ok 1st order prediction estimate
  dataloader->AddSpectator("woscdumb",  'F');
  // All necessary for determining signal / cosmic bkg
  dataloader->AddSpectator("isGENIE",   'F');
  dataloader->AddSpectator("pdg",       'F');
  dataloader->AddSpectator("pdgorig",   'F');
  dataloader->AddSpectator("isCC",      'F');


  // Define a preselection here, may well influence results
  // Also add in a truth cut to reject beam bkg from training
  TCut presel = "cvnsse>0.5 && prongl > 0";
  //dataloader->AddCut(presel);

  //TCut sigCut = "pdg==-12 && pdgorig==-14";
  //TCut bkgCut = "pdg!=-12 || pdgorig!=-14"

  dataloader->AddSignalTree    (sigTree);
  dataloader->AddBackgroundTree(bkgTree);

  dataloader->AddCut(presel);

  // Apply osc weight to beam events, scale cosmics by lt
  std::string weightExp = "woscdumb*" + std::to_string(nuScale) + "*isGENIE*" + "(pdg==12&&pdgorig==14)+" + std::to_string(cosScale) + "*!isGENIE";
  std::cout << weightExp << std::endl;
  dataloader->SetWeightExpression(weightExp.c_str());

  dataloader->PrepareTrainingAndTestTree("", "SplitMode=random");

  // A few methods here, to easily swap out and try
  //factory->BookMethod(TMVA::Types::kCuts, "Cuts","H:!V:FitMethod=MC:EffSel:SampleSize=200000:VarProp=FSmart" );
  //factory->BookMethod(TMVA::Types::kLikelihood, "Likelihood", "V");
  //factory->BookMethod(TMVA::Types::kLD, "LD", "H:!V:VarTransform=None:CreateMVAPdfs:PDFInterpolMVAPdf=Spline2:NbinsMVAPdf=50:NsmoothMVAPdf=10");
  // This was the method we settled on in 2016, and retrained with in 2017
  factory->BookMethod(dataloader, TMVA::Types::kBDT, "test_current_bdt","H:V:NTrees=500:InverseBoostNegWeights::BoostType=RealAdaBoost:Shrinkage=0.1:nCuts=40:MaxDepth=5:MinNodeSize=2%:NodePurityLimit=0.9");
  //factory->BookMethod( TMVA::Types::kMLP, "MLP", "H:!V:NeuronType=tanh:VarTransform=N:NCycles=600:HiddenLayers=N+5:TestRate=5:!UseRegulator" );


  //  factory->BookMethod( TMVA::Types::kBDT, "BDTG",
  //                       "H:V:NTrees=4000:InverseBoostNegWeights::BoostType=Grad:Shrinkage=0.15:nCuts=40:MaxDepth=30:MinNodeSize=2%:GradBaggingFraction=0.5:NodePurityLimit=0.9");

  factory->TrainAllMethods();
  factory->TestAllMethods();
  factory->EvaluateAllMethods();

  outputFile->Close();
  delete factory;
  delete dataloader;
}