From fbb0f7061ca9013cf43f4167905004a869b95404 Mon Sep 17 00:00:00 2001
From: Chao Peng <cpeng@anl.gov>
Date: Fri, 11 Jun 2021 04:15:19 +0000
Subject: [PATCH] Add a generic plugin for homogeneous calorimeter
---
src/GeometryHelpers.cpp | 5 +
src/HomogeneousCalorimeter_geo.cpp | 310 +++++++++++++++++++++++++++++
2 files changed, 315 insertions(+)
create mode 100644 src/HomogeneousCalorimeter_geo.cpp
diff --git a/src/GeometryHelpers.cpp b/src/GeometryHelpers.cpp
index c60c239b..201b88d6 100644
--- a/src/GeometryHelpers.cpp
+++ b/src/GeometryHelpers.cpp
@@ -60,6 +60,11 @@ void add_square(Point p, std::vector<Point> &res, double lside, double rmin, dou
// fill squares
std::vector<Point> fillSquares(Point ref, double lside, double rmin, double rmax, double phmin, double phmax)
{
+ // convert (0, 2pi) to (-pi, pi)
+ if (phmax > M_PI) {
+ phmin -= M_PI;
+ phmax -= M_PI;
+ }
// start with a seed square and find one in the ring
// move to center
ref = ref - Point(int(ref.x()/lside)*lside, int(ref.y()/lside)*lside);
diff --git a/src/HomogeneousCalorimeter_geo.cpp b/src/HomogeneousCalorimeter_geo.cpp
new file mode 100644
index 00000000..f3d74c7f
--- /dev/null
+++ b/src/HomogeneousCalorimeter_geo.cpp
@@ -0,0 +1,310 @@
+//==========================================================================
+// A general implementation for homogeneous calorimeter
+// it supports three types of placements
+// 1. Module placement with module dimensions and positions
+// 2. Array placement with module dimensions and numbers of row and columns
+// 3. Disk placement with module dimensions and (Rmin, Rmax), and (Phimin, Phimax)
+//--------------------------------------------------------------------------
+// Author: Chao Peng (ANL)
+// Date: 06/09/2021
+//==========================================================================
+#include "GeometryHelpers.h"
+#include "DD4hep/DetFactoryHelper.h"
+#include <XML/Helper.h>
+#include <iostream>
+#include <algorithm>
+#include <math.h>
+
+using namespace dd4hep;
+using namespace dd4hep::detail;
+
+/** \addtogroup calorimeters Calorimeters
+ */
+
+/** \addtogroup Homogeneous Calorimeter
+ * \brief Type: **HomogeneousCalorimeter**.
+ * \author C. Peng
+ * \ingroup calorimeters
+ *
+ *
+ * \code
+ * <detector id="1" name="HyCal" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis">
+ * <dimensions shape="box" sizex="120*cm" sizey="120*cm" sizez="46*cm"/>
+ * <position x="0" y="0" z="0"/>
+ * <rotation x="0" y="0" z="0"/>
+ * <placements>
+ * <array nrow="34" ncol="34" sector="1">
+ * <position x="0" y="0" z="-9.73*cm"/>
+ * <module sizex="2.05*cm" sizey="2.05*cm" sizez="18*cm" vis="GreenVis" material="PbWO4"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * <removal row="16" col="16"/>
+ * <removal row="16" col="17"/>
+ * <removal row="17" col="16"/>
+ * <removal row="17" col="17"/>
+ * </array>
+ * <array nrow="6" ncol="24" sector="2">
+ * <position x="-17*(2.05+0.015)*cm+12*(3.8+0.015)*cm" y="17*(2.05+0.015)*cm+3*(3.8+0.015)*cm" z="0"/>
+ * <module sizex="3.8*cm" sizey="3.8*cm" sizez="45*cm" vis="BlueVis" material="PbGlass"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * </array>
+ * <array nrow="24" ncol="6" sector="3">
+ * <position x="17*(2.05+0.015)*cm+3*(3.8+0.015)*cm" y="17*(2.05+0.015)*cm-12*(3.8+0.015)*cm" z="0"/>
+ * <module sizex="3.8*cm" sizey="3.8*cm" sizez="45*cm" vis="BlueVis" material="PbGlass"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * </array>
+ * <array nrow="6" ncol="24" sector="4">
+ * <position x="17*(2.05+0.015)*cm-12*(3.8+0.015)*cm" y="-17*(2.05+0.015)*cm-3*(3.8+0.015)*cm" z="0"/>
+ * <module sizex="3.8*cm" sizey="3.8*cm" sizez="45*cm" vis="BlueVis" material="PbGlass"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * </array>
+ * <array nrow="24" ncol="6" sector="3">
+ * <position x="-17*(2.05+0.015)*cm-3*(3.8+0.015)*cm" y="-17*(2.05+0.015)*cm+12*(3.8+0.015)*cm" z="0"/>
+ * <module sizex="3.8*cm" sizey="3.8*cm" sizez="45*cm" vis="BlueVis" material="PbGlass"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * </array>
+ * </placements>
+ * </detector>
+ *
+ * <detector id="2" name="SomeBlocks" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis">
+ * <dimensions shape="box" sizex="100*cm" sizey="100*cm" sizez="20.5*cm"/>
+ * <position x="0" y="0" z="30*cm"/>
+ * <rotation x="0" y="0" z="0"/>
+ * <placements>
+ * <blocks sector="1"/>
+ * <module sizex="2.05*cm" sizey="2.05*cm" sizez="20*cm" vis="GreenVis" material="PbWO4"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * <placement x="1*cm" y="1*cm" z="0"/>
+ * <placement x="-1*cm" y="1*cm" z="0"/>
+ * <placement x="1*cm" y="-1*cm" z="0"/>
+ * <placement x="-1*cm" y="-1*cm" z="0"/>
+ * </blocks>
+ * </placements>
+ * </detector>
+ *
+ * <detector id="2" name="DiskShapeCalorimeter" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis">
+ * <dimensions shape="disk" rmin="25*cm" rmax="125*cm" length="20.5*cm" phimin="0" phimax="360*degree"/>
+ * <position x="0" y="0" z="-30*cm"/>
+ * <rotation x="0" y="0" z="0"/>
+ * <placements>
+ * <disk rmin="25*cm" rmax="125*cm" length="20.5*cm" phimin="0" phimax="360*degree" sector="1"/>
+ * <module sizex="2.05*cm" sizey="2.05*cm" sizez="20*cm" vis="GreenVis" material="PbWO4"/>
+ * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/>
+ * </placements>
+ * </detector>
+ * \endcode
+ *
+ * @{
+ */
+
+// headers
+static void add_blocks(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id);
+static void add_array(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id);
+static void add_disk(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id);
+
+// helper function to get x, y, z if defined in a xml component
+template<class XmlComp>
+Position get_xml_xyz(XmlComp &comp, dd4hep::xml::Strng_t name)
+{
+ Position pos(0., 0., 0.);
+ if (comp.hasChild(name)) {
+ auto child = comp.child(name);
+ pos.SetX(dd4hep::getAttrOrDefault<double>(child, _Unicode(x), 0.));
+ pos.SetY(dd4hep::getAttrOrDefault<double>(child, _Unicode(y), 0.));
+ pos.SetZ(dd4hep::getAttrOrDefault<double>(child, _Unicode(z), 0.));
+ }
+ return pos;
+}
+
+// main
+static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDetector sens)
+{
+ static const std::string func = "HomogeneousCalorimeter";
+ xml::DetElement detElem = handle;
+ std::string detName = detElem.nameStr();
+ int detID = detElem.id();
+ DetElement det(detName, detID);
+ sens.setType("calorimeter");
+
+ // top-level children
+ xml::Component dims = detElem.dimensions();
+
+ // build envelop from dimensions
+ std::string shape = dd4hep::getAttrOrDefault(dims, _Unicode(shape), "");
+ // no shape input, try to determine shape from dimension variables
+ if (shape.empty()) {
+ if (dims.hasAttr(_Unicode(rmin)) && dims.hasAttr(_Unicode(rmax)) && dims.hasAttr(_Unicode(length))) {
+ shape = "disk";
+ } else if (dims.hasAttr(_Unicode(sizex)) && dims.hasAttr(_Unicode(sizey)) && dims.hasAttr(_Unicode(sizez))) {
+ shape = "box";
+ } else {
+ std::cerr << func << " Error: Cannot determine shape of the calorimeter. "
+ "Add shape (box, or disk) into dimensions\n";
+ return det;
+ }
+ }
+
+ Volume envVol(detName + "_envelope");
+ std::string filler = dd4hep::getAttrOrDefault(detElem, _Unicode(filler), "Air");
+ envVol.setMaterial(desc.material(filler));
+ envVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
+
+ // convert to lower case
+ std::transform(shape.begin(), shape.end(), shape.begin(), [](unsigned char c){ return std::tolower(c); });
+
+ if (shape == "box") {
+ double sx = dims.attr<double>(_Unicode(sizex));
+ double sy = dims.attr<double>(_Unicode(sizey));
+ double sz = dims.attr<double>(_Unicode(sizez));
+ envVol.setSolid(Box(sx/2., sy/2., sz/2.));
+ } else if (shape == "disk") {
+ double rmin = dims.rmin();
+ double rmax = dims.rmax();
+ double length = dims.length();
+ double phimin = dd4hep::getAttrOrDefault<double>(dims, _Unicode(phimin), 0.);
+ double phimax = dd4hep::getAttrOrDefault<double>(dims, _Unicode(phimax), 2*M_PI);
+ envVol.setSolid(Tube(rmin, rmax, length/2., phimin, phimax));
+ } else {
+ std::cerr << func << " Error: Unsupported shape " << shape << ", use box or disk\n";
+ return det;
+ }
+
+ // module placement
+ xml::Component plm = detElem.child(_Unicode(placements));
+ int sector = 1;
+ for (xml::Collection_t arr(plm, _Unicode(array)); arr; ++arr) {
+ add_array(desc, envVol, arr, sens, sector++);
+ }
+ for (xml::Collection_t mod(plm, _Unicode(blocks)); mod; ++mod) {
+ add_blocks(desc, envVol, mod, sens, sector++);
+ }
+ for (xml::Collection_t disk(plm, _Unicode(disk)); disk; ++disk) {
+ add_disk(desc, envVol, disk, sens, sector++);
+ }
+
+ // detector position and rotation
+ auto pos = get_xml_xyz(detElem, _Unicode(position));
+ auto rot = get_xml_xyz(detElem, _Unicode(rotation));
+ Volume motherVol = desc.pickMotherVolume(det);
+ Transform3D tr = Translation3D(pos.x(), pos.y(), pos.z()) * RotationZYX(rot.z(), rot.y(), rot.x());
+ PlacedVolume envPV = motherVol.placeVolume(envVol, tr);
+ envPV.addPhysVolID("system", detID);
+ det.setPlacement(envPV);
+ return det;
+}
+
+// helper function to build module with or w/o wrapper
+Volume build_module(Detector &desc, xml::Collection_t &plm, SensitiveDetector &sens, Position &dim)
+{
+ auto mod = plm.child(_Unicode(module));
+ auto sx = mod.attr<double>(_Unicode(sizex));
+ auto sy = mod.attr<double>(_Unicode(sizey));
+ auto sz = mod.attr<double>(_Unicode(sizez));
+ dim = Position{sx, sy, sz};
+ Box modShape(sx/2., sy/2., sz/2.);
+ auto modMat = desc.material(mod.attr<std::string>(_Unicode(material)));
+ Volume modVol("module_vol", modShape, modMat);
+ modVol.setSensitiveDetector(sens);
+ modVol.setVisAttributes(desc.visAttributes(mod.attr<std::string>(_Unicode(vis))));
+
+ // no wrapper
+ if (!plm.hasChild(_Unicode(wrapper))) {
+ return modVol;
+ // build wrapper
+ } else {
+ auto wrp = plm.child(_Unicode(wrapper));
+ auto thickness = wrp.attr<double>(_Unicode(thickness));
+ auto wrpMat = desc.material(wrp.attr<std::string>(_Unicode(material)));
+ Box wrpShape((sx + thickness)/2., (sy + thickness)/2., sz/2.);
+ Volume wrpVol("wrapper_vol", wrpShape, wrpMat);
+ wrpVol.placeVolume(modVol, Position(0., 0., 0.));
+ wrpVol.setVisAttributes(desc.visAttributes(wrp.attr<std::string>(_Unicode(vis))));
+ dim = Position{sx + thickness, sy + thickness, sz};
+ return wrpVol;
+ }
+}
+
+// place array of modules
+static void add_array(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid)
+{
+ Position modSize;
+ auto modVol = build_module(desc, plm, sens, modSize);
+ int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid);
+ int nrow = plm.attr<int>(_Unicode(nrow));
+ int ncol = plm.attr<int>(_Unicode(ncol));
+
+ // compute array position
+ double begx = -modSize.x()*ncol/2. + modSize.x()/2.;
+ double begy = modSize.y()*nrow/2. - modSize.y()/2.;
+
+ std::vector<std::pair<int, int>> removals;
+ // get the removal list
+ for (xml::Collection_t rm(plm, _Unicode(removal)); rm; ++rm) {
+ removals.push_back({rm.attr<int>(_Unicode(row)), rm.attr<int>(_Unicode(col))});
+ }
+
+ // placement to mother
+ auto pos = get_xml_xyz(plm, _Unicode(position));
+ auto rot = get_xml_xyz(plm, _Unicode(rotation));
+ for (int i = 0; i < nrow; ++i) {
+ for (int j = 0; j < ncol; ++j) {
+ if (std::find(removals.begin(), removals.end(), std::pair<int, int>(i, j)) != removals.end()) {
+ continue;
+ }
+ double px = begx + modSize.x()*j;
+ double py = begy - modSize.y()*i;
+ Transform3D tr = Translation3D(pos.x() + px, pos.y() + py, pos.z())
+ * RotationZYX(rot.z(), rot.y(), rot.x());
+ auto modPV = env.placeVolume(modVol, tr);
+ modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", i*ncol + j);
+ }
+ }
+}
+
+// place modules
+static void add_blocks(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid)
+{
+ Position modSize;
+ auto modVol = build_module(desc, plm, sens, modSize);
+ int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid);
+
+ int mid = 1;
+ for (xml::Collection_t pl(plm, _Unicode(placement)); pl; ++pl, ++mid) {
+ Position pos(dd4hep::getAttrOrDefault<double>(pl, _Unicode(x), 0.),
+ dd4hep::getAttrOrDefault<double>(pl, _Unicode(y), 0.),
+ dd4hep::getAttrOrDefault<double>(pl, _Unicode(z), 0.));
+ Position rot(dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotx), 0.),
+ dd4hep::getAttrOrDefault<double>(pl, _Unicode(roty), 0.),
+ dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotz), 0.));
+ Transform3D tr = Translation3D(pos.x(), pos.y(), pos.z())
+ * RotationZYX(rot.z(), rot.y(), rot.x());
+ auto modPV = env.placeVolume(modVol, tr);
+ modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", mid);
+ }
+}
+
+// place disk of modules
+static void add_disk(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid)
+{
+ Position modSize;
+ auto modVol = build_module(desc, plm, sens, modSize);
+ int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid);
+ double rmin = plm.attr<double>(_Unicode(rmin));
+ double rmax = plm.attr<double>(_Unicode(rmax));
+ double phimin = dd4hep::getAttrOrDefault<double>(plm, _Unicode(phimin), 0.);
+ double phimax = dd4hep::getAttrOrDefault<double>(plm, _Unicode(phimax), 2.*M_PI);
+
+ auto points = ref::utils::fillSquares({0., 0.}, modSize.x(), rmin, rmax, phimin, phimax);
+ // placement to mother
+ auto pos = get_xml_xyz(plm, _Unicode(position));
+ auto rot = get_xml_xyz(plm, _Unicode(rotation));
+ int mid = 1;
+ for (auto &p : points) {
+ Transform3D tr = Translation3D(pos.x() + p.x(), pos.y() + p.y(), pos.z())
+ * RotationZYX(rot.z(), rot.y(), rot.x());
+ auto modPV = env.placeVolume(modVol, tr);
+ modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", mid++);
+ }
+}
+//@}
+DECLARE_DETELEMENT(HomogeneousCalorimeter, create_detector)
+
--
GitLab