Egamma/macros/DoOptimizationOldVars.C

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