diff --git a/Simulation/Tools/HitAnalysis/CMakeLists.txt b/Simulation/Tools/HitAnalysis/CMakeLists.txt index 4c38d7c1c047a4da322fb7613859b2236ffb0bed..adacc0f9cafdb7c0670c24fe00fff63b5e8cddae 100644 --- a/Simulation/Tools/HitAnalysis/CMakeLists.txt +++ b/Simulation/Tools/HitAnalysis/CMakeLists.txt @@ -7,20 +7,29 @@ atlas_subdir( HitAnalysis ) # Declare the package's dependencies: atlas_depends_on_subdirs( PUBLIC - Control/AthenaBaseComps - Control/StoreGate GaudiKernel PRIVATE Calorimeter/CaloDetDescr Calorimeter/CaloIdentifier + Calorimeter/CaloSimEvent + Control/AthenaBaseComps DetectorDescription/GeoModel/GeoAdaptors Event/EventInfo + ForwardDetectors/AFP/AFP_SimEv + ForwardDetectors/ALFA/ALFA_SimEv + ForwardDetectors/LUCID/LUCID_SimUtils/LUCID_SimEvent + ForwardDetectors/ZDC/ZDC_SimEvent + Generators/GeneratorObjects + InnerDetector/InDetSimEvent LArCalorimeter/LArSimEvent + MuonSpectrometer/MuonSimEvent + Simulation/G4Sim/TrackRecord TileCalorimeter/TileDetDescr TileCalorimeter/TileSimEvent ) # External dependencies: find_package( CLHEP ) +find_package( HepMC ) find_package( ROOT COMPONENTS Core Tree MathCore Hist RIO pthread Table MathMore Minuit Minuit2 Matrix Physics HistPainter Rint Graf Graf3d Gpad Html Postscript Gui GX11TTF GX11 ) # tag ROOTBasicLibs was not recognized in automatic conversion in cmt2cmake @@ -31,8 +40,8 @@ find_package( ROOT COMPONENTS Core Tree MathCore Hist RIO pthread Table MathMore atlas_add_component( HitAnalysis src/*.cxx src/components/*.cxx - INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${CLHEP_INCLUDE_DIRS} - LINK_LIBRARIES ${ROOT_LIBRARIES} ${CLHEP_LIBRARIES} AthenaBaseComps StoreGateLib SGtests GaudiKernel CaloDetDescrLib CaloIdentifier GeoAdaptors EventInfo LArSimEvent TileDetDescr TileSimEvent ) + INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${CLHEP_INCLUDE_DIRS} ${HEPMC_INCLUDE_DIRS} + LINK_LIBRARIES ${ROOT_LIBRARIES} ${CLHEP_LIBRARIES} ${HEPMC_LIBRARIES} GaudiKernel CaloDetDescrLib CaloIdentifier CaloSimEvent AthenaBaseComps GeoAdaptors EventInfo AFP_SimEv ALFA_SimEv LUCID_SimEvent ZDC_SimEvent GeneratorObjects InDetSimEvent LArSimEvent MuonSimEvent TileDetDescr TileSimEvent ) # Install files from the package: atlas_install_headers( HitAnalysis ) diff --git a/Simulation/Tools/HitAnalysis/HitAnalysis/CaloHitAnalysis.h b/Simulation/Tools/HitAnalysis/HitAnalysis/CaloHitAnalysis.h deleted file mode 100755 index 946093b3948c937c81a75bc65850f4257723027b..0000000000000000000000000000000000000000 --- a/Simulation/Tools/HitAnalysis/HitAnalysis/CaloHitAnalysis.h +++ /dev/null @@ -1,59 +0,0 @@ -/* - Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration -*/ - -#ifndef CALO_HIT_ANALYSIS_H -#define CALO_HIT_ANALYSIS_H - -#include "GaudiKernel/ToolHandle.h" -#include "GaudiKernel/Algorithm.h" -#include "GaudiKernel/ObjectVector.h" -#include "CLHEP/Units/SystemOfUnits.h" -#include "StoreGate/StoreGateSvc.h" - -#include "AthenaBaseComps/AthAlgorithm.h" - -#include <string> - -#include "TH1.h" - -class TileID; -class TileDetDescrManager; -class TTree; -class ITHistSvc; - -class CaloHitAnalysis : public AthAlgorithm { - - public: - - CaloHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); - ~CaloHitAnalysis(); - - virtual StatusCode initialize(); - virtual StatusCode finalize(); - virtual StatusCode execute(); - - private: - - /** a handle on Store Gate for access to the Event Store */ - StoreGateSvc* m_storeGate; - - /** Simple variables by Ketevi */ - std::vector<float>* m_cell_eta; - std::vector<float>* m_cell_phi; - std::vector<float>* m_cell_e; - std::vector<float>* m_cell_radius; - - TTree * m_tree; - std::string m_ntupleFileName; - std::string m_ntupleDirName; - std::string m_ntupleTreeName; - - ITHistSvc * m_thistSvc; - - const TileID * m_tileID; - const TileDetDescrManager * m_tileMgr; -}; - -#endif // CALO_HIT_ANALYSIS_H - diff --git a/Simulation/Tools/HitAnalysis/cmt/requirements b/Simulation/Tools/HitAnalysis/cmt/requirements index df7ce418dc34978a5164af242a057dbd6cdb8a03..1d8154fcff756ad12a02f8090a24ad1f80259391 100755 --- a/Simulation/Tools/HitAnalysis/cmt/requirements +++ b/Simulation/Tools/HitAnalysis/cmt/requirements @@ -7,21 +7,29 @@ use AtlasPolicy AtlasPolicy-* branches run use GaudiInterface GaudiInterface-* External -use AtlasCLHEP AtlasCLHEP-* External use AtlasROOT AtlasROOT-* External -use StoreGate StoreGate-* Control - -use AthenaBaseComps AthenaBaseComps-* Control - private -use TileSimEvent TileSimEvent-* TileCalorimeter -use TileDetDescr TileDetDescr-* TileCalorimeter -use LArSimEvent LArSimEvent-* LArCalorimeter -use GeoAdaptors GeoAdaptors-* DetectorDescription/GeoModel -use CaloIdentifier CaloIdentifier-* Calorimeter -use CaloDetDescr CaloDetDescr-* Calorimeter -use EventInfo EventInfo-* Event +use AthenaBaseComps AthenaBaseComps-* Control +use InDetSimEvent InDetSimEvent-* InnerDetector +use TileSimEvent TileSimEvent-* TileCalorimeter +use TileDetDescr TileDetDescr-* TileCalorimeter +use LArSimEvent LArSimEvent-* LArCalorimeter +use GeoAdaptors GeoAdaptors-* DetectorDescription/GeoModel +use CaloIdentifier CaloIdentifier-* Calorimeter +use CaloDetDescr CaloDetDescr-* Calorimeter +use CaloSimEvent CaloSimEvent-* Calorimeter +use InDetSimEvent InDetSimEvent-* InnerDetector +use MuonSimEvent MuonSimEvent-* MuonSpectrometer +use AtlasCLHEP AtlasCLHEP-* External +use EventInfo EventInfo-* Event +use AtlasHepMC AtlasHepMC-* External +use GeneratorObjects GeneratorObjects-* Generators +use TrackRecord TrackRecord-* Simulation/G4Sim +use ALFA_SimEv ALFA_SimEv-* ForwardDetectors/ALFA +use LUCID_SimEvent LUCID_SimEvent-* ForwardDetectors/LUCID/LUCID_SimUtils +use ZDC_SimEvent ZDC_SimEvent-* ForwardDetectors/ZDC +use AFP_SimEv AFP_SimEv-* ForwardDetectors/AFP end_private use AtlasSimulationRunTime AtlasSimulationRunTime-* diff --git a/Simulation/Tools/HitAnalysis/share/AFPHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/AFPHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..2afbcca4cdb77defa0e2bd4de9df3fd27e423b30 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/AFPHitAnalysis_topOptions.py @@ -0,0 +1,44 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import AFPHitAnalysis +topSequence += AFPHitAnalysis() +AFPHitAnalysis = AFPHitAnalysis() +AFPHitAnalysis.NtupleFileName = '/AFPHitAnalysis/ntuples/' +AFPHitAnalysis.HistPath = '/AFPHitAnalysis/histos/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "AFPHitAnalysis DATAFILE='AFPHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/ALFAHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/ALFAHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..05ece5efe46658e466076b68e58401a4c7622553 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/ALFAHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import ALFAHitAnalysis +topSequence += ALFAHitAnalysis() +ALFAHitAnalysis = ALFAHitAnalysis() +ALFAHitAnalysis.HistPath = '/ALFAHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "ALFAHitAnalysis DATAFILE='ALFAHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/AllHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/AllHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..54b3b85785444b952715d591f944b1a53fb96924 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/AllHitAnalysis_topOptions.py @@ -0,0 +1,104 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/afs/cern.ch/work/a/ancuetog/HitsAreas/20.3.X.Y-VAL-r3/atlasG4_10ttbar_20.3.X.Y-VAL-r3-modtime.NoFrozenShower.DeadLAr.hits.pool.root" ) +#athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() +from HitAnalysis.HitAnalysisConf import * + +topSequence += SiHitAnalysis('PixelHitAnalysis') +topSequence.PixelHitAnalysis.CollectionName='PixelHits' +topSequence += SiHitAnalysis('SCTHitAnalysis') +topSequence.SCTHitAnalysis.CollectionName='SCT_Hits' +topSequence += SiHitAnalysis('BCMHitAnalysis') +topSequence.BCMHitAnalysis.CollectionName='BCMHits' +topSequence += SiHitAnalysis('BLMHitAnalysis') +topSequence.BLMHitAnalysis.CollectionName='BLMHits' +topSequence += TRTHitAnalysis('TRTHitAnalysis') +topSequence += RPCHitAnalysis('RPCHitAnalysis') +topSequence += MDTHitAnalysis('MDTHitAnalysis') +topSequence += CSCHitAnalysis('CSCHitAnalysis') +topSequence += TGCHitAnalysis('TGCHitAnalysis') +#topSequence += ALFAHitAnalysis('ALFAHitAnalysis') +#topSequence += LucidHitAnalysis('LucidHitAnalysis') +#topSequence += ZDCHitAnalysis('ZDCHitAnalysis') +topSequence += TrackRecordAnalysis('TrackRecordAnalysis') +topSequence += TruthHitAnalysis('TruthHitAnalysis') +topSequence += CaloHitAnalysis('CaloHitAnalysis') + + + +topSequence.PixelHitAnalysis.HistPath='/HitAnalysis/Pixel/histos/' +topSequence.SCTHitAnalysis.HistPath='/HitAnalysis/SCT/histos/' +topSequence.BCMHitAnalysis.HistPath='/HitAnalysis/BCM/histos/' +topSequence.BLMHitAnalysis.HistPath='/HitAnalysis/BLM/histos/' +topSequence.TRTHitAnalysis.HistPath='/HitAnalysis/TRT/histos/' +topSequence.RPCHitAnalysis.HistPath='/HitAnalysis/RPC/histos/' +topSequence.MDTHitAnalysis.HistPath='/HitAnalysis/MDT/histos/' +topSequence.CSCHitAnalysis.HistPath='/HitAnalysis/CSC/histos/' +topSequence.TGCHitAnalysis.HistPath='/HitAnalysis/TGC/histos/' +#topSequence.ALFAHitAnalysis.HistPath='/HitAnalysis/' +#topSequence.LucidHitAnalysis.HistPath='/HitAnalysis/' +#topSequence.ZDCHitAnalysis.HistPath='/HitAnalysis/' +topSequence.TrackRecordAnalysis.HistPath='/HitAnalysis/Track/histos/' +topSequence.TruthHitAnalysis.HistPath = '/HitAnalysis/Truth/histos/' +topSequence.CaloHitAnalysis.HistPath = '/HitAnalysis/Calo/histos/' +topSequence.PixelHitAnalysis.NtupleFileName='/HitAnalysis/Pixel/ntuple/' +topSequence.SCTHitAnalysis.NtupleFileName='/HitAnalysis/SCT/ntuple/' +topSequence.BCMHitAnalysis.NtupleFileName='/HitAnalysis/BCM/ntuple/' +topSequence.BLMHitAnalysis.NtupleFileName='/HitAnalysis/BLM/ntuple/' +topSequence.TRTHitAnalysis.NtupleFileName='/HitAnalysis/TRT/ntuple/' +topSequence.RPCHitAnalysis.NtupleFileName='/HitAnalysis/RPC/ntuple/' +topSequence.MDTHitAnalysis.NtupleFileName='/HitAnalysis/MDT/ntuple/' +topSequence.CSCHitAnalysis.NtupleFileName='/HitAnalysis/CSC/ntuple/' +topSequence.TGCHitAnalysis.NtupleFileName='/HitAnalysis/TGC/ntuple/' +#topSequence.ALFAHitAnalysis.NtupleFileName='/HitAnalysis/' +#topSequence.LucidHitAnalysis.NtupleFileName='/HitAnalysis/' +#topSequence.ZDCHitAnalysis.NtupleFileName='/HitAnalysis/' +topSequence.TrackRecordAnalysis.NtupleFileName='/HitAnalysis/Track/ntuple/' +topSequence.TruthHitAnalysis.NtupleFileName = '/HitAnalysis/Truth/ntuple/' +topSequence.CaloHitAnalysis.NtupleFileName = '/HitAnalysis/Calo/ntuple/' +#Add some more TH2 histograms + + +topSequence.PixelHitAnalysis.ExpertMode= "off" +topSequence.SCTHitAnalysis.ExpertMode= "off" +topSequence.BCMHitAnalysis.ExpertMode= "off" +topSequence.BLMHitAnalysis.ExpertMode= "off" +topSequence.CaloHitAnalysis.ExpertMode = "off" +topSequence.CaloHitAnalysis.CalibHits = "off" + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += ["HitAnalysis DATAFILE='AllHitAnalysis.root' OPT='RECREATE'" ] + + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = TRUE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include( "TrkDetDescrSvc/AtlasTrackingGeometrySvc.py" ) # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/CSCHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/CSCHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..3cc733d76a21145f721bd06370c054a5b4ca8ce1 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/CSCHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +q#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import CSCHitAnalysis +topSequence += CSCHitAnalysis() +CSCHitAnalysis = CSCHitAnalysis() +CSCHitAnalysis.HistPath = '/CSCHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "CSCHitAnalysis DATAFILE='CSCHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/CaloHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/CaloHitAnalysis_topOptions.py index 6c11ce957d6c01658de90a88560e7a02729dbb5a..95306b508b4560a8e68ed3738c9a6459dc8a4ba3 100755 --- a/Simulation/Tools/HitAnalysis/share/CaloHitAnalysis_topOptions.py +++ b/Simulation/Tools/HitAnalysis/share/CaloHitAnalysis_topOptions.py @@ -5,7 +5,7 @@ import AthenaPoolCnvSvc.ReadAthenaPool from PartPropSvc.PartPropSvcConf import PartPropSvc include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") -include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) +#include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) import os from glob import glob @@ -17,13 +17,19 @@ from AthenaCommon.AlgSequence import AlgSequence topSequence = AlgSequence() from HitAnalysis.HitAnalysisConf import CaloHitAnalysis -topSequence += CaloHitAnalysis() CaloHitAnalysis = CaloHitAnalysis() -CaloHitAnalysis.NtupleFileName = 'CaloHitAnalysis' +topSequence += CaloHitAnalysis +CaloHitAnalysis.HistPath = '/CaloHitAnalysis/histos/' +CaloHitAnalysis.NtupleFileName = '/CaloHitAnalysis/ntuple/' +#ExpertMode adds more histograms to the output. Default mode is off +CaloHitAnalysis.ExpertMode = "off" +#CalibHits adds Calibrated hits histograms to the output. Default mode is off +CaloHitAnalysis.CalibHits = "off" from GaudiSvc.GaudiSvcConf import THistSvc ServiceMgr += THistSvc() -ServiceMgr.THistSvc.Output += [ "CaloHitAnalysis DATAFILE='CaloHitAnalysis.root' OPT='RECREATE'" ] +ServiceMgr.THistSvc.Output += ["CaloHitAnalysis DATAFILE='CaloHitAnalysis.root' OPT='RECREATE'" ] + ServiceMgr.MessageSvc.OutputLevel = INFO ServiceMgr.MessageSvc.defaultLimit = 9999999 diff --git a/Simulation/Tools/HitAnalysis/share/HitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/HitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..2719aae9383ac07bcffb4fe2bfe496fb8db01a3a --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/HitAnalysis_topOptions.py @@ -0,0 +1,83 @@ +#from AthenaCommon.AppMgr import ToolSvc +#from AthenaCommon.AppMgr import ServiceMgr +#import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() +from HitAnalysis.HitAnalysisConf import * + +topSequence += SiHitAnalysis('PixelHitAnalysis') +topSequence.PixelHitAnalysis.CollectionName='PixelHits' +topSequence += SiHitAnalysis('SCTHitAnalysis') +topSequence.SCTHitAnalysis.CollectionName='SCT_Hits' +topSequence += SiHitAnalysis('BCMHitAnalysis') +topSequence.BCMHitAnalysis.CollectionName='BCMHits' +topSequence += SiHitAnalysis('BLMHitAnalysis') +topSequence.BLMHitAnalysis.CollectionName='BLMHits' +topSequence += TRTHitAnalysis('TRTHitAnalysis') +topSequence += RPCHitAnalysis('RPCHitAnalysis') +topSequence += MDTHitAnalysis('MDTHitAnalysis') +topSequence += CSCHitAnalysis('CSCHitAnalysis') +topSequence += TGCHitAnalysis('TGCHitAnalysis') +#topSequence += ALFAHitAnalysis('ALFAHitAnalysis') +#topSequence += LucidHitAnalysis('LucidHitAnalysis') +#topSequence += ZDCHitAnalysis('ZDCHitAnalysis') +topSequence += TrackRecordAnalysis('TrackRecordAnalysis') +topSequence += TruthHitAnalysis('TruthHitAnalysis') +topSequence += CaloHitAnalysis('CaloHitAnalysis') + + + +topSequence.PixelHitAnalysis.HistPath='/HitAnalysis/' +topSequence.SCTHitAnalysis.HistPath='/HitAnalysis/' +topSequence.BCMHitAnalysis.HistPath='/HitAnalysis/' +topSequence.BLMHitAnalysis.HistPath='/HitAnalysis/' +topSequence.TRTHitAnalysis.HistPath='/HitAnalysis/' +topSequence.RPCHitAnalysis.HistPath='/HitAnalysis/' +topSequence.MDTHitAnalysis.HistPath='/HitAnalysis/' +topSequence.CSCHitAnalysis.HistPath='/HitAnalysis/' +topSequence.TGCHitAnalysis.HistPath='/HitAnalysis/' +#topSequence.ALFAHitAnalysis.HistPath='/HitAnalysis/' +#topSequence.LucidHitAnalysis.HistPath='/HitAnalysis/' +#topSequence.ZDCHitAnalysis.HistPath='/HitAnalysis/' +topSequence.TrackRecordAnalysis.HistPath='/HitAnalysis/' +topSequence.TruthHitAnalysis.HistPath = '/HitAnalysis/' +topSequence.CaloHitAnalysis.HistPath = '/HitAnalysis/' +#Add some more TH2 histograms + + +topSequence.PixelHitAnalysis.ExpertMode= "off" +topSequence.SCTHitAnalysis.ExpertMode= "off" +topSequence.BCMHitAnalysis.ExpertMode= "off" +topSequence.BLMHitAnalysis.ExpertMode= "off" +topSequence.CaloHitAnalysis.ExpertMode = "off" +topSequence.CaloHitAnalysis.CalibHits = "off" + +#from GaudiSvc.GaudiSvcConf import THistSvc +#ServiceMgr += THistSvc() +#ServiceMgr.THistSvc.Output += ["CaloHitAnalysis DATAFILE='CaloHitAnalysis_v2.root' OPT='RECREATE'" ] + + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +#ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = TRUE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include( "TrkDetDescrSvc/AtlasTrackingGeometrySvc.py" ) # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/LucidHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/LucidHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..b4404d517df75c2370517461d5c8b9b189b22fb9 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/LucidHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import LucidHitAnalysis +topSequence += LucidHitAnalysis() +LucidHitAnalysis = LucidHitAnalysis() +LucidHitAnalysis.HistPath = '/LucidHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "LucidHitAnalysis DATAFILE='LucidHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/MDTHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/MDTHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..8cac6fc0153a0037b429515aef9e68ce77ff7977 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/MDTHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import MDTHitAnalysis +topSequence += MDTHitAnalysis() +MDTHitAnalysis = MDTHitAnalysis() +MDTHitAnalysis.NtupleFileName = 'MDTHitAnalysis' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "MDTHitAnalysis DATAFILE='MDTHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/MuonHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/MuonHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..1d887bb4d46208b112e8ac80ab6265457a4b1bb5 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/MuonHitAnalysis_topOptions.py @@ -0,0 +1,51 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import RPCHitAnalysis +topSequence += RPCHitAnalysis('RPCHitAnalysis') +topSequence.RPCHitAnalysis.HistPath = '/MuonHitAnalysis/' +from HitAnalysis.HitAnalysisConf import MDTHitAnalysis +topSequence += MDTHitAnalysis('MDTHitAnalysis') +topSequence.MDTHitAnalysis.HistPath = '/MuonHitAnalysis/' +from HitAnalysis.HitAnalysisConf import CSCHitAnalysis +topSequence += CSCHitAnalysis('CSCHitAnalysis') +topSequence.CSCHitAnalysis.HistPath = '/MuonHitAnalysis/' +from HitAnalysis.HitAnalysisConf import TGCHitAnalysis +topSequence += TGCHitAnalysis('TGCHitAnalysis') +topSequence.TGCHitAnalysis.HistPath = '/MuonHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "MuonHitAnalysis DATAFILE='MuonHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/RPCHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/RPCHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..17160c375da896c878b09f565799b6ccea8f3ac1 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/RPCHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import RPCHitAnalysis +topSequence += RPCHitAnalysis() +RPCHitAnalysis = RPCHitAnalysis() +RPCHitAnalysis.NtupleFileName = 'RPCHitAnalysis' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "RPCHitAnalysis DATAFILE='RPCHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/SiHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/SiHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..e274b8fb193a21712ce067a50f7577251a529a34 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/SiHitAnalysis_topOptions.py @@ -0,0 +1,65 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import SiHitAnalysis +topSequence += SiHitAnalysis('PixelHitAnalysis') +topSequence.PixelHitAnalysis.CollectionName='PixelHits' +topSequence += SiHitAnalysis('SCTHitAnalysis') +topSequence.SCTHitAnalysis.CollectionName='SCT_Hits' +topSequence += SiHitAnalysis('BCMHitAnalysis') +topSequence.BCMHitAnalysis.CollectionName='BCMHits' +topSequence += SiHitAnalysis('BLMHitAnalysis') +topSequence.BLMHitAnalysis.CollectionName='BLMHits' + +topSequence.PixelHitAnalysis.HistPath='/SiHitAnalysis/histos/' +topSequence.SCTHitAnalysis.HistPath='/SiHitAnalysis/histos/' +topSequence.BCMHitAnalysis.HistPath='/SiHitAnalysis/histos/' +topSequence.BLMHitAnalysis.HistPath='/SiHitAnalysis/histos/' + +topSequence.PixelHitAnalysis.NtupleFileName='/SiHitAnalysis/ntuples/' +topSequence.SCTHitAnalysis.NtupleFileName='/SiHitAnalysis/ntuples/' +topSequence.BCMHitAnalysis.NtupleFileName='/SiHitAnalysis/ntuples/' +topSequence.BLMHitAnalysis.NtupleFileName='/SiHitAnalysis/ntuples/' + + +#Add some more TH2 histograms +topSequence.PixelHitAnalysis.ExpertMode= "off" +topSequence.SCTHitAnalysis.ExpertMode= "off" +topSequence.BCMHitAnalysis.ExpertMode= "off" +topSequence.BLMHitAnalysis.ExpertMode= "off" + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "SiHitAnalysis DATAFILE='SiHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = TRUE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include( "TrkDetDescrSvc/AtlasTrackingGeometrySvc.py" ) # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/TGCHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/TGCHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..ec86239c82d51d66d7391370a48fd040bb55f902 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/TGCHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import TGCHitAnalysis +topSequence += TGCHitAnalysis() +TGCHitAnalysis = TGCHitAnalysis() +TGCHitAnalysis.HistPath = '/TGCHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "TGCHitAnalysis DATAFILE='TGCHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/TRTHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/TRTHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..799b65c1b5a9d0bb709b2c1f0ff9ed0b89c4d956 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/TRTHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import TRTHitAnalysis +topSequence += TRTHitAnalysis() +TRTHitAnalysis = TRTHitAnalysis() +TRTHitAnalysis.HistPath = '/TRTHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "TRTHitAnalysis DATAFILE='TRTHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = TRUE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include( "TrkDetDescrSvc/AtlasTrackingGeometrySvc.py" ) # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/TrackRecordAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/TrackRecordAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..eb6c8b3d313cbd161f5cb9e6616ad9dbf7b9152e --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/TrackRecordAnalysis_topOptions.py @@ -0,0 +1,51 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import TrackRecordAnalysis +topSequence += TrackRecordAnalysis('TrackRecordCaloEntry') +topSequence.TrackRecordCaloEntry.CollectionName='CaloEntryLayer' +topSequence.TrackRecordCaloEntry.HistPath='/TrackRecordAnalysis/' +topSequence += TrackRecordAnalysis('TrackRecordMuonEntry') +topSequence.TrackRecordMuonEntry.CollectionName='MuonEntryLayer' +topSequence.TrackRecordMuonEntry.HistPath='/TrackRecordAnalysis/' +topSequence += TrackRecordAnalysis('TrackRecordMuonExit') +topSequence.TrackRecordMuonExit.CollectionName='MuonExitLayer' +topSequence.TrackRecordMuonExit.HistPath='/TrackRecordAnalysis/' + + + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "TrackRecordAnalysis DATAFILE='TrackRecordAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/TruthHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/TruthHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..8a872418ae88dfc76cde9daed85f38f7a64b0ffd --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/TruthHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import TruthHitAnalysis +topSequence += TruthHitAnalysis() +TruthHitAnalysis = TruthHitAnalysis() +TruthHitAnalysis.HistPath = '/TruthHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "TruthHitAnalysis DATAFILE='TruthHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/share/ZDCHitAnalysis_topOptions.py b/Simulation/Tools/HitAnalysis/share/ZDCHitAnalysis_topOptions.py new file mode 100755 index 0000000000000000000000000000000000000000..5e0ef00009e0d667a5008bb3fb4e165867b3074c --- /dev/null +++ b/Simulation/Tools/HitAnalysis/share/ZDCHitAnalysis_topOptions.py @@ -0,0 +1,43 @@ +#from AthenaCommon.AppMgr import ToolSvc +from AthenaCommon.AppMgr import ServiceMgr +import AthenaPoolCnvSvc.ReadAthenaPool + +from PartPropSvc.PartPropSvcConf import PartPropSvc + +include( "ParticleBuilderOptions/McAOD_PoolCnv_jobOptions.py") +include( "EventAthenaPool/EventAthenaPool_joboptions.py" ) + +import os +from glob import glob +from AthenaCommon.AthenaCommonFlags import athenaCommonFlags +athenaCommonFlags.FilesInput = glob( "/tmp/"+os.environ['USER']+"HITS*root*" ) +ServiceMgr.EventSelector.InputCollections = athenaCommonFlags.FilesInput() # This is stupid and redundant, but necessary + +from AthenaCommon.AlgSequence import AlgSequence +topSequence = AlgSequence() + +from HitAnalysis.HitAnalysisConf import ZDCHitAnalysis +topSequence += ZDCHitAnalysis() +ZDCHitAnalysis = ZDCHitAnalysis() +ZDCHitAnalysis.HistPath = '/ZDCHitAnalysis/' + +from GaudiSvc.GaudiSvcConf import THistSvc +ServiceMgr += THistSvc() +ServiceMgr.THistSvc.Output += [ "ZDCHitAnalysis DATAFILE='ZDCHitAnalysis.root' OPT='RECREATE'" ] + +ServiceMgr.MessageSvc.OutputLevel = INFO +ServiceMgr.MessageSvc.defaultLimit = 9999999 + +theApp.EvtMax = -1 + +ServiceMgr.AuditorSvc.Auditors += [ "ChronoAuditor"] + +AthenaPoolCnvSvc = Service("AthenaPoolCnvSvc") +AthenaPoolCnvSvc.UseDetailChronoStat = FALSE + +# To set up a geometry +from RecExConfig.AutoConfiguration import * +ConfigureFieldAndGeo() # Configure the settings for the geometry +include("RecExCond/AllDet_detDescr.py") # Actually load the geometry +#include("TrkDetDescrSvc/AtlasTrackingGeometrySvc.py") # Tracking geometry, handy for ID work + diff --git a/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..8e2dac96fca9f38ba29c7886299275fd21f7447c --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.cxx @@ -0,0 +1,186 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "AFPHitAnalysis.h" + +#include "AFP_SimEv/AFP_SIDSimHit.h" +#include "AFP_SimEv/AFP_SIDSimHitCollection.h" + + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +AFPHitAnalysis::AFPHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hitID(0) + , h_pdgID(0) + , h_trackID(0) + , h_kine(0) + , h_edep(0) + , h_stepX(0) + , h_stepY(0) + , h_stepZ(0) + , h_time(0) + , h_stationID(0) + , h_detID(0) + , h_pixelRow(0) + , h_pixelCol(0) + , m_tree(0) + , m_ntupleFileName("/AFPHitAnalysis/ntuples/") + , m_path("/AFPHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("HistPath", m_path); + declareProperty("NtupleFileName", m_ntupleFileName); +} + +StatusCode AFPHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing AFPHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + /** now define the histograms**/ + h_hitID = new TH1D("h_hitID", "hitID",100, 0., 100.); + h_hitID->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_hitID->GetName(), h_hitID)); + + h_pdgID = new TH1D("h_pdgID", "pdgID", 200, -100,100); + h_pdgID->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_pdgID->GetName(), h_pdgID)); + + h_trackID = new TH1D("h_trackID", "trackID", 100, 0,100); + h_trackID->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_trackID->GetName(), h_trackID)); + + h_kine = new TH1D("h_kine", "kine", 100, 0,1000); + h_kine->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_kine->GetName(), h_kine)); + + h_edep = new TH1D("h_edep", "edep", 100, 0,1000); + h_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_edep->GetName(), h_edep)); + + h_stepX = new TH1D("h_stepX", "stepX", 100, 0,1000); + h_stepX->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_stepX->GetName(), h_stepX)); + + h_stepY = new TH1D("h_stepY", "stepY", 100, 0,1000); + h_stepY->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_stepY->GetName(), h_stepY)); + + h_stepZ = new TH1D("h_stepZ", "stepZ", 100, 0,1000); + h_stepZ->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_stepZ->GetName(), h_stepZ)); + + h_stationID = new TH1D("h_stationID", "stationID", 50, 0,50); + h_stationID->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_stationID->GetName(), h_stationID)); + + h_detID = new TH1D("h_detID", "detID", 50, 0,50); + h_detID->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_detID->GetName(), h_detID)); + + h_pixelRow = new TH1D("h_pixelRow", "pixelRow", 20, 0,20); + h_pixelRow->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_pixelRow->GetName(), h_pixelRow)); + + + h_pixelCol = new TH1D("h_pixelCol", "pixelCol", 20, 0,20); + h_pixelCol->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_pixelCol->GetName(), h_pixelCol)); + + /** now add branches and leaves to the tree */ + m_tree= new TTree("NtupleAFPHitAnanalysis","AFPHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + if (m_tree){ + m_tree->Branch("hitID", &m_hitID); + m_tree->Branch("pdgID", &m_pdgID); + m_tree->Branch("trackID", &m_trackID); + m_tree->Branch("kine", &m_kine); + m_tree->Branch("edep", &m_edep); + m_tree->Branch("stepX", &m_stepX); + m_tree->Branch("stepY", &m_stepY); + m_tree->Branch("stepZ", &m_stepZ); + m_tree->Branch("stationID", &m_stationID); + m_tree->Branch("detID", &m_detID); + m_tree->Branch("pixelRow", &m_pixelRow); + m_tree->Branch("pixelCol", &m_pixelCol); + + + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + return StatusCode::SUCCESS; +} + + + +StatusCode AFPHitAnalysis::execute() { + ATH_MSG_DEBUG( "In AFPHitAnalysis::execute()" ); + + m_hitID->clear(); + m_pdgID->clear(); + m_trackID->clear(); + m_kine->clear(); + m_edep->clear(); + m_stepX->clear(); + m_stepY->clear(); + m_stepZ->clear(); + m_time->clear(); + m_stationID->clear(); + m_detID->clear(); + m_pixelRow->clear(); + m_pixelCol->clear(); + + AFP_SIDSimHitConstIter hi; + const DataHandle<AFP_SIDSimHitCollection> iter; + CHECK(evtStore()->retrieve(iter,"AFP_SIDSimHitCollection")); + for(hi=(*iter).begin(); hi != (*iter).end();++hi){ + AFP_SIDSimHit ghit(*hi); + h_hitID->Fill(ghit.m_nHitID); + h_pdgID->Fill(ghit.m_nParticleEncoding); + h_trackID->Fill(ghit.m_nTrackID); + h_kine->Fill(ghit.m_fKineticEnergy); + h_edep->Fill(ghit.m_fEnergyDeposit); + h_stepX->Fill(ghit.m_fPostStepX-ghit.m_fPreStepX); + h_stepY->Fill(ghit.m_fPostStepY-ghit.m_fPreStepY); + h_stepX->Fill(ghit.m_fPostStepZ-ghit.m_fPreStepZ); + h_time->Fill(ghit.m_fGlobalTime); + h_stationID->Fill(ghit.m_nStationID); + h_detID->Fill(ghit.m_nDetectorID); + h_pixelRow->Fill(ghit.m_nPixelRow); + h_pixelCol->Fill(ghit.m_nPixelCol); + + + m_hitID->push_back(ghit.m_nHitID); + m_pdgID->push_back(ghit.m_nParticleEncoding); + m_trackID->push_back(ghit.m_nTrackID); + m_kine->push_back(ghit.m_fKineticEnergy); + m_edep->push_back(ghit.m_fEnergyDeposit); + m_stepX->push_back(ghit.m_fPostStepX-ghit.m_fPreStepX); + m_stepY->push_back(ghit.m_fPostStepY-ghit.m_fPreStepY); + m_stepX->push_back(ghit.m_fPostStepZ-ghit.m_fPreStepZ); + m_time->push_back(ghit.m_fGlobalTime); + m_stationID->push_back(ghit.m_nStationID); + m_detID->push_back(ghit.m_nDetectorID); + m_pixelRow->push_back(ghit.m_nPixelRow); + m_pixelCol->push_back(ghit.m_nPixelCol); + } + + if (m_tree) m_tree->Fill(); + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..7dfd58e0e680fd1f5f66404d7b06b290852ce774 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/AFPHitAnalysis.h @@ -0,0 +1,70 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef AFP_HIT_ANALYSIS_H +#define AFP_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TTree.h" + + +class TH1; +class TTree; + +class AFPHitAnalysis : public AthAlgorithm { + + public: + + AFPHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~AFPHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some histograms**/ + TH1* h_hitID; + TH1* h_pdgID; + TH1* h_trackID; + TH1* h_kine; + TH1* h_edep; + TH1* h_stepX; + TH1* h_stepY; + TH1* h_stepZ; + TH1* h_time; + TH1* h_stationID; + TH1* h_detID; + TH1* h_pixelRow; + TH1* h_pixelCol; + + std::vector<float>* m_hitID; + std::vector<float>* m_pdgID; + std::vector<float>* m_trackID; + std::vector<float>* m_kine; + std::vector<float>* m_edep; + std::vector<float>* m_stepX; + std::vector<float>* m_stepY; + std::vector<float>* m_stepZ; + std::vector<float>* m_time; + std::vector<int>* m_stationID; + std::vector<int>* m_detID; + std::vector<int>* m_pixelRow; + std::vector<int>* m_pixelCol; + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // AFP_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.cxx new file mode 100644 index 0000000000000000000000000000000000000000..b1074bce7cd5d450b89eae56ecb73011be75cf41 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.cxx @@ -0,0 +1,169 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "ALFAHitAnalysis.h" + +// Section of includes for Pixel and SCT tests +#include "ALFA_SimEv/ALFA_HitCollection.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> + +ALFAHitAnalysis::ALFAHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , m_station(0) + , m_plate(0) + , m_fiber(0) + , m_sign(0) + , m_energy(0) + , m_tree(0) + , m_ntupleFileName("/ALFAHitsAnalysis/ntuples/") + , m_path("/ALFAHitsAnalysis/histos") + , m_thistSvc("THistSvc", name) +{ + for(int i(0); i<8; i++){ + h_E_full_sum_h[i]=0; + h_E_layer_sum_h[i]=0; + h_hit_layer[i]=0; + h_hit_fiber[i]=0; + + } + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); +} + +StatusCode ALFAHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing ALFAHitAnalysis" ); + + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + + /** now add branches and leaves to the tree */ + std::stringstream s; + for(unsigned int j=0; j<8;j++){ + s.str(""); + s <<"edep_in_det_no."<< j+1; + float Emax = 5; + if(j==3) Emax=150; + h_E_full_sum_h[j] = new TH1D(s.str().c_str(),s.str().c_str(), 100,0,Emax); + h_E_full_sum_h[j]->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_E_full_sum_h[j]->GetName(), h_E_full_sum_h[j])); + + s.str(""); + s <<"edep_per_layer_det_no."<< j+1; + h_E_layer_sum_h[j] = new TH1D(s.str().c_str(),s.str().c_str(), 100,0,15); + h_E_layer_sum_h[j]->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_E_layer_sum_h[j]->GetName(), h_E_layer_sum_h[j])); + + s.str(""); + s <<"hit_layer_det_no."<< j+1; + h_hit_layer[j] = new TH1D(s.str().c_str(),s.str().c_str(),50,0,50); + h_hit_layer[j]->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_hit_layer[j]->GetName(), h_hit_layer[j])); + + s.str(""); + s <<"hit_fiber_det_no."<< j+1; + h_hit_fiber[j] = new TH1D(s.str().c_str(),s.str().c_str(),100,0,60); + h_hit_fiber[j]->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path+h_hit_fiber[j]->GetName(), h_hit_fiber[j])); + + } + + m_tree= new TTree("NtupleALFAHitAnalysis", "ALFAHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("station", &m_station); + m_tree->Branch("plate", &m_plate); + m_tree->Branch("fiber", &m_fiber); + m_tree->Branch("energy", &m_energy); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + return StatusCode::SUCCESS; +} + + + +StatusCode ALFAHitAnalysis::execute() { + ATH_MSG_DEBUG( "In ALFAHitAnalysis::execute()" ); + + m_station->clear(); + m_plate->clear(); + m_fiber->clear(); + m_sign->clear(); + m_energy->clear(); + + //cleaning + int fiber,plate,sign,station; + double E_fiber_sum[8][10][64][2], E_full_sum[8], E_layer_sum[8][20]; + for (int l= 0;l<8;l++){ + E_full_sum[l] = 0.; + for ( int i = 0; i < 10; i++ ){ + E_layer_sum[l][i] = 0.; + E_layer_sum[l][i+10] = 0.; + for ( int j = 0; j < 64; j++ ){ + for ( int k = 0; k < 2; k++ ){ + E_fiber_sum[l][i][j][k] = 0.; + } + } + } + } + + ALFA_HitConstIter iter; + const DataHandle<ALFA_HitCollection> col_alfa; + CHECK(evtStore()->retrieve(col_alfa,"ALFA_HitCollection")); + for(iter=(*col_alfa).begin(); iter!=(*col_alfa).end();++iter){ + station = (*iter).GetStationNumber(); + plate = (*iter).GetPlateNumber(); + fiber = (*iter).GetFiberNumber(); + sign = (*iter).GetSignFiber(); + E_fiber_sum[station-1][plate-1][fiber-1][(1-sign)/2]+=((*iter).GetEnergyDeposit()); + + m_station->push_back(station); + m_plate->push_back(plate); + m_fiber->push_back(fiber); + m_sign->push_back(sign); + m_energy->push_back((*iter).GetEnergyDeposit()); + } // End iteration + + for(int l=0;l<8;l++){ + for(int i=0;i<10; i++){ + for(int j=0;j<64;j++){ + for(int k=0;k<2;k++){ + E_full_sum[l]+=E_fiber_sum[l][i][j][k]; + E_layer_sum[l][2*i+k]+=E_fiber_sum[l][i][j][k]; + if(E_fiber_sum[l][i][j][k]>0.){ + h_hit_layer[l]->Fill(2*i+k+1); + h_hit_fiber[l]->Fill(j+1); + } + } + } + } + } + for(int l = 0;l<8;l++){ + h_E_full_sum_h[l]->Fill(E_full_sum[l]); + for(int i = 0; i< 20;i++){ + h_E_layer_sum_h[l]->Fill(E_layer_sum[l][i]); + } + } + + if (m_tree) m_tree->Fill(); + return StatusCode::SUCCESS; +} + diff --git a/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..8d833d12f703fef1138de93d749e410991476558 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/ALFAHitAnalysis.h @@ -0,0 +1,58 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef ALFA_HIT_ANALYSIS_H +#define ALFA_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + + +#include <string> +#include <vector> +#include "TH1.h" +#include "TTree.h" + +class TH1; +class TTree; + +class ALFAHitAnalysis : public AthAlgorithm { + + public: + + ALFAHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~ALFAHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + std::string m_collection; + /** Some variables**/ + TH1* h_E_full_sum_h[8]; + TH1* h_E_layer_sum_h[8]; + TH1* h_hit_layer[8]; + TH1* h_hit_fiber[8]; + + std::vector<int>* m_station; + std::vector<int>* m_plate; + std::vector<int>* m_fiber; + std::vector<int>* m_sign; + std::vector<double>* m_energy; + + + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + + +}; + +#endif // ALFA_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..a03e6f7fad3582a399d1a77acf02478b6fa8cbd7 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.cxx @@ -0,0 +1,248 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "CSCHitAnalysis.h" + +// Section of includes for CSC of the Muon Spectrometer tests +#include "GeoAdaptors/GeoMuonHits.h" + +#include "MuonSimEvent/CSCSimHitCollection.h" +#include "MuonSimEvent/CSCSimHit.h" +#include "CLHEP/Vector/LorentzVector.h" + +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +CSCHitAnalysis::CSCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_zr(0) + , h_hits_eta(0) + , h_hits_phi(0) + , h_hits_sx(0) + , h_hits_sy(0) + , h_hits_sz(0) + , h_hits_ex(0) + , h_hits_ey(0) + , h_hits_ez(0) + , h_hits_time(0) + , h_hits_edep(0) + , h_hits_kine(0) + , m_hits_x(0) + , m_hits_y(0) + , m_hits_z(0) + , m_hits_r(0) + , m_hits_eta(0) + , m_hits_phi(0) + , m_hits_start_x(0) + , m_hits_start_y(0) + , m_hits_start_z(0) + , m_hits_end_x(0) + , m_hits_end_y(0) + , m_hits_end_z(0) + , m_hits_time(0) + , m_hits_edep(0) + , m_hits_kine(0) + , m_tree(0) + , m_ntupleFileName("/CSCHitAnalysis/ntuples/") + , m_path("/CSCHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); +} + +StatusCode CSCHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing CSCHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + /** Histograms**/ + h_hits_x = new TH1D("h_csc_hits_x","hits_x", 100,-2000, 2000); + h_hits_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + + h_hits_y = new TH1D("h_csc_hits_y", "hits_y", 100,-2000,2000); + h_hits_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + + h_hits_z = new TH1D("h_csc_hits_z", "hits_z", 100,-10000,10000); + h_hits_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + + h_hits_r = new TH1D("h_csc_hits_r", "hits_r", 100,500,2500); + h_hits_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + + h_xy = new TH2D("h_csc_xy", "xy", 100,-2000.,2000.,100, -2000., 2000.); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D("h_csc_zr", "zr", 100,-10000.,10000.,100, 500., 2500.); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_hits_eta = new TH1D("h_csc_hits_eta", "hits_eta", 100,-3.0,3.0); + h_hits_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_eta->GetName(), h_hits_eta)); + + h_hits_phi = new TH1D("h_csc_hits_phi", "hits_phi", 100,-3.2,3.2); + h_hits_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_phi->GetName(), h_hits_phi)); + + h_hits_sx = new TH1D("h_csc_hits_sx","hits_sx", 100,-10, 10); + h_hits_sx->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_sx->GetName(), h_hits_sx)); + + h_hits_sy = new TH1D("h_csc_hits_sy", "hits_sy", 100,-500,500); + h_hits_sy->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_sy->GetName(), h_hits_sy)); + + h_hits_sz = new TH1D("h_csc_hits_sz", "hits_sz", 100,-1000,1000); + h_hits_sz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_sz->GetName(), h_hits_sz)); + + + h_hits_ex = new TH1D("h_csc_hits_ex","hits_ex", 100,-10, 10); + h_hits_ex->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ex->GetName(), h_hits_ex)); + + h_hits_ey = new TH1D("h_csc_hits_ey", "hits_ey", 100,-500,500); + h_hits_ey->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ey->GetName(), h_hits_ey)); + + h_hits_ez = new TH1D("h_csc_hits_ez", "hits_ez", 100,-1000,1000); + h_hits_ez->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ez->GetName(), h_hits_ez)); + + h_hits_time = new TH1D("h_csc_hits_time","hits_time", 100,20, 40); + h_hits_time->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_time->GetName(), h_hits_time)); + + h_hits_edep = new TH1D("h_csc_hits_edep", "hits_edep", 100,0,0.1); + h_hits_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_edep->GetName(), h_hits_edep)); + + h_hits_kine = new TH1D("h_csc_hits_kine", "hits_kine", 100,0,3000); + h_hits_kine->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_kine->GetName(), h_hits_kine)); + + m_tree= new TTree("NtupleCSCHitAnalysis","CSCHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("x", &m_hits_x); + m_tree->Branch("y", &m_hits_y); + m_tree->Branch("z", &m_hits_z); + m_tree->Branch("r", &m_hits_r); + m_tree->Branch("eta", &m_hits_eta); + m_tree->Branch("phi", &m_hits_phi); + m_tree->Branch("startX", &m_hits_start_x); + m_tree->Branch("startY", &m_hits_start_y); + m_tree->Branch("startZ", &m_hits_start_z); + m_tree->Branch("endX", &m_hits_end_x); + m_tree->Branch("endY", &m_hits_end_y); + m_tree->Branch("endZ", &m_hits_end_z); + m_tree->Branch("time", &m_hits_time); + m_tree->Branch("edep", &m_hits_edep); + m_tree->Branch("kine", &m_hits_kine); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + return StatusCode::SUCCESS; +} + + + +StatusCode CSCHitAnalysis::execute() { + ATH_MSG_DEBUG( "In CSCHitAnalysis::execute()" ); + + m_hits_x->clear(); + m_hits_y->clear(); + m_hits_z->clear(); + m_hits_r->clear(); + m_hits_eta->clear(); + m_hits_phi->clear(); + m_hits_start_x->clear(); + m_hits_start_y->clear(); + m_hits_start_z->clear(); + m_hits_end_x->clear(); + m_hits_end_y->clear(); + m_hits_end_z->clear(); + m_hits_time->clear(); + m_hits_edep->clear(); + m_hits_kine->clear(); + + + + const DataHandle<CSCSimHitCollection> csc_container; + if (evtStore()->retrieve(csc_container, "CSC_Hits" )==StatusCode::SUCCESS) { + for(CSCSimHitCollection::const_iterator i_hit = csc_container->begin(); + i_hit != csc_container->end();++i_hit){ + //CSCSimHitCollection::const_iterator i_hit; + //for(auto i_hit : *csc_container){ + GeoCSCHit ghit(*i_hit); + if(!ghit) continue; + Amg::Vector3D p = ghit.getGlobalPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_zr->Fill(p.z(), p.perp()); + h_hits_eta->Fill(p.eta()); + h_hits_phi->Fill(p.phi()); + h_hits_sx->Fill( (*i_hit).getHitStart().x()); + h_hits_sy->Fill( (*i_hit).getHitStart().y()); + h_hits_sz->Fill( (*i_hit).getHitStart().z()); + h_hits_ex->Fill( (*i_hit).getHitEnd().x()); + h_hits_ey->Fill( (*i_hit).getHitEnd().y()); + h_hits_ez->Fill( (*i_hit).getHitEnd().z()); + h_hits_edep->Fill((*i_hit).energyDeposit()); + h_hits_time->Fill((*i_hit).globalTime()); + h_hits_kine->Fill((*i_hit).kineticEnergy()); + + m_hits_x->push_back(p.x()); + m_hits_y->push_back(p.y()); + m_hits_z->push_back(p.z()); + m_hits_r->push_back(p.perp()); + m_hits_eta->push_back(p.eta()); + m_hits_phi->push_back(p.phi()); + m_hits_start_x->push_back( (*i_hit).getHitStart().x()); + m_hits_start_y->push_back( (*i_hit).getHitStart().y()); + m_hits_start_z->push_back( (*i_hit).getHitStart().z()); + m_hits_end_x->push_back( (*i_hit).getHitEnd().x()); + m_hits_end_y->push_back( (*i_hit).getHitEnd().y()); + m_hits_end_z->push_back( (*i_hit).getHitEnd().z()); + m_hits_edep->push_back((*i_hit).energyDeposit()); + m_hits_time->push_back((*i_hit).globalTime()); + m_hits_kine->push_back((*i_hit).kineticEnergy()); + + } + } // End while hits + if (m_tree) m_tree->Fill(); + + + + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..c85d3e63a672d986c28635058f09125ae7c9415e --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/CSCHitAnalysis.h @@ -0,0 +1,80 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef CSC_HIT_ANALYSIS_H +#define CSC_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + + +class TH1; +class TH2; +class TTree; + +class CSCHitAnalysis : public AthAlgorithm { + + public: + + CSCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~CSCHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_zr; + TH1* h_hits_eta; + TH1* h_hits_phi; + TH1* h_hits_sx; + TH1* h_hits_sy; + TH1* h_hits_sz; + TH1* h_hits_ex; + TH1* h_hits_ey; + TH1* h_hits_ez; + TH1* h_hits_time; + TH1* h_hits_edep; + TH1* h_hits_kine; + + std::vector<float>* m_hits_x; + std::vector<float>* m_hits_y; + std::vector<float>* m_hits_z; + std::vector<float>* m_hits_r; + std::vector<float>* m_hits_eta; + std::vector<float>* m_hits_phi; + std::vector<float>* m_hits_start_x; + std::vector<float>* m_hits_start_y; + std::vector<float>* m_hits_start_z; + std::vector<float>* m_hits_end_x; + std::vector<float>* m_hits_end_y; + std::vector<float>* m_hits_end_z; + std::vector<float>* m_hits_time; + std::vector<float>* m_hits_edep; + std::vector<float>* m_hits_kine; + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // CSC_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.cxx index 9bf5fd2cdca5e9c62fce14e7fe5365f7292404b2..ae15f65605266565fcca81e86df14fb00676fe6d 100755 --- a/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.cxx +++ b/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.cxx @@ -2,6 +2,9 @@ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration */ +// Base class +#include "CaloHitAnalysis.h" + // Section of includes for LAr calo tests #include "LArSimEvent/LArHitContainer.h" #include "CaloDetDescr/CaloDetDescrElement.h" @@ -14,20 +17,14 @@ #include "TileSimEvent/TileHit.h" #include "TileSimEvent/TileHitVector.h" -#include "GaudiKernel/MsgStream.h" -#include "GaudiKernel/AlgFactory.h" -#include "GaudiKernel/IToolSvc.h" -#include "GaudiKernel/ITHistSvc.h" - -#include "StoreGate/StoreGateSvc.h" +//Section of includes for Calibrated Calo hits +#include "GeoAdaptors/GeoCaloCalibHit.h" +#include "CaloSimEvent/CaloCalibrationHitContainer.h" -#include "EventInfo/EventInfo.h" -#include "EventInfo/EventID.h" - -#include "TTree.h" #include "TString.h" - -#include "HitAnalysis/CaloHitAnalysis.h" +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" #include <algorithm> #include <math.h> @@ -36,131 +33,354 @@ CaloHitAnalysis::CaloHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) : AthAlgorithm(name, pSvcLocator) - , m_storeGate(0) + , h_cell_eta(0) + , h_cell_phi(0) + , h_cell_e(0) + , h_cell_radius(0) + , h_xy(0) + , h_zr(0) + , h_etaphi(0) + , h_time_e(0) + , h_eta_e(0) + , h_phi_e(0) + , h_r_e(0) + , h_calib_eta(0) + , h_calib_phi(0) + , h_calib_rz(0) + , h_calib_etaphi(0) + , h_calib_eEM(0) + , h_calib_eNonEM(0) + , h_calib_eInv(0) + , h_calib_eEsc(0) + , h_calib_eTot(0) + , h_calib_eTotpartID(0) , m_cell_eta(0) , m_cell_phi(0) + , m_cell_x(0) + , m_cell_y(0) + , m_cell_z(0) , m_cell_e(0) , m_cell_radius(0) + , m_time(0) + , m_calib_eta(0) + , m_calib_phi(0) + , m_calib_radius(0) + , m_calib_z(0) + , m_calib_eEM(0) + , m_calib_eNonEM(0) + , m_calib_eInv(0) + , m_calib_eEsc(0) + , m_calib_eTot(0) + , m_calib_partID(0) + , m_expert("off") + , m_calib("off") + , m_thistSvc("THistSvc",name) + , m_path("/CaloHitAnalysis/") , m_tree(0) - , m_ntupleFileName("/ntuples/file1") - , m_ntupleDirName("CaloHitAnalysis") - , m_ntupleTreeName("CaloHitAna") - , m_thistSvc(0) + , m_ntupleFileName("/ntuples/") , m_tileID(0) , m_tileMgr(0) { - declareProperty("NtupleFileName", m_ntupleFileName); - declareProperty("NtupleDirectoryName", m_ntupleDirName); - declareProperty("NtupleTreeName", m_ntupleTreeName); + declareProperty("HistPath", m_path); + declareProperty("ExpertMode", m_expert="off"); + declareProperty("CalibHits", m_calib="off"); + declareProperty("NtupleFileName", m_ntupleFileName); } - -CaloHitAnalysis::~CaloHitAnalysis() -{;} - StatusCode CaloHitAnalysis::initialize() { ATH_MSG_DEBUG( "Initializing CaloHitAnalysis" ); - StoreGateSvc* detStore; - StatusCode sc = service("DetectorStore",detStore); - if (sc.isFailure()) { - ATH_MSG_ERROR( "Unable to get pointer to Detector Store Service" ); - return sc; - } + CHECK( detStore()->retrieve(m_tileMgr) ); - sc = detStore->retrieve(m_tileMgr); - if (sc.isFailure()) { - ATH_MSG_ERROR( "Unable to retrieve TileDetDescrManager from DetectorStore" ); - m_tileMgr=0; - } - - sc = detStore->retrieve(m_tileID); - if (sc.isFailure()) { - ATH_MSG_ERROR( "Unable to retrieve TileID helper from DetectorStore" ); - m_tileID=0; - } - - /** get a handle of StoreGate for access to the Event Store */ - sc = service("StoreGateSvc", m_storeGate); - if (sc.isFailure()) { - ATH_MSG_ERROR( "Unable to retrieve pointer to StoreGateSvc" ); - return sc; - } + CHECK( detStore()->retrieve(m_tileID) ); // Grab the Ntuple and histogramming service for the tree - sc = service("THistSvc",m_thistSvc); - if (sc.isFailure()) { - ATH_MSG_ERROR( "Unable to retrieve pointer to THistSvc" ); - return sc; - } + CHECK( m_thistSvc.retrieve() ); - m_tree = new TTree( TString(m_ntupleTreeName), "CaloHitAna" ); - std::string fullNtupleName = "/"+m_ntupleFileName+"/"+m_ntupleDirName+"/"+m_ntupleTreeName; - sc = m_thistSvc->regTree(fullNtupleName, m_tree); - if (sc.isFailure()) { - ATH_MSG_ERROR("Unable to register TTree: " << fullNtupleName); - return sc; - } - - /** now add branches and leaves to the tree */ - if (m_tree){ - m_tree->Branch("CellEta", &m_cell_eta); - m_tree->Branch("CellPhi", &m_cell_phi); - m_tree->Branch("CellE", &m_cell_e); - m_tree->Branch("CellRadius", &m_cell_radius); + h_cell_e = new TH1D("h_Calo_cell_e", "cell_e", 100,0.,500.); + h_cell_e->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_cell_e->GetName(), h_cell_e)); + + h_cell_eta = new TH1D("h_Calo_cell_eta", "cell_eta", 50,-5.,5.); + h_cell_eta->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_cell_eta->GetName(), h_cell_eta)); + + h_cell_phi = new TH1D("h_Calo_cell_phi", "cell_phi", 50,-3.1416,3.1416); + h_cell_phi->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_cell_phi->GetName(), h_cell_phi)); + + h_cell_radius = new TH1D("h_Calo_cell_radius", "cell_radius", 100, 0., 6000.); + h_cell_radius->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_cell_radius->GetName(), h_cell_radius)); + + h_xy = new TH2F("h_Calo_xy", "xy", 100,-4000,4000,100, -4000, 4000); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D("h_Calo_zr", "zr", 100,-7000.,7000.,100, 0., 6000.); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_etaphi = new TH2D("h_Calo_etaphi", "eta_phi", 50,-5.,5.,50, -3.1416, 3.1416); + h_etaphi->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_etaphi->GetName(), h_etaphi)); + + //These histograms will be filled only if expert mode is set on + h_time_e = new TH2D("h_Calo_time_e", "energy vs time", 100, 0,50, 100,0,500); + h_time_e->StatOverflows(); + + h_eta_e = new TH2D("h_Calo_eta_e", "energy vs eta", 50, -5,5, 100,0,500); + h_eta_e->StatOverflows(); + + h_phi_e = new TH2D("h_Calo_phi_e", "energy vs phi", 50, -3.1416,3.1416, 100,0,500); + h_phi_e->StatOverflows(); + + h_r_e = new TH2D("h_Calo_r_e", "energy vs radius", 100, 0,6000, 100,0,500); + h_r_e->StatOverflows(); + + if( m_expert == "on") + { + CHECK(m_thistSvc->regHist( m_path+h_time_e->GetName(), h_time_e)); + CHECK(m_thistSvc->regHist( m_path+h_eta_e->GetName(), h_eta_e)); + CHECK(m_thistSvc->regHist( m_path+h_phi_e->GetName(), h_phi_e)); + CHECK(m_thistSvc->regHist( m_path+h_r_e->GetName(), h_r_e)); + } + + //Histograms for calibrated hits + h_calib_eta = new TH1D("h_calib_eta", "calib. hits eta", 50,-5,5); + h_calib_eta->StatOverflows(); + + h_calib_phi = new TH1D("h_calib_phi", "calib. hits phi", 50,-3.1416,3.1416); + h_calib_phi->StatOverflows(); + + h_calib_rz = new TH2D("h_calib_rz", "calib. hits r vs z", 100,-7000,7000,1000, 0,6000); + h_calib_rz->StatOverflows(); + + h_calib_etaphi = new TH2D("h_calib_etaphi", "calib. hits eta vs phi",50,-5.,5., 50,-3.1416,3.1416); + h_calib_etaphi->StatOverflows(); + + h_calib_eEM = new TH1D("h_calib_eEM", "calib. hits EM energy", 100,0,100); + h_calib_eEM->StatOverflows(); + + h_calib_eNonEM = new TH1D("h_calib_nonEM", "calib. hits non EM energy", 100,0,100); + h_calib_eNonEM->StatOverflows(); + + h_calib_eInv = new TH1D("h_calib_eInv", "calib. hits invisible energy", 100,0,100); + h_calib_eInv->StatOverflows(); + + h_calib_eEsc = new TH1D("h_calib_eEsc", "calib. hits escaped energy", 100,0,100); + h_calib_eEsc->StatOverflows(); + + h_calib_eTot = new TH1D("h_calib_eTot", "calib. hits energy", 100,0,100); + h_calib_eTot->StatOverflows(); + + h_calib_eTotpartID = new TH1D("h_calib_eTotpartID", "calib. hits partID weighted with energy",600,0,300000); + h_calib_eTotpartID->StatOverflows(); + + if( m_calib == "on") + { + CHECK(m_thistSvc->regHist( m_path+h_calib_eta->GetName(), h_calib_eta)); + CHECK(m_thistSvc->regHist( m_path+h_calib_phi->GetName(), h_calib_phi)); + CHECK(m_thistSvc->regHist( m_path+h_calib_rz->GetName(), h_calib_rz)); + CHECK(m_thistSvc->regHist( m_path+h_calib_etaphi->GetName(), h_calib_etaphi)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eEM->GetName(), h_calib_eEM)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eNonEM->GetName(), h_calib_eNonEM)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eInv->GetName(), h_calib_eInv)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eEsc->GetName(), h_calib_eEsc)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eTot->GetName(), h_calib_eTot)); + CHECK(m_thistSvc->regHist( m_path+h_calib_eTotpartID->GetName(), h_calib_eTotpartID)); + } + + + m_tree = new TTree( "CaloHitNtuple", "CaloHitAna" ); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK( m_thistSvc->regTree(fullNtupleName, m_tree) ); + + if(m_tree) + { + m_tree->Branch("CellEta", &m_cell_eta); + m_tree->Branch("CellPhi", &m_cell_phi); + m_tree->Branch("CellX", &m_cell_x); + m_tree->Branch("CellY", &m_cell_y); + m_tree->Branch("CellZ", &m_cell_z); + m_tree->Branch("CellE", &m_cell_e); + m_tree->Branch("CellRadius", &m_cell_radius); + m_tree->Branch("Time", &m_time); + m_tree->Branch("CalibEta", &m_calib_eta); + m_tree->Branch("CalibPhi", &m_calib_phi); + m_tree->Branch("CalibRadius", &m_calib_radius); + m_tree->Branch("CalibZ", &m_calib_z); + m_tree->Branch("Calib_eEM", &m_calib_eEM); + m_tree->Branch("Calib_eNonEM", &m_calib_eNonEM); + m_tree->Branch("Calib_eInv", &m_calib_eInv); + m_tree->Branch("Calib_eEsc", &m_calib_eEsc); + m_tree->Branch("Calib_eTot", &m_calib_eTot); + m_tree->Branch("Calib_partID", &m_calib_partID); + + }else{ + ATH_MSG_ERROR( "No tree found!" ); } return StatusCode::SUCCESS; } -StatusCode CaloHitAnalysis::finalize() { - return StatusCode::SUCCESS; -} - StatusCode CaloHitAnalysis::execute() { ATH_MSG_DEBUG( "In CaloHitAnalysis::execute()" ); m_cell_eta->clear(); m_cell_phi->clear(); m_cell_e->clear(); + m_cell_x->clear(); + m_cell_y->clear(); + m_cell_z->clear(); m_cell_radius->clear(); + m_time->clear(); + m_calib_eta->clear(); + m_calib_phi->clear(); + m_calib_radius->clear(); + m_calib_z->clear(); + m_calib_eEM->clear(); + m_calib_eNonEM->clear(); + m_calib_eInv->clear(); + m_calib_eEsc->clear(); + m_calib_eTot->clear(); + m_calib_partID->clear(); std::string lArKey [4] = {"LArHitEMB", "LArHitEMEC", "LArHitFCAL", "LArHitHEC"}; for (unsigned int i=0;i<4;i++){ const DataHandle<LArHitContainer> iter; - if (m_storeGate->retrieve(iter,lArKey[i])==StatusCode::SUCCESS) { + if (evtStore()->retrieve(iter,lArKey[i])==StatusCode::SUCCESS) { LArHitContainer::const_iterator hi; - for (hi=(*iter).begin();hi!=(*iter).end();hi++){ - GeoLArHit ghit(**hi); + for (auto hi : *iter ){ + GeoLArHit ghit(*hi); if (!ghit) continue; - const CaloDetDescrElement *hitElement = ghit.getDetDescrElement(); - m_cell_e->push_back( ghit.Energy() ); - m_cell_eta->push_back( hitElement->eta() ); - m_cell_phi->push_back( hitElement->phi() ); - m_cell_radius->push_back( hitElement->r() ); + const CaloDetDescrElement *hitElement = ghit.getDetDescrElement(); + double energy = ghit.Energy(); + double time = ghit.Time(); + double eta = hitElement->eta(); + double phi = hitElement->phi(); + double radius = hitElement->r(); + float x = hitElement->x(); + float y = hitElement->y(); + double z = hitElement->z(); + + h_cell_e->Fill( energy ); + h_cell_eta->Fill( eta ); + h_cell_phi->Fill( phi ); + h_cell_radius->Fill( radius ); + h_xy->Fill(x,y); + h_zr->Fill(z,radius); + h_etaphi->Fill(eta, phi); + if( m_expert == "on"){ + h_time_e->Fill(time, energy); + h_eta_e->Fill(eta, energy); + h_phi_e->Fill(phi, energy); + h_r_e->Fill(radius, energy); + + } + m_cell_eta->push_back(eta); + m_cell_phi->push_back(phi); + m_cell_e->push_back(energy); + m_cell_x->push_back(x); + m_cell_y->push_back(y); + m_cell_z->push_back(z); + m_cell_radius->push_back(radius); + m_time->push_back(time); + } // End while hits } // End statuscode success upon retrieval of hits } // End detector type loop - - const TileHitVector * hitVec; - if (m_storeGate->retrieve(hitVec,"TileHitVec")==StatusCode::SUCCESS && + + const DataHandle<TileHitVector> hitVec; + //const TileHitVector* hitVec; + if (evtStore()->retrieve(hitVec,"TileHitVec")==StatusCode::SUCCESS && m_tileMgr && m_tileID ){ - for(TileHitVecConstIterator i_hit=hitVec->begin() ; i_hit!=hitVec->end() ; ++i_hit){ - Identifier pmt_id = (*i_hit).identify(); + for(auto i_hit : *hitVec ){ + Identifier pmt_id = (i_hit).identify(); Identifier cell_id = m_tileID->cell_id(pmt_id); - const CaloDetDescrElement* ddElement = m_tileMgr->get_cell_element(cell_id); - double tot_e=0.; - for (int t=0;t<(*i_hit).size();++t) tot_e += (*i_hit).energy(t); - m_cell_e->push_back( tot_e ); - m_cell_eta->push_back( ddElement->eta() ); - m_cell_phi->push_back( ddElement->phi() ); - m_cell_radius->push_back( ddElement->r() ); + const CaloDetDescrElement* ddElement = (m_tileID->is_tile_aux(cell_id)) ? 0 : m_tileMgr->get_cell_element(cell_id); + if(ddElement){ + double tot_e = 0.; + double tot_time = 0.; + for (int t=0;t<(i_hit).size();++t) tot_e += (i_hit).energy(t); + for (int t=0;t<(i_hit).size();++t) tot_time += (i_hit).time(t); + h_cell_e->Fill( tot_e ); + h_cell_eta->Fill( ddElement->eta() ); + h_cell_phi->Fill( ddElement->phi() ); + h_cell_radius->Fill( ddElement->z() ); + h_xy->Fill(ddElement->x(),ddElement->y()); + h_zr->Fill(ddElement->r(), ddElement->r()); + h_etaphi->Fill(ddElement->eta(), ddElement->phi()); + + if( m_expert == "on"){ + h_time_e->Fill(tot_time, tot_e); + h_eta_e->Fill(ddElement->eta(), tot_e); + h_phi_e->Fill(ddElement->phi(), tot_e); + h_r_e->Fill(ddElement->r(), tot_e); + } + m_cell_eta->push_back(ddElement->eta()); + m_cell_phi->push_back(ddElement->phi()); + m_cell_e->push_back(tot_e); + m_cell_x->push_back(ddElement->x()); + m_cell_y->push_back(ddElement->y()); + m_cell_z->push_back(ddElement->z()); + m_cell_radius->push_back(ddElement->r()); + m_time->push_back(tot_time); + } } } - if (m_tree) m_tree->Fill(); + //For calibrated hits + std::string LArCalibKey [3] = {"LArCalibrationHitActive", "LArCalibrationHitInactive","LArCalibrationHitDeadMaterial"}; + for(unsigned int j=0; j<3;j++){ + if(m_calib != "on") continue; + const DataHandle<CaloCalibrationHitContainer> iterator; + CaloCalibrationHitContainer::const_iterator hit_i; + if(evtStore()->retrieve(iterator, LArCalibKey[j])==StatusCode::SUCCESS){ + //Not tested + for(auto hit_i : *iterator){ + GeoCaloCalibHit geoHit(*hit_i, LArCalibKey[j]); + if(!geoHit) continue; + const CaloDetDescrElement* Element = geoHit.getDetDescrElement(); + double eta = Element->eta(); + double phi = Element->phi(); + double radius = Element->r(); + double z =Element->z(); + double emEnergy = geoHit.energyEM(); + double nonEmEnergy = geoHit.energyNonEM(); + double invEnergy = geoHit.energyInvisible(); + double escEnergy = geoHit.energyEscaped(); + double totEnergy = geoHit.energyTotal(); + double particleID = (*hit_i).particleID(); + + h_calib_eta->Fill(eta); + h_calib_phi->Fill(phi); + h_calib_rz->Fill(z,radius); + h_calib_etaphi->Fill(eta,phi); + h_calib_eEM->Fill(emEnergy); + h_calib_eNonEM->Fill(nonEmEnergy); + h_calib_eInv->Fill(invEnergy); + h_calib_eEsc->Fill(escEnergy); + h_calib_eTot->Fill(totEnergy); + h_calib_eTotpartID->Fill(particleID, totEnergy); + + m_calib_eta->push_back(eta); + m_calib_phi->push_back(phi); + m_calib_radius->push_back(radius); + m_calib_z->push_back(z); + m_calib_eEM->push_back(emEnergy); + m_calib_eNonEM->push_back(nonEmEnergy); + m_calib_eInv->push_back(invEnergy); + m_calib_eEsc->push_back(escEnergy); + m_calib_eTot->push_back(totEnergy); + m_calib_partID->push_back(particleID); + + } + } + } + if(m_tree) m_tree->Fill(); return StatusCode::SUCCESS; } diff --git a/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..a0e04559e5ab75fabb3002f417ca6da6356b5eec --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/CaloHitAnalysis.h @@ -0,0 +1,93 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef CALO_HIT_ANALYSIS_H +#define CALO_HIT_ANALYSIS_H + +// Base class +#include "AthenaBaseComps/AthAlgorithm.h" + +// Members +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +// STL includes +#include <string> +//#include <vector> + +class TileID; +class TileDetDescrManager; +class TH1; +class TH2; +class TTree; + +class CaloHitAnalysis : public AthAlgorithm { + + public: + + CaloHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~CaloHitAnalysis() {} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Simple variables by Ketevi */ + TH1* h_cell_eta; + TH1* h_cell_phi; + TH1* h_cell_e; + TH1* h_cell_radius; + TH2* h_xy; + TH2* h_zr; + TH2* h_etaphi; + TH2* h_time_e; + TH2* h_eta_e; + TH2* h_phi_e; + TH2* h_r_e; + TH1* h_calib_eta; + TH1* h_calib_phi; + TH2* h_calib_rz; + TH2* h_calib_etaphi; + TH1* h_calib_eEM; + TH1* h_calib_eNonEM; + TH1* h_calib_eInv; + TH1* h_calib_eEsc; + TH1* h_calib_eTot; + TH1* h_calib_eTotpartID; + + std::vector<float>* m_cell_eta; + std::vector<float>* m_cell_phi; + std::vector<float>* m_cell_x; + std::vector<float>* m_cell_y; + std::vector<float>* m_cell_z; + std::vector<float>* m_cell_e; + std::vector<float>* m_cell_radius; + std::vector<float>* m_time; + std::vector<float>* m_calib_eta; + std::vector<float>* m_calib_phi; + std::vector<float>* m_calib_radius; + std::vector<float>* m_calib_z; + std::vector<float>* m_calib_eEM; + std::vector<float>* m_calib_eNonEM; + std::vector<float>* m_calib_eInv; + std::vector<float>* m_calib_eEsc; + std::vector<float>* m_calib_eTot; + std::vector<float>* m_calib_partID; + + + std::string m_expert; + std::string m_calib; + ServiceHandle<ITHistSvc> m_thistSvc; + std::string m_path; + + TTree* m_tree; + std::string m_ntupleFileName; + const TileID * m_tileID; + const TileDetDescrManager * m_tileMgr; + +}; + +#endif // CALO_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..57b312d3454ad69f889c3a27749002e15d845798 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.cxx @@ -0,0 +1,212 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "LucidHitAnalysis.h" + + +#include "LUCID_SimEvent/LUCID_SimHitCollection.h" + + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +LucidHitAnalysis::LucidHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hit_x(0) + , h_hit_y(0) + , h_hit_z(0) + , h_xy(0) + , h_zr(0) + , h_hit_post_x(0) + , h_hit_post_y(0) + , h_hit_post_z(0) + , h_hit_edep(0) + , h_hit_pdgid(0) + , h_hit_pretime(0) + , h_hit_posttime(0) + , h_genvolume(0) + , h_wavelength(0) + , m_hit_x(0) + , m_hit_y(0) + , m_hit_z(0) + , m_hit_post_x(0) + , m_hit_post_y(0) + , m_hit_post_z(0) + , m_hit_edep(0) + , m_hit_pdgid(0) + , m_hit_pretime(0) + , m_hit_posttime(0) + , m_gen_volume(0) + , m_wavelength(0) + , m_tree(0) + , m_ntupleFileName("/LucidHitAnalysis/ntuple/") + , m_path("/LucidHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); + +} + +StatusCode LucidHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing LucidHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + /** Histograms**/ + h_hit_x = new TH1D("h_hit_x", "hit_x", 100,-150.,150.); + h_hit_x->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_x->GetName(), h_hit_x)); + + h_hit_y = new TH1D("h_hit_y", "hit_y", 100,-150.,150.); + h_hit_y->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_y->GetName(), h_hit_y)); + + h_hit_z = new TH1D("h_hit_z", "hit_z", 100,-20000.,20000.); + h_hit_z->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_z->GetName(), h_hit_z)); + + h_xy = new TH2D("h_xy", "hit_xy", 100,-150.,150.,100,-150,150); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D("h_zr", "hit_zr", 100,-20000.,20000.,100,0,250); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_hit_post_x = new TH1D("h_hit_post_x", "hit_post_x", 100,-150.,150.); + h_hit_post_x->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_post_x->GetName(), h_hit_post_x)); + + h_hit_post_y = new TH1D("h_hit_post_y", "hit_post_y", 100,-150,150.); + h_hit_post_y->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_post_y->GetName(), h_hit_post_y)); + + h_hit_post_z = new TH1D("h_hit_post_z", "hit_post_z", 100,-15000,15000.); + h_hit_post_z->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_post_z->GetName(), h_hit_post_z)); + + h_hit_edep = new TH1D("h_hit_edep", "hit_edep", 100,0.,20.); + h_hit_edep->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_edep->GetName(), h_hit_edep)); + + h_hit_pdgid = new TH1D("h_hit_pdgid", "hit_pdgid", 100,0.,7e6); + h_hit_pdgid->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_pdgid->GetName(), h_hit_pdgid)); + + h_hit_pretime = new TH1D("h_hit_pretime", "hit_pretime", 100,0.,100.); + h_hit_pretime->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_pretime->GetName(), h_hit_pretime)); + + h_hit_posttime = new TH1D("h_hit_posttime", "hit_posttime", 100,0.,100.); + h_hit_posttime->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_hit_posttime->GetName(), h_hit_posttime)); + + h_genvolume = new TH1D("h_genvolume", "genvolume", 20,0.,5.); + h_genvolume->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_genvolume->GetName(), h_genvolume)); + + h_wavelength = new TH1D("m_wavelength", "wavelength", 150,0.,800.); + h_wavelength->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_wavelength->GetName(), h_wavelength)); + + + m_tree= new TTree("NtupleLucidHitAnalysis","LucidHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + if (m_tree){ + m_tree->Branch("hit_x", &m_hit_x); + m_tree->Branch("hit_y", &m_hit_y); + m_tree->Branch("hit_z", &m_hit_z); + m_tree->Branch("hit_post_x", &m_hit_post_x); + m_tree->Branch("hit_post_y", &m_hit_post_y); + m_tree->Branch("hit_post_z", &m_hit_post_z); + m_tree->Branch("edep", &m_hit_edep); + m_tree->Branch("pdgid", &m_hit_pdgid); + m_tree->Branch("pretime", &m_hit_pretime); + m_tree->Branch("posttime", &m_hit_posttime); + m_tree->Branch("gen_volume", &m_gen_volume); + m_tree->Branch("wavelength", &m_wavelength); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + return StatusCode::SUCCESS; +} + + + +StatusCode LucidHitAnalysis::execute() { + ATH_MSG_DEBUG( "In LucidHitAnalysis::execute()" ); + + m_hit_x->clear(); + m_hit_y->clear(); + m_hit_z->clear(); + m_hit_post_x->clear(); + m_hit_post_y->clear(); + m_hit_post_z->clear(); + m_hit_edep->clear(); + m_hit_pdgid->clear(); + m_hit_pretime->clear(); + m_hit_posttime->clear(); + m_gen_volume->clear(); + m_wavelength->clear(); + + + const DataHandle<LUCID_SimHitCollection> iter; + if (evtStore()->retrieve(iter)==StatusCode::SUCCESS) { + for(LUCID_SimHitCollection::const_iterator i_hit = (*iter).begin(); + i_hit != (*iter).end();++i_hit){ + double x = i_hit->GetX(); + double y = i_hit->GetY(); + double z = i_hit->GetZ(); + double r = sqrt(x*x+y*y); + h_hit_x->Fill(i_hit->GetX()); + h_hit_y->Fill(i_hit->GetY()); + h_hit_z->Fill(i_hit->GetZ()); + h_xy->Fill(x, y); + h_zr->Fill(z, r); + h_hit_post_x->Fill(i_hit->GetEPX()); + h_hit_post_y->Fill(i_hit->GetEPY()); + h_hit_post_z->Fill(i_hit->GetEPZ()); + h_hit_edep->Fill(i_hit->GetEnergy()); + h_hit_pdgid->Fill(i_hit->GetPdgCode()); + h_hit_pretime->Fill(i_hit->GetPreStepTime()); + h_hit_posttime->Fill(i_hit->GetPostStepTime()); + h_genvolume->Fill(i_hit->GetGenVolume()); + h_wavelength->Fill(i_hit->GetWavelength()); + + + m_hit_x->push_back(i_hit->GetX()); + m_hit_y->push_back(i_hit->GetY()); + m_hit_z->push_back(i_hit->GetZ()); + m_hit_post_x->push_back(i_hit->GetEPX()); + m_hit_post_y->push_back(i_hit->GetEPY()); + m_hit_post_z->push_back(i_hit->GetEPZ()); + m_hit_edep->push_back(i_hit->GetEnergy()); + m_hit_pdgid->push_back(i_hit->GetPdgCode()); + m_hit_pretime->push_back(i_hit->GetPreStepTime()); + m_hit_posttime->push_back(i_hit->GetPostStepTime()); + m_gen_volume->push_back(i_hit->GetGenVolume()); + m_wavelength->push_back(i_hit->GetWavelength()); + } + } // End while hits + + if (m_tree) m_tree->Fill(); + + + + + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..056751fbcbc68d8f623c38711cd83336c9772dc3 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/LucidHitAnalysis.h @@ -0,0 +1,77 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef LUCID_HIT_ANALYSIS_H +#define LUCID_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + +class TH1; +class TH2; +class TTree; + +class LucidHitAnalysis : public AthAlgorithm { + + public: + + LucidHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~LucidHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + std::string m_collection; + /** Some histograms**/ + TH1* h_hit_x; + TH1* h_hit_y; + TH1* h_hit_z; + TH2* h_xy; + TH2* h_zr; + TH1* h_hit_post_x; + TH1* h_hit_post_y; + TH1* h_hit_post_z; + TH1* h_hit_edep; + TH1* h_hit_pdgid; + TH1* h_hit_pretime; + TH1* h_hit_posttime; + TH1* h_genvolume; + TH1* h_wavelength; + + + std::vector<float>* m_hit_x; + std::vector<float>* m_hit_y; + std::vector<float>* m_hit_z; + std::vector<float>* m_hit_post_x; + std::vector<float>* m_hit_post_y; + std::vector<float>* m_hit_post_z; + std::vector<float>* m_hit_edep; + std::vector<float>* m_hit_pdgid; + std::vector<float>* m_hit_pretime; + std::vector<float>* m_hit_posttime; + std::vector<float>* m_gen_volume; + std::vector<float>* m_wavelength; + + + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // LUCID_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..1f7ddd9d3fe849ae48c5ecee1303cded275b5e90 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.cxx @@ -0,0 +1,236 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "MDTHitAnalysis.h" + +// Section of includes for MDT of the Muon Spectrometer tests +#include "GeoAdaptors/GeoMuonHits.h" + +#include "MuonSimEvent/MDTSimHitCollection.h" +#include "MuonSimEvent/MDTSimHit.h" +#include "CLHEP/Vector/LorentzVector.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +MDTHitAnalysis::MDTHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_zr(0) + , h_hits_eta(0) + , h_hits_phi(0) + , h_hits_lx(0) + , h_hits_ly(0) + , h_hits_lz(0) + , h_hits_driftR(0) + , h_hits_time(0) + , h_hits_edep(0) + , h_hits_kine(0) + , h_hits_step(0) + , m_hits_x(0) + , m_hits_y(0) + , m_hits_z(0) + , m_hits_r(0) + , m_hits_eta(0) + , m_hits_phi(0) + , m_hits_lx(0) + , m_hits_ly(0) + , m_hits_lz(0) + , m_hits_driftR(0) + , m_hits_time(0) + , m_hits_edep(0) + , m_hits_kine(0) + , m_hits_step(0) + , m_tree(0) + , m_ntupleFileName("MDTHitAnalysis/ntuple/") + , m_path("/MDTHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); + +} + +StatusCode MDTHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing MDTHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + /** Histograms */ + + h_hits_x = new TH1D("h_hits_mdt_x","hits_x", 100,-25000, 25000); + h_hits_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + + h_hits_y = new TH1D("h_hits_mdt_y", "hits_y", 100,-25000,25000); + h_hits_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + + h_hits_z = new TH1D("h_hits_mdt_z", "hits_z", 100,-45000,45000); + h_hits_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + + h_hits_r = new TH1D("h_hits_mdt_r", "hits_r", 100,4000,26000); + h_hits_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + + h_xy = new TH2D("h_mdt_xy", "xy", 100,-25000.,25000.,100, -25000., 25000.); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D("h_mdt_zr", "zr", 100,-45000.,45000.,100, 4000., 26000.); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_hits_eta = new TH1D("h_hits_mdt_eta", "hits_eta", 100,-3.0,3.0); + h_hits_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_eta->GetName(), h_hits_eta)); + + h_hits_phi = new TH1D("h_hits_mdt_phi", "hits_phi", 100,-3.2,3.2); + h_hits_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_phi->GetName(), h_hits_phi)); + + h_hits_lx = new TH1D("h_hits_mdt_lx","hits_lx", 100,-20, 20); + h_hits_lx->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lx->GetName(), h_hits_lx)); + + h_hits_ly = new TH1D("h_hits_mdt_ly", "hits_ly", 100,-20,20); + h_hits_ly->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ly->GetName(), h_hits_ly)); + + h_hits_lz = new TH1D("h_hits_mdt_lz", "hits_lz", 100,-2000,2000); + h_hits_lz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lz->GetName(), h_hits_lz)); + + h_hits_driftR = new TH1D("h_hits_mdt_driftR", "hits_driftR", 100,0,15); + h_hits_driftR->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_driftR->GetName(), h_hits_driftR)); + + h_hits_time = new TH1D("h_hits_mdt_time","hits_time", 100,0, 150); + h_hits_time->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_time->GetName(), h_hits_time)); + + h_hits_edep = new TH1D("h_hits_mdt_edep", "hits_edep", 100,0,0.2); + h_hits_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_edep->GetName(), h_hits_edep)); + + h_hits_kine = new TH1D("h_hits_mdt_kine", "hits_kine", 100,0,20000); + h_hits_kine->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_kine->GetName(), h_hits_kine)); + + h_hits_step = new TH1D("h_hits_mdt_step", "hits_step", 100,0,100); + h_hits_step->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_step->GetName(), h_hits_step)); + + m_tree= new TTree("NtupleMDTHitAnalysis","MDTHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("x", &m_hits_x); + m_tree->Branch("y", &m_hits_y); + m_tree->Branch("z", &m_hits_z); + m_tree->Branch("r", &m_hits_r); + m_tree->Branch("eta", &m_hits_eta); + m_tree->Branch("phi", &m_hits_phi); + m_tree->Branch("lx", &m_hits_lx); + m_tree->Branch("ly", &m_hits_ly); + m_tree->Branch("lz", &m_hits_lz); + m_tree->Branch("driftR", &m_hits_driftR); + m_tree->Branch("time", &m_hits_time); + m_tree->Branch("edep", &m_hits_edep); + m_tree->Branch("kine", &m_hits_kine); + m_tree->Branch("step", &m_hits_step); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + return StatusCode::SUCCESS; +} + + + +StatusCode MDTHitAnalysis::execute() { + ATH_MSG_DEBUG( "In MDTHitAnalysis::execute()" ); + + m_hits_x->clear(); + m_hits_y->clear(); + m_hits_z->clear(); + m_hits_r->clear(); + m_hits_eta->clear(); + m_hits_phi->clear(); + m_hits_lx->clear(); + m_hits_ly->clear(); + m_hits_lz->clear(); + m_hits_driftR->clear(); + m_hits_time->clear(); + m_hits_edep->clear(); + m_hits_kine->clear(); + m_hits_step->clear(); + + const DataHandle<MDTSimHitCollection> mdt_container; + if (evtStore()->retrieve(mdt_container, "MDT_Hits" )==StatusCode::SUCCESS) { + for(MDTSimHitCollection::const_iterator i_hit = mdt_container->begin(); + i_hit != mdt_container->end();++i_hit){ + //MDTSimHitCollection::const_iterator i_hit; + //for(auto i_hit : *mdt_container){ + GeoMDTHit ghit(*i_hit); + if(!ghit) continue; + Amg::Vector3D p = ghit.getGlobalPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_zr->Fill(p.z(),p.perp()); + h_hits_eta->Fill(p.eta()); + h_hits_phi->Fill(p.phi()); + h_hits_driftR->Fill((*i_hit).driftRadius()); + h_hits_lx->Fill( (*i_hit).localPosition().x()); + h_hits_ly->Fill( (*i_hit).localPosition().y()); + h_hits_lz->Fill( (*i_hit).localPosition().z()); + h_hits_edep->Fill((*i_hit).energyDeposit()); + h_hits_time->Fill((*i_hit).globalTime()); + h_hits_step->Fill((*i_hit).stepLength()); + h_hits_kine->Fill((*i_hit).kineticEnergy()); + + + m_hits_x->push_back(p.x()); + m_hits_y->push_back(p.y()); + m_hits_z->push_back(p.z()); + m_hits_r->push_back(p.perp()); + m_hits_eta->push_back(p.eta()); + m_hits_phi->push_back(p.phi()); + m_hits_driftR->push_back((*i_hit).driftRadius()); + m_hits_lx->push_back( (*i_hit).localPosition().x()); + m_hits_ly->push_back( (*i_hit).localPosition().y()); + m_hits_lz->push_back( (*i_hit).localPosition().z()); + m_hits_edep->push_back((*i_hit).energyDeposit()); + m_hits_time->push_back((*i_hit).globalTime()); + m_hits_step->push_back((*i_hit).stepLength()); + m_hits_kine->push_back((*i_hit).kineticEnergy()); + } + } // End while hits + + if (m_tree) m_tree->Fill(); + + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..85b154dad3711c524a9d9549bb93727fc1a1b0e1 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/MDTHitAnalysis.h @@ -0,0 +1,79 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef MDT_HIT_ANALYSIS_H +#define MDT_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + + +class TH1; +class TH2; +class TTree; + + +class MDTHitAnalysis : public AthAlgorithm { + + public: + + MDTHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~MDTHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_zr; + TH1* h_hits_eta; + TH1* h_hits_phi; + TH1* h_hits_lx; + TH1* h_hits_ly; + TH1* h_hits_lz; + TH1* h_hits_driftR; + TH1* h_hits_time; + TH1* h_hits_edep; + TH1* h_hits_kine; + TH1* h_hits_step; + std::vector<float>* m_hits_x; + std::vector<float>* m_hits_y; + std::vector<float>* m_hits_z; + std::vector<float>* m_hits_r; + std::vector<float>* m_hits_eta; + std::vector<float>* m_hits_phi; + std::vector<float>* m_hits_lx; + std::vector<float>* m_hits_ly; + std::vector<float>* m_hits_lz; + std::vector<float>* m_hits_driftR; + std::vector<float>* m_hits_time; + std::vector<float>* m_hits_edep; + std::vector<float>* m_hits_kine; + std::vector<float>* m_hits_step; + + TTree * m_tree; + std::string m_ntupleFileName; + + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // MDT_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..b83c04c183da51f37ba9d459e6794163912e1399 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.cxx @@ -0,0 +1,225 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "RPCHitAnalysis.h" + +// Section of includes for the RPC of the Muon Spectrometer tests +#include "GeoAdaptors/GeoMuonHits.h" + +#include "MuonSimEvent/RPCSimHitCollection.h" +#include "MuonSimEvent/RPCSimHit.h" +#include "CLHEP/Vector/LorentzVector.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +RPCHitAnalysis::RPCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_zr(0) + , h_hits_eta(0) + , h_hits_phi(0) + , h_hits_lx(0) + , h_hits_ly(0) + , h_hits_lz(0) + , h_hits_time(0) + , h_hits_edep(0) + , h_hits_kine(0) + , h_hits_step(0) + , m_hits_x(0) + , m_hits_y(0) + , m_hits_z(0) + , m_hits_r(0) + , m_hits_eta(0) + , m_hits_phi(0) + , m_hits_lx(0) + , m_hits_ly(0) + , m_hits_lz(0) + , m_hits_time(0) + , m_hits_edep(0) + , m_hits_kine(0) + , m_hits_step(0) + , m_tree(0) + , m_ntupleFileName("/RPCHitAnalysis/ntuple/") + , m_path("/RPCHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); + +} + +StatusCode RPCHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing RPCHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + /** Histograms**/ + h_hits_x = new TH1D("h_hits_rpc_x","hits_x", 100,-11000, 11000); + h_hits_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + + h_hits_y = new TH1D("h_hits_rpc_y", "hits_y", 100,-11000,11000); + h_hits_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + + h_hits_z = new TH1D("h_hits_rpc_z", "hits_z", 100,-12500, 12500); + h_hits_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + + h_hits_r = new TH1D("h_hits_rpc_r", "hits_r", 100,6000,14000); + h_hits_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + + h_xy = new TH2D("h_rpc_xy", "xy", 100,-11000.,11000.,100, -11000., 11000.); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D("m_rpc_zr", "zr", 100,-12500.,12500.,100, 6000., 14000.); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_hits_eta = new TH1D("h_hits_rpc_eta", "hits_eta", 100,-1.5,1.5); + h_hits_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_eta->GetName(), h_hits_eta)); + + h_hits_phi = new TH1D("h_hits_rpc_phi", "hits_phi", 100,-3.2,3.2); + h_hits_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_phi->GetName(), h_hits_phi)); + + h_hits_lx = new TH1D("h_hits_rpc_lx","hits_lx", 100,-10, 10); + h_hits_lx->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lx->GetName(), h_hits_lx)); + + h_hits_ly = new TH1D("h_hits_rpc_ly", "hits_ly", 100,-1500,1500); + h_hits_ly->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ly->GetName(), h_hits_ly)); + + h_hits_lz = new TH1D("h_hits_rpc_lz", "hits_lz", 100,-600,600); + h_hits_lz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lz->GetName(), h_hits_lz)); + + + h_hits_time = new TH1D("h_hits_rpc_time","hits_time", 100,0, 120); + h_hits_time->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_time->GetName(), h_hits_time)); + + h_hits_edep = new TH1D("h_hits_rpc_edep", "hits_edep", 100,0,0.15); + h_hits_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_edep->GetName(), h_hits_edep)); + + h_hits_kine = new TH1D("h_hits_rpc_kine", "hits_kine", 100,0,500); + h_hits_kine->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_kine->GetName(), h_hits_kine)); + + h_hits_step = new TH1D("h_hits_rpc_step", "hits_step", 100,0,25); + h_hits_step->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_step->GetName(), h_hits_step)); + + + m_tree= new TTree("NtupleRPCHitAnalysis","RPCHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("x", &m_hits_x); + m_tree->Branch("y", &m_hits_y); + m_tree->Branch("z", &m_hits_z); + m_tree->Branch("r", &m_hits_r); + m_tree->Branch("eta", &m_hits_eta); + m_tree->Branch("phi", &m_hits_phi); + m_tree->Branch("lx", &m_hits_lx); + m_tree->Branch("ly", &m_hits_ly); + m_tree->Branch("lz", &m_hits_lz); + m_tree->Branch("time", &m_hits_time); + m_tree->Branch("edep", &m_hits_edep); + m_tree->Branch("kine", &m_hits_kine); + m_tree->Branch("step", &m_hits_step); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + return StatusCode::SUCCESS; +} + + + +StatusCode RPCHitAnalysis::execute() { + ATH_MSG_DEBUG( "In RPCHitAnalysis::execute()" ); + + m_hits_x->clear(); + m_hits_y->clear(); + m_hits_z->clear(); + m_hits_r->clear(); + m_hits_eta->clear(); + m_hits_phi->clear(); + m_hits_lx->clear(); + m_hits_ly->clear(); + m_hits_lz->clear(); + m_hits_time->clear(); + m_hits_edep->clear(); + m_hits_kine->clear(); + m_hits_step->clear(); + + const DataHandle<RPCSimHitCollection> rpc_container; + if (evtStore()->retrieve(rpc_container, "RPC_Hits" )==StatusCode::SUCCESS) { + for(RPCSimHitCollection::const_iterator i_hit = rpc_container->begin(); + i_hit != rpc_container->end();++i_hit){ + //RPCSimHitCollection::const_iterator i_hit; + //for(auto i_hit : *rpc_container){ + GeoRPCHit ghit(*i_hit); + if(!ghit) continue; + Amg::Vector3D p = ghit.getGlobalPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_zr->Fill(p.z(),p.perp()); + h_hits_eta->Fill(p.eta()); + h_hits_phi->Fill(p.phi()); + h_hits_lx->Fill( (*i_hit).localPosition().x()); + h_hits_ly->Fill( (*i_hit).localPosition().y()); + h_hits_lz->Fill( (*i_hit).localPosition().z()); + h_hits_edep->Fill((*i_hit).energyDeposit()); + h_hits_time->Fill((*i_hit).globalTime()); + h_hits_step->Fill((*i_hit).stepLength()); + h_hits_kine->Fill((*i_hit).kineticEnergy()); + + m_hits_x->push_back(p.x()); + m_hits_y->push_back(p.y()); + m_hits_z->push_back(p.z()); + m_hits_r->push_back(p.perp()); + m_hits_eta->push_back(p.eta()); + m_hits_phi->push_back(p.phi()); + m_hits_lx->push_back( (*i_hit).localPosition().x()); + m_hits_ly->push_back( (*i_hit).localPosition().y()); + m_hits_lz->push_back( (*i_hit).localPosition().z()); + m_hits_edep->push_back((*i_hit).energyDeposit()); + m_hits_time->push_back((*i_hit).globalTime()); + m_hits_step->push_back((*i_hit).stepLength()); + m_hits_kine->push_back((*i_hit).kineticEnergy()); + } + } // End while hits + + + if (m_tree) m_tree->Fill(); + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..fadb854c2127018a69936086aba16a5841417d52 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/RPCHitAnalysis.h @@ -0,0 +1,77 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef RPC_HIT_ANALYSIS_H +#define RPC_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + + +class TH1; +class TH2; +class TTree; + +class RPCHitAnalysis : public AthAlgorithm { + + public: + + RPCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~RPCHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_zr; + TH1* h_hits_eta; + TH1* h_hits_phi; + TH1* h_hits_lx; + TH1* h_hits_ly; + TH1* h_hits_lz; + TH1* h_hits_time; + TH1* h_hits_edep; + TH1* h_hits_kine; + TH1* h_hits_step; + + + std::vector<float>* m_hits_x; + std::vector<float>* m_hits_y; + std::vector<float>* m_hits_z; + std::vector<float>* m_hits_r; + std::vector<float>* m_hits_eta; + std::vector<float>* m_hits_phi; + std::vector<float>* m_hits_lx; + std::vector<float>* m_hits_ly; + std::vector<float>* m_hits_lz; + std::vector<float>* m_hits_time; + std::vector<float>* m_hits_edep; + std::vector<float>* m_hits_kine; + std::vector<float>* m_hits_step; + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // RPC_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx new file mode 100644 index 0000000000000000000000000000000000000000..8c757c7893afc2c1ee7a44037b905f510cebe616 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx @@ -0,0 +1,232 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "SiHitAnalysis.h" + +// Section of includes for Pixel and SCT tests +#include "InDetSimEvent/SiHitCollection.h" +#include "GeoAdaptors/GeoSiHit.h" + +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> + +//m_collection: "PixelHits", "SCT_Hits", "BCMHits" and "BLMHits" +SiHitAnalysis::SiHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , m_collection("BCMHits") + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_zr(0) + , h_hits_time(0) + , h_hits_eloss(0) + , h_hits_step(0) + , h_hits_barcode(0) + , h_time_eloss(0) + , h_z_eloss(0) + , h_r_eloss(0) + , m_hits_x(0) + , m_hits_y(0) + , m_hits_z(0) + , m_hits_r(0) + , m_hits_time(0) + , m_hits_eloss(0) + , m_hits_step(0) + , m_hits_barcode(0) + , m_tree(0) + , m_ntupleFileName("/ntuples/") + , m_expert("off") + , m_path("/SiHitAnalysis/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("CollectionName", m_collection="BCMHits"); + declareProperty("HistPath", m_path); + declareProperty("ExpertMode", m_expert = "off"); + declareProperty("NtupleFileName", m_ntupleFileName); +} + +StatusCode SiHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing SiHitAnalysis" ); + + std::string detName("Pixel"); + //initialise pixel or SCT variables + if(m_collection=="PixelHits"){ + detName = "Pixel"; + }else if (m_collection=="SCT_Hits"){ + detName = "SCT"; + }else if(m_collection=="BCMHits"){ + detName = "BCM"; + }else if(m_collection=="BLMHits"){ + detName = "BLM"; + }else{ + ATH_MSG_ERROR("SiHitsAnalysis for "<< name()<<"not supported!!! \n"); + return StatusCode::FAILURE; + } + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + /** Histograms**/ + float bin_down = -600; + float bin_up = 600; + float radius_up = 600; + float radius_down = 200; + if(detName=="Pixel"){ + bin_down = -170; + bin_up = 170; + radius_up = 170; + radius_down = 0; + } + + h_hits_x = new TH1D(("h_"+detName+"_x").c_str(),("h_"+detName+"_x").c_str(), 100,bin_down, bin_up); + h_hits_x->StatOverflows(); + + h_hits_y = new TH1D(("h_"+detName+"_y").c_str(), ("h_"+detName+"_y").c_str(), 100,bin_down,bin_up); + h_hits_y->StatOverflows(); + + h_hits_z = new TH1D(("h_"+detName+"_z").c_str(), ("h_"+detName+"_z").c_str(), 200,-1200,1200); + h_hits_z->StatOverflows(); + + h_hits_r = new TH1D(("h_"+detName+"_r").c_str(), ("h_"+detName+"_r").c_str(), 100,radius_down,radius_up); + h_hits_r->StatOverflows(); + + h_xy = new TH2D(("h_"+detName+"_xy").c_str(), ("h_"+detName+"_xy").c_str(), 100,bin_down,bin_up,100, bin_down, bin_up); + h_xy->StatOverflows(); + + h_zr = new TH2D(("h_"+detName+"_zr").c_str(), ("h_"+detName+"_zr").c_str(), 100,-1200,1200.,100, radius_down, radius_up); + h_zr->StatOverflows(); + + h_hits_time = new TH1D(("h_"+detName+"_time").c_str(), ("h_"+detName+"_time").c_str(), 100,0,500); + h_hits_time->StatOverflows(); + + h_hits_eloss = new TH1D(("h_"+detName+"_eloss").c_str(), ("h_"+detName+"_eloss").c_str(), 100,0,50); + h_hits_eloss->StatOverflows(); + + h_hits_step = new TH1D(("h_"+detName+"_step").c_str(), ("h_"+detName+"_step").c_str(), 100,0,50); + h_hits_step->StatOverflows(); + + h_hits_barcode = new TH1D(("h_"+detName+"_barcode").c_str(), ("h_"+detName+"_barcode").c_str(), 200,0,250000); + h_hits_barcode->StatOverflows(); + + h_time_eloss = new TH2D(("h_"+detName+"_time_eloss").c_str(), ("h_"+detName+" Eloss vs. time").c_str(),100, 0,500,100,0,50); + h_time_eloss->StatOverflows(); + + h_z_eloss = new TH2D(("h_"+detName+"_z_eloss").c_str(), ("h_"+detName+" Eloss vs. z").c_str(),100, -1200,1200,100,0,50); + h_z_eloss->StatOverflows(); + + h_r_eloss = new TH2D(("h_"+detName+"_r_eloss").c_str(), ("h_"+detName+ " Eloss vs. r").c_str(),100, radius_down,radius_down,100,0,50); + h_r_eloss->StatOverflows(); + + + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + CHECK(m_thistSvc->regHist(m_path + h_hits_time->GetName(), h_hits_time)); + CHECK(m_thistSvc->regHist(m_path + h_hits_eloss->GetName(), h_hits_eloss)); + CHECK(m_thistSvc->regHist(m_path + h_hits_step->GetName(), h_hits_step)); + CHECK(m_thistSvc->regHist(m_path + h_hits_barcode->GetName(), h_hits_barcode)); + + + //To be filled only when the expert mode is on. + if(m_expert == "on") + { + CHECK(m_thistSvc->regHist(m_path + h_time_eloss->GetName(), h_time_eloss)); + CHECK(m_thistSvc->regHist(m_path + h_z_eloss->GetName(), h_z_eloss)); + CHECK(m_thistSvc->regHist(m_path + h_r_eloss->GetName(), h_r_eloss)); + } + + /**ntuple**/ + CHECK(m_thistSvc.retrieve()); + m_tree= new TTree(detName.c_str(), detName.c_str()); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"+detName; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch((detName+"_x").c_str(), &m_hits_x); + m_tree->Branch((detName+"_y").c_str(), &m_hits_y); + m_tree->Branch((detName+"_z").c_str(), &m_hits_z); + m_tree->Branch((detName+"_r").c_str(), &m_hits_r); + m_tree->Branch((detName+"_time").c_str(), &m_hits_time); + m_tree->Branch((detName+"_eloss").c_str(), &m_hits_eloss); + m_tree->Branch((detName+"_step").c_str(), &m_hits_step); + m_tree->Branch((detName+"_barcode").c_str(), &m_hits_barcode); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + return StatusCode::SUCCESS; +} + + + +StatusCode SiHitAnalysis::execute() { + ATH_MSG_DEBUG( "In SiHitAnalysis::execute()" ); + m_hits_x->clear(); + m_hits_y->clear(); + m_hits_z->clear(); + m_hits_r->clear(); + m_hits_time->clear(); + m_hits_eloss->clear(); + m_hits_step->clear(); + m_hits_barcode->clear(); + + const DataHandle<SiHitCollection> p_collection; + if (evtStore()->retrieve(p_collection, m_collection )==StatusCode::SUCCESS) { + for(SiHitConstIterator i_hit = p_collection->begin(); + i_hit != p_collection->end();++i_hit){ + GeoSiHit ghit(*i_hit); + HepGeom::Point3D<double> p = ghit.getGlobalPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_zr->Fill(p.z(),p.perp()); + h_hits_eloss->Fill(i_hit->energyLoss()); + h_hits_time->Fill(i_hit->meanTime()); + double step_length=(i_hit->localStartPosition()-i_hit->localEndPosition()).mag(); + h_hits_step->Fill(step_length); + h_hits_barcode->Fill(i_hit->particleLink().barcode()); + + if(m_expert == "on") + { + h_time_eloss->Fill(i_hit->meanTime(), i_hit->energyLoss()); + if(i_hit->getBarrelEndcap()==0) { + h_z_eloss->Fill(p.z(), i_hit->energyLoss()); + }else{ + h_r_eloss->Fill(p.perp(), i_hit->energyLoss()); + } + } + + m_hits_x->push_back(p.x()); + m_hits_y->push_back(p.y()); + m_hits_z->push_back(p.z()); + m_hits_r->push_back(p.perp()); + m_hits_eloss->push_back(i_hit->energyLoss()); + m_hits_time->push_back(i_hit->meanTime()); + m_hits_step->push_back(step_length); + m_hits_barcode->push_back(i_hit->particleLink().barcode()); + } // End while hits + } // End statuscode success upon retrieval of hits + if (m_tree) m_tree->Fill(); + + + return StatusCode::SUCCESS; +} + diff --git a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..29fdc72de8e1852d763dfefe9ed2c2d22d8ef289 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h @@ -0,0 +1,70 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef SI_HIT_ANALYSIS_H +#define SI_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" + +class TH1; +class TH2; +class TTree; + +class SiHitAnalysis : public AthAlgorithm { + + public: + + SiHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~SiHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + std::string m_collection; + /** Some variables**/ + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_zr; + TH1* h_hits_time; + TH1* h_hits_eloss; + TH1* h_hits_step; + TH1* h_hits_barcode; + TH2* h_time_eloss; + TH2* h_z_eloss; + TH2* h_r_eloss; + + std::vector<float>* m_hits_x; + std::vector<float>* m_hits_y; + std::vector<float>* m_hits_z; + std::vector<float>* m_hits_r; + std::vector<float>* m_hits_time; + std::vector<float>* m_hits_eloss; + std::vector<float>* m_hits_step; + std::vector<float>* m_hits_barcode; + + TTree* m_tree; + std::string m_ntupleFileName; + + std::string m_expert; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // SI_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..437baef10b8bbd0e5440bb93b24a4177142f00a5 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.cxx @@ -0,0 +1,255 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "TGCHitAnalysis.h" + +// Section of includes for TGC of the Muon Spectrometer tests +#include "GeoAdaptors/GeoMuonHits.h" + +#include "MuonSimEvent/TGCSimHitCollection.h" +#include "MuonSimEvent/TGCSimHit.h" +#include "CLHEP/Vector/LorentzVector.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +TGCHitAnalysis::TGCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_rz(0) + , h_hits_eta(0) + , h_hits_phi(0) + , h_hits_lx(0) + , h_hits_ly(0) + , h_hits_lz(0) + , h_hits_dcx(0) + , h_hits_dcy(0) + , h_hits_dcz(0) + , h_hits_time(0) + , h_hits_edep(0) + , h_hits_kine(0) + , h_hits_step(0) + , m_hits_x(0) + , m_hits_y(0) + , m_hits_z(0) + , m_hits_r(0) + , m_hits_eta(0) + , m_hits_phi(0) + , m_hits_lx(0) + , m_hits_ly(0) + , m_hits_lz(0) + , m_hits_dcx(0) + , m_hits_dcy(0) + , m_hits_dcz(0) + , m_hits_time(0) + , m_hits_edep(0) + , m_hits_kine(0) + , m_hits_step(0) + , m_tree(0) + , m_ntupleFileName("/TGCHitAnalysis/ntuple/") + , m_path("/TGCHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); +} + +StatusCode TGCHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing TGCHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + /** Histograms**/ + h_hits_x = new TH1D("h_hits_tgc_x","hits_x", 100,-5000, 5000); + h_hits_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + + h_hits_y = new TH1D("h_hits_tgc_y", "hits_y", 100,-5000,5000); + h_hits_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + + h_hits_z = new TH1D("h_hits_tgc_z", "hits_z", 100,-12000,12000); + h_hits_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + + h_hits_r = new TH1D("h_hits_tgc_r", "hits_r", 100,2000,10000); + h_hits_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + + h_xy = new TH2D("h_tgc_xy", "xy", 100,-5000.,5000.,100, -5000., 5000.); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_rz = new TH2D("h_tgc_rz", "rz", 100,2000.,10000.,100, -12000., 12000.); + h_rz->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_rz->GetName(), h_rz)); + + + h_hits_eta = new TH1D("h_hits_tgc_eta", "hits_eta", 100,-10.0,10.0); + h_hits_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_eta->GetName(), h_hits_eta)); + + h_hits_phi = new TH1D("h_hits_tgc_phi", "hits_phi", 100,-3.2,3.2); + h_hits_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_phi->GetName(), h_hits_phi)); + + h_hits_lx = new TH1D("h_hits_tgc_lx","hits_lx", 100,-800, 800); + h_hits_lx->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lx->GetName(), h_hits_lx)); + + h_hits_ly = new TH1D("h_hits_tgc_ly", "hits_ly", 100,-800,800); + h_hits_ly->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_ly->GetName(), h_hits_ly)); + + h_hits_lz = new TH1D("h_hits_tgc_lz", "hits_lz", 100,-800,800); + h_hits_lz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_lz->GetName(), h_hits_lz)); + + h_hits_dcx = new TH1D("h_hits_tgc_dcx","hits_dcx", 100,-1, 1); + h_hits_dcx->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_dcx->GetName(), h_hits_dcx)); + + h_hits_dcy = new TH1D("h_hits_tgc_dcy", "hits_dcy", 100,-1,1); + h_hits_dcy->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_dcy->GetName(), h_hits_dcy)); + + h_hits_dcz = new TH1D("h_hits_tgc_dcz", "hits_dcz", 100,-1,1); + h_hits_dcz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_dcz->GetName(), h_hits_dcz)); + + h_hits_time = new TH1D("h_hits_tgc_time","hits_time", 100,0, 250); + h_hits_time->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_time->GetName(), h_hits_time)); + + h_hits_edep = new TH1D("h_hits_tgc_edep", "hits_edep", 100,0,0.5); + h_hits_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_edep->GetName(), h_hits_edep)); + + h_hits_kine = new TH1D("h_hits_tgc_kine", "hits_kine", 100,0,1000); + h_hits_kine->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_kine->GetName(), h_hits_kine)); + + h_hits_step = new TH1D("h_hits_tgc_step", "hits_step", 100,0,50); + h_hits_step->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_step->GetName(), h_hits_step)); + + m_tree= new TTree("NtupleTGCHitAnalysis","TGCHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("x", &m_hits_x); + m_tree->Branch("y", &m_hits_y); + m_tree->Branch("z", &m_hits_z); + m_tree->Branch("r", &m_hits_r); + m_tree->Branch("eta", &m_hits_eta); + m_tree->Branch("phi", &m_hits_phi); + m_tree->Branch("lx", &m_hits_lx); + m_tree->Branch("ly", &m_hits_ly); + m_tree->Branch("lz", &m_hits_lz); + m_tree->Branch("dcx", &m_hits_dcx); + m_tree->Branch("dcy", &m_hits_dcy); + m_tree->Branch("dcz", &m_hits_dcz); + m_tree->Branch("time", &m_hits_time); + m_tree->Branch("edep", &m_hits_edep); + m_tree->Branch("kine", &m_hits_kine); + m_tree->Branch("step", &m_hits_step); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + + return StatusCode::SUCCESS; +} + + + +StatusCode TGCHitAnalysis::execute() { + ATH_MSG_DEBUG( "In TGCHitAnalysis::execute()" ); + + m_hits_x->clear(); + m_hits_y->clear(); + m_hits_z->clear(); + m_hits_r->clear(); + m_hits_eta->clear(); + m_hits_phi->clear(); + m_hits_lx->clear(); + m_hits_ly->clear(); + m_hits_lz->clear(); + m_hits_dcx->clear(); + m_hits_dcy->clear(); + m_hits_dcz->clear(); + m_hits_time->clear(); + m_hits_edep->clear(); + m_hits_kine->clear(); + m_hits_step->clear(); + + const DataHandle<TGCSimHitCollection> tgc_container; + if (evtStore()->retrieve(tgc_container, "TGC_Hits" )==StatusCode::SUCCESS) { + for(TGCSimHitCollection::const_iterator i_hit = tgc_container->begin(); + i_hit != tgc_container->end();++i_hit){ + //TGCSimHitCollection::const_iterator i_hit; + //for(auto i_hit : *tgc_container){ + GeoTGCHit ghit(*i_hit); + if(!ghit) continue; + Amg::Vector3D p = ghit.getGlobalPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_rz->Fill(p.perp(), p.z()); + h_hits_eta->Fill(p.eta()); + h_hits_phi->Fill(p.phi()); + h_hits_lx->Fill( (*i_hit).localPosition().x()); + h_hits_ly->Fill( (*i_hit).localPosition().y()); + h_hits_lz->Fill( (*i_hit).localPosition().z()); + h_hits_dcx->Fill((*i_hit).localDireCos().x()); + h_hits_dcy->Fill((*i_hit).localDireCos().y()); + h_hits_dcz->Fill((*i_hit).localDireCos().z()); + h_hits_edep->Fill((*i_hit).energyDeposit()); + h_hits_time->Fill((*i_hit).globalTime()); + h_hits_step->Fill((*i_hit).stepLength()); + h_hits_kine->Fill((*i_hit).kineticEnergy()); + + m_hits_x->push_back(p.x()); + m_hits_y->push_back(p.y()); + m_hits_z->push_back(p.z()); + m_hits_r->push_back(p.perp()); + m_hits_eta->push_back(p.eta()); + m_hits_phi->push_back(p.phi()); + m_hits_lx->push_back( (*i_hit).localPosition().x()); + m_hits_ly->push_back( (*i_hit).localPosition().y()); + m_hits_lz->push_back( (*i_hit).localPosition().z()); + m_hits_dcx->push_back((*i_hit).localDireCos().x()); + m_hits_dcy->push_back((*i_hit).localDireCos().y()); + m_hits_dcz->push_back((*i_hit).localDireCos().z()); + m_hits_edep->push_back((*i_hit).energyDeposit()); + m_hits_time->push_back((*i_hit).globalTime()); + m_hits_step->push_back((*i_hit).stepLength()); + m_hits_kine->push_back((*i_hit).kineticEnergy()); + } + } // End while hits + + if (m_tree) m_tree->Fill(); + + + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..74a73d450312bf0dc6c073df970cd3f3b8f80ef1 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TGCHitAnalysis.h @@ -0,0 +1,85 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef TGC_HIT_ANALYSIS_H +#define TGC_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + +class TH1; +class TH2; +class TTree; + +class TGCHitAnalysis : public AthAlgorithm { + + public: + + TGCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~TGCHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + + + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_rz; + TH1* h_hits_eta; + TH1* h_hits_phi; + TH1* h_hits_lx; + TH1* h_hits_ly; + TH1* h_hits_lz; + TH1* h_hits_dcx; + TH1* h_hits_dcy; + TH1* h_hits_dcz; + TH1* h_hits_time; + TH1* h_hits_edep; + TH1* h_hits_kine; + TH1* h_hits_step; + + std::vector<float>* m_hits_x; + std::vector<float>* m_hits_y; + std::vector<float>* m_hits_z; + std::vector<float>* m_hits_r; + std::vector<float>* m_hits_eta; + std::vector<float>* m_hits_phi; + std::vector<float>* m_hits_lx; + std::vector<float>* m_hits_ly; + std::vector<float>* m_hits_lz; + std::vector<float>* m_hits_dcx; + std::vector<float>* m_hits_dcy; + std::vector<float>* m_hits_dcz; + std::vector<float>* m_hits_time; + std::vector<float>* m_hits_edep; + std::vector<float>* m_hits_kine; + std::vector<float>* m_hits_step; + + + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // TGC_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..e5fdacefa9374b5726c2b48cd0260f4011ffe2d5 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.cxx @@ -0,0 +1,206 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "TRTHitAnalysis.h" + +// Section of includes for TRT tests +#include "InDetSimEvent/TRTUncompressedHitCollection.h" +#include "GeoAdaptors/GeoTRTUncompressedHit.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> + +TRTHitAnalysis::TRTHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_TRT_x(0) + , h_TRT_y(0) + , h_TRT_z(0) + , h_TRT_r(0) + , h_TRT_xy(0) + , h_TRT_zr(0) + , h_TRT_time_photons(0) + , h_TRT_time_nonphotons(0) + , h_TRT_edep_photons(0) + , h_TRT_edep_nonphotons(0) + , h_TRT_kine_photons(0) + , h_TRT_kine_nonphotons(0) + , h_TRT_barcode(0) + , m_TRT_x(0) + , m_TRT_y(0) + , m_TRT_z(0) + , m_TRT_r(0) + , m_TRT_time_photons(0) + , m_TRT_time_nonphotons(0) + , m_TRT_edep_photons(0) + , m_TRT_edep_nonphotons(0) + , m_TRT_kine_photons(0) + , m_TRT_kine_nonphotons(0) + , m_TRT_barcode(0) + , m_tree(0) + , m_path("/histos/TRTHitAnalysis/") + , m_ntupleFileName("/ntuples/TRTHitAnalysis/") + , m_thistSvc("THistSvc", name) + +{ + declareProperty("HistPath", m_path); + declareProperty("NtupleFileName", m_ntupleFileName); +} + +StatusCode TRTHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing TRTHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + /** Histograms **/ + + h_TRT_x = new TH1D("h_TRT_x","hits_x", 100,-1100, 1100); + h_TRT_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_x->GetName(), h_TRT_x)); + + h_TRT_y = new TH1D("h_TRT_y", "hits_y", 100,-1100,1100); + h_TRT_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_y->GetName(), h_TRT_y)); + + h_TRT_z = new TH1D("h_TRT_z", "hits_z", 100,-3000,3000); + h_TRT_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_z->GetName(), h_TRT_z)); + + h_TRT_r = new TH1D("h_TRT_r", "hits_r", 100,500,1100); + h_TRT_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_r->GetName(), h_TRT_r)); + + h_TRT_xy = new TH2D("h_TRT_xy", "xy", 100,-1100.,1100.,100, -1100., 1100.); + h_TRT_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_TRT_xy->GetName(), h_TRT_xy)); + + h_TRT_zr = new TH2D("h_TRT_zr", "zr", 100,-3000,3000,100,500.,1100.); + h_TRT_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_TRT_zr->GetName(), h_TRT_zr)); + + h_TRT_time_photons = new TH1D("h_TRT_time_photons", "hits_time_photons", 100,0,500); + h_TRT_time_photons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_time_photons->GetName(), h_TRT_time_photons)); + + h_TRT_time_nonphotons = new TH1D("h_TRT_time_nonphotons", "hits_time_nonphotons", 100,0,500); + h_TRT_time_nonphotons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_time_nonphotons->GetName(), h_TRT_time_nonphotons)); + + h_TRT_edep_photons = new TH1D("h_TRT_edep_photons", "hits_edep_photons", 100,0,500); + h_TRT_edep_photons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_edep_photons->GetName(), h_TRT_edep_photons)); + + h_TRT_edep_nonphotons = new TH1D("m_TRT_edep_nonphotons", "hits_edep_nonphotons", 100,0,500); + h_TRT_edep_nonphotons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_edep_nonphotons->GetName(), h_TRT_edep_nonphotons)); + + h_TRT_kine_photons = new TH1D("h_TRT_kine_photons", "hits_kine_photons", 100,0,2); + h_TRT_kine_photons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_kine_photons->GetName(), h_TRT_kine_photons)); + + h_TRT_kine_nonphotons = new TH1D("h_TRT_kine_nonphotons", "hits_kine_nonphotons", 100,0,1000); + h_TRT_kine_nonphotons->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_kine_nonphotons->GetName(), h_TRT_kine_nonphotons)); + + h_TRT_barcode = new TH1D("h_TRT_barcode", "hits_barcode", 100,-500,300000); + h_TRT_barcode->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_TRT_barcode->GetName(), h_TRT_barcode)); + + + /*ntuples*/ + + m_tree= new TTree("TRTHitNtuple","TRTHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + if (m_tree){ + m_tree->Branch("x", &m_TRT_x); + m_tree->Branch("y", &m_TRT_y); + m_tree->Branch("z", &m_TRT_z); + m_tree->Branch("r", &m_TRT_r); + m_tree->Branch("time_photons", &m_TRT_time_photons); + m_tree->Branch("time_nonphotons", &m_TRT_time_nonphotons); + m_tree->Branch("EnergyDeposit_photons", &m_TRT_edep_photons); + m_tree->Branch("EnergyDeposit_nonphotons", &m_TRT_edep_nonphotons); + m_tree->Branch("KineticEnergy_photons", &m_TRT_kine_photons); + m_tree->Branch("KineticEnergy_nonphotons", &m_TRT_kine_nonphotons); + m_tree->Branch("barcode", &m_TRT_barcode); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + + return StatusCode::SUCCESS; +} + + + +StatusCode TRTHitAnalysis::execute() { + ATH_MSG_DEBUG( "In TRTHitAnalysis::execute()" ); + + m_TRT_x->clear(); + m_TRT_y->clear(); + m_TRT_z->clear(); + m_TRT_r->clear(); + m_TRT_time_photons->clear(); + m_TRT_time_nonphotons->clear(); + m_TRT_edep_photons->clear(); + m_TRT_edep_nonphotons->clear(); + m_TRT_kine_photons->clear(); + m_TRT_kine_nonphotons->clear(); + m_TRT_barcode->clear(); + + const DataHandle<TRTUncompressedHitCollection> p_collection; + if (evtStore()->retrieve(p_collection, "TRTUncompressedHits" )==StatusCode::SUCCESS) { + for(TRTUncompressedHitConstIter i_hit = p_collection->begin(); + i_hit != p_collection->end();++i_hit){ + GeoTRTUncompressedHit ghit(*i_hit); + HepGeom::Point3D<double> p = ghit.getGlobalPosition(); + + h_TRT_x->Fill(p.x()); + h_TRT_y->Fill(p.y()); + h_TRT_z->Fill(p.z()); + h_TRT_r->Fill(p.perp()); + h_TRT_xy->Fill(p.x(), p.y()); + h_TRT_zr->Fill(p.z(),sqrt(pow(p.x(),2)+pow(p.y(),2))); + h_TRT_barcode->Fill(i_hit->particleLink().barcode()); + + m_TRT_x->push_back(p.x()); + m_TRT_y->push_back(p.y()); + m_TRT_z->push_back(p.z()); + m_TRT_r->push_back(p.perp()); + m_TRT_barcode->push_back(i_hit->particleLink().barcode()); + int particleId(i_hit->GetParticleEncoding()); + if(particleId == 22 || static_cast<int>(abs(particleId)/100000)==41 || static_cast<int>(abs(particleId)/10000000) == 1){ + h_TRT_time_photons->Fill(i_hit->GetGlobalTime()); + h_TRT_edep_photons->Fill(i_hit->GetEnergyDeposit()); + h_TRT_kine_photons->Fill(i_hit->GetKineticEnergy()); + m_TRT_time_photons->push_back(i_hit->GetGlobalTime()); + m_TRT_edep_photons->push_back(i_hit->GetEnergyDeposit()); + m_TRT_kine_photons->push_back(i_hit->GetKineticEnergy()); + }else{ + h_TRT_time_nonphotons->Fill(i_hit->GetGlobalTime()); + h_TRT_edep_nonphotons->Fill(i_hit->GetEnergyDeposit()); + h_TRT_kine_nonphotons->Fill(i_hit->GetKineticEnergy()); + + m_TRT_time_nonphotons->push_back(i_hit->GetGlobalTime()); + m_TRT_edep_nonphotons->push_back(i_hit->GetEnergyDeposit()); + m_TRT_kine_nonphotons->push_back(i_hit->GetKineticEnergy()); + } + } // End while hits + } + if (m_tree) m_tree->Fill(); + + + + return StatusCode::SUCCESS; +} + diff --git a/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..7c34cfe043cdb34c2f6f8f843b906f72ca32ff59 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TRTHitAnalysis.h @@ -0,0 +1,73 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef TRT_HIT_ANALYSIS_H +#define TRT_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + +class TH1; +class TH2; +class TTree; + +class TRTHitAnalysis : public AthAlgorithm { + + public: + + TRTHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~TRTHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_TRT_x; + TH1* h_TRT_y; + TH1* h_TRT_z; + TH1* h_TRT_r; + TH2* h_TRT_xy; + TH2* h_TRT_zr; + TH1* h_TRT_time_photons; + TH1* h_TRT_time_nonphotons; + TH1* h_TRT_edep_photons; + TH1* h_TRT_edep_nonphotons; + TH1* h_TRT_kine_photons; + TH1* h_TRT_kine_nonphotons; + TH1* h_TRT_barcode; + + + std::vector<float>* m_TRT_x; + std::vector<float>* m_TRT_y; + std::vector<float>* m_TRT_z; + std::vector<float>* m_TRT_r; + std::vector<float>* m_TRT_time_photons; + std::vector<float>* m_TRT_time_nonphotons; + std::vector<float>* m_TRT_edep_photons; + std::vector<float>* m_TRT_edep_nonphotons; + std::vector<float>* m_TRT_kine_photons; + std::vector<float>* m_TRT_kine_nonphotons; + std::vector<float>* m_TRT_barcode; + + TTree * m_tree; + std::string m_path; + std::string m_ntupleFileName; + ServiceHandle<ITHistSvc> m_thistSvc; + + +}; + +#endif // TRT_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..92e687af1d599a7859e1de4445ddf53f13b7ffe3 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.cxx @@ -0,0 +1,264 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "TrackRecordAnalysis.h" + +// Section of includes for TrackRecord tests + +#include "TrackRecord/TrackRecord.h" +#include "TrackRecord/TrackRecordCollection.h" +#include "CLHEP/Vector/LorentzVector.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +TrackRecordAnalysis::TrackRecordAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_hits_x(0) + , h_hits_y(0) + , h_hits_z(0) + , h_hits_r(0) + , h_xy(0) + , h_zr(0) + , h_hits_eta(0) + , h_hits_phi(0) + , h_hits_px(0) + , h_hits_py(0) + , h_hits_pz(0) + , h_hits_pt(0) + , h_time(0) + , h_edep(0) + , h_pdg(0) + , m_x(0) + , m_y(0) + , m_z(0) + , m_r(0) + , m_eta(0) + , m_phi(0) + , m_px(0) + , m_py(0) + , m_pz(0) + , m_pt(0) + , m_time(0) + , m_edep(0) + , m_pdg(0) + , m_tree(0) + , m_ntupleFileName("/TrackRecordAnalysis/ntuple/") + , m_path("/TrackRecordAnalysis/histos/") + , m_thistSvc("THistSvc", name) + , m_collection("CaloEntryLayer") +{ + declareProperty("NtupleFileName", m_ntupleFileName); + declareProperty("HistPath", m_path); + declareProperty("CollectionName", m_collection="CaloEntryLayer"); +} + +StatusCode TrackRecordAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing TrackRecordAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + std::string detName("CaloEntry"); + std::string ntupName("TrackRecordCaloEntry"); + if(m_collection=="CaloEntryLayer"){ + detName="CaloEntry"; + ntupName="TrackRecordCaloEntry"; + } else if(m_collection=="MuonEntryLayer"){ + detName="MuonEntry"; + ntupName="TrackRecordMuonEntry"; + } else if(m_collection=="MuonExitLayer"){ + detName="MuonExit"; + ntupName="TrackRecordMuonExit"; + }else{ + ATH_MSG_ERROR("TrackRecordAnalysis for "<< name() << "not supported !!! \n"); + return StatusCode::FAILURE; + + } + std::cout<<"Name "<<name()<<std::endl; + //Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + float x_down = -5000; + float x_up = 5000; + float radius = 4500; + float eta_down = -5.8; + float eta_up = 5.8; + float z_down = -7000; + float z_up = 7000; + if(detName =="CaloEntry"){ + x_down = -1180; + x_up = 1180; + radius = 1200; + eta_down = - 5.6; + eta_up = 5.6; + z_down = -3700; + z_up = 3700; + + } + + /** Histograms**/ + h_hits_x = new TH1D((detName+"_x").c_str(),"hits_x", 100,x_down, x_up); + h_hits_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_x->GetName(), h_hits_x)); + + h_hits_y = new TH1D((detName+"_y").c_str(), "hits_y", 100,x_down,x_up); + h_hits_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_y->GetName(), h_hits_y)); + + h_hits_z = new TH1D((detName+"_z").c_str(), "hits_z", 100,z_down,z_up); + h_hits_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_z->GetName(), h_hits_z)); + + h_hits_r = new TH1D((detName+"_r").c_str(), "hits_r", 100,0,radius); + h_hits_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_r->GetName(), h_hits_r)); + + h_xy = new TH2D((detName+"_xy").c_str(), "xy", 100,x_down,x_up,100, x_down, x_up); + h_xy->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_xy->GetName(), h_xy)); + + h_zr = new TH2D((detName+"_zr").c_str(), "zr", 100,z_down,z_up,100, 0., radius); + h_zr->StatOverflows(); + CHECK(m_thistSvc->regHist( m_path+h_zr->GetName(), h_zr)); + + h_hits_eta = new TH1D((detName+"_eta").c_str(), "hits_eta", 100,eta_down,eta_up); + h_hits_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_eta->GetName(), h_hits_eta)); + + h_hits_phi = new TH1D((detName+"_phi").c_str(), "hits_phi", 100,-3.1416,3.1416); + h_hits_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_phi->GetName(), h_hits_phi)); + + h_hits_px = new TH1D((detName+"_px").c_str(),"Px", 100, -1500, 1500); + h_hits_px->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_px->GetName(), h_hits_px)); + + h_hits_py = new TH1D((detName+"_py").c_str(), "Py", 100,-1500,1500); + h_hits_py->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_py->GetName(), h_hits_py)); + + h_hits_pz = new TH1D((detName+"_pz").c_str(), "Pz", 100,-1500,1500); + h_hits_pz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_pz->GetName(), h_hits_pz)); + + h_hits_pt = new TH1D((detName+"_pt").c_str(), "hits_pt", 100,0,2500); + h_hits_pt->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_hits_pt->GetName(), h_hits_pt)); + + h_time = new TH1D((detName+"_time").c_str(),"time", 100,0, 250); + h_time->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_time->GetName(), h_time)); + + h_edep = new TH1D((detName+"_edep").c_str(), "energy", 100,0,7500); + h_edep->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_edep->GetName(), h_edep)); + + h_pdg = new TH1D((detName+"_pdg").c_str(), "pdg", 100,-1000,1000); + h_pdg->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_pdg->GetName(), h_pdg)); + + + m_tree= new TTree(ntupName.c_str(),ntupName.c_str()); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"+detName; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("x", &m_x); + m_tree->Branch("y", &m_y); + m_tree->Branch("z", &m_z); + m_tree->Branch("r", &m_r); + m_tree->Branch("eta", &m_eta); + m_tree->Branch("phi", &m_phi); + m_tree->Branch("px", &m_px); + m_tree->Branch("py", &m_py); + m_tree->Branch("pz", &m_pz); + m_tree->Branch("pt", &m_pt); + m_tree->Branch("time", &m_time); + m_tree->Branch("energy", &m_edep); + m_tree->Branch("pdg", &m_pdg); + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + + return StatusCode::SUCCESS; +} + + + +StatusCode TrackRecordAnalysis::execute() { + ATH_MSG_DEBUG( "In TrackRecordAnalysis::execute()" ); + + m_x->clear(); + m_y->clear(); + m_z->clear(); + m_r->clear(); + m_eta->clear(); + m_phi->clear(); + m_px->clear(); + m_py->clear(); + m_pz->clear(); + m_pt->clear(); + m_time->clear(); + m_edep->clear(); + m_pdg->clear(); + + const DataHandle<TrackRecordCollection> TRcoll; + if (evtStore()->retrieve(TRcoll, m_collection )==StatusCode::SUCCESS) { + for(TrackRecordCollection::const_iterator track = TRcoll->begin(); + track != TRcoll->end();++track){ + //TrackRecordCollection::const_iterator track; + //for(auto track : *TRcoll){ + std::cout<<"Entra en el loop"<<std::endl; + CLHEP::Hep3Vector p =(*track).GetPosition(); + h_hits_x->Fill(p.x()); + h_hits_y->Fill(p.y()); + h_hits_z->Fill(p.z()); + h_hits_r->Fill(p.perp()); + h_xy->Fill(p.x(), p.y()); + h_zr->Fill(p.z(),p.perp()); + h_hits_eta->Fill(p.eta()); + h_hits_phi->Fill(p.phi()); + CLHEP::Hep3Vector mom = (*track).GetMomentum(); + h_hits_px->Fill( mom.x()); + h_hits_py->Fill( mom.y()); + h_hits_pz->Fill( mom.z()); + h_hits_pt->Fill( mom.perp()); + h_edep->Fill((*track).GetEnergy()); + h_time->Fill((*track).GetTime()); + h_pdg->Fill((*track).GetPDGCode()); + + m_x->push_back(p.x()); + m_y->push_back(p.y()); + m_z->push_back(p.z()); + m_r->push_back(p.perp()); + m_eta->push_back(p.eta()); + m_phi->push_back(p.phi()); + m_px->push_back( mom.x()); + m_py->push_back( mom.y()); + m_pz->push_back( mom.z()); + m_pt->push_back( mom.perp()); + m_edep->push_back((*track).GetEnergy()); + m_time->push_back((*track).GetTime()); + m_pdg->push_back((*track).GetPDGCode()); + } + } // End while hits + + + if (m_tree) m_tree->Fill(); + + + + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.h b/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..9f530143fde1b4695f4e4139edb67d72c2183d0b --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TrackRecordAnalysis.h @@ -0,0 +1,80 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef TRACK_RECORD_ANALYSIS_H +#define TRACK_RECORD_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + + + +class TH1; +class TH2; +class TTree; + + +class TrackRecordAnalysis : public AthAlgorithm { + + public: + + TrackRecordAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~TrackRecordAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + + /** Some variables**/ + TH1* h_hits_x; + TH1* h_hits_y; + TH1* h_hits_z; + TH1* h_hits_r; + TH2* h_xy; + TH2* h_zr; + TH1* h_hits_eta; + TH1* h_hits_phi; + TH1* h_hits_px; + TH1* h_hits_py; + TH1* h_hits_pz; + TH1* h_hits_pt; + TH1* h_time; + TH1* h_edep; + TH1* h_pdg; + + + + std::vector<float>* m_x; + std::vector<float>* m_y; + std::vector<float>* m_z; + std::vector<float>* m_r; + std::vector<float>* m_eta; + std::vector<float>* m_phi; + std::vector<float>* m_px; + std::vector<float>* m_py; + std::vector<float>* m_pz; + std::vector<float>* m_pt; + std::vector<float>* m_time; + std::vector<float>* m_edep; + std::vector<float>* m_pdg; + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + std::string m_collection; +}; + +#endif // TRACK_RECORD_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..36025ca049e38c4eca204129de9f74d55a1b5d5b --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.cxx @@ -0,0 +1,358 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "TruthHitAnalysis.h" + +// Section of includes for Truth tests +#include "HepMC/GenEvent.h" +#include "GeneratorObjects/McEventCollection.h" +#include "EventInfo/EventInfo.h" +#include "EventInfo/EventID.h" + +#include "TH1.h" +#include "TTree.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> + +TruthHitAnalysis::TruthHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_n_vert(0) + , h_n_part(0) + , h_n_vert_prim(0) + , h_n_part_prim(0) + , h_n_vert_sec(0) + , h_n_part_sec(0) + , h_vtx_x(0) + , h_vtx_y(0) + , h_vtx_z(0) + , h_vtx_r(0) + , h_vtx_prim_xy(0) + , h_vtx_prim_zr(0) + , h_vtx_sec_xy(0) + , h_vtx_sec_zr(0) + , h_n_generations(0) + , h_truth_px(0) + , h_truth_py(0) + , h_truth_pz(0) + , h_truth_pt(0) + , h_truth_eta(0) + , h_truth_phi(0) + , h_barcode(0) + , h_part_status(0) + , h_part_pdgid(0) + , h_part_pdgid_sec(0) + , h_part_eta(0) + , h_part_phi(0) + , h_part_p(0) + , m_vtx_x(0) + , m_vtx_y(0) + , m_vtx_z(0) + , m_vtx_r(0) + , m_vtx_barcode(0) + , m_truth_px(0) + , m_truth_py(0) + , m_truth_pz(0) + , m_truth_pt(0) + , m_truth_eta(0) + , m_truth_phi(0) + , m_barcode(0) + , m_status(0) + , m_pdgid(0) + , m_tree(0) + , m_ntupleFileName("/TruthHitAnalysis/ntuples/") + , m_path("/TruthHitAnalysis/histos/") + , m_thistSvc("THistSvc", name) + +{ + declareProperty("HistPath", m_path); + declareProperty("NtupleFileName", m_ntupleFileName); +} + +StatusCode TruthHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing TruthHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + + /** histograms declaration */ + h_n_vert = new TH1D("h_n_vert","n_vert", 100,200, 1500); + h_n_vert->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_vert->GetName(), h_n_vert)); + + h_n_part = new TH1D("h_n_part","n_part", 100,1000, 10000); + h_n_part->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_part->GetName(), h_n_part)); + + h_n_vert_prim = new TH1D("h_n_vert_prim","n_vert prim", 100,0, 1000); + h_n_vert_prim->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_vert_prim->GetName(), h_n_vert_prim)); + + + h_n_part_prim = new TH1D("h_n_part_prim","n_part prim", 100,200, 1500); + h_n_part_prim->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_part_prim->GetName(), h_n_part_prim)); + + h_n_vert_sec = new TH1D("h_n_vert_sec","n_vert sec", 100,0, 1000); + h_n_vert_sec->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_vert_sec->GetName(), h_n_vert_sec)); + + h_n_part_sec = new TH1D("h_n_part_sec","n_part sec", 100,0, 5000); + h_n_part_sec->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_part_sec->GetName(), h_n_part_sec)); + + h_vtx_x = new TH1D("h_vtx_x","vtx_x", 100,-1300, 1300); + h_vtx_x->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_x->GetName(), h_vtx_x)); + + h_vtx_y = new TH1D("h_vtx_y","vtx_y", 100,-1200, 1200); + h_vtx_y->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_y->GetName(), h_vtx_y)); + + h_vtx_z = new TH1D("h_vtx_z","vtx_z", 100,-5000, 5000); + h_vtx_z->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_z->GetName(), h_vtx_z)); + + h_vtx_r = new TH1D("h_vtx_r","vtx_r", 100,0, 1160); + h_vtx_r->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_r->GetName(), h_vtx_r)); + + h_vtx_prim_xy = new TH2D("h_vtx_prim_xy","vtx_prim_xy", 100,-100, 100, 100,-100, 100); + h_vtx_prim_xy->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_prim_xy->GetName(), h_vtx_prim_xy)); + + h_vtx_prim_zr = new TH2D("h_vtx_prim_zr","vtx_prim_zr", 100,-1500, 1500, 100,0, 150); + h_vtx_prim_zr->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_prim_zr->GetName(), h_vtx_prim_zr)); + + h_vtx_sec_xy = new TH2D("h_vtx_sec_xy","vtx_sec_xy", 100,-1200, 1200, 100,-1200, 1200); + h_vtx_sec_xy->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_sec_xy->GetName(), h_vtx_sec_xy)); + + h_vtx_sec_zr = new TH2D("h_vtx_sec_zr","vtx_sec_zr", 100,-6000, 6000, 100,0, 1160); + h_vtx_sec_zr->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_vtx_sec_zr->GetName(), h_vtx_sec_zr)); + + + h_n_generations = new TH1D("h_n_generations","h_generations", 100,0, 25); + h_n_generations->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_n_generations->GetName(), h_n_generations)); + + h_truth_px = new TH1D("h_turht_px","truth_px", 100,0, 4000); + h_truth_px->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_px->GetName(), h_truth_px)); + + h_truth_py = new TH1D("h_turht_py","truth_py", 100,0, 4000); + h_truth_py->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_py->GetName(), h_truth_py)); + + h_truth_pz = new TH1D("h_truth_pz","truth_pz", 100,0, 4000); + h_truth_pz->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_pz->GetName(), h_truth_pz)); + + h_truth_pt = new TH1D("h_truth_pt","truth_pt", 100,0, 4000); + h_truth_pt->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_pt->GetName(), h_truth_pt)); + + h_truth_eta = new TH1D("h_truth_eta","truth_eta", 50,-10, 10); + h_truth_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_eta->GetName(), h_truth_eta)); + + h_truth_phi = new TH1D("h_truth_phi","truth_phi", 25,-3.1416, 3.1416); + h_truth_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_truth_phi->GetName(), h_truth_phi)); + + h_barcode = new TH1D("h_truth_barcode","truth_barcode", 100,0, 300000); + h_barcode->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_barcode->GetName(), h_barcode)); + + h_part_status = new TH1D("h_part_status","part status", 100,0,50); + h_part_status->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_status->GetName(), h_part_status)); + + h_part_pdgid = new TH1D("h_part_pdgid","part pdgid", 100,-5000, 5000); + h_part_pdgid->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_pdgid->GetName(), h_part_pdgid)); + + h_part_pdgid_sec = new TH1D("h_part_pdgid_sec","part pdgid sec", 100,-5000, 5000); + h_part_pdgid_sec->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_pdgid_sec->GetName(), h_part_pdgid_sec)); + + h_part_eta = new TH1D("h_part_eta","part eta", 100,-10, 10); + h_part_eta->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_eta->GetName(), h_part_eta)); + + h_part_phi = new TH1D("h_part_phi","part phi", 100,-3.2, 3.2); + h_part_phi->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_phi->GetName(), h_part_phi)); + + h_part_p = new TH1D("h_part_p","part p", 100,0, 5000); + h_part_p->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_part_p->GetName(), h_part_p)); + + + m_tree= new TTree("TruthHitNtuple ","TruthHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + + + /** now add branches and leaves to the tree */ + if (m_tree){ + + m_tree->Branch("vtx_x", &m_vtx_x); + m_tree->Branch("vtx_y", &m_vtx_y); + m_tree->Branch("vtx_z", &m_vtx_z); + m_tree->Branch("vtx_r", &m_vtx_r); + m_tree->Branch("vtx_barcode", &m_vtx_barcode); + m_tree->Branch("truth_px", &m_truth_px); + m_tree->Branch("truth_py", &m_truth_py); + m_tree->Branch("truth_pz", &m_truth_pz); + m_tree->Branch("truth_pt", &m_truth_pt); + m_tree->Branch("truth_eta", &m_truth_eta); + m_tree->Branch("truth_phi", &m_truth_phi); + m_tree->Branch("barcode", &m_barcode); + m_tree->Branch("status", &m_status); + m_tree->Branch("pdg_id", &m_pdgid); + + }else{ + ATH_MSG_ERROR("No tree found!"); + } + + + return StatusCode::SUCCESS; +} + + + +StatusCode TruthHitAnalysis::execute() { + ATH_MSG_DEBUG( "In TruthHitAnalysis::execute()" ); + + m_vtx_x->clear(); + m_vtx_y->clear(); + m_vtx_z->clear(); + m_vtx_r->clear(); + m_vtx_barcode->clear(); + m_truth_px->clear(); + m_truth_py->clear(); + m_truth_pz->clear(); + m_truth_pt->clear(); + m_truth_eta->clear(); + m_truth_phi->clear(); + m_barcode->clear(); + m_status->clear(); + m_pdgid->clear(); + + + const DataHandle<EventInfo> event; + if (!evtStore()->retrieve(event, "McEventInfo" ).isSuccess()) + return StatusCode::FAILURE; + const DataHandle<McEventCollection> mcCollection; + if(evtStore()->retrieve(mcCollection,"TruthEvent")==StatusCode::SUCCESS){ + McEventCollection::const_iterator currentGenEventIter = mcCollection->begin(); + if(currentGenEventIter != mcCollection->end()){ + int nvtx = 0; + int nvtx_sec=0; + for(HepMC::GenEvent::vertex_const_iterator vtx=(*currentGenEventIter)->vertices_begin(); vtx!=(*currentGenEventIter)->vertices_end();++vtx){ + + + double x = (*vtx)->position().x(); + double y = (*vtx)->position().y(); + double z = (*vtx)->position().z(); + double r = sqrt(x*x+y*y); + h_vtx_x->Fill(x); + h_vtx_y->Fill(y); + h_vtx_r->Fill(r); + h_vtx_z->Fill(z); + + + int bcode = (*vtx)->barcode(); + m_vtx_x->push_back(x); + m_vtx_y->push_back(y); + m_vtx_r->push_back(r); + m_vtx_z->push_back(z); + m_vtx_barcode->push_back(bcode); + + if((*vtx)->barcode()>-20000){ + h_vtx_prim_xy->Fill(x,y); + h_vtx_prim_zr->Fill(z,r); + ++nvtx; + } else{ + h_vtx_sec_xy->Fill(x,y); + h_vtx_sec_zr->Fill(z,r); + ++nvtx_sec; + } + + } //End iteration over vertices + + h_n_vert->Fill(nvtx+nvtx_sec); + h_n_vert_prim->Fill(nvtx); + h_n_vert_sec->Fill(nvtx_sec); + int npart_prim=0; + int npart_sec=0; + HepMC::GenEvent::particle_const_iterator currentGenParticleIter; + for(currentGenParticleIter=(*currentGenEventIter)->particles_begin(); currentGenParticleIter!=(*currentGenEventIter)->particles_end(); ++currentGenParticleIter){ + + + const HepMC::FourVector mom=(*currentGenParticleIter)->momentum(); + + + h_truth_px->Fill(mom.x()); + h_truth_py->Fill(mom.y()); + h_truth_pz->Fill(mom.z()); + h_truth_pt->Fill(mom.perp()); + h_truth_eta->Fill(mom.eta()); + h_truth_phi->Fill(mom.phi()); + h_barcode->Fill((*currentGenParticleIter)->barcode()); + + h_part_status->Fill((*currentGenParticleIter)->status()); + + m_truth_px->push_back(mom.x()); + m_truth_py->push_back(mom.y()); + m_truth_pz->push_back(mom.z()); + m_truth_pt->push_back(mom.perp()); + m_truth_eta->push_back(mom.eta()); + m_truth_phi->push_back(mom.phi()); + m_barcode->push_back((*currentGenParticleIter)->barcode()); + m_status->push_back((*currentGenParticleIter)->status()); + int pdg = (*currentGenParticleIter)->pdg_id(); + m_pdgid->push_back(pdg); + + if((*currentGenParticleIter)->barcode()<200000){ + h_part_pdgid->Fill(pdg); + h_part_p->Fill(mom.rho()); + h_part_eta->Fill(mom.eta()); + h_part_phi->Fill(mom.phi()); + ++npart_prim; + if((*currentGenParticleIter)->barcode()<10000){ + h_n_generations->Fill(0); + }else{ + h_n_generations->Fill(1); + } + } //End barcode <200000 + else{ + h_part_pdgid_sec->Fill(pdg); + ++npart_sec; + const int gen = (*currentGenParticleIter)->barcode()/1000000 +2; + h_n_generations ->Fill(gen); + } } // End iteration over particles + + h_n_part_prim->Fill(npart_prim); + h_n_part_sec->Fill(npart_sec); + h_n_part->Fill(npart_prim+npart_sec); + + } // End mcCollection + } // End statuscode success upon retrieval of events + + + if (m_tree) m_tree->Fill(); + + return StatusCode::SUCCESS; +} + diff --git a/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..e868d555a3d09d1d57ae541141adbe7c23d875d3 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/TruthHitAnalysis.h @@ -0,0 +1,92 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef TRUTH_HIT_ANALYSIS_H +#define TRUTH_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + + +#include <string> +#include <vector> +#include "TH1.h" +#include "TH2.h" +#include "TTree.h" + +class TH1; +class TH2; +class TTree; + +class TruthHitAnalysis : public AthAlgorithm { + + public: + + TruthHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~TruthHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_n_vert; + TH1* h_n_part; + TH1* h_n_vert_prim; + TH1* h_n_part_prim; + TH1* h_n_vert_sec; + TH1* h_n_part_sec; + TH1* h_vtx_x; + TH1* h_vtx_y; + TH1* h_vtx_z; + TH1* h_vtx_r; + TH2* h_vtx_prim_xy; + TH2* h_vtx_prim_zr; + TH2* h_vtx_sec_xy; + TH2* h_vtx_sec_zr; + TH1* h_n_generations; + TH1* h_truth_px; + TH1* h_truth_py; + TH1* h_truth_pz; + TH1* h_truth_pt; + TH1* h_truth_eta; + TH1* h_truth_phi; + TH1* h_barcode; + TH1* h_part_status; + TH1* h_part_pdgid; + TH1* h_part_pdgid_sec; + TH1* h_part_eta; + TH1* h_part_phi; + TH1* h_part_p; + + + std::vector<float>* m_vtx_x; + std::vector<float>* m_vtx_y; + std::vector<float>* m_vtx_z; + std::vector<float>* m_vtx_r; + std::vector<float>* m_vtx_barcode; + std::vector<float>* m_truth_px; + std::vector<float>* m_truth_py; + std::vector<float>* m_truth_pz; + std::vector<float>* m_truth_pt; + std::vector<float>* m_truth_eta; + std::vector<float>* m_truth_phi; + std::vector<float>* m_barcode; + std::vector<float>* m_status; + std::vector<float>* m_pdgid; + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + + + +}; + +#endif // TRUTH_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.cxx new file mode 100755 index 0000000000000000000000000000000000000000..7b838ddb1266acf312544ec72cbc17520d4545ab --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.cxx @@ -0,0 +1,218 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "ZDCHitAnalysis.h" + + +#include "ZDC_SimEvent/ZDC_SimStripHit_Collection.h" +#include "ZDC_SimEvent/ZDC_SimStripHit.h" +#include "ZDC_SimEvent/ZDC_SimPixelHit_Collection.h" +#include "ZDC_SimEvent/ZDC_SimPixelHit.h" + + +#include "TH1.h" +#include "TString.h" + + +#include <algorithm> +#include <math.h> +#include <functional> +#include <iostream> +#include <stdio.h> + +ZDCHitAnalysis::ZDCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator) + : AthAlgorithm(name, pSvcLocator) + , h_zdc_sidea_0(0) + , h_zdc_sidea_1(0) + , h_zdc_sidea_2(0) + , h_zdc_sidea_3(0) + , h_zdc_sidec_0(0) + , h_zdc_sidec_1(0) + , h_zdc_sidec_2(0) + , h_zdc_sidec_3(0) + , m_zdc_strip_side(0) + , m_zdc_strip_mod(0) + , m_zdc_strip_energy(0) + , m_zdc_pix_side(0) + , m_zdc_pix_mod(0) + , m_zdc_pix_energy(0) + , m_path("/ZDCHitAnalysis/") + , m_thistSvc("THistSvc", name) +{ + declareProperty("HistPath", m_path); +} + +StatusCode ZDCHitAnalysis::initialize() { + ATH_MSG_DEBUG( "Initializing ZDCHitAnalysis" ); + + // Grab the Ntuple and histogramming service for the tree + CHECK(m_thistSvc.retrieve()); + + /** Histograms*/ + + h_zdc_sidea_0 = new TH1D("m_edep_side_a0","edep_side_a0", 100,0, 1000); + h_zdc_sidea_0->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidea_0->GetName(), h_zdc_sidea_0)); + + h_zdc_sidea_1 = new TH1D("m_edep_side_a1","edep_side_a1", 100,0, 1000); + h_zdc_sidea_1->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidea_1->GetName(), h_zdc_sidea_1)); + + h_zdc_sidea_2 = new TH1D("m_edep_side_a2","edep_side_a2", 100,0, 1000); + h_zdc_sidea_2->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidea_2->GetName(), h_zdc_sidea_2)); + + h_zdc_sidea_3 = new TH1D("m_edep_side_a3","edep_side_a3", 100,0, 1000); + h_zdc_sidea_3->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidea_3->GetName(), h_zdc_sidea_3)); + + h_zdc_sidec_0 = new TH1D("m_edep_side_c0","edep_side_c0", 100,0, 1000); + h_zdc_sidec_0->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidec_0->GetName(), h_zdc_sidec_0)); + + h_zdc_sidec_1 = new TH1D("m_edep_side_c1","edep_side_c1", 100,0, 1000); + h_zdc_sidec_1->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidec_1->GetName(), h_zdc_sidec_1)); + + h_zdc_sidec_2 = new TH1D("m_edep_side_c2","edep_side_c2", 100,0, 1000); + h_zdc_sidec_2->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidec_2->GetName(), h_zdc_sidec_2)); + + h_zdc_sidec_3 = new TH1D("m_edep_side_c3","edep_side_c3", 100,0, 1000); + h_zdc_sidec_3->StatOverflows(); + CHECK(m_thistSvc->regHist(m_path + h_zdc_sidec_3->GetName(), h_zdc_sidec_3)); + + + m_tree= new TTree("NtupleZDCHitAnalysis","ZDCHitAna"); + std::string fullNtupleName = "/"+m_ntupleFileName+"/"; + CHECK(m_thistSvc->regTree(fullNtupleName,m_tree)); + + /** now add branches and leaves to the tree */ + if (m_tree){ + m_tree->Branch("strip_side", &m_zdc_strip_side); + m_tree->Branch("strip_mode", &m_zdc_strip_mod); + m_tree->Branch("strip_energy", &m_zdc_strip_energy); + + m_tree->Branch("pix_side", &m_zdc_pix_side); + m_tree->Branch("pix_mode", &m_zdc_pix_mod); + m_tree->Branch("pix_energy", &m_zdc_pix_energy); + + }else{ + ATH_MSG_ERROR("No tree found!"); + } + return StatusCode::SUCCESS; +} + + + +StatusCode ZDCHitAnalysis::execute() { + ATH_MSG_DEBUG( "In ZDCHitAnalysis::execute()" ); + + m_zdc_strip_side->clear(); + m_zdc_strip_mod->clear(); + m_zdc_strip_energy->clear(); + m_zdc_pix_side->clear(); + m_zdc_pix_mod->clear(); + m_zdc_pix_energy->clear(); + double ene_strip = -999.; + int side_strip = -1; + int mod_strip = -1; + double ene_pix = -999.; + int side_pix = -1; + int mod_pix = -1; + + ZDC_SimStripHit_ConstIterator striphi; + const DataHandle<ZDC_SimStripHit_Collection> stripiter; + CHECK(evtStore()->retrieve(stripiter,"ZDC_SimStripHit_Collection")); + for(striphi=(*stripiter).begin(); striphi != (*stripiter).end();++striphi){ + ZDC_SimStripHit ghit(*striphi); + ene_strip=ghit.GetEdep(); + side_strip = ghit.GetSide(); + mod_strip = ghit.GetMod(); + + if(side_pix==1){ + switch(mod_strip){ + case 0: + h_zdc_sidea_0->Fill(ene_strip); + break; + case 1: + h_zdc_sidea_1->Fill(ene_strip); + break; + case 2: + h_zdc_sidea_2->Fill(ene_strip); + break; + case 3: + h_zdc_sidea_3->Fill(ene_strip); + break; + } + }else{ + switch (mod_strip){ + case 0: + h_zdc_sidec_0->Fill(ene_strip); + break; + case 1: + h_zdc_sidec_1->Fill(ene_strip); + break; + case 2: + h_zdc_sidec_2->Fill(ene_strip); + break; + case 3: + h_zdc_sidec_3->Fill(ene_strip); + break; + } + } + } + + ZDC_SimPixelHit_ConstIterator pixelhi; + const DataHandle<ZDC_SimPixelHit_Collection> pixeliter; + CHECK(evtStore()->retrieve(pixeliter,"ZDC_SimPixelHit_Collection")); + for(pixelhi=(*pixeliter).begin(); pixelhi != (*pixeliter).end();++pixelhi){ + ZDC_SimPixelHit ghit(*pixelhi); + ene_pix=ghit.GetEdep(); + side_pix = ghit.GetSide(); + mod_pix = ghit.GetMod(); + + if(side_pix==1){ + switch(mod_pix){ + case 0: + h_zdc_sidea_0->Fill(ene_pix); + break; + case 1: + h_zdc_sidea_1->Fill(ene_pix); + break; + case 2: + h_zdc_sidea_2->Fill(ene_pix); + break; + case 3: + h_zdc_sidea_3->Fill(ene_pix); + break; + } + }else{ + switch (mod_pix){ + case 0: + h_zdc_sidec_0->Fill(ene_pix); + break; + case 1: + h_zdc_sidec_1->Fill(ene_pix); + break; + case 2: + h_zdc_sidec_2->Fill(ene_pix); + break; + case 3: + h_zdc_sidec_3->Fill(ene_pix); + break; + } + } + + m_zdc_strip_side->push_back(side_strip); + m_zdc_strip_mod->push_back(mod_strip); + m_zdc_strip_energy->push_back(ene_strip); + m_zdc_pix_side->push_back(side_pix); + m_zdc_pix_mod->push_back(mod_pix); + m_zdc_pix_energy->push_back(ene_pix); + } + + if (m_tree) m_tree->Fill(); + return StatusCode::SUCCESS; +} diff --git a/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.h new file mode 100755 index 0000000000000000000000000000000000000000..cf7678b4e77f7e72e20a8687084d36fd68b6afd2 --- /dev/null +++ b/Simulation/Tools/HitAnalysis/src/ZDCHitAnalysis.h @@ -0,0 +1,61 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef ZDC_HIT_ANALYSIS_H +#define ZDC_HIT_ANALYSIS_H + +#include "AthenaBaseComps/AthAlgorithm.h" + +#include "GaudiKernel/ServiceHandle.h" +#include "GaudiKernel/ITHistSvc.h" + +#include <string> +#include <vector> +#include "TH1.h" +#include "TTree.h" + + +class TH1; +class TTree; + + +class ZDCHitAnalysis : public AthAlgorithm { + + public: + + ZDCHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator); + ~ZDCHitAnalysis(){} + + virtual StatusCode initialize(); + virtual StatusCode execute(); + + private: + + /** Some variables**/ + TH1* h_zdc_sidea_0; + TH1* h_zdc_sidea_1; + TH1* h_zdc_sidea_2; + TH1* h_zdc_sidea_3; + TH1* h_zdc_sidec_0; + TH1* h_zdc_sidec_1; + TH1* h_zdc_sidec_2; + TH1* h_zdc_sidec_3; + + std::vector<int>* m_zdc_strip_side; + std::vector<int>* m_zdc_strip_mod; + std::vector<double>* m_zdc_strip_energy; + + std::vector<int>* m_zdc_pix_side; + std::vector<int>* m_zdc_pix_mod; + std::vector<double>* m_zdc_pix_energy; + + TTree * m_tree; + std::string m_ntupleFileName; + std::string m_path; + ServiceHandle<ITHistSvc> m_thistSvc; + +}; + +#endif // ZDC_HIT_ANALYSIS_H + diff --git a/Simulation/Tools/HitAnalysis/src/components/HitAnalysis_entries.cxx b/Simulation/Tools/HitAnalysis/src/components/HitAnalysis_entries.cxx index 0b28107422c7765358dbdc36425f3c87caed8bc7..ee89af53ce7dcd141d5329f44748de814851eefc 100755 --- a/Simulation/Tools/HitAnalysis/src/components/HitAnalysis_entries.cxx +++ b/Simulation/Tools/HitAnalysis/src/components/HitAnalysis_entries.cxx @@ -1,10 +1,48 @@ -#include "HitAnalysis/CaloHitAnalysis.h" - #include "GaudiKernel/DeclareFactoryEntries.h" +#include "../SiHitAnalysis.h" +#include "../CaloHitAnalysis.h" +#include "../TRTHitAnalysis.h" +#include "../TGCHitAnalysis.h" +#include "../RPCHitAnalysis.h" +#include "../MDTHitAnalysis.h" +#include "../CSCHitAnalysis.h" +#include "../TruthHitAnalysis.h" +#include "../TrackRecordAnalysis.h" +#include "../ALFAHitAnalysis.h" +#include "../LucidHitAnalysis.h" +#include "../ZDCHitAnalysis.h" +#include "../AFPHitAnalysis.h" + DECLARE_ALGORITHM_FACTORY( CaloHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( SiHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( TRTHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( TGCHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( RPCHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( MDTHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( CSCHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( TruthHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( TrackRecordAnalysis ) +DECLARE_ALGORITHM_FACTORY( ALFAHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( LucidHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( ZDCHitAnalysis ) +DECLARE_ALGORITHM_FACTORY( AFPHitAnalysis ) DECLARE_FACTORY_ENTRIES( HitAnalysis ) { DECLARE_ALGORITHM( CaloHitAnalysis ) + DECLARE_ALGORITHM( SiHitAnalysis ) + DECLARE_ALGORITHM( TRTHitAnalysis ) + DECLARE_ALGORITHM( TGCHitAnalysis ) + DECLARE_ALGORITHM( RPCHitAnalysis ) + DECLARE_ALGORITHM( MDTHitAnalysis ) + DECLARE_ALGORITHM( CSCHitAnalysis ) + DECLARE_ALGORITHM( TruthHitAnalysis ) + DECLARE_ALGORITHM( TrackRecordAnalysis ) + DECLARE_ALGORITHM( ALFAHitAnalysis ) + DECLARE_ALGORITHM( LucidHitAnalysis ) + DECLARE_ALGORITHM( ZDCHitAnalysis ) + DECLARE_ALGORITHM( AFPHitAnalysis ) + + }