50 """Initialize the module"""
55 self.
rootfile = ROOT.TFile(
'cosmicAnalysis.root',
'recreate')
62 ROOT.gStyle.SetCanvasColor(0)
63 ROOT.gStyle.SetStatBorderSize(1)
64 ROOT.gStyle.SetStatColor(0)
65 ROOT.gStyle.SetTitleColor(1)
66 ROOT.gStyle.SetCanvasColor(0)
67 ROOT.gStyle.SetPadColor(0)
68 ROOT.gStyle.SetPadBorderMode(0)
69 ROOT.gStyle.SetCanvasBorderMode(0)
70 ROOT.gStyle.SetFrameBorderMode(0)
71 ROOT.gStyle.SetOptStat(0)
76 for key
in TrackData.__dict__:
79 if isinstance(self.
trackData.__getattribute__(key), int):
85 for key
in DEDXData.__dict__:
88 if isinstance(self.
dedxData.__getattribute__(key), int):
94 self.
HitsVsLayer = ROOT.TH2F(
'HitsVsLayer',
'', 6, 0.5, 6.5, 6, 0.5, 6.5)
96 self.
HitsVsSensor = ROOT.TH2F(
'HitsVsSensor',
'', 6, 0.5, 6.5, 5, 0.5, 5.5)
98 self.
LayerVsSensor = ROOT.TH2F(
'LayerVsSensor',
'', 6, 0.5, 6.5, 5, 0.5, 5.5)
102 self.
Chi2 = ROOT.TH1F(
'Chi2',
'', 300, 0.0, 500)
104 self.
NDF = ROOT.TH1F(
'NDF',
'', 200, 0.0, 200)
108 self.
Momentum = ROOT.TH1F(
'Momentum',
'', 500, 0.0, 1000)
113 'ADCCountOverNumberOfHitsInCDCVsMomentum',
'', 200, 0.0, 300, 0.0, 1000)
116 'MomentumVsADCCountOverNumberOfHitsInCDC',
'', 100, 0.0, 100, 0.0, 1000)
122 """ Fill histograms """
126 nCDCDedxTracks = cdcDedxTracks.getEntries()
127 for DEDXTrack_index
in range(nCDCDedxTracks):
128 dedxTrack = cdcDedxTracks[DEDXTrack_index]
129 self.dedxData.momentumCDC = dedxTrack.getMomentum()
130 self.dedxData.acdCount = 0
132 self.dedxData.path = 0
133 self.dedxData.numberOfHits = dedxTrack.size()
135 nKLMClusters = KLMClusters.getEntries()
136 if nKLMClusters == 2:
137 for nHit
in range(dedxTrack.size()):
138 self.dedxData.acdCount = self.dedxData.acdCount + dedxTrack.getADCCount(nHit)
139 self.dedxData.dE = self.dedxData.dE + dedxTrack.getDE(nHit)
140 self.dedxData.path = self.dedxData.path + dedxTrack.getPath(nHit)
141 self.dedxData.distanceKLM = math.sqrt(
143 KLMClusters[0].getClusterPosition().Z() -
144 KLMClusters[1].getClusterPosition().Z(),
147 self.ADCCountOverNumberOfHitsInCDC.Fill(self.dedxData.acdCount / dedxTrack.size())
148 self.MomentumVsADCCountOverNumberOfHitsInCDC.Fill(
149 self.dedxData.momentumCDC, self.dedxData.acdCount / dedxTrack.size())
150 self.ADCCountOverNumberOfHitsInCDCVsMomentum.Fill(
151 self.dedxData.acdCount / dedxTrack.size(), self.dedxData.momentumCDC)
154 self.tree_DEDX.Fill()
158 nRecoTracks = RecoTracks.getEntries()
160 for track_index
in range(nRecoTracks):
161 track = RecoTracks[track_index]
163 if track.wasFitSuccessful():
166 nTrackFitResults = TrackFitResults.getEntries()
167 if nTrackFitResults == 1:
168 self.trackData.momentum = TrackFitResults[0].getMomentum().Mag()
169 self.Momentum.Fill(TrackFitResults[0].getMomentum().Mag())
171 if track.hasPXDHits()
or track.hasSVDHits():
173 if track.getNumberOfSVDHits() % 2 == 0:
174 totalNumberOfHits = track.getNumberOfPXDHits() + track.getNumberOfSVDHits() / 2
176 totalNumberOfHits = track.getNumberOfPXDHits() + (track.getNumberOfSVDHits() - 1) / 2
178 self.TotalNumberOfHits.Fill(totalNumberOfHits)
181 self.Chi2.Fill(track.getTrackFitStatus().getChi2())
182 self.NDF.Fill(track.getTrackFitStatus().getNdf())
183 self.Chi2OverNDF.Fill(track.getTrackFitStatus().getChi2() / track.getTrackFitStatus().getNdf())
185 self.trackData.chi2 = track.getTrackFitStatus().getChi2()
186 self.trackData.ndf = track.getTrackFitStatus().getNdf()
187 self.trackData.chiSquaredOverNdf = track.getTrackFitStatus().getChi2() / track.getTrackFitStatus().getNdf()
190 self.tree_track.Fill()
193 if track.hasPXDHits():
202 for n
in range(0, len(track.getPXDHitList())):
203 pxdHit = track.getPXDHitList()[n]
205 info = geoCache.get(sensorID)
206 layer = sensorID.getLayerNumber()
207 sensor = sensorID.getSensorNumber()
209 self.HitsVsLayer.Fill(totalNumberOfHits, layer)
210 self.HitsVsSensor.Fill(totalNumberOfHits, sensor)
211 self.LayerVsSensor.Fill(layer, sensor)
212 self.PXDClusterSize.Fill(track.getPXDHitList()[n].getSize())
215 if track.hasSVDHits():
224 for n
in range(0, len(track.getSVDHitList())):
225 svdHit = track.getSVDHitList()[n]
227 info = geoCache.get(sensorID)
228 layer = sensorID.getLayerNumber()
229 sensor = sensorID.getSensorNumber()
231 if svdHit.isUCluster():
233 for lst
in range(0, len(track.getSVDHitList())):
234 svdHitA = track.getSVDHitList()[lst]
236 if info == geoCache.get(sensorIDA):
237 if svdHitA.isUCluster() == 0:
239 self.HitsVsLayer.Fill(totalNumberOfHits, layer)
240 self.HitsVsSensor.Fill(totalNumberOfHits, sensor)
241 self.LayerVsSensor.Fill(layer, sensor)
245 print(
'Fit was not successful')
248 """Close & write output files"""
251 labelHits = (
'',
'One',
'Two',
'Three',
'Four',
'Five',
'Six')
252 labelLayers = (
'',
'First',
'Second',
'Third',
'Fourth',
'Fifth',
'Sixth')
253 labelSensors = (
'',
'First',
'Second',
'Third',
'Fourth',
'Fifth')
255 for i
in range(1, 7):
270 for i
in range(1, 7):
271 self.
HitsVsLayer.GetXaxis().SetBinLabel(i, labelHits[i])
272 self.
HitsVsLayer.GetYaxis().SetBinLabel(i, labelLayers[i])
275 self.
HitsVsLayer.GetXaxis().SetTitle(
"Number hits in VXD")
281 ROOT.gPad.SetMargin(0.11, 0.11, 0.1, 0.1)
285 for i
in range(1, 7):
286 self.
HitsVsSensor.GetXaxis().SetBinLabel(i, labelHits[i])
287 for i
in range(1, 6):
288 self.
HitsVsSensor.GetYaxis().SetBinLabel(i, labelSensors[i])
291 self.
HitsVsSensor.GetXaxis().SetTitle(
"Number hits in VXD")
297 ROOT.gPad.SetMargin(0.11, 0.11, 0.1, 0.1)
301 for i
in range(1, 7):
303 for i
in range(1, 6):
304 self.
LayerVsSensor.GetYaxis().SetBinLabel(i, labelSensors[i])
313 ROOT.gPad.SetMargin(0.11, 0.11, 0.1, 0.1)
317 self.
Chi2.GetXaxis().SetTitle(
"#Chi^{2}")
318 self.
Chi2.GetYaxis().SetTitle(
"Number of tracks")
319 self.
Chi2.GetYaxis().CenterTitle()
320 self.
Chi2.GetXaxis().CenterTitle()
321 self.
Chi2.GetYaxis().SetTitleOffset(1.3)
322 self.
Chi2.GetXaxis().SetTitleOffset(1.3)
323 self.
Chi2.SetFillStyle(3365)
324 self.
Chi2.SetFillColor(9)
325 self.
Chi2.SetLineColor(9)
329 self.
NDF.GetXaxis().SetTitle(
"Degrees of freedom")
330 self.
NDF.GetYaxis().SetTitle(
"Number of tracks")
331 self.
NDF.GetYaxis().CenterTitle()
332 self.
NDF.GetXaxis().CenterTitle()
333 self.
NDF.GetYaxis().SetTitleOffset(1.3)
334 self.
NDF.GetXaxis().SetTitleOffset(1.3)
335 self.
NDF.SetFillStyle(3365)
336 self.
NDF.SetFillColor(9)
337 self.
NDF.SetLineColor(9)
341 self.
Chi2OverNDF.SetTitle(
"#Chi^{2}/Degrees of freedom")
342 self.
Chi2OverNDF.GetYaxis().SetTitle(
"Number of tracks")
353 self.
Momentum.GetXaxis().SetTitle(
"Momentum [GeVc^{-1}]")
354 self.
Momentum.GetYaxis().SetTitle(
"Number of tracks")
355 self.
Momentum.GetYaxis().CenterTitle()
356 self.
Momentum.GetXaxis().CenterTitle()
357 self.
Momentum.GetYaxis().SetTitleOffset(1.3)
358 self.
Momentum.GetXaxis().SetTitleOffset(1.3)
377 function = ROOT.TF1(
"function",
"[0]+[1]*log(x)")
378 function.SetParameters(0.0, 0.0)
389 inverse_function = ROOT.TF1(
"inverse_function",
"exp((-[0]+x)/[1])", 0, 115)
390 inverse_function.SetParameter(0, function.GetParameter(0))
391 inverse_function.SetParameter(1, function.GetParameter(1))
392 constant = ROOT.TF1(
"constant",
"[0]", 115, 300)
393 constant.SetParameter(0, 20)
405 print(
'Parameters of fitted function y = exp((-p0+x)/p1) + p2, where x ~ ADCCounts/NumberOfCDCHits, y ~ momentum')
406 print(
'p0 =', inverse_function.GetParameter(0),
'p1 =', inverse_function.GetParameter(1),
'p2 =', constant.GetParameter(0))