Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:26:57

0001 from __future__ import print_function
0002 from builtins import range
0003 import sys
0004 import numpy as np
0005 
0006 import ROOT
0007 ROOT.PyConfig.IgnoreCommandLineOptions = True
0008 from optparse import OptionParser
0009 from RecoTauTag.Configuration.tools.DisplayManager import DisplayManager
0010 
0011 ROOT.gROOT.SetBatch(True)
0012 ROOT.gStyle.SetOptStat(0)
0013 
0014 colours = [1, 2, 3, 6, 8]
0015 styles = [1, 2, 3, 4, 5]
0016 
0017 def findTree(f):
0018     for key in f.GetListOfKeys():
0019         if key.GetName() == "Events":
0020             tree = f.Get(key.GetName())
0021             if isinstance(tree, ROOT.TTree):
0022                 return tree
0023             elif isinstance(tree, ROOT.TDirectory):
0024                 return findTree(tree)
0025     print('Failed to find TTree Events in file', f)
0026     return None
0027 
0028 def applyHistStyle(h, i):
0029     h.SetLineColor(colours[i])
0030     h.SetLineStyle(styles[i])
0031     h.SetLineWidth(3)
0032     h.SetStats(False)
0033 
0034 def comparisonPlots(u_names, trees, titles, pname='compareTauVariables.pdf', ratio=True):
0035     display = DisplayManager(pname, ratio)
0036     for branch in u_names:
0037         nbins = 50
0038         min_x = min(t.GetMinimum(branch) for t in trees)
0039         max_x = max(t.GetMaximum(branch) for t in trees)
0040         title_x = branch
0041         if min_x == max_x or all(t.GetMinimum(branch) == t.GetMaximum(branch) for t in trees):
0042             continue
0043         if min_x < -900 and max_x < -min_x * 1.5:
0044             min_x = - max_x
0045         min_x = min(0., min_x)
0046         hists = []
0047         for i, t in enumerate(trees):
0048             h_name = branch+t.GetName()+str(i)
0049             h = ROOT.TH1F(h_name, branch, nbins, min_x, max_x + (max_x - min_x) * 0.01)
0050             h.Sumw2()
0051             h.GetXaxis().SetTitle(title_x)
0052             h.GetYaxis().SetTitle('Entries')
0053             applyHistStyle(h, i)
0054             t.Project(h_name, branch, '1') # Should introduce weight...
0055             hists.append(h)
0056         display.Draw(hists, titles)
0057 
0058 def interSect(tree1, tree2, var='ull_dumpTauVariables_EventNumber_DUMP.obj', common=False, save=False,  titles=[]):
0059     tlist1 = ROOT.TEntryList()
0060     tlist2 = ROOT.TEntryList()
0061     tree1.Draw(var)
0062     r_evt1 = tree1.GetV1()
0063     if len(titles) > 0 and titles[0] == 'Reference':
0064         evt1 = np.array([int(r_evt1[i]) & 0xffffffff for i in range(tree2.GetEntries())], dtype=int)
0065     else:
0066         evt1 = np.array([r_evt1[i] for i in range(tree1.GetEntries())], dtype=int)
0067     tree2.Draw(var)
0068     r_evt2 = tree2.GetV1()
0069     if len(titles) > 1 and titles[1] == 'Reference':
0070         evt2 = np.array([int(r_evt2[i]) & 0xffffffff for i in range(tree2.GetEntries())], dtype=int)
0071     else:
0072         evt2 = np.array([int(r_evt2[i]) for i in range(tree2.GetEntries())], dtype=int)
0073     if common:
0074         indices1 = np.nonzero(np.in1d(evt1, evt2))
0075         indices2 = np.nonzero(np.in1d(evt2, evt1))
0076     else:
0077         indices1 = np.nonzero(np.in1d(evt1, evt2) == 0)
0078         indices2 = np.nonzero(np.in1d(evt2, evt1) == 0)
0079     if save:
0080         if len(titles) < 2:
0081             titles = ['tree1', 'tree2']
0082         evt1[indices1].tofile(titles[0]+'.csv', sep=',', format='%d')
0083         evt2[indices2].tofile(titles[1]+'.csv', sep=',', format='%d')
0084     for ind1 in indices1[0]:
0085         tlist1.Enter(ind1)
0086     for ind2 in indices2[0]:
0087         tlist2.Enter(ind2)
0088     return tlist1, tlist2
0089 
0090 def scanForDiff(tree1, tree2, branch_names, scan_var='floats_dumpTauVariables_pt_DUMP.obj', index_var='ull_dumpTauVariables_EventNumber_DUMP.obj'):
0091     tree2.BuildIndex(index_var)
0092     diff_events = []
0093     for entry_1 in tree1:
0094         ind = int(getattr(tree1, index_var))
0095         tree2.GetEntryWithIndex(ind)
0096         var1 = getattr(tree1, scan_var)
0097         var2 = getattr(tree2, scan_var)
0098         if tree1.evt != tree2.evt:
0099             continue
0100         if round(var1, 6) != round(var2, 6): 
0101             diff_events.append(ind)
0102             print('Event', ind)
0103             for branch in branch_names:
0104                 v1 = getattr(tree1, branch)
0105                 v2 = getattr(tree2, branch)
0106                 if round(v1, 6) != round(v2, 6) and v1 > -99.:
0107                     print('{b:>43}: {v1:>8.4f}, {v2:>8.4f}'.format(b=branch, v1=v1, v2=v2))
0108             print()
0109     print('Found', len(diff_events), 'events with differences in', scan_var)
0110     print(diff_events)
0111 
0112 
0113 if __name__ == '__main__':
0114         
0115     usage = '''
0116 %prog [options] arg1 arg2 ... argN
0117     Compares first found trees in N different root files; 
0118     in case of two root files, additional information about the event overlap
0119     can be calculated.
0120     Example run commands:
0121 > python compareTauVariables.py dumpTauVariables_CMSSW_8_0_X.root dumpTauVariables_CMSSW_8_1_X.root -t CMSSW_8_0_X,CMSSW_8_1_X
0122 '''
0123 
0124     parser = OptionParser(usage=usage)
0125 
0126     parser.add_option('-t', '--titles', type='string', dest='titles', default='Reference,Test', help='Comma-separated list of titles for the N input files (e.g. CMSSW_8_0_X,CMSSW_8_1_X)')
0127     parser.add_option('-i', '--do-intersection', dest='do_intersect', action='store_true', default=False, help='Produce plots for events not in common')
0128     parser.add_option('-c', '--do-common', dest='do_common', action='store_true', default=False, help='Produce plots for events in common')
0129     parser.add_option('-r', '--do-ratio', dest='do_ratio', action='store_true', default=False, help='Show ratio plots')
0130     parser.add_option('-d', '--diff', dest='do_diff', action='store_true', default=False, help='Print events where single variable differs')
0131     parser.add_option('-v', '--var-diff', dest='var_diff', default='floats_dumpTauVariables_pt_DUMP.obj', help='Variable for printing single event diffs')
0132 
0133     (options,args) = parser.parse_args()
0134 
0135     if len(args) < 2:
0136         print('provide at least 2 input root files')
0137         sys.exit(1)
0138 
0139     titles = options.titles.split(',')    
0140     if len(titles) < len(args):
0141         print('Provide at least as many titles as input files')
0142         sys.exit(1)
0143 
0144     for i, arg in enumerate(args):
0145         if arg.endswith('.txt'):
0146             f_txt = open(arg)
0147             for line in f_txt.read().splitlines():
0148                 line.strip()
0149                 if line.startswith('/afs'):
0150                     args[i] = line
0151                     break
0152 
0153     tfiles = [ROOT.TFile(arg) for arg in args]    
0154     trees = [findTree(f) for f in tfiles]
0155     b_names = []
0156     b_event_number = ''
0157     for t in trees:
0158         bs = []
0159         for b in t.GetListOfBranches():
0160             if "DUMP" in b.GetName():
0161                 bs.append(b.GetName().replace("'","\'")+"obj")
0162         b_names.append(set(bs))
0163         if "EventNumber" in bs :
0164             b_event_number = bs
0165     u_names = set.intersection(*b_names)
0166     u_names = sorted(u_names)
0167 
0168     print('Making plots for all common branches (', len(u_names), ')')
0169     comparisonPlots(u_names, trees, titles, 'compareTauVariables.pdf', options.do_ratio)
0170 
0171     if len(trees) == 2 and options.do_intersect:
0172         intersect = interSect(trees[0], trees[1], var=b_event_number, save=True, titles=titles)
0173         trees[0].SetEntryList(intersect[0])
0174         trees[1].SetEntryList(intersect[1])
0175         if not all(l.GetN() == 0 for l in intersect):
0176             comparisonPlots(u_names, trees, titles, 'compareTauVariables_intersect.pdf', options.do_ratio)
0177 
0178     if len(trees) == 2 and options.do_common:
0179         intersect = interSect(trees[0], trees[1], var=b_event_number, common=True)
0180         trees[0].SetEntryList(intersect[0])
0181         trees[1].SetEntryList(intersect[1])
0182         comparisonPlots(u_names, trees, titles, 'compareTauVariables_common.pdf', options.do_ratio)
0183 
0184     if len(trees) == 2 and options.do_diff:
0185         scanForDiff(trees[0], trees[1], u_names, scan_var=options.var_diff, index_var=b_event_number)