Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

MIT tutorial first draft ready #144

Merged
merged 1 commit into from
Mar 21, 2024
Merged

MIT tutorial first draft ready #144

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
39 changes: 31 additions & 8 deletions full-detector-simulations/FCCeeGeneralOverview/FCCeeGeneralOverview.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@ So, let's start playing with Full Sim!
ssh -X [email protected]
# set-up the Key4hep environment
source /cvmfs/sw-nightlies.hsf.org/key4hep/setup.sh
# create a repository with the tutorial
# create a repository for the tutorial
mkdir FCC_Full_Sim_Tutorial
cd FCC_Full_Sim_Tutorial
git clone https://github.com/HEP-FCC/fcc-tutorials
Expand Down Expand Up @@ -56,7 +56,7 @@ k4run CLDReconstruction.py --inputFiles ../../fcc-tutorials/full-detector-simula
cd ../../fcc-tutorials/full-detector-simulations/FCCeeGeneralOverview/
```

This created an edm4hep ROOT file with a bunch of new RECO level collections, including `edm4hep::ReconstructedParticle` from Particle Flow (PandoraPFA). You can inspect the ROOT file content with
This created an edm4hep ROOT file with a bunch of new DIGI/RECO level collections, including `edm4hep::ReconstructedParticle` from Particle Flow (PandoraPFA). You can inspect the ROOT file content with

```
podio-dump wzp6_ee_mumuH_ecm240_CLD_RECO_edm4hep.root
Expand All @@ -69,7 +69,7 @@ NB: this step also produces a file named `wzp6_ee_mumuH_ecm240_CLD_RECO_aida.roo

### Plotting the Higgs recoil mass

Let's now use the reconstructed sample to produce some physics quantities. As an example, we will plot the Higgs recoil mass using the Python bindings of edm4hep. The following very simple script already does a decent job:
Let's now use the reconstructed sample to produce some physics quantities. As an example, we will plot the Higgs recoil mass using the Python bindings of [PODIO](https://github.com/key4hep/key4hep-tutorials/blob/main/edm4hep_analysis/edm4hep_api_intro.md#reading-edm4hep-files). The following very simple script already does a decent job:

<!-- Explain a bit the script -->
```
Expand Down Expand Up @@ -169,10 +169,11 @@ Now let's plot the energy resolution for raw clusters and MVA calibrated cluster

```
python plot_calo_energy_resolution.py electron_gun_10GeV_ALLEGRO_RECO.root
# display with your favorite png renderer
display electron_gun_10GeV_ALLEGRO_RECO_clusterEnergyResolution.png

```

Look at both distributions to see how the MVA calibration modifies the cluster energy. NB: we are in the middle of a transition for the Geant4 interface and we are here using the new version which is not yet fully validated, this explains the drop in energy.
Look at both distributions to see how the MVA calibration improves the performance (both in terms of response and in terms of resolution). NB: we are in the middle of a transition for the Geant4 interface and we are here using the new workflow (ddsim) which is not yet fully validated. For instance, the low efficiency is under investigation.

### Changing Liquid Argon to Liquid Krypton

Expand All @@ -187,7 +188,7 @@ And to accomodate this change, we have to update the sampling fraction in the st

```
vim run_ALLEGRO_RECO.py
# comment out line 105 uncomment line 106
# comment out line 105 and uncomment line 106
```
The derivation of new sampling fractions upon geometry change is done by ECAL experts and is out of scope for this tutorial.

Expand All @@ -196,15 +197,38 @@ Let's now re-run the simulation with the modified detector:
ddsim --enableGun --gun.distribution uniform --gun.energy "10*GeV" --gun.particle e- --numberOfEvents 100 --outputFile electron_gun_10GeV_ALLEGRO_LKr_SIM.root --random.enableEventSeed --random.seed 42 --compactFile $K4GEO/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml
k4run run_ALLEGRO_RECO.py --EventDataSvc.input="electron_gun_10GeV_ALLEGRO_LKr_SIM.root" --out.filename="electron_gun_10GeV_ALLEGRO_LKr_RECO.root"
python plot_calo_energy_resolution.py electron_gun_10GeV_ALLEGRO_LKr_RECO.root
display electron_gun_10GeV_ALLEGRO_LKr_RECO_calibratedClusterEnergyResolution.png
```

See how the resolution changed.
You can see that the energy resolution did not improve. This can be due to multiple reasons: the derivation of the sampling fractions was done with little statistics for the LKr option, the MVA calibration was not re-trained, we plot only 100 events, ... But at least now you know how to change some of the detector free parameters :-)


## Towards IDEA tracking with the detailed Drift Chamber

The IDEA tracking system Full Sim description is getting complete (we are just missing the Silicon Wrapper which will arrive soon) and we now have to design tracking algorithms. In this exercise we will produce a dataset containing Vertex and Drift Chamber digitized hits and which can therefore be used to devlop tracking algorithms.

Let's run the IDEA simulation and digitization:

```
ddsim --enableGun --gun.distribution uniform --gun.energy "10*GeV" --gun.particle e- --numberOfEvents 100 --outputFile electron_gun_10GeV_IDEA_SIM.root --random.enableEventSeed --random.seed 42 --compactFile $K4GEO/FCCee/IDEA/compact/IDEA_o1_v02/IDEA_o1_v02.xml
k4run run_IDEA_DIGI.py --EventDataSvc.input="electron_gun_10GeV_IDEA_SIM.root" --out.filename="electron_gun_10GeV_IDEA_DIGI.root"
```

The following collections will be useful for tracking:
```
CDCHDigis extension::DriftChamberDigi
VTXDDigis edm4hep::TrackerHit3D
VTXIBDigis edm4hep::TrackerHit3D
VTXOBDigis edm4hep::TrackerHit3D
```
Note that the drift chamber digitized hit has two positions (on the left or on the right of the wire) because the only accessible information from a drift chamber is the distance to the wire which is degenerated. This left/right ambiguity is alleviated in the tracking algorithm.

You can now play with the digitized hits:
```
root electron_gun_10GeV_IDEA_DIGI.root
events->Draw("leftHitSimHitDeltaDistToWire") // this is the difference of the distance to the wire between the left digi and the sim hit (only produced in debug mode)
events->Draw("CDCHDigis.rightPosition.x:CDCHDigis.rightPosition.y:CDCHDigis.rightPosition.z") // you can see that there was no magnetic field in this simulation and you can see the secondaries created inside the tracking volume
events->Draw("CDCHHits.eDep/CDCHHits.pathLength", "CDCHHits.eDep/CDCHHits.pathLength < 4e-7") // this shows the de/dx from the simHits
```


Expand Down Expand Up @@ -246,4 +270,3 @@ This tool is useful but not perfect and will not meet all the needs (especially
- [FCC Full Sim webpage](https://fcc-ee-detector-full-sim.docs.cern.ch/)
- Bi-weekly FCC Full Sim [working meeting](https://indico.cern.ch/category/16938/)


78 changes: 64 additions & 14 deletions full-detector-simulations/FCCeeGeneralOverview/run_IDEA_DIGI.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,13 +1,15 @@
import os

from Gaudi.Configuration import INFO, DEBUG
from Gaudi.Configuration import *

# input
# Loading the input SIM file
from Configurables import k4DataSvc, PodioInput
evtsvc = k4DataSvc('EventDataSvc')
evtsvc.input = "electron_gun_10GeV_IDEA_SIM.root"
inp = PodioInput('InputReader')

# Detector geometry (needed for the digitizer because we need to go in local coordinates)
################## Simulation setup
# Detector geometry
from Configurables import GeoSvc
geoservice = GeoSvc("GeoSvc")
path_to_detector = os.environ.get("K4GEO", "")
Expand All @@ -19,34 +21,78 @@
geoservice.detectors = [os.path.join(path_to_detector, _det) for _det in detectors_to_use]
geoservice.OutputLevel = INFO

# Digitize drift chamber sim hits
# digitize vertex hits
from Configurables import VTXdigitizer
import math
innerVertexResolution_x = 0.003 # [mm], assume 5 µm resolution for ARCADIA sensor
innerVertexResolution_y = 0.003 # [mm], assume 5 µm resolution for ARCADIA sensor
innerVertexResolution_t = 1000 # [ns]
outerVertexResolution_x = 0.050/math.sqrt(12) # [mm], assume ATLASPix3 sensor with 50 µm pitch
outerVertexResolution_y = 0.150/math.sqrt(12) # [mm], assume ATLASPix3 sensor with 150 µm pitch
outerVertexResolution_t = 1000 # [ns]

vtxib_digitizer = VTXdigitizer("VTXIBdigitizer",
inputSimHits = "VTXIBCollection",
outputDigiHits = "VTXIBDigis",
detectorName = "Vertex",
readoutName = "VTXIBCollection",
xResolution = innerVertexResolution_x, # mm, r-phi direction
yResolution = innerVertexResolution_y, # mm, z direction
tResolution = innerVertexResolution_t,
forceHitsOntoSurface = False,
OutputLevel = INFO
)

vtxob_digitizer = VTXdigitizer("VTXOBdigitizer",
inputSimHits = "VTXOBCollection",
outputDigiHits = "VTXOBDigis",
detectorName = "Vertex",
readoutName = "VTXOBCollection",
xResolution = outerVertexResolution_x, # mm, r-phi direction
yResolution = outerVertexResolution_y, # mm, z direction
tResolution = outerVertexResolution_t, # ns
forceHitsOntoSurface = False,
OutputLevel = INFO
)

vtxd_digitizer = VTXdigitizer("VTXDdigitizer",
inputSimHits = "VTXDCollection",
outputDigiHits = "VTXDDigis",
detectorName = "Vertex",
readoutName = "VTXDCollection",
xResolution = outerVertexResolution_x, # mm, r direction
yResolution = outerVertexResolution_y, # mm, phi direction
tResolution = outerVertexResolution_t, # ns
forceHitsOntoSurface = False,
OutputLevel = INFO
)

# digitize drift chamber hits
from Configurables import DCHsimpleDigitizerExtendedEdm
dch_digitizer = DCHsimpleDigitizerExtendedEdm("DCHsimpleDigitizerExtendedEdm",
inputSimHits = "",
outputDigiHits = savetrackertool.SimTrackHits.Path.replace("sim", "digi"),
inputSimHits = "CDCHHits",
outputDigiHits = "CDCHDigis",
outputSimDigiAssociation = "DC_simDigiAssociation",
readoutName = "CDCHHits",
xyResolution = 0.1, # mm
zResolution = 1, # mm
debugMode = False,
debugMode = True,
OutputLevel = INFO
)


# Output
################ Output
from Configurables import PodioOutput
out = PodioOutput("out",
OutputLevel=INFO)
out.outputCommands = ["keep *"]
out.filename = "output_simplifiedDriftChamber_MagneticField_"+str(magneticField)+"_pMin_"+str(momentum*1000)+"_MeV"+"_ThetaMinMax_"+str(thetaMin)+"_"+str(thetaMax)+"_pdgId_"+str(pdgCode)+"_stepLength_default.root"

out.filename = "electron_gun_10GeV_IDEA_DIGI.root"

#CPU information
from Configurables import AuditorSvc, ChronoAuditor
chra = ChronoAuditor()
audsvc = AuditorSvc()
audsvc.Auditors = [chra]
genAlg.AuditExecute = True
hepmc_converter.AuditExecute = True
out.AuditExecute = True

from Configurables import EventCounter
Expand All @@ -56,12 +102,16 @@
from Configurables import ApplicationMgr
ApplicationMgr(
TopAlg = [
inp,
event_counter,
vtxib_digitizer,
vtxob_digitizer,
vtxd_digitizer,
dch_digitizer,
out
],
EvtSel = 'NONE',
EvtMax = 100,
ExtSvc = [geoservice, podioevent, , audsvc],
EvtMax = -1,
ExtSvc = [geoservice, evtsvc, audsvc],
StopOnSignal = True,
)
Loading