Add parametrized fake DIRC with support

CMakeLists.txt

@@ -21,6 +21,7 @@ dd4hep_configure_output()
 
 dd4hep_add_plugin(${a_lib_name} SOURCES 
   src/BackwardGEMTracker_geo.cpp
+  src/BarrelBarDetectorWithSideFrame_geo.cpp
   src/BarrelCalorimeter_geo.cpp
   src/BarrelCalorimeterHybrid_geo.cpp
   src/BarrelCalorimeterInterlayers_geo.cpp

athena.xml

@@ -165,17 +165,14 @@
     ## Central tracking detectors
   </documentation>
 
-  <include ref="compact/tracking_config_accadia.xml"/>
+  <include ref="compact/tracking_config_acadia.xml"/>
 
   <!--include ref="compact/tracking_config_mgpds.xml"/-->
 
   <documentation level="10">
   ### PID detectors
   </documentation>
-  <!--include ref="compact/dirc.xml"/-->
-  <include ref="compact/mrich.xml"/>
-  <!--include ref="compact/forward_trd.xml"/-->
-  <include ref="compact/drich.xml"/>
+  <include ref="compact/pid_config_acadia.xml" />
 
   <documentation level="10">
   ## Central calorimetry

compact/definitions.xml

@@ -335,7 +335,7 @@ Examples:
     <constant name="BackwardTOF_length"           value="0.0*cm"/>
 
     <comment> Note: PID has space for DIRC and optional uRWEL layers</comment>
-    <constant name="BarrelPIDThickness"           value="10.0 * cm"/>
+    <constant name="BarrelPIDThickness"           value="15.0 * cm"/>
 
 
     <documentation>
@@ -454,7 +454,7 @@ end of the solenoid coils.
     <constant name="SolenoidYokeEndcapP_zmin"  value="Solenoid_length/2.0 + Solenoid_offset + EndcapPExtra_length"/>
     <constant name="SolenoidYokeEndcapN_zmin"  value="Solenoid_length/2.0 - Solenoid_offset + EndcapNExtra_length"/>
 
-    <constant name="EcalBarrel_rmin"            value="BarrelTrackingAndPID_rmax + 35.0 * mm"/>
+    <constant name="EcalBarrel_rmin"            value="BarrelTrackingAndPID_rmax + 10.0 * mm"/>
     <constant name="Barrel_rmax"                value="Solenoid_rmin - 1.0 *cm "/>
     <constant name="Barrel_TotalCalThickness"   value="Barrel_rmax - EcalBarrel_rmin"/>
 
@@ -614,7 +614,7 @@ The logic goes like this:
     <constant name="HcalEndcapN_thickness" value="HcalEndcapNLayer_NRepeat * HcalSingleLayerThickness"/>
 
     <comment> Special DIRC parameters </comment>
-    <constant name="DIRCFull_length"            value="max(CentralTracking_length + BackwardInnerEndcap_length + EcalEndcapN_length, EcalBarrel_length) + 50*cm"/>
+    <constant name="DIRCFull_length"            value="max(CentralTracking_length + BackwardInnerEndcap_length + EcalEndcapN_length, EcalBarrel_length) + 10*cm"/>
 		<constant name="DIRCExtra_offset"           value="-5*cm"/>
 
 

compact/display.xml

@@ -57,7 +57,14 @@
     <comment>
       PID Detectors
     </comment>
-    <vis name="RICHVis"      ref="AnlViolet"     showDaughters="true"/>
+    <vis name="RICHVis"         ref="AnlViolet"     showDaughters="true"/>
+
+    <vis name="DIRCVis"         ref="AnlTeal"       showDaughters="true" visible="false"/>
+    <vis name="DIRCLayerVis"    ref="AnlTeal"       showDaughters="true" visible="false"/>
+    <vis name="DIRCModuleVis"   ref="AnlBlue"       showDaughters="true" visible="false"/>
+    <vis name="DIRCSupportVis"  ref="VertexSupportVis" showDaughters="true" visible="true"/>
+    <vis name="DIRCFrameVis"    ref="AnlOrange"     showDaughters="true" visible="true"/>
+    <vis name="DIRCBarVis"      ref="AnlTeal"       showDaughters="true" visible="true"/>
 
     <comment>
       Beam line

compact/fake_dirc.xml

+<lccdd>
+
+  <comment>
+    DIRC consists of 16 (default) identical modules making a barrel.
+  </comment>
+
+  <define>
+    <constant name="DIRC_rmin"                  value="BarrelTracking_rmax"/>
+    <constant name="DIRC_rmax"                  value="DIRC_rmin + BarrelPIDThickness"/>
+    <constant name="DIRC_length"                value="DIRCFull_length"/>
+    <constant name="DIRC_offset"                value="(CentralTracking_length - DIRC_length)/2.0 - DIRCExtra_offset"/>
+
+    <constant name="DIRCModule_count"           value="16"  comment="Number of DIRC boxes per... DIRC. 16 - default"/>
+
+    <constant name="DIRCBar_thickness"          value="17*mm" />
+    <constant name="DIRCBar_length"             value="DIRC_length-8*cm" />
+
+    <constant name="DIRCSkinFront_thickness"    value="2.5*mm" />
+    <constant name="DIRCSkinBack_thickness"     value="2.5*mm" />
+
+    <constant name="DIRCFrame_width"            value="20*mm"/>
+    <constant name="DIRCFrame_length"           value="DIRC_length"/>
+    <constant name="DIRCFrame_thickness"        value="40*mm"/>
+
+    <constant name="DIRCModule_rmax"            value="DIRC_rmin + DIRCBar_thickness + DIRCSkinFront_thickness + DIRCSkinBack_thickness"/>
+    <constant name="DIRCModule_thickness"       value="DIRCModule_rmax - DIRC_rmin" />
+    <constant name="DIRCFrame_rmax"             value="DIRC_rmin + DIRCFrame_thickness" />
+
+    <constant name="DIRCModule_halfangle" value="180*degree/DIRCModule_count" />
+    <constant name="DIRCModule_width"           value="2*DIRC_rmin * tan(DIRCModule_halfangle)"/>
+
+    <constant name="DIRCModule_rmin"            value="DIRC_rmin + 1*cm"/>
+
+  </define>
+
+
+  <materials>
+  </materials>
+
+  <limits>
+  </limits>
+
+  <regions>
+  </regions>
+
+  <display>
+  </display>
+
+  <detectors>
+    <detector id="BarrelDIRC_ID" name="FakeDIRC" type="athena_FakeDIRC" readout="DIRCBarHits" vis="DIRCVis">
+      <dimensions rmin="DIRC_rmin" rmax="DIRC_rmax" length="DIRC_length" />
+      <position x="0" y="0" z="DIRC_offset" />
+      <comment> Fake DIRC modules </comment>
+      <module name="DIRCModule" vis="DIRCModuleVis">
+        <module_component name="FrontSkin"
+               material="CarbonFiber"
+               sensitive="false"
+               width="DIRCModule_width"
+               thickness="DIRCSkinFront_thickness"
+               vis="VertexSupportVis"
+               length="DIRCBar_length" />
+        <module_component name="QuartzBar"
+               material="Quartz"
+               sensitive="true"
+               width="DIRCModule_width"
+               thickness="DIRCBar_thickness"
+               vis="DIRCBarVis"
+               length="DIRCBar_length" />
+        <module_component name="BackSkin"
+               material="CarbonFiber"
+               sensitive="false"
+               width="DIRCModule_width"
+               thickness="DIRCSkinBack_thickness"
+               vis="DIRCSupportVis"
+               length="DIRCBar_length" />
+        <frame material="StainlessSteel"
+               width="DIRCFrame_width"
+               length="DIRCFrame_length"
+               vis="DIRCFrameVis"
+               thickness="DIRCFrame_thickness" />
+      </module>
+      <comment> Fake DIRC layers </comment>
+      <layer module="DIRCModule" id="1" vis="DIRCLayerVis">
+        <barrel_envelope
+          inner_r="DIRC_rmin"
+          outer_r="DIRC_rmax"
+          z_length="DIRC_length" />
+        <rphi_layout
+          phi_tilt="0"
+          nphi="DIRCModule_count"
+          phi0="0"
+          rc="0.5*(DIRCModule_rmin+DIRCModule_rmax)"
+          dr="0" />
+        <z_layout 
+          dr="0.0*mm"
+          z0="0.0*mm"
+          nz="1" />
+      </layer>
+    </detector>
+  </detectors>
+
+  <readouts>
+    <readout name="DIRCBarHits">
+      <segmentation type="CartesianGridXY" grid_size_x="3.0*mm" grid_size_y="3.0*mm" />
+      <id>system:8,barrel:2,layer:4,module:8,section:4,x:32:-16,y:-16</id>
+    </readout>
+  </readouts>
+
+  <plugins>
+  </plugins>
+
+  <fields>
+  </fields>
+</lccdd>

compact/pid_config_acadia.xml

+<lccdd>
+
+<include ref="fake_dirc.xml"/>
+<include ref="mrich.xml"/>
+<include ref="drich.xml"/>
+
+</lccdd>
+

compact/subsystem_views/README.mds

+# subsytem_views
+
+These compact detector files are copies of the main `athena.xml` compact description file, modified to only show certain subsystems.
+
+To add or update a subsystem:
+
+1. Copy `athena.xml` to this directory (`compact/subsystem_views`) with the appropriate name.
+2. Remove all the undesired detectors from this xml file. **Do not modify any of the included xml files** -- only modify the top level xml file here.
+3. Check that that you have the desired color scheme (eg. `colors.xml`) and display attributes (eg. `display_detailed.xml`) included.
+
+These compact detector files should not be used for any real simulations or studies.

compact/subsystem_views/inner_detector.xml

@@ -162,17 +162,14 @@
     ## Central tracking detectors
   </documentation>
 
-  <include ref="compact/tracking_config_accadia.xml"/>
+  <include ref="compact/tracking_config_acadia.xml"/>
 
   <!--include ref="compact/tracking_config_mgpds.xml"/-->
 
   <documentation level="10">
   ### PID detectors
   </documentation>
-  <!--include ref="compact/dirc.xml"/-->
-  <include ref="compact/mrich.xml"/>
-  <!--include ref="compact/forward_trd.xml"/-->
-  <include ref="compact/drich.xml"/>
+  <include ref="compact/pid_config_acadia.xml"/>
 
   <documentation level="10">
   ## Central calorimetry

compact/subsystem_views/pid_only.xml

@@ -163,9 +163,7 @@
       PID detectors
       -------------
   </comment>
-  <include ref="compact/dirc.xml"/>
-  <include ref="compact/mrich.xml"/>
-  <include ref="compact/drich.xml"/>
+  <include ref="compact/pid_config_acadia.xml"/>
 
 
   <comment>

compact/subsystem_views/tracking_only.xml

@@ -162,7 +162,7 @@
       Central tracking detectors
       --------------------------
   </comment>
-  <include ref="compact/tracking_config_accadia.xml"/>
+  <include ref="compact/tracking_config_acadia.xml"/>
 
   <comment>
       PID detectors

src/BarrelBarDetectorWithSideFrame_geo.cpp

+/** \addtogroup PID
+ * \brief Type: **Barrel bar detector (think DIRC) with longitudinal frame**.
+ * \author S. Joosten
+ *
+ * \ingroup PID
+ *
+ *
+ * \code
+ * \endcode
+ *
+ * @{
+ */
+#include "DD4hep/DetFactoryHelper.h"
+#include "DD4hep/Printout.h"
+#include "DD4hep/Shapes.h"
+#include "DDRec/DetectorData.h"
+#include "DDRec/Surface.h"
+#include "XML/Layering.h"
+
+using namespace std;
+using namespace dd4hep;
+
+/** Barrel Bar detector with optional framek
+ *
+ * - Optional "frame" tag within the module element (frame eats part of the module width
+ *   and is longitudinal at both sides)
+ * - Detector is setup as a "tracker" so we can use the hits to perform fast MC
+ *   for e.g. the DIRC
+ *
+ */
+static Ref_t create_BarrelBarDetectorWithSideFrame(Detector& description, xml_h e, SensitiveDetector sens)
+{
+  typedef vector<PlacedVolume>            Placements;
+  xml_det_t                               x_det    = e;
+  Material                                air      = description.air();
+  int                                     det_id   = x_det.id();
+  string                                  det_name = x_det.nameStr();
+  DetElement                              sdet(det_name, det_id);
+  map<string, Volume>                     volumes;
+  map<string, Placements>                 sensitives;
+  map<string, std::vector<rec::VolPlane>> volplane_surfaces;
+  PlacedVolume                            pv;
+  dd4hep::xml::Dimension                  dimensions(x_det.dimensions());
+  xml_dim_t                               dirc_pos = x_det.position();
+
+  map<string, std::array<double, 2>> module_thicknesses;
+
+  Tube   topVolumeShape(dimensions.rmin(), dimensions.rmax(), dimensions.length() * 0.5);
+  Volume assembly(det_name, topVolumeShape, air);
+
+  sens.setType("tracker");
+
+  // loop over the modules
+  for (xml_coll_t mi(x_det, _U(module)); mi; ++mi) {
+    xml_comp_t x_mod = mi;
+    string     m_nam = x_mod.nameStr();
+
+    if (volumes.find(m_nam) != volumes.end()) {
+      printout(ERROR, "BarrelBarDetectorWithSideFrame",
+               string((string("Module with named ") + m_nam + string(" already exists."))).c_str());
+      throw runtime_error("Logics error in building modules.");
+    }
+
+    int    ncomponents     = 0;
+    int    sensor_number   = 1;
+    double total_thickness = 0;
+
+    // Compute module total thickness from components
+    xml_coll_t ci(x_mod, _U(module_component));
+    for (ci.reset(), total_thickness = 0.0; ci; ++ci) {
+      total_thickness += xml_comp_t(ci).thickness();
+    }
+    // the module assembly volume
+    Assembly m_vol(m_nam);
+    volumes[m_nam] = m_vol;
+    m_vol.setVisAttributes(description, x_mod.visStr());
+
+    // Optional module frame.
+    // frame is 2 longitudinal bars at the side of the module
+    //
+    // |___| <-- example module cross section (x-y plane), frame is flush with the
+    //           bottom of the module and protrudes on the top if needed
+    //
+
+    // Get frame width, as it impacts the main module for being built. We
+    // construct the actual frame structure later (once we know the module width)
+    double frame_width = 0;
+    if (x_mod.hasChild("frame")) {
+      xml_comp_t m_frame = x_mod.child(_U(frame));
+      frame_width        = m_frame.width();
+    }
+
+    double thickness_so_far    = 0.0;
+    double thickness_sum       = -total_thickness / 2.0;
+    double max_component_width = 0;
+    for (xml_coll_t ci(x_mod, _U(module_component)); ci; ++ci, ++ncomponents) {
+      xml_comp_t x_comp = ci;
+      xml_comp_t x_pos  = x_comp.position(false);
+      string     c_nam  = _toString(ncomponents, "component%d");
+
+      double box_width    = x_comp.width() - 2 * frame_width;
+      max_component_width = fmax(max_component_width, box_width);
+
+      Box    c_box{box_width / 2, x_comp.length() / 2, x_comp.thickness() / 2};
+      Volume c_vol{c_nam, c_box, description.material(x_comp.materialStr())};
+
+      pv = m_vol.placeVolume(c_vol, Position(0, 0, thickness_sum + x_comp.thickness() / 2.0));
+
+      c_vol.setRegion(description, x_comp.regionStr());
+      c_vol.setLimitSet(description, x_comp.limitsStr());
+      c_vol.setVisAttributes(description, x_comp.visStr());
+      if (x_comp.isSensitive()) {
+        c_vol.setSensitiveDetector(sens);
+        sensitives[m_nam].push_back(pv);
+        module_thicknesses[m_nam] = {thickness_so_far + x_comp.thickness() / 2.0,
+                                     total_thickness - thickness_so_far - x_comp.thickness() / 2.0};
+        // -------- create a measurement plane for the tracking surface attched to the sensitive volume -----
+        rec::Vector3D u(0., 1., 0.);
+        rec::Vector3D v(0., 0., 1.);
+        rec::Vector3D n(1., 0., 0.);
+
+        // compute the inner and outer thicknesses that need to be assigned to the tracking surface
+        // depending on wether the support is above or below the sensor
+        double inner_thickness = module_thicknesses[m_nam][0];
+        double outer_thickness = module_thicknesses[m_nam][1];
+
+        rec::SurfaceType type(rec::SurfaceType::Sensitive);
+
+        rec::VolPlane surf(c_vol, type, inner_thickness, outer_thickness, u, v, n); //,o ) ;
+        volplane_surfaces[m_nam].push_back(surf);
+      }
+      thickness_sum += x_comp.thickness();
+      thickness_so_far += x_comp.thickness();
+    }
+    // Now add-on the frame
+    if (x_mod.hasChild("frame")) {
+      xml_comp_t m_frame         = x_mod.child(_U(frame));
+      double     frame_thickness = getAttrOrDefault<double>(m_frame, _U(thickness), total_thickness);
+
+      Box lframe_box{m_frame.width() / 2., m_frame.length() / 2., frame_thickness / 2.};
+      Box rframe_box{m_frame.width() / 2., m_frame.length() / 2., frame_thickness / 2.};
+
+      // Keep track of frame with so we can adjust the module bars appropriately
+
+      Volume lframe_vol{"left_frame", lframe_box, description.material(m_frame.materialStr())};
+      Volume rframe_vol{"right_frame", rframe_box, description.material(m_frame.materialStr())};
+      lframe_vol.setVisAttributes(description, m_frame.visStr());
+      rframe_vol.setVisAttributes(description, m_frame.visStr());
+
+      m_vol.placeVolume(lframe_vol, Position(frame_width / 2. + max_component_width / 2, 0.,
+                                             frame_thickness / 2. - total_thickness / 2.0));
+      m_vol.placeVolume(rframe_vol, Position(-frame_width / 2. - max_component_width / 2, 0.,
+                                             frame_thickness / 2. - total_thickness / 2.0));
+    }
+  }
+
+  // now build the layers
+  for (xml_coll_t li(x_det, _U(layer)); li; ++li) {
+    xml_comp_t x_layer  = li;
+    xml_comp_t x_barrel = x_layer.child(_U(barrel_envelope));
+    xml_comp_t x_layout = x_layer.child(_U(rphi_layout));
+    xml_comp_t z_layout = x_layer.child(_U(z_layout)); // Get the <z_layout> element.
+    int        lay_id   = x_layer.id();
+    string     m_nam    = x_layer.moduleStr();
+    string     lay_nam  = _toString(x_layer.id(), "layer%d");
+    Tube       lay_tub(x_barrel.inner_r(), x_barrel.outer_r(), x_barrel.z_length() / 2.0);
+    Volume     lay_vol(lay_nam, lay_tub, air); // Create the layer envelope volume.
+    lay_vol.setVisAttributes(description, x_layer.visStr());
+
+    double phi0     = x_layout.phi0();     // Starting phi of first module.
+    double phi_tilt = x_layout.phi_tilt(); // Phi tilt of a module.
+    double rc       = x_layout.rc();       // Radius of the module center.
+    int    nphi     = x_layout.nphi();     // Number of modules in phi.
+    double rphi_dr  = x_layout.dr();       // The delta radius of every other module.
+    double phi_incr = (M_PI * 2) / nphi;   // Phi increment for one module.
+    double phic     = phi0;                // Phi of the module center.
+    double z0       = z_layout.z0();       // Z position of first module in phi.
+    double nz       = z_layout.nz();       // Number of modules to place in z.
+    double z_dr     = z_layout.dr();       // Radial displacement parameter, of every other module.
+
+    Volume      module_env = volumes[m_nam];
+    DetElement  lay_elt(sdet, _toString(x_layer.id(), "layer%d"), lay_id);
+    Placements& sensVols = sensitives[m_nam];
+
+    // Z increment for module placement along Z axis.
+    // Adjust for z0 at center of module rather than
+    // the end of cylindrical envelope.
+    double z_incr = nz > 1 ? (2.0 * z0) / (nz - 1) : 0.0;
+    // Starting z for module placement along Z axis.
+    double module_z = -z0;
+    int    module   = 1;
+
+    // Loop over the number of modules in phi.
+    for (int ii = 0; ii < nphi; ii++) {
+      double dx = z_dr * std::cos(phic + phi_tilt); // Delta x of module position.
+      double dy = z_dr * std::sin(phic + phi_tilt); // Delta y of module position.
+      double x  = rc * std::cos(phic);              // Basic x module position.
+      double y  = rc * std::sin(phic);              // Basic y module position.
+
+      // Loop over the number of modules in z.
+      for (int j = 0; j < nz; j++) {
+        string     module_name = _toString(module, "module%d");
+        DetElement mod_elt(lay_elt, module_name, module);
+
+        Transform3D tr(RotationZYX(0, ((M_PI / 2) - phic - phi_tilt), -M_PI / 2), Position(x, y, module_z));
+
+        pv = lay_vol.placeVolume(module_env, tr);
+        pv.addPhysVolID("module", module);
+        pv.addPhysVolID("section", j);
+        mod_elt.setPlacement(pv);
+        for (size_t ic = 0; ic < sensVols.size(); ++ic) {
+          PlacedVolume sens_pv = sensVols[ic];
+          DetElement   comp_de(mod_elt, std::string("de_") + sens_pv.volume().name(), module);
+          comp_de.setPlacement(sens_pv);
+          rec::volSurfaceList(comp_de)->push_back(volplane_surfaces[m_nam][ic]);
+        }
+
+        /// Increase counters etc.
+        module++;
+        // Adjust the x and y coordinates of the module.
+        x += dx;
+        y += dy;
+        // Flip sign of x and y adjustments.
+        dx *= -1;
+        dy *= -1;
+        // Add z increment to get next z placement pos.
+        module_z += z_incr;
+      }
+      phic += phi_incr; // Increment the phi placement of module.
+      rc += rphi_dr;    // Increment the center radius according to dr parameter.
+      rphi_dr *= -1;    // Flip sign of dr parameter.
+      module_z = -z0;   // Reset the Z placement parameter for module.
+    }
+    // Create the PhysicalVolume for the layer.
+    pv = assembly.placeVolume(lay_vol); // Place layer in mother
+    pv.addPhysVolID("layer", lay_id);   // Set the layer ID.
+    lay_elt.setAttributes(description, lay_vol, x_layer.regionStr(), x_layer.limitsStr(), x_layer.visStr());
+    lay_elt.setPlacement(pv);
+  }
+  sdet.setAttributes(description, assembly, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());
+  assembly.setVisAttributes(description.invisible());
+  pv = description.pickMotherVolume(sdet).placeVolume(assembly, Position(0, 0, dirc_pos.z()));
+  pv.addPhysVolID("system", det_id); // Set the subdetector system ID.
+  pv.addPhysVolID("barrel", 1);      // Flag this as a barrel subdetector.
+  sdet.setPlacement(pv);
+  return sdet;
+}
+
+//@}
+// clang-format off
+DECLARE_DETELEMENT(BarrelBarDetectorWithSideFrame, create_BarrelBarDetectorWithSideFrame)
+DECLARE_DETELEMENT(athena_FakeDIRC, create_BarrelBarDetectorWithSideFrame)