ATLAS Offline Software
LayerMaterialUtils.py
Go to the documentation of this file.
1 # Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
2 
3 import math
4 from ROOT import TCanvas, TH2F, TH3F, TProfile, gDirectory, gROOT, kCanDelete, kMustCleanup, kRed, kWhite
5 
6 def drawCylinderLayer(g4CylTree,cylLayer,tgCylTree,saveLocation) :
7 
8  gROOT.cd()
9  # create a canvas first
10  clname = cylLayer.GetName()
11 
12  # now fill it
13  for evt in cylLayer :
14 
15  radius = evt.cylRadius
16  position = evt.cylPosZ
17  halfz = evt.cylHalfZ
18  binsz = evt.cylMatBinsZ
19  binsphi = evt.cylMatBinsPhi
20 
21  numbinsz = 100
22  numbinsphi = 100
23  if binsz > numbinsz : numbinsz = binsz
24  if binsphi > numbinsphi : numbinsphi = binsphi
25 
26 
27  print ('----------------------------------------------------------------')
28  print (' Layer: ', clname, ' ( LayerIndex : ',evt.cylLayerIndex, ')')
29  print (' r, halflength : ', radius, halfz)
30  print (' bins (z/phi) : ', evt.cylMatBinsZ, evt.cylMatBinsPhi)
31 
32 
33  # create the x0 histogram
34  clmap_x0 = TH2F('clmap_x0_'+clname,clname,binsz,-halfz,halfz,binsphi,-math.pi,math.pi)
35  x0max = 0.
36  # loop over bins and fill it
37  for ibin in range(evt.cylMatBinsTotal) :
38  # first get the ibin phi
39  ibinz = ibin/numbinsphi
40  ibinphi = ibin - ibinz*numbinsphi
41  x0 = evt.cylMatX0[ibin]
42  if x0 > x0max : x0max = x0
43  clmap_x0.SetBinContent(ibinz+1,ibinphi+1,x0)
44  pass
45 
46  # and draw it
47  canvasX0map = TCanvas(clname+'_x0',clname,100,100,800,1000)
48  canvasX0map.Divide(2,3)
49 
50  # ---------------------------------------------------------
51  # draw profiles first to get a good estimate for the maximum
52  canvasX0map.cd(5)
53  x0profZG4 = TProfile(clname+'x0profZG4',clname+'x0profZG4',numbinsz,-halfz+position,halfz+position) # noqa: F841
54  x0profZTG = TProfile(clname+'x0profZTG',clname+'x0profZTG',numbinsz,-halfz+position,halfz+position) # noqa: F841
55  g4CylTree.Draw('PathInX0:DensedHitZ>>'+clname+'x0profZG4','','prof')
56  tgCylTree.Draw('PathInX0:DensedHitZ>>'+clname+'x0profZTG','','prof')
57  x0profG4z = gDirectory.Get(clname+'x0profZG4')
58  x0profTGz = gDirectory.Get(clname+'x0profZTG')
59  # style
60  x0profTGz.SetMarkerStyle(8)
61  x0profTGz.SetMarkerSize(0.5)
62  x0profTGz.SetMarkerColor(kRed)
63  x0profTGz.SetLineColor(kRed)
64  x0profG4z.GetYaxis().SetTitle('t [X_{0}]')
65  x0profG4z.GetXaxis().SetTitle('z [mm]')
66  x0maxProf = 1.5*x0profG4z.GetMaximum()
67  x0profG4z.GetYaxis().SetRangeUser(0.,x0maxProf)
68  x0profG4z.Draw('hist')
69  x0profTGz.SetMarkerColor(kRed)
70  x0profTGz.SetMarkerSize(0.8)
71  x0profTGz.Draw('pe,same')
72  canvasX0map.cd(6)
73  x0profPhiG4 = TProfile(clname+'x0profPhiG4',clname+'x0profPhiG4',numbinsphi,-math.pi,math.pi) # noqa: F841
74  x0profPhiTG = TProfile(clname+'x0profPhiTG',clname+'x0profPhiTG',numbinsphi,-math.pi,math.pi) # noqa: F841
75  g4CylTree.Draw('PathInX0:DensedHitPhi>>'+clname+'x0profPhiG4','','prof')
76  tgCylTree.Draw('PathInX0:DensedHitPhi>>'+clname+'x0profPhiTG','','prof')
77  x0profG4phi = gDirectory.Get(clname+'x0profPhiG4')
78  x0profTGphi = gDirectory.Get(clname+'x0profPhiTG')
79  # style
80  x0profTGphi.SetMarkerStyle(8)
81  x0profTGphi.SetMarkerSize(0.5)
82  x0profTGphi.SetMarkerColor(kRed)
83  x0profTGphi.SetLineColor(kRed)
84  x0profG4phi.GetXaxis().SetTitle('#phi')
85  x0profG4phi.GetYaxis().SetTitle('t [X_{0}]')
86  if 1.5*x0profG4phi.GetMaximum() > x0maxProf : x0maxProf = 1.5*x0profG4phi.GetMaximum()
87  x0profG4phi.GetYaxis().SetRangeUser(0.,1.5*x0profG4phi.GetMaximum())
88  x0profG4phi.Draw('hist')
89  x0profTGphi.SetMarkerColor(kRed)
90  x0profTGphi.SetMarkerSize(0.8)
91  x0profTGphi.Draw('pe,same')
92  # ---- now the maps
93  canvasX0map.cd(1)
94  try :
95  g4CylTree.Draw('PathInX0>>'+clname+'range')
96  x0range = gDirectory.Get(clname+'range')
97  x0maxProf = 1.2*x0range.GetXaxis().GetXmax()
98  except TypeError :
99  print (' -> TTree does not exist in TrackingGeometry ')
100  x0mapG4 = TH3F(clname+'x0mapG4',clname+'x0mapG4', numbinsphi,-math.pi,math.pi,numbinsz,-halfz+position,halfz+position,50,0.,x0maxProf)
101  try :
102  g4CylTree.Draw('PathInX0:DensedHitZ:DensedHitPhi>>'+clname+'x0mapG4','','')
103  x0mapG4 = gDirectory.Get(clname+'x0mapG4')
104  x0mapG4Prof = x0mapG4.Project3DProfile('xy')
105  x0mapG4Prof.GetXaxis().SetTitle('z [mm]')
106  x0mapG4Prof.GetYaxis().SetTitle('#phi')
107  x0mapG4Prof.GetZaxis().SetTitle('t [X_{0}]')
108  x0mapG4Prof.GetZaxis().SetTitleOffset(1.5*x0mapG4Prof.GetZaxis().GetTitleOffset())
109  x0mapG4Prof.GetZaxis().SetRangeUser(0.,x0maxProf)
110  x0mapG4Prof.SetBit(kCanDelete)
111  x0mapG4Prof.SetBit(kMustCleanup)
112  x0mapG4Prof.Draw('colz')
113  except TypeError :
114  print (' -> TTree does not exist in TrackingGeometry ')
115  canvasX0map.cd(2)
116  x0mapTG = TH3F(clname+'x0mapTG',clname+'x0mapTG',numbinsphi,-math.pi,math.pi,numbinsz,-halfz+position,halfz+position,50,0.,x0maxProf)
117  try :
118  tgCylTree.Draw('PathInX0:DensedHitZ:DensedHitPhi>>'+clname+'x0mapTG','','')
119  x0mapTG = gDirectory.Get(clname+'x0mapTG')
120  x0mapTGProf = x0mapTG.Project3DProfile('xy')
121  x0mapTGProf.GetXaxis().SetTitle('z [mm]')
122  x0mapTGProf.GetYaxis().SetTitle('#phi')
123  x0mapTGProf.GetZaxis().SetTitle('t [X_{0}]')
124  x0mapTGProf.GetZaxis().SetTitleOffset(1.5*x0mapTGProf.GetZaxis().GetTitleOffset())
125  x0mapTGProf.GetZaxis().SetRangeUser(0.,x0maxProf)
126  x0mapTGProf.SetBit(kCanDelete)
127  x0mapTGProf.SetBit(kMustCleanup)
128  x0mapTGProf.Draw('colz')
129  except TypeError :
130  print (' -> TTree does not exist in Geant4 ')
131  canvasX0map.cd(3)
132  clmap_x0.GetZaxis().SetRangeUser(0.,x0max)
133  clmap_x0.GetXaxis().SetTitle('z [mm]')
134  clmap_x0.GetYaxis().SetTitle('#phi')
135  clmap_x0.GetZaxis().SetTitle('t [X_{0}]')
136  clmap_x0.GetZaxis().SetTitleOffset(1.2*clmap_x0.GetZaxis().GetTitleOffset())
137  clmap_x0.Draw('colz')
138  # ---------------------------------------------------------
139  canvasX0map.cd(4)
140  try :
141  x0mapDiff = x0mapTGProf.Clone()
142  x0mapDiff.Add(x0mapG4Prof,-1)
143  x0mapDiff.GetXaxis().SetTitle('z [mm]')
144  x0mapDiff.GetYaxis().SetTitle('#phi')
145  x0mapDiff.GetZaxis().SetTitle('#Delta(t) [X_{0}]')
146  x0mapDiff.GetZaxis().SetRangeUser(-1.,1.)
147  x0mapDiff.GetZaxis().SetTitleOffset(1.5*x0mapDiff.GetZaxis().GetTitleOffset())
148  x0mapDiff.Divide(x0mapG4Prof)
149  x0mapDiff.Draw('colz')
150  except UnboundLocalError :
151  print (' -> No TrackingGeometry/G4 histograms to divide')
152 
153 
154  canvasX0map.SaveAs(saveLocation+'/LayerMaterial_CylinderLayer_'+str(evt.cylLayerIndex)+'.png')
155 
156  canvasAnalysis = TCanvas(clname+'_analyse',clname,100,100,600,1200)
157  canvasAnalysis.Divide(1,3)
158  canvasAnalysis.cd(1)
159  cylLayer.Draw('LayerMaterialL0:LayerMaterialX0','','box')
160  canvasAnalysis.cd(2)
161  cylLayer.Draw('LayerMaterialZARho:LayerMaterialX0','','box')
162  canvasAnalysis.cd(3)
163  cylLayer.Draw('LayerMaterialZARho:LayerMaterialL0','','box')
164  canvasAnalysis.SaveAs(saveLocation+'/'+clname+'_analysis.png')
165 
166  # memory cleanup
167  gDirectory.Delete(clname+'x0mapTG')
168  gDirectory.Delete(clname+'x0mapG4')
169  gDirectory.Delete('clmap_x0_'+clname)
170  gDirectory.Delete('clmap_x0_pz_'+clname)
171  gDirectory.Delete('clmap_x0_pphi_'+clname)
172 
173  #clmap_x0_profZ.Draw('hist')
174  #canvasX0map.cd(3)
175  #clmap_x0_profPhi.Draw('hist')
176 
177 
178 def drawDiscLayer(g4DiscTree,discLayer,tgDiscTree,saveLocation) :
179 
180  gROOT.cd()
181  # create a canvas first
182  dlname = discLayer.GetName()
183  # now fill it
184  for evt in discLayer :
185 
186  innerR = evt.discInnerR
187  outerR = evt.discOuterR
188  position = evt.discPositionZ
189  binsR = evt.discMatBinsR
190  binsPhi = evt.discMatBinsPhi
191  scalor = math.pi/outerR
192 
193  numbinsr = 100
194  numbinsphi = 100
195  if binsR > numbinsr : numbinsr = binsR
196  if binsPhi > numbinsphi : numbinsphi = binsPhi
197 
198  print ('----------------------------------------------------------------')
199  print (' Layer: ', dlname, ' ( LayerIndex : ',evt.discLayerIndex, ')')
200  print (' position : ', position)
201  print (' innerR/outerR : ', innerR, outerR)
202  print (' bins (r/phi) : ', binsR, binsPhi)
203 
204  # POLAR plot needs to be in -pi,pi,-pi,pi, NOT yet fixed
205 
206  # create the x0 histogram
207  dlmap_x0 = TH2F('dlmap_x0_'+dlname,dlname,binsR,-math.pi,math.pi,binsPhi,-math.pi,math.pi)
208  # the steps
209  stepR = (outerR-innerR)/binsR
210  stepPhi = 2*math.pi/binsPhi
211  x0max = 0.
212  # loop over bins and fill it
213  for ibin in range(evt.discMatBinsTotal) :
214  ibinphi = ibin/binsR
215  ibinr = ibin - ibinphi*binsR
216  #ibinz = ibin/binsphi
217  #ibinphi = ibin - ibinz*binsphi
218  x0 = evt.discMatX0[ibin]
219  if x0 > x0max : x0max = x0
220  iphival = -math.pi+(ibinphi+0.5)*stepPhi
221  irvalscal = (innerR + (ibinr+0.5)*stepR)*scalor
222  dlmap_x0.Fill(irvalscal,iphival,x0)
223  # z / phi value for the profile plots
224 
225  pass
226 
227  canvasX0map = TCanvas(dlname,dlname,100,100,800,1100)
228  canvasX0map.Divide(2,3)
229  # draw profiles tog et good estimate for x0max first
230  x0maxProf = x0max
231  canvasX0map.cd(5)
232  x0profRG4 = TProfile(dlname+'x0profRG4',dlname+'x0profRG4',numbinsr,innerR,outerR) # noqa: F841
233  x0profRTG = TProfile(dlname+'x0profRTG',dlname+'x0profRTG',numbinsr,innerR,outerR) # noqa: F841
234  g4DiscTree.Draw('PathInX0:DensedHitR>>'+dlname+'x0profRG4','','prof')
235  tgDiscTree.Draw('PathInX0:DensedHitR>>'+dlname+'x0profRTG','','prof')
236  x0profG4r = gDirectory.Get(dlname+'x0profRG4')
237  x0profTGr = gDirectory.Get(dlname+'x0profRTG')
238  # style
239  x0profTGr.SetMarkerStyle(8)
240  x0profTGr.SetMarkerSize(0.5)
241  x0profTGr.SetMarkerColor(kRed)
242  x0profTGr.SetLineColor(kRed)
243  x0profG4r.GetXaxis().SetTitle('r [mm]')
244  x0profG4r.GetYaxis().SetTitle('t [X_{0}]')
245  x0maxProf = 1.5*x0profG4r.GetMaximum()
246  x0profG4r.GetYaxis().SetRangeUser(0.,x0maxProf)
247  x0profG4r.Draw('hist')
248  x0profTGr.SetMarkerColor(kRed)
249  x0profTGr.SetMarkerSize(0.8)
250  x0profTGr.Draw('pe,same')
251  canvasX0map.cd(6)
252  x0profPhiG4 = TProfile(dlname+'x0profPhiG4',dlname+'x0profPhiG4',numbinsphi,-math.pi,math.pi) # noqa: F841
253  x0profPhiTG = TProfile(dlname+'x0profPhiTG',dlname+'x0profPhiTG',numbinsphi,-math.pi,math.pi) # noqa: F841
254  g4DiscTree.Draw('PathInX0:DensedHitPhi>>'+dlname+'x0profPhiG4','','prof')
255  tgDiscTree.Draw('PathInX0:DensedHitPhi>>'+dlname+'x0profPhiTG','','prof')
256  x0profG4phi = gDirectory.Get(dlname+'x0profPhiG4')
257  x0profTGphi = gDirectory.Get(dlname+'x0profPhiTG')
258  # style
259  x0profTGphi.SetMarkerStyle(8)
260  x0profTGphi.SetMarkerSize(0.5)
261  x0profTGphi.SetMarkerColor(kRed)
262  x0profTGphi.SetLineColor(kRed)
263  x0profG4phi.GetXaxis().SetTitle('m [mm]')
264  x0profG4phi.GetYaxis().SetTitle('t [X_{0}]')
265  if 1.5*x0profG4phi.GetMaximum() > x0maxProf : x0maxProf = 1.5*x0profG4phi.GetMaximum()
266  x0profG4phi.GetYaxis().SetRangeUser(0.,1.5*x0profG4phi.GetMaximum())
267  x0profG4phi.Draw('hist')
268  x0profTGphi.SetMarkerColor(kRed)
269  x0profTGphi.SetMarkerSize(0.8)
270  x0profTGphi.Draw('pe,same')
271 
272  canvasX0map.cd(1)
273  try :
274  g4DiscTree.Draw('PathInX0>>'+dlname+'range')
275  x0range = gDirectory.Get(dlname+'range')
276  x0maxProf = 1.2*x0range.GetXaxis().GetXmax()
277  except TypeError :
278  print (' -> TTree does not exist in TrackingGeometry ')
279  x0mapG4 = TH3F(dlname+'x0mapG4',dlname+'x0mapG4',numbinsr,-1.1*outerR,1.1*outerR,numbinsr,-1.1*outerR,1.1*outerR,50,0.,x0maxProf)
280  try :
281  g4DiscTree.Draw('PathInX0:DensedHitY:DensedHitX>>'+dlname+'x0mapG4','','')
282  x0mapG4 = gDirectory.Get(dlname+'x0mapG4')
283  x0mapG4Prof = x0mapG4.Project3DProfile('xy')
284  x0mapG4Prof.GetXaxis().SetTitle('x [mm]')
285  x0mapG4Prof.GetYaxis().SetTitle('y [mm]')
286  x0mapG4Prof.GetZaxis().SetTitle('t [X_{0}]')
287  x0mapG4Prof.GetZaxis().SetRangeUser(0.,x0maxProf)
288  x0mapG4Prof.GetZaxis().SetTitleOffset(1.5*x0mapG4Prof.GetZaxis().GetTitleOffset())
289  x0mapG4Prof.SetBit(kCanDelete)
290  x0mapG4Prof.SetBit(kMustCleanup)
291  x0mapG4Prof.Draw('colz')
292  except TypeError :
293  print (' -> TTree does not exist in TrackingGeometry ')
294  canvasX0map.cd(2)
295  try :
296  x0mapTG = TH3F(dlname+'x0mapTG',dlname+'x0mapTG',numbinsr,-1.1*outerR,1.1*outerR,numbinsr,-1.1*outerR,1.1*outerR,50,0.,x0maxProf)
297  tgDiscTree.Draw('PathInX0:DensedHitY:DensedHitX>>'+dlname+'x0mapTG','','')
298  x0mapTG = gDirectory.Get(dlname+'x0mapTG')
299  x0mapTGProf = x0mapTG.Project3DProfile('xy')
300  x0mapTGProf.GetXaxis().SetTitle('x [mm]')
301  x0mapTGProf.GetYaxis().SetTitle('y [mm]')
302  x0mapTGProf.GetZaxis().SetTitle('t [X_{0}]')
303  x0mapTGProf.GetZaxis().SetRangeUser(0.,x0maxProf)
304  x0mapTGProf.GetZaxis().SetTitleOffset(1.5*x0mapTGProf.GetZaxis().GetTitleOffset())
305  x0mapTGProf.SetBit(kCanDelete)
306  x0mapTGProf.SetBit(kMustCleanup)
307  x0mapTGProf.Draw('colz')
308  except TypeError :
309  print (' -> TTree does not exist in Geant4 ')
310  canvasX0map.cd(3)
311  dlmap_x0.SetLineColor(kWhite)
312  dlmap_x0.GetZaxis().SetRangeUser(0.,1.1*x0max)
313  dlmap_x0.Draw('polcolz')
314  canvasX0map.cd(4)
315  try :
316  x0mapDiff = x0mapTGProf.Clone()
317  x0mapDiff.Add(x0mapG4Prof,-1)
318  x0mapDiff.GetXaxis().SetTitle('x [mm]')
319  x0mapDiff.GetYaxis().SetTitle('y [mm]')
320  x0mapDiff.GetZaxis().SetTitle('t [X_{0}]')
321  x0mapDiff.GetZaxis().SetRangeUser(-1.,1.)
322  x0mapDiff.GetZaxis().SetTitleOffset(1.5*x0mapDiff.GetZaxis().GetTitleOffset())
323  x0mapDiff.Divide(x0mapG4Prof)
324  x0mapDiff.Draw('colz')
325  except UnboundLocalError :
326  print (' -> No TrackingGeometry/G4 histograms to divide')
327 
328  canvasX0map.SaveAs(saveLocation+'/LayerMaterial_DiscLayer_'+str(evt.discLayerIndex)+'.png')
329 
330  # memory cleanup
331  gDirectory.Delete(dlname+'x0mapTG')
332  gDirectory.Delete(dlname+'x0mapG4')
python.TrigEgammaMonitorHelper.TH2F
def TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:45
python.TrigEgammaMonitorHelper.TProfile
def TProfile(*args, **kwargs)
Definition: TrigEgammaMonitorHelper.py:81
LayerMaterialUtils.drawCylinderLayer
def drawCylinderLayer(g4CylTree, cylLayer, tgCylTree, saveLocation)
Definition: LayerMaterialUtils.py:6
plotBeamSpotVxVal.range
range
Definition: plotBeamSpotVxVal.py:195
LayerMaterialUtils.drawDiscLayer
def drawDiscLayer(g4DiscTree, discLayer, tgDiscTree, saveLocation)
Definition: LayerMaterialUtils.py:178
str
Definition: BTagTrackIpAccessor.cxx:11
Generated on Sun Dec 22 2024 21:13:51 for ATLAS Offline Software by   doxygen 1.8.18