Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • EIC/detectors/athena
  • zwzhao/athena
  • FernandoTA/athena
  • palspeic/athena
4 results
Show changes
Showing
with 1032 additions and 640 deletions
......@@ -33,34 +33,45 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector s
for(xml_coll_t i(x_det,_U(layer)); i; ++i, ++n) {
xml_comp_t x_layer = i;
string l_name = det_name+_toString(n,"_layer%d");
double z = x_layer.outer_z();
double outer_z = x_layer.outer_z();
double rmin = x_layer.inner_r();
double r = rmin;
DetElement layer(sdet,_toString(n,"layer%d"),x_layer.id());
Tube l_tub (rmin,2*rmin,z);
Tube l_tub (rmin,2*rmin,outer_z);
Volume l_vol(l_name,l_tub,air);
int im = 0;
for(xml_coll_t j(x_layer,_U(slice)); j; ++j, ++im) {
// If slices are only given a thickness attribute, they are radially concentric slices
// If slices are given an inner_z attribute, they are longitudinal slices with equal rmin
xml_comp_t x_slice = j;
Material mat = description.material(x_slice.materialStr());
string s_name= l_name+_toString(im,"_slice%d");
double thickness = x_slice.thickness();
Tube s_tub(r,r+thickness,z,2*M_PI);
double s_outer_z = dd4hep::getAttrOrDefault(x_slice, _Unicode(outer_z), outer_z);
double s_inner_z = dd4hep::getAttrOrDefault(x_slice, _Unicode(inner_z), 0.0*dd4hep::cm);
Tube s_tub(r,r+thickness,(s_inner_z > 0? 0.5*(s_outer_z-s_inner_z): s_outer_z),2*M_PI);
Volume s_vol(s_name, s_tub, mat);
r += thickness;
if ( x_slice.isSensitive() ) {
sens.setType("tracker");
s_vol.setSensitiveDetector(sens);
}
// Set Attributes
s_vol.setAttributes(description,x_slice.regionStr(),x_slice.limitsStr(),x_slice.visStr());
pv = l_vol.placeVolume(s_vol);
if (s_inner_z > 0) {
// Place off-center volumes twice
Position s_pos(0, 0, 0.5*(s_outer_z+s_inner_z));
pv = l_vol.placeVolume(s_vol, -s_pos);
pv = l_vol.placeVolume(s_vol, +s_pos);
} else {
r += thickness;
pv = l_vol.placeVolume(s_vol);
}
// Slices have no extra id. Take the ID of the layer!
pv.addPhysVolID("slice",im);
}
l_tub.setDimensions(rmin,r,z);
l_tub.setDimensions(rmin,r,outer_z);
//cout << l_name << " " << rmin << " " << r << " " << z << endl;
l_vol.setVisAttributes(description,x_layer.visStr());
......@@ -78,4 +89,4 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector s
return sdet;
}
DECLARE_DETELEMENT(refdet_SolenoidCoil,create_detector)
DECLARE_DETELEMENT(athena_SolenoidCoil,create_detector)
/** \addtogroup Trackers Trackers
* \brief Type: **BarrelTrackerWithFrame**.
* \author W. Armstrong
*
* \ingroup trackers
*
* @{
*/
#include <array>
#include <map>
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include "DD4hep/Shapes.h"
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
#include "XML/Utilities.h"
#include "XML/Layering.h"
#if defined(USE_ACTSDD4HEP)
#include "ActsDD4hep/ActsExtension.hpp"
#include "ActsDD4hep/ConvertMaterial.hpp"
#else
#include "Acts/Plugins/DD4hep/ActsExtension.hpp"
#include "Acts/Plugins/DD4hep/ConvertDD4hepMaterial.hpp"
#endif
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
using namespace dd4hep::detail;
/** Endcap Trapezoidal Tracker.
*
* @author Whitney Armstrong
*
*/
static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens)
{
typedef vector<PlacedVolume> Placements;
xml_det_t x_det = e;
Material vacuum = description.vacuum();
int det_id = x_det.id();
string det_name = x_det.nameStr();
bool reflect = x_det.reflect(false);
DetElement sdet(det_name, det_id);
Assembly assembly(det_name);
Material air = description.material("Air");
// Volume assembly (det_name,Box(10000,10000,10000),vacuum);
Volume motherVol = description.pickMotherVolume(sdet);
int m_id = 0, c_id = 0, n_sensor = 0;
map<string, Volume> modules;
map<string, Placements> sensitives;
map<string, std::vector<VolPlane>> volplane_surfaces;
map<string, std::array<double, 2>> module_thicknesses;
PlacedVolume pv;
// ACTS extension
{
Acts::ActsExtension* detWorldExt = new Acts::ActsExtension();
detWorldExt->addType("endcap", "detector");
// SJJ probably need to set the envelope here, as ACTS can't figure
// that out for Assembly volumes. May also need binning to properly pick up
// on the support material @TODO
//
// Add the volume boundary material if configured
for (xml_coll_t bmat(x_det, _Unicode(boundary_material)); bmat; ++bmat) {
xml_comp_t x_boundary_material = bmat;
Acts::xmlToProtoSurfaceMaterial(x_boundary_material, *detWorldExt, "boundary_material");
}
sdet.addExtension<Acts::ActsExtension>(detWorldExt);
}
assembly.setVisAttributes(description.invisible());
sens.setType("tracker");
for (xml_coll_t su(x_det, _U(support)); su; ++su) {
xml_comp_t x_support = su;
double support_thickness = getAttrOrDefault(x_support, _U(thickness), 2.0 * mm);
double support_length = getAttrOrDefault(x_support, _U(length), 2.0 * mm);
double support_rmin = getAttrOrDefault(x_support, _U(rmin), 2.0 * mm);
double support_zstart = getAttrOrDefault(x_support, _U(zstart), 2.0 * mm);
std::string support_name = getAttrOrDefault<std::string>(x_support, _Unicode(name), "support_tube");
std::string support_vis = getAttrOrDefault<std::string>(x_support, _Unicode(vis), "AnlRed");
xml_dim_t pos (x_support.child(_U(position), false));
xml_dim_t rot (x_support.child(_U(rotation), false));
Solid support_solid;
if(x_support.hasChild("shape")){
xml_comp_t shape(x_support.child(_U(shape)));
string shape_type = shape.typeStr();
support_solid = xml::createShape(description, shape_type, shape);
} else {
support_solid = Tube(support_rmin, support_rmin + support_thickness, support_length / 2);
}
Transform3D tr = Transform3D(Rotation3D(),Position(0,0,(reflect?-1.0:1.0) * (support_zstart + support_length / 2)));
if ( pos.ptr() && rot.ptr() ) {
Rotation3D rot3D(RotationZYX(rot.z(0),rot.y(0),rot.x(0)));
Position pos3D(pos.x(0),pos.y(0),pos.z(0));
tr = Transform3D(rot3D, pos3D);
}
else if ( pos.ptr() ) {
tr = Transform3D(Rotation3D(),Position(pos.x(0),pos.y(0),pos.z(0)));
}
else if ( rot.ptr() ) {
Rotation3D rot3D(RotationZYX(rot.z(0),rot.y(0),rot.x(0)));
tr = Transform3D(rot3D,Position());
}
Material support_mat = description.material(x_support.materialStr());
Volume support_vol(support_name, support_solid, support_mat);
support_vol.setVisAttributes(description.visAttributes(support_vis));
pv = assembly.placeVolume(support_vol, tr);
// pv = assembly.placeVolume(support_vol, Position(0, 0, support_zstart + support_length / 2));
}
for (xml_coll_t mi(x_det, _U(module)); mi; ++mi, ++m_id) {
xml_comp_t x_mod = mi;
string m_nam = x_mod.nameStr();
xml_comp_t trd = x_mod.trd();
double posY;
double x1 = trd.x1();
double x2 = trd.x2();
double z = trd.z();
double total_thickness = 0.;
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();
double thickness_so_far = 0.0;
double thickness_sum = -total_thickness / 2.0;
double y1 = total_thickness / 2;
double y2 = total_thickness / 2;
Trapezoid m_solid(x1, x2, y1, y2, z);
Volume m_volume(m_nam, m_solid, vacuum);
m_volume.setVisAttributes(description.visAttributes(x_mod.visStr()));
Solid frame_s;
if(x_mod.hasChild("frame")){
// build frame from trd (assumed to be smaller)
xml_comp_t m_frame = x_mod.child(_U(frame));
xml_comp_t f_pos = m_frame.child(_U(position));
xml_comp_t frame_trd = m_frame.trd();
double frame_thickness = getAttrOrDefault(m_frame, _U(thickness), total_thickness);
double frame_x1 = frame_trd.x1();
double frame_x2 = frame_trd.x2();
double frame_z = frame_trd.z();
// make the frame match the total thickness if thickness attribute is not given
Trapezoid f_solid1(x1, x2,frame_thickness / 2.0, frame_thickness / 2.0, z);
Trapezoid f_solid(frame_x1, frame_x2, frame_thickness / 2.0, frame_thickness / 2.0, frame_z) ;
SubtractionSolid frame_shape(f_solid1, f_solid);
frame_s = frame_shape;
Material f_mat = description.material(m_frame.materialStr());
Volume f_vol(m_nam + "_frame", frame_shape, f_mat);
f_vol.setVisAttributes(description.visAttributes(m_frame.visStr()));
// figure out how to best place
pv = m_volume.placeVolume(f_vol, Position(f_pos.x(), f_pos.y(), f_pos.z()));
}
for (ci.reset(), n_sensor = 1, c_id = 0, posY = -y1; ci; ++ci, ++c_id) {
xml_comp_t c = ci;
double c_thick = c.thickness();
auto comp_x1 = getAttrOrDefault(c, _Unicode(x1), x1);
auto comp_x2 = getAttrOrDefault(c, _Unicode(x2), x2);
auto comp_height = getAttrOrDefault(c, _Unicode(height), z);
Material c_mat = description.material(c.materialStr());
string c_name = _toString(c_id, "component%d");
Trapezoid comp_s1(comp_x1, comp_x2, c_thick / 2e0, c_thick / 2e0, comp_height);
Solid comp_shape = comp_s1;
if(frame_s.isValid()) {
comp_shape = SubtractionSolid( comp_s1, frame_s);
}
Volume c_vol(c_name, comp_shape, c_mat);
c_vol.setVisAttributes(description.visAttributes(c.visStr()));
pv = m_volume.placeVolume(c_vol, Position(0, posY + c_thick / 2, 0));
if (c.isSensitive()) {
module_thicknesses[m_nam] = {thickness_so_far + c_thick/2.0, total_thickness-thickness_so_far - c_thick/2.0};
//std::cout << " adding sensitive volume" << c_name << "\n";
sdet.check(n_sensor > 2, "SiTrackerEndcap2::fromCompact: " + c_name + " Max of 2 modules allowed!");
pv.addPhysVolID("sensor", n_sensor);
sens.setType("tracker");
c_vol.setSensitiveDetector(sens);
sensitives[m_nam].push_back(pv);
++n_sensor;
// -------- create a measurement plane for the tracking surface attched to the sensitive volume -----
Vector3D u(0., 0., -1.);
Vector3D v(-1., 0., 0.);
Vector3D n(0., 1., 0.);
// Vector3D o( 0. , 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];
SurfaceType type(SurfaceType::Sensitive);
// if( isStripDetector )
// type.setProperty( SurfaceType::Measurement1D , true ) ;
VolPlane surf(c_vol, type, inner_thickness, outer_thickness, u, v, n); //,o ) ;
volplane_surfaces[m_nam].push_back(surf);
//--------------------------------------------
}
posY += c_thick;
thickness_sum += c_thick;
thickness_so_far += c_thick;
}
modules[m_nam] = m_volume;
}
for (xml_coll_t li(x_det, _U(layer)); li; ++li) {
xml_comp_t x_layer(li);
int l_id = x_layer.id();
int mod_num = 1;
xml_comp_t l_env = x_layer.child(_U(envelope));
string layer_name = det_name + std::string("_layer") + std::to_string(l_id);
std::string layer_vis = l_env.attr<std::string>(_Unicode(vis));
double layer_rmin = l_env.attr<double>(_Unicode(rmin));
double layer_rmax = l_env.attr<double>(_Unicode(rmax));
double layer_length = l_env.attr<double>(_Unicode(length));
double layer_zstart = l_env.attr<double>(_Unicode(zstart));
double layer_center_z = layer_zstart + layer_length/2.0;
//printout(INFO,"ROOTGDMLParse","+++ Read geometry from GDML file file:%s",input.c_str());
//std::cout << "SiTracker Endcap layer " << l_id << " zstart = " << layer_zstart/dd4hep::mm << "mm ( " << layer_length/dd4hep::mm << " mm thick )\n";
//Assembly layer_assembly(layer_name);
//assembly.placeVolume(layer_assembly);
Tube layer_tub(layer_rmin, layer_rmax, layer_length / 2);
Volume layer_vol(layer_name, layer_tub, air); // Create the layer envelope volume.
layer_vol.setVisAttributes(description.visAttributes(layer_vis));
PlacedVolume layer_pv;
if (reflect) {
layer_pv =
assembly.placeVolume(layer_vol, Transform3D(RotationZYX(0.0, -M_PI, 0.0), Position(0, 0, -layer_center_z)));
layer_pv.addPhysVolID("layer", l_id);
layer_name += "_N";
} else {
layer_pv = assembly.placeVolume(layer_vol, Position(0, 0, layer_center_z));
layer_pv.addPhysVolID("layer", l_id);
layer_name += "_P";
}
DetElement layer_element(sdet, layer_name, l_id);
layer_element.setPlacement(layer_pv);
Acts::ActsExtension* layerExtension = new Acts::ActsExtension();
layerExtension->addType("sensitive disk", "layer");
//layerExtension->addType("axes", "definitions", "XZY");
//layerExtension->addType("sensitive disk", "layer");
//layerExtension->addType("axes", "definitions", "XZY");
for (xml_coll_t lmat(x_layer, _Unicode(layer_material)); lmat; ++lmat) {
xml_comp_t x_layer_material = lmat;
xmlToProtoSurfaceMaterial(x_layer_material, *layerExtension, "layer_material");
}
layer_element.addExtension<Acts::ActsExtension>(layerExtension);
for (xml_coll_t ri(x_layer, _U(ring)); ri; ++ri) {
xml_comp_t x_ring = ri;
double r = x_ring.r();
double phi0 = x_ring.phi0(0);
double zstart = x_ring.zstart();
double dz = x_ring.dz(0);
int nmodules = x_ring.nmodules();
string m_nam = x_ring.moduleStr();
Volume m_vol = modules[m_nam];
double iphi = 2 * M_PI / nmodules;
double phi = phi0;
Placements& sensVols = sensitives[m_nam];
for (int k = 0; k < nmodules; ++k) {
string m_base = _toString(l_id, "layer%d") + _toString(mod_num, "_module%d");
double x = -r * std::cos(phi);
double y = -r * std::sin(phi);
if (!reflect) {
DetElement module(layer_element, m_base + "_pos", det_id);
pv = layer_vol.placeVolume(
m_vol, Transform3D(RotationZYX(0, -M_PI / 2 - phi, -M_PI / 2), Position(x, y, zstart + dz)));
pv.addPhysVolID("module", mod_num);
module.setPlacement(pv);
for (size_t ic = 0; ic < sensVols.size(); ++ic) {
PlacedVolume sens_pv = sensVols[ic];
DetElement comp_elt(module, sens_pv.volume().name(), mod_num);
comp_elt.setPlacement(sens_pv);
//std::cout << " adding ACTS extension" << "\n";
Acts::ActsExtension* moduleExtension = new Acts::ActsExtension("XZY");
comp_elt.addExtension<Acts::ActsExtension>(moduleExtension);
volSurfaceList(comp_elt)->push_back(volplane_surfaces[m_nam][ic]);
}
} else {
pv = layer_vol.placeVolume(
m_vol, Transform3D(RotationZYX(0, -M_PI / 2 - phi, -M_PI / 2), Position(x, y, -zstart - dz)));
pv.addPhysVolID("module", mod_num);
DetElement r_module(layer_element, m_base + "_neg", det_id);
r_module.setPlacement(pv);
for (size_t ic = 0; ic < sensVols.size(); ++ic) {
PlacedVolume sens_pv = sensVols[ic];
DetElement comp_elt(r_module, sens_pv.volume().name(), mod_num);
comp_elt.setPlacement(sens_pv);
//std::cout << " adding ACTS extension" << "\n";
Acts::ActsExtension* moduleExtension = new Acts::ActsExtension("XZY");
comp_elt.addExtension<Acts::ActsExtension>(moduleExtension);
volSurfaceList(comp_elt)->push_back(volplane_surfaces[m_nam][ic]);
}
}
dz = -dz;
phi += iphi;
++mod_num;
}
}
}
pv = motherVol.placeVolume(assembly,Position(0,0,(reflect?-1.0e-9:1.0e-9)) );
pv.addPhysVolID("system", det_id);
sdet.setPlacement(pv);
return sdet;
}
//@}
// clang-format off
DECLARE_DETELEMENT(athena_TrapEndcapTracker, create_detector)
DECLARE_DETELEMENT(athena_GEMTrackerEndcap, create_detector)
#include <XML/Helper.h>
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
#include "DD4hep/OpticalSurfaces.h"
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
//////////////////////////////////
// Central Barrel Tracker Silicon
//////////////////////////////////
using namespace std;
using namespace dd4hep;
static Ref_t createDetector(Detector& desc, xml_h e, SensitiveDetector sens)
{
xml_det_t x_det = e;
string detName = x_det.nameStr();
int detID = x_det.id();
xml_dim_t dim = x_det.dimensions();
double RIn = dim.rmin();
double ROut = dim.rmax();
double SizeZ = dim.length();
double SizeZCut = dim.zmax();
double SiLayerGap = dim.gap();
Material Vacuum = desc.material("Vacuum");
// Create Global Volume
Tube cb_CTD_GVol_Solid(RIn, ROut, SizeZ / 2.0, 0., 360.0 * deg);
Volume detVol("cb_CTD_GVol_Logic", cb_CTD_GVol_Solid, Vacuum);
detVol.setVisAttributes(desc.visAttributes(x_det.visStr()));
// Construct Silicon Layers
xml_comp_t x_layer = x_det.child(_U(layer));
const int repeat = x_layer.repeat();
xml_comp_t x_slice = x_layer.child(_U(slice));
Material slice_mat = desc.material(x_slice.materialStr());
double layerRIn[100];
double layerROut[100];
// Loop over layers
for(int i = 0; i < repeat; i++) {
layerRIn[i] = RIn + (SiLayerGap * i);
layerROut[i] = RIn + (0.01 + SiLayerGap * i);
if (layerROut[i] > ROut)
continue;
string logic_layer_name = detName + _toString(i, "_Logic_lay_%d");
if (i==7){logic_layer_name = detName + _toString(20, "_Logic_lay_%d");}
Volume layerVol(logic_layer_name,Tube(layerRIn[i], layerROut[i], SizeZ / 2.0, 0.0, 360.0 * deg), slice_mat);
layerVol.setVisAttributes(desc,x_layer.visStr());
sens.setType("tracker");
layerVol.setSensitiveDetector(sens);
Position layer_pos = Position(0.0, 0.0, 0.0);
PlacedVolume layerPV = detVol.placeVolume(layerVol, layer_pos);
layerPV.addPhysVolID("layer", i+1);
}
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
Transform3D tr(RotationZYX(0.0, 0.0, 0.0), Position(0.0, 0.0, 0.0));
PlacedVolume detPV = motherVol.placeVolume(detVol, tr);
detPV.addPhysVolID("system", detID);
detPV.addPhysVolID("barrel", 1);
det.setPlacement(detPV);
return det;
}
DECLARE_DETELEMENT(cb_CTD_Si, createDetector)
#include <XML/Helper.h>
#include "TMath.h"
#include "TString.h"
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
#include "DD4hep/OpticalSurfaces.h"
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include "ref_utils.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
void addModules(Volume &mother, xml::DetElement &detElem, Detector &desc, SensitiveDetector &sens);
// create the detector
static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetector sens)
{
xml::DetElement detElem = handle;
std::string detName = detElem.nameStr();
int detID = detElem.id();
DetElement det(detName, detID);
xml::Component dims = detElem.dimensions();
// xml::Component rads = detElem.child(_Unicode(radiator));
auto rmin = dims.rmin();
auto rmax = dims.rmax();
auto length = dims.length();
auto z0 = dims.z();
auto gasMat = desc.material("AirOptical");
// detector envelope
Tube envShape(rmin, rmax, length/2., 0., 2*M_PI);
Volume envVol("ce_MRICH_GVol", envShape, gasMat);
envVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
// modules
addModules(envVol, detElem, desc, sens);
// place envelope
Volume motherVol = desc.pickMotherVolume(det);
PlacedVolume envPV = motherVol.placeVolume(envVol, Position(0, 0, z0));
envPV.addPhysVolID("system", detID);
det.setPlacement(envPV);
return det;
}
void addModules(Volume &mother, xml::DetElement &detElem, Detector &desc, SensitiveDetector &sens)
{
xml::Component dims = detElem.dimensions();
xml::Component mods = detElem.child(_Unicode(modules));
auto rmin = dims.rmin();
auto rmax = dims.rmax();
auto mThick = mods.attr<double>(_Unicode(thickness));
auto mWidth = mods.attr<double>(_Unicode(width));
auto mGap = mods.attr<double>(_Unicode(gap));
auto modMat = desc.material(mods.materialStr());
auto gasMat = desc.material("AirOptical");
// single module
Box mShape(mWidth/2., mWidth/2., mThick/2. - 0.1*mm);
Volume mVol("ce_MRICH_mod_Solid", mShape, modMat);
// a thin gas layer to detect optical photons
Box modShape(mWidth/2., mWidth/2., mThick/2.);
Volume modVol("ce_MRICH_mod_Solid_v", modShape, gasMat);
// thin gas layer is on top (+z) of the material
modVol.placeVolume(mVol, Position(0., 0., -0.1*mm));
modVol.setVisAttributes(desc.visAttributes(mods.visStr()));
sens.setType("photoncounter");
modVol.setSensitiveDetector(sens);
// place modules in the sectors (disk)
auto points = ref::utils::fillSquares({0., 0.}, mWidth + mGap, rmin - mGap, rmax + mGap);
// determine module direction, always facing z = 0
double roty = dims.z() > 0. ? M_PI/2. : -M_PI/2.;
int imod = 1;
for (auto &p : points) {
// operations are inversely ordered
Transform3D tr = Translation3D(p.x(), p.y(), 0.) // move to position
* RotationY(roty); // facing z = 0.
auto modPV = mother.placeVolume(modVol, tr);
modPV.addPhysVolID("sector", 1).addPhysVolID("module", imod ++);
}
}
// clang-format off
DECLARE_DETELEMENT(refdet_ce_MRICH, createDetector)
#include <XML/Helper.h>
///////////////////////////
// Central Ion GEM
///////////////////////////
using namespace dd4hep;
static Ref_t createDetector(Detector& desc, xml_h handle, SensitiveDetector sens)
{
xml::DetElement detElem = handle;
std::string detName = detElem.nameStr();
int detID = detElem.id();
xml::Component dims = detElem.dimensions();
xml::Component pos = detElem.position();
double SizeZ = dims.z_length(); // Size in Z direction
double ROut = dims.rmax(); // Outer radius
double RIn = dims.rmin(); // Inner radius
double ShiftZ = dims.z_offset();
double X = dims.x();
double Z = dims.z();
int Nlayers = dims.number();
double HCAL_rmin = dims.rmax1(); // Maximum radius that the layer can be
Material mat = desc.material(detElem.materialStr());
Material vac = desc.material("Vacuum");
// Outer Volume
Tube ci_GEM_GVol_Solid(RIn, ROut, SizeZ / 2., 0., 360 * deg);
Volume detVol("ci_GEM_GVol_Logic", ci_GEM_GVol_Solid, vac);
// Adding layers to placed volume
for (xml_coll_t li(detElem, _U(layer)); li; ++li) {
xml_comp_t x_layer = li;
std::string layer_name = detName + _toString(x_layer.id(), "_layer%d");
double outer_r = x_layer.outer_r();
if (outer_r > HCAL_rmin) {
outer_r = HCAL_rmin;
}
Volume layer_vol(layer_name, Tube(x_layer.inner_r(), outer_r, x_layer.dz()), mat);
layer_vol.setVisAttributes(desc.visAttributes(detElem.visStr()));
sens.setType("tracker");
layer_vol.setSensitiveDetector(sens);
Position layer_pos(0, 0, x_layer.z());
PlacedVolume layer_phv = detVol.placeVolume(layer_vol, layer_pos);
layer_phv.addPhysVolID("layer", x_layer.id());
}
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
Transform3D tr(RotationZYX(0.0, 0.0, 0.0), Position(pos.x(), pos.x(), pos.z() + SizeZ / 2.0));
PlacedVolume detPV = motherVol.placeVolume(detVol, tr);
detPV.addPhysVolID("system", detID);
det.setPlacement(detPV);
return det;
}
// clang-format off
DECLARE_DETELEMENT(ci_GEM, createDetector)
#include <XML/Helper.h>
///////////////////////////
// Central Ion GEM
///////////////////////////
using namespace dd4hep;
static Ref_t createDetector(Detector& desc, xml_h handle, SensitiveDetector sens)
{
xml::DetElement detElem = handle;
std::string detName = detElem.nameStr();
int detID = detElem.id();
xml::Component dims = detElem.dimensions();
double RIn = dims.rmin();
double ROut = dims.rmax();
double SizeZ = dims.z_length();
double ShiftZ = dims.z_offset();
double PosZ = dims.z();
double lay_RIn = dims.rmin1();
double lay_ROut = dims.rmax1();
double lay_dz = dims.dz();
Material mat_iron = desc.material("Iron");
Material mat_vac = desc.material("Vacuum");
//Outer volume
Tube ci_Hcal_Solid(RIn, ROut, SizeZ / 2., 0., 360 * dd4hep::deg);
Volume envelopeVol("ci_Hcal_Logic", ci_Hcal_Solid, mat_vac);
//Iron tube for the layers
Tube ci_Hcal_detSolid(lay_RIn, lay_ROut, lay_dz / 2., 0., 360 * dd4hep::deg);
//Adding layers to placed detector volume
for (xml_coll_t li(detElem,_U(layer)); li; ++li){
xml_comp_t x_layer = li;
std::string layer_name = detName + _toString(x_layer.id(), "_layer%d");
Volume layer_vol(layer_name, ci_Hcal_detSolid, mat_iron);
layer_vol.setVisAttributes(detElem.visStr());
sens.setType("calorimeter");
layer_vol.setSensitiveDetector(sens);
Position layer_pos(0, 0, x_layer.z());
PlacedVolume layer_phv = envelopeVol.placeVolume(layer_vol, layer_pos);
layer_phv.addPhysVolID("layer", x_layer.id());
}
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
Transform3D tr(RotationZYX(0,0,0), Position(0, 0, ShiftZ));
PlacedVolume detPV = motherVol.placeVolume(envelopeVol, tr);
detPV.addPhysVolID("system", detID);
det.setPlacement(detPV);
return det;
}
// clang-format off
DECLARE_DETELEMENT(ci_HCAL, createDetector)
......@@ -44,7 +44,7 @@ static Ref_t createDetector(Detector& desc, xml_h e, SensitiveDetector sens)
// Loop over layers
for(int i = 0; i < repeat; i++) {
layerRIn = RIn + 1.0 * cm + ((double)i * 0.5) * cm;
layerROut = ROut - 25.0 * cm + ((double)i * 2.0) * cm;
layerROut = ROut;//RIn + ((double)i * 0.5) * cm;
layerPosZ = SizeZ / 2.0 - 5.0 * cm - ((double)i * 3.0) * cm;
layerSizeZ = 1.0 * cm;
......
//==========================================================================
// AIDA Detector description implementation
//--------------------------------------------------------------------------
// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
// All rights reserved.
//
// For the licensing terms see $DD4hepINSTALL/LICENSE.
// For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
//
// Author : M.Frank
//
//==========================================================================
//
// Specialized generic detector constructor
//
//==========================================================================
//
// Implementation of the Sci Fiber geometry: M. Żurek 07/19/2021
#include "DD4hep/DetFactoryHelper.h"
#include "XML/Layering.h"
#include "Math/Point2D.h"
#include "TGeoPolygon.h"
#include "TMath.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::detail;
typedef ROOT::Math::XYPoint Point;
// Fill fiber lattice into trapezoid starting from position (0,0) in x-z coordinate system
vector<vector<Point>> fiberPositions(double radius, double x_spacing, double z_spacing, double x, double z, double phi, double spacing_tol = 1e-2) {
// z_spacing - distance between fiber layers in z
// x_spacing - distance between fiber centers in x
// x - half-length of the shorter (bottom) base of the trapezoid
// z - height of the trapezoid
// phi - angle between z and trapezoid arm
vector<vector<Point>> positions;
int z_layers = floor((z/2-radius-spacing_tol)/z_spacing); // number of layers that fit in z/2
double z_pos = 0.;
double x_pos = 0.;
for(int l = -z_layers; l < z_layers+1; l++) {
vector<Point> xline;
z_pos = l*z_spacing;
double x_max = x + (z/2. + z_pos)*tan(phi) - spacing_tol; // calculate max x at particular z_pos
(l % 2 == 0) ? x_pos = 0. : x_pos = x_spacing/2; // account for spacing/2 shift
while(x_pos < (x_max - radius)) {
xline.push_back(Point(x_pos, z_pos));
if(x_pos != 0.) xline.push_back(Point(-x_pos, z_pos)); // using symmetry around x=0
x_pos += x_spacing;
}
// Sort fiber IDs for a better organization
sort(xline.begin(), xline.end(), [](const Point &p1, const Point &p2) {
return p1.x() < p2.x();
});
positions.emplace_back(std::move(xline));
}
return positions;
}
// Calculate number of divisions for the readout grid for the fiber layers
std::pair<int, int> getNdivisions(double x, double z, double dx, double dz) {
// x and z defined as in vector<Point> fiberPositions
// dx, dz - size of the grid in x and z we want to get close to with the polygons
// See also descripltion when the function is called
double SiPMsize = 13.0*mm;
double grid_min = SiPMsize + 3.0*mm;
if(dz < grid_min) {
dz = grid_min;
}
if(dx < grid_min) {
dx = grid_min;
}
int nfit_cells_z = floor(z/dz);
int n_cells_z = nfit_cells_z;
if(nfit_cells_z == 0) n_cells_z++;
int nfit_cells_x = floor((2*x)/dx);
int n_cells_x = nfit_cells_x;
if(nfit_cells_x == 0) n_cells_x++;
return std::make_pair(n_cells_x, n_cells_z);
}
// Calculate dimensions of the polygonal grid in the cartesian coordinate system x-z
vector< tuple<int, Point, Point, Point, Point> > gridPoints(int div_x, int div_z, double x, double z, double phi) {
// x, z and phi defined as in vector<Point> fiberPositions
// div_x, div_z - number of divisions in x and z
double dz = z/div_z;
std::vector<std::tuple<int, Point, Point, Point, Point>> points;
for(int iz = 0; iz < div_z + 1; iz++){
for(int ix = 0; ix < div_x + 1; ix++){
double A_z = -z/2 + iz*dz;
double B_z = -z/2 + (iz+1)*dz;
double len_x_for_z = 2*(x+iz*dz*tan(phi));
double len_x_for_z_plus_1 = 2*(x + (iz+1)*dz*tan(phi));
double dx_for_z = len_x_for_z/div_x;
double dx_for_z_plus_1 = len_x_for_z_plus_1/div_x;
double A_x = -len_x_for_z/2. + ix*dx_for_z;
double B_x = -len_x_for_z_plus_1/2. + ix*dx_for_z_plus_1;
double C_z = B_z;
double D_z = A_z;
double C_x = B_x + dx_for_z_plus_1;
double D_x = A_x + dx_for_z;
int id = ix + div_x * iz;
auto A = Point(A_x, A_z);
auto B = Point(B_x, B_z);
auto C = Point(C_x, C_z);
auto D = Point(D_x, D_z);
// vertex points filled in the clock-wise direction
points.push_back(make_tuple(id, A, B, C, D));
}
}
return points;
}
// Create detector
static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) {
static double tolerance = 0e0;
Layering layering (e);
xml_det_t x_det = e;
Material air = description.air();
int det_id = x_det.id();
string det_name = x_det.nameStr();
xml_comp_t x_staves = x_det.staves();
xml_comp_t x_dim = x_det.dimensions();
int nsides = x_dim.numsides();
double inner_r = x_dim.rmin();
double dphi = (2*M_PI/nsides);
double hphi = dphi/2;
double support_thickness = 0.0;
if(x_staves.hasChild("support")){
support_thickness = getAttrOrDefault(x_staves.child(_U(support)), _U(thickness), 5.0 * cm);
}
double mod_z = layering.totalThickness() + support_thickness;
double outer_r = inner_r + mod_z;
double totThick = mod_z;
double offset = x_det.attr<double>(_Unicode(offset));
DetElement sdet (det_name,det_id);
Volume motherVol = description.pickMotherVolume(sdet);
PolyhedraRegular hedra (nsides,inner_r,inner_r+totThick+tolerance*2e0,x_dim.z());
Volume envelope (det_name,hedra,air);
PlacedVolume env_phv = motherVol.placeVolume(envelope,Transform3D(Translation3D(0,0,offset)*RotationZ(M_PI/nsides)));
env_phv.addPhysVolID("system",det_id);
sdet.setPlacement(env_phv);
DetElement stave_det("stave0",det_id);
double dx = 0.0; //mod_z / std::sin(dphi); // dx per layer
// Compute the top and bottom face measurements.
double trd_x2 = (2 * std::tan(hphi) * outer_r - dx)/2 - tolerance;
double trd_x1 = (2 * std::tan(hphi) * inner_r + dx)/2 - tolerance;
double trd_y1 = x_dim.z()/2 - tolerance;
double trd_y2 = trd_y1;
double trd_z = mod_z/2 - tolerance;
// Create the trapezoid for the stave.
Trapezoid trd(trd_x1, // Outer side, i.e. the "long" X side.
trd_x2, // Inner side, i.e. the "short" X side.
trd_y1, // Corresponds to subdetector (or module) Z.
trd_y2, //
trd_z); // Thickness, in Y for top stave, when rotated.
Volume mod_vol("stave",trd,air);
double l_pos_z = -(layering.totalThickness() / 2) - support_thickness/2.0;
//double trd_x2_support = trd_x1;
double trd_x1_support = (2 * std::tan(hphi) * outer_r - dx- support_thickness)/2 - tolerance;
Solid support_frame_s;
// optional stave support
if(x_staves.hasChild("support")){
xml_comp_t x_support = x_staves.child(_U(support));
// is the support on the inside surface?
bool is_inside_support = getAttrOrDefault<bool>(x_support, _Unicode(inside), true);
// number of "beams" running the length of the stave.
int n_beams = getAttrOrDefault<int>(x_support, _Unicode(n_beams), 3);
double beam_thickness = support_thickness / 4.0; // maybe a parameter later...
trd_x1_support = (2 * std::tan(hphi) * (outer_r - support_thickness + beam_thickness)) / 2 - tolerance;
double grid_size = getAttrOrDefault(x_support, _Unicode(grid_size), 25.0 * cm);
double beam_width = 2.0 * trd_x1_support / (n_beams + 1); // quick hack to make some gap between T beams
double cross_beam_thickness = support_thickness/4.0;
//double trd_x1_support = (2 * std::tan(hphi) * (inner_r + beam_thickness)) / 2 - tolerance;
double trd_x2_support = trd_x2;
int n_cross_supports = std::floor((trd_y1-cross_beam_thickness)/grid_size);
Box beam_vert_s(beam_thickness / 2.0 - tolerance, trd_y1, support_thickness / 2.0 - tolerance);
Box beam_hori_s(beam_width / 2.0 - tolerance, trd_y1, beam_thickness / 2.0 - tolerance);
UnionSolid T_beam_s(beam_vert_s, beam_hori_s, Position(0, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
// cross supports
Trapezoid trd_support(trd_x1_support,trd_x2_support,
beam_thickness / 2.0 - tolerance, beam_thickness / 2.0 - tolerance,
support_thickness / 2.0 - tolerance - cross_beam_thickness/2.0);
UnionSolid support_array_start_s(T_beam_s,trd_support,Position(0,0,cross_beam_thickness/2.0));
for (int isup = 0; isup < n_cross_supports; isup++) {
support_array_start_s = UnionSolid(support_array_start_s, trd_support, Position(0, -1.0 * isup * grid_size, cross_beam_thickness/2.0));
support_array_start_s = UnionSolid(support_array_start_s, trd_support, Position(0, 1.0 * isup * grid_size, cross_beam_thickness/2.0));
}
support_array_start_s =
UnionSolid(support_array_start_s, beam_hori_s,
Position(-1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
support_array_start_s =
UnionSolid(support_array_start_s, beam_hori_s,
Position(1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, -support_thickness / 2.0 + beam_thickness / 2.0));
support_array_start_s =
UnionSolid(support_array_start_s, beam_vert_s, Position(-1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, 0));
support_array_start_s =
UnionSolid(support_array_start_s, beam_vert_s, Position(1.8 * 0.5*(trd_x1+trd_x2_support) / n_beams, 0, 0));
support_frame_s = support_array_start_s;
Material support_mat = description.material(x_support.materialStr());
Volume support_vol("support_frame_v", support_frame_s, support_mat);
support_vol.setVisAttributes(description,x_support.visStr());
// figure out how to best place
//auto pv = mod_vol.placeVolume(support_vol, Position(0.0, 0.0, l_pos_z + support_thickness / 2.0));
auto pv = mod_vol.placeVolume(support_vol, Position(0.0, 0.0, -l_pos_z - support_thickness / 2.0));
}
//l_pos_z += support_thickness;
sens.setType("calorimeter");
{ // ===== buildBarrelStave(description, sens, module_volume) =====
// Parameters for computing the layer X dimension:
double stave_z = trd_y1;
double tan_hphi = std::tan(hphi);
double l_dim_x = trd_x1; // Starting X dimension for the layer.
// Loop over the sets of layer elements in the detector.
int l_num = 1;
for(xml_coll_t li(x_det,_U(layer)); li; ++li) {
xml_comp_t x_layer = li;
int repeat = x_layer.repeat();
// Loop over number of repeats for this layer.
for (int j=0; j<repeat; j++) {
string l_name = _toString(l_num,"layer%d");
double l_thickness = layering.layer(l_num-1)->thickness(); // Layer's thickness.
Position l_pos(0,0,l_pos_z+l_thickness/2); // Position of the layer.
double l_trd_x1 = l_dim_x - tolerance;
double l_trd_x2 = l_dim_x + l_thickness*tan_hphi - tolerance;
double l_trd_y1 = stave_z-tolerance;
double l_trd_y2 = l_trd_y1;
double l_trd_z = l_thickness/2-tolerance;
Trapezoid l_trd(l_trd_x1,l_trd_x2,l_trd_y1,l_trd_y2,l_trd_z);
Volume l_vol(l_name,l_trd,air);
DetElement layer(stave_det, l_name, det_id);
// Loop over the sublayers or slices for this layer.
int s_num = 1;
double s_pos_z = -(l_thickness / 2);
for(xml_coll_t si(x_layer,_U(slice)); si; ++si) {
xml_comp_t x_slice = si;
string s_name = _toString(s_num,"slice%d");
double s_thick = x_slice.thickness();
Volume s_vol(s_name);
DetElement slice(layer,s_name,det_id);
double s_trd_x1 = l_dim_x + (s_pos_z+l_thickness/2)*tan_hphi - tolerance;
double s_trd_x2 = l_dim_x + (s_pos_z+l_thickness/2+s_thick)*tan_hphi - tolerance;
double s_trd_y1 = stave_z-tolerance;
double s_trd_y2 = s_trd_y1;
double s_trd_z = s_thick/2-tolerance;
Trapezoid s_trd(s_trd_x1, s_trd_x2, s_trd_y1, s_trd_y2, s_trd_z);
s_vol.setSolid(s_trd);
s_vol.setMaterial(description.material(x_slice.materialStr()));
if (x_slice.hasChild("fiber")) {
xml_comp_t x_fiber = x_slice.child(_Unicode(fiber));
double f_radius = getAttrOrDefault(x_fiber, _U(radius), 0.1 * cm);
double f_spacing_x = getAttrOrDefault(x_fiber, _Unicode(spacing_x), 0.122 * cm);
double f_spacing_z = getAttrOrDefault(x_fiber, _Unicode(spacing_z), 0.134 * cm);
std::string f_id_grid = getAttrOrDefault(x_fiber, _Unicode(identifier_grid), "grid");
std::string f_id_fiber = getAttrOrDefault(x_fiber, _Unicode(identifier_fiber), "fiber");
// Calculate fiber positions inside the slice
Tube f_tube(0, f_radius, stave_z-tolerance);
// Set up the readout grid for the fiber layers
// Trapezoid is divided into segments with equal dz and equal number of divisions in x
// Every segment is a polygon that can be attached later to the lightguide
// The grid size is assumed to be ~2x2 cm (starting values). This is to be larger than
// SiPM chip (for GlueX 13mmx13mm: 4x4 grid 3mmx3mm with 3600 50×50 μm pixels each)
// See, e.g., https://arxiv.org/abs/1801.03088 Fig. 2d
// Calculate number of divisions
pair<int, int> grid_div = getNdivisions(s_trd_x1, s_thick-tolerance, 2.0*cm, 2.0*cm);
// Calculate polygonal grid coordinates (vertices)
vector<tuple<int, Point, Point, Point, Point>> grid_vtx = gridPoints(grid_div.first, grid_div.second, s_trd_x1, s_thick-tolerance, hphi);
vector<int> f_id_count(grid_div.first*grid_div.second,0);
auto f_pos = fiberPositions(f_radius, f_spacing_x, f_spacing_z, s_trd_x1, s_thick-tolerance, hphi);
for (auto &line : f_pos) {
for (auto &p : line) {
int f_grid_id = -1;
int f_id = -1;
// Check to which grid fiber belongs to
for (auto &poly_vtx : grid_vtx) {
auto [grid_id, vtx_a, vtx_b, vtx_c, vtx_d] = poly_vtx;
double poly_x[4] = {vtx_a.x(), vtx_b.x(), vtx_c.x(), vtx_d.x()};
double poly_y[4] = {vtx_a.y(), vtx_b.y(), vtx_c.y(), vtx_d.y()};
double f_xy[2] = {p.x(), p.y()};
TGeoPolygon poly(4);
poly.SetXY(poly_x,poly_y);
poly.FinishPolygon();
if(poly.Contains(f_xy)) {
f_grid_id = grid_id;
f_id = f_id_count[grid_id];
f_id_count[grid_id]++;
}
}
string f_name = "fiber" + to_string(f_grid_id) + "_" + to_string(f_id);
Volume f_vol(f_name, f_tube, description.material(x_fiber.materialStr()));
DetElement fiber(slice, f_name, det_id);
if ( x_fiber.isSensitive() ) {
f_vol.setSensitiveDetector(sens);
}
fiber.setAttributes(description,f_vol,x_fiber.regionStr(),x_fiber.limitsStr(),x_fiber.visStr());
// Fiber placement
Transform3D f_tr(RotationZYX(0,0,M_PI*0.5),Position(p.x(), 0 ,p.y()));
PlacedVolume fiber_phv = s_vol.placeVolume(f_vol, f_tr);
fiber_phv.addPhysVolID(f_id_grid, f_grid_id + 1).addPhysVolID(f_id_fiber, f_id + 1);
fiber.setPlacement(fiber_phv);
}
}
}
if ( x_slice.isSensitive() ) {
s_vol.setSensitiveDetector(sens);
}
slice.setAttributes(description,s_vol,x_slice.regionStr(),x_slice.limitsStr(),x_slice.visStr());
// Slice placement.
PlacedVolume slice_phv = l_vol.placeVolume(s_vol,Position(0,0,s_pos_z+s_thick/2));
slice_phv.addPhysVolID("slice", s_num);
slice.setPlacement(slice_phv);
// Increment Z position of slice.
s_pos_z += s_thick;
// Increment slice number.
++s_num;
}
// Set region, limitset, and vis of layer.
layer.setAttributes(description,l_vol,x_layer.regionStr(),x_layer.limitsStr(),x_layer.visStr());
PlacedVolume layer_phv = mod_vol.placeVolume(l_vol,l_pos);
layer_phv.addPhysVolID("layer", l_num);
layer.setPlacement(layer_phv);
// Increment to next layer Z position.
double xcut = l_thickness * tan_hphi;
l_dim_x += xcut;
l_pos_z += l_thickness;
++l_num;
}
}
}
// Set stave visualization.
if ( x_staves ) {
mod_vol.setVisAttributes(description.visAttributes(x_staves.visStr()));
}
// Phi start for a stave.
double phi = M_PI / nsides;
double mod_x_off = dx / 2; // Stave X offset, derived from the dx.
double mod_y_off = inner_r + mod_z/2; // Stave Y offset
// Create nsides staves.
for (int i = 0; i < nsides; i++, phi -= dphi) { // i is module number
// Compute the stave position
double m_pos_x = mod_x_off * std::cos(phi) - mod_y_off * std::sin(phi);
double m_pos_y = mod_x_off * std::sin(phi) + mod_y_off * std::cos(phi);
Transform3D tr(RotationZYX(0,phi,M_PI*0.5),Translation3D(-m_pos_x,-m_pos_y,0));
PlacedVolume pv = envelope.placeVolume(mod_vol,tr);
pv.addPhysVolID("system",det_id);
pv.addPhysVolID("module",i+1);
DetElement sd = i==0 ? stave_det : stave_det.clone(_toString(i,"stave%d"));
sd.setPlacement(pv);
sdet.add(sd);
}
// Set envelope volume attributes.
envelope.setAttributes(description,x_det.regionStr(),x_det.limitsStr(),x_det.visStr());
return sdet;
}
DECLARE_DETELEMENT(athena_EcalBarrelHybrid,create_detector)
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
#include "DD4hep/OpticalSurfaces.h"
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include <XML/Helper.h>
///////////////////////////////////////////
// Far Forward Ion Zero Degree Calorimeter
///////////////////////////////////////////
using namespace std;
using namespace dd4hep;
static Ref_t createDetector(Detector& desc, xml_h e, SensitiveDetector sens)
{
xml_det_t x_det = e;
string detName = x_det.nameStr();
int detID = x_det.id();
xml_dim_t dim = x_det.dimensions();
double Width = dim.x();
double Thickness = dim.z();
xml_dim_t pos = x_det.position();
xml_dim_t rot = x_det.rotation();
Material Vacuum = desc.material("Vacuum");
xml_comp_t mod = x_det.child(_Unicode(module));
string modName = mod.nameStr();
Material mPbWO4 = desc.material(mod.materialStr());
double mThickness = mod.attr<double>(_Unicode(thickness));
double mRmin = mod.attr<double>(_Unicode(rmin));
double mRmax = mod.attr<double>(_Unicode(rmax));
double mWidth = mod.attr<double>(_Unicode(width));
double mGap = mod.attr<double>(_Unicode(gap));
int mNTowers = mod.attr<double>(_Unicode(ntower));
// Create Global Volume
Box ffi_ZDC_GVol_Solid(Width * 0.5, Width * 0.5, Thickness * 0.5);
Volume detVol("ffi_ZDC_GVol_Logic", ffi_ZDC_GVol_Solid, Vacuum);
detVol.setVisAttributes(desc.visAttributes(x_det.visStr()));
// Construct Tower
// Single Module
Box ffi_ZDC_ECAL_Solid_Tower(mWidth * 0.5, mWidth * 0.5, mThickness * 0.5);
Volume modVol("ffi_ZDC_ECAL_Logic_Tower", ffi_ZDC_ECAL_Solid_Tower, mPbWO4);
modVol.setVisAttributes(desc.visAttributes(mod.visStr()));
sens.setType("calorimeter");
modVol.setSensitiveDetector(sens);
// Module Position
double mod_x = 0.0 * mm;
double mod_y = 0.0 * mm;
double mod_z = -1.0 * Thickness / 2.0 + mThickness / 2.0 + 2.0 * mm;
int k = -1;
// Place Modules
for (int j = 0; j < mNTowers; j++) {
if (j == 0)
mod_y = Width / 2.0 - mWidth / 2.0 - mGap;
else
mod_y -= (mWidth + mGap);
if (abs(mod_y + mWidth / 2.0) > Width / 2.0)
continue;
mod_x = Width / 2.0 - (mWidth + mGap) * 0.5;
for (int i = 0; i < mNTowers; i++) {
if (i > 0)
mod_x -= (mWidth + mGap);
if (abs(mod_x + mWidth / 2.0) > Width / 2.0)
continue;
k++;
string module_name = detName + _toString(k,"_ECAL_Phys_%d");
PlacedVolume pv_mod = detVol.placeVolume(modVol, Position(mod_x,mod_y,mod_z));
pv_mod.addPhysVolID("sector", 1).addPhysVolID("module",k+1);
}
}
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
Transform3D tr(RotationZYX(rot.z(), -rot.y(), rot.x()), Position(pos.x(), pos.y(), pos.z()));
PlacedVolume detPV = motherVol.placeVolume(detVol, tr);
detPV.addPhysVolID("system", detID);
det.setPlacement(detPV);
return det;
}
DECLARE_DETELEMENT(ffi_ZDC, createDetector)
//==========================================================================
// AIDA Detector description implementation
//--------------------------------------------------------------------------
// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
// All rights reserved.
//
// For the licensing terms see $DD4hepINSTALL/LICENSE.
// For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
//
// Author : M.Frank
//
//==========================================================================
//
// Specialized generic detector constructor
//
//==========================================================================
#include "DD4hep/DetFactoryHelper.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::detail;
static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector sens) {
xml_det_t x_det = e;
Material air = description.air();
string det_name = x_det.nameStr();
bool reflect = x_det.reflect();
DetElement sdet(det_name,x_det.id());
Assembly assembly(det_name);
PlacedVolume pv;
int l_num = 0;
xml::Component pos = x_det.position();
for(xml_coll_t i(x_det,_U(layer)); i; ++i, ++l_num) {
xml_comp_t x_layer = i;
string l_nam = det_name+_toString(l_num,"_layer%d");
double zmin = x_layer.inner_z();
double rmin = x_layer.inner_r();
double rmax = x_layer.outer_r();
double z = zmin, layerWidth = 0.;
int s_num = 0;
for(xml_coll_t j(x_layer,_U(slice)); j; ++j) {
double thickness = xml_comp_t(j).thickness();
layerWidth += thickness;
}
Tube l_tub(rmin,rmax,layerWidth,2*M_PI);
Volume l_vol(l_nam,l_tub,air);
l_vol.setVisAttributes(description,x_layer.visStr());
for(xml_coll_t j(x_layer,_U(slice)); j; ++j, ++s_num) {
xml_comp_t x_slice = j;
double thick = x_slice.thickness();
Material mat = description.material(x_slice.materialStr());
string s_nam = l_nam+_toString(s_num,"_slice%d");
Volume s_vol(s_nam, Tube(rmin,rmax,thick), mat);
if ( x_slice.isSensitive() ) {
sens.setType("tracker");
s_vol.setSensitiveDetector(sens);
}
s_vol.setAttributes(description,x_slice.regionStr(),x_slice.limitsStr(),x_slice.visStr());
pv = l_vol.placeVolume(s_vol,Position(0,0,z-zmin-layerWidth/2+thick/2));
pv.addPhysVolID("slice",s_num);
}
DetElement layer(sdet,l_nam+"_pos",l_num);
pv = assembly.placeVolume(l_vol,Position(0,0,zmin+layerWidth/2.));
pv.addPhysVolID("layer",l_num);
pv.addPhysVolID("barrel",1);
layer.setPlacement(pv);
if ( reflect ) {
pv = assembly.placeVolume(l_vol,Transform3D(RotationY(M_PI),Position(0,0,-zmin-layerWidth/2)));
pv.addPhysVolID("layer",l_num);
pv.addPhysVolID("barrel",2);
DetElement layerR = layer.clone(l_nam+"_neg");
sdet.add(layerR.setPlacement(pv));
}
}
if ( x_det.hasAttr(_U(combineHits)) ) {
sdet.setCombineHits(x_det.attr<bool>(_U(combineHits)),sens);
}
pv = description.pickMotherVolume(sdet).placeVolume(assembly,Position(pos.x(),pos.y(),pos.z()));
pv.addPhysVolID("system", x_det.id()); // Set the subdetector system ID.
sdet.setPlacement(pv);
return sdet;
}
DECLARE_DETELEMENT(ref_DiskTracker,create_detector)
DECLARE_DETELEMENT(ref_SolenoidEndcap,create_detector)
#include "ref_utils.h"
// some utility functions that can be shared
namespace ref::utils {
typedef ROOT::Math::XYPoint Point;
// check if a square in a ring
inline bool in_ring(const Point &pt, double side, double rmin, double rmax, double phmin, double phmax)
{
if (pt.r() > rmax || pt.r() < rmin) {
return false;
}
// check four corners
std::vector<Point> pts {
Point(pt.x() - side/2., pt.y() - side/2.),
Point(pt.x() - side/2., pt.y() + side/2.),
Point(pt.x() + side/2., pt.y() - side/2.),
Point(pt.x() + side/2., pt.y() + side/2.),
};
for (auto &p : pts) {
if (p.r() > rmax || p.r() < rmin || p.phi() > phmax || p.phi() < phmin) {
return false;
}
}
return true;
}
// check if a square is overlapped with the others
inline bool overlap(const Point &pt, double side, const std::vector<Point> &pts)
{
for (auto &p : pts) {
auto pn = (p - pt)/side;
if ((std::abs(pn.x()) < 1. - 1e-6) && (std::abs(pn.y()) < 1. - 1e-6)) {
return true;
}
}
return false;
}
// a helper function to recursively fill square in a ring
void add_square(Point p, std::vector<Point> &res, double lside, double rmin, double rmax,
double phmin, double phmax)
{
// outside of the ring or overlapping
if (!in_ring(p, lside, rmin, rmax, phmin, phmax) || overlap(p, lside, res)) {
return;
}
res.emplace_back(p);
// check adjacent squares
add_square(Point(p.x() + lside, p.y()), res, lside, rmin, rmax, phmin, phmax);
add_square(Point(p.x() - lside, p.y()), res, lside, rmin, rmax, phmin, phmax);
add_square(Point(p.x(), p.y() + lside), res, lside, rmin, rmax, phmin, phmax);
add_square(Point(p.x(), p.y() - lside), res, lside, rmin, rmax, phmin, phmax);
}
// fill squares
std::vector<Point> fillSquares(Point ref, double lside, double rmin, double rmax, double phmin, double phmax)
{
// 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);
auto find_seed = [] (const Point &ref, int n, double side, double rmin, double rmax, double phmin, double phmax) {
for (int ix = -n; ix < n; ++ix) {
for (int iy = -n; iy < n; ++iy) {
Point pt(ref.x() + ix*side, ref.y() + iy*side);
if (in_ring(pt, side, rmin, rmax, phmin, phmax)) {
return pt;
}
}
}
return ref;
};
std::vector<Point> res;
ref = find_seed(ref, int(rmax/lside) + 2, lside, rmin, rmax, phmin, phmax);
add_square(ref, res, lside, rmin, rmax, phmin, phmax);
return res;
}
} // ref::utils
#pragma once
#include <vector>
#include "Math/Point2D.h"
// some utility functions that can be shared
namespace ref::utils {
typedef ROOT::Math::XYPoint Point;
// fill squares in a ring
std::vector<Point> fillSquares(Point ref, double lside, double rmin, double rmax,
double phmin = -M_PI, double phmax = M_PI);
} // ref::utils
#dawn_view_03:detector:
#dawn_view_13:detector:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view03 -d scripts/view3 -D
#dawn_view_03:ev001:
# - ./bin/make_dawn_views -t view13 -d scripts/view13 -D
#dawn_view_14:detector:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view03_ev001 -d scripts/view3 -s 2
#dawn_view_03:ev003:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view03_ev100 -d scripts/view3 -s 3
#dawn_view_03:ev004:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view03_ev001 -d scripts/view3 -s 4
#dawn_view_06:detector:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view06 -d scripts/view6 -D
#dawn_view_07:detector:
# extends: .views
# script:
# - ./bin/make_dawn_views -t view07 -d scripts/view7 -D
dawn_view_11:detector:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11 -D
dawn_view_11:ev000:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11
dawn_view_11:ev001:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11 -s 1
dawn_view_11:ev002:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11 -s 2
dawn_view_11:ev003:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11 -s 3
dawn_view_11:ev004:
extends: .views
script:
- ./bin/make_dawn_views -t view11 -d scripts/view11 -s 4
dawn_view_12:detector:
extends: .views
script:
- ./bin/make_dawn_views -t view12 -d scripts/view12 -D
dawn_view_13:detector:
extends: .views
script:
- ./bin/make_dawn_views -t view13 -d scripts/view13 -D
dawn_view_14:detector:
extends: .views
script:
- ./bin/make_dawn_views -t view14 -d scripts/view14 -D
# - ./bin/make_dawn_views -t view14 -d scripts/view14 -D
view_prim:detector_only:
extends: .views
stage: test
script:
- ./bin/generate_prim_file -o ${LOCAL_DATA_PATH} -D -t detector_view
- ls -lrth && ls -lrth ${LOCAL_DATA_PATH}
view_prim:ev001:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- ./bin/generate_prim_file -o ${LOCAL_DATA_PATH} -t view_ev001 -s 1
view_prim:ev002:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- ./bin/generate_prim_file -o ${LOCAL_DATA_PATH} -t view_ev002 -s 2
view_prim:ev003:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- ./bin/generate_prim_file -o ${LOCAL_DATA_PATH} -t view_ev003 -s 3
view_prim:ev004:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- ./bin/generate_prim_file -o ${LOCAL_DATA_PATH} -t view_ev004 -s 4
view_prim:calorimeters:
extends: .views
stage: test
script:
- cp "compact/subsystem_views/calorimeters.xml" "${DETECTOR_PATH}/."
- ./bin/generate_prim_file -c ${DETECTOR_PATH}/calorimeters.xml -o ${LOCAL_DATA_PATH} -D -t calorimeters_view
- ls -lrth && ls -lrth ${LOCAL_DATA_PATH}
view_prim:calorimeters_ev001:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- cp "compact/subsystem_views/calorimeters.xml" "${DETECTOR_PATH}/."
- ./bin/generate_prim_file -c ${DETECTOR_PATH}/calorimeters.xml -o ${LOCAL_DATA_PATH} -t calorimeters_view_ev001 -s 1
view_prim:calorimeters_ev002:
extends: .views
stage: test
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
script:
- cp "compact/subsystem_views/calorimeters.xml" "${DETECTOR_PATH}/."
- ./bin/generate_prim_file -c ${DETECTOR_PATH}/calorimeters.xml -o ${LOCAL_DATA_PATH} -t calorimeters_view_ev002 -s 2
dawn_view_01:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -t view01 -d scripts/view1 -D
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view01 -d scripts/view1 -D
dawn_view_01:ev001:
extends: .views
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
needs:
- job: view_prim:ev001
optional: true
script:
- ./bin/make_dawn_views -t view01_ev001 -d scripts/view1 -s 1
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/view_ev001.prim -t view01_ev001 -d scripts/view1 -s 1
dawn_view_01:ev002:
extends: .views
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
needs:
- job: view_prim:ev002
optional: true
script:
- ./bin/make_dawn_views -t view01_ev001 -d scripts/view1 -s 2
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/view_ev002.prim -t view01_ev002 -d scripts/view1 -s 2
view_01:
stage: test
stage: collect
needs:
- ["dawn_view_01:detector", "dawn_view_01:ev001","dawn_view_01:ev002"]
- job: dawn_view_01:detector
optional: false
- job: dawn_view_01:ev001
optional: true
- job: dawn_view_01:ev002
optional: true
script:
- ls -lrth *
- ls -lrth images/*
......
dawn_view_11:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view11 -d scripts/view11 -D
dawn_view_11:ev000:
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
extends: .views
needs:
- job: view_prim:ev001
optional: true
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/view_ev001.prim -t view11 -d scripts/view11
dawn_view_11:ev001:
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
extends: .views
needs:
- job: view_prim:ev001
optional: true
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/view_ev002.prim -t view11 -d scripts/view11 -s 1
dawn_view_11:ev002:
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
extends: .views
needs:
- job: view_prim:ev002
optional: true
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/view_ev003.prim -t view11 -d scripts/view11 -s 2
dawn_view_11:ev003:
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
extends: .views
needs:
- job: view_prim:ev003
optional: true
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH} -t view11 -d scripts/view11 -s 3
dawn_view_11:ev004:
rules:
- if: '$DETECTOR_EVENT_VIEWS == "ON"'
extends: .views
needs:
- job: view_prim:ev004
optional: true
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH} -t view11 -d scripts/view11 -s 4
view_11:
stage: collect
rules:
- if: '$CI_SERVER_HOST == "eicweb.phy.anl.gov"'
needs:
- job: dawn_view_11:detector
optional: false
- job: dawn_view_11:ev001
optional: true
- job: dawn_view_11:ev002
optional: true
- job: dawn_view_11:ev003
optional: true
- job: dawn_view_11:ev004
optional: true
script:
- ls -lrth *
- ls -lrth images/*
dawn_view_12:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view12 -d scripts/view12 -D -- ${SLICE}
- ls -lrth *
- ls -lrth images/*
parallel:
matrix:
- SLICE: ["100", "300", "500", "700", "900", "1100", "1300", "1500", "1700", "1900"]
view_12:
stage: collect
rules:
- if: '$CI_SERVER_HOST == "eicweb.phy.anl.gov"'
needs:
- ["dawn_view_12:detector"]
script:
- ls -lrth *
- ls -lrth images/*
dawn_view_13:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view13 -d scripts/view13 -D
view_13:
stage: collect
rules:
- if: '$CI_SERVER_HOST == "eicweb.phy.anl.gov"'
needs:
- ["dawn_view_13:detector"]
script:
- ls -lrth *
- ls -lrth images/*
dawn_view_14:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view14 -d scripts/view14 -D -- ${SLICE}
- ls -lrth *
- ls -lrth images/*
parallel:
matrix:
- SLICE: ["100", "300", "500", "700", "900", "1100", "1300", "1500", "1700", "1900"]
view_14:
stage: collect
rules:
- if: '$CI_SERVER_HOST == "eicweb.phy.anl.gov"'
needs:
- ["dawn_view_14:detector"]
script:
- ls -lrth *
- ls -lrth images/*
dawn_view_15:detector:
extends: .views
needs:
- job: view_prim:detector_only
optional: false
script:
- ./bin/make_dawn_views -i ${LOCAL_DATA_PATH}/detector_view.prim -t view15 -d scripts/view15 -D -- ${SLICE}
- ls -lrth *
- ls -lrth images/*
parallel:
matrix:
- SLICE: ["100", "300", "500", "700", "900", "1100", "1300", "1500", "1700", "1900"]
view_15:
stage: collect
rules:
- if: '$CI_SERVER_HOST == "eicweb.phy.anl.gov"'
needs:
- ["dawn_view_15:detector"]
script:
- ls -lrth *
- ls -lrth images/*