Egamma/macros/DoOptimizationNewVars.C

0001 // Revision 1.3  2008/01/21 02:25:22  wittich
0002 // cleanup, tweak of case where eff is the same
0003 //
0004 #ifndef __CINT__
0005 #include <map>
0006 #include <iostream>
0007 #include <vector>
0008 #include <fstream>
0009 using std::cout;
0010 using std::endl;
0011 using std::vector;
0012 #endif // __CINT__
0013
0014 #include "TFile.h"
0015 #include "TChain.h"
0016 #include "TH1F.h"
0017 #include "TH2F.h"
0018 #include "TCanvas.h"
0019 #include "TStyle.h"
0020 #include "TRegexp.h"
0021 #include "TGraphErrors.h"
0022 #include "TTreeFormula.h"
0023 #include "TROOT.h"
0024
0025
0026
0027 struct sample {
0028   TString name;
0029   TChain *chain;
0030   Double_t nEvents;
0031   Double_t xSection;
0032   Bool_t isSignal;
0033   Bool_t doPresel;
0034   Int_t nElecs;
0035   Double_t etPresel;
0036   Double_t effPresel;
0037 };
0038
0039 struct path {
0040   TString name;
0041   Int_t nCandsCut;
0042 };
0043
0044 void DoOptimization() {
0045   /* Function: DoOptimization() Performs the optimization algorithm as
0046      described in Jan. 10 TSG meeting.  Luminosity can be edited at
0047      the top.  The optimal cuts are written to myfile.
0048   */
0049
0050   /* Set Luminosity to use.  Conversion factor to get rate in Hz. */
0051   Double_t luminosity = 1.0E32; // in cm^-2 s^-1
0052   Double_t conversion = 1.0E-27;
0053
0054   /* Set output file name.  This file will contain a list of the
0055      optimal cuts on the various variables. */
0056   ofstream myfile;
0057   myfile.open("OptimalCutsNewVars.txt");
0058
0059   /* Set initial cuts on the variables */
0060 //   Double_t cutEt = 17.;
0061 //   Double_t cutHoE = 0.07;
0062 //   Double_t cutEpMatchBarrel = 0.9;
0063 //   Double_t cutEpMatchEndcap = 0.9;
0064 //   Double_t cutItrack = 0.15;
0065 //   Double_t cutEt = 20.;
0066 //   Double_t cutHoE = 0.125;
0067 //   Double_t cutEpMatchBarrel = 1.5;
0068 //   Double_t cutEpMatchEndcap = 2.5;
0069 //   Double_t cutItrack = 0.3;
0070   Double_t cutEt = 17.;
0071   Double_t cutHoE = 0.125;
0072   Double_t cutEpMatchBarrel = 2.1;
0073   Double_t cutEpMatchEndcap = 3.1;
0074   Double_t cutItrack = 0.1;
0075
0076   //17, 0.125, 2.1. 3.1, 0.1
0077
0078   /* IMPORTANT: the size of this array must be 3^N, where N is the
0079      number of cuts that are varying.  This should way of doing things
0080      should be changed.  However, it is relatively complicated an I
0081      have not had a chance yet.  Keep in mind that if you change N,
0082      you MUST change many other things in the code. */
0083   const unsigned int kNDim = 243;
0084   Double_t newEff[kNDim];
0085
0086   /* Set threshold on rate. */
0087   Double_t rateLimit = 20.;
0088
0089   /* Set the initial step sizes.  These sizes will be divided by 2
0090      every time an optimal point is reached until they reach some
0091      threshold on the first step (also below) */
0092   Double_t stepSizes[5];
0093   stepSizes[0] = 1.;
0094   stepSizes[1] = 0.02;
0095   stepSizes[2] = 0.2;
0096   stepSizes[3] = 0.2;
0097   stepSizes[4] = 0.2;
0098   Double_t stepThresh = 0.25;
0099
0100   /* Set the HLT paths to study.  Choices are (Relaxed)SingleElecs,
0101      (Relaxed)DoubleElecs.  The only difference between these is what
0102      L1 was used and whether L1 Non-Iso electrons are included For
0103      details on what cuts were used, see
0104      http://cmslxr.fnal.gov/lxr/source/HLTrigger/Egamma/data/EgammaHLTLocal_1e32.cff?v=CMSSW_1_6_7
0105   */
0106   vector<path> paths;
0107   path thisPath;
0108   thisPath.name = "SingleElecs";
0109   thisPath.nCandsCut = 1; // Set to 2 for double paths
0110   paths.push_back(thisPath);
0111
0112   /* More parameters (including samples) set below */
0113
0114   /* Variable definitions */
0115
0116   TChain *thisChain;
0117
0118   vector<sample> samples;
0119   sample thisSample;
0120
0121   Double_t pass = 0.;
0122   //  Double_t eff = 0.;
0123   Double_t bestEff = 0., bestEffRate = 0.;
0124   //  Double_t rate = 0., rateL1 = 0., rateTot = 0.;
0125   Double_t rate = 0.;
0126   Double_t thisCutEt = cutEt;
0127   Double_t thisCutHoE = cutHoE;
0128   Double_t thisCutEpMatchBarrel = cutEpMatchBarrel;
0129   Double_t thisCutEpMatchEndcap = cutEpMatchEndcap;
0130   Double_t thisCutItrack = cutItrack;
0131
0132   unsigned int pathNum, sampleNum;
0133
0134   unsigned int direction = 0, bestDir = 0;
0135
0136   TString cutText = "";
0137
0138   Bool_t isBetter = false;
0139
0140   //Int_t xSecNum = 0, fileNum = 0;
0141
0142   /* Set samples as follows:
0143      thisChain = new TChain("Events");
0144      thisChain->Add("fileN.root");
0145      where fileN.root is one of the root files contributing to the samples
0146      ...
0147      thisSample.name = name for display purposes if necessary;
0148      thisSample.chain = thisChain; always necessary
0149      thisSample.nEvents = number of events in sample BEFORE ANY CUTS; used to calculate rate
0150      thisSample.xSection = cross section for this sample in mb;
0151      thisSample.isSignal = set to true for signal sample(s);
0152      thisSample.doPresel = set to true to do MC preselection on
0153                            electrons (only used to mimic previous studies);
0154      thisSample.nElecs = number of electrons in MC sample (1 for W->ev,
0155                          2 for Z->ee, etc.) for preselection;
0156      thisSample.etPresel = MC Et cut for preselection;
0157      thisSample.effPresel = Currently, a pre-calculated efficiency for this
0158                             cut is needed to do preselection.  See Lorenzo's note
0159                 for theoretical values;
0160      samples.push_back(thisSample); always needed at the end.
0161   */
0162
0163   thisChain = new TChain("Events");
0164   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-0.root");
0165   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-1.root");
0166   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-2.root");
0167   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-3.root");
0168   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-4.root");
0169   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-5.root");
0170   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-6.root");
0171   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-7.root");
0172   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-8.root");
0173   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-9.root");
0174   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-10.root");
0175   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-11.root");
0176   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-12.root");
0177   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-13.root");
0178   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-14.root");
0179   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-15.root");
0180   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-16.root");
0181   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-17.root");
0182   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-18.root");
0183   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-19.root");
0184   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-20.root");
0185   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-21.root");
0186   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-22.root");
0187   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-23.root");
0188   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-24.root");
0189   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-25.root");
0190   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-26.root");
0191   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-27.root");
0192   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-28.root");
0193   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-29.root");
0194   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-30.root");
0195   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-31.root");
0196   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-32.root");
0197   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-33.root");
0198   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-34.root");
0199   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-35.root");
0200   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-36.root");
0201   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-37.root");
0202   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-38.root");
0203   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-39.root");
0204   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-40.root");
0205   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-41.root");
0206   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-42.root");
0207   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-43.root");
0208   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-44.root");
0209   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-45.root");
0210   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-46.root");
0211   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-47.root");
0212   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-48.root");
0213   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-49.root");
0214   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-50.root");
0215   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-51.root");
0216   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-52.root");
0217   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-53.root");
0218   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-54.root");
0219   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-55.root");
0220   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-56.root");
0221   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-57.root");
0222   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-58.root");
0223   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-59.root");
0224   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-60.root");
0225   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-61.root");
0226   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-62.root");
0227   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-63.root");
0228   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-64.root");
0229   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-65.root");
0230   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-66.root");
0231   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-67.root");
0232   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-68.root");
0233   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-69.root");
0234   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-70.root");
0235   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-71.root");
0236   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-72.root");
0237   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-73.root");
0238   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-74.root");
0239   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-75.root");
0240   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-76.root");
0241   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-77.root");
0242   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-78.root");
0243   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-79.root");
0244   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-80.root");
0245   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-81.root");
0246   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-82.root");
0247   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-83.root");
0248   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-84.root");
0249   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-85.root");
0250   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-86.root");
0251   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-87.root");
0252   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-88.root");
0253   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-89.root");
0254   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-90.root");
0255   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-91.root");
0256   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-92.root");
0257   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-93.root");
0258   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-94.root");
0259   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-95.root");
0260   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-96.root");
0261   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-97.root");
0262   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-98.root");
0263   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-99.root");
0264   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-100.root");
0265   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-101.root");
0266   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-102.root");
0267   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-103.root");
0268   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-104.root");
0269   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-105.root");
0270   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-106.root");
0271   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-107.root");
0272   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-108.root");
0273   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-109.root");
0274   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-110.root");
0275   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-111.root");
0276   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-112.root");
0277   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-113.root");
0278   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-114.root");
0279   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-115.root");
0280   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-116.root");
0281   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-117.root");
0282   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-118.root");
0283   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-119.root");
0284   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-120.root");
0285   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-121.root");
0286   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-122.root");
0287   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-123.root");
0288   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-124.root");
0289   thisChain->Add("~/data/hlt/egamma/1_6_7/QCD-15-20/QCD-15-20-HLTVars-125.root");
0290   thisSample.name = "QCD 15-20";
0291   thisSample.chain = thisChain;
0292   thisSample.nEvents = 1260000;
0293   thisSample.xSection = 1.46;
0294   thisSample.isSignal = false;
0295   thisSample.doPresel = false;
0296   thisSample.nElecs = 0;
0297   thisSample.etPresel = 0.;
0298   thisSample.effPresel = 0.;
0299   samples.push_back(thisSample);
0300
0301   thisChain = new TChain("Events");
0302   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-20-30-HLTVars.root");
0303   thisSample.name = "QCD 20-30";
0304   thisSample.chain = thisChain;
0305   thisSample.nEvents = 100000;
0306   thisSample.xSection = 6.32E-1;
0307   thisSample.isSignal = false;
0308   thisSample.doPresel = false;
0309   thisSample.nElecs = 0;
0310   thisSample.etPresel = 0.;
0311   thisSample.effPresel = 0.;
0312   samples.push_back(thisSample);
0313
0314   thisChain = new TChain("Events");
0315   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-30-50-HLTVars.root");
0316   thisSample.name = "QCD 30-50";
0317   thisSample.chain = thisChain;
0318   thisSample.nEvents = 100000;
0319   thisSample.xSection = 1.63E-1;
0320   thisSample.isSignal = false;
0321   thisSample.doPresel = false;
0322   thisSample.nElecs = 0;
0323   thisSample.etPresel = 0.;
0324   thisSample.effPresel = 0.;
0325   samples.push_back(thisSample);
0326
0327   thisChain = new TChain("Events");
0328   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-1.root");
0329   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-2.root");
0330   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-3.root");
0331   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-4.root");
0332   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-5.root");
0333   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-6.root");
0334   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-7.root");
0335   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-8.root");
0336   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-9.root");
0337   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars-10.root");
0338   thisSample.name = "QCD 50-80";
0339   thisSample.chain = thisChain;
0340   thisSample.nEvents = 100000;
0341   thisSample.xSection = 2.16E-2;
0342   thisSample.isSignal = false;
0343   thisSample.doPresel = false;
0344   thisSample.nElecs = 0;
0345   thisSample.etPresel = 0.;
0346   thisSample.effPresel = 0.;
0347   samples.push_back(thisSample);
0348
0349 //   thisChain = new TChain("Events");;
0350 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-0.root");
0351 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-1.root");
0352 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-2.root");
0353 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-3.root");
0354 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-4.root");
0355 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-5.root");
0356 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-6.root");
0357 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-7.root");
0358 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-8.root");
0359 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-HLTVars80120-9.root");
0360 //   thisSample.name = "QCD 80-120";
0361 //   thisSample.chain = thisChain;
0362 //   thisSample.nEvents = 239993;
0363 //   thisSample.xSection = 3.08E-3;
0364 //   thisSample.isSignal = false;
0365 //   thisSample.doPresel = false;
0366 //   thisSample.nElecs = 0;
0367 //   thisSample.etPresel = 0.;
0368 //   thisSample.effPresel = 0.;
0369 //   samples.push_back(thisSample);
0370
0371 //   thisChain = new TChain("Events");;
0372 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-120-170-HLTVars-2.root");
0373 //   thisSample.name = "QCD 120-170";
0374 //   thisSample.chain = thisChain;
0375 //   thisSample.nEvents = 11092;
0376 //   thisSample.xSection = 4.94E-4;
0377 //   thisSample.isSignal = false;
0378 //   thisSample.doPresel = false;
0379 //   thisSample.nElecs = 0;
0380 //   thisSample.etPresel = 0.;
0381 //   thisSample.effPresel = 0.;
0382 //   samples.push_back(thisSample);
0383
0384 //   thisChain = new TChain("Events");;
0385 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-170-230-HLTVars-0.root");
0386 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-170-230-HLTVars-1.root");
0387 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-170-230-HLTVars-2.root");
0388 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-170-230-HLTVars-3.root");
0389 //   thisSample.name = "QCD 170-230";
0390 //   thisSample.chain = thisChain;
0391 //   thisSample.nEvents = 400000;
0392 //   thisSample.xSection = 1.01E-4;
0393 //   thisSample.isSignal = false;
0394 //   thisSample.doPresel = false;
0395 //   thisSample.nElecs = 0;
0396 //   thisSample.etPresel = 0.;
0397 //   thisSample.effPresel = 0.;
0398 //   samples.push_back(thisSample);
0399
0400 //   thisChain = new TChain("Events");;
0401 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-230-300-HLTVars-0.root");
0402 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-230-300-HLTVars-1.root");
0403 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-230-300-HLTVars-2.root");
0404 //   thisChain->Add("~/data/hlt/egamma/1_6_0/QCD-230-300-HLTVars-3.root");
0405 //   thisSample.name = "QCD 230-300";
0406 //   thisSample.chain = thisChain;
0407 //   thisSample.nEvents = 400000;
0408 //   thisSample.xSection = 2.45E-5;
0409 //   thisSample.isSignal = false;
0410 //   thisSample.doPresel = false;
0411 //   thisSample.nElecs = 0;
0412 //   thisSample.etPresel = 0.;
0413 //   thisSample.effPresel = 0.;
0414 //   samples.push_back(thisSample);
0415
0416   thisChain = new TChain("Events");;
0417   thisChain->Add("../test/HLTStudyData/ZEE-HLTVars-NoPresel.root");
0418   thisSample.name = "Z->ee";
0419   thisSample.chain = thisChain;
0420   thisSample.nEvents = 3800;
0421   thisSample.xSection = 1.62E-6;
0422   thisSample.isSignal = false;
0423   thisSample.doPresel = true;
0424   thisSample.nElecs = 2;
0425   thisSample.etPresel = 5.;
0426   thisSample.effPresel = 0.504;
0427   samples.push_back(thisSample);
0428
0429   thisChain = new TChain("Events");
0430   thisChain->Add("../test/HLTStudyData/WENU-HLTVars-NoPresel.root");
0431   thisSample.name = "W->ev";
0432   thisSample.chain = thisChain;
0433   thisSample.nEvents = 2000;
0434   thisSample.xSection = 1.72E-5;
0435   thisSample.isSignal = true;
0436   thisSample.doPresel = true;
0437   thisSample.nElecs = 1;
0438   thisSample.etPresel = 20.;
0439   thisSample.effPresel = 0.463;
0440   samples.push_back(thisSample);
0441
0442   thisChain = new TChain("Events");;
0443   thisChain->Add("~/data/hlt/egamma/1_6_0/TTbar-HLTVars-1e.root");
0444   thisSample.name = "t-tbar";
0445   thisSample.chain = thisChain;
0446   thisSample.nEvents = 100000;
0447   thisSample.xSection = 8.33E-7;
0448   thisSample.isSignal = false;
0449   thisSample.doPresel = false;
0450   thisSample.nElecs = 1;
0451   thisSample.etPresel = 0.;
0452   thisSample.effPresel = 0.;
0453   samples.push_back(thisSample);
0454
0455   /* Now do the initial cut to set a base rate and efficiency */
0456
0457   /* For each path, look at each sample and apply the cuts to get
0458      the efficiency and rate.  DOES NOT WORK WITH MULTIPLE PATHS
0459      AT ONCE YET */
0460   for (pathNum = 0; pathNum < paths.size(); pathNum++) {
0461     rate = 0.;
0462     for (sampleNum = 0; sampleNum < samples.size(); sampleNum++) {
0463       /* Set the cut... Messy for now until update to add H/E, |1/E-1/p| */
0464       cutText = "Sum$(";
0465       cutText += paths[pathNum].name;
0466       cutText += ".l1Match && ";
0467       cutText += paths[pathNum].name;
0468       cutText += ".Et > ";
0469       cutText += thisCutEt;
0470       cutText += " && ";
0471
0472       /*
0473       cutText += paths[pathNum].name;
0474       cutText += ".HoE < ";
0475       cutText += thisCutHoE;
0476       To be uncommented when all samples have HOE */
0477
0478       cutText += paths[pathNum].name;
0479       cutText += ".IHcal / ";
0480       cutText += paths[pathNum].name;
0481       cutText += ".Et < ";
0482       cutText += thisCutHoE;
0483       cutText += " && ";
0484       cutText += paths[pathNum].name;
0485
0486       /* From here, we must do EpMatch in barrel and endcap */
0487       /*
0488       cutText += ".pixMatch > 0 && ((";
0489       cutText += paths[pathNum].name;
0490       cutText += ".EpMatch < ";
0491       cutText += thisCutEpMatchBarrel;
0492       cutText += " && fabs(";
0493       cutText += paths[pathNum].name;
0494       cutText += ".eta) < 1.5) || (";
0495       cutText += paths[pathNum].name;
0496       cutText += ".EpMatch < ";
0497       cutText += thisCutEpMatchEndcap;
0498       cutText += " && fabs(";
0499       cutText += paths[pathNum].name;
0500       cutText += ".eta) > 1.5 && fabs(";
0501       cutText += paths[pathNum].name;
0502       cutText += ".eta) < 2.5)) && ";
0503       To be uncommented when all samples have EpMatch */
0504
0505       cutText += ".pixMatch > 0 && ((fabs((1. - ";
0506       cutText += paths[pathNum].name;
0507       cutText += ".Eoverp) * sin(2.*fabs(exp(-";
0508       cutText += paths[pathNum].name;
0509       cutText += ".eta))) / ";
0510       cutText += paths[pathNum].name;
0511       cutText += ".Et) < ";
0512       cutText += thisCutEpMatchBarrel;
0513       cutText += " && fabs(";
0514       cutText += paths[pathNum].name;
0515       cutText += ".eta) < 1.5) || (fabs((1. - ";
0516       cutText += paths[pathNum].name;
0517       cutText += ".Eoverp) * sin(2.*fabs(exp(-";
0518       cutText += paths[pathNum].name;
0519       cutText += ".eta))) / ";
0520       cutText += paths[pathNum].name;
0521       cutText += ".Et) < ";
0522       cutText += thisCutEpMatchEndcap;
0523       cutText += " && fabs(";
0524       cutText += paths[pathNum].name;
0525       cutText += ".eta) > 1.5 && fabs(";
0526       cutText += paths[pathNum].name;
0527       cutText += ".eta) < 2.5)) && ";
0528       cutText += paths[pathNum].name;
0529       cutText += ".Itrack < ";
0530       cutText += thisCutItrack;
0531       cutText += ") >= ";
0532       cutText += paths[pathNum].nCandsCut;
0533       if (samples[sampleNum].doPresel) {
0534     cutText += " && Sum$(CaloVarss_hltCutVars_mcSingleElecs_EGAMMAHLT.obj.Et > ";
0535     cutText += samples[pathNum].etPresel;
0536     cutText += " && fabs(CaloVarss_hltCutVars_mcSingleElecs_EGAMMAHLT.obj.eta) < 2.7) >= ";
0537     cutText += samples[sampleNum].nElecs;
0538       }
0539
0540       pass = (Double_t)samples[sampleNum].chain->Draw("", cutText);
0541
0542       rate += pass / samples[sampleNum].nEvents * samples[sampleNum].xSection
0543     * luminosity * conversion;
0544
0545       if (samples[sampleNum].isSignal) {
0546     bestEff = pass;
0547     bestEffRate = rate;
0548       }
0549     }
0550     if (rate > rateLimit ) {
0551       std::cout<<"Invalid initial choice!"<< std::endl;
0552       std::cout<<"Rate is "<<rate<< std::endl;
0553       return;
0554     }
0555     int cnt = 0; // how many steps
0556     while (stepSizes[0] > stepThresh) {
0557       isBetter = false;
0558       std::cout << cutEt <<  ", "
0559        << cutHoE << ", "
0560        << cutEpMatchBarrel << ". "
0561        << cutEpMatchEndcap << ", "
0562        << cutItrack
0563        << endl;
0564       for (direction = 0; direction < kNDim; direction++) {
0565     rate = 0.;
0566     if (direction % 3 == 0) thisCutEt = cutEt + stepSizes[0];
0567     else if  (direction % 3 == 1) thisCutEt = cutEt;
0568     else thisCutEt = cutEt - stepSizes[0];
0569     if ((direction / 3) % 3 == 0) thisCutHoE = cutHoE + stepSizes[1];
0570     else if ((direction / 3) % 3 == 1) thisCutHoE = cutHoE;
0571     else thisCutHoE = cutHoE - stepSizes[1];
0572     if ((direction / 9) % 3 == 0) thisCutEpMatchBarrel = cutEpMatchBarrel + stepSizes[2];
0573     else if ((direction / 9) % 3 == 1) thisCutEpMatchBarrel = cutEpMatchBarrel;
0574     else thisCutEpMatchBarrel = cutEpMatchBarrel - stepSizes[2];
0575     if ((direction / 27) % 3 == 0) thisCutEpMatchEndcap = cutEpMatchEndcap + stepSizes[3];
0576     else if ((direction / 27) % 3 == 1) thisCutEpMatchEndcap = cutEpMatchEndcap;
0577     else thisCutEpMatchEndcap = cutEpMatchEndcap - stepSizes[3];
0578     if ((direction / 81) % 3 == 0)
0579       thisCutItrack = cutItrack + stepSizes[4];
0580     else
0581       if ((direction / 81) % 3 == 1)
0582         thisCutItrack = cutItrack;
0583     else
0584       thisCutItrack = std::max(0.05, cutItrack - stepSizes[4]);
0585     std::cout << thisCutEt << ", "
0586          <<thisCutHoE << ", "
0587          <<thisCutEpMatchBarrel << ", "
0588          <<thisCutEpMatchEndcap << ", "
0589          <<thisCutItrack
0590          << std::endl;
0591     for (sampleNum = 0; sampleNum < samples.size(); sampleNum++) {
0592       cutText = "Sum$(";
0593       cutText += paths[pathNum].name;
0594       cutText += ".l1Match && ";
0595       cutText += paths[pathNum].name;
0596       cutText += ".Et > ";
0597       cutText += thisCutEt;
0598       cutText += " && ";
0599
0600       /*
0601       cutText += paths[pathNum].name;
0602       cutText += ".HoE < ";
0603           cutText += thisCutHoE;
0604       To be uncommented when all samples have HOE */
0605
0606       cutText += paths[pathNum].name;
0607       cutText += ".IHcal / ";
0608       cutText += paths[pathNum].name;
0609       cutText += ".Et < ";
0610       cutText += thisCutHoE;
0611       cutText += " && ";
0612       cutText += paths[pathNum].name;
0613
0614       /* From here, we must do EpMatch in barrel and endcap */
0615       /*
0616       cutText += ".pixMatch > 0 && ((";
0617       cutText += paths[pathNum].name;
0618       cutText += ".EpMatch < ";
0619       cutText += thisCutEpMatchBarrel;
0620       cutText += " && fabs(";
0621       cutText += paths[pathNum].name;
0622       cutText += ".eta) < 1.5) || (";
0623       cutText += paths[pathNum].name;
0624       cutText += ".EpMatch < ";
0625       cutText += thisCutEpMatchEndcap;
0626       cutText += " && fabs(";
0627       cutText += paths[pathNum].name;
0628       cutText += ".eta) > 1.5 && fabs(";
0629       cutText += paths[pathNum].name;
0630       cutText += ".eta) < 2.5)) && ";
0631       To be uncommented when all samples have EpMatch */
0632
0633       cutText += ".pixMatch > 0 && ((fabs((1. - ";
0634       cutText += paths[pathNum].name;
0635       cutText += ".Eoverp) * sin(2.*fabs(exp(-";
0636       cutText += paths[pathNum].name;
0637       cutText += ".eta))) / ";
0638       cutText += paths[pathNum].name;
0639       cutText += ".Et) < ";
0640       cutText += thisCutEpMatchBarrel;
0641       cutText += " && fabs(";
0642       cutText += paths[pathNum].name;
0643       cutText += ".eta) < 1.5) || (fabs((1. - ";
0644       cutText += paths[pathNum].name;
0645       cutText += ".Eoverp) * sin(2.*fabs(exp(-";
0646       cutText += paths[pathNum].name;
0647       cutText += ".eta))) / ";
0648       cutText += paths[pathNum].name;
0649       cutText += ".Et) < ";
0650       cutText += thisCutEpMatchEndcap;
0651       cutText += " && fabs(";
0652       cutText += paths[pathNum].name;
0653       cutText += ".eta) > 1.5 && fabs(";
0654       cutText += paths[pathNum].name;
0655       cutText += ".eta) < 2.5)) && ";
0656       cutText += paths[pathNum].name;
0657       cutText += ".Itrack < ";
0658       cutText += thisCutItrack;
0659       cutText += ") >= ";
0660       cutText += paths[pathNum].nCandsCut;
0661       if (samples[sampleNum].doPresel) {
0662         cutText += " && Sum$(CaloVarss_hltCutVars_mcSingleElecs_EGAMMAHLT.obj.Et > ";
0663         cutText += samples[pathNum].etPresel;
0664         cutText += " && fabs(CaloVarss_hltCutVars_mcSingleElecs_EGAMMAHLT.obj.eta) < 2.7) >= ";
0665         cutText += samples[sampleNum].nElecs;
0666       }
0667
0668       pass = (Double_t)samples[sampleNum].chain->Draw("", cutText);
0669
0670       rate += pass / samples[sampleNum].nEvents * samples[sampleNum].xSection * luminosity * conversion;
0671
0672       if (samples[sampleNum].isSignal) {
0673         newEff[direction] = pass;
0674       }
0675     }
0676     // To keep track of progress during output
0677     std::cout <<"Pass " << cnt << ": "<<newEff[direction]<<", Rate: "<<rate
0678           << ", direction: " << direction
0679           << std::endl;
0680
0681     if (rate <= rateLimit ) {
0682       if ( (newEff[direction] > bestEff)
0683            || ((newEff[direction] == bestEff) && (rate < bestEffRate))) {
0684         bestEff = newEff[direction];
0685         bestDir = direction;
0686         bestEffRate = rate;
0687         isBetter = true;
0688       }
0689     }
0690       }
0691       /* THIS PART SHOULD BE CHANGED.  Right now, if you change the number
0692      of cuts varying, you must edit this part non-trivially.  The Nth
0693      cut should have:
0694      if ((bestDir / 3^(N-1)) % 3 == 0) cut = cut + stepSizes[N-1];
0695      else if ((bestDir / 3^(N-1)) % 3 == 1) cut = cut;
0696      else cut = cut - stepSizes[N-1];
0697       */
0698       if (isBetter) {
0699     std::cout << "There is a better cut!"<< std::endl;
0700
0701     if (bestDir % 3 == 0)
0702       cutEt = cutEt + stepSizes[0];
0703     else if  (bestDir % 3 == 1)
0704       cutEt = cutEt;
0705     else
0706       cutEt = cutEt - stepSizes[0];
0707
0708     if ((bestDir / 3) % 3 == 0)
0709       cutHoE = cutHoE + stepSizes[1];
0710     else if ((bestDir / 3) % 3 == 1)
0711       cutHoE = cutHoE;
0712     else
0713       cutHoE = cutHoE - stepSizes[1];
0714
0715     if ((bestDir / 9) % 3 == 0)
0716       cutEpMatchBarrel = cutEpMatchBarrel + stepSizes[2];
0717     else if ((bestDir / 9) % 3 == 1)
0718       cutEpMatchBarrel = cutEpMatchBarrel;
0719     else
0720       cutEpMatchBarrel = cutEpMatchBarrel - stepSizes[2];
0721
0722     if ((bestDir / 27) % 3 == 0)
0723       cutEpMatchEndcap = cutEpMatchEndcap + stepSizes[3];
0724     else if ((bestDir / 27) % 3 == 1)
0725       cutEpMatchEndcap = cutEpMatchEndcap;
0726     else
0727       cutEpMatchEndcap = cutEpMatchEndcap - stepSizes[3];
0728
0729     if ((bestDir / 81) % 3 == 0)
0730       cutItrack = cutItrack + stepSizes[4];
0731     else if ((bestDir / 81) % 3 == 1)
0732       cutItrack = cutItrack;
0733     else
0734       cutItrack = cutItrack - stepSizes[4];
0735       }
0736       else {
0737     std::cout << "Nothing better... Shrinking window"<< std::endl;
0738     stepSizes[0] = stepSizes[0] / 2.;
0739     stepSizes[1] = stepSizes[1] / 2.;
0740     stepSizes[2] = stepSizes[2] / 2.;
0741     stepSizes[3] = stepSizes[3] / 2.;
0742     stepSizes[4] = stepSizes[4] / 2.;
0743       }
0744       std::cout << " CURRENT BEST below (direction = "
0745         << bestDir << ", best eff = " << bestEff
0746         << ", rate = " << bestEffRate
0747         << ")"
0748         << std::endl;
0749       std::cout <<"Optimal cuts: "<< endl;
0750       std::cout <<"Et = "<< cutEt << endl;
0751       std::cout <<"HoE = "<< cutHoE<< endl;
0752       std::cout <<"EpMatchBarrel = "<< cutEpMatchBarrel<< endl;
0753       std::cout <<"EpMatchEndcap = "<< cutEpMatchEndcap<< endl;
0754       std::cout <<"Itrack = "<< cutItrack << endl;
0755       ++cnt;
0756     }
0757
0758     myfile << "Optimal cuts: " << std::endl;
0759     myfile << "eff = " << bestEff << std::endl;
0760     myfile << "Et = " << cutEt << std::endl;
0761     myfile << "HoE = " << cutHoE << std::endl;
0762     myfile << "EpMatchBarrel = " << cutEpMatchBarrel << std::endl;
0763     myfile << "EpMatchEndcap = " << cutEpMatchEndcap << std::endl;
0764     myfile << "Itrack = " << cutItrack << std::endl;
0765
0766   }
0767
0768   myfile.close();
0769
0770 }
0771
0772 #ifndef __CINT__
0773 TROOT theroot("theroot", "of all evil");
0774
0775 int main()
0776 {
0777   DoOptimization();
0778   return 0;
0779 }
0780 #endif // __CINT__