Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:31:15

0001 #include <iostream>
0002 #include <cassert>
0003 #include <TROOT.h>
0004 #include <TSystem.h>
0005 #include <TFile.h>
0006 #include <TTree.h>
0007 #include <TBranch.h>
0008 #include <TCanvas.h>
0009 #include <TH1.h>
0010 #include <TH2.h>
0011 #include <TGraph.h>
0012 #include <TF1.h>
0013 #include <TFormula.h>
0014 #include <TStyle.h>
0015 #include <TKey.h>
0016 #include <vector>
0017 #include "FWCore/FWLite/interface/FWLiteEnabler.h"
0018 #include "AnalysisDataFormats/TopObjects/interface/TtSemiEvtSolution.h"
0019 #include "TopQuarkAnalysis/TopTools/interface/LRHelpFunctions.h"
0020 
0021 ///////////////////////
0022 // Constants         //
0023 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
0024 
0025 //input files
0026 const int signal_nrDir = 5;
0027 const int signal_nrFiles[signal_nrDir] = {30, 25, 20, 15, 10};
0028 const TString signal_path[signal_nrDir] = {
0029     "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt0j/Alpgen_tt0j_TtSemiMuEvents_",
0030     "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt1j/Alpgen_tt1j_TtSemiMuEvents_",
0031     "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt2j/Alpgen_tt2j_TtSemiMuEvents_",
0032     "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt3j/Alpgen_tt3j_TtSemiMuEvents_",
0033     "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt4j/Alpgen_tt4j_TtSemiMuEvents_"};
0034 const int signal_NrEv[signal_nrDir] = {76250, 68000, 40000, 16000, 24425};  //nb of events you want to process
0035 
0036 const int bckgd_nrDir = 8;
0037 const int bckgd_nrFiles[bckgd_nrDir] = {30, 25, 20, 15, 10, 10, 10, 5};
0038 const TString bckgd_path[bckgd_nrDir] = {"dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt0j/Alpgen_tt0j_TtOtherEvents_",
0039                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt1j/Alpgen_tt1j_TtOtherEvents_",
0040                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt2j/Alpgen_tt2j_TtOtherEvents_",
0041                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt3j/Alpgen_tt3j_TtOtherEvents_",
0042                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/tt4j/Alpgen_tt4j_TtOtherEvents_",
0043                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/W4j/Alpgen_W4j_TtSemiMuEvents_",
0044                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/W5j/Alpgen_W5j_TtSemiMuEvents_",
0045                                          "dcap:///pnfs/iihe/becms/ghammad/TQAF136Final/W6j/Alpgen_W6j_TtSemiMuEvents_"};
0046 
0047 const int bckgd_NrEv[bckgd_nrDir] = {
0048     76250, 68000, 40000, 16000, 24425, 69700, 18000, 12500};  //nb of events you want to process
0049 
0050 //observable histogram variables
0051 const int nrSignalSelObs = 38;
0052 const int SignalSelObs[nrSignalSelObs] = {1,  2,  3,  4,  5,  6,  7,  8,  9,  10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
0053                                           20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38};
0054 const int nrSignalSelHistBins = 35;
0055 const double SignalSelObsMin[nrSignalSelObs] = {10, 0,   0, 0, 0, 0.5, 0, -5, -2, -10, 0, 0,    0, 0.15, 0, 0, 0, 0, 0,
0056                                                 0,  0.1, 0, 0, 0, 0,   0, 0,  0,  0,   0, 0.55, 0, 0,    0, 0, 0, 0, 0};
0057 const double SignalSelObsMax[nrSignalSelObs] = {200, 200, 1000, 300,  0.5, 1,    1600, 5,   20,  15,   50, 1,  1,
0058                                                 1,   1,   0.3,  0.9,  0.3, 0.35, 1000, 500, 0.8, 70,   2,  1,  1,
0059                                                 0.8, 1,   0.65, 0.75, 1,   0.7,  0.65, 0.6, 1,   0.35, 1,  0.4};
0060 
0061 //observable fit functions
0062 
0063 TFormula gauss("gauss", "gaus");
0064 TFormula symgauss("symgauss", "[0]*(exp(-0.5*(x/[1])**2))");
0065 TFormula dblgauss("dblgauss", "[0]*(exp(-0.5*((x-[1])/[2])**2)+exp(-0.5*((x+[3])/[4])**2))");
0066 TFormula symdblgauss("symdblgauss", "[0]*(exp(-0.5*((x-[1])/[2])**2)+exp(-0.5*((x+[1])/[2])**2))");
0067 TFormula sigm("sigm", "[0]/(1 + 1/exp([1]*([2] - x)))");
0068 TFormula sigmc("sigmc", "[0]/(1 + 1/exp([1]*([2] - x)))+[3]");
0069 TFormula dblsigm("dblsigm", "[0]/(1 + 1/exp([1]**2*([2] - x)))/(1 + 1/exp([3]**2*(x - [4])))");
0070 TFormula symdblsigm("symdblsigm", "[0]/(1 + 1/exp([1]**2*([2] - x)))/(1 + 1/exp([1]**2*([2] + x)))");
0071 
0072 const char *SignalSelObsFits[nrSignalSelObs] = {
0073 
0074     "[0]/(1 + 1/exp([1]*([2] - x)))",  //"[0]*(1-exp(-[1]*x))-[2]*(1-exp(-[3]*x))",//"[0]*exp(-pow((x-[1])/[2],2))+[3]*(1-exp(-[4]*x))", //obs1
0075     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //"[0]/(1 + 1/exp([1]*([2] - x)))*(exp(-[3]*x))", //obs2
0076     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4",                //obs3
0077     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //obs4
0078     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs5
0079     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]+pol3(4)+pol5(8)",     //obs6
0080     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4",                //obs7
0081     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol3+[4]",            //obs8
0082     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol3+[4]",            //obs9
0083     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //obs10
0084     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //obs11
0085     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //obs12
0086     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]",                     //obs13
0087     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs14
0088     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol3",                //obs15
0089     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4",                //obs16
0090     "pol9",                                                   //"[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4", //obs17
0091     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4",                //obs18
0092     "pol6",                                                   //"[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4", //obs19
0093     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*exp([4]*([5] - x))",  //obs20
0094     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*pol4",                //obs21
0095     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs22
0096     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs23
0097     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs24
0098     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs25
0099     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs26
0100     "pol6",                                                   //"[0]/(1 + 1/exp([1]*([2] - x)))", //obs27
0101     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs28
0102     "[0]/(1 + 1/exp([1]*([2] - x)))+[3]*exp([4]*([5] - x))",  //obs29
0103     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs30
0104     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs31
0105     "pol6",                                                   //"[0]/(1 + 1/exp([1]*([2] - x)))", //obs32
0106     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs33
0107     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs34
0108     "pol3+[3]*exp([4]*([5] - x))",                            //"[0]/(1 + 1/exp([1]*([2] - x)))", //obs35
0109     "[0]/(1 + 1/exp([1]*([2] - x)))",                         //obs36
0110     "pol3+[3]*exp([4]*([5] - x))",                            //"[0]/(1 + 1/exp([1]*([2] - x)))", //obs37
0111     "pol5+[3]*exp([4]*([5] - x))",                            //"[0]/(1 + 1/exp([1]*([2] - x)))", //obs38
0112                                                               //                             "pol3+[3]*exp([4]*([5] - x))",    //obs39
0113                                                               //                             "[0]/(1 + 1/exp([1]*([2] - x)))", //obs40
0114                                                               //                             "[0]/(1 + 1/exp([1]*([2] - x)))", //obs41
0115                                                               //                             "[0]/(1 + 1/exp([1]*([2] - x)))", //obs42
0116                                                               //                             "[0]/(1 + 1/exp([1]*([2] - x)))"  //obs43
0117 };
0118 
0119 //output files ps/root
0120 const TString SignalSelOutfileName = "./TtSemiLRSignalSelAllObs.root";
0121 const TString SignalSelPSfile = "./TtSemiLRSignalSelAllObs.ps";
0122 
0123 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
0124 
0125 //
0126 // Global variables
0127 //
0128 LRHelpFunctions *myLRhelper;
0129 void doEventloop();
0130 std::vector<int> obsNrs;
0131 std::vector<double> obsMin, obsMax;
0132 std::vector<const char *> obsFits;
0133 
0134 bool MuonIso = true;
0135 
0136 //
0137 // Main analysis
0138 //
0139 
0140 int main() {
0141   gSystem->Load("libFWCoreFWLite");
0142   FWLiteEnabler::enable();
0143 
0144   // define all histograms & fit functions
0145   //to replace with something more elegant
0146   for (int j = 0; j < nrSignalSelObs; j++) {
0147     obsNrs.push_back(SignalSelObs[j]);
0148     obsMin.push_back(SignalSelObsMin[j]);
0149     obsMax.push_back(SignalSelObsMax[j]);
0150     obsFits.push_back(SignalSelObsFits[j]);
0151   }
0152   myLRhelper = new LRHelpFunctions(obsNrs, nrSignalSelHistBins, obsMin, obsMax, obsFits);
0153 
0154   // manually set some initial values for fit function parameters
0155   std::vector<double> parsFobs1;
0156   parsFobs1.push_back(20);
0157   parsFobs1.push_back(0.04);
0158   parsFobs1.push_back(21);
0159   parsFobs1.push_back(0.04);
0160   myLRhelper->setObsFitParameters(1, parsFobs1);
0161 
0162   std::vector<double> parsFobs2;
0163   parsFobs2.push_back(0.495);
0164   parsFobs2.push_back(-0.148);
0165   parsFobs2.push_back(60.33);
0166   parsFobs2.push_back(0.396);  //parsFobs2.push_back(0.03);
0167   myLRhelper->setObsFitParameters(2, parsFobs2);
0168 
0169   std::vector<double> parsFobs3;
0170   parsFobs3.push_back(7.60418e-01);
0171   parsFobs3.push_back(-3.31635e-02);
0172   parsFobs3.push_back(1.57387e+02);
0173   parsFobs3.push_back(-1.23931e-08);
0174   parsFobs3.push_back(-1.90918e-04);
0175   myLRhelper->setObsFitParameters(3, parsFobs3);
0176 
0177   std::vector<double> parsFobs4;
0178   parsFobs4.push_back(1.087);
0179   parsFobs4.push_back(-0.1978);
0180   parsFobs4.push_back(22.803);
0181   parsFobs4.push_back(-0.126);  //parsFobs4.push_back(0.04);
0182   myLRhelper->setObsFitParameters(4, parsFobs4);
0183 
0184   std::vector<double> parsFobs7;
0185   parsFobs7.push_back(0.606878);
0186   parsFobs7.push_back(-1.52796e-02);
0187   parsFobs7.push_back(2.50574e+02);
0188   parsFobs7.push_back(-4.46936e-10);
0189   parsFobs7.push_back(-1.48804e-04);
0190   myLRhelper->setObsFitParameters(7, parsFobs7);
0191 
0192   std::vector<double> parsFobs8;
0193   parsFobs8.push_back(3.30611e-01);
0194   parsFobs8.push_back(-8.34406e+00);
0195   parsFobs8.push_back(1.04307e+00);
0196   parsFobs8.push_back(-1.75190e-03);
0197   parsFobs8.push_back(5.66972e-01);
0198   myLRhelper->setObsFitParameters(8, parsFobs8);
0199 
0200   std::vector<double> parsFobs9;
0201   parsFobs9.push_back(6.37793e-01);
0202   parsFobs9.push_back(-1.71768e+00);
0203   parsFobs9.push_back(1.88952e+00);
0204   parsFobs9.push_back(-1.03833e-03);
0205   parsFobs9.push_back(3.30284e-01);
0206   myLRhelper->setObsFitParameters(9, parsFobs9);
0207 
0208   std::vector<double> parsFobs10;
0209   parsFobs10.push_back(0.618);
0210   parsFobs10.push_back(-1.579);
0211   parsFobs10.push_back(-0.10812);
0212   parsFobs10.push_back(0.342);
0213   myLRhelper->setObsFitParameters(10, parsFobs10);
0214 
0215   std::vector<double> parsFobs11;
0216   parsFobs11.push_back(0.7624);
0217   parsFobs11.push_back(-0.64975);
0218   parsFobs11.push_back(3.1225);
0219   parsFobs11.push_back(0.218675);
0220   myLRhelper->setObsFitParameters(11, parsFobs11);
0221 
0222   std::vector<double> parsFobs12;
0223   parsFobs12.push_back(1.57736e-01);
0224   parsFobs12.push_back(-2.01467e+01);
0225   parsFobs12.push_back(5.97867e-01);
0226   parsFobs12.push_back(3.81101e-01);
0227   myLRhelper->setObsFitParameters(12, parsFobs12);
0228 
0229   std::vector<double> parsFobs13;
0230   parsFobs13.push_back(1.57736e-01);
0231   parsFobs13.push_back(-2.01467e+01);
0232   parsFobs13.push_back(5.97867e-01);
0233   parsFobs13.push_back(3.81101e-01);
0234   myLRhelper->setObsFitParameters(13, parsFobs13);
0235 
0236   std::vector<double> parsFobs15;
0237   parsFobs15.push_back(0.6672);
0238   parsFobs15.push_back(-9.3022);
0239   parsFobs15.push_back(0.03384);
0240   parsFobs15.push_back(0.00014967);
0241   parsFobs15.push_back(-4315.96);
0242   myLRhelper->setObsFitParameters(15, parsFobs15);
0243 
0244   std::vector<double> parsFobs16;
0245   parsFobs16.push_back(0.56855);
0246   parsFobs16.push_back(-165.768);
0247   parsFobs16.push_back(0.0021543);
0248   parsFobs16.push_back(0.0148839);
0249   parsFobs16.push_back(4391.8);
0250   myLRhelper->setObsFitParameters(16, parsFobs16);
0251 
0252   std::vector<double> parsFobs17;
0253   parsFobs17.push_back(0.45862);
0254   parsFobs17.push_back(-42.3119);
0255   parsFobs17.push_back(0.0024431);
0256   parsFobs17.push_back(-0.0082168);
0257   parsFobs17.push_back(-41.3239);
0258   myLRhelper->setObsFitParameters(17, parsFobs17);
0259 
0260   std::vector<double> parsFobs18;
0261   parsFobs18.push_back(0.57713);
0262   parsFobs18.push_back(-88.4547);
0263   parsFobs18.push_back(-0.0079014);
0264   parsFobs18.push_back(-0.025394);
0265   parsFobs18.push_back(4512.33);
0266   myLRhelper->setObsFitParameters(18, parsFobs18);
0267 
0268   std::vector<double> parsFobs33;
0269   parsFobs33.push_back(5.99882e-01);
0270   parsFobs33.push_back(-1.33575e+01);
0271   parsFobs33.push_back(1.24161e-01);
0272   myLRhelper->setObsFitParameters(33, parsFobs33);
0273 
0274   std::vector<double> parsFobs35;
0275   parsFobs35.push_back(2.49026e-01);
0276   parsFobs35.push_back(1.08819e+00);
0277   parsFobs35.push_back(-7.26373e-01);
0278   parsFobs35.push_back(1.26367e-07);
0279   parsFobs35.push_back(5.51754e+02);
0280   parsFobs35.push_back(3.94562e-02);
0281   myLRhelper->setObsFitParameters(35, parsFobs35);
0282 
0283   std::vector<double> parsFobs37;
0284   parsFobs37.push_back(1.43676e-01);
0285   parsFobs37.push_back(2.44475e+00);
0286   parsFobs37.push_back(-4.56374e+00);
0287   parsFobs37.push_back(3.01449e+00);
0288   parsFobs37.push_back(4.65671e+01);
0289   parsFobs37.push_back(-4.40296e-02);
0290   myLRhelper->setObsFitParameters(37, parsFobs37);
0291 
0292   // fill signal and background contributions to S and B histograms
0293   doEventloop();
0294 
0295   // normalize the S and B histograms to construct the pdf's
0296   //myLRhelper -> normalizeSandBhists();
0297   myLRhelper->normalizeSandBhists();
0298 
0299   // produce and fit the S/S+N histograms
0300   myLRhelper->makeAndFitSoverSplusBHists();
0301 
0302   // store histograms and fits in root-file
0303   myLRhelper->storeToROOTfile(SignalSelOutfileName);
0304 
0305   // make some control plots and put them in a .ps file
0306   myLRhelper->storeControlPlots(SignalSelPSfile);
0307 }
0308 
0309 //
0310 // Loop over the events (with the definition of what is considered signal and background)
0311 //
0312 
0313 void doEventloop() {
0314   std::cout << std::endl << std::endl << "**** STARTING EVENT LOOP FOR SIGNAL ****" << std::endl;
0315 
0316   /********************************************** for the signal **********************************************/
0317 
0318   int okEvents = 0, totNrEv = 0;
0319   for (int nrDir = 0; nrDir < signal_nrDir; nrDir++) {
0320     std::cout << " Signal : " << signal_path[nrDir] << std::endl;
0321 
0322     int Signal_totNrEv = 0, Signal_okEvents = 0;
0323     for (int nr = 1; nr <= signal_nrFiles[nrDir]; nr++) {
0324       TString signal_ft = signal_path[nrDir];
0325       signal_ft += nr;
0326       signal_ft += ".root";
0327       if (!gSystem->AccessPathName(signal_ft)) {
0328         TFile *signal_file = TFile::Open(signal_ft);
0329         TTree *signal_events = dynamic_cast<TTree *>(signal_file->Get("Events"));
0330         assert(signal_events != nullptr);
0331 
0332         TBranch *signal_solsbranch = signal_events->GetBranch("TtSemiEvtSolutions_solutions__TEST.obj");
0333         assert(signal_solsbranch != nullptr);
0334         std::vector<TtSemiEvtSolution> signal_sols;
0335         //signal_solsbranch->SetAddress( & signal_sols );
0336 
0337         //loop over all events in a file
0338         for (int ev = 0; ev < signal_events->GetEntries(); ++ev) {
0339           if (Signal_totNrEv > signal_NrEv[nrDir] && signal_NrEv[nrDir] != -1)
0340             continue;
0341           ++Signal_totNrEv;
0342           ++totNrEv;
0343           if ((double)((totNrEv * 1.) / 5000.) == (double)(totNrEv / 5000))
0344             std::cout << "  Processing signal event " << totNrEv << std::endl;
0345           signal_solsbranch->SetAddress(&signal_sols);
0346           signal_solsbranch->GetEntry(ev);
0347           signal_events->GetEntry(ev, 0);
0348           if (signal_sols.size() == 12) {
0349             // get observable values
0350             std::vector<double> signal_obsVals;
0351             for (int j = 0; j < nrSignalSelObs; j++) {
0352               if (myLRhelper->obsFitIncluded(obsNrs[j]))
0353                 signal_obsVals.push_back(signal_sols[0].getLRSignalEvtObsVal(obsNrs[j]));
0354             }
0355 
0356             // Fill the observables
0357             // signal: semileptonic top event
0358             myLRhelper->fillToSignalHists(signal_obsVals);
0359             ++Signal_okEvents;
0360             ++okEvents;
0361           }
0362         }
0363         signal_file->Close();
0364 
0365       } else {
0366         std::cout << signal_ft << " doesn't exist" << std::endl;
0367       }
0368     }
0369     std::cout << std::endl
0370               << "********************  STATISTICS FOR SIGNAL " << signal_path[nrDir] << " ***********************"
0371               << std::endl;
0372     std::cout << std::endl << " Nb of requested events  :" << (signal_NrEv[nrDir]) << std::endl;
0373     std::cout << std::endl << " Nb of processed events  :" << (Signal_totNrEv) << std::endl;
0374     std::cout << std::endl << " Nb of events filled in the histo :" << (Signal_okEvents) << std::endl;
0375     std::cout << std::endl << "******************************************************************" << std::endl;
0376   }
0377   std::cout << std::endl << "********************  STATISTICS FOR SIGNAL ***********************" << std::endl;
0378   std::cout << std::endl << " Nb of processed events  :" << (totNrEv) << std::endl;
0379   std::cout << std::endl << " Nb of events filled in the histo :" << (okEvents) << std::endl;
0380   std::cout << std::endl << "******************************************************************" << std::endl;
0381 
0382   std::cout << std::endl << std::endl << "**** STARTING EVENT LOOP FOR BCKGD ****" << std::endl;
0383 
0384   /********************************************** for the background **********************************************/
0385 
0386   okEvents = 0, totNrEv = 0;
0387   for (int nrDir = 0; nrDir < bckgd_nrDir; nrDir++) {
0388     std::cout << " Background : " << bckgd_path[nrDir] << std::endl;
0389 
0390     int Bckgd_totNrEv = 0, Bckgd_okEvents = 0;
0391     for (int nr = 1; nr <= bckgd_nrFiles[nrDir]; nr++) {
0392       TString bckgd_ft = bckgd_path[nrDir];
0393       bckgd_ft += nr;
0394       bckgd_ft += ".root";
0395       if (!gSystem->AccessPathName(bckgd_ft)) {
0396         TFile *bckgd_file = TFile::Open(bckgd_ft);
0397         TTree *bckgd_events = dynamic_cast<TTree *>(bckgd_file->Get("Events"));
0398         assert(bckgd_events != nullptr);
0399 
0400         TBranch *bckgd_solsbranch = bckgd_events->GetBranch("TtSemiEvtSolutions_solutions__TEST.obj");
0401         assert(bckgd_solsbranch != nullptr);
0402         std::vector<TtSemiEvtSolution> bckgd_sols;
0403         //bckgd_solsbranch->SetAddress( & bckgd_sols );
0404 
0405         //loop over all events in a file
0406         for (int ev = 0; ev < bckgd_events->GetEntries(); ++ev) {
0407           if (Bckgd_totNrEv > bckgd_NrEv[nrDir] && bckgd_NrEv[nrDir] != -1)
0408             continue;
0409           ++Bckgd_totNrEv;
0410           ++totNrEv;
0411           if ((double)((totNrEv * 1.) / 5000.) == (double)(totNrEv / 5000))
0412             std::cout << "  Processing bckgd event " << totNrEv << std::endl;
0413           bckgd_solsbranch->SetAddress(&bckgd_sols);
0414           bckgd_solsbranch->GetEntry(ev);
0415           bckgd_events->GetEntry(ev, 0);
0416           if (bckgd_sols.size() == 12) {
0417             // get observable values
0418             std::vector<double> bckgd_obsVals;
0419             for (int j = 0; j < nrSignalSelObs; j++) {
0420               if (myLRhelper->obsFitIncluded(obsNrs[j]))
0421                 bckgd_obsVals.push_back(bckgd_sols[0].getLRSignalEvtObsVal(obsNrs[j]));
0422             }
0423             // Fill the observables
0424             // bckgd: semileptonic top event
0425             myLRhelper->fillToBackgroundHists(bckgd_obsVals);
0426             ++okEvents;
0427             ++Bckgd_okEvents;
0428           }
0429         }
0430         bckgd_file->Close();
0431 
0432       } else {
0433         std::cout << bckgd_ft << " doesn't exist" << std::endl;
0434       }
0435     }
0436     std::cout << std::endl
0437               << "********************  STATISTICS FOR BCKGD " << bckgd_path[nrDir] << " ***********************"
0438               << std::endl;
0439     std::cout << std::endl << " Nb of requested events  :" << (bckgd_NrEv[nrDir]) << std::endl;
0440     std::cout << std::endl << " Nb of processed events  :" << (Bckgd_totNrEv) << std::endl;
0441     std::cout << std::endl << " Nb of events filled in the histo :" << (Bckgd_okEvents) << std::endl;
0442     std::cout << std::endl << "******************************************************************" << std::endl;
0443   }
0444   std::cout << std::endl << "********************  STATISTICS FOR BCKGD ***********************" << std::endl;
0445   std::cout << std::endl << " Nb of processed events  :" << (totNrEv) << std::endl;
0446   std::cout << std::endl << " Nb of events filled in the histo :" << (okEvents) << std::endl;
0447   std::cout << std::endl << "******************************************************************" << std::endl;
0448 }