Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:17:06

0001 from __future__ import print_function
0002 import ROOT
0003 from DataFormats.FWLite import Handle, Events
0004 import numpy as np
0005 import pickle
0006 
0007 events_c = Events('output_test_DDX.root')
0008 
0009 handleJ  = Handle("std::vector<pat::Jet>")
0010 labelJ = ("selectedUpdatedPatJets","","PATtest")
0011 
0012 h_probQ_ddb = ROOT.TH1F('h_probQ_ddb', ';prob Q;', 40, 0., 1.)
0013 h_probH_ddb = ROOT.TH1F('h_probH_ddb', ';prob H;', 40, 0., 1.)
0014 
0015 h_probQ_ddc = ROOT.TH1F('h_probQ_ddc', ';prob Q;', 40, 0., 1.)
0016 h_probH_ddc = ROOT.TH1F('h_probH_ddc', ';prob H;', 40, 0., 1.)
0017 
0018 info = {}
0019 info['pt'] = []
0020 info['eta'] = []
0021 info['mass'] = []
0022 info['BvL'] = []
0023 info['CvL'] = []
0024 info['CvB'] = []
0025 
0026 
0027 for iev, event in enumerate(events_c):
0028     event.getByLabel (labelJ, handleJ)
0029     jets = handleJ.product()
0030     print(iev)
0031     for jet in jets:
0032         #if jet.pt() < 300 or jet.pt() > 2000: continue
0033         #if jet.mass() < 40 or jet.mass() > 200: continue
0034 
0035         print(jet.pt(), jet.mass(), jet.eta())
0036         print("DDB", jet.bDiscriminator("pfDeepDoubleBvLJetTags:probQCD"), jet.bDiscriminator("pfDeepDoubleBvLJetTags:probHbb"))
0037         # print("DDB", jet.bDiscriminator("pfMassIndependentDeepDoubleBvLJetTags:probQCD"), jet.bDiscriminator("pfMassIndependentDeepDoubleBvLJetTags:probHbb"))
0038         print("DDCvL", jet.bDiscriminator("pfDeepDoubleCvLJetTags:probQCD"), jet.bDiscriminator("pfDeepDoubleCvLJetTags:probHcc"))
0039         # print("DDCvL", jet.bDiscriminator("pfMassIndependentDeepDoubleCvLJetTags:probQCD"), jet.bDiscriminator("pfMassIndependentDeepDoubleCvLJetTags:probHcc"))
0040         print("DDCvB", jet.bDiscriminator("pfDeepDoubleCvBJetTags:probHbb"), jet.bDiscriminator("pfDeepDoubleCvBJetTags:probHcc") )
0041         # print("DDCvB", jet.bDiscriminator("pfMassIndependentDeepDoubleCvBJetTags:probHbb"), jet.bDiscriminator("pfMassIndependentDeepDoubleCvBJetTags:probHcc"))
0042         h_probQ_ddb.Fill(jet.bDiscriminator("pfDeepDoubleBvLJetTags:probQCD"))
0043         h_probH_ddb.Fill(jet.bDiscriminator("pfDeepDoubleBvLJetTags:probHbb"))
0044         h_probQ_ddc.Fill(jet.bDiscriminator("pfDeepDoubleCvLJetTags:probQCD"))
0045         h_probH_ddc.Fill(jet.bDiscriminator("pfDeepDoubleCvLJetTags:probHcc"))
0046 
0047         info['mass'].append(jet.mass())
0048         info['pt'].append(jet.pt())
0049         info['eta'].append(jet.eta())
0050         info['BvL'].append(jet.bDiscriminator("pfDeepDoubleBvLJetTags:probHbb"))
0051         info['CvL'].append(jet.bDiscriminator("pfDeepDoubleCvLJetTags:probHcc"))
0052         info['CvB'].append(jet.bDiscriminator("pfDeepDoubleCvBJetTags:probHcc"))
0053 
0054 with open('outputs.pkl', 'wb') as handle:
0055     pickle.dump(info, handle)
0056 
0057 c1a = ROOT.TCanvas()
0058 h_probH_ddb.Draw("HISTO")
0059 h_probH_ddb.SetLineColor(632)
0060 h_probH_ddb.SetLineStyle(10)
0061 h_probQ_ddb.Draw("SAME")
0062 c1a.Draw()
0063 c1a.SaveAs("ProbQ_vc_vb.png")
0064 
0065 c1b = ROOT.TCanvas()
0066 h_probH_ddc.Draw("HISTO")
0067 h_probH_ddc.SetLineColor(632)
0068 h_probH_ddc.SetLineStyle(10)
0069 h_probQ_ddc.Draw("SAME")
0070 c1b.Draw()
0071 c1b.SaveAs("ProbH_vc_vb.png")