Add benchmark ex8.py

TIMBER/Analyzer.py

@@ -459,7 +459,8 @@ def MergeCollections(self,name,collectionNames):
                         concat_str = 'Concatenate(%s,%s)'%(concat_str,collName+'_'+var)
 
                 self.Define(name+'_'+var,concat_str)
-
+
+        self.Define('n'+name,'+'.join(['n'+n for n in collectionNames]))
 
         # for i in range(1,len(collectionNames)):
         #     collName = collectionNames[i]

TIMBER/Framework/include/common.h

@@ -28,6 +28,19 @@ namespace analyzer {
         return v;
     }
 
+    RVec<ROOT::Math::PtEtaPhiMVector> TLvector(RVec<float> pt,RVec<float> eta,RVec<float> phi,RVec<float> m) {
+        RVec<ROOT::Math::PtEtaPhiMVector> vs;
+        for (int i = 0; i < pt.size(); i++) {
+            ROOT::Math::PtEtaPhiMVector v(pt[i],eta[i],phi[i],m[i]);
+            vs.push_back(v);
+        }
+        return vs;
+    }
+
+    float transverseMass(float MET_pt, float obj_pt, float MET_phi, float obj_phi) {
+        return sqrt(2.0*MET_pt*obj_pt-(1-cos(deltaPhi(MET_phi,obj_phi))));
+    }
+
     double invariantMass(ROOT::Math::PtEtaPhiMVector v1, ROOT::Math::PtEtaPhiMVector v2) {
         return (v1+v2).M();
     }

benchmark/ex.cc

@@ -44,4 +44,21 @@ RVec<float> CloseLepVeto (RVec<float> Lepton_pt, RVec<float> Lepton_eta, RVec<fl
         jet_bool_vect.push_back(found_lep);
     }
     return jet_bool_vect;
+}
+
+int NonZlep(RVec<ROOT::Math::PtEtaPhiMVector> Lepton_vect, RVec<int> Lepton_pdgId, RVec<int> Lepton_charge) {
+    RVec<RVec<int>> combos = Combinations(Lepton_vect,3); // build combinations where first two are the Z
+    int NonZlep_idx = -1;
+    float deltaMZ = 1000.; // start at large value for comparison
+    float dMZ;
+    bool sameFlavOppSign;
+    for (int i = 0; i < combos[0].size(); i++) { // loop over combinations
+        dMZ = abs(91.2-(Lepton_vect[combos[0][i]]+Lepton_vect[combos[1][i]]).M());
+        sameFlavOppSign = (Lepton_pdgId[combos[0][i]] == -1*Lepton_pdgId[combos[1][i]]);
+        if ((dMZ < deltaMZ) && sameFlavOppSign) {
+            deltaMZ = dMZ;
+            NonZlep_idx = combos[2][i];
+        }
+    }
+    return NonZlep_idx;
 }

benchmark/ex7.py

@@ -12,9 +12,7 @@
 a.Define('goodJet_pt','Jet_pt[!(%s)]'%nearLep)
 a.Define('goodJet_pt30','goodJet_pt[goodJet_pt > 30]')
 a.Define('PtSum','Sum(goodJet_pt30)')
-print a.DataFrame.GetColumnType('PtSum')
 h = a.DataFrame.Histo1D('PtSum')
 h.Draw('hist e')
 #---------------------
-print ('%s secs'%(time.time() - start))
-input('')
+print ('%s secs'%(time.time() - start))

benchmark/ex8.py

@@ -0,0 +1,24 @@
+'''Plot the transverse mass of the Missing ET and a lepton for events with at least
+3 leptons and two of the leptons make up a Z (same flavor, opp sign, closest to invar
+mass of 91.2 GeV), where the lepton being plotted is the third lepton.'''
+# Lucas' Note: There could be four leptons in the event where two satisfy the Z condition
+# and then either of the remaining two could be the one used in the transverse mass
+# calculation. I do not consider that case (the first that comes up will be picked here).
+import time
+start = time.time()
+#---------------------
+from TIMBER.Analyzer import analyzer
+from TIMBER.Tools.Common import CompileCpp
+CompileCpp('benchmark/ex.cc')
+
+a = analyzer('examples/ttbar16_sample.root')
+a.MergeCollections("Lepton",["Electron","Muon"])
+a.Cut('nLepton2','nLepton>2')
+a.Define('Lepton_vect','analyzer::TLvector(Lepton_pt, Lepton_eta, Lepton_phi, Lepton_mass)')
+a.Define('NonZlep_idx','NonZlep(Lepton_vect,Lepton_pdgId,Lepton_charge)')
+a.Define('MT','NonZlep_idx == -1 ? -1 : analyzer::transverseMass(MET_pt, Lepton_pt[NonZlep_idx], MET_phi, Lepton_phi[NonZlep_idx])')
+a.Cut('MT_cut','MT>=0')
+h = a.DataFrame.Histo1D('MT')
+h.Draw('hist e')
+#---------------------
+print ('%s secs'%(time.time() - start))