Update version 2 and version 3 of LAr ECal (ECalBarrel_NobelLiquid_InclinedTrapezoid_o1_v01_geo.cpp) by including layered calorimeter data.

FCCee/CLD/compact/CLD_o4_v05/LAr_ECalBarrel.xml

@@ -45,13 +45,13 @@
     <constant name="LArBathThicknessBack" value="40*mm"/>
 
     <!-- air margin around calorimeter -->
-    <constant name="BarCryoECal_rmin" value="ECalBarrel_inner_radius+AirMarginThickness"/>
-    <constant name="BarCryoECal_rmax" value="ECalBarrel_outer_radius-AirMarginThickness"/>
-    <constant name="BarCryoECal_dz" value="ECalBarrel_half_length"/>
+    <constant name="BarCryoECal_rmin" value="BarECal_rmin+AirMarginThickness"/>
+    <constant name="BarCryoECal_rmax" value="BarECal_rmax-AirMarginThickness"/>
+    <constant name="BarCryoECal_dz" value="BarECal_dz"/>
     <!-- calorimeter active volume -->
     <constant name="EMBarrel_rmin" value="BarCryoECal_rmin+CryoBarrelFront+LArBathThicknessFront"/>
     <constant name="EMBarrel_rmax" value="BarCryoECal_rmax-CryoBarrelBack-LArBathThicknessBack"/>
-    <constant name="EMBarrel_dz" value="(ECalBarrel_half_length)-CryoBarrelSide"/>
+    <constant name="EMBarrel_dz" value="BarECal_dz-CryoBarrelSide"/>
     <!-- thickness of active volume between two absorber plates at EMBarrel_rmin, measuring perpendicular to the readout plate -->
     <constant name="LAr_thickness" value="LArGapThickness"/>
     <!-- passive layer consists of lead in the middle and steel on the outside, glued -->
@@ -71,7 +71,6 @@
     <constant name="Glue_thickness" value="0.1*mm"/>
     <!-- readout in between two absorber plates -->
     <constant name="readout_thickness" value="1.2*mm"/>
-    <constant name="ECal_cell_size" value="5.1*mm"/>
   </define>
 
   <display>
@@ -82,21 +81,21 @@
     <!-- readout for the simulation -->
     <!-- offset in eta is eta max value including cryostat -->
     <readout name="ECalBarrelEta">
-        <!-- segmentation type="GridTheta" grid_size_theta="0.5625" offset_theta="-0.83"/ -->
-	<segmentation type="GridEta" grid_size_eta="0.01" offset_eta="-1.2"/>
-      <id>system:5,cryo:1,type:3,subtype:3,layer:8,module:11,eta:9</id>
+        <!-- segmentation type="GridTheta_k4geo" grid_size_theta="0.5625" offset_theta="-0.83"/ -->
+        <segmentation type="GridEta_k4geo" grid_size_eta="0.01" offset_eta="-1.2"/>
+      <id>system:4,cryo:1,type:3,subtype:3,layer:8,module:11,eta:9</id>
     </readout>
     <!-- readout for the reconstruction -->
     <!-- phi position is calculated based on the centre of volume (hence it cannot be done in the simulation from energy deposits position) -->
     <readout name="ECalBarrelPhiEta">
-        <!-- segmentation type="FCCSWGridPhiTheta" grid_size_theta="0.5625" phi_bins="768" offset_theta="-0.83" offset_phi="-pi+(pi/768.)"/ -->
-        <segmentation type="FCCSWGridPhiEta" grid_size_eta="0.01" phi_bins="768" offset_eta="-1.2" offset_phi="-pi+(pi/768.)"/>
-      <id>system:5,cryo:1,type:3,subtype:3,layer:8,eta:9,phi:10</id>
+        <!-- segmentation type="FCCSWGridPhiTheta_k4geo" grid_size_theta="0.5625" phi_bins="768" offset_theta="-0.83" offset_phi="-pi+(pi/768.)"/ -->
+        <segmentation type="FCCSWGridPhiEta_k4geo" grid_size_eta="0.01" phi_bins="768" offset_eta="-1.2" offset_phi="-pi+(pi/768.)"/>
+      <id>system:4,cryo:1,type:3,subtype:3,layer:8,eta:9,phi:10</id>
     </readout>
   </readouts>
 
   <detectors>
-    <detector id="DetID_ECal_Barrel" name="ECalBarrel" type="ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01" readout="ECalBarrelEta">
+    <detector id="BarECal_id" name="ECalBarrel" type="ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01" readout="ECalBarrelEta">
       <type_flags type=" DetType_CALORIMETER + DetType_ELECTROMAGNETIC + DetType_BARREL"/>
       <sensitive type="SimpleCalorimeterSD"/>
       <dimensions rmin="BarCryoECal_rmin" rmax="BarCryoECal_rmax" dz="BarCryoECal_dz" vis="ecal_envelope"/>

detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v01_geo.cpp

@@ -624,8 +624,8 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd,
     caloLayer.outer_thickness           = difference_bet_r1r2 / 2;
 
     caloLayer.absorberThickness         = absorberThickness;
-    caloLayer.cellSize0 = 2 * dd4hep::mm;
-    caloLayer.cellSize1 = 2 * dd4hep::mm;
+    caloLayer.cellSize0 = 20 * dd4hep::mm;
+    caloLayer.cellSize1 = 20 * dd4hep::mm;
 
     caloData->layers.push_back(caloLayer);
   }

detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v02_geo.cpp

@@ -1,6 +1,8 @@
 #include "DD4hep/DetFactoryHelper.h"
 #include "DD4hep/Handle.h"
 #include "XML/Utilities.h"
+#include "DDRec/MaterialManager.h"
+#include "DDRec/Vector3D.h"
 
 #include <DDRec/DetectorData.h>
 
@@ -656,9 +658,70 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd,
   caloData->layoutType = dd4hep::rec::LayeredCalorimeterData::BarrelLayout;
   caloDetElem.addExtension<dd4hep::rec::LayeredCalorimeterData>(caloData);
 
+  caloData->extent[0] = Rmin;
+  caloData->extent[1] = Rmax; // or r_max ?
+  caloData->extent[2] = 0.;      // NN: for barrel detectors this is 0
+  caloData->extent[3] = caloDim.dz();
+
   // Set type flags
   dd4hep::xml::setDetectorTypeFlag(xmlDetElem, caloDetElem);
-
+  dd4hep::rec::MaterialManager matMgr(envelopeVol);
+  dd4hep::rec::LayeredCalorimeterData::Layer caloLayer;
+
+  double rad_first = Rmin;
+  double rad_last = 0;
+  double scale_fact = dR / (-Rmin * cos(angle) + sqrt(pow(Rmax, 2) - pow(Rmin * sin(angle), 2)));
+  // since the layer height is given along the electrode and not along the radius it needs to be scaled to get the values of layer height radially
+  std::cout << "Scaling factor " << scale_fact << std::endl;
+  for (size_t il = 0; il < layerHeight.size(); il++) {
+    double thickness_sen = 0.;
+    double absorberThickness = 0.;
+
+    rad_last = rad_first + (layerHeight[il] * scale_fact);
+    dd4hep::rec::Vector3D ivr1 = dd4hep::rec::Vector3D(0., rad_first, 0); // defining starting vector points of the given layer
+    dd4hep::rec::Vector3D ivr2 = dd4hep::rec::Vector3D(0., rad_last, 0);  // defining end vector points of the given layer
+
+    std::cout << "radius first " << rad_first << " radius last " << rad_last << std::endl;
+    const dd4hep::rec::MaterialVec &materials = matMgr.materialsBetween(ivr1, ivr2); // calling material manager to get material info between two points
+    auto mat = matMgr.createAveragedMaterial(materials);                             // creating average of all the material between two points to calculate X0 and lambda of averaged material
+    const double nRadiationLengths = mat.radiationLength();
+    const double nInteractionLengths = mat.interactionLength();
+    const double difference_bet_r1r2 = (ivr1 - ivr2).r();
+    const double value_of_x0 = layerHeight[il] / nRadiationLengths;
+    const double value_of_lambda = layerHeight[il] / nInteractionLengths;
+    std::string str1("LAr");
+
+    for (size_t imat = 0; imat < materials.size(); imat++) {
+
+      std::string str2(materials.at(imat).first.name());
+      if (str1.compare(str2) == 0){
+	  thickness_sen += materials.at(imat).second;
+      }
+      else {
+	absorberThickness += materials.at(imat).second;
+      }
+    }
+    rad_first = rad_last;
+    std::cout << "The sensitive thickness is " << thickness_sen << std::endl;
+    std::cout << "The absorber thickness is " << absorberThickness << std::endl;
+    std::cout << "The radiation length is " << value_of_x0 << " and the interaction length is " << value_of_lambda << std::endl;
+
+    caloLayer.distance = rad_first;
+    caloLayer.sensitive_thickness       = thickness_sen;
+    caloLayer.inner_nRadiationLengths   = value_of_x0 / 2.0;
+    caloLayer.inner_nInteractionLengths = value_of_lambda / 2.0;
+    caloLayer.inner_thickness           = difference_bet_r1r2 / 2.0;
+
+    caloLayer.outer_nRadiationLengths   = value_of_x0 / 2.0;
+    caloLayer.outer_nInteractionLengths = value_of_lambda / 2.0;
+    caloLayer.outer_thickness           = difference_bet_r1r2 / 2;
+
+    caloLayer.absorberThickness         = absorberThickness;
+    caloLayer.cellSize0 = 20 * dd4hep::mm;
+    caloLayer.cellSize1 = 20 * dd4hep::mm;
+
+    caloData->layers.push_back(caloLayer);
+  }
   return caloDetElem;
 }
 }  // namespace det

detector/calorimeter/ECalBarrel_NobleLiquid_InclinedTrapezoids_o1_v03_geo.cpp

@@ -1,6 +1,8 @@
 #include "DD4hep/DetFactoryHelper.h"
 #include "DD4hep/Handle.h"
 #include "XML/Utilities.h"
+#include "DDRec/MaterialManager.h"
+#include "DDRec/Vector3D.h"
 
 #include <DDRec/DetectorData.h>
 
@@ -687,8 +689,72 @@ static dd4hep::detail::Ref_t createECalBarrelInclined(dd4hep::Detector& aLcdd,
   caloData->layoutType = dd4hep::rec::LayeredCalorimeterData::BarrelLayout;
   caloDetElem.addExtension<dd4hep::rec::LayeredCalorimeterData>(caloData);
 
+  caloData->extent[0] = Rmin;
+  caloData->extent[1] = Rmax; // or r_max ?
+  caloData->extent[2] = 0.;      // NN: for barrel detectors this is 0
+  caloData->extent[3] = caloDim.dz();
+
+
   // Set type flags
   dd4hep::xml::setDetectorTypeFlag(xmlDetElem, caloDetElem);
+  dd4hep::rec::MaterialManager matMgr(envelopeVol);
+  dd4hep::rec::LayeredCalorimeterData::Layer caloLayer;
+
+  double rad_first = Rmin;
+  double rad_last = 0;
+  double scale_fact = dR / (-Rmin * cos(angle) + sqrt(pow(Rmax, 2) - pow(Rmin * sin(angle), 2)));
+  // since the layer height is given along the electrode and not along the radius it needs to be scaled to get the values of layer height radially
+  std::cout << "Scaling factor " << scale_fact << std::endl;
+  for (size_t il = 0; il < layerHeight.size(); il++) {
+    double thickness_sen = 0.;
+    double absorberThickness = 0.;
+
+    rad_last = rad_first + (layerHeight[il] * scale_fact);
+    dd4hep::rec::Vector3D ivr1 = dd4hep::rec::Vector3D(0., rad_first, 0); // defining starting vector points of the given layer
+    dd4hep::rec::Vector3D ivr2 = dd4hep::rec::Vector3D(0., rad_last, 0);  // defining end vector points of the given layer
+
+    std::cout << "radius first " << rad_first << " radius last " << rad_last << std::endl;
+    const dd4hep::rec::MaterialVec &materials = matMgr.materialsBetween(ivr1, ivr2); // calling material manager to get material info between two points
+    auto mat = matMgr.createAveragedMaterial(materials);                             // creating average of all the material between two points to calculate X0 and lambda of averaged material
+    const double nRadiationLengths = mat.radiationLength();
+    const double nInteractionLengths = mat.interactionLength();
+    const double difference_bet_r1r2 = (ivr1 - ivr2).r();
+    const double value_of_x0 = layerHeight[il] / nRadiationLengths;
+    const double value_of_lambda = layerHeight[il] / nInteractionLengths;
+    std::string str1("LAr");
+
+    for (size_t imat = 0; imat < materials.size(); imat++) {
+
+      std::string str2(materials.at(imat).first.name());
+      if (str1.compare(str2) == 0){
+	  thickness_sen += materials.at(imat).second;
+      }
+      else {
+	absorberThickness += materials.at(imat).second;
+      }
+    }
+    rad_first = rad_last;
+    std::cout << "The sensitive thickness is " << thickness_sen << std::endl;
+    std::cout << "The absorber thickness is " << absorberThickness << std::endl;
+    std::cout << "The radiation length is " << value_of_x0 << " and the interaction length is " << value_of_lambda << std::endl;
+
+    caloLayer.distance = rad_first;
+    caloLayer.sensitive_thickness       = thickness_sen;
+    caloLayer.inner_nRadiationLengths   = value_of_x0 / 2.0;
+    caloLayer.inner_nInteractionLengths = value_of_lambda / 2.0;
+    caloLayer.inner_thickness           = difference_bet_r1r2 / 2.0;
+
+    caloLayer.outer_nRadiationLengths   = value_of_x0 / 2.0;
+    caloLayer.outer_nInteractionLengths = value_of_lambda / 2.0;
+    caloLayer.outer_thickness           = difference_bet_r1r2 / 2;
+
+    caloLayer.absorberThickness         = absorberThickness;
+    caloLayer.cellSize0 = 20 * dd4hep::mm;
+    caloLayer.cellSize1 = 20 * dd4hep::mm;
+
+    caloData->layers.push_back(caloLayer);
+  }
+
 
   return caloDetElem;
 }