Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:30:45

0001 from ROOT import TFile, TCanvas, TPad, TH1F, TProfile, TProfile2D, TLine, gStyle, kOrange, kOpenCircle, kFullCircle
0002 from math import sqrt
0003 
0004 def isTrackerMuon(chain, index):
0005   if chain.muon_type[index] & 4 == 4:
0006     return True
0007 
0008 def inBarrel(chain, index):
0009   """
0010   Establish if the outer hit of a muon is in the barrel region.
0011   """
0012   if abs(chain.muon_outerPositionz[index]) < 108:
0013     return True
0014 
0015 def inForward(chain, index):
0016   """
0017   Establish if the outer hit of a muon is in the last layer of TEC in
0018   the forward region.
0019   """
0020   if abs(chain.muon_outerPositionz[index]) > 260:
0021     return True
0022 
0023 def isWithinZ(chain, index, value):
0024   """
0025   Establish if a muon is close to 0 in Z by value @value.
0026   """
0027   if abs(chain.muon_dz[index]) <= value:
0028     return True
0029 
0030 def overlayAndRatio(canvas, min_ratio, max_ratio, h1, h2):
0031   canvas.ResetAttPad()
0032   canvas.Clear()
0033   pad = TPad("main","main", 0, 0.3, 1, 1)
0034   pad.SetBottomMargin(0.05);
0035   pad.Draw()
0036   pad.cd()
0037   h1.Draw()
0038   h2.SetLineColor(kOrange+10)
0039   h2.Draw('SAME')
0040   canvas.cd()
0041   ratio = TPad("ratio","ratio", 0, 0.05, 1, 0.3);
0042   ratio.SetTopMargin(0.05);
0043   ratio.Draw()
0044   ratio.cd()
0045   if isinstance(h1, TProfile) and isinstance(h2, TProfile):
0046     h1_p = h1.ProjectionX()
0047     h2_p = h2.ProjectionX()
0048     h1_p.Divide(h2_p)
0049     h1_p.SetMaximum(max_ratio)
0050     h1_p.SetMinimum(min_ratio)
0051     h1_p.SetMarkerStyle(kFullCircle)
0052     h1_p.SetMarkerSize(0.6)
0053     h1_p.SetTitle('')
0054     h1_p.GetXaxis().SetLabelFont(42)
0055     h1_p.GetXaxis().SetLabelSize(0.1)
0056     h1_p.GetYaxis().SetLabelFont(42)
0057     h1_p.GetYaxis().SetLabelSize(0.1)
0058     h1_p.GetYaxis().SetNdivisions(505)
0059     equality = TLine()
0060     h1_p.Draw("SAME HIST P")
0061     equality.SetLineColor(kOrange+10)
0062     equality.DrawLine(h1_p.GetXaxis().GetXmin(), 1, h1_p.GetXaxis().GetXmax(), 1)
0063   
0064 
0065 
0066 gStyle.SetOptStat(0)
0067 
0068 files = ['trackTupla.root', 'trackTuplaNewMaterial.root']
0069 destinations = ['/afs/cern.ch/work/r/rovere/public/temporary/materialEffects/CurrentGeometry',
0070                 '/afs/cern.ch/work/r/rovere/public/temporary/materialEffects/NewMaterialGeometry']
0071 
0072 delta_p_xy = []
0073 delta_p_xy.append(TProfile2D('Delta_p_Old', 'Delta_p', 240, -120, 120, 240, -120, 120))
0074 delta_p_xy.append(TProfile2D('Delta_p_New', 'Delta_p', 240, -120, 120, 240, -120, 120))
0075 outermost_z = []
0076 outermost_z.append(TH1F("outermost_z_Old", "outermost_z", 560, -280, 280))
0077 outermost_z.append(TH1F("outermost_z_New", "outermost_z", 560, -280, 280))
0078 delta_p_outermost_z = []
0079 delta_p_outermost_z.append(TProfile("Delta_p_outermost_z_Old", "Delta_p_outermost_z", 560, -280, 280))
0080 delta_p_outermost_z.append(TProfile("Delta_p_outermost_z_New", "Delta_p_outermost_z", 560, -280, 280))
0081 delta_p_rz = []
0082 delta_p_rz.append(TProfile2D("Delta_p_rz_Old", "Delta_p_rz", 600, -300, 300, 120, 0, 120))
0083 delta_p_rz.append(TProfile2D("Delta_p_rz_New", "Delta_p_rz", 600, -300, 300, 120, 0, 120))
0084 delta_p_eta = []
0085 delta_p_eta.append(TProfile("Delta_p_eta_Old", "Delta_p_eta", 100, -2.5, 2.5))
0086 delta_p_eta.append(TProfile("Delta_p_eta_New", "Delta_p_eta", 100, -2.5, 2.5))
0087 delta_pt_eta = []
0088 delta_pt_eta.append(TProfile("Delta_pt_eta_Old", "Delta_pt_eta", 100, -2.5, 2.5))
0089 delta_pt_eta.append(TProfile("Delta_pt_eta_New", "Delta_pt_eta", 100, -2.5, 2.5))
0090 delta_pl_eta = []
0091 delta_pl_eta.append(TProfile("Delta_pl_eta_Old", "Delta_pl_eta", 100, -2.5, 2.5))
0092 delta_pl_eta.append(TProfile("Delta_pl_eta_New", "Delta_pl_eta", 100, -2.5, 2.5))
0093 
0094 c = TCanvas("c", "c", 1024, 1024)
0095 counter = 0
0096 for f in files:
0097   fh = TFile.Open(f)
0098   chain = fh.Get('Tracks')
0099   entries = chain.GetEntriesFast()
0100 
0101   for e in xrange(entries):
0102     nb = chain.GetEntry(e)
0103     num_muons = chain.nmuon
0104     for m in xrange(num_muons):
0105       delta_p_value = abs(chain.muon_innerMom[m] - chain.muon_outerMom[m])
0106       delta_pt_value = sqrt((chain.muon_innerMomx[m] - chain.muon_outerMomx[m])**2 +
0107                             (chain.muon_innerMomy[m] - chain.muon_outerMomy[m])**2)
0108       delta_pl_value = abs(chain.muon_innerMomz[m] - chain.muon_outerMomz[m])
0109       if isTrackerMuon(chain, m):
0110         outermost_z[counter].Fill(chain.muon_outerPositionz[m])
0111         delta_p_outermost_z[counter].Fill(chain.muon_outerPositionz[m], delta_p_value)
0112         delta_p_xy[counter].Fill(chain.muon_outerPositionx[m],
0113                                  chain.muon_outerPositiony[m],
0114                                  delta_p_value)
0115         delta_p_rz[counter].Fill(chain.muon_outerPositionz[m],
0116                                  sqrt(chain.muon_outerPositionx[m]**2+chain.muon_outerPositiony[m]**2),
0117                                  delta_p_value)
0118         delta_p_eta[counter].Fill(chain.muon_eta[m], delta_p_value)
0119         delta_pt_eta[counter].Fill(chain.muon_eta[m], delta_pt_value)
0120         delta_pl_eta[counter].Fill(chain.muon_eta[m], delta_pl_value)
0121 
0122   outermost_z[counter].Draw()
0123   c.SaveAs("%s/OuterMostMuonHitInTracker.png" % destinations[counter])
0124 
0125   delta_p_outermost_z[counter].SetMaximum(0.2)
0126   delta_p_outermost_z[counter].Draw()
0127   c.SaveAs("%s/DeltaP_OuterMostMuonHitInTracker.png" % destinations[counter])
0128 
0129   delta_p_rz[counter].SetMaximum(4)
0130   delta_p_rz[counter].Draw('COLZ')
0131   c.SaveAs("%s/DeltaP_RZ.png" % destinations[counter])
0132   
0133   delta_p_eta[counter].SetMaximum(0.2)
0134   delta_p_eta[counter].Draw('COLZ')
0135   c.SaveAs("%s/DeltaP_ETA.png" % destinations[counter])
0136 
0137   for val in [1, 0.5, 0.1]:
0138     delta_p_xy[counter].SetMaximum(val)
0139     delta_p_xy[counter].Draw("COLZ")
0140     c.SaveAs("%s/DeltaP_XY_Max%f.png" % (destinations[counter],val))
0141   counter +=1
0142 
0143 overlayAndRatio(c, 0.85, 1.15, *delta_p_outermost_z)
0144 c.SaveAs("%s/DeltaP_OuterMostMuonHitInTracker_Comparison.png" % destinations[0])
0145 c.SaveAs("%s/DeltaP_OuterMostMuonHitInTracker_Comparison.png" % destinations[1])
0146 
0147 overlayAndRatio(c, 0.85, 1.15, *delta_p_eta)
0148 c.SaveAs("%s/DeltaP_Eta_Comparison.png" % destinations[0])
0149 c.SaveAs("%s/DeltaP_Eta_Comparison.png" % destinations[1])
0150 
0151 overlayAndRatio(c, 0.999, 1.001, *delta_pt_eta)
0152 c.SaveAs("%s/DeltaPt_Eta_Comparison.png" % destinations[0])
0153 c.SaveAs("%s/DeltaPt_Eta_Comparison.png" % destinations[1])
0154 
0155 delta_pl_eta[0].SetMaximum(0.2)
0156 delta_pl_eta[1].SetMaximum(0.2)
0157 overlayAndRatio(c, 0.85, 1.15, *delta_pl_eta)
0158 c.SaveAs("%s/DeltaPl_Eta_Comparison.png" % destinations[0])
0159 c.SaveAs("%s/DeltaPl_Eta_Comparison.png" % destinations[1])
0160 
0161