import ROOT import sys import numpy as np import os #The Purpose of this code is to produce histograms from both data and MC root files, and superimpose the created #histograms onto each other. This code can work in two ways. The first way is by putting in system arguments for the #files that you want to produce the histograms for.The order of the system arguments matters. Put all the data files #first and then put the MC simulation files after.The files need to be ordered the same way in order to superimpose #them correctly. For example the first data file is paired with the first MC file, the second data with the second MC #etc. The second way that the code works is by going through every file in two different directories given in the #code. It produces the histograms for these files and then superimposes them onto each other. #This first directory is for the data. Change the directory to where your data files are stored directory_1='data_exotics_1' #This second directory is for the MC. Change the directory to where you MC files are stored directory_2='mc_exotics_1' #This is a function to create the histograms based off the system arguments given. def hist_func_select(FileName,i): #This is used to create the output file with the correct numbering. Since the first half of the files should be data #and the second half should be MC after going through the first half of the files the histogram names need to change if i<=(len(sys.argv)-1)/2: hist='data_hist_%s.root'%i else: j=i-int((len(sys.argv)-1)/2) hist='mc_hist_%s.root'%j #This is to read the tree from the input file inFile = ROOT.TFile.Open(FileName, "READ" ) #This is the tree in the root file that the histograms will be produced from. This should be changed depending on the #root file tree = inFile.Get ( "trees_DV_" ) #This is to produce the histograms. This should be adjusted depending on the root file jet_et_hist = ROOT.TH1D ( "jet_et" ," HLT jet et " ,20 ,0 ,200) MSVtx_r_hist = ROOT.TH1D ( "MSVtx_r" ," MSVtx_r " ,100 ,0 ,7000) muon_pt_hist = ROOT.TH1D ( "muon_pT" ," muon_pT " ,100 ,-10 ,100) met_phi_hist = ROOT.TH1D ( "met_phi" ," met_phi " ,100 ,-4 ,4) track_eta_hist = ROOT.TH1D ( "track_eta" ," track_eta " ,100 ,-4 ,4) MSVtx_ntrks_hist = ROOT.TH1D ( "MSVtx_nTracks" ," MSVtx_nTracks " ,15 ,0 ,20) #Sumw2 is used to deal with the uncertainty in the histogram jet_et_hist.Sumw2() MSVtx_r_hist.Sumw2() muon_pt_hist.Sumw2() met_phi_hist.Sumw2() track_eta_hist.Sumw2() MSVtx_ntrks_hist.Sumw2() #This loop is used to fill the hisograms for entryNum in range (0,tree.GetEntries()): tree.GetEntry(entryNum) jet = getattr(tree, "HLT_jet_ET") vtx_r = getattr(tree, "MSVtx_R") muon_pt = getattr(tree, "muon_pT") met_phi = getattr(tree,"met_phi") track_eta = getattr(tree,"track_eta") MSVtx_trks = getattr(tree, "MSVtx_nTrks") if len(jet)!=0: #jet_et_hist.Fill (jet[0]) for i in range (len(jet)): jet_et_hist.Fill (jet[i]) if len(vtx_r)!=0: #jet_et_hist.Fill (vtx_r[0]) for i in range (len(vtx_r)): MSVtx_r_hist.Fill (vtx_r[i]) if len(muon_pt)!=0: for i in range(len(muon_pt)): muon_pt_hist.Fill (muon_pt[i]) if len(track_eta)!=0: for i in range(len(track_eta)): track_eta_hist.Fill(track_eta[i]) if len(MSVtx_trks)!=0: for i in range (len(MSVtx_trks)): MSVtx_ntrks_hist.Fill (MSVtx_trks[i]) met_phi_hist.Fill (met_phi) jet_et_hist . SetDirectory (0) MSVtx_r_hist . SetDirectory (0) muon_pt_hist . SetDirectory(0) met_phi_hist . SetDirectory(0) track_eta_hist . SetDirectory(0) MSVtx_ntrks_hist . SetDirectory(0) inFile.Close () outHistFile = ROOT . TFile . Open ( hist , "RECREATE" ) outHistFile . cd () jet_et_hist . Write () MSVtx_r_hist . Write() muon_pt_hist . Write() met_phi_hist . Write() track_eta_hist . Write() MSVtx_ntrks_hist . Write() outHistFile . Close () #This is a function to create the histograms by going through all the files in the directories given at the beginning #of the code def hist_func_all(FileName,i,hist_title): hist='%s_hist_%s.root'%(hist_title,i) #This is to read the tree from the input file. This should be changed depending on the root file tree = FileName.Get ( "trees_DV_" ) #This is to produce the histograms jet_et_hist = ROOT.TH1D ( "jet_et" ," HLT jet et "+str(i) ,20 ,0 ,200) MSVtx_r_hist = ROOT.TH1D ( "MSVtx_r" ," MSVtx_r " ,100 ,0 ,7000) muon_pt_hist = ROOT.TH1D ( "muon_pT" ," muon_pT " ,100 ,-10 ,100) met_phi_hist = ROOT.TH1D ( "met_phi" ," met_phi " ,100 ,-4 ,4) track_eta_hist = ROOT.TH1D ( "track_eta" ," track_eta " ,100 ,-4 ,4) MSVtx_ntrks_hist = ROOT.TH1D ( "MSVtx_nTracks" ," MSVtx_nTracks " ,15 ,0 ,20) #Sumw2 is used to deal with the uncertainty in the histogram jet_et_hist.Sumw2() MSVtx_r_hist.Sumw2() muon_pt_hist.Sumw2() met_phi_hist.Sumw2() track_eta_hist.Sumw2() MSVtx_ntrks_hist.Sumw2() #This loop is used to fill the hisograms for entryNum in range (0,tree.GetEntries()): tree.GetEntry(entryNum) jet = getattr(tree, "HLT_jet_ET") vtx_r = getattr(tree, "MSVtx_R") muon_pt = getattr(tree, "muon_pT") met_phi = getattr(tree,"met_phi") track_eta = getattr(tree,"track_eta") MSVtx_trks = getattr(tree, "MSVtx_nTrks") if len(jet)!=0: #jet_et_hist.Fill (jet[0]) for i in range (len(jet)): jet_et_hist.Fill (jet[i]) if len(vtx_r)!=0: #jet_et_hist.Fill (vtx_r[0]) for i in range (len(vtx_r)): MSVtx_r_hist.Fill (vtx_r[i]) if len(muon_pt)!=0: for i in range(len(muon_pt)): muon_pt_hist.Fill (muon_pt[i]) if len(track_eta)!=0: for i in range(len(track_eta)): track_eta_hist.Fill(track_eta[i]) if len(MSVtx_trks)!=0: for i in range (len(MSVtx_trks)): MSVtx_ntrks_hist.Fill (MSVtx_trks[i]) met_phi_hist.Fill (met_phi) jet_et_hist . SetDirectory (0) MSVtx_r_hist . SetDirectory (0) muon_pt_hist . SetDirectory(0) met_phi_hist . SetDirectory(0) track_eta_hist . SetDirectory(0) MSVtx_ntrks_hist . SetDirectory(0) inFile.Close () outHistFile = ROOT.TFile.Open(hist,"RECREATE" ) outHistFile . cd () #Normalize the histograms if (jet_et_hist.Integral())!=0: jet_et_hist.Scale(1/jet_et_hist.Integral()) if (MSVtx_r_hist.Integral())!=0: MSVtx_r_hist.Scale(1/MSVtx_r_hist.Integral()) if (muon_pt_hist.Integral())!=0: muon_pt_hist.Scale(1/muon_pt_hist.Integral()) if (met_phi_hist.Integral())!=0: met_phi_hist.Scale(1/met_phi_hist.Integral()) if (track_eta_hist.Integral())!=0: track_eta_hist.Scale(1/track_eta_hist.Integral()) if (MSVtx_ntrks_hist.Integral())!=0: MSVtx_ntrks_hist.Scale(1/MSVtx_ntrks_hist.Integral()) #Writing the histograms jet_et_hist . Write () MSVtx_r_hist . Write() muon_pt_hist . Write() met_phi_hist . Write() track_eta_hist . Write() MSVtx_ntrks_hist . Write() #closing the file outHistFile . Close () #This function is used to superimpose the data and mc simulations into a pdf def plot_func(histFileName_1, histFileName_2,i): #plot_name contains the names of all the histograms that were created plot_name=['MSVtx_nTracks','jet_et','met_phi','muon_pT','track_eta','MSVtx_r'] #This loops through all the histograms in order to create all the plots for k in range(len(plot_name)): #Creates the name for the pdf file plotFileName = '%s_plot_%s.pdf'%(plot_name[k],i) #Opens the histograms that were made previously histFile_1 = ROOT.TFile.Open(histFileName_1, " READ " ) histFile_2 = ROOT.TFile.Open(histFileName_2, " READ " ) #Gets the specific plot out of the hist files dataHisto = histFile_1.Get(plot_name[k]) mcHisto = histFile_2.Get(plot_name[k]) #These if statements stop the code from crashing if there is an issue with the hist files if not dataHisto: print("Failed to get data histo") if not mcHisto: print ("Failed to get MC histogram") #closes hist files dataHisto.SetDirectory(0) mcHisto.SetDirectory(0) histFile_1.Close() histFile_2.Close() #For root plots you need a canvas. Here we are making a canvis canvas = ROOT.TCanvas("canvas") canvas.cd() #This line is nessesary to make several plots in a single pdf canvas.Print(plotFileName+"[") #Removes statitics box. Remove these two lines if you want the statistics boxes mcHisto.SetStats(0) dataHisto.SetStats(0) #Here is the line width and color being set mcHisto.SetLineColor(ROOT.kRed) dataHisto.SetLineColor(ROOT.kBlack) mcHisto.SetLineWidth(2) dataHisto.SetLineWidth(2) #Here is the Yaxis being set mcHisto.GetYaxis().SetTitle(" Number of events " ) dataHisto.GetYaxis().SetTitle( " Number of events " ) #Makes the y axis logrithmic canvas.SetLogy(True) #Creates the MC plot mcHisto.Draw('h') canvas.Print(plotFileName) #Creates the data plot dataHisto.Draw('h') canvas.Print(plotFileName) #This line scales the MC to the data mcHisto.Scale(dataHisto.Integral()/mcHisto.Integral()) #Creates the super imposed plot mcHisto.Draw('h') dataHisto.Draw('h,same') #Here is the legend being set legend = ROOT.TLegend (0.7,0.6,0.85,0.75) legend.AddEntry(mcHisto ,"MC") legend.AddEntry(dataHisto ,"Data") legend.Draw("same") canvas.Print(plotFileName) canvas.Print(plotFileName+"]") #The purpose of this if statement is to determine if you put in system arguments or not. if len(sys.argv)>1: #This loop is to loop all the input files through the function for i in range(1,len(sys.argv)): FileName = sys.argv [i] hist_func_select(FileName,i) for j in range(1,int((len(sys.argv)-1)/2)+1): hist_1='data_hist_%s.root'%j hist_2='mc_hist_%s.root'%j plot_func(hist_1,hist_2,j) else: #This loop is to loop all the input data files through the function and produce histograms from them #This i=1 is to keep track of the file number i=1 for filename in os.listdir(directory_1): filename = os.path.join(directory_1, filename) print(filename) inFile = ROOT.TFile.Open(filename, "READ" ) hist_func_all(inFile,i,'data') #The i=i+1 increases the file number by 1 for each time it goes through the loop i=i+1 #This loop is to loop all the input MC files through the function and produce histograms from them #The file number is reset with the i=1 again i=1 for filename in os.listdir(directory_2): filename = os.path.join(directory_2, filename) print(filename) inFile = ROOT.TFile.Open(filename, "READ" ) hist_func_all(inFile,i,'mc') i=i+1 #This loop is used to produce the superimposed plots for each of the histograms created for j in range(1,i): hist_1='data_hist_%s.root'%j hist_2='mc_hist_%s.root'%j plot_func(hist_1,hist_2,j)