From a07e90a928857746546ffa9524995b242bee30a4 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 26 Oct 2023 11:37:44 +0200 Subject: [PATCH 01/15] Migration of new ALLEGRO barrel ECAL from FCCDetectors, https://github.com/HEP-FCC/FCCDetectors/pull/56 --- .../compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml | 46 + .../ALLEGRO_o1_v02/BeamInstrumentation.xml | 36 + .../compact/ALLEGRO_o1_v02/Beampipe.xml | 150 +++ .../ALLEGRO_o1_v02/FCCee_DectDimensions.xml | 230 +++++ .../ALLEGRO_o1_v02/FCCee_DectEmptyMaster.xml | 30 + .../ALLEGRO_o1_v02/FCCee_ECalBarrel.xml | 143 +++ .../FCCee_ECalBarrel_thetamodulemerged.xml | 144 +++ .../FCCee_EcalEndcaps_coneCryo.xml | 167 ++++ .../FCCee_HCalBarrel_TileCal.xml | 133 +++ .../FCCee_HCalEndcaps_ThreeParts_TileCal.xml | 141 +++ .../compact/ALLEGRO_o1_v02/HOMAbsorber.xml | 64 ++ .../compact/ALLEGRO_o1_v02/LumiCal.xml | 191 ++++ .../compact/ALLEGRO_o1_v02/MuonTagger.xml | 50 + .../ALLEGRO_o1_v02/SimplifiedDriftChamber.xml | 58 ++ .../ALLEGRO_o1_v02/Solenoid_o1_v01_02.xml | 57 ++ .../ALLEGRO/compact/ALLEGRO_o1_v02/Vertex.xml | 234 +++++ .../compact/ALLEGRO_o1_v02/elements.xml | 884 ++++++++++++++++++ .../compact/ALLEGRO_o1_v02/materials.xml | 261 ++++++ FCCee/ALLEGRO/compact/README.md | 4 +- .../calorimeter/ECalBarrelInclined_geo.cpp | 614 ++++++++++++ .../FCCSWGridModuleThetaMerged.h | 117 +++ .../FCCSWGridModuleThetaMergedHandle.h | 124 +++ .../include/detectorSegmentations/GridTheta.h | 102 ++ .../src/FCCSWGridModuleThetaMerged.cpp | 156 ++++ .../src/FCCSWGridModuleThetaMergedHandle.cpp | 4 + detectorSegmentations/src/GridEta.cpp | 57 ++ .../src/plugins/SegmentationFactories.cpp | 27 +- 27 files changed, 4211 insertions(+), 13 deletions(-) create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/BeamInstrumentation.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Beampipe.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectDimensions.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectEmptyMaster.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_EcalEndcaps_coneCryo.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalBarrel_TileCal.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalEndcaps_ThreeParts_TileCal.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/HOMAbsorber.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/LumiCal.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/MuonTagger.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/SimplifiedDriftChamber.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Solenoid_o1_v01_02.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Vertex.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/elements.xml create mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/materials.xml create mode 100644 detector/calorimeter/ECalBarrelInclined_geo.cpp create mode 100644 detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h create mode 100644 detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h create mode 100644 detectorSegmentations/include/detectorSegmentations/GridTheta.h create mode 100644 detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp create mode 100644 detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp create mode 100644 detectorSegmentations/src/GridEta.cpp diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml new file mode 100644 index 00000000..49c76c53 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/ALLEGRO_o1_v02.xml @@ -0,0 +1,46 @@ + + + + + + Master compact file describing the latest developments of the FCCee IDEA detector concept with a LAr calorimeter. + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/BeamInstrumentation.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/BeamInstrumentation.xml new file mode 100644 index 00000000..49647b2b --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/BeamInstrumentation.xml @@ -0,0 +1,36 @@ + + + + COmpensating and screening solenoids for FCCee + + + + + Beampipe Instrumentation + + + + + + +
+ + + + + + + + +
+ +
+ + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Beampipe.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Beampipe.xml new file mode 100644 index 00000000..b86af521 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Beampipe.xml @@ -0,0 +1,150 @@ + + + + A beampipe for FCCee, R(central) = 1.5 cm + + + + + + + + + + + + + + + + + + + + Part of beampipe made of Beryllium + + + + + + + +
+ +
+ + + + + + + + + + + + + Golden foil in the inner part of the Be beampipe + +
+ +
+ + Part of beampipe made of Copper + +
+ + + + + +
+ + + +
+ + +
+ + +Full Cone Tungsten Shield + + + + Before HOM space +
+ + After HOM space (1197.5*m - 1298.7*mm) +18 cm as solenoid is now closer to IP +
+ + +Asymmetric Tungsten Shield no Rotation + + + + +
+ + was 370 +
+ + one degree less, to fit lumical window +
+ +
+ + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectDimensions.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectDimensions.xml new file mode 100644 index 00000000..439e20cd --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectDimensions.xml @@ -0,0 +1,230 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectEmptyMaster.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectEmptyMaster.xml new file mode 100644 index 00000000..a5ad3e23 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_DectEmptyMaster.xml @@ -0,0 +1,30 @@ + + + + + + + + + + + + Use this one if you want to use official dimensions but only place one detector inside + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml new file mode 100644 index 00000000..b1e0e9f2 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml @@ -0,0 +1,143 @@ + + + + + + Settings for the inclined EM calorimeter. + The barrel is filled with liquid argon. Passive material includes lead in the middle and steal on the outside, glued together. + Passive plates are inclined by a certain angle from the radial direction. + In between of two passive plates there is a readout. + Space between the plate and readout is of trapezoidal shape and filled with liquid argon. + Definition of sizes, visualization settings, readout and longitudinal segmentation are specified. + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + system:4,cryo:1,type:3,subtype:3,layer:8,module:11,eta:9 + + + + + + + system:4,cryo:1,type:3,subtype:3,layer:8,eta:9,phi:10 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml new file mode 100644 index 00000000..9411b7dc --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml @@ -0,0 +1,144 @@ + + + + + + Settings for the inclined EM calorimeter. + The barrel is filled with liquid argon. Passive material includes lead in the middle and steal on the outside, glued together. + Passive plates are inclined by a certain angle from the radial direction. + In between of two passive plates there is a readout. + Space between the plate and readout is of trapezoidal shape and filled with liquid argon. + Definition of sizes, visualization settings, readout and longitudinal segmentation are specified. + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + system:4,cryo:1,type:3,subtype:3,layer:8,module:11,theta:10 + + + + + + system:4,cryo:1,type:3,subtype:3,layer:8,module:11,theta:10 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_EcalEndcaps_coneCryo.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_EcalEndcaps_coneCryo.xml new file mode 100644 index 00000000..33abcf03 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_EcalEndcaps_coneCryo.xml @@ -0,0 +1,167 @@ + + + + + + Liquid argon EM calorimeter endcap design. + Electromagnetic part (EMEC) includes lead+steel absorber. + Absorbers have shape of simple discs. + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + system:4,subsystem:1,type:3,subtype:3,layer:8,sublayer:8,eta:10,phi:10 + + + + system:4,subsystem:1,type:3,subtype:3,layer:8,eta:10,phi:10 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalBarrel_TileCal.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalBarrel_TileCal.xml new file mode 100644 index 00000000..0acecdb1 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalBarrel_TileCal.xml @@ -0,0 +1,133 @@ + + + + The first FCCee HCal layout based on ATLAS HCal, not optimised yet + 1. Nov 2022, J. Faltova: update material and radial segmentation for FCCee + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + system:4,layer:5,row:9,eta:9,phi:10 + + + + system:4,layer:5,eta:9,phi:10 + + + + + system:4,layer:5,row:9,eta:0,phi:10 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalEndcaps_ThreeParts_TileCal.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalEndcaps_ThreeParts_TileCal.xml new file mode 100644 index 00000000..c5cbe72a --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_HCalEndcaps_ThreeParts_TileCal.xml @@ -0,0 +1,141 @@ + + + + HCal layout based on ATLAS HCal, with realistic longitudinal segmentation and steel support + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + system:4,type:2,layer:5,row:9,eta:11,phi:10 + + + + system:4,type:2,layer:4,eta:11,phi:10 + + + + + system:4,type:2,layer:5,row:9,eta:10,phi:10 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/HOMAbsorber.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/HOMAbsorber.xml new file mode 100644 index 00000000..862207ea --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/HOMAbsorber.xml @@ -0,0 +1,64 @@ + + + + Higher mode absorber for FCCee + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
+ + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/LumiCal.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/LumiCal.xml new file mode 100644 index 00000000..15f0434e --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/LumiCal.xml @@ -0,0 +1,191 @@ + + + + LumiCal for FCCee detector based on CLD + + + + + + + + + + system:8,barrel:3,layer:8,slice:8,r:32:-16,phi:-16 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/MuonTagger.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/MuonTagger.xml new file mode 100644 index 00000000..79b0ed18 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/MuonTagger.xml @@ -0,0 +1,50 @@ + + + + + + Simple muon tagger - barrel and endcaps + + + + + + + + + + + system:4,subsystem:1,type:3,subtype:3,layer:8,sublayer:8,eta:10,phi:10 + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/SimplifiedDriftChamber.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/SimplifiedDriftChamber.xml new file mode 100644 index 00000000..11586a5d --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/SimplifiedDriftChamber.xml @@ -0,0 +1,58 @@ + + + + + + A simplified implementation of the drift chamber for the FCCee-IDEA concept + + + + + + + + + + + + + + + + + + + + + + + + + + + + + ${GlobalTrackerReadoutID_DCH} + + + + + + + + Dimensions for the drift chamber + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Solenoid_o1_v01_02.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Solenoid_o1_v01_02.xml new file mode 100644 index 00000000..69a117f8 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Solenoid_o1_v01_02.xml @@ -0,0 +1,57 @@ + + + + + + + + + + + + + + + + + + + Solenoid + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Vertex.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Vertex.xml new file mode 100644 index 00000000..a28449d6 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/Vertex.xml @@ -0,0 +1,234 @@ + + + + + + CLD Vertex Detector for FCCee + + + Tracking detectors + + + + + + + + + + + + Vertex Assembly + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + ${GlobalTrackerReadoutID} + + + ${GlobalTrackerReadoutID} + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + Vertex Detector Endcaps + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/elements.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/elements.xml new file mode 100644 index 00000000..f35eb345 --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/elements.xml @@ -0,0 +1,884 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/materials.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/materials.xml new file mode 100644 index 00000000..a899266c --- /dev/null +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/materials.xml @@ -0,0 +1,261 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +

+ + + + + + + + + + + + + + + + + + + + + + + + diff --git a/FCCee/ALLEGRO/compact/README.md b/FCCee/ALLEGRO/compact/README.md index ba8b5c28..9a6ded62 100644 --- a/FCCee/ALLEGRO/compact/README.md +++ b/FCCee/ALLEGRO/compact/README.md @@ -1,3 +1,5 @@ ALLEGRO ======================== -ALLEGRO_o1_v01: it is a liquid Noble gas based detector. This version picked from the latest version in FCCDetectors repo. +ALLEGRO_o1_v01: it is a liquid Noble gas based detector. This version picked from the latest version in FCCDetectors repo. + +ALLEGRO_o1_v02: evolves from o1_v01, replacing the barrel ECAL. This version has a constant cell size in theta for the ECAL barrel (instead of eta as in o1_v01) and now it is possible to have a different number of cells merged for each longitudinal layer. diff --git a/detector/calorimeter/ECalBarrelInclined_geo.cpp b/detector/calorimeter/ECalBarrelInclined_geo.cpp new file mode 100644 index 00000000..51178159 --- /dev/null +++ b/detector/calorimeter/ECalBarrelInclined_geo.cpp @@ -0,0 +1,614 @@ +#include "DD4hep/DetFactoryHelper.h" +#include "DD4hep/Handle.h" +#include "XML/Utilities.h" + +#include + +// todo: remove gaudi logging and properly capture output +#define endmsg std::endl +#define lLog std::cout +namespace MSG { +const std::string ERROR = " Error: "; +const std::string DEBUG = " Debug: "; +const std::string INFO = " Info: "; +} + +namespace det { +static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, + dd4hep::xml::Handle_t aXmlElement, + dd4hep::SensitiveDetector aSensDet) { + + dd4hep::xml::DetElement xmlDetElem = aXmlElement; + std::string nameDet = xmlDetElem.nameStr(); + dd4hep::xml::Dimension dim(xmlDetElem.dimensions()); + dd4hep::DetElement caloDetElem(nameDet, xmlDetElem.id()); + + // Create air envelope for the whole barrel + dd4hep::Volume envelopeVol(nameDet + "_vol", dd4hep::Tube(dim.rmin(), dim.rmax(), dim.dz()), + aLcdd.material("Air")); + envelopeVol.setVisAttributes(aLcdd, dim.visStr()); + + // Retrieve cryostat data + dd4hep::xml::DetElement cryostat = aXmlElement.child(_Unicode(cryostat)); + dd4hep::xml::Dimension cryoDim(cryostat.dimensions()); + double cryoThicknessFront = cryoDim.rmin2() - cryoDim.rmin1(); + dd4hep::xml::DetElement cryoFront = cryostat.child(_Unicode(front)); + dd4hep::xml::DetElement cryoBack = cryostat.child(_Unicode(back)); + dd4hep::xml::DetElement cryoSide = cryostat.child(_Unicode(side)); + bool cryoFrontSensitive = cryoFront.isSensitive(); + bool cryoBackSensitive = cryoBack.isSensitive(); + bool cryoSideSensitive = cryoSide.isSensitive(); + + // Retrieve active and passive material data + dd4hep::xml::DetElement calo = aXmlElement.child(_Unicode(calorimeter)); + dd4hep::xml::Dimension caloDim(calo.dimensions()); + dd4hep::xml::DetElement active = calo.child(_Unicode(active)); + std::string activeMaterial = active.materialStr(); + double activeThickness = active.thickness(); + + dd4hep::xml::DetElement overlap = active.child(_Unicode(overlap)); + double activePassiveOverlap = overlap.offset(); + if (activePassiveOverlap < 0 || activePassiveOverlap > 0.5) { + // todo: ServiceHandle incidentSvc("IncidentSvc", "ECalConstruction"); + lLog << MSG::ERROR << "Overlap between active and passive cannot be more than half of passive plane!" << endmsg; + //todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + } + dd4hep::xml::DetElement layers = calo.child(_Unicode(layers)); + uint numLayers = 0; + std::vector layerHeight; + double layersTotalHeight = 0; + for (dd4hep::xml::Collection_t layer_coll(layers, _Unicode(layer)); layer_coll; ++layer_coll) { + dd4hep::xml::Component layer = layer_coll; + numLayers += layer.repeat(); + for (int iLay = 0; iLay < layer.repeat(); iLay++) { + layerHeight.push_back(layer.thickness()); + } + layersTotalHeight += layer.repeat() * layer.thickness(); + } + lLog << MSG::DEBUG << "Number of layers: " << numLayers << " total thickness " << layersTotalHeight << endmsg; + // The following code checks if the xml geometry file contains a constant defining + // the number of layers the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of layers is needed + // in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). + int nLayers = -1; + try { + nLayers = aLcdd.constant("ECalBarrelNumLayers"); + } + catch(...) { + ; + } + if (nLayers > 0 && nLayers != numLayers) { + lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + assert(nLayers == numLayers); + } + + dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); + std::string readoutMaterial = readout.materialStr(); + double readoutThickness = readout.thickness(); + + dd4hep::xml::DetElement passive = calo.child(_Unicode(passive)); + dd4hep::xml::DetElement passiveInner = passive.child(_Unicode(inner)); + dd4hep::xml::DetElement passiveInnerMax = passive.child(_Unicode(innerMax)); + dd4hep::xml::DetElement passiveOuter = passive.child(_Unicode(outer)); + dd4hep::xml::DetElement passiveGlue = passive.child(_Unicode(glue)); + std::string passiveInnerMaterial = passiveInner.materialStr(); + std::string passiveOuterMaterial = passiveOuter.materialStr(); + std::string passiveGlueMaterial = passiveGlue.materialStr(); + double passiveInnerThicknessMin = passiveInner.thickness(); + double passiveInnerThicknessMax = passiveInnerMax.thickness(); + double passiveOuterThickness = passiveOuter.thickness(); + double passiveGlueThickness = passiveGlue.thickness(); + double passiveThickness = passiveInnerThicknessMin + passiveOuterThickness + passiveGlueThickness; + double angle = passive.rotation().angle(); + + double bathRmin = caloDim.rmin(); // - margin for inclination + double bathRmax = caloDim.rmax(); // + margin for inclination + dd4hep::Tube bathOuterShape(bathRmin, bathRmax, caloDim.dz()); // make it 4 volumes + 5th for detector envelope + dd4hep::Tube bathAndServicesOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), caloDim.dz()); // make it 4 volumes + 5th for detector envelope + if (cryoThicknessFront > 0) { + // 1. Create cryostat + dd4hep::Tube cryoFrontShape(cryoDim.rmin1(), cryoDim.rmin2(), cryoDim.dz()); + dd4hep::Tube cryoBackShape(cryoDim.rmax1(), cryoDim.rmax2(), cryoDim.dz()); + dd4hep::Tube cryoSideOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), cryoDim.dz()); + dd4hep::SubtractionSolid cryoSideShape(cryoSideOuterShape, bathAndServicesOuterShape); + lLog << MSG::INFO << "ECAL cryostat: front: rmin (cm) = " << cryoDim.rmin1() << " rmax (cm) = " << cryoDim.rmin2() << " dz (cm) = " << cryoDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL cryostat: back: rmin (cm) = " << cryoDim.rmax1() << " rmax (cm) = " << cryoDim.rmax2() << " dz (cm) = " << cryoDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL cryostat: side: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << cryoDim.dz() - caloDim.dz() << endmsg; + dd4hep::Volume cryoFrontVol(cryostat.nameStr()+"_front", cryoFrontShape, aLcdd.material(cryostat.materialStr())); + dd4hep::Volume cryoBackVol(cryostat.nameStr()+"_back", cryoBackShape, aLcdd.material(cryostat.materialStr())); + dd4hep::Volume cryoSideVol(cryostat.nameStr()+"_side", cryoSideShape, aLcdd.material(cryostat.materialStr())); + dd4hep::PlacedVolume cryoFrontPhysVol = envelopeVol.placeVolume(cryoFrontVol); + dd4hep::PlacedVolume cryoBackPhysVol = envelopeVol.placeVolume(cryoBackVol); + dd4hep::PlacedVolume cryoSidePhysVol = envelopeVol.placeVolume(cryoSideVol); + if (cryoFrontSensitive) { + cryoFrontVol.setSensitiveDetector(aSensDet); + cryoFrontPhysVol.addPhysVolID("cryo", 1); + cryoFrontPhysVol.addPhysVolID("type", 1); + lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; + } + if (cryoBackSensitive) { + cryoBackVol.setSensitiveDetector(aSensDet); + cryoBackPhysVol.addPhysVolID("cryo", 1); + cryoBackPhysVol.addPhysVolID("type", 2); + lLog << MSG::INFO << "Cryostat back volume set as sensitive" << endmsg; + } + if (cryoSideSensitive) { + cryoSideVol.setSensitiveDetector(aSensDet); + cryoSidePhysVol.addPhysVolID("cryo", 1); + cryoSidePhysVol.addPhysVolID("type", 3); + lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; + } + dd4hep::DetElement cryoFrontDetElem(caloDetElem, "cryo_front", 0); + cryoFrontDetElem.setPlacement(cryoFrontPhysVol); + dd4hep::DetElement cryoBackDetElem(caloDetElem, "cryo_back", 0); + cryoBackDetElem.setPlacement(cryoBackPhysVol); + dd4hep::DetElement cryoSideDetElem(caloDetElem, "cryo_side", 0); + cryoSideDetElem.setPlacement(cryoSidePhysVol); + // 1.2. Create place-holder for services + dd4hep::Tube servicesFrontShape(cryoDim.rmin2(), bathRmin, caloDim.dz()); + dd4hep::Tube servicesBackShape(bathRmax, cryoDim.rmax1(), caloDim.dz()); + lLog << MSG::INFO << "ECAL services: front: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << bathRmin << " dz (cm) = " << caloDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL services: back: rmin (cm) = " << bathRmax << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << caloDim.dz() << endmsg; + dd4hep::Volume servicesFrontVol("services_front", servicesFrontShape, aLcdd.material(activeMaterial)); + dd4hep::Volume servicesBackVol("services_back", servicesBackShape, aLcdd.material(activeMaterial)); + dd4hep::PlacedVolume servicesFrontPhysVol = envelopeVol.placeVolume(servicesFrontVol); + dd4hep::PlacedVolume servicesBackPhysVol = envelopeVol.placeVolume(servicesBackVol); + if (cryoFrontSensitive) { + servicesFrontVol.setSensitiveDetector(aSensDet); + servicesFrontPhysVol.addPhysVolID("cryo", 1); + servicesFrontPhysVol.addPhysVolID("type", 4); + lLog << MSG::INFO << "Services front volume set as sensitive" << endmsg; + } + if (cryoBackSensitive) { + servicesBackVol.setSensitiveDetector(aSensDet); + servicesBackPhysVol.addPhysVolID("cryo", 1); + servicesBackPhysVol.addPhysVolID("type", 5); + lLog << MSG::INFO << "Services back volume set as sensitive" << endmsg; + } + dd4hep::DetElement servicesFrontDetElem(caloDetElem, "services_front", 0); + servicesFrontDetElem.setPlacement(servicesFrontPhysVol); + dd4hep::DetElement servicesBackDetElem(caloDetElem, "services_back", 0); + servicesBackDetElem.setPlacement(servicesBackPhysVol); + } + // 2. Create bath that is inside the cryostat and surrounds the detector + // Bath is filled with active material -> but not sensitive + dd4hep::Volume bathVol(activeMaterial + "_bath", bathOuterShape, aLcdd.material(activeMaterial)); + lLog << MSG::INFO << "ECAL bath: material = " << activeMaterial << " rmin (cm) = " << bathRmin + << " rmax (cm) = " << bathRmax << " thickness in front of ECal (cm) = " << caloDim.rmin() - cryoDim.rmin2() + << " thickness behind ECal (cm) = " << cryoDim.rmax1() - caloDim.rmax() << endmsg; + + // 3. Create the calorimeter by placing the passive material, trapezoid active layers, readout and again trapezoid + // active layers in the bath. + // sensitive detector for the layers + dd4hep::SensitiveDetector sd = aSensDet; + dd4hep::xml::Dimension sdType = xmlDetElem.child(_U(sensitive)); + sd.setType(sdType.typeStr()); + lLog << MSG::INFO << "ECAL calorimeter volume rmin (cm) = " << caloDim.rmin() << " rmax (cm) = " << caloDim.rmax() + << endmsg; + + // 3.a. Create the passive planes, readout in between of 2 passive planes and the remaining space fill with active + // material + ////////////////////////////// + // PASSIVE PLANES + ////////////////////////////// + lLog << MSG::INFO << "passive inner material = " << passiveInnerMaterial << "\n" + << " and outer material = " << passiveOuterMaterial << "\n" + << " thickness of inner part at inner radius (cm) = " << passiveInnerThicknessMin << "\n" + << " thickness of inner part at outer radius (cm) = " << passiveInnerThicknessMax << "\n" + << " thickness of outer part (cm) = " << passiveOuterThickness << "\n" + << " thickness of total (cm) = " << passiveThickness << "\n" + << " rotation angle = " << angle << endmsg; + uint numPlanes = + round(M_PI / asin((passiveThickness + activeThickness + readoutThickness) / (2. * caloDim.rmin() * cos(angle)))); + + double dPhi = 2. * M_PI / numPlanes; + lLog << MSG::INFO << "number of passive plates = " << numPlanes << " azim. angle difference = " << dPhi << endmsg; + lLog << MSG::INFO << " distance at inner radius (cm) = " << 2. * M_PI * caloDim.rmin() / numPlanes << "\n" + << " distance at outer radius (cm) = " << 2. * M_PI * caloDim.rmax() / numPlanes << endmsg; + + // The following code checks if the xml geometry file contains a constant defining + // the number of planes in the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of plane information (retrieved from the + // xml) is used in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). In principle the code above should be refactored so that the number + // of planes is one of the inputs of the calculation and other geometrical parameters + // are adjusted accordingly. This is left for the future, and we use the workaround + // below to enforce for the time being that the number of planes is "correct" + int nModules = -1; + try { + nModules = aLcdd.constant("ECalBarrelNumPlanes"); + } + catch(...) { + ; + } + if (nModules > 0 && nModules != numPlanes) { + lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + assert(nModules == numPlanes); + } + // Readout is in the middle between two passive planes + double offsetPassivePhi = caloDim.offset() + dPhi / 2.; + double offsetReadoutPhi = caloDim.offset() + 0; + lLog << MSG::INFO << "readout material = " << readoutMaterial << "\n" + << " thickness of readout planes (cm) = " << readoutThickness << "\n number of readout layers = " << numLayers + << endmsg; + double Rmin = caloDim.rmin(); + double Rmax = caloDim.rmax(); + double dR = Rmax - Rmin; + double planeLength = -Rmin * cos(angle) + sqrt(pow(Rmax, 2) - pow(Rmin * sin(angle), 2)); + lLog << MSG::INFO << "thickness of calorimeter (cm) = " << dR << "\n" + << " length of passive or readout planes (cm) = " << planeLength << endmsg; + + + // fill the thickness in the boundary of each layer + std::vector passiveInnerThicknessLayer(numLayers+1); + double runningHeight = 0; + for (uint iLay = 0; iLay < numLayers; iLay++) { + passiveInnerThicknessLayer[iLay] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * + (runningHeight) / (Rmax - Rmin); + runningHeight += layerHeight[iLay]; + } + passiveInnerThicknessLayer[numLayers] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * + (runningHeight) / (Rmax - Rmin); + + double passiveAngle = atan2((passiveInnerThicknessMax - passiveInnerThicknessMin) / 2., planeLength); + double cosPassiveAngle = cos(passiveAngle); + double rotatedOuterThickness = passiveOuterThickness / cosPassiveAngle; + double rotatedGlueThickness = passiveGlueThickness / cosPassiveAngle; + + // rescale layer thicknesses + double scaleLayerThickness = planeLength / layersTotalHeight; + layersTotalHeight = 0; + for (uint iLay = 0; iLay < numLayers; iLay++) { + layerHeight[iLay] *= scaleLayerThickness; + + layersTotalHeight += layerHeight[iLay]; + lLog << MSG::DEBUG << "Thickness of layer " << iLay << " : " << layerHeight[iLay] << endmsg; + } + double layerFirstOffset = -planeLength / 2. + layerHeight[0] / 2.; + + //dd4hep::Box passiveShape(passiveThickness / 2., caloDim.dz(), planeLength / 2.); + dd4hep::Trd1 passiveShape(passiveInnerThicknessMin / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., + passiveInnerThicknessMax / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., + caloDim.dz(), planeLength / 2.); + // inner layer is not in the first calo layer (to sample more uniformly in the layer where upstream correction is + // applied) + //dd4hep::Box passiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); + dd4hep::Trd1 passiveInnerShape(passiveInnerThicknessLayer[1] / 2., passiveInnerThicknessMax / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); + //dd4hep::Box passiveInnerShapeFirstLayer(passiveInnerThickness / 2., caloDim.dz(), layerHeight[0] / 2.); + dd4hep::Trd1 passiveInnerShapeFirstLayer(passiveInnerThicknessMin / 2., passiveInnerThicknessLayer[1] / 2., caloDim.dz(), layerHeight[0] / 2.); + dd4hep::Box passiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); + dd4hep::Box passiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); + // passive volume consists of inner part and two outer, joind by glue + dd4hep::Volume passiveVol("passive", passiveShape, aLcdd.material("Air")); + dd4hep::Volume passiveInnerVol(passiveInnerMaterial + "_passive", passiveInnerShape, + aLcdd.material(passiveInnerMaterial)); + dd4hep::Volume passiveInnerVolFirstLayer(activeMaterial + "_passive", passiveInnerShapeFirstLayer, + aLcdd.material(activeMaterial)); + dd4hep::Volume passiveOuterVol(passiveOuterMaterial + "_passive", passiveOuterShape, + aLcdd.material(passiveOuterMaterial)); + dd4hep::Volume passiveGlueVol(passiveGlueMaterial + "_passive", passiveGlueShape, + aLcdd.material(passiveGlueMaterial)); + + if (passiveInner.isSensitive()) { + lLog << MSG::DEBUG << "Passive inner volume set as sensitive" << endmsg; + // inner part starts at second layer + double layerOffset = layerFirstOffset + layerHeight[1] / 2.; + for (uint iLayer = 1; iLayer < numLayers; iLayer++) { + //dd4hep::Box layerPassiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Trd1 layerPassiveInnerShape(passiveInnerThicknessLayer[iLayer] / 2., passiveInnerThicknessLayer[iLayer+1] / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Volume layerPassiveInnerVol(passiveInnerMaterial, layerPassiveInnerShape, + aLcdd.material(passiveInnerMaterial)); + layerPassiveInnerVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveInnerPhysVol = + passiveInnerVol.placeVolume(layerPassiveInnerVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveInnerPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveInnerDetElem("layer", iLayer); + layerPassiveInnerDetElem.setPlacement(layerPassiveInnerPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + } + if (passiveOuter.isSensitive()) { + lLog << MSG::DEBUG << "Passive outer volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset / cosPassiveAngle; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerPassiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); + dd4hep::Volume layerPassiveOuterVol(passiveOuterMaterial, layerPassiveOuterShape, + aLcdd.material(passiveOuterMaterial)); + layerPassiveOuterVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveOuterPhysVol = + passiveOuterVol.placeVolume(layerPassiveOuterVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveOuterPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveOuterDetElem("layer", iLayer); + layerPassiveOuterDetElem.setPlacement(layerPassiveOuterPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; + } + } + } + if (passiveGlue.isSensitive()) { + lLog << MSG::DEBUG << "Passive glue volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset / cosPassiveAngle; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerPassiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); + dd4hep::Volume layerPassiveGlueVol(passiveGlueMaterial, layerPassiveGlueShape, + aLcdd.material(passiveGlueMaterial)); + layerPassiveGlueVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveGluePhysVol = + passiveGlueVol.placeVolume(layerPassiveGlueVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveGluePhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveGlueDetElem("layer", iLayer); + layerPassiveGlueDetElem.setPlacement(layerPassiveGluePhysVol); + if (iLayer != numLayers - 1) { + layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; + } + } + } + + dd4hep::PlacedVolume passiveInnerPhysVol = + passiveVol.placeVolume(passiveInnerVol, dd4hep::Position(0, 0, layerHeight[0] / 2.)); + dd4hep::PlacedVolume passiveInnerPhysVolFirstLayer = + passiveVol.placeVolume(passiveInnerVolFirstLayer, dd4hep::Position(0, 0, layerFirstOffset)); + dd4hep::PlacedVolume passiveOuterPhysVolBelow = passiveVol.placeVolume( + passiveOuterVol, + dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), + dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - + rotatedGlueThickness / 2. - rotatedOuterThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveOuterPhysVolAbove = passiveVol.placeVolume( + passiveOuterVol, + dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), + dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + + rotatedGlueThickness / 2. + rotatedOuterThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveGluePhysVolBelow = passiveVol.placeVolume( + passiveGlueVol, + dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), + dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - + rotatedGlueThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveGluePhysVolAbove = passiveVol.placeVolume( + passiveGlueVol, + dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), + dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + + rotatedGlueThickness / 4., 0, 0))); + passiveInnerPhysVol.addPhysVolID("subtype", 0); + passiveInnerPhysVolFirstLayer.addPhysVolID("subtype", 0); + passiveOuterPhysVolBelow.addPhysVolID("subtype", 1); + passiveOuterPhysVolAbove.addPhysVolID("subtype", 2); + passiveGluePhysVolBelow.addPhysVolID("subtype", 3); + passiveGluePhysVolAbove.addPhysVolID("subtype", 4); + if (passiveInner.isSensitive()) { + passiveInnerVolFirstLayer.setSensitiveDetector(aSensDet); + passiveInnerPhysVolFirstLayer.addPhysVolID("layer", 0); + dd4hep::DetElement passiveInnerDetElemFirstLayer("layer", 0); + passiveInnerDetElemFirstLayer.setPlacement(passiveInnerPhysVolFirstLayer); + } + + ////////////////////////////// + // READOUT PLANES + ////////////////////////////// + dd4hep::Box readoutShape(readoutThickness / 2., caloDim.dz(), planeLength / 2.); + dd4hep::Volume readoutVol(readoutMaterial, readoutShape, aLcdd.material(readoutMaterial)); + if (readout.isSensitive()) { + lLog << MSG::INFO << "Readout volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerReadoutShape(readoutThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Volume layerReadoutVol(readoutMaterial, layerReadoutShape, aLcdd.material(readoutMaterial)); + layerReadoutVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerReadoutPhysVol = + readoutVol.placeVolume(layerReadoutVol, dd4hep::Position(0, 0, layerOffset)); + layerReadoutPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerReadoutDetElem("layer", iLayer); + layerReadoutDetElem.setPlacement(layerReadoutPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + } + + ////////////////////////////// + // ACTIVE + ////////////////////////////// + // thickness of active layers at inner radius and outer ( = distance between passive plane and readout plane) + // at inner radius: distance projected at plane perpendicular to readout plane + double activeInThickness = Rmin * sin(dPhi / 2.) * cos(angle); + activeInThickness -= passiveThickness * (0.5 - activePassiveOverlap); + // at outer radius: distance projected at plane perpendicular to readout plane + double activeOutThickness = (Rmin + planeLength) * sin(dPhi / 2.) * cos(angle); + // make correction for outer readius caused by inclination angle + // first calculate intersection of readout plane and plane parallel to shifted passive plane + double xIntersect = (Rmin * (tan(angle) - cos(dPhi / 2.) * tan(angle + dPhi / 2.)) - planeLength * sin(dPhi / 2.)) / + (tan(angle) - tan(angle + dPhi / 2.)); + double yIntersect = tan(angle) * xIntersect + Rmin * (sin(dPhi / 2.) - tan(angle)) + planeLength * sin(dPhi / 2.); + // distance from inner radius to intersection + double correction = + planeLength - sqrt(pow(xIntersect - Rmin * cos(dPhi / 2), 2) + pow(yIntersect - Rmin * sin(dPhi / 2), 2)); + // correction to the active thickness + activeOutThickness += 2. * correction * sin(dPhi / 4.); + activeOutThickness -= passiveThickness * (0.5 - activePassiveOverlap); + // print the active layer dimensions + double activeInThicknessAfterSubtraction = + 2. * activeInThickness - readoutThickness - 2. * activePassiveOverlap * passiveThickness; + double activeOutThicknessAfterSubtraction = + 2. * activeOutThickness - readoutThickness - 2. * activePassiveOverlap * + (passiveThickness + passiveInnerThicknessMax - passiveInnerThicknessMin); // correct thickness for trapezoid + lLog << MSG::INFO << "active material = " << activeMaterial + << " active layers thickness at inner radius (cm) = " << activeInThicknessAfterSubtraction + << " thickness at outer radious (cm) = " << activeOutThicknessAfterSubtraction << " making " + << (activeOutThicknessAfterSubtraction - activeInThicknessAfterSubtraction) * 100 / + activeInThicknessAfterSubtraction + << " % increase." << endmsg; + lLog << MSG::INFO + << "active passive initial overlap (before subtraction) (cm) = " << passiveThickness * activePassiveOverlap + << " = " << activePassiveOverlap * 100 << " %" << endmsg; + + // creating shape for rows of layers (active material between two passive planes, with readout in the middle) + // first define area between two passive planes, area can reach up to the symmetry axis of passive plane + dd4hep::Trd1 activeOuterShape(activeInThickness, activeOutThickness, caloDim.dz(), planeLength / 2.); + // subtract readout shape from the middle + dd4hep::SubtractionSolid activeShapeNoReadout(activeOuterShape, readoutShape); + + // make calculation for active plane that is inclined with 0 deg (= offset + angle) + double Cx = Rmin * cos(-angle) + planeLength / 2.; + double Cy = Rmin * sin(-angle); + double Ax = Rmin * cos(-angle + dPhi / 2.) + planeLength / 2. * cos(dPhi / 2.); + double Ay = Rmin * sin(-angle + dPhi / 2.) + planeLength / 2. * sin(dPhi / 2.); + double CAx = fabs(Ax - Cx); + double CAy = fabs(Ay - Cy); + double zprim, xprim; + zprim = CAx; + xprim = CAy; + + double Bx = Rmin * cos(-angle - dPhi / 2.) + planeLength / 2. * cos(-dPhi / 2.); + double By = Rmin * sin(-angle - dPhi / 2.) + planeLength / 2. * sin(-dPhi / 2.); + double CBx = fabs(Bx - Cx); + double CBy = fabs(By - Cy); + double zprimB, xprimB; + zprimB = CBx; + xprimB = CBy; + + // subtract passive volume above + dd4hep::SubtractionSolid activeShapeNoPassiveAbove( + activeShapeNoReadout, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), + dd4hep::Position(-fabs(xprim), 0, fabs(zprim)))); + // subtract passive volume below + dd4hep::SubtractionSolid activeShape( + activeShapeNoPassiveAbove, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), + dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB)))); + dd4hep::Volume activeVol("active", activeShape, aLcdd.material("Air")); + + std::vector layerPhysVols; + // place layers within active volume + std::vector layerInThickness; + std::vector layerOutThickness; + double layerIncreasePerUnitThickness = (activeOutThickness - activeInThickness) / layersTotalHeight; + for (uint iLay = 0; iLay < numLayers; iLay++) { + if (iLay == 0) { + layerInThickness.push_back(activeInThickness); + } else { + layerInThickness.push_back(layerOutThickness[iLay - 1]); + } + layerOutThickness.push_back(layerInThickness[iLay] + layerIncreasePerUnitThickness * layerHeight[iLay]); + } + double layerOffset = layerFirstOffset; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Trd1 layerOuterShape(layerInThickness[iLayer], layerOutThickness[iLayer], caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::SubtractionSolid layerShapeNoReadout(layerOuterShape, readoutShape); + dd4hep::SubtractionSolid layerShapeNoPassiveAbove( + layerShapeNoReadout, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), + dd4hep::Position(-fabs(xprim), 0, fabs(zprim) - layerOffset))); + // subtract passive volume below + dd4hep::SubtractionSolid layerShape( + layerShapeNoPassiveAbove, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), + dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB) - layerOffset))); + dd4hep::Volume layerVol("layer", layerShape, aLcdd.material(activeMaterial)); + layerVol.setSensitiveDetector(aSensDet); + layerPhysVols.push_back(activeVol.placeVolume(layerVol, dd4hep::Position(0, 0, layerOffset))); + layerPhysVols.back().addPhysVolID("layer", iLayer); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + + dd4hep::DetElement bathDetElem(caloDetElem, "bath", 1); + std::vector activePhysVols; + // Next place elements: passive planes, readout planes and rows of layers + for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { + // first calculate positions of passive and readout planes + // PASSIVE + // calculate centre position of the plane without plane rotation + double phi = offsetPassivePhi + iPlane * dPhi; + double xRadial = (Rmin + planeLength / 2.) * cos(phi); + double yRadial = (Rmin + planeLength / 2.) * sin(phi); + // calculate position of the beginning of plane + double xRmin = Rmin * cos(phi); + double yRmin = Rmin * sin(phi); + // rotate centre by angle wrt beginning of plane + double xRotated = xRmin + (xRadial - xRmin) * cos(angle) - (yRadial - yRmin) * sin(angle); + double yRotated = yRmin + (xRadial - xRmin) * sin(angle) + (yRadial - yRmin) * cos(angle); + dd4hep::Transform3D transform(dd4hep::RotationX(-M_PI / 2.) // to get in XY plane + * + dd4hep::RotationY(M_PI / 2. // to get pointed towards centre + - + phi - angle), + dd4hep::Position(xRotated, yRotated, 0)); + dd4hep::PlacedVolume passivePhysVol = bathVol.placeVolume(passiveVol, transform); + passivePhysVol.addPhysVolID("module", iPlane); + passivePhysVol.addPhysVolID("type", 1); // 0 = active, 1 = passive, 2 = readout + dd4hep::DetElement passiveDetElem(bathDetElem, "passive" + std::to_string(iPlane), iPlane); + passiveDetElem.setPlacement(passivePhysVol); + + // READOUT + // calculate centre position of the plane without plane rotation + double phiRead = offsetReadoutPhi + iPlane * dPhi; + double xRadialRead = (Rmin + planeLength / 2.) * cos(phiRead); + double yRadialRead = (Rmin + planeLength / 2.) * sin(phiRead); + // calculate position of the beginning of plane + double xRminRead = Rmin * cos(phiRead); + double yRminRead = Rmin * sin(phiRead); + // rotate centre by angle wrt beginning of plane + double xRotatedRead = xRminRead + (xRadialRead - xRminRead) * cos(angle) - (yRadialRead - yRminRead) * sin(angle); + double yRotatedRead = yRminRead + (xRadialRead - xRminRead) * sin(angle) + (yRadialRead - yRminRead) * cos(angle); + dd4hep::Transform3D transformRead( + dd4hep::RotationX(-M_PI / 2.) // to get in XY plane + * + dd4hep::RotationY(M_PI / 2. // to get pointed towards centre + - + phiRead - angle), + dd4hep::Position(xRotatedRead, yRotatedRead, 0)); + dd4hep::PlacedVolume readoutPhysVol = bathVol.placeVolume(readoutVol, transformRead); + readoutPhysVol.addPhysVolID("module", iPlane); + readoutPhysVol.addPhysVolID("type", 2); // 0 = active, 1 = passive, 2 = readout + dd4hep::DetElement readoutDetElem(bathDetElem, "readout" + std::to_string(iPlane), iPlane); + readoutDetElem.setPlacement(readoutPhysVol); + + // ACTIVE + dd4hep::Rotation3D rotationActive(dd4hep::RotationX(-M_PI / 2) * + dd4hep::RotationY(M_PI / 2 - phiRead - angle)); + activePhysVols.push_back(bathVol.placeVolume( + activeVol, + dd4hep::Transform3D(rotationActive, dd4hep::Position(xRotatedRead, yRotatedRead, 0)))); + activePhysVols.back().addPhysVolID("module", iPlane); + activePhysVols.back().addPhysVolID("type", 0); // 0 = active, 1 = passive, 2 = readout + } + dd4hep::PlacedVolume bathPhysVol = envelopeVol.placeVolume(bathVol); + bathDetElem.setPlacement(bathPhysVol); + for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { + dd4hep::DetElement activeDetElem(bathDetElem, "active" + std::to_string(iPlane), iPlane); + activeDetElem.setPlacement(activePhysVols[iPlane]); + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::DetElement layerDetElem(activeDetElem, "layer" + std::to_string(iLayer), iLayer); + layerDetElem.setPlacement(layerPhysVols[iLayer]); + } + } + + // Place the envelope + dd4hep::Volume motherVol = aLcdd.pickMotherVolume(caloDetElem); + dd4hep::PlacedVolume envelopePhysVol = motherVol.placeVolume(envelopeVol); + envelopePhysVol.addPhysVolID("system", xmlDetElem.id()); + caloDetElem.setPlacement(envelopePhysVol); + + // Create caloData object + auto caloData = new dd4hep::rec::LayeredCalorimeterData; + caloData->layoutType = dd4hep::rec::LayeredCalorimeterData::BarrelLayout; + caloDetElem.addExtension(caloData); + + // Set type flags + dd4hep::xml::setDetectorTypeFlag(xmlDetElem, caloDetElem); + + return caloDetElem; +} +} // namespace det + +DECLARE_DETELEMENT(EmCaloBarrelInclined, det::createECalBarrelInclined) diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h new file mode 100644 index 00000000..12cbafe9 --- /dev/null +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h @@ -0,0 +1,117 @@ +#ifndef DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H +#define DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H + +// FCCSW +#include "DetSegmentation/GridTheta.h" +#include "DD4hep/VolumeManager.h" + +/** FCCSWGridModuleThetaMerged Detector/DetSegmentation/DetSegmentation/FCCSWGridModuleThetaMerged.h FCCSWGridModuleThetaMerged.h + * + * Segmentation in theta and module. + * Based on GridTheta, merges modules and theta cells based on layer number + * + */ + +namespace dd4hep { +namespace DDSegmentation { +class FCCSWGridModuleThetaMerged : public GridTheta { +public: + /// default constructor using an arbitrary type + FCCSWGridModuleThetaMerged(const std::string& aCellEncoding); + /// Default constructor used by derived classes passing an existing decoder + FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder); + + /// destructor + virtual ~FCCSWGridModuleThetaMerged() = default; + + /// read n(modules) from detector metadata + void GetNModulesFromGeom(); + + /// read n(layers) from detector metadata + void GetNLayersFromGeom(); + + /** Determine the local position based on the cell ID. + * @param[in] aCellId ID of a cell. + * return Position (relative to R, phi of Geant4 volume it belongs to, scaled for R=1). + */ + virtual Vector3D position(const CellID& aCellID) const; + /** Determine the cell ID based on the position. + * @param[in] aLocalPosition (not used). + * @param[in] aGlobalPosition + * @param[in] aVolumeId ID of the Geant4 volume + * return Cell ID. + */ + virtual CellID cellID(const Vector3D& aLocalPosition, const Vector3D& aGlobalPosition, + const VolumeID& aVolumeID) const; + /** Determine the azimuthal angle (relative to the G4 volume) based on the cell ID. + * @param[in] aCellId ID of a cell. + * return Phi. + */ + double phi(const CellID& aCellID) const; + /** Determine the polar angle (relative to the G4 volume) based on the cell ID. + * @param[in] aCellId ID of a cell. + * return Theta. + */ + double theta(const CellID& aCellID) const; + /** Determine the radius based on the cell ID. + * @param[in] aCellId ID of a cell. + * return Radius. + */ + // double radius(const CellID& aCellID) const; + /** Get the number of merged cells in theta for given layer + * @param[in] layer + * return The number of merged cells in theta + */ + inline int mergedThetaCells(const int layer) const { + if (layer < (int) m_mergedCellsTheta.size()) + return m_mergedCellsTheta[layer]; + else + return 1; + } + /** Get the number of merged modules (inclined in phi) + * @param[in] layer + * return The number of merged modules + */ + inline int mergedModules(const int layer) const { + if (layer < (int) m_mergedModules.size()) + return m_mergedModules[layer]; + else + return 1; + } + /** Get the total number of modules of detector + * return The number of modules (as it was set by the user in the xml file..) + */ + inline int nModules() const { return m_nModules; } + /** Get the number of layers + * return The number of layers + */ + inline int nLayers() const { return m_nLayers; } + /** Get the field name used for the layer + * return The field name for the layer. + */ + inline const std::string& fieldNameLayer() const { return m_layerID; } + /** Get the field name used for the module number + * return The field name for the module. + */ + inline const std::string& fieldNameModule() const { return m_moduleID; } + +protected: + /// the field name used for layer + std::string m_layerID; + /// the field name used for the read-out module (can differ from module due to merging) + std::string m_moduleID; + /// vector of number of cells to be merged along theta for each layer (typically between 1 and 4) + std::vector m_mergedCellsTheta; + /// vector of number of modules to be merged for each layer (typically 1 or 2) + std::vector m_mergedModules; + + /// number of modules (or, equivalently, the deltaPhi between adjacent modules) + int m_nModules; + + /// number of layers (from the geometry) + int m_nLayers; + +}; +} +} +#endif /* DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H */ diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h new file mode 100644 index 00000000..bf6b6bcb --- /dev/null +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h @@ -0,0 +1,124 @@ +#ifndef DD4HEP_DDCORE_GRIDMODULETHETAMERGED_H +#define DD4HEP_DDCORE_GRIDMODULETHETAMERGED_H 1 + +// FCCSW +#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" + +// DD4hep +#include "DD4hep/Segmentations.h" +#include "DD4hep/detail/SegmentationsInterna.h" + +/// Namespace for the AIDA detector description toolkit +namespace dd4hep { + +/// Namespace for base segmentations + + +// Forward declarations +class Segmentation; +template +class SegmentationWrapper; + +/// We need some abbreviation to make the code more readable. +typedef Handle> FCCSWGridModuleThetaMergedHandle; + +/// Implementation class for the module-theta (with per-layer merging) segmentation. +/** + * Concrete user handle to serve specific needs of client code + * which requires access to the base functionality not served + * by the super-class Segmentation. + * + * Note: + * We only check the validity of the underlying handle. + * If for whatever reason the implementation object is not valid + * This is not checked. + * In principle this CANNOT happen unless some brain-dead has + * fiddled with the handled object directly..... + * + * Note: + * The handle base corrsponding to this object in for + * conveniance reasons instantiated in DD4hep/src/Segmentations.cpp. + * + */ +class FCCSWGridModuleThetaMerged : public FCCSWGridModuleThetaMergedHandle { +public: + /// Definition of the basic handled object + typedef FCCSWGridModuleThetaMergedHandle::Object Object; + +public: + /// Default constructor + FCCSWGridModuleThetaMerged() = default; + /// Copy constructor + FCCSWGridModuleThetaMerged(const FCCSWGridModuleThetaMerged& e) = default; + /// Copy Constructor from segmentation base object + FCCSWGridModuleThetaMerged(const Segmentation& e) : Handle(e) {} + /// Copy constructor from handle + FCCSWGridModuleThetaMerged(const Handle& e) : Handle(e) {} + /// Copy constructor from other polymorph/equivalent handle + template + FCCSWGridModuleThetaMerged(const Handle& e) : Handle(e) {} + /// Assignment operator + FCCSWGridModuleThetaMerged& operator=(const FCCSWGridModuleThetaMerged& seg) = default; + /// Equality operator + bool operator==(const FCCSWGridModuleThetaMerged& seg) const { return m_element == seg.m_element; } + + /// determine the position based on the cell ID + inline Position position(const CellID& id) const { return Position(access()->implementation->position(id)); } + + /// determine the cell ID based on the position + inline dd4hep::CellID cellID(const Position& local, const Position& global, const VolumeID& volID) const { + return access()->implementation->cellID(local, global, volID); + } + + /// access the grid size in theta + inline double gridSizeTheta() const { return access()->implementation->gridSizeTheta(); } + + /// access the coordinate offset in theta + inline double offsetTheta() const { return access()->implementation->offsetTheta(); } + + /// access the number of theta cells merged for given layer + inline int mergedThetaCells(const int layer) const { return access()->implementation->mergedThetaCells(layer); } + + /// access the number of modules + inline int nModules() const { return access()->implementation->nModules(); } + + /// access the number of layers + inline int nLayers() const { return access()->implementation->nLayers(); } + + /// access the number of merged modules for given layer + inline int mergedModules(const int layer) const { return access()->implementation->mergedModules(layer); } + + /// set the coordinate offset in theta + inline void setOffsetTheta(double offset) const { access()->implementation->setOffsetTheta(offset); } + + /// set the grid size in theta + inline void setGridSizeTheta(double cellSize) const { access()->implementation->setGridSizeTheta(cellSize); } + + /// access the field name used for theta + inline const std::string& fieldNameTheta() const { return access()->implementation->fieldNameTheta(); } + + /// access the field name used for module + inline const std::string& fieldNameModule() const { return access()->implementation->fieldNameModule(); } + + /// access the field name used for layer + inline const std::string& fieldNameLayer() const { return access()->implementation->fieldNameLayer(); } + + + /** \brief Returns a std::vector of the cellDimensions of the given cell ID + in natural order of dimensions (nModules, dTheta) + + Returns a std::vector of the cellDimensions of the given cell ID + \param cellID + \return std::vector size 2: + -# size in module + -# size in theta + */ + inline std::vector cellDimensions(const CellID& /*id*/) const { + //return {access()->implementation->gridSizePhi(), access()->implementation->gridSizeTheta()}; + // not implemented + return {0., 0.}; + } +}; + +} /* End namespace dd4hep */ +#endif // DD4HEP_DDCORE_GRIDMODULETHETAMERGED_H diff --git a/detectorSegmentations/include/detectorSegmentations/GridTheta.h b/detectorSegmentations/include/detectorSegmentations/GridTheta.h new file mode 100644 index 00000000..5707d2e7 --- /dev/null +++ b/detectorSegmentations/include/detectorSegmentations/GridTheta.h @@ -0,0 +1,102 @@ +#ifndef DETSEGMENTATION_GRIDTHETA_H +#define DETSEGMENTATION_GRIDTHETA_H + +#include "DDSegmentation/Segmentation.h" + +/* #include "DDSegmentation/SegmentationUtil.h" */ +#include "TVector3.h" +#include + +/** GridTheta Detector/DetSegmentation/DetSegmentation/GridTheta.h GridTheta.h + * + * Segmentation in theta. + * + */ + +namespace dd4hep { +namespace DDSegmentation { +class GridTheta : public Segmentation { +public: + /// default constructor using an arbitrary type + GridTheta(const std::string& aCellEncoding); + /// Default constructor used by derived classes passing an existing decoder + GridTheta(const BitFieldCoder* decoder); + /// destructor + virtual ~GridTheta() = default; + + /** Determine the global position based on the cell ID. + * @warning This segmentation has no knowledge of radius, so radius = 1 is taken into calculations. + * @param[in] aCellId ID of a cell. + * return Position (radius = 1). + */ + virtual Vector3D position(const CellID& aCellID) const; + /** Determine the cell ID based on the position. + * @param[in] aLocalPosition (not used). + * @param[in] aGlobalPosition position in the global coordinates. + * @param[in] aVolumeId ID of a volume. + * return Cell ID. + */ + virtual CellID cellID(const Vector3D& aLocalPosition, const Vector3D& aGlobalPosition, + const VolumeID& aVolumeID) const; + /** Determine the theta angle based on the cell ID. + * @param[in] aCellId ID of a cell. + * return Pseudorapidity. + */ + double theta(const CellID& aCellID) const; + /** Get the grid size in theta angle. + * return Grid size in theta. + */ + inline double gridSizeTheta() const { return m_gridSizeTheta; } + /** Get the coordinate offset in theta angle. + * return The offset in theta. + */ + inline double offsetTheta() const { return m_offsetTheta; } + /** Get the field name used for theta angle + * return The field name for theta. + */ + inline const std::string& fieldNameTheta() const { return m_thetaID; } + /** Set the grid size in theta angle. + * @param[in] aCellSize Cell size in theta. + */ + inline void setGridSizeTheta(double aCellSize) { m_gridSizeTheta = aCellSize; } + /** Set the coordinate offset in theta angle. + * @param[in] aOffset Offset in theta. + */ + inline void setOffsetTheta(double offset) { m_offsetTheta = offset; } + /** Set the field name used for theta angle. + * @param[in] aFieldName Field name for theta. + */ + inline void setFieldNameTheta(const std::string& fieldName) { m_thetaID = fieldName; } + /// calculates the Cartesian position from cylindrical coordinates (r, phi, theta) + inline Vector3D positionFromRThetaPhi(double ar, double atheta, double aphi) const { + return Vector3D(ar * std::cos(aphi), ar * std::sin(aphi), ar * std::cos(atheta)/std::sin(atheta)); + } + /// calculates the theta angle from Cartesian coordinates + inline double thetaFromXYZ(const Vector3D& aposition) const { + TVector3 vec(aposition.X, aposition.Y, aposition.Z); + return vec.Theta(); + } + /// from SegmentationUtil + /// to be removed once SegmentationUtil can be included w/o linker error + /// calculates the azimuthal angle phi from Cartesian coordinates + inline double phiFromXYZ(const Vector3D& aposition) const { return std::atan2(aposition.Y, aposition.X); } + /// from SegmentationUtil + /// to be removed once SegmentationUtil can be included w/o linker error + /// calculates the radius in the xy-plane from Cartesian coordinates + inline double radiusFromXYZ(const Vector3D& aposition) const { + return std::sqrt(aposition.X * aposition.X + aposition.Y * aposition.Y); + } + +protected: + /// determine the theta angle based on the current cell ID + double theta() const; + /// the grid size in theta + double m_gridSizeTheta; + /// the coordinate offset in theta + double m_offsetTheta; + /// the field name used for theta + std::string m_thetaID; +}; +} +} +#endif /* DETSEGMENTATION_GRIDTHETA_H */ diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp new file mode 100644 index 00000000..8feb2681 --- /dev/null +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp @@ -0,0 +1,156 @@ +#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" + +#include +#include "DD4hep/Detector.h" + +namespace dd4hep { +namespace DDSegmentation { + +/// default constructor using an encoding string +FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEncoding) : GridTheta(cellEncoding) { + // define type and description + _type = "FCCSWGridModuleThetaMerged"; + _description = "Module-theta segmentation with per-layer merging along theta and/or module"; + + // register all necessary parameters (additional to those registered in GridTheta) + registerIdentifier("identifier_layer", "Cell ID identifier for layer", m_layerID, "layer"); + registerIdentifier("identifier_module", "Cell ID identifier for readout module", m_moduleID, "module"); + registerParameter("mergedCells_Theta", "Numbers of merged cells in theta per layer", m_mergedCellsTheta, std::vector()); + registerParameter("mergedModules", "Numbers of merged modules per layer", m_mergedModules, std::vector()); + GetNModulesFromGeom(); + GetNLayersFromGeom(); +} + +FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder) : GridTheta(decoder) { + // define type and description + _type = "FCCSWGridModuleThetaMerged"; + _description = "Module-theta segmentation with per-layer merging along theta and/or module"; + + // register all necessary parameters (additional to those registered in GridTheta) + registerIdentifier("identifier_layer", "Cell ID identifier for layer", m_layerID, "layer"); + registerIdentifier("identifier_module", "Cell ID identifier for module", m_moduleID, "module"); + registerParameter("mergedCells_Theta", "Numbers of merged cells in theta per layer", m_mergedCellsTheta, std::vector()); + registerParameter("mergedModules", "Numbers of merged modules per layer", m_mergedModules, std::vector()); + GetNModulesFromGeom(); + GetNLayersFromGeom(); +} + +void FCCSWGridModuleThetaMerged::GetNModulesFromGeom() { + dd4hep::Detector* dd4hepgeo = &(dd4hep::Detector::getInstance()); + try { + m_nModules = dd4hepgeo->constant("ECalBarrelNumPlanes"); + } + catch(...) { + std::cout << "Number of modules not found in detector metadata, exiting..." << std::endl; + exit(1); + } + std::cout << "Number of modules read from detector metadata and used in readout class: " << m_nModules << std::endl; +} + +void FCCSWGridModuleThetaMerged::GetNLayersFromGeom() { + dd4hep::Detector* dd4hepgeo = &(dd4hep::Detector::getInstance()); + try { + m_nLayers = dd4hepgeo->constant("ECalBarrelNumLayers"); + } + catch(...) { + std::cout << "Number of layers not found in detector metadata, exiting..." << std::endl; + exit(1); + } + std::cout << "Number of layers read from detector metadata and used in readout class: " << m_nLayers << std::endl; +} + +/// determine the local position based on the cell ID +Vector3D FCCSWGridModuleThetaMerged::position(const CellID& cID) const { + + // debug + // std::cout << "cellID: " << cID << std::endl; + + // cannot return for R=0 otherwise will lose phi info, return for R=1 + return positionFromRThetaPhi(1.0, theta(cID), phi(cID)); +} + +/// determine the cell ID based on the global position +CellID FCCSWGridModuleThetaMerged::cellID(const Vector3D& /* localPosition */, const Vector3D& globalPosition, + const VolumeID& vID) const { + CellID cID = vID; + + // retrieve layer (since merging depends on layer) + int layer = _decoder->get(vID, m_layerID); + + // retrieve theta + double lTheta = thetaFromXYZ(globalPosition); + + // calculate theta bin with original segmentation + int thetaBin = positionToBin(lTheta, m_gridSizeTheta, m_offsetTheta); + + // adjust theta bin if cells are merged along theta in this layer + // assume that m_mergedCellsTheta[layer]>=1 + thetaBin -= (thetaBin % m_mergedCellsTheta[layer]); + + // set theta field of cellID + _decoder->set(cID, m_thetaID, thetaBin); + + // retrieve module number + int module = _decoder->get(vID, m_moduleID); + + // adjust module number if modules are merged in this layer + // assume that m_mergedModules[layer]>=1 + module -= (module % m_mergedModules[layer]); + + // set module field of cellID + _decoder->set(cID, m_moduleID, module); + + return cID; +} + +/// determine the azimuth based on the cell ID +/// the value returned is the relative shift in phi +/// with respect to the first module in the group of +/// merged ones - which will be then added on top of +/// the phi of the volume containing the first cell +/// by the positioning tool +double FCCSWGridModuleThetaMerged::phi(const CellID& cID) const { + + // retrieve layer + int layer = _decoder->get(cID, m_layerID); + + // calculate phi offset due to merging + // assume that m_mergedModules[layer]>=1 + double phi = (m_mergedModules[layer]-1) * M_PI / m_nModules ; + + // debug + // std::cout << "layer: " << layer << std::endl; + // std::cout << "merged modules: " << m_mergedModules[layer] << std::endl; + // std::cout << "phi: " << phi << std::endl; + + return phi; +} + +/// determine the polar angle based on the cell ID and the +/// number of merged theta cells +double FCCSWGridModuleThetaMerged::theta(const CellID& cID) const { + + // retrieve layer + int layer = _decoder->get(cID, m_layerID); + + // retrieve theta bin from cellID and determine theta position + CellID thetaValue = _decoder->get(cID, m_thetaID); + double _theta = binToPosition(thetaValue, m_gridSizeTheta, m_offsetTheta); + + // adjust return value if cells are merged along theta in this layer + // shift by (N-1)*half theta grid size + // assume that m_mergedCellsTheta[layer]>=1 + _theta += (m_mergedCellsTheta[layer]-1) * m_gridSizeTheta / 2.0 ; + + // debug + // std::cout << "layer: " << layer << std::endl; + // std::cout << "theta bin: " << thetaValue << std::endl; + // std::cout << "gridSizeTheta, offsetTheta: " << m_gridSizeTheta << " , " << m_offsetTheta << std::endl; + // std::cout << "merged cells: " << m_mergedCellsTheta[layer] << std::endl; + // std::cout << "theta: " << _theta << std::endl; + + return _theta; +} + +} +} diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp new file mode 100644 index 00000000..e4b49941 --- /dev/null +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp @@ -0,0 +1,4 @@ +#include "DetSegmentation/FCCSWGridModuleThetaMergedHandle.h" +#include "DD4hep/detail/Handle.inl" + +DD4HEP_INSTANTIATE_HANDLE_UNNAMED(dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged); diff --git a/detectorSegmentations/src/GridEta.cpp b/detectorSegmentations/src/GridEta.cpp new file mode 100644 index 00000000..232d0b50 --- /dev/null +++ b/detectorSegmentations/src/GridEta.cpp @@ -0,0 +1,57 @@ +#include "DetSegmentation/GridEta.h" +#include "DD4hep/Factories.h" + +namespace dd4hep { +namespace DDSegmentation { + +/// default constructor using an encoding string +GridEta::GridEta(const std::string& cellEncoding) : Segmentation(cellEncoding) { + // define type and description + _type = "GridEta"; + _description = "Eeta segmentation in the global coordinates"; + + // register all necessary parameters + registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); + registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); + registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); +} + +GridEta::GridEta(const BitFieldCoder* decoder) : Segmentation(decoder) { + // define type and description + _type = "GridEta"; + _description = "Eta segmentation in the global coordinates"; + + // register all necessary parameters + registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); + registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); + registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); +} + +/// determine the local based on the cell ID +Vector3D GridEta::position(const CellID& cID) const { + return positionFromREtaPhi(1.0, eta(cID), 0.); +} + +/// determine the cell ID based on the position +CellID GridEta::cellID(const Vector3D& /* localPosition */, const Vector3D& globalPosition, const VolumeID& vID) const { + CellID cID = vID; + double lEta = etaFromXYZ(globalPosition); + _decoder->set(cID, m_etaID, positionToBin(lEta, m_gridSizeEta, m_offsetEta)); + return cID; +} + +/// determine the pseudorapidity based on the current cell ID +//double GridEta::eta() const { +// CellID etaValue = (*_decoder)[m_etaID].value(); +// return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); +//} + +/// determine the pseudorapidity based on the cell ID +double GridEta::eta(const CellID& cID) const { + CellID etaValue = _decoder->get(cID, m_etaID); + return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); +} + +} // namespace DDSegmentation +} // namespace dd4hep + diff --git a/detectorSegmentations/src/plugins/SegmentationFactories.cpp b/detectorSegmentations/src/plugins/SegmentationFactories.cpp index 5591f66c..7a1f8c15 100644 --- a/detectorSegmentations/src/plugins/SegmentationFactories.cpp +++ b/detectorSegmentations/src/plugins/SegmentationFactories.cpp @@ -8,21 +8,24 @@ dd4hep::SegmentationObject* create_segmentation(const dd4hep::BitFieldCoder* dec } } -#include "detectorSegmentations/GridEta_k4geo.h" -DECLARE_SEGMENTATION(GridEta_k4geo, create_segmentation) +#include "DetSegmentation/GridEta.h" +DECLARE_SEGMENTATION(GridEta, create_segmentation) -#include "detectorSegmentations/GridTheta_k4geo.h" -DECLARE_SEGMENTATION(GridTheta_k4geo, create_segmentation) +#include "DetSegmentation/GridTheta.h" +DECLARE_SEGMENTATION(GridTheta, create_segmentation) -#include "detectorSegmentations/FCCSWGridPhiTheta_k4geo.h" -DECLARE_SEGMENTATION(FCCSWGridPhiTheta_k4geo, create_segmentation) +#include "DetSegmentation/FCCSWGridPhiTheta.h" +DECLARE_SEGMENTATION(FCCSWGridPhiTheta, create_segmentation) -#include "detectorSegmentations/FCCSWGridPhiEta_k4geo.h" -DECLARE_SEGMENTATION(FCCSWGridPhiEta_k4geo, create_segmentation) +#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" +DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged, create_segmentation) -#include "detectorSegmentations/GridRPhiEta_k4geo.h" -DECLARE_SEGMENTATION(GridRPhiEta_k4geo, create_segmentation) +#include "DetSegmentation/FCCSWGridPhiEta.h" +DECLARE_SEGMENTATION(FCCSWGridPhiEta, create_segmentation) -#include "detectorSegmentations/GridSimplifiedDriftChamber_k4geo.h" -DECLARE_SEGMENTATION(GridSimplifiedDriftChamber_k4geo, create_segmentation) +#include "DetSegmentation/GridRPhiEta.h" +DECLARE_SEGMENTATION(GridRPhiEta, create_segmentation) + +#include "DetSegmentation/GridSimplifiedDriftChamber.h" +DECLARE_SEGMENTATION(GridSimplifiedDriftChamber, create_segmentation) From 0649533f53d3b1893859c4fd17a9fbba2ce72dd7 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 26 Oct 2023 11:50:45 +0200 Subject: [PATCH 02/15] path for include fixed, file included in previous commit removed --- .../FCCSWGridModuleThetaMerged.h | 2 +- .../FCCSWGridModuleThetaMergedHandle.h | 2 +- .../src/FCCSWGridModuleThetaMerged.cpp | 2 +- .../src/FCCSWGridModuleThetaMergedHandle.cpp | 2 +- detectorSegmentations/src/GridEta.cpp | 57 ------------------- detectorSegmentations/src/GridEta_k4geo.cpp | 8 +-- 6 files changed, 8 insertions(+), 65 deletions(-) delete mode 100644 detectorSegmentations/src/GridEta.cpp diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h index 12cbafe9..09b955b6 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h @@ -2,7 +2,7 @@ #define DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H // FCCSW -#include "DetSegmentation/GridTheta.h" +#include "detectorSegmentations/GridTheta.h" #include "DD4hep/VolumeManager.h" /** FCCSWGridModuleThetaMerged Detector/DetSegmentation/DetSegmentation/FCCSWGridModuleThetaMerged.h FCCSWGridModuleThetaMerged.h diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h index bf6b6bcb..f7bd05ea 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h @@ -2,7 +2,7 @@ #define DD4HEP_DDCORE_GRIDMODULETHETAMERGED_H 1 // FCCSW -#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" // DD4hep #include "DD4hep/Segmentations.h" diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp index 8feb2681..6c394589 100644 --- a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp @@ -1,4 +1,4 @@ -#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" #include #include "DD4hep/Detector.h" diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp index e4b49941..0573a2b9 100644 --- a/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp @@ -1,4 +1,4 @@ -#include "DetSegmentation/FCCSWGridModuleThetaMergedHandle.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h" #include "DD4hep/detail/Handle.inl" DD4HEP_INSTANTIATE_HANDLE_UNNAMED(dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged); diff --git a/detectorSegmentations/src/GridEta.cpp b/detectorSegmentations/src/GridEta.cpp deleted file mode 100644 index 232d0b50..00000000 --- a/detectorSegmentations/src/GridEta.cpp +++ /dev/null @@ -1,57 +0,0 @@ -#include "DetSegmentation/GridEta.h" -#include "DD4hep/Factories.h" - -namespace dd4hep { -namespace DDSegmentation { - -/// default constructor using an encoding string -GridEta::GridEta(const std::string& cellEncoding) : Segmentation(cellEncoding) { - // define type and description - _type = "GridEta"; - _description = "Eeta segmentation in the global coordinates"; - - // register all necessary parameters - registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); - registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); - registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); -} - -GridEta::GridEta(const BitFieldCoder* decoder) : Segmentation(decoder) { - // define type and description - _type = "GridEta"; - _description = "Eta segmentation in the global coordinates"; - - // register all necessary parameters - registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); - registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); - registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); -} - -/// determine the local based on the cell ID -Vector3D GridEta::position(const CellID& cID) const { - return positionFromREtaPhi(1.0, eta(cID), 0.); -} - -/// determine the cell ID based on the position -CellID GridEta::cellID(const Vector3D& /* localPosition */, const Vector3D& globalPosition, const VolumeID& vID) const { - CellID cID = vID; - double lEta = etaFromXYZ(globalPosition); - _decoder->set(cID, m_etaID, positionToBin(lEta, m_gridSizeEta, m_offsetEta)); - return cID; -} - -/// determine the pseudorapidity based on the current cell ID -//double GridEta::eta() const { -// CellID etaValue = (*_decoder)[m_etaID].value(); -// return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); -//} - -/// determine the pseudorapidity based on the cell ID -double GridEta::eta(const CellID& cID) const { - CellID etaValue = _decoder->get(cID, m_etaID); - return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); -} - -} // namespace DDSegmentation -} // namespace dd4hep - diff --git a/detectorSegmentations/src/GridEta_k4geo.cpp b/detectorSegmentations/src/GridEta_k4geo.cpp index 064613a0..5c387e5e 100644 --- a/detectorSegmentations/src/GridEta_k4geo.cpp +++ b/detectorSegmentations/src/GridEta_k4geo.cpp @@ -11,7 +11,7 @@ GridEta_k4geo::GridEta_k4geo(const std::string& cellEncoding) : Segmentation(cel _description = "Eeta segmentation in the global coordinates"; // register all necessary parameters - registerParameter("grid_size_eta", "Cell size in Eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); + registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); } @@ -19,10 +19,10 @@ GridEta_k4geo::GridEta_k4geo(const std::string& cellEncoding) : Segmentation(cel GridEta_k4geo::GridEta_k4geo(const BitFieldCoder* decoder) : Segmentation(decoder) { // define type and description _type = "GridEta_k4geo"; - _description = "Eeta segmentation in the global coordinates"; + _description = "Eta segmentation in the global coordinates"; // register all necessary parameters - registerParameter("grid_size_eta", "Cell size in Eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); + registerParameter("grid_size_eta", "Cell size in eta", m_gridSizeEta, 1., SegmentationParameter::LengthUnit); registerParameter("offset_eta", "Angular offset in eta", m_offsetEta, 0., SegmentationParameter::AngleUnit, true); registerIdentifier("identifier_eta", "Cell ID identifier for eta", m_etaID, "eta"); } @@ -46,7 +46,7 @@ CellID GridEta_k4geo::cellID(const Vector3D& /* localPosition */, const Vector3D // return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); //} -/// determine the polar angle theta based on the cell ID +/// determine the seudorapidity based on the cell ID double GridEta_k4geo::eta(const CellID& cID) const { CellID etaValue = _decoder->get(cID, m_etaID); return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); From 40df04825018a2764af0c8cefa13e80c7fd4a269 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 26 Oct 2023 12:14:27 +0200 Subject: [PATCH 03/15] it compiles --- .../FCCSWGridModuleThetaMerged.h | 4 +- .../include/detectorSegmentations/GridTheta.h | 102 ------------------ .../detectorSegmentations/GridTheta_k4geo.h | 2 +- .../src/FCCSWGridModuleThetaMerged.cpp | 8 +- .../src/plugins/SegmentationFactories.cpp | 27 ++--- 5 files changed, 21 insertions(+), 122 deletions(-) delete mode 100644 detectorSegmentations/include/detectorSegmentations/GridTheta.h diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h index 09b955b6..2b3fccb7 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h @@ -2,7 +2,7 @@ #define DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H // FCCSW -#include "detectorSegmentations/GridTheta.h" +#include "detectorSegmentations/GridTheta_k4geo.h" #include "DD4hep/VolumeManager.h" /** FCCSWGridModuleThetaMerged Detector/DetSegmentation/DetSegmentation/FCCSWGridModuleThetaMerged.h FCCSWGridModuleThetaMerged.h @@ -14,7 +14,7 @@ namespace dd4hep { namespace DDSegmentation { -class FCCSWGridModuleThetaMerged : public GridTheta { +class FCCSWGridModuleThetaMerged : public GridTheta_k4geo { public: /// default constructor using an arbitrary type FCCSWGridModuleThetaMerged(const std::string& aCellEncoding); diff --git a/detectorSegmentations/include/detectorSegmentations/GridTheta.h b/detectorSegmentations/include/detectorSegmentations/GridTheta.h deleted file mode 100644 index 5707d2e7..00000000 --- a/detectorSegmentations/include/detectorSegmentations/GridTheta.h +++ /dev/null @@ -1,102 +0,0 @@ -#ifndef DETSEGMENTATION_GRIDTHETA_H -#define DETSEGMENTATION_GRIDTHETA_H - -#include "DDSegmentation/Segmentation.h" - -/* #include "DDSegmentation/SegmentationUtil.h" */ -#include "TVector3.h" -#include - -/** GridTheta Detector/DetSegmentation/DetSegmentation/GridTheta.h GridTheta.h - * - * Segmentation in theta. - * - */ - -namespace dd4hep { -namespace DDSegmentation { -class GridTheta : public Segmentation { -public: - /// default constructor using an arbitrary type - GridTheta(const std::string& aCellEncoding); - /// Default constructor used by derived classes passing an existing decoder - GridTheta(const BitFieldCoder* decoder); - /// destructor - virtual ~GridTheta() = default; - - /** Determine the global position based on the cell ID. - * @warning This segmentation has no knowledge of radius, so radius = 1 is taken into calculations. - * @param[in] aCellId ID of a cell. - * return Position (radius = 1). - */ - virtual Vector3D position(const CellID& aCellID) const; - /** Determine the cell ID based on the position. - * @param[in] aLocalPosition (not used). - * @param[in] aGlobalPosition position in the global coordinates. - * @param[in] aVolumeId ID of a volume. - * return Cell ID. - */ - virtual CellID cellID(const Vector3D& aLocalPosition, const Vector3D& aGlobalPosition, - const VolumeID& aVolumeID) const; - /** Determine the theta angle based on the cell ID. - * @param[in] aCellId ID of a cell. - * return Pseudorapidity. - */ - double theta(const CellID& aCellID) const; - /** Get the grid size in theta angle. - * return Grid size in theta. - */ - inline double gridSizeTheta() const { return m_gridSizeTheta; } - /** Get the coordinate offset in theta angle. - * return The offset in theta. - */ - inline double offsetTheta() const { return m_offsetTheta; } - /** Get the field name used for theta angle - * return The field name for theta. - */ - inline const std::string& fieldNameTheta() const { return m_thetaID; } - /** Set the grid size in theta angle. - * @param[in] aCellSize Cell size in theta. - */ - inline void setGridSizeTheta(double aCellSize) { m_gridSizeTheta = aCellSize; } - /** Set the coordinate offset in theta angle. - * @param[in] aOffset Offset in theta. - */ - inline void setOffsetTheta(double offset) { m_offsetTheta = offset; } - /** Set the field name used for theta angle. - * @param[in] aFieldName Field name for theta. - */ - inline void setFieldNameTheta(const std::string& fieldName) { m_thetaID = fieldName; } - /// calculates the Cartesian position from cylindrical coordinates (r, phi, theta) - inline Vector3D positionFromRThetaPhi(double ar, double atheta, double aphi) const { - return Vector3D(ar * std::cos(aphi), ar * std::sin(aphi), ar * std::cos(atheta)/std::sin(atheta)); - } - /// calculates the theta angle from Cartesian coordinates - inline double thetaFromXYZ(const Vector3D& aposition) const { - TVector3 vec(aposition.X, aposition.Y, aposition.Z); - return vec.Theta(); - } - /// from SegmentationUtil - /// to be removed once SegmentationUtil can be included w/o linker error - /// calculates the azimuthal angle phi from Cartesian coordinates - inline double phiFromXYZ(const Vector3D& aposition) const { return std::atan2(aposition.Y, aposition.X); } - /// from SegmentationUtil - /// to be removed once SegmentationUtil can be included w/o linker error - /// calculates the radius in the xy-plane from Cartesian coordinates - inline double radiusFromXYZ(const Vector3D& aposition) const { - return std::sqrt(aposition.X * aposition.X + aposition.Y * aposition.Y); - } - -protected: - /// determine the theta angle based on the current cell ID - double theta() const; - /// the grid size in theta - double m_gridSizeTheta; - /// the coordinate offset in theta - double m_offsetTheta; - /// the field name used for theta - std::string m_thetaID; -}; -} -} -#endif /* DETSEGMENTATION_GRIDTHETA_H */ diff --git a/detectorSegmentations/include/detectorSegmentations/GridTheta_k4geo.h b/detectorSegmentations/include/detectorSegmentations/GridTheta_k4geo.h index 038fdf06..84236c9d 100644 --- a/detectorSegmentations/include/detectorSegmentations/GridTheta_k4geo.h +++ b/detectorSegmentations/include/detectorSegmentations/GridTheta_k4geo.h @@ -69,7 +69,7 @@ class GridTheta_k4geo : public Segmentation { inline void setFieldNameTheta(const std::string& fieldName) { m_thetaID = fieldName; } /// calculates the Cartesian position from cylindrical coordinates (r, phi, theta) inline Vector3D positionFromRThetaPhi(double ar, double atheta, double aphi) const { - return Vector3D(ar * std::cos(aphi), ar * std::sin(aphi), ar / std::tan(atheta)); + return Vector3D(ar * std::cos(aphi), ar * std::sin(aphi), ar * std::cos(atheta)/std::sin(atheta)); } /// calculates the theta angle from Cartesian coordinates inline double thetaFromXYZ(const Vector3D& aposition) const { diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp index 6c394589..1f1763ad 100644 --- a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp @@ -7,12 +7,12 @@ namespace dd4hep { namespace DDSegmentation { /// default constructor using an encoding string -FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEncoding) : GridTheta(cellEncoding) { +FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEncoding) : GridTheta_k4geo(cellEncoding) { // define type and description _type = "FCCSWGridModuleThetaMerged"; _description = "Module-theta segmentation with per-layer merging along theta and/or module"; - // register all necessary parameters (additional to those registered in GridTheta) + // register all necessary parameters (additional to those registered in GridTheta_k4geo) registerIdentifier("identifier_layer", "Cell ID identifier for layer", m_layerID, "layer"); registerIdentifier("identifier_module", "Cell ID identifier for readout module", m_moduleID, "module"); registerParameter("mergedCells_Theta", "Numbers of merged cells in theta per layer", m_mergedCellsTheta, std::vector()); @@ -21,12 +21,12 @@ FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEn GetNLayersFromGeom(); } -FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder) : GridTheta(decoder) { +FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder) : GridTheta_k4geo(decoder) { // define type and description _type = "FCCSWGridModuleThetaMerged"; _description = "Module-theta segmentation with per-layer merging along theta and/or module"; - // register all necessary parameters (additional to those registered in GridTheta) + // register all necessary parameters (additional to those registered in GridTheta_k4geo) registerIdentifier("identifier_layer", "Cell ID identifier for layer", m_layerID, "layer"); registerIdentifier("identifier_module", "Cell ID identifier for module", m_moduleID, "module"); registerParameter("mergedCells_Theta", "Numbers of merged cells in theta per layer", m_mergedCellsTheta, std::vector()); diff --git a/detectorSegmentations/src/plugins/SegmentationFactories.cpp b/detectorSegmentations/src/plugins/SegmentationFactories.cpp index 7a1f8c15..cc308c42 100644 --- a/detectorSegmentations/src/plugins/SegmentationFactories.cpp +++ b/detectorSegmentations/src/plugins/SegmentationFactories.cpp @@ -8,24 +8,25 @@ dd4hep::SegmentationObject* create_segmentation(const dd4hep::BitFieldCoder* dec } } -#include "DetSegmentation/GridEta.h" -DECLARE_SEGMENTATION(GridEta, create_segmentation) +#include "detectorSegmentations/GridEta_k4geo.h" +DECLARE_SEGMENTATION(GridEta_k4geo, create_segmentation) -#include "DetSegmentation/GridTheta.h" -DECLARE_SEGMENTATION(GridTheta, create_segmentation) +#include "detectorSegmentations/GridTheta_k4geo.h" +DECLARE_SEGMENTATION(GridTheta_k4geo, create_segmentation) -#include "DetSegmentation/FCCSWGridPhiTheta.h" -DECLARE_SEGMENTATION(FCCSWGridPhiTheta, create_segmentation) +#include "detectorSegmentations/FCCSWGridPhiTheta_k4geo.h" +DECLARE_SEGMENTATION(FCCSWGridPhiTheta_k4geo, create_segmentation) -#include "DetSegmentation/FCCSWGridModuleThetaMerged.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged, create_segmentation) -#include "DetSegmentation/FCCSWGridPhiEta.h" -DECLARE_SEGMENTATION(FCCSWGridPhiEta, create_segmentation) +#include "detectorSegmentations/FCCSWGridPhiEta_k4geo.h" +DECLARE_SEGMENTATION(FCCSWGridPhiEta_k4geo, create_segmentation) -#include "DetSegmentation/GridRPhiEta.h" -DECLARE_SEGMENTATION(GridRPhiEta, create_segmentation) +#include "detectorSegmentations/GridRPhiEta_k4geo.h" +DECLARE_SEGMENTATION(GridRPhiEta_k4geo, create_segmentation) + +#include "detectorSegmentations/GridSimplifiedDriftChamber_k4geo.h" +DECLARE_SEGMENTATION(GridSimplifiedDriftChamber_k4geo, create_segmentation) -#include "DetSegmentation/GridSimplifiedDriftChamber.h" -DECLARE_SEGMENTATION(GridSimplifiedDriftChamber, create_segmentation) From f34cb12c6339ef5ba618d06121defcbb3b0b33f6 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 26 Oct 2023 12:24:26 +0200 Subject: [PATCH 04/15] new segmentation class refactored by appending _k4geo to name --- .../FCCee_ECalBarrel_thetamodulemerged.xml | 4 ++-- ... FCCSWGridModuleThetaMergedHandle_k4geo.h} | 22 +++++++++---------- ...d.h => FCCSWGridModuleThetaMerged_k4geo.h} | 16 +++++++------- ...CCSWGridModuleThetaMergedHandle_k4geo.cpp} | 4 ++-- ...p => FCCSWGridModuleThetaMerged_k4geo.cpp} | 22 +++++++++---------- .../src/plugins/SegmentationFactories.cpp | 4 ++-- 6 files changed, 36 insertions(+), 36 deletions(-) rename detectorSegmentations/include/detectorSegmentations/{FCCSWGridModuleThetaMergedHandle.h => FCCSWGridModuleThetaMergedHandle_k4geo.h} (82%) rename detectorSegmentations/include/detectorSegmentations/{FCCSWGridModuleThetaMerged.h => FCCSWGridModuleThetaMerged_k4geo.h} (86%) rename detectorSegmentations/src/{FCCSWGridModuleThetaMergedHandle.cpp => FCCSWGridModuleThetaMergedHandle_k4geo.cpp} (80%) rename detectorSegmentations/src/{FCCSWGridModuleThetaMerged.cpp => FCCSWGridModuleThetaMerged_k4geo.cpp} (85%) diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml index 9411b7dc..a454914d 100644 --- a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml @@ -84,13 +84,13 @@ - + system:4,cryo:1,type:3,subtype:3,layer:8,module:11,theta:10 - + system:4,cryo:1,type:3,subtype:3,layer:8,module:11,theta:10 diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h similarity index 82% rename from detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h rename to detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h index f7bd05ea..43d17740 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h @@ -2,7 +2,7 @@ #define DD4HEP_DDCORE_GRIDMODULETHETAMERGED_H 1 // FCCSW -#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h" // DD4hep #include "DD4hep/Segmentations.h" @@ -20,7 +20,7 @@ template class SegmentationWrapper; /// We need some abbreviation to make the code more readable. -typedef Handle> FCCSWGridModuleThetaMergedHandle; +typedef Handle> FCCSWGridModuleThetaMergedHandle_k4geo; /// Implementation class for the module-theta (with per-layer merging) segmentation. /** @@ -40,27 +40,27 @@ typedef Handle> * conveniance reasons instantiated in DD4hep/src/Segmentations.cpp. * */ -class FCCSWGridModuleThetaMerged : public FCCSWGridModuleThetaMergedHandle { +class FCCSWGridModuleThetaMerged_k4geo : public FCCSWGridModuleThetaMergedHandle_k4geo { public: /// Definition of the basic handled object - typedef FCCSWGridModuleThetaMergedHandle::Object Object; + typedef FCCSWGridModuleThetaMergedHandle_k4geo::Object Object; public: /// Default constructor - FCCSWGridModuleThetaMerged() = default; + FCCSWGridModuleThetaMerged_k4geo() = default; /// Copy constructor - FCCSWGridModuleThetaMerged(const FCCSWGridModuleThetaMerged& e) = default; + FCCSWGridModuleThetaMerged_k4geo(const FCCSWGridModuleThetaMerged_k4geo& e) = default; /// Copy Constructor from segmentation base object - FCCSWGridModuleThetaMerged(const Segmentation& e) : Handle(e) {} + FCCSWGridModuleThetaMerged_k4geo(const Segmentation& e) : Handle(e) {} /// Copy constructor from handle - FCCSWGridModuleThetaMerged(const Handle& e) : Handle(e) {} + FCCSWGridModuleThetaMerged_k4geo(const Handle& e) : Handle(e) {} /// Copy constructor from other polymorph/equivalent handle template - FCCSWGridModuleThetaMerged(const Handle& e) : Handle(e) {} + FCCSWGridModuleThetaMerged_k4geo(const Handle& e) : Handle(e) {} /// Assignment operator - FCCSWGridModuleThetaMerged& operator=(const FCCSWGridModuleThetaMerged& seg) = default; + FCCSWGridModuleThetaMerged_k4geo& operator=(const FCCSWGridModuleThetaMerged_k4geo& seg) = default; /// Equality operator - bool operator==(const FCCSWGridModuleThetaMerged& seg) const { return m_element == seg.m_element; } + bool operator==(const FCCSWGridModuleThetaMerged_k4geo& seg) const { return m_element == seg.m_element; } /// determine the position based on the cell ID inline Position position(const CellID& id) const { return Position(access()->implementation->position(id)); } diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h similarity index 86% rename from detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h rename to detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h index 2b3fccb7..796ba400 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h @@ -1,11 +1,11 @@ -#ifndef DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H -#define DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H +#ifndef DETSEGMENTATION_FCCSWGridModuleThetaMerged_k4geo_H +#define DETSEGMENTATION_FCCSWGridModuleThetaMerged_k4geo_H // FCCSW #include "detectorSegmentations/GridTheta_k4geo.h" #include "DD4hep/VolumeManager.h" -/** FCCSWGridModuleThetaMerged Detector/DetSegmentation/DetSegmentation/FCCSWGridModuleThetaMerged.h FCCSWGridModuleThetaMerged.h +/** FCCSWGridModuleThetaMerged_k4geo Detector/DetSegmentation/DetSegmentation/FCCSWGridModuleThetaMerged_k4geo.h FCCSWGridModuleThetaMerged_k4geo.h * * Segmentation in theta and module. * Based on GridTheta, merges modules and theta cells based on layer number @@ -14,15 +14,15 @@ namespace dd4hep { namespace DDSegmentation { -class FCCSWGridModuleThetaMerged : public GridTheta_k4geo { +class FCCSWGridModuleThetaMerged_k4geo : public GridTheta_k4geo { public: /// default constructor using an arbitrary type - FCCSWGridModuleThetaMerged(const std::string& aCellEncoding); + FCCSWGridModuleThetaMerged_k4geo(const std::string& aCellEncoding); /// Default constructor used by derived classes passing an existing decoder - FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder); + FCCSWGridModuleThetaMerged_k4geo(const BitFieldCoder* decoder); /// destructor - virtual ~FCCSWGridModuleThetaMerged() = default; + virtual ~FCCSWGridModuleThetaMerged_k4geo() = default; /// read n(modules) from detector metadata void GetNModulesFromGeom(); @@ -114,4 +114,4 @@ class FCCSWGridModuleThetaMerged : public GridTheta_k4geo { }; } } -#endif /* DETSEGMENTATION_FCCSWGRIDMODULETHETAMERGED_H */ +#endif /* DETSEGMENTATION_FCCSWGridModuleThetaMerged_k4geo_H */ diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle_k4geo.cpp similarity index 80% rename from detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp rename to detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle_k4geo.cpp index 0573a2b9..c70004a6 100644 --- a/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle.cpp +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMergedHandle_k4geo.cpp @@ -1,4 +1,4 @@ -#include "detectorSegmentations/FCCSWGridModuleThetaMergedHandle.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h" #include "DD4hep/detail/Handle.inl" -DD4HEP_INSTANTIATE_HANDLE_UNNAMED(dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged); +DD4HEP_INSTANTIATE_HANDLE_UNNAMED(dd4hep::DDSegmentation::FCCSWGridModuleThetaMerged_k4geo); diff --git a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp b/detectorSegmentations/src/FCCSWGridModuleThetaMerged_k4geo.cpp similarity index 85% rename from detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp rename to detectorSegmentations/src/FCCSWGridModuleThetaMerged_k4geo.cpp index 1f1763ad..ac54d86a 100644 --- a/detectorSegmentations/src/FCCSWGridModuleThetaMerged.cpp +++ b/detectorSegmentations/src/FCCSWGridModuleThetaMerged_k4geo.cpp @@ -1,4 +1,4 @@ -#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" +#include "detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h" #include #include "DD4hep/Detector.h" @@ -7,9 +7,9 @@ namespace dd4hep { namespace DDSegmentation { /// default constructor using an encoding string -FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEncoding) : GridTheta_k4geo(cellEncoding) { +FCCSWGridModuleThetaMerged_k4geo::FCCSWGridModuleThetaMerged_k4geo(const std::string& cellEncoding) : GridTheta_k4geo(cellEncoding) { // define type and description - _type = "FCCSWGridModuleThetaMerged"; + _type = "FCCSWGridModuleThetaMerged_k4geo"; _description = "Module-theta segmentation with per-layer merging along theta and/or module"; // register all necessary parameters (additional to those registered in GridTheta_k4geo) @@ -21,9 +21,9 @@ FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const std::string& cellEn GetNLayersFromGeom(); } -FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const BitFieldCoder* decoder) : GridTheta_k4geo(decoder) { +FCCSWGridModuleThetaMerged_k4geo::FCCSWGridModuleThetaMerged_k4geo(const BitFieldCoder* decoder) : GridTheta_k4geo(decoder) { // define type and description - _type = "FCCSWGridModuleThetaMerged"; + _type = "FCCSWGridModuleThetaMerged_k4geo"; _description = "Module-theta segmentation with per-layer merging along theta and/or module"; // register all necessary parameters (additional to those registered in GridTheta_k4geo) @@ -35,7 +35,7 @@ FCCSWGridModuleThetaMerged::FCCSWGridModuleThetaMerged(const BitFieldCoder* deco GetNLayersFromGeom(); } -void FCCSWGridModuleThetaMerged::GetNModulesFromGeom() { +void FCCSWGridModuleThetaMerged_k4geo::GetNModulesFromGeom() { dd4hep::Detector* dd4hepgeo = &(dd4hep::Detector::getInstance()); try { m_nModules = dd4hepgeo->constant("ECalBarrelNumPlanes"); @@ -47,7 +47,7 @@ void FCCSWGridModuleThetaMerged::GetNModulesFromGeom() { std::cout << "Number of modules read from detector metadata and used in readout class: " << m_nModules << std::endl; } -void FCCSWGridModuleThetaMerged::GetNLayersFromGeom() { +void FCCSWGridModuleThetaMerged_k4geo::GetNLayersFromGeom() { dd4hep::Detector* dd4hepgeo = &(dd4hep::Detector::getInstance()); try { m_nLayers = dd4hepgeo->constant("ECalBarrelNumLayers"); @@ -60,7 +60,7 @@ void FCCSWGridModuleThetaMerged::GetNLayersFromGeom() { } /// determine the local position based on the cell ID -Vector3D FCCSWGridModuleThetaMerged::position(const CellID& cID) const { +Vector3D FCCSWGridModuleThetaMerged_k4geo::position(const CellID& cID) const { // debug // std::cout << "cellID: " << cID << std::endl; @@ -70,7 +70,7 @@ Vector3D FCCSWGridModuleThetaMerged::position(const CellID& cID) const { } /// determine the cell ID based on the global position -CellID FCCSWGridModuleThetaMerged::cellID(const Vector3D& /* localPosition */, const Vector3D& globalPosition, +CellID FCCSWGridModuleThetaMerged_k4geo::cellID(const Vector3D& /* localPosition */, const Vector3D& globalPosition, const VolumeID& vID) const { CellID cID = vID; @@ -109,7 +109,7 @@ CellID FCCSWGridModuleThetaMerged::cellID(const Vector3D& /* localPosition */, c /// merged ones - which will be then added on top of /// the phi of the volume containing the first cell /// by the positioning tool -double FCCSWGridModuleThetaMerged::phi(const CellID& cID) const { +double FCCSWGridModuleThetaMerged_k4geo::phi(const CellID& cID) const { // retrieve layer int layer = _decoder->get(cID, m_layerID); @@ -128,7 +128,7 @@ double FCCSWGridModuleThetaMerged::phi(const CellID& cID) const { /// determine the polar angle based on the cell ID and the /// number of merged theta cells -double FCCSWGridModuleThetaMerged::theta(const CellID& cID) const { +double FCCSWGridModuleThetaMerged_k4geo::theta(const CellID& cID) const { // retrieve layer int layer = _decoder->get(cID, m_layerID); diff --git a/detectorSegmentations/src/plugins/SegmentationFactories.cpp b/detectorSegmentations/src/plugins/SegmentationFactories.cpp index cc308c42..3e65e7ba 100644 --- a/detectorSegmentations/src/plugins/SegmentationFactories.cpp +++ b/detectorSegmentations/src/plugins/SegmentationFactories.cpp @@ -17,8 +17,8 @@ DECLARE_SEGMENTATION(GridTheta_k4geo, create_segmentation) -#include "detectorSegmentations/FCCSWGridModuleThetaMerged.h" -DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged, create_segmentation) +#include "detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h" +DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged, create_segmentation) #include "detectorSegmentations/FCCSWGridPhiEta_k4geo.h" DECLARE_SEGMENTATION(FCCSWGridPhiEta_k4geo, create_segmentation) From d7663aafd5369821516b6a7dc310629c792a5f6c Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Mon, 30 Oct 2023 18:15:33 +0100 Subject: [PATCH 05/15] Merging Giovanni changes of ECal barrel inclined constructor into the refactored new constructor ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp --- .../FCCee_ECalBarrel_thetamodulemerged.xml | 2 +- .../calorimeter/ECalBarrelInclined_geo.cpp | 614 ------------------ ...leLiquid_InclinedTrapezoids_o1_v01_geo.cpp | 43 ++ 3 files changed, 44 insertions(+), 615 deletions(-) delete mode 100644 detector/calorimeter/ECalBarrelInclined_geo.cpp diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml index a454914d..c8655ce7 100644 --- a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml @@ -96,7 +96,7 @@ - + diff --git a/detector/calorimeter/ECalBarrelInclined_geo.cpp b/detector/calorimeter/ECalBarrelInclined_geo.cpp deleted file mode 100644 index 51178159..00000000 --- a/detector/calorimeter/ECalBarrelInclined_geo.cpp +++ /dev/null @@ -1,614 +0,0 @@ -#include "DD4hep/DetFactoryHelper.h" -#include "DD4hep/Handle.h" -#include "XML/Utilities.h" - -#include - -// todo: remove gaudi logging and properly capture output -#define endmsg std::endl -#define lLog std::cout -namespace MSG { -const std::string ERROR = " Error: "; -const std::string DEBUG = " Debug: "; -const std::string INFO = " Info: "; -} - -namespace det { -static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, - dd4hep::xml::Handle_t aXmlElement, - dd4hep::SensitiveDetector aSensDet) { - - dd4hep::xml::DetElement xmlDetElem = aXmlElement; - std::string nameDet = xmlDetElem.nameStr(); - dd4hep::xml::Dimension dim(xmlDetElem.dimensions()); - dd4hep::DetElement caloDetElem(nameDet, xmlDetElem.id()); - - // Create air envelope for the whole barrel - dd4hep::Volume envelopeVol(nameDet + "_vol", dd4hep::Tube(dim.rmin(), dim.rmax(), dim.dz()), - aLcdd.material("Air")); - envelopeVol.setVisAttributes(aLcdd, dim.visStr()); - - // Retrieve cryostat data - dd4hep::xml::DetElement cryostat = aXmlElement.child(_Unicode(cryostat)); - dd4hep::xml::Dimension cryoDim(cryostat.dimensions()); - double cryoThicknessFront = cryoDim.rmin2() - cryoDim.rmin1(); - dd4hep::xml::DetElement cryoFront = cryostat.child(_Unicode(front)); - dd4hep::xml::DetElement cryoBack = cryostat.child(_Unicode(back)); - dd4hep::xml::DetElement cryoSide = cryostat.child(_Unicode(side)); - bool cryoFrontSensitive = cryoFront.isSensitive(); - bool cryoBackSensitive = cryoBack.isSensitive(); - bool cryoSideSensitive = cryoSide.isSensitive(); - - // Retrieve active and passive material data - dd4hep::xml::DetElement calo = aXmlElement.child(_Unicode(calorimeter)); - dd4hep::xml::Dimension caloDim(calo.dimensions()); - dd4hep::xml::DetElement active = calo.child(_Unicode(active)); - std::string activeMaterial = active.materialStr(); - double activeThickness = active.thickness(); - - dd4hep::xml::DetElement overlap = active.child(_Unicode(overlap)); - double activePassiveOverlap = overlap.offset(); - if (activePassiveOverlap < 0 || activePassiveOverlap > 0.5) { - // todo: ServiceHandle incidentSvc("IncidentSvc", "ECalConstruction"); - lLog << MSG::ERROR << "Overlap between active and passive cannot be more than half of passive plane!" << endmsg; - //todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); - } - dd4hep::xml::DetElement layers = calo.child(_Unicode(layers)); - uint numLayers = 0; - std::vector layerHeight; - double layersTotalHeight = 0; - for (dd4hep::xml::Collection_t layer_coll(layers, _Unicode(layer)); layer_coll; ++layer_coll) { - dd4hep::xml::Component layer = layer_coll; - numLayers += layer.repeat(); - for (int iLay = 0; iLay < layer.repeat(); iLay++) { - layerHeight.push_back(layer.thickness()); - } - layersTotalHeight += layer.repeat() * layer.thickness(); - } - lLog << MSG::DEBUG << "Number of layers: " << numLayers << " total thickness " << layersTotalHeight << endmsg; - // The following code checks if the xml geometry file contains a constant defining - // the number of layers the barrel. In that case, it makes the program abort - // if the number of planes in the xml is different from the one calculated from - // the geometry. This is because the number of layers is needed - // in other parts of the code (the readout for the FCC-ee ECAL with - // inclined modules). - int nLayers = -1; - try { - nLayers = aLcdd.constant("ECalBarrelNumLayers"); - } - catch(...) { - ; - } - if (nLayers > 0 && nLayers != numLayers) { - lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; - // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); - // make the code crash (incidentSvc does not work) - assert(nLayers == numLayers); - } - - dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); - std::string readoutMaterial = readout.materialStr(); - double readoutThickness = readout.thickness(); - - dd4hep::xml::DetElement passive = calo.child(_Unicode(passive)); - dd4hep::xml::DetElement passiveInner = passive.child(_Unicode(inner)); - dd4hep::xml::DetElement passiveInnerMax = passive.child(_Unicode(innerMax)); - dd4hep::xml::DetElement passiveOuter = passive.child(_Unicode(outer)); - dd4hep::xml::DetElement passiveGlue = passive.child(_Unicode(glue)); - std::string passiveInnerMaterial = passiveInner.materialStr(); - std::string passiveOuterMaterial = passiveOuter.materialStr(); - std::string passiveGlueMaterial = passiveGlue.materialStr(); - double passiveInnerThicknessMin = passiveInner.thickness(); - double passiveInnerThicknessMax = passiveInnerMax.thickness(); - double passiveOuterThickness = passiveOuter.thickness(); - double passiveGlueThickness = passiveGlue.thickness(); - double passiveThickness = passiveInnerThicknessMin + passiveOuterThickness + passiveGlueThickness; - double angle = passive.rotation().angle(); - - double bathRmin = caloDim.rmin(); // - margin for inclination - double bathRmax = caloDim.rmax(); // + margin for inclination - dd4hep::Tube bathOuterShape(bathRmin, bathRmax, caloDim.dz()); // make it 4 volumes + 5th for detector envelope - dd4hep::Tube bathAndServicesOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), caloDim.dz()); // make it 4 volumes + 5th for detector envelope - if (cryoThicknessFront > 0) { - // 1. Create cryostat - dd4hep::Tube cryoFrontShape(cryoDim.rmin1(), cryoDim.rmin2(), cryoDim.dz()); - dd4hep::Tube cryoBackShape(cryoDim.rmax1(), cryoDim.rmax2(), cryoDim.dz()); - dd4hep::Tube cryoSideOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), cryoDim.dz()); - dd4hep::SubtractionSolid cryoSideShape(cryoSideOuterShape, bathAndServicesOuterShape); - lLog << MSG::INFO << "ECAL cryostat: front: rmin (cm) = " << cryoDim.rmin1() << " rmax (cm) = " << cryoDim.rmin2() << " dz (cm) = " << cryoDim.dz() << endmsg; - lLog << MSG::INFO << "ECAL cryostat: back: rmin (cm) = " << cryoDim.rmax1() << " rmax (cm) = " << cryoDim.rmax2() << " dz (cm) = " << cryoDim.dz() << endmsg; - lLog << MSG::INFO << "ECAL cryostat: side: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << cryoDim.dz() - caloDim.dz() << endmsg; - dd4hep::Volume cryoFrontVol(cryostat.nameStr()+"_front", cryoFrontShape, aLcdd.material(cryostat.materialStr())); - dd4hep::Volume cryoBackVol(cryostat.nameStr()+"_back", cryoBackShape, aLcdd.material(cryostat.materialStr())); - dd4hep::Volume cryoSideVol(cryostat.nameStr()+"_side", cryoSideShape, aLcdd.material(cryostat.materialStr())); - dd4hep::PlacedVolume cryoFrontPhysVol = envelopeVol.placeVolume(cryoFrontVol); - dd4hep::PlacedVolume cryoBackPhysVol = envelopeVol.placeVolume(cryoBackVol); - dd4hep::PlacedVolume cryoSidePhysVol = envelopeVol.placeVolume(cryoSideVol); - if (cryoFrontSensitive) { - cryoFrontVol.setSensitiveDetector(aSensDet); - cryoFrontPhysVol.addPhysVolID("cryo", 1); - cryoFrontPhysVol.addPhysVolID("type", 1); - lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; - } - if (cryoBackSensitive) { - cryoBackVol.setSensitiveDetector(aSensDet); - cryoBackPhysVol.addPhysVolID("cryo", 1); - cryoBackPhysVol.addPhysVolID("type", 2); - lLog << MSG::INFO << "Cryostat back volume set as sensitive" << endmsg; - } - if (cryoSideSensitive) { - cryoSideVol.setSensitiveDetector(aSensDet); - cryoSidePhysVol.addPhysVolID("cryo", 1); - cryoSidePhysVol.addPhysVolID("type", 3); - lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; - } - dd4hep::DetElement cryoFrontDetElem(caloDetElem, "cryo_front", 0); - cryoFrontDetElem.setPlacement(cryoFrontPhysVol); - dd4hep::DetElement cryoBackDetElem(caloDetElem, "cryo_back", 0); - cryoBackDetElem.setPlacement(cryoBackPhysVol); - dd4hep::DetElement cryoSideDetElem(caloDetElem, "cryo_side", 0); - cryoSideDetElem.setPlacement(cryoSidePhysVol); - // 1.2. Create place-holder for services - dd4hep::Tube servicesFrontShape(cryoDim.rmin2(), bathRmin, caloDim.dz()); - dd4hep::Tube servicesBackShape(bathRmax, cryoDim.rmax1(), caloDim.dz()); - lLog << MSG::INFO << "ECAL services: front: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << bathRmin << " dz (cm) = " << caloDim.dz() << endmsg; - lLog << MSG::INFO << "ECAL services: back: rmin (cm) = " << bathRmax << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << caloDim.dz() << endmsg; - dd4hep::Volume servicesFrontVol("services_front", servicesFrontShape, aLcdd.material(activeMaterial)); - dd4hep::Volume servicesBackVol("services_back", servicesBackShape, aLcdd.material(activeMaterial)); - dd4hep::PlacedVolume servicesFrontPhysVol = envelopeVol.placeVolume(servicesFrontVol); - dd4hep::PlacedVolume servicesBackPhysVol = envelopeVol.placeVolume(servicesBackVol); - if (cryoFrontSensitive) { - servicesFrontVol.setSensitiveDetector(aSensDet); - servicesFrontPhysVol.addPhysVolID("cryo", 1); - servicesFrontPhysVol.addPhysVolID("type", 4); - lLog << MSG::INFO << "Services front volume set as sensitive" << endmsg; - } - if (cryoBackSensitive) { - servicesBackVol.setSensitiveDetector(aSensDet); - servicesBackPhysVol.addPhysVolID("cryo", 1); - servicesBackPhysVol.addPhysVolID("type", 5); - lLog << MSG::INFO << "Services back volume set as sensitive" << endmsg; - } - dd4hep::DetElement servicesFrontDetElem(caloDetElem, "services_front", 0); - servicesFrontDetElem.setPlacement(servicesFrontPhysVol); - dd4hep::DetElement servicesBackDetElem(caloDetElem, "services_back", 0); - servicesBackDetElem.setPlacement(servicesBackPhysVol); - } - // 2. Create bath that is inside the cryostat and surrounds the detector - // Bath is filled with active material -> but not sensitive - dd4hep::Volume bathVol(activeMaterial + "_bath", bathOuterShape, aLcdd.material(activeMaterial)); - lLog << MSG::INFO << "ECAL bath: material = " << activeMaterial << " rmin (cm) = " << bathRmin - << " rmax (cm) = " << bathRmax << " thickness in front of ECal (cm) = " << caloDim.rmin() - cryoDim.rmin2() - << " thickness behind ECal (cm) = " << cryoDim.rmax1() - caloDim.rmax() << endmsg; - - // 3. Create the calorimeter by placing the passive material, trapezoid active layers, readout and again trapezoid - // active layers in the bath. - // sensitive detector for the layers - dd4hep::SensitiveDetector sd = aSensDet; - dd4hep::xml::Dimension sdType = xmlDetElem.child(_U(sensitive)); - sd.setType(sdType.typeStr()); - lLog << MSG::INFO << "ECAL calorimeter volume rmin (cm) = " << caloDim.rmin() << " rmax (cm) = " << caloDim.rmax() - << endmsg; - - // 3.a. Create the passive planes, readout in between of 2 passive planes and the remaining space fill with active - // material - ////////////////////////////// - // PASSIVE PLANES - ////////////////////////////// - lLog << MSG::INFO << "passive inner material = " << passiveInnerMaterial << "\n" - << " and outer material = " << passiveOuterMaterial << "\n" - << " thickness of inner part at inner radius (cm) = " << passiveInnerThicknessMin << "\n" - << " thickness of inner part at outer radius (cm) = " << passiveInnerThicknessMax << "\n" - << " thickness of outer part (cm) = " << passiveOuterThickness << "\n" - << " thickness of total (cm) = " << passiveThickness << "\n" - << " rotation angle = " << angle << endmsg; - uint numPlanes = - round(M_PI / asin((passiveThickness + activeThickness + readoutThickness) / (2. * caloDim.rmin() * cos(angle)))); - - double dPhi = 2. * M_PI / numPlanes; - lLog << MSG::INFO << "number of passive plates = " << numPlanes << " azim. angle difference = " << dPhi << endmsg; - lLog << MSG::INFO << " distance at inner radius (cm) = " << 2. * M_PI * caloDim.rmin() / numPlanes << "\n" - << " distance at outer radius (cm) = " << 2. * M_PI * caloDim.rmax() / numPlanes << endmsg; - - // The following code checks if the xml geometry file contains a constant defining - // the number of planes in the barrel. In that case, it makes the program abort - // if the number of planes in the xml is different from the one calculated from - // the geometry. This is because the number of plane information (retrieved from the - // xml) is used in other parts of the code (the readout for the FCC-ee ECAL with - // inclined modules). In principle the code above should be refactored so that the number - // of planes is one of the inputs of the calculation and other geometrical parameters - // are adjusted accordingly. This is left for the future, and we use the workaround - // below to enforce for the time being that the number of planes is "correct" - int nModules = -1; - try { - nModules = aLcdd.constant("ECalBarrelNumPlanes"); - } - catch(...) { - ; - } - if (nModules > 0 && nModules != numPlanes) { - lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; - // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); - // make the code crash (incidentSvc does not work) - assert(nModules == numPlanes); - } - // Readout is in the middle between two passive planes - double offsetPassivePhi = caloDim.offset() + dPhi / 2.; - double offsetReadoutPhi = caloDim.offset() + 0; - lLog << MSG::INFO << "readout material = " << readoutMaterial << "\n" - << " thickness of readout planes (cm) = " << readoutThickness << "\n number of readout layers = " << numLayers - << endmsg; - double Rmin = caloDim.rmin(); - double Rmax = caloDim.rmax(); - double dR = Rmax - Rmin; - double planeLength = -Rmin * cos(angle) + sqrt(pow(Rmax, 2) - pow(Rmin * sin(angle), 2)); - lLog << MSG::INFO << "thickness of calorimeter (cm) = " << dR << "\n" - << " length of passive or readout planes (cm) = " << planeLength << endmsg; - - - // fill the thickness in the boundary of each layer - std::vector passiveInnerThicknessLayer(numLayers+1); - double runningHeight = 0; - for (uint iLay = 0; iLay < numLayers; iLay++) { - passiveInnerThicknessLayer[iLay] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * - (runningHeight) / (Rmax - Rmin); - runningHeight += layerHeight[iLay]; - } - passiveInnerThicknessLayer[numLayers] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * - (runningHeight) / (Rmax - Rmin); - - double passiveAngle = atan2((passiveInnerThicknessMax - passiveInnerThicknessMin) / 2., planeLength); - double cosPassiveAngle = cos(passiveAngle); - double rotatedOuterThickness = passiveOuterThickness / cosPassiveAngle; - double rotatedGlueThickness = passiveGlueThickness / cosPassiveAngle; - - // rescale layer thicknesses - double scaleLayerThickness = planeLength / layersTotalHeight; - layersTotalHeight = 0; - for (uint iLay = 0; iLay < numLayers; iLay++) { - layerHeight[iLay] *= scaleLayerThickness; - - layersTotalHeight += layerHeight[iLay]; - lLog << MSG::DEBUG << "Thickness of layer " << iLay << " : " << layerHeight[iLay] << endmsg; - } - double layerFirstOffset = -planeLength / 2. + layerHeight[0] / 2.; - - //dd4hep::Box passiveShape(passiveThickness / 2., caloDim.dz(), planeLength / 2.); - dd4hep::Trd1 passiveShape(passiveInnerThicknessMin / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., - passiveInnerThicknessMax / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., - caloDim.dz(), planeLength / 2.); - // inner layer is not in the first calo layer (to sample more uniformly in the layer where upstream correction is - // applied) - //dd4hep::Box passiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); - dd4hep::Trd1 passiveInnerShape(passiveInnerThicknessLayer[1] / 2., passiveInnerThicknessMax / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); - //dd4hep::Box passiveInnerShapeFirstLayer(passiveInnerThickness / 2., caloDim.dz(), layerHeight[0] / 2.); - dd4hep::Trd1 passiveInnerShapeFirstLayer(passiveInnerThicknessMin / 2., passiveInnerThicknessLayer[1] / 2., caloDim.dz(), layerHeight[0] / 2.); - dd4hep::Box passiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); - dd4hep::Box passiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); - // passive volume consists of inner part and two outer, joind by glue - dd4hep::Volume passiveVol("passive", passiveShape, aLcdd.material("Air")); - dd4hep::Volume passiveInnerVol(passiveInnerMaterial + "_passive", passiveInnerShape, - aLcdd.material(passiveInnerMaterial)); - dd4hep::Volume passiveInnerVolFirstLayer(activeMaterial + "_passive", passiveInnerShapeFirstLayer, - aLcdd.material(activeMaterial)); - dd4hep::Volume passiveOuterVol(passiveOuterMaterial + "_passive", passiveOuterShape, - aLcdd.material(passiveOuterMaterial)); - dd4hep::Volume passiveGlueVol(passiveGlueMaterial + "_passive", passiveGlueShape, - aLcdd.material(passiveGlueMaterial)); - - if (passiveInner.isSensitive()) { - lLog << MSG::DEBUG << "Passive inner volume set as sensitive" << endmsg; - // inner part starts at second layer - double layerOffset = layerFirstOffset + layerHeight[1] / 2.; - for (uint iLayer = 1; iLayer < numLayers; iLayer++) { - //dd4hep::Box layerPassiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); - dd4hep::Trd1 layerPassiveInnerShape(passiveInnerThicknessLayer[iLayer] / 2., passiveInnerThicknessLayer[iLayer+1] / 2., caloDim.dz(), layerHeight[iLayer] / 2.); - dd4hep::Volume layerPassiveInnerVol(passiveInnerMaterial, layerPassiveInnerShape, - aLcdd.material(passiveInnerMaterial)); - layerPassiveInnerVol.setSensitiveDetector(aSensDet); - dd4hep::PlacedVolume layerPassiveInnerPhysVol = - passiveInnerVol.placeVolume(layerPassiveInnerVol, dd4hep::Position(0, 0, layerOffset)); - layerPassiveInnerPhysVol.addPhysVolID("layer", iLayer); - dd4hep::DetElement layerPassiveInnerDetElem("layer", iLayer); - layerPassiveInnerDetElem.setPlacement(layerPassiveInnerPhysVol); - if (iLayer != numLayers - 1) { - layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; - } - } - } - if (passiveOuter.isSensitive()) { - lLog << MSG::DEBUG << "Passive outer volume set as sensitive" << endmsg; - double layerOffset = layerFirstOffset / cosPassiveAngle; - for (uint iLayer = 0; iLayer < numLayers; iLayer++) { - dd4hep::Box layerPassiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); - dd4hep::Volume layerPassiveOuterVol(passiveOuterMaterial, layerPassiveOuterShape, - aLcdd.material(passiveOuterMaterial)); - layerPassiveOuterVol.setSensitiveDetector(aSensDet); - dd4hep::PlacedVolume layerPassiveOuterPhysVol = - passiveOuterVol.placeVolume(layerPassiveOuterVol, dd4hep::Position(0, 0, layerOffset)); - layerPassiveOuterPhysVol.addPhysVolID("layer", iLayer); - dd4hep::DetElement layerPassiveOuterDetElem("layer", iLayer); - layerPassiveOuterDetElem.setPlacement(layerPassiveOuterPhysVol); - if (iLayer != numLayers - 1) { - layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; - } - } - } - if (passiveGlue.isSensitive()) { - lLog << MSG::DEBUG << "Passive glue volume set as sensitive" << endmsg; - double layerOffset = layerFirstOffset / cosPassiveAngle; - for (uint iLayer = 0; iLayer < numLayers; iLayer++) { - dd4hep::Box layerPassiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); - dd4hep::Volume layerPassiveGlueVol(passiveGlueMaterial, layerPassiveGlueShape, - aLcdd.material(passiveGlueMaterial)); - layerPassiveGlueVol.setSensitiveDetector(aSensDet); - dd4hep::PlacedVolume layerPassiveGluePhysVol = - passiveGlueVol.placeVolume(layerPassiveGlueVol, dd4hep::Position(0, 0, layerOffset)); - layerPassiveGluePhysVol.addPhysVolID("layer", iLayer); - dd4hep::DetElement layerPassiveGlueDetElem("layer", iLayer); - layerPassiveGlueDetElem.setPlacement(layerPassiveGluePhysVol); - if (iLayer != numLayers - 1) { - layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; - } - } - } - - dd4hep::PlacedVolume passiveInnerPhysVol = - passiveVol.placeVolume(passiveInnerVol, dd4hep::Position(0, 0, layerHeight[0] / 2.)); - dd4hep::PlacedVolume passiveInnerPhysVolFirstLayer = - passiveVol.placeVolume(passiveInnerVolFirstLayer, dd4hep::Position(0, 0, layerFirstOffset)); - dd4hep::PlacedVolume passiveOuterPhysVolBelow = passiveVol.placeVolume( - passiveOuterVol, - dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), - dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - - rotatedGlueThickness / 2. - rotatedOuterThickness / 4., 0, 0))); - dd4hep::PlacedVolume passiveOuterPhysVolAbove = passiveVol.placeVolume( - passiveOuterVol, - dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), - dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + - rotatedGlueThickness / 2. + rotatedOuterThickness / 4., 0, 0))); - dd4hep::PlacedVolume passiveGluePhysVolBelow = passiveVol.placeVolume( - passiveGlueVol, - dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), - dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - - rotatedGlueThickness / 4., 0, 0))); - dd4hep::PlacedVolume passiveGluePhysVolAbove = passiveVol.placeVolume( - passiveGlueVol, - dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), - dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + - rotatedGlueThickness / 4., 0, 0))); - passiveInnerPhysVol.addPhysVolID("subtype", 0); - passiveInnerPhysVolFirstLayer.addPhysVolID("subtype", 0); - passiveOuterPhysVolBelow.addPhysVolID("subtype", 1); - passiveOuterPhysVolAbove.addPhysVolID("subtype", 2); - passiveGluePhysVolBelow.addPhysVolID("subtype", 3); - passiveGluePhysVolAbove.addPhysVolID("subtype", 4); - if (passiveInner.isSensitive()) { - passiveInnerVolFirstLayer.setSensitiveDetector(aSensDet); - passiveInnerPhysVolFirstLayer.addPhysVolID("layer", 0); - dd4hep::DetElement passiveInnerDetElemFirstLayer("layer", 0); - passiveInnerDetElemFirstLayer.setPlacement(passiveInnerPhysVolFirstLayer); - } - - ////////////////////////////// - // READOUT PLANES - ////////////////////////////// - dd4hep::Box readoutShape(readoutThickness / 2., caloDim.dz(), planeLength / 2.); - dd4hep::Volume readoutVol(readoutMaterial, readoutShape, aLcdd.material(readoutMaterial)); - if (readout.isSensitive()) { - lLog << MSG::INFO << "Readout volume set as sensitive" << endmsg; - double layerOffset = layerFirstOffset; - for (uint iLayer = 0; iLayer < numLayers; iLayer++) { - dd4hep::Box layerReadoutShape(readoutThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); - dd4hep::Volume layerReadoutVol(readoutMaterial, layerReadoutShape, aLcdd.material(readoutMaterial)); - layerReadoutVol.setSensitiveDetector(aSensDet); - dd4hep::PlacedVolume layerReadoutPhysVol = - readoutVol.placeVolume(layerReadoutVol, dd4hep::Position(0, 0, layerOffset)); - layerReadoutPhysVol.addPhysVolID("layer", iLayer); - dd4hep::DetElement layerReadoutDetElem("layer", iLayer); - layerReadoutDetElem.setPlacement(layerReadoutPhysVol); - if (iLayer != numLayers - 1) { - layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; - } - } - } - - ////////////////////////////// - // ACTIVE - ////////////////////////////// - // thickness of active layers at inner radius and outer ( = distance between passive plane and readout plane) - // at inner radius: distance projected at plane perpendicular to readout plane - double activeInThickness = Rmin * sin(dPhi / 2.) * cos(angle); - activeInThickness -= passiveThickness * (0.5 - activePassiveOverlap); - // at outer radius: distance projected at plane perpendicular to readout plane - double activeOutThickness = (Rmin + planeLength) * sin(dPhi / 2.) * cos(angle); - // make correction for outer readius caused by inclination angle - // first calculate intersection of readout plane and plane parallel to shifted passive plane - double xIntersect = (Rmin * (tan(angle) - cos(dPhi / 2.) * tan(angle + dPhi / 2.)) - planeLength * sin(dPhi / 2.)) / - (tan(angle) - tan(angle + dPhi / 2.)); - double yIntersect = tan(angle) * xIntersect + Rmin * (sin(dPhi / 2.) - tan(angle)) + planeLength * sin(dPhi / 2.); - // distance from inner radius to intersection - double correction = - planeLength - sqrt(pow(xIntersect - Rmin * cos(dPhi / 2), 2) + pow(yIntersect - Rmin * sin(dPhi / 2), 2)); - // correction to the active thickness - activeOutThickness += 2. * correction * sin(dPhi / 4.); - activeOutThickness -= passiveThickness * (0.5 - activePassiveOverlap); - // print the active layer dimensions - double activeInThicknessAfterSubtraction = - 2. * activeInThickness - readoutThickness - 2. * activePassiveOverlap * passiveThickness; - double activeOutThicknessAfterSubtraction = - 2. * activeOutThickness - readoutThickness - 2. * activePassiveOverlap * - (passiveThickness + passiveInnerThicknessMax - passiveInnerThicknessMin); // correct thickness for trapezoid - lLog << MSG::INFO << "active material = " << activeMaterial - << " active layers thickness at inner radius (cm) = " << activeInThicknessAfterSubtraction - << " thickness at outer radious (cm) = " << activeOutThicknessAfterSubtraction << " making " - << (activeOutThicknessAfterSubtraction - activeInThicknessAfterSubtraction) * 100 / - activeInThicknessAfterSubtraction - << " % increase." << endmsg; - lLog << MSG::INFO - << "active passive initial overlap (before subtraction) (cm) = " << passiveThickness * activePassiveOverlap - << " = " << activePassiveOverlap * 100 << " %" << endmsg; - - // creating shape for rows of layers (active material between two passive planes, with readout in the middle) - // first define area between two passive planes, area can reach up to the symmetry axis of passive plane - dd4hep::Trd1 activeOuterShape(activeInThickness, activeOutThickness, caloDim.dz(), planeLength / 2.); - // subtract readout shape from the middle - dd4hep::SubtractionSolid activeShapeNoReadout(activeOuterShape, readoutShape); - - // make calculation for active plane that is inclined with 0 deg (= offset + angle) - double Cx = Rmin * cos(-angle) + planeLength / 2.; - double Cy = Rmin * sin(-angle); - double Ax = Rmin * cos(-angle + dPhi / 2.) + planeLength / 2. * cos(dPhi / 2.); - double Ay = Rmin * sin(-angle + dPhi / 2.) + planeLength / 2. * sin(dPhi / 2.); - double CAx = fabs(Ax - Cx); - double CAy = fabs(Ay - Cy); - double zprim, xprim; - zprim = CAx; - xprim = CAy; - - double Bx = Rmin * cos(-angle - dPhi / 2.) + planeLength / 2. * cos(-dPhi / 2.); - double By = Rmin * sin(-angle - dPhi / 2.) + planeLength / 2. * sin(-dPhi / 2.); - double CBx = fabs(Bx - Cx); - double CBy = fabs(By - Cy); - double zprimB, xprimB; - zprimB = CBx; - xprimB = CBy; - - // subtract passive volume above - dd4hep::SubtractionSolid activeShapeNoPassiveAbove( - activeShapeNoReadout, passiveShape, - dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), - dd4hep::Position(-fabs(xprim), 0, fabs(zprim)))); - // subtract passive volume below - dd4hep::SubtractionSolid activeShape( - activeShapeNoPassiveAbove, passiveShape, - dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), - dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB)))); - dd4hep::Volume activeVol("active", activeShape, aLcdd.material("Air")); - - std::vector layerPhysVols; - // place layers within active volume - std::vector layerInThickness; - std::vector layerOutThickness; - double layerIncreasePerUnitThickness = (activeOutThickness - activeInThickness) / layersTotalHeight; - for (uint iLay = 0; iLay < numLayers; iLay++) { - if (iLay == 0) { - layerInThickness.push_back(activeInThickness); - } else { - layerInThickness.push_back(layerOutThickness[iLay - 1]); - } - layerOutThickness.push_back(layerInThickness[iLay] + layerIncreasePerUnitThickness * layerHeight[iLay]); - } - double layerOffset = layerFirstOffset; - for (uint iLayer = 0; iLayer < numLayers; iLayer++) { - dd4hep::Trd1 layerOuterShape(layerInThickness[iLayer], layerOutThickness[iLayer], caloDim.dz(), layerHeight[iLayer] / 2.); - dd4hep::SubtractionSolid layerShapeNoReadout(layerOuterShape, readoutShape); - dd4hep::SubtractionSolid layerShapeNoPassiveAbove( - layerShapeNoReadout, passiveShape, - dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), - dd4hep::Position(-fabs(xprim), 0, fabs(zprim) - layerOffset))); - // subtract passive volume below - dd4hep::SubtractionSolid layerShape( - layerShapeNoPassiveAbove, passiveShape, - dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), - dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB) - layerOffset))); - dd4hep::Volume layerVol("layer", layerShape, aLcdd.material(activeMaterial)); - layerVol.setSensitiveDetector(aSensDet); - layerPhysVols.push_back(activeVol.placeVolume(layerVol, dd4hep::Position(0, 0, layerOffset))); - layerPhysVols.back().addPhysVolID("layer", iLayer); - if (iLayer != numLayers - 1) { - layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; - } - } - - dd4hep::DetElement bathDetElem(caloDetElem, "bath", 1); - std::vector activePhysVols; - // Next place elements: passive planes, readout planes and rows of layers - for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { - // first calculate positions of passive and readout planes - // PASSIVE - // calculate centre position of the plane without plane rotation - double phi = offsetPassivePhi + iPlane * dPhi; - double xRadial = (Rmin + planeLength / 2.) * cos(phi); - double yRadial = (Rmin + planeLength / 2.) * sin(phi); - // calculate position of the beginning of plane - double xRmin = Rmin * cos(phi); - double yRmin = Rmin * sin(phi); - // rotate centre by angle wrt beginning of plane - double xRotated = xRmin + (xRadial - xRmin) * cos(angle) - (yRadial - yRmin) * sin(angle); - double yRotated = yRmin + (xRadial - xRmin) * sin(angle) + (yRadial - yRmin) * cos(angle); - dd4hep::Transform3D transform(dd4hep::RotationX(-M_PI / 2.) // to get in XY plane - * - dd4hep::RotationY(M_PI / 2. // to get pointed towards centre - - - phi - angle), - dd4hep::Position(xRotated, yRotated, 0)); - dd4hep::PlacedVolume passivePhysVol = bathVol.placeVolume(passiveVol, transform); - passivePhysVol.addPhysVolID("module", iPlane); - passivePhysVol.addPhysVolID("type", 1); // 0 = active, 1 = passive, 2 = readout - dd4hep::DetElement passiveDetElem(bathDetElem, "passive" + std::to_string(iPlane), iPlane); - passiveDetElem.setPlacement(passivePhysVol); - - // READOUT - // calculate centre position of the plane without plane rotation - double phiRead = offsetReadoutPhi + iPlane * dPhi; - double xRadialRead = (Rmin + planeLength / 2.) * cos(phiRead); - double yRadialRead = (Rmin + planeLength / 2.) * sin(phiRead); - // calculate position of the beginning of plane - double xRminRead = Rmin * cos(phiRead); - double yRminRead = Rmin * sin(phiRead); - // rotate centre by angle wrt beginning of plane - double xRotatedRead = xRminRead + (xRadialRead - xRminRead) * cos(angle) - (yRadialRead - yRminRead) * sin(angle); - double yRotatedRead = yRminRead + (xRadialRead - xRminRead) * sin(angle) + (yRadialRead - yRminRead) * cos(angle); - dd4hep::Transform3D transformRead( - dd4hep::RotationX(-M_PI / 2.) // to get in XY plane - * - dd4hep::RotationY(M_PI / 2. // to get pointed towards centre - - - phiRead - angle), - dd4hep::Position(xRotatedRead, yRotatedRead, 0)); - dd4hep::PlacedVolume readoutPhysVol = bathVol.placeVolume(readoutVol, transformRead); - readoutPhysVol.addPhysVolID("module", iPlane); - readoutPhysVol.addPhysVolID("type", 2); // 0 = active, 1 = passive, 2 = readout - dd4hep::DetElement readoutDetElem(bathDetElem, "readout" + std::to_string(iPlane), iPlane); - readoutDetElem.setPlacement(readoutPhysVol); - - // ACTIVE - dd4hep::Rotation3D rotationActive(dd4hep::RotationX(-M_PI / 2) * - dd4hep::RotationY(M_PI / 2 - phiRead - angle)); - activePhysVols.push_back(bathVol.placeVolume( - activeVol, - dd4hep::Transform3D(rotationActive, dd4hep::Position(xRotatedRead, yRotatedRead, 0)))); - activePhysVols.back().addPhysVolID("module", iPlane); - activePhysVols.back().addPhysVolID("type", 0); // 0 = active, 1 = passive, 2 = readout - } - dd4hep::PlacedVolume bathPhysVol = envelopeVol.placeVolume(bathVol); - bathDetElem.setPlacement(bathPhysVol); - for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { - dd4hep::DetElement activeDetElem(bathDetElem, "active" + std::to_string(iPlane), iPlane); - activeDetElem.setPlacement(activePhysVols[iPlane]); - for (uint iLayer = 0; iLayer < numLayers; iLayer++) { - dd4hep::DetElement layerDetElem(activeDetElem, "layer" + std::to_string(iLayer), iLayer); - layerDetElem.setPlacement(layerPhysVols[iLayer]); - } - } - - // Place the envelope - dd4hep::Volume motherVol = aLcdd.pickMotherVolume(caloDetElem); - dd4hep::PlacedVolume envelopePhysVol = motherVol.placeVolume(envelopeVol); - envelopePhysVol.addPhysVolID("system", xmlDetElem.id()); - caloDetElem.setPlacement(envelopePhysVol); - - // Create caloData object - auto caloData = new dd4hep::rec::LayeredCalorimeterData; - caloData->layoutType = dd4hep::rec::LayeredCalorimeterData::BarrelLayout; - caloDetElem.addExtension(caloData); - - // Set type flags - dd4hep::xml::setDetectorTypeFlag(xmlDetElem, caloDetElem); - - return caloDetElem; -} -} // namespace det - -DECLARE_DETELEMENT(EmCaloBarrelInclined, det::createECalBarrelInclined) diff --git a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp index 44e43f27..770b60e4 100644 --- a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp +++ b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp @@ -68,6 +68,25 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, layersTotalHeight += layer.repeat() * layer.thickness(); } lLog << MSG::DEBUG << "Number of layers: " << numLayers << " total thickness " << layersTotalHeight << endmsg; + // The following code checks if the xml geometry file contains a constant defining + // the number of layers the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of layers is needed + // in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). + int nLayers = -1; + try { + nLayers = aLcdd.constant("ECalBarrelNumLayers"); + } + catch(...) { + ; + } + if (nLayers > 0 && nLayers != int(numLayers)) { + lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + assert(nLayers == int(numLayers)); + } dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); std::string readoutMaterial = readout.materialStr(); @@ -187,10 +206,34 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, << " rotation angle = " << angle << endmsg; uint numPlanes = round(M_PI / asin((passiveThickness + activeThickness + readoutThickness) / (2. * caloDim.rmin() * cos(angle)))); + double dPhi = 2. * M_PI / numPlanes; lLog << MSG::INFO << "number of passive plates = " << numPlanes << " azim. angle difference = " << dPhi << endmsg; lLog << MSG::INFO << " distance at inner radius (cm) = " << 2. * M_PI * caloDim.rmin() / numPlanes << "\n" << " distance at outer radius (cm) = " << 2. * M_PI * caloDim.rmax() / numPlanes << endmsg; + + // The following code checks if the xml geometry file contains a constant defining + // the number of planes in the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of plane information (retrieved from the + // xml) is used in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). In principle the code above should be refactored so that the number + // of planes is one of the inputs of the calculation and other geometrical parameters + // are adjusted accordingly. This is left for the future, and we use the workaround + // below to enforce for the time being that the number of planes is "correct" + int nModules = -1; + try { + nModules = aLcdd.constant("ECalBarrelNumPlanes"); + } + catch(...) { + ; + } + if (nModules > 0 && nModules != int(numPlanes)) { + lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + assert(nModules == int(numPlanes)); + } // Readout is in the middle between two passive planes double offsetPassivePhi = caloDim.offset() + dPhi / 2.; double offsetReadoutPhi = caloDim.offset() + 0; From 84b54bee2cf475ced774d529901c67bf0ca1e2ec Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Mon, 30 Oct 2023 18:15:42 +0100 Subject: [PATCH 06/15] Typo in comment fixed --- detectorSegmentations/src/GridEta_k4geo.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/detectorSegmentations/src/GridEta_k4geo.cpp b/detectorSegmentations/src/GridEta_k4geo.cpp index 5c387e5e..fada906e 100644 --- a/detectorSegmentations/src/GridEta_k4geo.cpp +++ b/detectorSegmentations/src/GridEta_k4geo.cpp @@ -46,7 +46,7 @@ CellID GridEta_k4geo::cellID(const Vector3D& /* localPosition */, const Vector3D // return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); //} -/// determine the seudorapidity based on the cell ID +/// determine the pseudorapidity based on the cell ID double GridEta_k4geo::eta(const CellID& cID) const { CellID etaValue = _decoder->get(cID, m_etaID); return binToPosition(etaValue, m_gridSizeEta, m_offsetEta); From c87ad44e0c0300c42cc5f541f25ddf7c647acf53 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Mon, 30 Oct 2023 18:36:47 +0100 Subject: [PATCH 07/15] assert replaced by exception --- ...ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp index 770b60e4..60d7e85b 100644 --- a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp +++ b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp @@ -85,7 +85,8 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); // make the code crash (incidentSvc does not work) - assert(nLayers == int(numLayers)); + // Andre, Alvaro, assert replaced by exception + throw std::runtime_error("Incorrect number of layers (ECalBarrelNumLayers) in xml file!"); } dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); @@ -232,7 +233,8 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); // make the code crash (incidentSvc does not work) - assert(nModules == int(numPlanes)); + // Andre, Alvaro, assert replaced by exception + throw std::runtime_error("Incorrect number of planes (ECalBarrelNumPlanes) in xml file!"); } // Readout is in the middle between two passive planes double offsetPassivePhi = caloDim.offset() + dPhi / 2.; From d787908c24aa55cd03ccf911e8f78bc074c9c37a Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Tue, 31 Oct 2023 13:25:12 +0100 Subject: [PATCH 08/15] Function not implemented throws exception https://github.com/key4hep/k4geo/pull/293#discussion_r1377229794 --- .../FCCSWGridModuleThetaMergedHandle_k4geo.h | 2 ++ 1 file changed, 2 insertions(+) diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h index 43d17740..f30ea245 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h @@ -116,6 +116,8 @@ class FCCSWGridModuleThetaMerged_k4geo : public FCCSWGridModuleThetaMergedHandle inline std::vector cellDimensions(const CellID& /*id*/) const { //return {access()->implementation->gridSizePhi(), access()->implementation->gridSizeTheta()}; // not implemented + throw std::runtime_error( Form("Function %s not implemented", __PRETTY_FUNCTION__) ); + return {0., 0.}; } }; From 7c299efa0f19b50beab1777ee5692b623bc80480 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Wed, 1 Nov 2023 10:24:31 +0100 Subject: [PATCH 09/15] a phrase added in the docstring that this is not implemented --- .../FCCSWGridModuleThetaMergedHandle_k4geo.h | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h index f30ea245..c0b79ea8 100644 --- a/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h +++ b/detectorSegmentations/include/detectorSegmentations/FCCSWGridModuleThetaMergedHandle_k4geo.h @@ -105,13 +105,14 @@ class FCCSWGridModuleThetaMerged_k4geo : public FCCSWGridModuleThetaMergedHandle /** \brief Returns a std::vector of the cellDimensions of the given cell ID - in natural order of dimensions (nModules, dTheta) + in natural order of dimensions (nModules, dTheta). Not implemented yet. Returns a std::vector of the cellDimensions of the given cell ID \param cellID \return std::vector size 2: -# size in module -# size in theta + Not implemented yet. */ inline std::vector cellDimensions(const CellID& /*id*/) const { //return {access()->implementation->gridSizePhi(), access()->implementation->gridSizeTheta()}; From 50e1c0887e382b7695f26998a04e3312b04fb90a Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Wed, 1 Nov 2023 10:26:01 +0100 Subject: [PATCH 10/15] rm file that is not used --- .../ALLEGRO_o1_v02/FCCee_ECalBarrel.xml | 143 ------------------ 1 file changed, 143 deletions(-) delete mode 100644 FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml deleted file mode 100644 index b1e0e9f2..00000000 --- a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel.xml +++ /dev/null @@ -1,143 +0,0 @@ - - - - - - Settings for the inclined EM calorimeter. - The barrel is filled with liquid argon. Passive material includes lead in the middle and steal on the outside, glued together. - Passive plates are inclined by a certain angle from the radial direction. - In between of two passive plates there is a readout. - Space between the plate and readout is of trapezoidal shape and filled with liquid argon. - Definition of sizes, visualization settings, readout and longitudinal segmentation are specified. - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - system:4,cryo:1,type:3,subtype:3,layer:8,module:11,eta:9 - - - - - - - system:4,cryo:1,type:3,subtype:3,layer:8,eta:9,phi:10 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - From 2f530af7716958aeac40b6e53031898db82512c4 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Wed, 1 Nov 2023 15:48:05 +0100 Subject: [PATCH 11/15] new version of ECalBarrel_NobleLiquid_InclinedTrapezoids https://github.com/key4hep/k4geo/pull/293#discussion_r1378882726 --- ...leLiquid_InclinedTrapezoids_o1_v01_geo.cpp | 45 -- ...leLiquid_InclinedTrapezoids_o1_v02_geo.cpp | 618 ++++++++++++++++++ 2 files changed, 618 insertions(+), 45 deletions(-) create mode 100644 detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02_geo.cpp diff --git a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp index 60d7e85b..44e43f27 100644 --- a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp +++ b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp @@ -68,26 +68,6 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, layersTotalHeight += layer.repeat() * layer.thickness(); } lLog << MSG::DEBUG << "Number of layers: " << numLayers << " total thickness " << layersTotalHeight << endmsg; - // The following code checks if the xml geometry file contains a constant defining - // the number of layers the barrel. In that case, it makes the program abort - // if the number of planes in the xml is different from the one calculated from - // the geometry. This is because the number of layers is needed - // in other parts of the code (the readout for the FCC-ee ECAL with - // inclined modules). - int nLayers = -1; - try { - nLayers = aLcdd.constant("ECalBarrelNumLayers"); - } - catch(...) { - ; - } - if (nLayers > 0 && nLayers != int(numLayers)) { - lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; - // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); - // make the code crash (incidentSvc does not work) - // Andre, Alvaro, assert replaced by exception - throw std::runtime_error("Incorrect number of layers (ECalBarrelNumLayers) in xml file!"); - } dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); std::string readoutMaterial = readout.materialStr(); @@ -207,35 +187,10 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, << " rotation angle = " << angle << endmsg; uint numPlanes = round(M_PI / asin((passiveThickness + activeThickness + readoutThickness) / (2. * caloDim.rmin() * cos(angle)))); - double dPhi = 2. * M_PI / numPlanes; lLog << MSG::INFO << "number of passive plates = " << numPlanes << " azim. angle difference = " << dPhi << endmsg; lLog << MSG::INFO << " distance at inner radius (cm) = " << 2. * M_PI * caloDim.rmin() / numPlanes << "\n" << " distance at outer radius (cm) = " << 2. * M_PI * caloDim.rmax() / numPlanes << endmsg; - - // The following code checks if the xml geometry file contains a constant defining - // the number of planes in the barrel. In that case, it makes the program abort - // if the number of planes in the xml is different from the one calculated from - // the geometry. This is because the number of plane information (retrieved from the - // xml) is used in other parts of the code (the readout for the FCC-ee ECAL with - // inclined modules). In principle the code above should be refactored so that the number - // of planes is one of the inputs of the calculation and other geometrical parameters - // are adjusted accordingly. This is left for the future, and we use the workaround - // below to enforce for the time being that the number of planes is "correct" - int nModules = -1; - try { - nModules = aLcdd.constant("ECalBarrelNumPlanes"); - } - catch(...) { - ; - } - if (nModules > 0 && nModules != int(numPlanes)) { - lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; - // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); - // make the code crash (incidentSvc does not work) - // Andre, Alvaro, assert replaced by exception - throw std::runtime_error("Incorrect number of planes (ECalBarrelNumPlanes) in xml file!"); - } // Readout is in the middle between two passive planes double offsetPassivePhi = caloDim.offset() + dPhi / 2.; double offsetReadoutPhi = caloDim.offset() + 0; diff --git a/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02_geo.cpp b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02_geo.cpp new file mode 100644 index 00000000..9f62dec8 --- /dev/null +++ b/detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02_geo.cpp @@ -0,0 +1,618 @@ +#include "DD4hep/DetFactoryHelper.h" +#include "DD4hep/Handle.h" +#include "XML/Utilities.h" + +#include + +// Taken from https://github.com/HEP-FCC/FCCDetectors/blob/main/Detector/DetFCChhECalInclined/src/ECalBarrelInclined_geo.cpp (go there to see older commit history) + +// todo: remove gaudi logging and properly capture output +#define endmsg std::endl +#define lLog std::cout +namespace MSG { +const std::string ERROR = " Error: "; +const std::string DEBUG = " Debug: "; +const std::string INFO = " Info: "; +} + +namespace det { +static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd, + dd4hep::xml::Handle_t aXmlElement, + dd4hep::SensitiveDetector aSensDet) { + + dd4hep::xml::DetElement xmlDetElem = aXmlElement; + std::string nameDet = xmlDetElem.nameStr(); + dd4hep::xml::Dimension dim(xmlDetElem.dimensions()); + dd4hep::DetElement caloDetElem(nameDet, xmlDetElem.id()); + + // Create air envelope for the whole barrel + dd4hep::Volume envelopeVol(nameDet + "_vol", dd4hep::Tube(dim.rmin(), dim.rmax(), dim.dz()), + aLcdd.material("Air")); + envelopeVol.setVisAttributes(aLcdd, dim.visStr()); + + // Retrieve cryostat data + dd4hep::xml::DetElement cryostat = aXmlElement.child(_Unicode(cryostat)); + dd4hep::xml::Dimension cryoDim(cryostat.dimensions()); + double cryoThicknessFront = cryoDim.rmin2() - cryoDim.rmin1(); + dd4hep::xml::DetElement cryoFront = cryostat.child(_Unicode(front)); + dd4hep::xml::DetElement cryoBack = cryostat.child(_Unicode(back)); + dd4hep::xml::DetElement cryoSide = cryostat.child(_Unicode(side)); + bool cryoFrontSensitive = cryoFront.isSensitive(); + bool cryoBackSensitive = cryoBack.isSensitive(); + bool cryoSideSensitive = cryoSide.isSensitive(); + + // Retrieve active and passive material data + dd4hep::xml::DetElement calo = aXmlElement.child(_Unicode(calorimeter)); + dd4hep::xml::Dimension caloDim(calo.dimensions()); + dd4hep::xml::DetElement active = calo.child(_Unicode(active)); + std::string activeMaterial = active.materialStr(); + double activeThickness = active.thickness(); + + dd4hep::xml::DetElement overlap = active.child(_Unicode(overlap)); + double activePassiveOverlap = overlap.offset(); + if (activePassiveOverlap < 0 || activePassiveOverlap > 0.5) { + // todo: ServiceHandle incidentSvc("IncidentSvc", "ECalConstruction"); + lLog << MSG::ERROR << "Overlap between active and passive cannot be more than half of passive plane!" << endmsg; + //todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + } + dd4hep::xml::DetElement layers = calo.child(_Unicode(layers)); + uint numLayers = 0; + std::vector layerHeight; + double layersTotalHeight = 0; + for (dd4hep::xml::Collection_t layer_coll(layers, _Unicode(layer)); layer_coll; ++layer_coll) { + dd4hep::xml::Component layer = layer_coll; + numLayers += layer.repeat(); + for (int iLay = 0; iLay < layer.repeat(); iLay++) { + layerHeight.push_back(layer.thickness()); + } + layersTotalHeight += layer.repeat() * layer.thickness(); + } + lLog << MSG::DEBUG << "Number of layers: " << numLayers << " total thickness " << layersTotalHeight << endmsg; + // The following code checks if the xml geometry file contains a constant defining + // the number of layers the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of layers is needed + // in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). + int nLayers = -1; + try { + nLayers = aLcdd.constant("ECalBarrelNumLayers"); + } + catch(...) { + ; + } + if (nLayers > 0 && nLayers != int(numLayers)) { + lLog << MSG::ERROR << "Incorrect number of layers (ECalBarrelNumLayers) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + // Andre, Alvaro, assert replaced by exception + throw std::runtime_error("Incorrect number of layers (ECalBarrelNumLayers) in xml file!"); + } + + dd4hep::xml::DetElement readout = calo.child(_Unicode(readout)); + std::string readoutMaterial = readout.materialStr(); + double readoutThickness = readout.thickness(); + + dd4hep::xml::DetElement passive = calo.child(_Unicode(passive)); + dd4hep::xml::DetElement passiveInner = passive.child(_Unicode(inner)); + dd4hep::xml::DetElement passiveInnerMax = passive.child(_Unicode(innerMax)); + dd4hep::xml::DetElement passiveOuter = passive.child(_Unicode(outer)); + dd4hep::xml::DetElement passiveGlue = passive.child(_Unicode(glue)); + std::string passiveInnerMaterial = passiveInner.materialStr(); + std::string passiveOuterMaterial = passiveOuter.materialStr(); + std::string passiveGlueMaterial = passiveGlue.materialStr(); + double passiveInnerThicknessMin = passiveInner.thickness(); + double passiveInnerThicknessMax = passiveInnerMax.thickness(); + double passiveOuterThickness = passiveOuter.thickness(); + double passiveGlueThickness = passiveGlue.thickness(); + double passiveThickness = passiveInnerThicknessMin + passiveOuterThickness + passiveGlueThickness; + double angle = passive.rotation().angle(); + + double bathRmin = caloDim.rmin(); // - margin for inclination + double bathRmax = caloDim.rmax(); // + margin for inclination + dd4hep::Tube bathOuterShape(bathRmin, bathRmax, caloDim.dz()); // make it 4 volumes + 5th for detector envelope + dd4hep::Tube bathAndServicesOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), caloDim.dz()); // make it 4 volumes + 5th for detector envelope + if (cryoThicknessFront > 0) { + // 1. Create cryostat + dd4hep::Tube cryoFrontShape(cryoDim.rmin1(), cryoDim.rmin2(), cryoDim.dz()); + dd4hep::Tube cryoBackShape(cryoDim.rmax1(), cryoDim.rmax2(), cryoDim.dz()); + dd4hep::Tube cryoSideOuterShape(cryoDim.rmin2(), cryoDim.rmax1(), cryoDim.dz()); + dd4hep::SubtractionSolid cryoSideShape(cryoSideOuterShape, bathAndServicesOuterShape); + lLog << MSG::INFO << "ECAL cryostat: front: rmin (cm) = " << cryoDim.rmin1() << " rmax (cm) = " << cryoDim.rmin2() << " dz (cm) = " << cryoDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL cryostat: back: rmin (cm) = " << cryoDim.rmax1() << " rmax (cm) = " << cryoDim.rmax2() << " dz (cm) = " << cryoDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL cryostat: side: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << cryoDim.dz() - caloDim.dz() << endmsg; + dd4hep::Volume cryoFrontVol(cryostat.nameStr()+"_front", cryoFrontShape, aLcdd.material(cryostat.materialStr())); + dd4hep::Volume cryoBackVol(cryostat.nameStr()+"_back", cryoBackShape, aLcdd.material(cryostat.materialStr())); + dd4hep::Volume cryoSideVol(cryostat.nameStr()+"_side", cryoSideShape, aLcdd.material(cryostat.materialStr())); + dd4hep::PlacedVolume cryoFrontPhysVol = envelopeVol.placeVolume(cryoFrontVol); + dd4hep::PlacedVolume cryoBackPhysVol = envelopeVol.placeVolume(cryoBackVol); + dd4hep::PlacedVolume cryoSidePhysVol = envelopeVol.placeVolume(cryoSideVol); + if (cryoFrontSensitive) { + cryoFrontVol.setSensitiveDetector(aSensDet); + cryoFrontPhysVol.addPhysVolID("cryo", 1); + cryoFrontPhysVol.addPhysVolID("type", 1); + lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; + } + if (cryoBackSensitive) { + cryoBackVol.setSensitiveDetector(aSensDet); + cryoBackPhysVol.addPhysVolID("cryo", 1); + cryoBackPhysVol.addPhysVolID("type", 2); + lLog << MSG::INFO << "Cryostat back volume set as sensitive" << endmsg; + } + if (cryoSideSensitive) { + cryoSideVol.setSensitiveDetector(aSensDet); + cryoSidePhysVol.addPhysVolID("cryo", 1); + cryoSidePhysVol.addPhysVolID("type", 3); + lLog << MSG::INFO << "Cryostat front volume set as sensitive" << endmsg; + } + dd4hep::DetElement cryoFrontDetElem(caloDetElem, "cryo_front", 0); + cryoFrontDetElem.setPlacement(cryoFrontPhysVol); + dd4hep::DetElement cryoBackDetElem(caloDetElem, "cryo_back", 0); + cryoBackDetElem.setPlacement(cryoBackPhysVol); + dd4hep::DetElement cryoSideDetElem(caloDetElem, "cryo_side", 0); + cryoSideDetElem.setPlacement(cryoSidePhysVol); + // 1.2. Create place-holder for services + dd4hep::Tube servicesFrontShape(cryoDim.rmin2(), bathRmin, caloDim.dz()); + dd4hep::Tube servicesBackShape(bathRmax, cryoDim.rmax1(), caloDim.dz()); + lLog << MSG::INFO << "ECAL services: front: rmin (cm) = " << cryoDim.rmin2() << " rmax (cm) = " << bathRmin << " dz (cm) = " << caloDim.dz() << endmsg; + lLog << MSG::INFO << "ECAL services: back: rmin (cm) = " << bathRmax << " rmax (cm) = " << cryoDim.rmax1() << " dz (cm) = " << caloDim.dz() << endmsg; + dd4hep::Volume servicesFrontVol("services_front", servicesFrontShape, aLcdd.material(activeMaterial)); + dd4hep::Volume servicesBackVol("services_back", servicesBackShape, aLcdd.material(activeMaterial)); + dd4hep::PlacedVolume servicesFrontPhysVol = envelopeVol.placeVolume(servicesFrontVol); + dd4hep::PlacedVolume servicesBackPhysVol = envelopeVol.placeVolume(servicesBackVol); + if (cryoFrontSensitive) { + servicesFrontVol.setSensitiveDetector(aSensDet); + servicesFrontPhysVol.addPhysVolID("cryo", 1); + servicesFrontPhysVol.addPhysVolID("type", 4); + lLog << MSG::INFO << "Services front volume set as sensitive" << endmsg; + } + if (cryoBackSensitive) { + servicesBackVol.setSensitiveDetector(aSensDet); + servicesBackPhysVol.addPhysVolID("cryo", 1); + servicesBackPhysVol.addPhysVolID("type", 5); + lLog << MSG::INFO << "Services back volume set as sensitive" << endmsg; + } + dd4hep::DetElement servicesFrontDetElem(caloDetElem, "services_front", 0); + servicesFrontDetElem.setPlacement(servicesFrontPhysVol); + dd4hep::DetElement servicesBackDetElem(caloDetElem, "services_back", 0); + servicesBackDetElem.setPlacement(servicesBackPhysVol); + } + // 2. Create bath that is inside the cryostat and surrounds the detector + // Bath is filled with active material -> but not sensitive + dd4hep::Volume bathVol(activeMaterial + "_bath", bathOuterShape, aLcdd.material(activeMaterial)); + lLog << MSG::INFO << "ECAL bath: material = " << activeMaterial << " rmin (cm) = " << bathRmin + << " rmax (cm) = " << bathRmax << " thickness in front of ECal (cm) = " << caloDim.rmin() - cryoDim.rmin2() + << " thickness behind ECal (cm) = " << cryoDim.rmax1() - caloDim.rmax() << endmsg; + + // 3. Create the calorimeter by placing the passive material, trapezoid active layers, readout and again trapezoid + // active layers in the bath. + // sensitive detector for the layers + dd4hep::SensitiveDetector sd = aSensDet; + dd4hep::xml::Dimension sdType = xmlDetElem.child(_U(sensitive)); + sd.setType(sdType.typeStr()); + lLog << MSG::INFO << "ECAL calorimeter volume rmin (cm) = " << caloDim.rmin() << " rmax (cm) = " << caloDim.rmax() + << endmsg; + + // 3.a. Create the passive planes, readout in between of 2 passive planes and the remaining space fill with active + // material + ////////////////////////////// + // PASSIVE PLANES + ////////////////////////////// + lLog << MSG::INFO << "passive inner material = " << passiveInnerMaterial << "\n" + << " and outer material = " << passiveOuterMaterial << "\n" + << " thickness of inner part at inner radius (cm) = " << passiveInnerThicknessMin << "\n" + << " thickness of inner part at outer radius (cm) = " << passiveInnerThicknessMax << "\n" + << " thickness of outer part (cm) = " << passiveOuterThickness << "\n" + << " thickness of total (cm) = " << passiveThickness << "\n" + << " rotation angle = " << angle << endmsg; + uint numPlanes = + round(M_PI / asin((passiveThickness + activeThickness + readoutThickness) / (2. * caloDim.rmin() * cos(angle)))); + + double dPhi = 2. * M_PI / numPlanes; + lLog << MSG::INFO << "number of passive plates = " << numPlanes << " azim. angle difference = " << dPhi << endmsg; + lLog << MSG::INFO << " distance at inner radius (cm) = " << 2. * M_PI * caloDim.rmin() / numPlanes << "\n" + << " distance at outer radius (cm) = " << 2. * M_PI * caloDim.rmax() / numPlanes << endmsg; + + // The following code checks if the xml geometry file contains a constant defining + // the number of planes in the barrel. In that case, it makes the program abort + // if the number of planes in the xml is different from the one calculated from + // the geometry. This is because the number of plane information (retrieved from the + // xml) is used in other parts of the code (the readout for the FCC-ee ECAL with + // inclined modules). In principle the code above should be refactored so that the number + // of planes is one of the inputs of the calculation and other geometrical parameters + // are adjusted accordingly. This is left for the future, and we use the workaround + // below to enforce for the time being that the number of planes is "correct" + int nModules = -1; + try { + nModules = aLcdd.constant("ECalBarrelNumPlanes"); + } + catch(...) { + ; + } + if (nModules > 0 && nModules != int(numPlanes)) { + lLog << MSG::ERROR << "Incorrect number of planes (ECalBarrelNumPlanes) in xml file!" << endmsg; + // todo: incidentSvc->fireIncident(Incident("ECalConstruction", "GeometryFailure")); + // make the code crash (incidentSvc does not work) + // Andre, Alvaro, assert replaced by exception + throw std::runtime_error("Incorrect number of planes (ECalBarrelNumPlanes) in xml file!"); + } + // Readout is in the middle between two passive planes + double offsetPassivePhi = caloDim.offset() + dPhi / 2.; + double offsetReadoutPhi = caloDim.offset() + 0; + lLog << MSG::INFO << "readout material = " << readoutMaterial << "\n" + << " thickness of readout planes (cm) = " << readoutThickness << "\n number of readout layers = " << numLayers + << endmsg; + double Rmin = caloDim.rmin(); + double Rmax = caloDim.rmax(); + double dR = Rmax - Rmin; + double planeLength = -Rmin * cos(angle) + sqrt(pow(Rmax, 2) - pow(Rmin * sin(angle), 2)); + lLog << MSG::INFO << "thickness of calorimeter (cm) = " << dR << "\n" + << " length of passive or readout planes (cm) = " << planeLength << endmsg; + + + // fill the thickness in the boundary of each layer + std::vector passiveInnerThicknessLayer(numLayers+1); + double runningHeight = 0; + for (uint iLay = 0; iLay < numLayers; iLay++) { + passiveInnerThicknessLayer[iLay] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * + (runningHeight) / (Rmax - Rmin); + runningHeight += layerHeight[iLay]; + } + passiveInnerThicknessLayer[numLayers] = passiveInnerThicknessMin + (passiveInnerThicknessMax - passiveInnerThicknessMin) * + (runningHeight) / (Rmax - Rmin); + + double passiveAngle = atan2((passiveInnerThicknessMax - passiveInnerThicknessMin) / 2., planeLength); + double cosPassiveAngle = cos(passiveAngle); + double rotatedOuterThickness = passiveOuterThickness / cosPassiveAngle; + double rotatedGlueThickness = passiveGlueThickness / cosPassiveAngle; + + // rescale layer thicknesses + double scaleLayerThickness = planeLength / layersTotalHeight; + layersTotalHeight = 0; + for (uint iLay = 0; iLay < numLayers; iLay++) { + layerHeight[iLay] *= scaleLayerThickness; + + layersTotalHeight += layerHeight[iLay]; + lLog << MSG::DEBUG << "Thickness of layer " << iLay << " : " << layerHeight[iLay] << endmsg; + } + double layerFirstOffset = -planeLength / 2. + layerHeight[0] / 2.; + + //dd4hep::Box passiveShape(passiveThickness / 2., caloDim.dz(), planeLength / 2.); + dd4hep::Trd1 passiveShape(passiveInnerThicknessMin / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., + passiveInnerThicknessMax / 2. + rotatedOuterThickness / 2. + rotatedGlueThickness / 2., + caloDim.dz(), planeLength / 2.); + // inner layer is not in the first calo layer (to sample more uniformly in the layer where upstream correction is + // applied) + //dd4hep::Box passiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); + dd4hep::Trd1 passiveInnerShape(passiveInnerThicknessLayer[1] / 2., passiveInnerThicknessMax / 2., caloDim.dz(), planeLength / 2. - layerHeight[0] / 2.); + //dd4hep::Box passiveInnerShapeFirstLayer(passiveInnerThickness / 2., caloDim.dz(), layerHeight[0] / 2.); + dd4hep::Trd1 passiveInnerShapeFirstLayer(passiveInnerThicknessMin / 2., passiveInnerThicknessLayer[1] / 2., caloDim.dz(), layerHeight[0] / 2.); + dd4hep::Box passiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); + dd4hep::Box passiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), planeLength / 2. / cosPassiveAngle); + // passive volume consists of inner part and two outer, joind by glue + dd4hep::Volume passiveVol("passive", passiveShape, aLcdd.material("Air")); + dd4hep::Volume passiveInnerVol(passiveInnerMaterial + "_passive", passiveInnerShape, + aLcdd.material(passiveInnerMaterial)); + dd4hep::Volume passiveInnerVolFirstLayer(activeMaterial + "_passive", passiveInnerShapeFirstLayer, + aLcdd.material(activeMaterial)); + dd4hep::Volume passiveOuterVol(passiveOuterMaterial + "_passive", passiveOuterShape, + aLcdd.material(passiveOuterMaterial)); + dd4hep::Volume passiveGlueVol(passiveGlueMaterial + "_passive", passiveGlueShape, + aLcdd.material(passiveGlueMaterial)); + + if (passiveInner.isSensitive()) { + lLog << MSG::DEBUG << "Passive inner volume set as sensitive" << endmsg; + // inner part starts at second layer + double layerOffset = layerFirstOffset + layerHeight[1] / 2.; + for (uint iLayer = 1; iLayer < numLayers; iLayer++) { + //dd4hep::Box layerPassiveInnerShape(passiveInnerThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Trd1 layerPassiveInnerShape(passiveInnerThicknessLayer[iLayer] / 2., passiveInnerThicknessLayer[iLayer+1] / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Volume layerPassiveInnerVol(passiveInnerMaterial, layerPassiveInnerShape, + aLcdd.material(passiveInnerMaterial)); + layerPassiveInnerVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveInnerPhysVol = + passiveInnerVol.placeVolume(layerPassiveInnerVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveInnerPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveInnerDetElem("layer", iLayer); + layerPassiveInnerDetElem.setPlacement(layerPassiveInnerPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + } + if (passiveOuter.isSensitive()) { + lLog << MSG::DEBUG << "Passive outer volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset / cosPassiveAngle; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerPassiveOuterShape(passiveOuterThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); + dd4hep::Volume layerPassiveOuterVol(passiveOuterMaterial, layerPassiveOuterShape, + aLcdd.material(passiveOuterMaterial)); + layerPassiveOuterVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveOuterPhysVol = + passiveOuterVol.placeVolume(layerPassiveOuterVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveOuterPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveOuterDetElem("layer", iLayer); + layerPassiveOuterDetElem.setPlacement(layerPassiveOuterPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; + } + } + } + if (passiveGlue.isSensitive()) { + lLog << MSG::DEBUG << "Passive glue volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset / cosPassiveAngle; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerPassiveGlueShape(passiveGlueThickness / 4., caloDim.dz(), layerHeight[iLayer] / 2. / cosPassiveAngle); + dd4hep::Volume layerPassiveGlueVol(passiveGlueMaterial, layerPassiveGlueShape, + aLcdd.material(passiveGlueMaterial)); + layerPassiveGlueVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerPassiveGluePhysVol = + passiveGlueVol.placeVolume(layerPassiveGlueVol, dd4hep::Position(0, 0, layerOffset)); + layerPassiveGluePhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerPassiveGlueDetElem("layer", iLayer); + layerPassiveGlueDetElem.setPlacement(layerPassiveGluePhysVol); + if (iLayer != numLayers - 1) { + layerOffset += (layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.) / cosPassiveAngle; + } + } + } + + dd4hep::PlacedVolume passiveInnerPhysVol = + passiveVol.placeVolume(passiveInnerVol, dd4hep::Position(0, 0, layerHeight[0] / 2.)); + dd4hep::PlacedVolume passiveInnerPhysVolFirstLayer = + passiveVol.placeVolume(passiveInnerVolFirstLayer, dd4hep::Position(0, 0, layerFirstOffset)); + dd4hep::PlacedVolume passiveOuterPhysVolBelow = passiveVol.placeVolume( + passiveOuterVol, + dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), + dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - + rotatedGlueThickness / 2. - rotatedOuterThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveOuterPhysVolAbove = passiveVol.placeVolume( + passiveOuterVol, + dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), + dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + + rotatedGlueThickness / 2. + rotatedOuterThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveGluePhysVolBelow = passiveVol.placeVolume( + passiveGlueVol, + dd4hep::Transform3D(dd4hep::RotationY(-passiveAngle), + dd4hep::Position(-(passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. - + rotatedGlueThickness / 4., 0, 0))); + dd4hep::PlacedVolume passiveGluePhysVolAbove = passiveVol.placeVolume( + passiveGlueVol, + dd4hep::Transform3D(dd4hep::RotationY(passiveAngle), + dd4hep::Position((passiveInnerThicknessMin + passiveInnerThicknessMax) / 4. + + rotatedGlueThickness / 4., 0, 0))); + passiveInnerPhysVol.addPhysVolID("subtype", 0); + passiveInnerPhysVolFirstLayer.addPhysVolID("subtype", 0); + passiveOuterPhysVolBelow.addPhysVolID("subtype", 1); + passiveOuterPhysVolAbove.addPhysVolID("subtype", 2); + passiveGluePhysVolBelow.addPhysVolID("subtype", 3); + passiveGluePhysVolAbove.addPhysVolID("subtype", 4); + if (passiveInner.isSensitive()) { + passiveInnerVolFirstLayer.setSensitiveDetector(aSensDet); + passiveInnerPhysVolFirstLayer.addPhysVolID("layer", 0); + dd4hep::DetElement passiveInnerDetElemFirstLayer("layer", 0); + passiveInnerDetElemFirstLayer.setPlacement(passiveInnerPhysVolFirstLayer); + } + + ////////////////////////////// + // READOUT PLANES + ////////////////////////////// + dd4hep::Box readoutShape(readoutThickness / 2., caloDim.dz(), planeLength / 2.); + dd4hep::Volume readoutVol(readoutMaterial, readoutShape, aLcdd.material(readoutMaterial)); + if (readout.isSensitive()) { + lLog << MSG::INFO << "Readout volume set as sensitive" << endmsg; + double layerOffset = layerFirstOffset; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Box layerReadoutShape(readoutThickness / 2., caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::Volume layerReadoutVol(readoutMaterial, layerReadoutShape, aLcdd.material(readoutMaterial)); + layerReadoutVol.setSensitiveDetector(aSensDet); + dd4hep::PlacedVolume layerReadoutPhysVol = + readoutVol.placeVolume(layerReadoutVol, dd4hep::Position(0, 0, layerOffset)); + layerReadoutPhysVol.addPhysVolID("layer", iLayer); + dd4hep::DetElement layerReadoutDetElem("layer", iLayer); + layerReadoutDetElem.setPlacement(layerReadoutPhysVol); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + } + + ////////////////////////////// + // ACTIVE + ////////////////////////////// + // thickness of active layers at inner radius and outer ( = distance between passive plane and readout plane) + // at inner radius: distance projected at plane perpendicular to readout plane + double activeInThickness = Rmin * sin(dPhi / 2.) * cos(angle); + activeInThickness -= passiveThickness * (0.5 - activePassiveOverlap); + // at outer radius: distance projected at plane perpendicular to readout plane + double activeOutThickness = (Rmin + planeLength) * sin(dPhi / 2.) * cos(angle); + // make correction for outer readius caused by inclination angle + // first calculate intersection of readout plane and plane parallel to shifted passive plane + double xIntersect = (Rmin * (tan(angle) - cos(dPhi / 2.) * tan(angle + dPhi / 2.)) - planeLength * sin(dPhi / 2.)) / + (tan(angle) - tan(angle + dPhi / 2.)); + double yIntersect = tan(angle) * xIntersect + Rmin * (sin(dPhi / 2.) - tan(angle)) + planeLength * sin(dPhi / 2.); + // distance from inner radius to intersection + double correction = + planeLength - sqrt(pow(xIntersect - Rmin * cos(dPhi / 2), 2) + pow(yIntersect - Rmin * sin(dPhi / 2), 2)); + // correction to the active thickness + activeOutThickness += 2. * correction * sin(dPhi / 4.); + activeOutThickness -= passiveThickness * (0.5 - activePassiveOverlap); + // print the active layer dimensions + double activeInThicknessAfterSubtraction = + 2. * activeInThickness - readoutThickness - 2. * activePassiveOverlap * passiveThickness; + double activeOutThicknessAfterSubtraction = + 2. * activeOutThickness - readoutThickness - 2. * activePassiveOverlap * + (passiveThickness + passiveInnerThicknessMax - passiveInnerThicknessMin); // correct thickness for trapezoid + lLog << MSG::INFO << "active material = " << activeMaterial + << " active layers thickness at inner radius (cm) = " << activeInThicknessAfterSubtraction + << " thickness at outer radious (cm) = " << activeOutThicknessAfterSubtraction << " making " + << (activeOutThicknessAfterSubtraction - activeInThicknessAfterSubtraction) * 100 / + activeInThicknessAfterSubtraction + << " % increase." << endmsg; + lLog << MSG::INFO + << "active passive initial overlap (before subtraction) (cm) = " << passiveThickness * activePassiveOverlap + << " = " << activePassiveOverlap * 100 << " %" << endmsg; + + // creating shape for rows of layers (active material between two passive planes, with readout in the middle) + // first define area between two passive planes, area can reach up to the symmetry axis of passive plane + dd4hep::Trd1 activeOuterShape(activeInThickness, activeOutThickness, caloDim.dz(), planeLength / 2.); + // subtract readout shape from the middle + dd4hep::SubtractionSolid activeShapeNoReadout(activeOuterShape, readoutShape); + + // make calculation for active plane that is inclined with 0 deg (= offset + angle) + double Cx = Rmin * cos(-angle) + planeLength / 2.; + double Cy = Rmin * sin(-angle); + double Ax = Rmin * cos(-angle + dPhi / 2.) + planeLength / 2. * cos(dPhi / 2.); + double Ay = Rmin * sin(-angle + dPhi / 2.) + planeLength / 2. * sin(dPhi / 2.); + double CAx = fabs(Ax - Cx); + double CAy = fabs(Ay - Cy); + double zprim, xprim; + zprim = CAx; + xprim = CAy; + + double Bx = Rmin * cos(-angle - dPhi / 2.) + planeLength / 2. * cos(-dPhi / 2.); + double By = Rmin * sin(-angle - dPhi / 2.) + planeLength / 2. * sin(-dPhi / 2.); + double CBx = fabs(Bx - Cx); + double CBy = fabs(By - Cy); + double zprimB, xprimB; + zprimB = CBx; + xprimB = CBy; + + // subtract passive volume above + dd4hep::SubtractionSolid activeShapeNoPassiveAbove( + activeShapeNoReadout, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), + dd4hep::Position(-fabs(xprim), 0, fabs(zprim)))); + // subtract passive volume below + dd4hep::SubtractionSolid activeShape( + activeShapeNoPassiveAbove, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), + dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB)))); + dd4hep::Volume activeVol("active", activeShape, aLcdd.material("Air")); + + std::vector layerPhysVols; + // place layers within active volume + std::vector layerInThickness; + std::vector layerOutThickness; + double layerIncreasePerUnitThickness = (activeOutThickness - activeInThickness) / layersTotalHeight; + for (uint iLay = 0; iLay < numLayers; iLay++) { + if (iLay == 0) { + layerInThickness.push_back(activeInThickness); + } else { + layerInThickness.push_back(layerOutThickness[iLay - 1]); + } + layerOutThickness.push_back(layerInThickness[iLay] + layerIncreasePerUnitThickness * layerHeight[iLay]); + } + double layerOffset = layerFirstOffset; + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::Trd1 layerOuterShape(layerInThickness[iLayer], layerOutThickness[iLayer], caloDim.dz(), layerHeight[iLayer] / 2.); + dd4hep::SubtractionSolid layerShapeNoReadout(layerOuterShape, readoutShape); + dd4hep::SubtractionSolid layerShapeNoPassiveAbove( + layerShapeNoReadout, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(-dPhi / 2.), + dd4hep::Position(-fabs(xprim), 0, fabs(zprim) - layerOffset))); + // subtract passive volume below + dd4hep::SubtractionSolid layerShape( + layerShapeNoPassiveAbove, passiveShape, + dd4hep::Transform3D(dd4hep::RotationY(dPhi / 2.), + dd4hep::Position(fabs(xprimB), 0, -fabs(zprimB) - layerOffset))); + dd4hep::Volume layerVol("layer", layerShape, aLcdd.material(activeMaterial)); + layerVol.setSensitiveDetector(aSensDet); + layerPhysVols.push_back(activeVol.placeVolume(layerVol, dd4hep::Position(0, 0, layerOffset))); + layerPhysVols.back().addPhysVolID("layer", iLayer); + if (iLayer != numLayers - 1) { + layerOffset += layerHeight[iLayer] / 2. + layerHeight[iLayer + 1] / 2.; + } + } + + dd4hep::DetElement bathDetElem(caloDetElem, "bath", 1); + std::vector activePhysVols; + // Next place elements: passive planes, readout planes and rows of layers + for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { + // first calculate positions of passive and readout planes + // PASSIVE + // calculate centre position of the plane without plane rotation + double phi = offsetPassivePhi + iPlane * dPhi; + double xRadial = (Rmin + planeLength / 2.) * cos(phi); + double yRadial = (Rmin + planeLength / 2.) * sin(phi); + // calculate position of the beginning of plane + double xRmin = Rmin * cos(phi); + double yRmin = Rmin * sin(phi); + // rotate centre by angle wrt beginning of plane + double xRotated = xRmin + (xRadial - xRmin) * cos(angle) - (yRadial - yRmin) * sin(angle); + double yRotated = yRmin + (xRadial - xRmin) * sin(angle) + (yRadial - yRmin) * cos(angle); + dd4hep::Transform3D transform(dd4hep::RotationX(-M_PI / 2.) // to get in XY plane + * + dd4hep::RotationY(M_PI / 2. // to get pointed towards centre + - + phi - angle), + dd4hep::Position(xRotated, yRotated, 0)); + dd4hep::PlacedVolume passivePhysVol = bathVol.placeVolume(passiveVol, transform); + passivePhysVol.addPhysVolID("module", iPlane); + passivePhysVol.addPhysVolID("type", 1); // 0 = active, 1 = passive, 2 = readout + dd4hep::DetElement passiveDetElem(bathDetElem, "passive" + std::to_string(iPlane), iPlane); + passiveDetElem.setPlacement(passivePhysVol); + + // READOUT + // calculate centre position of the plane without plane rotation + double phiRead = offsetReadoutPhi + iPlane * dPhi; + double xRadialRead = (Rmin + planeLength / 2.) * cos(phiRead); + double yRadialRead = (Rmin + planeLength / 2.) * sin(phiRead); + // calculate position of the beginning of plane + double xRminRead = Rmin * cos(phiRead); + double yRminRead = Rmin * sin(phiRead); + // rotate centre by angle wrt beginning of plane + double xRotatedRead = xRminRead + (xRadialRead - xRminRead) * cos(angle) - (yRadialRead - yRminRead) * sin(angle); + double yRotatedRead = yRminRead + (xRadialRead - xRminRead) * sin(angle) + (yRadialRead - yRminRead) * cos(angle); + dd4hep::Transform3D transformRead( + dd4hep::RotationX(-M_PI / 2.) // to get in XY plane + * + dd4hep::RotationY(M_PI / 2. // to get pointed towards centre + - + phiRead - angle), + dd4hep::Position(xRotatedRead, yRotatedRead, 0)); + dd4hep::PlacedVolume readoutPhysVol = bathVol.placeVolume(readoutVol, transformRead); + readoutPhysVol.addPhysVolID("module", iPlane); + readoutPhysVol.addPhysVolID("type", 2); // 0 = active, 1 = passive, 2 = readout + dd4hep::DetElement readoutDetElem(bathDetElem, "readout" + std::to_string(iPlane), iPlane); + readoutDetElem.setPlacement(readoutPhysVol); + + // ACTIVE + dd4hep::Rotation3D rotationActive(dd4hep::RotationX(-M_PI / 2) * + dd4hep::RotationY(M_PI / 2 - phiRead - angle)); + activePhysVols.push_back(bathVol.placeVolume( + activeVol, + dd4hep::Transform3D(rotationActive, dd4hep::Position(xRotatedRead, yRotatedRead, 0)))); + activePhysVols.back().addPhysVolID("module", iPlane); + activePhysVols.back().addPhysVolID("type", 0); // 0 = active, 1 = passive, 2 = readout + } + dd4hep::PlacedVolume bathPhysVol = envelopeVol.placeVolume(bathVol); + bathDetElem.setPlacement(bathPhysVol); + for (uint iPlane = 0; iPlane < numPlanes; iPlane++) { + dd4hep::DetElement activeDetElem(bathDetElem, "active" + std::to_string(iPlane), iPlane); + activeDetElem.setPlacement(activePhysVols[iPlane]); + for (uint iLayer = 0; iLayer < numLayers; iLayer++) { + dd4hep::DetElement layerDetElem(activeDetElem, "layer" + std::to_string(iLayer), iLayer); + layerDetElem.setPlacement(layerPhysVols[iLayer]); + } + } + + // Place the envelope + dd4hep::Volume motherVol = aLcdd.pickMotherVolume(caloDetElem); + dd4hep::PlacedVolume envelopePhysVol = motherVol.placeVolume(envelopeVol); + envelopePhysVol.addPhysVolID("system", xmlDetElem.id()); + caloDetElem.setPlacement(envelopePhysVol); + + // Create caloData object + auto caloData = new dd4hep::rec::LayeredCalorimeterData; + caloData->layoutType = dd4hep::rec::LayeredCalorimeterData::BarrelLayout; + caloDetElem.addExtension(caloData); + + // Set type flags + dd4hep::xml::setDetectorTypeFlag(xmlDetElem, caloDetElem); + + return caloDetElem; +} +} // namespace det + +DECLARE_DETELEMENT(ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02, det::createECalBarrelInclined) From 1199569da09e4ace4ddaf4376257f0d2db4ae0ad Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 2 Nov 2023 10:43:19 +0100 Subject: [PATCH 12/15] call to proper detector constructor version --- .../ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml index c8655ce7..4fd78287 100644 --- a/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml +++ b/FCCee/ALLEGRO/compact/ALLEGRO_o1_v02/FCCee_ECalBarrel_thetamodulemerged.xml @@ -96,7 +96,7 @@ - + From ec1947fa046c175e7dbe3517328ef06db2aa4774 Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 2 Nov 2023 10:43:48 +0100 Subject: [PATCH 13/15] include new segmentation in the segmentation factory --- detectorSegmentations/src/plugins/SegmentationFactories.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/detectorSegmentations/src/plugins/SegmentationFactories.cpp b/detectorSegmentations/src/plugins/SegmentationFactories.cpp index 3e65e7ba..4dd15369 100644 --- a/detectorSegmentations/src/plugins/SegmentationFactories.cpp +++ b/detectorSegmentations/src/plugins/SegmentationFactories.cpp @@ -18,7 +18,7 @@ DECLARE_SEGMENTATION(GridTheta_k4geo, create_segmentation) #include "detectorSegmentations/FCCSWGridModuleThetaMerged_k4geo.h" -DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged, create_segmentation) +DECLARE_SEGMENTATION(FCCSWGridModuleThetaMerged_k4geo, create_segmentation) #include "detectorSegmentations/FCCSWGridPhiEta_k4geo.h" DECLARE_SEGMENTATION(FCCSWGridPhiEta_k4geo, create_segmentation) From 30163eb1f5397bb1a97151a39ab1d2435932710a Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado Date: Thu, 2 Nov 2023 10:45:32 +0100 Subject: [PATCH 14/15] README updated --- detector/calorimeter/README.md | 3 +++ 1 file changed, 3 insertions(+) diff --git a/detector/calorimeter/README.md b/detector/calorimeter/README.md index d4596f46..1175c34a 100644 --- a/detector/calorimeter/README.md +++ b/detector/calorimeter/README.md @@ -11,6 +11,9 @@ The documentation about its usage is [here](../../doc/detector/calorimeter/ECalB ### o1_v01 Original version taken from [FCCDetectors](https://github.com/HEP-FCC/FCCDetectors/blob/main/Detector/DetFCChhECalInclined/src/ECalBarrelInclined_geo.cpp). +### o1_v02 +New version adapted to the theta segmentation with the possibility to have different number of cell merged per layer. The main difference is that now one has to set ECalBarrelNumLayers in the xml while before it was just dynamically computed based on other xml parameters (also, to avoid silent mistakes, the number from the xml and the one computed dynamically must match). + ## CaloDisks This sub-detector makes calorimeter endcaps (original and reflected). It is used in ALLEGRO detector concept. From a1e9dab7b7d9db476521abbffabc68a7d992818b Mon Sep 17 00:00:00 2001 From: Alvaro Tolosa Delgado <114166109+atolosadelgado@users.noreply.github.com> Date: Thu, 2 Nov 2023 18:26:05 +0100 Subject: [PATCH 15/15] Update detector/calorimeter/README.md Co-authored-by: Brieuc Francois --- detector/calorimeter/README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/detector/calorimeter/README.md b/detector/calorimeter/README.md index 1175c34a..e99fa66f 100644 --- a/detector/calorimeter/README.md +++ b/detector/calorimeter/README.md @@ -12,7 +12,7 @@ The documentation about its usage is [here](../../doc/detector/calorimeter/ECalB Original version taken from [FCCDetectors](https://github.com/HEP-FCC/FCCDetectors/blob/main/Detector/DetFCChhECalInclined/src/ECalBarrelInclined_geo.cpp). ### o1_v02 -New version adapted to the theta segmentation with the possibility to have different number of cell merged per layer. The main difference is that now one has to set ECalBarrelNumLayers in the xml while before it was just dynamically computed based on other xml parameters (also, to avoid silent mistakes, the number from the xml and the one computed dynamically must match). +New version adapted to the theta segmentation with the possibility to have different number of cell merged per layer. The main difference is that now one has to set `ECalBarrelNumLayers` and `ECalBarrelNumPlanes` in the xml while before it was just dynamically computed based on other xml parameters (also, to avoid silent mistakes, the number from the xml and the one computed dynamically must match). ## CaloDisks This sub-detector makes calorimeter endcaps (original and reflected). It is used in ALLEGRO detector concept.