Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-10-08 05:12:16

0001 #!/usr/bin/env python3
0002 
0003 import ROOT
0004 from array import array
0005 
0006 from collections import OrderedDict
0007 import pickle
0008 
0009 from Validation.RecoTrack.python.plotting.ntuple import *
0010 
0011 import analysis as an
0012 import graphics as gr
0013 #import graphics_vertical as gr # Use this to plot helixes better, this, however, supports only vertical z-axis and no filters except fake filter.
0014 
0015 
0016 def main():
0017     '''
0018     The main function, below are commented templates of function calls to perform analysis.
0019     '''
0020 
0021     ntuple = TrackingNtuple("/afs/cern.ch/work/j/jmaanpaa/public/trackingNtuple_810pre6_v5.root")
0022     #ntuple = TrackingNtuple("trackingNtuple_pre9_1000ev.root") # Uncomment this for 1000 event data
0023     pl = gr.EventPlotter() 
0024 
0025     ###### USEFUL COMMANDS, UNCOMMENT TO USE ######
0026 
0027     PlotFakes(ntuple) # Plots the fakes in 3D one by one
0028 
0029     #Create_Histograms_MaxMatched_FakeTracks(ntuple)
0030     #Create_Histograms_MaxMatched_RealTracks(ntuple)
0031 
0032     #Create_Histograms_Classification(ntuple) # CREATES tHE MAIN CLASSIFICATION HISTOGRAMS
0033     #Create_Histograms_IncludedDecayHitFrac(ntuple)
0034     #Create_Histograms_UnmatchedInfo(ntuple)
0035     #Create_Histograms_PixInfo(ntuple)
0036 
0037     #Create_EndOfTrackingHistograms(ntuple)
0038     #Create_EndOfTrackingHistogramsReal(ntuple)
0039     #Create_EOT_DistanceHistogram(ntuple)
0040 
0041     #Create_ResolutionHistograms(ntuple)
0042 
0043     #print an.Calculate_AmountOfFakes(ntuple)
0044 
0045 
0046     ###### NOT THAT USEFUL COMMANDS, NOT TESTED AFTER UPDATES ######
0047      
0048     #ntuple_iter = iter(ntuple)#__iter__()
0049     #event = ntuple_iter.next()
0050     #tracks = event.tracks()
0051 
0052     #Create_PtHistograms(ntuple)
0053      
0054     #an.Save_Normalisation_Coefficients(ntuple)
0055     
0056     #Create_EndOfTrackingPointsPlot(ntuple)
0057     #an.EndOfTrackingDetectorInfo(Load_File("End_list_file_layermiss.dmp"),[0],True)
0058  
0059     #Create_Histograms_SharedHitsFracs_RealTracks(ntuple)
0060     #Create_Histograms_SharedHitsFracs(ntuple)  
0061 
0062     #PlotEvent3DHits(event, "PSG")
0063     #PlotXY(event, [-10,10])
0064     #PlotZY(event,flag="P")#, [-10,10])
0065     #PlotPixelGluedTracks3D(tracks)
0066     #PlotTracks3D(FindFakes(event))
0067     #PlotTracks3D(FindUntrackedParticles(event))
0068     #EfficiencyRate(ntuple, 100)
0069     #for track in an.FindFakes(event):
0070 #        pl.DrawTrackTest(track)
0071 
0072     #pl.PlotEvent3DHits(event)
0073 
0074     '''
0075     track = iter(tracks).next()
0076     pl.Plot3DHits(track)
0077     pl.Draw()
0078     
0079     for track in tracks:
0080         if track.hits():
0081             #pl.PlotTrack3D(track)
0082             pl.Plot3DHelixes([track])
0083             pl.Plot3DHits(track)
0084             #pl.Plot3DHelixes(an.FindFakes(event))
0085             pl.Draw()
0086     '''
0087     #pl.ParticleTest(event.trackingParticles(),True)
0088     #pl.PlotFakes(event,1,1,None)#[])
0089     #pl.PlotZY(event)
0090     #pl.Draw()
0091 
0092 
0093 def Save_To_File(items, name):
0094     f = open(name,'w')
0095     pickle.dump(items,f)
0096     f.close()
0097 
0098 def Load_File(name):
0099     f = open(name, 'r')
0100     items = pickle.load(f)
0101     f.close()
0102     return items
0103 
0104 def PlotFakes(ntuple):
0105     ''' Plots fakes with the PlotFakes command in the graphics module. '''
0106     pl = gr.EventPlotter()
0107     for event in ntuple:
0108         # Fill the filters with corresponding class indexes or detector layers to filter everything else from the plotting
0109         pl.PlotFakes(event, True, fake_filter = [], end_filter = [], last_filter = "", particle_end_filter = "")
0110 
0111 def Create_PtHistograms(ntuple):
0112     ''' Creates pt histograms for fake and true tracks. '''
0113     pl = gr.EventPlotter()
0114     
0115     fake_pts, true_pts = an.TrackPt(ntuple)
0116     fake_pts = [min(pt, 4.99) for pt in fake_pts]
0117     true_pts = [min(pt, 4.99) for pt in true_pts]
0118     pl.ClassHistogram(fake_pts, "Fake p_{t}", False, False)
0119     pl.ClassHistogram(true_pts, "True p_{t}", False, False)
0120 
0121 def Create_ResolutionHistograms(ntuple):
0122     '''
0123     Performs resolution analysis and plots histograms from the results.
0124     Comment and uncomment lines depending on if you have saved the files.
0125     '''
0126     pl = gr.EventPlotter()
0127 
0128     res_data_fakes = an.Resolution_Analysis_Fakes(ntuple, 100)
0129     #Save_To_File(res_data_fakes, "resolutions_fake_matches1000ev.dmp")
0130     #res_data_fakes = Load_File("resolutions_fake_matches1000ev.dmp")
0131 
0132     #pl.ResolutionHistograms(res_data_fakes)
0133     #pl.HistogramSigmaPt(res_data_fakes)
0134 
0135     res_data_trues = an.Resolution_Analysis_Trues(ntuple, 100)
0136     #Save_To_File(res_data_trues, "resolutions_trues1000ev.dmp")
0137     #res_data_trues = Load_File("resolutions_trues1000ev.dmp")
0138 
0139     res_data_fakes_NLay = an.Resolution_Analysis_Fakes(ntuple, 100, real_criterion = ["NLay",6])
0140     #Save_To_File(res_data_fakes_NLay, "resolutions_fake_NLay_matches1000ev.dmp")
0141     #res_data_fakes_NLay = Load_File("resolutions_fake_NLay_matches1000ev.dmp")
0142 
0143     #pl.ResolutionHistograms(res_data_fakes)
0144     #pl.HistogramSigmaPt(res_data_fakes)
0145 
0146     #pl.ResolutionHistograms(res_data_trues)
0147     #pl.HistogramSigmaPt(res_data_trues)
0148     
0149     #pl.ResolutionHistogramsCombine(res_data_trues, res_data_fakes)
0150     pl.ResolutionHistogramsNormalised(res_data_trues, res_data_fakes)
0151     pl.ResolutionHistogramsNormalisedCumulative(res_data_trues, res_data_fakes)
0152 
0153     #pl.HistogramSigmaPtCombined(res_data_trues, res_data_fakes)
0154     #pl.HistogramSigmaPtDual(res_data_trues, res_data_fakes)
0155     pl.ResolutionsSigmaVsPt(res_data_trues, res_data_fakes)
0156     
0157     pl.ResolutionHistogramsNormalised(res_data_trues, res_data_fakes_NLay)
0158     pl.ResolutionHistogramsNormalisedCumulative(res_data_trues, res_data_fakes_NLay)
0159 
0160     pl.ResolutionsSigmaVsPt(res_data_trues, res_data_fakes_NLay)
0161 
0162 
0163 
0164 def Create_EndOfTrackingHistograms(ntuple):
0165     '''
0166     Creates various histograms related to end of tracking
0167     '''
0168     pl = gr.EventPlotter()
0169 
0170     ''' Comment/uncomment these lines depending on if you have saved the end of tracking file'''
0171     end_list = an.Get_EndOfTrackingPoints(ntuple, 100, [])
0172     #Save_To_File(end_list,"End_list_file_layermiss_invalid4_fix_4class.dmp")
0173 
0174     #end_list = Load_File("End_list_file_layermiss_invalid4_fix_4class.dmp") #an.Get_EndOfTrackingPoints(ntuple, 10, [])
0175 
0176     #end_list = an.Get_EndOfTrackingPointsReal(ntuple, 100)
0177     #Save_To_File(end_list,"End_list_file_real_invalid.dmp")
0178 
0179     #end_list = Load_File("End_list_file_real_invalid.dmp")
0180     '''
0181     for end in end_list:
0182         if end.end_class == 0 and end.last_str == end.particle_end_str and end.last_str == "BPix3":
0183             print end.last
0184             print end.particle_end
0185             print end.fake_end
0186             print end.last_str, end.particle_end_str, end.fake_end_str
0187             print end.end_class
0188             print "*****"
0189     '''
0190     an.Analyse_EOT_ParticleDoubleHit(end_list)
0191     pl.EndOfTrackingHistogram(end_list, "end_class",[],[], False, False) #last false will adjust the bpix3 filter
0192     #pl.EndOfTrackingHistogram(end_list, "invalid_reason_last",[],[],False, True)
0193     pl.EndOfTrackingHistogram(end_list, "invalid_reason_fake",[0],[],False, False)
0194     pl.EndOfTrackingHistogram(end_list, "particle_pdgId",[0],[],False, False)
0195     #pl.EndOfTrackingHistogram(end_list, "invalid_reason_particle",[],[],False, False)
0196     #pl.EndOfTrackingHistogram(end_list, "particle_end_layer",[0],[],False, True)
0197 
0198     pl.EndOfTrackingHistogram(end_list, "last_layer",[0],[],False, False)# 11,13,21
0199     pl.EndOfTrackingHistogram(end_list, "fake_end_layer",[0],[],False, False)
0200     pl.EndOfTrackingHistogram(end_list, "particle_end_layer",[0],[],False, False)
0201 
0202     errors = an.Analyse_EOT_Error(end_list, True, "same", False)  # Same detectior id
0203     pl.ClassHistogram(errors,"Distance between EOT fake and particle hit", False, False)
0204     errors = an.Analyse_EOT_Error(end_list, True, "different", False)  # Different detectior id
0205     pl.ClassHistogram(errors,"Distance between EOT fake and particle hit", False, False)
0206     an.EndOfTrackingDetectorInfo(end_list,[0],True)
0207 
0208 
0209 def Create_EndOfTrackingHistogramsReal(ntuple):
0210     '''
0211     Create end of tracking histograms for true tracks.
0212     '''
0213     pl = gr.EventPlotter()
0214 
0215     ''' Comment/uncomment these lines depending on if you have saved the end of tracking file'''
0216     end_list = an.Get_EndOfTrackingPointsReal(ntuple, 100)
0217     #Save_To_File(end_list,"End_list_file_real_fixed.dmp")
0218 
0219     #end_list = Load_File("End_list_file_real_fixed.dmp") #an.Get_EndOfTrackingPoints(ntuple, 10, [])
0220 
0221     '''
0222     for end in end_list:
0223         if end.end_class == 0 and end.last_str == end.particle_end_str and end.last_str == "BPix3":
0224             print end.last
0225             print end.particle_end
0226             print end.fake_end
0227             print end.last_str, end.particle_end_str, end.fake_end_str
0228             print end.end_class
0229             print "*****"
0230     '''
0231     pl.EndOfTrackingHistogram(end_list, "end_class",[],[])
0232     pl.EndOfTrackingHistogram(end_list, "last_layer",[0],[],False)
0233     pl.EndOfTrackingHistogram(end_list, "fake_end_layer",[0],[],False)
0234     pl.EndOfTrackingHistogram(end_list, "particle_end_layer",[0],[],False)
0235 
0236 def Create_EOT_DistanceHistogram(ntuple):
0237     pl = gr.EventPlotter()
0238     end_list = an.Get_EndOfTrackingPoints(ntuple, 100, [])
0239     #Save_To_File(end_list,"End_list_file_layermiss.dmp")
0240     #end_list = Load_File("End_list_file_layermiss.dmp")
0241 
0242     errors = an.Analyse_EOT_Error(end_list, True)
0243     
0244     pl.ClassHistogram(errors,"Distance between EOT fake and particle hit", False, False)
0245     
0246 def Create_EndOfTrackingPointsPlot(ntuple):
0247     # NOT TESTED AFTER UPDATE
0248     pl = gr.EventPlotter()
0249 
0250     end_list = an.Get_EndOfTrackingPoints(ntuple, 100, [])
0251     last = [], fake_ends = [], particle_ends = [], end_classes = [], fake_classes = []
0252     for end in end_list:
0253         last.append(end.last)
0254         fake_ends.append(end.fake_end)
0255         particle_ends.append(end.particle_end)
0256         end_classes.append(end.end_class)
0257         fake_classes.append(end.fake_class)
0258 
0259     '''
0260     pl.PlotEndOfTrackingPoints3D(last,[])#, fail)
0261     pl.PlotEndOfTrackingPoints3D(fake_ends,[],[],2)
0262     pl.PlotEndOfTrackingPoints3D(particle_ends,[],[],3)
0263     #pl.PlotDetectorRange()
0264     pl.Draw()
0265     '''
0266     pl.PlotEndOfTracking3D(last, fake_ends, particle_ends, end_classes, fake_classes)
0267     pl.PlotDetectorRange()
0268     pl.Draw()
0269 
0270     pl.Reset()
0271     pl.PlotEndOfTracking3D(last, fake_ends, particle_ends, end_classes, fake_classes, [0],[])
0272     pl.PlotDetectorRange()
0273     pl.Draw()
0274 
0275     pl.Reset()
0276     pl.PlotEndOfTracking3D(last, fake_ends, particle_ends, end_classes, fake_classes, [],[10])
0277     pl.PlotDetectorRange()
0278     pl.Draw()
0279 
0280     pl.Reset()
0281     pl.PlotEndOfTracking3D(last, fake_ends, particle_ends, end_classes, fake_classes, [],[13])
0282     pl.PlotDetectorRange()
0283     pl.Draw()
0284 
0285 def Create_Histograms_IncludedDecayHitFrac(ntuple):
0286     '''
0287     Creates plots related to the "multiple matches including a decay interaction" -class.
0288     '''
0289     pl = gr.EventPlotter()
0290     
0291     daughter_frac, parent_frac, total_frac = an.Calculate_IndludedDecayHitFractions(ntuple, 100)
0292     pl.ClassHistogram(daughter_frac, "Daughter share of hits from fake", False, False)
0293     pl.ClassHistogram(parent_frac, "Parent share of hits from fake", False, False)
0294     pl.ClassHistogram(total_frac, "Total share of hits from fake", False, False)
0295 
0296 
0297 def Create_Histograms_UnmatchedInfo(ntuple):
0298     pl = gr.EventPlotter()
0299     
0300     results = an.Calculate_UnmatchedSharedHits(ntuple, 100)
0301     pl.ClassHistogram(results, "Maximum number of shared hits in unmatched and unclassified classes", False, False)
0302 
0303 def Create_Histograms_SharedHitsFracs(ntuple):
0304     pl = gr.EventPlotter()
0305     
0306     hits_shared, hit_fractions = an.Calculate_SharedHits(ntuple, 100, []) 
0307     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in all classes", False, False)
0308     pl.ClassHistogram(hit_fractions, "Shared hit fractions in all classes", False, False)
0309 
0310     hits_shared, hit_fractions = an.Calculate_SharedHits(ntuple, 100, [10,11,12,13]) 
0311     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in single match classes", False, False)
0312     pl.ClassHistogram(hit_fractions, "Shared hit fractions in single match classes", False, False)
0313 
0314     hits_shared, hit_fractions = an.Calculate_SharedHits(ntuple, 100, [20,21,22,23]) 
0315     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in multiple match classes", False, False)
0316     pl.ClassHistogram(hit_fractions, "Shared hit fractions in multiple match classes", False, False)
0317 
0318     hits_shared, hit_fractions = an.Calculate_SharedHits(ntuple, 100, [10]) 
0319     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in other single match class", False, False)
0320     pl.ClassHistogram(hit_fractions, "Shared hit fractions in other single match class", False, False)
0321 
0322     hits_shared, hit_fractions = an.Calculate_SharedHits(ntuple, 100, [13]) 
0323     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in single decay match class", False, False)
0324     pl.ClassHistogram(hit_fractions, "Shared hit fractions in single decay match class", False, False)
0325 
0326 def Create_Histograms_SharedHitsFracs_RealTracks(ntuple):
0327     pl = gr.EventPlotter()
0328     
0329     hits_shared, hit_fractions = an.Calculate_SharedHits_RealTracks(ntuple, 100) 
0330     pl.ClassHistogram(hits_shared, "Maximum number of shared hits in real reconstructed tracks", False, False)
0331     pl.ClassHistogram(hit_fractions, "Shared hit fractions in real reconstructed tracks", False, False)
0332 
0333 def Create_Histograms_MaxMatched_FakeTracks(ntuple):
0334     pl = gr.EventPlotter()
0335     #results = an.Calculate_MaxMatchedHits(ntuple, 100)
0336     #pl.ClassHistogram(results, "Max Matched Hits", False, False) 
0337      
0338     results = an.Calculate_MaxMatchedConsecutiveHits(ntuple, 100) 
0339     pl.ClassHistogram(results, "Matched consecutive hits of fake", False, False, True, 0, 20, 20, True)
0340 
0341     results = an.Calculate_MaxMatchedConsecutiveHits(ntuple, 100, True) 
0342     pl.ClassHistogram(results, "Consecutive hit fraction of fake", False, False, False, 0, 1.1, 11, True)
0343 
0344 def Create_Histograms_MaxMatched_RealTracks(ntuple):
0345     pl = gr.EventPlotter()
0346     results = an.Calculate_MaxMatchedHits_RealTracks(ntuple, 100)
0347     pl.ClassHistogram(results, "Matched Hits of real reconstructed tracks", False, False) 
0348     
0349     results = an.Calculate_MaxMatchedConsecutiveHits_RealTracks(ntuple, 100)
0350     pl.ClassHistogram(results, "Matched consecutive hits of real reconstructed track", False, False)
0351 
0352 def Create_Histograms_PixInfo(ntuple):
0353     pl = gr.EventPlotter()
0354     
0355     nPix_data, nLay_data = an.Calculate_MatchPixelHits(ntuple, 100)
0356     pl.ClassHistogram(nPix_data, "Number of pixelhits in matches", False, False)
0357     pl.ClassHistogram(nLay_data, "Number of pixel layers in matches", False, False)
0358 
0359 def Create_Histograms_Classification(ntuple):
0360     '''
0361     Creates histograms from classification of fakes.
0362     '''
0363     pl = gr.EventPlotter()
0364     '''
0365     results = an.ClassifyEventFakes(ntuple, 10, False)
0366     pl.ClassHistogram(results, "Classification of all fakes", an.classes, True) 
0367     
0368     results = an.Calculate_MaxMatchedConsecutiveHits(ntuple, 100)
0369     pl.ClassHistogram(results, "Max matched consecutive hits", False, False) 
0370     
0371     results = an.Calculate_MaxMatchedHits(ntuple, 100)
0372     pl.ClassHistogram(results, "Max Matched Hits", False, False) 
0373      
0374     results = an.Calculate_MaxMatchedConsecutiveHits(ntuple, 100, True)
0375     pl.ClassHistogram(results, "Hit fraction of fake", False, False) #"Max matched consecutive hits with respect to fake", False, False)
0376     '''
0377  
0378     results = an.ClassifyEventFakes(ntuple, 100, False) # ORIGINAL MATCH CRITERION
0379     #results = an.ClassifyEventFakes(ntuple, 100, False, real_criterion = ["NLay", 6]) # SHARED LAYERS (NLay) MATCH CRITERION
0380 
0381     #pl.ClassHistogram(results, "Classification of all fakes", an.classes, True) 
0382 
0383     all_classified = OrderedDict()
0384     all_classified["No matches"] = results[0]
0385     all_classified["One match"] = sum([results[i] for i in [10,11,12,13]])
0386     all_classified["Multiple matches"] = sum([results[i] for i in [20,21,22,23]])
0387     all_classified["Loose match"] = results[-1]
0388 
0389     pl.ClassHistogram(all_classified, "Classification of all fakes", False, True)
0390 
0391     single_classified = OrderedDict()
0392     single_classified["Decayed and reconstructed"] = results[11]
0393     single_classified["Only reconstructed"] = results[12]
0394     single_classified["Only decayed"] = results[13]
0395     single_classified["Neither"] = results[10]
0396     pl.ClassHistogram(single_classified, "Classification of fakes with single match", False, True)
0397 
0398     multi_classified = OrderedDict()
0399     multi_classified["Includes a decay interaction"] = results[21]
0400     multi_classified["All matches reconstructed"] = results[22]
0401     multi_classified["Some matches reconstructed"] = results[23]
0402     multi_classified["Other"] = results[20]
0403     pl.ClassHistogram(multi_classified, "Classification of fakes with multiple matches", False, True)
0404 
0405 
0406 
0407 if __name__ == "__main__":
0408     main()
0409