From 5a38bed43b7893b72b2c51198ad5bf0dd495f54f Mon Sep 17 00:00:00 2001
From: Zhan Li <lizhan@ihep.ac.cn>
Date: Wed, 13 Mar 2024 09:33:50 +0000
Subject: [PATCH] onlyVXD

---
 Examples/options/detsim_tracker.py   |  85 +++++++++++++++++++
 Examples/options/tracker_analysis.py | 120 +++++++++++++++++++++++++++
 2 files changed, 205 insertions(+)
 create mode 100644 Examples/options/detsim_tracker.py
 create mode 100644 Examples/options/tracker_analysis.py

diff --git a/Examples/options/detsim_tracker.py b/Examples/options/detsim_tracker.py
new file mode 100644
index 00000000..fc0ac912
--- /dev/null
+++ b/Examples/options/detsim_tracker.py
@@ -0,0 +1,85 @@
+#!/usr/bin/env python
+#Author: Zhan Li <lizhan@ihep.ac.cn>
+#Created [2024-03-07 Thu 14:53]
+
+import os
+import sys
+
+import Gaudi.Configuration
+from Configurables import RndmGenSvc, HepRndm__Engine_CLHEP__RanluxEngine_, k4DataSvc, GeomSvc
+from Configurables import TimeProjectionChamberSensDetTool
+from Configurables import GenAlgo
+from Configurables import GtGunTool
+from Configurables import StdHepRdr
+from Configurables import SLCIORdr
+from Configurables import HepMCRdr
+from Configurables import GenPrinter
+from Configurables import GtBeamBackgroundTool
+from Configurables import DetSimSvc
+from Configurables import DetSimAlg
+from Configurables import AnExampleDetElemTool
+from Configurables import PodioOutput
+from Configurables import ApplicationMgr
+
+seed = [42]
+
+rndmengine = Gaudi.Configuration.HepRndm__Engine_CLHEP__HepJamesRandom_("RndmGenSvc.Engine")
+rndmengine.SetSingleton = True
+rndmengine.Seeds = seed
+
+rndmgensvc = RndmGenSvc("RndmGenSvc")
+rndmgensvc.Engine = rndmengine.name()
+
+
+dsvc = k4DataSvc("EventDataSvc")
+#geometry_option = "CepC_v4-onlyVXD.xml"
+geometry_option = "CepC_v4_onlyTracker.xml"
+#geometry_option = "CepC_v4.xml"
+
+geometry_path = os.path.join(os.getenv("DETCEPCV4ROOT"), "compact", geometry_option)
+geosvc = GeomSvc("GeomSvc")
+geosvc.compact = geometry_path
+
+#Previously I do not have these 2 lines
+tpc_sensdettool = TimeProjectionChamberSensDetTool("TimeProjectionChamberSensDetTool")
+tpc_sensdettool.TypeOption = 1
+
+
+# Physics Generator
+bg = GtBeamBackgroundTool("GtBeamBackgroundTool")
+bg.InputFileMap = {"default":"/scratchfs/atlas/lizhan/cepc/CEPCSW/ToCEPCSWsingle.out"}
+#bg.InputFileMap = {"default":"/cefs/higgs/shihy/tools/CEPCSW/CEPCSW/Test/0/ToCEPCSW-1.out"}
+bg.InputBeamEnergyMap = {"default":120}
+bg.RotationAlongYMap = {"default":16.5e-3}
+
+
+genprinter = GenPrinter("GenPrinter")
+
+genalg = GenAlgo("GenAlgo")
+genalg.GenTools = ["GtBeamBackgroundTool"]
+
+detsimsvc = DetSimSvc("DetSimSvc")
+
+detsimalg = DetSimAlg("DetSimAlg")
+detsimalg.RandomSeeds = seed
+
+
+detsimalg.RunCmds = []
+detsimalg.AnaElems = [
+    "Edm4hepWriterAnaElemTool"
+]
+detsimalg.RootDetElem = "WorldDetElemTool"
+
+example_dettool = AnExampleDetElemTool("AnExampleDetElemTool")
+
+
+# POD I/O
+out = PodioOutput("outputalg")
+out.filename = "test-SIT.root"
+out.outputCommands = ["keep *"]
+
+ApplicationMgr( TopAlg = [genalg, detsimalg, out],
+                EvtSel = 'NONE',
+                EvtMax = 1,
+                ExtSvc = [rndmengine, rndmgensvc, dsvc, geosvc],
+)
\ No newline at end of file
diff --git a/Examples/options/tracker_analysis.py b/Examples/options/tracker_analysis.py
new file mode 100644
index 00000000..896a31ca
--- /dev/null
+++ b/Examples/options/tracker_analysis.py
@@ -0,0 +1,120 @@
+#!/usr/bin/env python
+#Author: Zhan Li <lizhan@ihep.ac.cn>
+#Created [2024-03-07 Thu 14:53]
+
+import ROOT
+import math
+
+#read
+#path1="/scratchfs/atlas/lizhan/cepc/CEPCSW/test-onlyVXD-single.root"
+#path1="/scratchfs/atlas/lizhan/cepc/CEPCSW/test-detsim10.root"
+
+def GetVXDHitMapFromFile(infile_name,out_file_name):
+    file = ROOT.TFile(infile_name)
+    tree = file.Get("events")
+    posx_VXD=[]
+    posy_VXD=[]
+    posz_VXD=[]
+    posr_VXD=[]
+    for entry in tree:
+        for elem in entry.VXDCollection:
+            x=elem.position.x
+            y=elem.position.y
+            z=elem.position.z
+            r = math.sqrt(x*x + y*y)
+            posx_VXD.append(x)
+            posy_VXD.append(y)
+            posz_VXD.append(z)
+            posr_VXD.append(r)
+    file.Close()
+    VXDHitMap = ROOT.TH2F("VXDHitMap","Hit Map of VXD",40,-125,125,20,0,80)
+    for i in range(0,len(posz_VXD)):
+        VXDHitMap.Fill(posz_VXD[i],posr_VXD[i])
+    ofile=ROOT.TFile(out_file_name,"Recreate")
+    VXDHitMap.Write()
+    ofile.Close()
+    
+    
+def GetSITHitMapFromFile(infile_path,out_file_name):    
+    file = ROOT.TFile(infile_path)
+    tree = file.Get("events")
+    posx_SIT=[]
+    posy_SIT=[]
+    posz_SIT=[]
+    posr_SIT=[]
+    for entry in tree:
+        for elem in entry.SITCollection:
+            x=elem.position.x
+            y=elem.position.y
+            z=elem.position.z
+            r = math.sqrt(x*x + y*y)
+            posx_SIT.append(x)
+            posy_SIT.append(y)
+            posz_SIT.append(z)
+            posr_SIT.append(r)
+    file.Close()
+    SITHitMap = ROOT.TH2F("SITHitMap","Hit Map of SIT",200,-2400,2400,20,120,320)
+    for i in range(0,len(posz_SIT)):
+        SITHitMap.Fill(posz_SIT[i],posr_SIT[i])
+    ofile=ROOT.TFile(out_file_name,"Recreate")
+    SITHitMap.Write()
+    ofile.Close()
+
+print("test")
+
+def GetSETHitMapFromFile(infile_path,out_file_name):    
+    file = ROOT.TFile(infile_path)
+    tree = file.Get("events")
+    posx_SET=[]
+    posy_SET=[]
+    posz_SET=[]
+    posr_SET=[]
+    for entry in tree:
+        for elem in entry.SETCollection:
+            x=elem.position.x
+            y=elem.position.y
+            z=elem.position.z
+            r = math.sqrt(x*x + y*y)
+            posx_SET.append(x)
+            posy_SET.append(y)
+            posz_SET.append(z)
+            posr_SET.append(r)
+    file.Close()
+    SETHitMap = ROOT.TH2F("SETHitMap","Hit Map of SET",200,-2400,2400,20,1700,2000)
+    for i in range(0,len(posz_SET)):
+        SETHitMap.Fill(posz_SET[i],posr_SET[i])
+    ofile=ROOT.TFile(out_file_name,"Recreate")
+    SETHitMap.Write()
+    ofile.Close()
+
+
+
+#draw
+def DrawHitMap(in_file_name,TH2name,out_file_name):
+    file = ROOT.TFile(in_file_name)
+    TH2_HitMap=file.Get(TH2name)
+    c = ROOT.TCanvas("c", "c", 800, 800)
+    ROOT.gStyle.SetOptStat(0)
+    TH2_HitMap.GetXaxis().SetTitle("Z [mm]")
+    TH2_HitMap.GetXaxis().SetTitleOffset(0.77)
+    TH2_HitMap.GetXaxis().SetTitleSize(0.05)
+    TH2_HitMap.GetYaxis().SetTitle("R [mm]")
+    TH2_HitMap.GetYaxis().SetTitleOffset(0.8)
+    TH2_HitMap.GetYaxis().SetTitleSize(0.045)
+    TH2_HitMap.Draw("COL Z CJUST SAME")
+
+    c.SaveAs(out_file_name+".root")
+    c.SaveAs(out_file_name+".pdf")
+
+path1="/scratchfs/atlas/lizhan/cepc/CEPCSW/test-SIT.root"
+
+VXDOutFile="/scratchfs/atlas/lizhan/cepc/CEPCSW/VXDHM.root"
+SITOutFile="/scratchfs/atlas/lizhan/cepc/CEPCSW/SITHM.root"
+SETOutFile="/scratchfs/atlas/lizhan/cepc/CEPCSW/SETHM.root"
+GetVXDHitMapFromFile(path1,VXDOutFile)
+GetSITHitMapFromFile(path1,SITOutFile)
+GetSETHitMapFromFile(path1,SETOutFile)
+
+DrawHitMap(VXDOutFile,"VXDHitMap","Full_VXDHM_1evt")
+DrawHitMap(SITOutFile,"SITHitMap","Full_SITHM_1evt")
+DrawHitMap(SETOutFile,"SETHitMap","Full_SETHM_1evt")
\ No newline at end of file
-- 
GitLab