Skip to content
Snippets Groups Projects
Commit be848cad authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Merge branch 'master' of eicweb.phy.anl.gov:jihee.kim/NPDet into jihee

parents ce8e86d0 8be3ad65
Branches
Tags
No related merge requests found
Showing
with 1080 additions and 21 deletions
...@@ -28,13 +28,16 @@ dd4hep_add_plugin(GenDetectors ...@@ -28,13 +28,16 @@ dd4hep_add_plugin(GenDetectors
trackers/src/SiliconTrackerBarrel_geo.cpp trackers/src/SiliconTrackerBarrel_geo.cpp
trackers/src/SiTrackerRomanPot_geo.cpp trackers/src/SiTrackerRomanPot_geo.cpp
trackers/src/SimpleRomanPot_geo.cpp trackers/src/SimpleRomanPot_geo.cpp
trackers/src/RomanPot_geo.cpp # WIP
calorimeters/src/CylindricalEndcapCalorimeter_geo.cpp calorimeters/src/CylindricalEndcapCalorimeter_geo.cpp
calorimeters/src/EcalBarrel_geo.cpp calorimeters/src/EcalBarrel_geo.cpp
calorimeters/src/PolyhedraEndcapCalorimeter3_geo.cpp calorimeters/src/PolyhedraEndcapCalorimeter3_geo.cpp
calorimeters/src/HexagonalShashlykSamplingECAL_geo.cpp
calorimeters/src/EndcapECAL_geo.cpp calorimeters/src/EndcapECAL_geo.cpp
calorimeters/src/CrystalEndcapECAL_geo.cpp calorimeters/src/CrystalEndcapECAL_geo.cpp
beamline/src/Beampipe_geo.cpp beamline/src/Beampipe_geo.cpp
pid/src/GenericRICH_geo.cpp pid/src/GenericRICH_geo.cpp
pid/src/HexagonalScintPreShower_geo.cpp
USES ActsCore ActsDD4hepPlugin USES ActsCore ActsDD4hepPlugin
NOINSTALL NOINSTALL
) )
......
<lccdd>
<info name="GenericRICH_example"
title="Generic RICH detector example"
author="Whitney Armstrong"
url="https://eicweb.phy.anl.gov/EIC/NPDet"
status="development"
version="1.0"
>
<comment>A generic RICH detector</comment>
</info>
<includes>
<gdmlFile ref="elements.xml"/>
<gdmlFile ref="materials.xml"/>
</includes>
<define>
<constant name="world_side" value="10*m"/>
<constant name="world_x" value="world_side"/>
<constant name="world_y" value="world_side"/>
<constant name="world_z" value="world_side"/>
<constant name="CrossingAngle" value="0.020*rad"/>
<constant name="tracker_region_zmax" value="5*m"/>
<constant name="tracker_region_rmax" value="5*m"/>
</define>
<limits>
<limitset name="cal_limits">
<limit name="step_length_max" particles="*" value="5.0" unit="mm" />
</limitset>
<limitset name="GenericRICHRegionLimitSet">
<limit name="step_length_max" particles="*" value="1.0" unit="mm" />
<limit name="track_length_max" particles="*" value="1.0" unit="mm" />
<limit name="time_max" particles="*" value="0.1" unit="ns" />
<limit name="ekin_min" particles="*" value="0.001" unit="MeV" />
<limit name="range_min" particles="*" value="0.1" unit="mm" />
</limitset>
</limits>
<regions>
<region name="GenericRICHRegion" eunit="MeV" lunit="mm" cut="0.0001" threshold="0.0001">
<limitsetref name="GenericRICHRegionLimitSet"/>
</region>
</regions>
<comment>Common Generic visualization attributes</comment>
<display>
<vis name="InvisibleNoDaughters" showDaughters="false" visible="false"/>
<vis name="InvisibleWithDaughters" showDaughters="true" visible="false"/>
<vis name="GreenVis" alpha="0.5" r= "0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedVis" alpha="1.0" r= "1.0" g="0.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueVis" alpha="1.0" r= "0.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="OrangeVis" alpha="0.5" r= "1.0" g="0.45" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedGreenVis" alpha="0.5" r= "1.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueGreenVis" alpha="0.5" r= "0.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="PurpleVis" alpha="1.0" r= "1.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="DoubleRedG" alpha="0.5" r= "2.0" g=".10" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG015" alpha="0.5" r= "0.0" g=".2" b="1.0" showDaughters="true" visible="true"/>
<vis name="RBG510" alpha="0.5" r= "1.0" g=".2" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG" alpha="0.5" r= "1.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="GrayVis" alpha="0.5" r= "0.75" g="0.75" b="0.75" showDaughters="true" visible="true"/>
</display>
<detectors>
<detector id="1" name="ShashlykECAL" type="HexagonalShashlykSamplingECAL" readout="ShashlykECALHits" vis="RedVis">
<module name="hex_module" side="5.0*cm" layers="194">
<slice name="Pb_slice" thickness="0.5*cm" number="0" material="Lead"/>
<slice name="Scint_slice" thickness="1.5*cm" number="1" material="PlasticScint"/>
<dimensions rmin="10*cm" rmax1="40*cm" rmax2="80*cm" zmin="20*cm" zmax="120*cm"/>
</module>
</detector>
</detectors>
<!-- Definition of the readout segmentation/definition -->
<readouts>
<readout name="ShashlykECALHits">
<id>system:5,sector:5,module:14,sensor:2,side:32:-2,strip:24</id>
</readout>
</readouts>
<plugins>
<!--
<plugin name="DD4hep_GenericSurfaceInstallerPlugin">
<argument value="ForwardRomanPot"/>
<argument value="dimension=2"/>
<argument value="u_x=-1."/>
<argument value="v_y=-1."/>
<argument value="n_z=1."/>
</plugin>
<plugin name="InstallSurfaceManager"/>
-->
</plugins>
<fields>
<field name="GlobalSolenoid" type="solenoid"
inner_field="4.0*tesla"
outer_field="-0.6*tesla"
zmax="3*m"
outer_radius="2*m">
</field>
</fields>
</lccdd>
...@@ -5,6 +5,15 @@ ...@@ -5,6 +5,15 @@
http://www.engineeringtoolbox.com/air-composition-24_212.html http://www.engineeringtoolbox.com/air-composition-24_212.html
--> -->
<material name="PlasticScint">
<D type="density" unit="g/cm3" value="1.032"/>
<composite n="9" ref="C"/>
<composite n="10" ref="H"/>
</material>
<material name="Lead">
<D type="density" value="11.34" unit="g/cm3"/>
<fraction n="1.0" ref="Pb"/>
</material>
<material name="Air"> <material name="Air">
<D type="density" unit="g/cm3" value="0.0012"/> <D type="density" unit="g/cm3" value="0.0012"/>
<fraction n="0.754" ref="N"/> <fraction n="0.754" ref="N"/>
......
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include <XML/Helper.h>
#include "TMath.h"
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
using namespace dd4hep;
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();
//xml::Component dims = detElem.dimensions();
//double rInner = dims.rmin();
//double rOuter1 = dims.rmax1();
//double rOuter2 = dims.rmax2();
//double zMin = dims.zmin();
//double zMax = dims.zmax();
Material air = desc.air();
Material PyrexGlass = desc.material("PyrexGlass");
Assembly assembly( detName+"_assembly" );
//Cone envShape(0.5 * (zMax - zMin), rInner, rOuter1, rInner, rOuter2);
//Volume envVol(detName + "_envelope", envShape, air);
double l_side = 6.25 * cm;
double r_nearest = l_side * std::sqrt(3.0) / 2.0;
double thickness = 2.0 * cm;
double l_side_spacing = 6.25 * cm + 1.0 * mm;
double r_nearest_spacing = l_side_spacing * std::sqrt(3.0) / 2.0;
double y_spacing = 1.5*l_side_spacing;
int nx = 2;
int ny = 2;
int nlayers = 194;
double layer0_thickness = 0.5*cm;
double layer1_thickness = 1.5*cm;
double total_layer_thickness = layer0_thickness + layer1_thickness;
std::vector<double> ptx = { r_nearest, 0.0, -r_nearest, -r_nearest, 0.0, r_nearest };
std::vector<double> pty = { l_side/2.0, l_side, l_side/2.0, -l_side/2.0, -l_side, -l_side/2.0 };
std::vector<double> sec_x = { 0.0, 0.0 };
std::vector<double> sec_y = { 0.0, 0.0 };
std::vector<double> sec_z0 = { -layer0_thickness/2.0, layer0_thickness/2.0 };
std::vector<double> sec_z1 = { -layer1_thickness/2.0, layer1_thickness/2.0 };
std::vector<double> z_scale = { 1.0,1.0};
std::vector<double> sec_r = {r_nearest, r_nearest};
double width_x = (2.0*r_nearest_spacing)*nx + r_nearest_spacing;
double width_y = (y_spacing)*ny + 0.5*l_side_spacing;
double offset_x = -(width_x / 2.0) + r_nearest_spacing;
double offset_y = -(width_y / 2.0) + l_side_spacing;
auto blue_vis = desc.visAttributes("PurpleVis");
auto red_vis = desc.visAttributes("RedVis");
//dd4hep::getAttrOrDefault(desc, _Unicode(vis), "BlueVis")
Assembly module_assembly( "module_assembly" );
//Box Pb_wall_shape(width_x/2.0, width_y/2.0,1.2*cm/2.0);
//ExtrudedPolygon lead_layer_shape("hex_lead",ptx, pty, sec_z0, sec_x, sec_y, z_scale);
Polyhedra lead_layer_shape("hex_lead", 6, 0.0, 2.0*M_PI, sec_z0, sec_r);
Volume lead_layer_Vol("lead_layer_Vol", lead_layer_shape, desc.material("Lead"));
lead_layer_Vol.setVisAttributes(blue_vis);
//ExtrudedPolygon scint_layer_shape("hex_scint",ptx, pty, sec_z1, sec_x, sec_y, z_scale);
Polyhedra scint_layer_shape("hex_scint", 6, 0.0, 2.0*M_PI, sec_z1, sec_r);
Volume scint_layer_Vol("hex_scint_vol", scint_layer_shape, desc.material("PlasticScint"));
scint_layer_Vol.setVisAttributes(red_vis);
for(int ilayer = 0; ilayer < nlayers; ilayer++) {
double z_layer = ilayer * total_layer_thickness + layer0_thickness / 2.0;
auto lead_PV = module_assembly.placeVolume(lead_layer_Vol, Translation3D(0.0, 0.0, z_layer)*RotationZ(M_PI*30.0/180.0));
lead_PV.addPhysVolID("layer", ilayer).addPhysVolID("slice", 1);
double dz_scint = layer0_thickness / 2.0 + layer1_thickness / 2.0;
auto scint_PV =
module_assembly.placeVolume(scint_layer_Vol, Translation3D(0.0, 0.0, z_layer + dz_scint)*RotationZ(M_PI*30.0/180.0));
scint_PV.addPhysVolID("layer", ilayer).addPhysVolID("slice", 2);
}
for(int ix =0; ix<nx; ix++){
for(int iy =0; iy<ny; iy++){
int imod = ix + iy*nx;
double extra_x = -(iy%2)*r_nearest_spacing + r_nearest_spacing;
PlacedVolume module_PV = assembly.placeVolume(
module_assembly, Position(offset_x + ix * (2.0 * r_nearest_spacing) + extra_x,
offset_y + iy * (y_spacing), 0.0));
module_PV.addPhysVolID("module", imod);
}
}
//mirrorVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
//assembly.setVisAttributes(desc.visAttributes(detElem.visStr()));
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
PlacedVolume envPV = motherVol.placeVolume(assembly, Position(0, 0, 0));
envPV.addPhysVolID("system", detID);
det.setPlacement(envPV);
return det;
}
// clang-format off
DECLARE_DETELEMENT(HexagonalShashlykSamplingECAL, createDetector)
...@@ -9,20 +9,75 @@ ...@@ -9,20 +9,75 @@
<comment>A generic RICH detector</comment> <comment>A generic RICH detector</comment>
</info> </info>
<includes>
<gdmlFile ref="elements.xml"/>
<gdmlFile ref="materials.xml"/>
</includes>
<define> <define>
<constant name="world_side" value="10*m"/> <constant name="world_side" value="10*m"/>
<constant name="world_x" value="world_side"/> <constant name="world_x" value="world_side"/>
<constant name="world_y" value="world_side"/> <constant name="world_y" value="world_side"/>
<constant name="world_z" value="world_side"/> <constant name="world_z" value="world_side"/>
<constant name="CrossingAngle" value="0.020*rad"/> <constant name="CrossingAngle" value="0.020*rad"/>
<constant name="tracker_region_zmax" value="5*m"/>
<constant name="tracker_region_rmax" value="5*m"/>
<constant name="PhotMomWaveConv" value="1243.125*eV"/>
</define> </define>
<properties>
<matrix name="RINDEX__N2" coldim="2" values="
1.0*eV 1.00033
4.0*eV 1.00033
5.1*eV 1.00033
"/>
<matrix name="RINDEX__Pyrex" coldim="2" values="
1.0*eV 1.5
4.0*eV 1.5
5.1*eV 1.5
"/>
<matrix name= "REFLECTIVITY_mirror" coldim="2" values="
1.0*eV 0.9
4.0*eV 0.9
5.1*eV 0.9
"/>
</properties>
<includes>
<gdmlFile ref="elements.xml"/>
<gdmlFile ref="materials.xml"/>
</includes>
<materials>
<material name="N2cherenkov">
<D type="density" value="0.00125" unit="g/cm3"/>
<composite n="1" ref="N"/>
<property name="RINDEX" ref="RINDEX__N2"/>
</material>
<material name="PyrexCherenkov">
<D type="density" value="2.23" unit="g/cm3"/>
<fraction n="0.806" ref="SiliconOxide"/>
<fraction n="0.130" ref="BoronOxide"/>
<fraction n="0.040" ref="SodiumOxide"/>
<fraction n="0.023" ref="AluminumOxide"/>
<property name="RINDEX" ref="RINDEX__Pyrex"/>
</material>
</materials>
<surfaces>
<comment> For the values of "finish", model and type, see TGeoOpticalSurface.h !
</comment>
<opticalsurface finish="polished" model="glisur" name="MirrorOpticalSurface" type="dielectric_metal" value="0">
<property name="REFLECTIVITY" ref="REFLECTIVITY_mirror"/>
<property name="RINDEX" coldim="2" values="1.034*eV 1.5 4.136*eV 1.5"/>
<!--<property name="EFFICIENCY" ref="EFFICIENCY0x8b77240"/>-->
</opticalsurface>
<opticalsurface name="mirror2" finish="polished" model="glisur" type="dielectric_dielectric">
<property name="REFLECTIVITY" coldim="2" values="1.034*eV 0.8 4.136*eV 0.9"/>
<property name="EFFICIENCY" coldim="2" values="2.034*eV 0.8 4.136*eV 1.0"/>
<property name="RINDEX" coldim="2" values="1.034*eV 1.5 4.136*eV 1.5"/>
</opticalsurface>
</surfaces>
<limits> <limits>
<limitset name="cal_limits"> <limitset name="cal_limits">
<limit name="step_length_max" particles="*" value="5.0" unit="mm" /> <limit name="step_length_max" particles="*" value="5.0" unit="mm" />
...@@ -60,7 +115,7 @@ ...@@ -60,7 +115,7 @@
</display> </display>
<detectors> <detectors>
<detector id="1" name="ForwardRICH" type="GenericRICH" readout="ForwardRICHHits" vis="RedVis"> <detector id="1" name="ForwardRICH" type="GenericRICH" readout="ForwardRICHHits" vis="RedVis" material="N2cherenkov">
<dimensions rmin="10*cm" rmax1="40*cm" rmax2="80*cm" zmin="20*cm" zmax="120*cm"/> <dimensions rmin="10*cm" rmax1="40*cm" rmax2="80*cm" zmin="20*cm" zmax="120*cm"/>
</detector> </detector>
</detectors> </detectors>
......
<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0"
xmlns:xs="http://www.w3.org/2001/XMLSchema"
xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd">
<info name="GenericRICH_example" title="Generic RICH detector example"
author="Whitney Armstrong"
url="https://eicweb.phy.anl.gov/EIC/NPDet"
status="development"
version="$Id: compact.xml v1.0 2016-12-21$">
<comment>A generic RICH detector</comment>
</info>
<includes>
<gdmlFile ref="elements.xml"/>
<gdmlFile ref="materials.xml"/>
</includes>
<define>
<constant name="world_side" value="10*m"/>
<constant name="world_x" value="world_side"/>
<constant name="world_y" value="world_side"/>
<constant name="world_z" value="world_side"/>
<constant name="CrossingAngle" value="0.020*rad"/>
</define>
<limits>
<limitset name="cal_limits">
<limit name="step_length_max" particles="*" value="5.0" unit="mm" />
</limitset>
<limitset name="GenericRICHRegionLimitSet">
<limit name="step_length_max" particles="*" value="1.0" unit="mm" />
<limit name="track_length_max" particles="*" value="1.0" unit="mm" />
<limit name="time_max" particles="*" value="0.1" unit="ns" />
<limit name="ekin_min" particles="*" value="0.001" unit="MeV" />
<limit name="range_min" particles="*" value="0.1" unit="mm" />
</limitset>
</limits>
<regions>
<region name="GenericRICHRegion" eunit="MeV" lunit="mm" cut="0.0001" threshold="0.0001">
<limitsetref name="GenericRICHRegionLimitSet"/>
</region>
</regions>
<comment>Common Generic visualization attributes</comment>
<display>
<vis name="InvisibleNoDaughters" showDaughters="false" visible="false"/>
<vis name="InvisibleWithDaughters" showDaughters="true" visible="false"/>
<vis name="GreenVis" alpha="0.5" r= "0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedVis" alpha="0.5" r= "1.0" g="0.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueVis" alpha="0.5" r= "0.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="OrangeVis" alpha="0.5" r= "1.0" g="0.45" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedGreenVis" alpha="0.5" r= "1.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueGreenVis" alpha="0.5" r= "0.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="PurpleVis" alpha="0.5" r= "1.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="DoubleRedG" alpha="0.5" r= "2.0" g=".10" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG015" alpha="0.5" r= "0.0" g=".2" b="1.0" showDaughters="true" visible="true"/>
<vis name="RBG510" alpha="0.5" r= "1.0" g=".2" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG" alpha="0.5" r= "1.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="GrayVis" alpha="0.5" r= "0.75" g="0.75" b="0.75" showDaughters="true" visible="true"/>
</display>
<detectors>
<detector id="1" name="ForwardRICH" type="HexagonalScintPreShower" readout="ForwardRICHHits" vis="RedVis">
<dimensions rmin="10*cm" rmax1="40*cm" rmax2="80*cm" zmin="20*cm" zmax="120*cm"/>
</detector>
</detectors>
<!-- Definition of the readout segmentation/definition -->
<readouts>
<readout name="ForwardRICHHits">
<id>system:5,layer:4,module:14,sensor:2,side:32:-2,strip:24</id>
</readout>
</readouts>
<plugins>
<!--
<plugin name="DD4hep_GenericSurfaceInstallerPlugin">
<argument value="ForwardRomanPot"/>
<argument value="dimension=2"/>
<argument value="u_x=-1."/>
<argument value="v_y=-1."/>
<argument value="n_z=1."/>
</plugin>
<plugin name="InstallSurfaceManager"/>
-->
</plugins>
<fields>
<field name="GlobalSolenoid" type="solenoid"
inner_field="4.0*tesla"
outer_field="-0.6*tesla"
zmax="3*m"
outer_radius="2*m">
</field>
</fields>
</lccdd>
/control/verbose 2
/run/initialize
/gps/verbose 2
/gps/particle e-
/gps/number 1
/gps/ene/type Gauss
/gps/ene/mono 2.5 GeV
/gps/ene/sigma 2.0 GeV
/gps/pos/type Volume
/gps/pos/shape Cylinder
/gps/pos/centre 0.0 0.0 0.0 cm
/gps/pos/radius 0.01 cm
/gps/pos/halfz 0.01 cm
/gps/position 0 0 0.0 m
/gps/direction 0 0.2 1.0
#/gps/ang/type iso
#/control/execute vis.mac
/run/beamOn 1
<?xml version="1.0" encoding="UTF-8"?> <?xml version="1.0" encoding="UTF-8"?>
<materials> <materials>
<!-- <!--
Air by weight from Air by weight from
...@@ -91,6 +93,11 @@ ...@@ -91,6 +93,11 @@
<fraction n="0.023" ref="AluminumOxide"/> <fraction n="0.023" ref="AluminumOxide"/>
</material> </material>
<material name="Lead">
<D type="density" value="11.34" unit="g/cm3"/>
<fraction n="1.0" ref="Pb"/>
</material>
<material name="CarbonFiber"> <material name="CarbonFiber">
<D type="density" value="1.5" unit="g/cm3"/> <D type="density" value="1.5" unit="g/cm3"/>
<fraction n="0.65" ref="C"/> <fraction n="0.65" ref="C"/>
......
#!/usr/bin/env python
#
from __future__ import absolute_import, unicode_literals
import os
import time
import logging
import DDG4
from DDG4 import OutputLevel as Output
from g4units import keV, GeV, mm, ns, MeV
#
#
logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.INFO)
logger = logging.getLogger(__name__)
"""
dd4hep simulation example setup using the python configuration
@author M.Frank
@version 1.0
"""
def run():
kernel = DDG4.Kernel()
description = kernel.detectorDescription()
install_dir = os.environ['DD4hepINSTALL']
#kernel.loadGeometry(str("file:" + install_dir + "/DDDetectors/compact/SiD.xml"))
kernel.loadGeometry(str("file:" + "GenericRICH_example.xml"))
DDG4.importConstants(description)
geant4 = DDG4.Geant4(kernel, tracker='Geant4TrackerCombineAction')
geant4.printDetectors()
logger.info("# Configure UI")
geant4.setupUI('qt',vis=True,macro='vis.mac')
logger.info("# Configure G4 magnetic field tracking")
geant4.setupTrackingField()
logger.info("# Setup random generator")
rndm = DDG4.Action(kernel, 'Geant4Random/Random')
rndm.Seed = 987654321
rndm.initialize()
# rndm.showStatus()
logger.info("# Configure Run actions")
run1 = DDG4.RunAction(kernel, 'Geant4TestRunAction/RunInit')
run1.Property_int = 12345
run1.Property_double = -5e15 * keV
run1.Property_string = 'Startrun: Hello_2'
logger.info("%s %s %s", run1.Property_string, str(run1.Property_double), str(run1.Property_int))
run1.enableUI()
kernel.registerGlobalAction(run1)
kernel.runAction().adopt(run1)
logger.info("# Configure Event actions")
prt = DDG4.EventAction(kernel, 'Geant4ParticlePrint/ParticlePrint')
prt.OutputLevel = Output.INFO
prt.OutputType = 3 # Print both: table and tree
kernel.eventAction().adopt(prt)
logger.info("""
Configure I/O
""")
# evt_lcio = geant4.setupLCIOOutput('LcioOutput','CLICSiD_'+time.strftime('%Y-%m-%d_%H-%M'))
# evt_lcio.OutputLevel = Output.ERROR
geant4.setupROOTOutput('RootOutput', 'derp_' + time.strftime('%Y-%m-%d_%H-%M'))
gen = DDG4.GeneratorAction(kernel, "Geant4GeneratorActionInit/GenerationInit")
kernel.generatorAction().adopt(gen)
gen = DDG4.GeneratorAction(kernel, "Geant4GeneratorWrapper/GPS")
gen.Uses = 'G4GeneralParticleSource'
gen.OutputLevel = Output.INFO
gen.Mask = 1
gen.enableUI()
kernel.generatorAction().adopt(gen)
#gen = DDG4.GeneratorAction(kernel, "Geant4IsotropeGenerator/IsotropPi+")
#gen.Mask = 1
#gen.Particle = 'e-'
#gen.Energy = 5 * GeV
#gen.Multiplicity = 2
#gen.Distribution = 'cos(theta)'
#geant4.setupGun('electron','e-',3.0*GeV)
#logger.info("# Install vertex smearing for this interaction")
#gen = DDG4.GeneratorAction(kernel, "Geant4InteractionVertexSmear/SmearPi+")
#gen.Mask = 1
#gen.Offset = (20 * mm, 10 * mm, 10 * mm, 0 * ns)
#gen.Sigma = (4 * mm, 1 * mm, 1 * mm, 0 * ns)
#kernel.generatorAction().adopt(gen)
#logger.info("# Second particle generator: e-")
#gen = DDG4.GeneratorAction(kernel, "Geant4IsotropeGenerator/IsotropE-")
#gen.Mask = 2
#gen.Particle = 'e-'
#gen.Energy = 25 * GeV
#gen.Multiplicity = 3
#gen.Distribution = 'uniform'
#kernel.generatorAction().adopt(gen)
#logger.info(" Install vertex smearing for this interaction")
#gen = DDG4.GeneratorAction(kernel, "Geant4InteractionVertexSmear/SmearE-")
#gen.Mask = 2
#gen.Offset = (-20 * mm, -10 * mm, -10 * mm, 0 * ns)
#gen.Sigma = (12 * mm, 8 * mm, 8 * mm, 0 * ns)
#kernel.generatorAction().adopt(gen)
# ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
#logger.info("# Merge all existing interaction records")
#gen = DDG4.GeneratorAction(kernel, "Geant4InteractionMerger/InteractionMerger")
#gen.OutputLevel = 4 # generator_output_level
#gen.enableUI()
#kernel.generatorAction().adopt(gen)
#logger.info("# Finally generate Geant4 primaries")
#gen = DDG4.GeneratorAction(kernel, "Geant4PrimaryHandler/PrimaryHandler")
#gen.OutputLevel = 4 # generator_output_level
#gen.enableUI()
#kernel.generatorAction().adopt(gen)
#logger.info("# ....and handle the simulation particles.")
#part = DDG4.GeneratorAction(kernel, "Geant4ParticleHandler/ParticleHandler")
#kernel.generatorAction().adopt(part)
## part.SaveProcesses = ['conv','Decay']
#part.SaveProcesses = ['Decay']
#part.MinimalKineticEnergy = 100 * MeV
#part.OutputLevel = 5 # generator_output_level
#part.enableUI()
#user = DDG4.Action(kernel, "Geant4TCUserParticleHandler/UserParticleHandler")
#user.TrackingVolume_Zmax = DDG4.EcalEndcap_zmin
#user.TrackingVolume_Rmax = DDG4.EcalBarrel_rmin
#user.enableUI()
#part.adopt(user)
logger.info("# Setup global filters fur use in sensitive detectors")
#f1 = DDG4.Filter(kernel, 'GeantinoRejectFilter/GeantinoRejector')
#f2 = DDG4.Filter(kernel, 'ParticleRejectFilter/OpticalPhotonRejector')
#f2.particle = 'opticalphoton'
#f3 = DDG4.Filter(kernel, 'ParticleSelectFilter/OpticalPhotonSelector')
#f3.particle = 'opticalphoton'
#f4 = DDG4.Filter(kernel, 'EnergyDepositMinimumCut')
#f4.Cut = 10 * MeV
#f4.enableUI()
#kernel.registerGlobalFilter(f1)
#kernel.registerGlobalFilter(f2)
#kernel.registerGlobalFilter(f3)
#kernel.registerGlobalFilter(f4)
logger.info("# First the tracking detectors")
#seq, act = geant4.setupTracker('SiVertexBarrel')
#seq.adopt(f1)
## seq.adopt(f4)
#act.adopt(f1)
#seq, act = geant4.setupTracker('SiVertexEndcap')
#seq.adopt(f1)
## seq.adopt(f4)
#seq, act = geant4.setupTracker('SiTrackerBarrel')
#seq, act = geant4.setupTracker('SiTrackerEndcap')
#seq, act = geant4.setupTracker('SiTrackerForward')
#logger.info("# Now setup the calorimeters")
#seq, act = geant4.setupCalorimeter('EcalBarrel')
#seq, act = geant4.setupCalorimeter('EcalEndcap')
#seq, act = geant4.setupCalorimeter('HcalBarrel')
#seq, act = geant4.setupCalorimeter('HcalEndcap')
#seq, act = geant4.setupCalorimeter('HcalPlug')
#seq, act = geant4.setupCalorimeter('MuonBarrel')
#seq, act = geant4.setupCalorimeter('MuonEndcap')
#seq, act = geant4.setupCalorimeter('LumiCal')
#seq, act = geant4.setupCalorimeter('BeamCal')
logger.info("# Now build the physics list:")
phys = geant4.setupPhysics('QGSP_BERT')
geant4.addPhysics(str('Geant4PhysicsList/Myphysics'))
ph = DDG4.PhysicsList(kernel, 'Geant4OpticalPhotonPhysics/OpticalPhotonPhys')
ph.VerboseLevel = 2
ph.enableUI()
phys.adopt(ph)
ph = DDG4.PhysicsList(kernel, 'Geant4CerenkovPhysics/CerenkovPhys')
ph.MaxNumPhotonsPerStep = 10
ph.MaxBetaChangePerStep = 10.0
ph.TrackSecondariesFirst = True
ph.VerboseLevel = 2
ph.enableUI()
phys.adopt(ph)
## Add special particle types from specialized physics constructor
#part = geant4.addPhysics('Geant4ExtraParticles/ExtraParticles')
#part.pdgfile = 'checkout/DDG4/examples/particle.tbl'
# Add global range cut
rg = geant4.addPhysics('Geant4DefaultRangeCut/GlobalRangeCut')
rg.RangeCut = 0.7 * mm
phys.dump()
#ui_action = dd4hep.sim.createAction(kernel, "Geant4UIManager/UI")
#ui_action.HaveVIS = True
#ui_action.HaveUI = True
#ui_action.SessionType = qt
#ui_action.SetupUI = macro
#kernel.registerGlobalAction(ui_action)
#ui = geant4.setupUI("qt",vis=True,macro="vis.mac")
#ui.enableUI()
kernel.configure()
kernel.initialize()
# DDG4.setPrintLevel(Output.DEBUG)
kernel.run()
kernel.terminate()
if __name__ == "__main__":
run()
/vis/open OGL 800x800-0+0
/vis/drawVolume
/vis/viewer/set/viewpointThetaPhi 30 30
#/vis/viewer/zoom 30.
#/vis/viewer/set/style wireframe
#/vis/scene/add/axes 0 0 0 1 m
/vis/scene/add/trajectories rich smooth
/vis/modeling/trajectories/create/drawByCharge
/vis/modeling/trajectories/drawByCharge-0/default/setDrawStepPts true
/vis/modeling/trajectories/drawByCharge-0/default/setStepPtsSize 2
/vis/scene/add/hits
/vis/scene/endOfEventAction accumulate 200
/vis/viewer/flush
/run/beamOn 1
/control/execute gps.mac
...@@ -3,6 +3,8 @@ ...@@ -3,6 +3,8 @@
#include <XML/Helper.h> #include <XML/Helper.h>
#include "TMath.h" #include "TMath.h"
#include "DDRec/Surface.h" #include "DDRec/Surface.h"
#include "DD4hep/OpticalSurfaces.h"
#include "DDRec/DetectorData.h" #include "DDRec/DetectorData.h"
using namespace std; using namespace std;
...@@ -17,6 +19,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec ...@@ -17,6 +19,7 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
std::string detName = detElem.nameStr(); std::string detName = detElem.nameStr();
int detID = detElem.id(); int detID = detElem.id();
DetElement det(detName, detID);
xml::Component dims = detElem.dimensions(); xml::Component dims = detElem.dimensions();
double rInner = dims.rmin(); double rInner = dims.rmin();
double rOuter1 = dims.rmax1(); double rOuter1 = dims.rmax1();
...@@ -26,22 +29,34 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec ...@@ -26,22 +29,34 @@ static Ref_t createDetector(Detector& desc, xml::Handle_t handle, SensitiveDetec
Material air = desc.air(); Material air = desc.air();
Material PyrexGlass = desc.material("PyrexGlass"); Material PyrexGlass = desc.material("PyrexGlass");
Material N2cherenkov = desc.material("N2cherenkov");
Cone envShape(0.5 * (zMax - zMin), rInner, rOuter1, rInner, rOuter2); Cone envShape(0.5 * (zMax - zMin), rInner, rOuter1, rInner, rOuter2);
Volume envVol(detName + "_envelope", envShape, air); Volume envVol(detName + "_envelope", envShape, N2cherenkov);
DetElement mirror_DE(det,"Mirror_DE",0);
Tube mirrorShape(rInner+1*cm, rOuter2-2*cm, 6*mm/2.0); Tube mirrorShape(rInner+1*cm, rOuter2-2*cm, 6*mm/2.0);
Volume mirrorVol("RICH_mirror_dummy", mirrorShape, PyrexGlass); Volume mirrorVol("RICH_mirror_dummy", mirrorShape, PyrexGlass);
PlacedVolume mirrorPV = envVol.placeVolume(mirrorVol, Position(0, 0, 0.5 * (zMax - zMin) - 1*cm)); PlacedVolume mirrorPV = envVol.placeVolume(mirrorVol, Position(0, 0, 0.5 * (zMax - zMin) - 1*cm));
mirror_DE.setPlacement(mirrorPV);
envVol.setVisAttributes(desc.visAttributes(detElem.visStr())); envVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det); Volume motherVol = desc.pickMotherVolume(det);
PlacedVolume envPV = motherVol.placeVolume(envVol, Position(0, 0, 0.5 * (zMin + zMax))); PlacedVolume envPV = motherVol.placeVolume(envVol, Position(0, 0, 0.5 * (zMin + zMax)));
envPV.addPhysVolID("system", detID); envPV.addPhysVolID("system", detID);
det.setPlacement(envPV); det.setPlacement(envPV);
OpticalSurfaceManager surfMgr = desc.surfaceManager();
OpticalSurface mirrorSurf = surfMgr.opticalSurface("MirrorOpticalSurface");
//OpticalSurface airSurf = surfMgr.opticalSurface("/world/"+det_name+"#AirSurface");
//BorderSurface mirrorBorder_Surf = BorderSurface(desc, det, "RICHmirror", mirrorSurf, mirrorPV, envPV);
SkinSurface mirrorBorder_Surf(desc,mirror_DE,"RICHmirror", mirrorSurf, mirrorVol);
//BorderSurface bubbleSurf = BorderSurface(description, sdet, "TankBubble", airSurf, bubblePlace, tankPlace);
mirrorBorder_Surf.isValid();
//tankSurf.isValid();
return det; return det;
} }
......
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include <XML/Helper.h>
#include "TMath.h"
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
using namespace dd4hep;
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();
xml::Component dims = detElem.dimensions();
double rInner = dims.rmin();
double rOuter1 = dims.rmax1();
double rOuter2 = dims.rmax2();
double zMin = dims.zmin();
double zMax = dims.zmax();
Material air = desc.air();
Material PyrexGlass = desc.material("PyrexGlass");
Assembly assembly( detName+"_assembly" );
//Cone envShape(0.5 * (zMax - zMin), rInner, rOuter1, rInner, rOuter2);
//Volume envVol(detName + "_envelope", envShape, air);
double l_side = 6.25 * cm;
double r_nearest = l_side * std::sqrt(3.0) / 2.0;
double thickness = 2.0 * cm;
double l_side_spacing = 6.25 * cm + 1.0 * mm;
double r_nearest_spacing = l_side_spacing * std::sqrt(3.0) / 2.0;
double y_spacing = 1.5*l_side_spacing;
std::vector<double> ptx = { r_nearest, 0.0, -r_nearest, -r_nearest, 0.0, r_nearest };
std::vector<double> pty = { l_side/2.0, l_side, l_side/2.0, -l_side/2.0, -l_side, -l_side/2.0 };
std::vector<double> sec_x = { 0.0, 0.0 };
std::vector<double> sec_y = { 0.0, 0.0 };
std::vector<double> sec_z = { -thickness/2.0, thickness/2.0 };
std::vector<double> z_scale = { 1.0, 1.0 };
std::vector<double> sec_r = {r_nearest, r_nearest};
int nx = 10;
int ny = 10;
double width_x = (2.0*r_nearest_spacing)*nx + r_nearest_spacing;
double width_y = (y_spacing)*ny + 0.5*l_side_spacing;
double offset_x = -(width_x / 2.0) + r_nearest_spacing;
double offset_y = -(width_y / 2.0) + l_side_spacing;
Box Pb_wall_shape(width_x/2.0, width_y/2.0,1.2*cm/2.0);
Volume Pb_wall_Vol("Pb_preshower_wall", Pb_wall_shape, desc.material("Lead"));
PlacedVolume Pb_wall_PV = assembly.placeVolume(Pb_wall_Vol, Position(0.0,0.0, -5.0*cm));
//ExtrudedPolygon mirrorShape("hex_scint",ptx, pty, sec_z, sec_x, sec_y, z_scale);
Polyhedra mirrorShape("hex_lead", 6, 0.0, 2.0*M_PI, sec_z, sec_r);
Volume mirrorVol("RICH_mirror_dummy", mirrorShape, PyrexGlass);
for(int ix =0; ix<nx; ix++){
for(int iy =0; iy<ny; iy++){
double extra_x = -(iy%2)*r_nearest_spacing + r_nearest_spacing;
PlacedVolume mirrorPV = assembly.placeVolume(
mirrorVol, Translation3D(offset_x + ix * (2.0 * r_nearest_spacing) + extra_x,
offset_y + iy * (y_spacing),
0.0)*RotationZ(M_PI*30.0/180.0));
}
}
mirrorVol.setVisAttributes(desc.visAttributes(detElem.visStr()));
assembly.setVisAttributes(desc.visAttributes(detElem.visStr()));
DetElement det(detName, detID);
Volume motherVol = desc.pickMotherVolume(det);
PlacedVolume envPV = motherVol.placeVolume(assembly, Position(0, 0, 0.5 * (zMin + zMax)));
envPV.addPhysVolID("system", detID);
det.setPlacement(envPV);
return det;
}
// clang-format off
DECLARE_DETELEMENT(HexagonalScintPreShower, createDetector)
<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0"
xmlns:xs="http://www.w3.org/2001/XMLSchema"
xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd">
<info name="RomanPot_example" title="Roman Pot detector example"
author="Tomas Polakovic"
url="https://eicweb.phy.anl.gov/EIC/NPDet"
status="development"
version="$Id: compact.xml v1.0 2016-12-21$">
<comment>A simplified Roman pot detector</comment>
</info>
<includes>
<gdmlFile ref="elements.xml"/>
<gdmlFile ref="materials.xml"/>
</includes>
<define>
<constant name="world_side" value="10*m"/>
<constant name="world_x" value="world_side"/>
<constant name="world_y" value="world_side"/>
<constant name="world_z" value="world_side"/>
<constant name="CrossingAngle" value="0.020*rad"/>
</define>
<limits>
<limitset name="cal_limits">
<limit name="step_length_max" particles="*" value="5.0" unit="mm" />
</limitset>
<limitset name="SimpleRomanPotRegionLimitSet">
<limit name="step_length_max" particles="*" value="1.0" unit="mm" />
<limit name="track_length_max" particles="*" value="1.0" unit="mm" />
<limit name="time_max" particles="*" value="0.1" unit="ns" />
<limit name="ekin_min" particles="*" value="0.001" unit="MeV" />
<limit name="range_min" particles="*" value="0.1" unit="mm" />
</limitset>
</limits>
<regions>
<region name="SimpleRomanPotRegion" eunit="MeV" lunit="mm" cut="0.0001" threshold="0.0001">
<limitsetref name="SimpleRomanPotRegionLimitSet"/>
</region>
</regions>
<comment>Common Generic visualization attributes</comment>
<display>
<vis name="InvisibleNoDaughters" showDaughters="false" visible="false"/>
<vis name="InvisibleWithDaughters" showDaughters="true" visible="false"/>
<vis name="GreenVis" alpha="0.5" r= "0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedVis" alpha="0.0" r= "1.0" g="0.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueVis" alpha="0.0" r= "0.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="OrangeVis" alpha="0.5" r= "1.0" g="0.45" b="0.0" showDaughters="true" visible="true"/>
<vis name="RedGreenVis" alpha="0.5" r= "1.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="BlueGreenVis" alpha="0.5" r= "0.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="PurpleVis" alpha="0.5" r= "1.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="DoubleRedG" alpha="0.5" r= "2.0" g=".10" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG015" alpha="0.5" r= "0.0" g=".2" b="1.0" showDaughters="true" visible="true"/>
<vis name="RBG510" alpha="0.5" r= "1.0" g=".2" b="0.0" showDaughters="true" visible="true"/>
<vis name="RBG" alpha="0.5" r= "1.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="GrayVis" alpha="0.5" r= "0.75" g="0.75" b="0.75" showDaughters="true" visible="true"/>
</display>
<detectors>
<detector id = "1" name = "MyRomanPot" type = "RomanPot" readout =
"ForwardRomanPotHits" vis = "RedVis">
<dimensions x = "3.0*cm" y = "3.0*cm" delta = "0.005*cm" />
<frame x = "10.0*cm" y = "5.0*cm" z = "2*cm" />
<position z_offset = "0.0*m" rotation = "false" vmax = "10*cm" v = "2.0*cm" />
<layer repeat = "5">
<slice material = "Carbon" thickness = "0.5*mm" vis = "BlueVis" />
<slice material = "Silicon" thickness = "0.03*cm" vis = "GreenVis" sensitive = "true" />
<slice material = "Carbon" thickness = "0.5*mm" vis = "BlueVis" />
<slice material = "Vacuum" thickness = "1.0*mm" vis = "InvisibleWithDaughters" />
</layer>
</detector>
</detectors>
<!-- Definition of the readout segmentation/definition -->
<readouts>
<readout name="ForwardRomanPotHits">
<!-- <id>system:5,element:5,layer:4,module:14,sensor:2,side:32:-2,strip:24</id> -->
<id>system:5,element:5,layer:4,slice:5</id>
</readout>
</readouts>
<plugins>
<!--
<plugin name="DD4hep_GenericSurfaceInstallerPlugin">
<argument value="MyRomanPot"/>
<argument value="dimension=2"/>
<argument value="u_x=-1."/>
<argument value="v_y=-1."/>
<argument value="n_z=1."/>
</plugin>
-->
<plugin name="InstallSurfaceManager"/>
</plugins>
<fields>
<field name="GlobalSolenoid" type="solenoid"
inner_field="4.0*tesla"
outer_field="-0.6*tesla"
zmax="3*m"
outer_radius="2*m">
</field>
</fields>
</lccdd>
#include "DD4hep/DetFactoryHelper.h"
#include "DD4hep/Printout.h"
#include "DD4hep/Shapes.h"
#include "TMath.h"
#include "DDRec/Surface.h"
#include "DDRec/DetectorData.h"
#include "XML/Layering.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
using namespace ROOT::Math;
static Ref_t build_detector(Detector& dtor, xml_h e, SensitiveDetector sens) {
xml_det_t x_det = e;
Layering layering (e);
int det_id = x_det.id();
string det_name = x_det.nameStr();
xml_dim_t dim = x_det.dimensions();
double pixel_x = dim.x();
double pixel_y = dim.y();
double delta = dim.delta();
xml_dim_t frame_dim = x_det.child(_U(frame));
double frame_x = frame_dim.x();
double frame_y = frame_dim.y();
double frame_z = frame_dim.z();
xml_dim_t pos = x_det.position();
double zoffset = pos.z_offset();
bool rotated = pos.attr<bool>(_Unicode(rotation));
double pos_out = pos.vmax();
double curr_pos = pos.v();
Material vacuum = dtor.material("Vacuum");
Material aluminum = dtor.material("Aluminum");
Material frame_mat = dtor.material("Aluminum");
PlacedVolume pv;
DetElement sdet(det_name, det_id);
Assembly assembly(det_name + "_assembly");
sens.setType("tracker");
string module_name = "RomanPot";
string vis0 = dd4hep::getAttrOrDefault(x_det, _Unicode(vis), "BlueVis");
sdet.setAttributes(dtor, assembly, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());
double rp_detector_tube_length = frame_x + pos_out; // TODO: Shortest possible tube given the RP dimensions. Will probably change in the future.
double rp_detector_tube_thickness = 5.0 * dd4hep::mm;
double rp_detector_tube_radius = 0.5*sqrt(frame_z * frame_z + frame_y * frame_y); // TODO: Tightest possible fit. Will probably change in the future.
Tube rp_detector_vacuum_tube(0.0, rp_detector_tube_radius + rp_detector_tube_thickness, rp_detector_tube_length);
Tube rp_detector_vacuum_tube2(0.0, rp_detector_tube_radius, rp_detector_tube_length);
// These are currently hardcoded. Will probably change in the future.
double rp_chamber_thickness = 5.0 * dd4hep::mm;
double rp_chamber_radius = 5.0 * dd4hep::cm;
double rp_chamber_length = rp_detector_tube_radius;
// All the rotation and translation shenanigans related to RP unit orientation happen here.
// `rot1` and `rot2` are used for in construction of the tubes, `shift_pos1` and `shift_pos2`
// properly transform the RP detector positions.
auto rot1 = rotated ? Rotation3D(RotationX(0.5*M_PI)) : Rotation3D(RotationY(0.5*M_PI));
auto rot2 = rotated ? Rotation3D(RotationX(-0.5*M_PI)) : Rotation3D(RotationY(-0.5*M_PI));
auto shift_pos1 = rotated ? Transform3D(RotationZ(0.5*M_PI) * Translation3D(curr_pos + 0.5*frame_x, 0.0, 0.0)) : Transform3D(Translation3D(curr_pos + 0.5*frame_x, 0, 0));
auto shift_pos2 = rotated ? Transform3D(RotationZ(-0.5*M_PI) * Translation3D(curr_pos + 0.5*frame_x, 0.0, 0.0)) : Transform3D(RotationZ(M_PI) * Translation3D(curr_pos + 0.5*frame_x, 0, 0));
auto det_offset = Position(0.5*(pixel_x + delta - frame_x), 0, 0);
// Construct the RP enclosure shell.
Tube rp_beam_pipe_tube(rp_chamber_radius, rp_chamber_radius + rp_chamber_thickness, rp_chamber_length/2.0);
Tube rp_beam_vacuum_tube(0.0, rp_chamber_radius + rp_chamber_thickness, rp_chamber_length);
Tube rp_beam_vacuum_tube2(0.0, rp_chamber_radius, rp_chamber_length);
UnionSolid rp_chamber_tee_outer1(rp_beam_vacuum_tube, rp_detector_vacuum_tube, rot1);
UnionSolid rp_chamber_tee_outer(rp_chamber_tee_outer1, rp_detector_vacuum_tube, rot2);
UnionSolid rp_chamber_tee_inner1(rp_beam_vacuum_tube2, rp_detector_vacuum_tube2, rot1);
UnionSolid rp_chamber_tee_inner(rp_chamber_tee_inner1, rp_detector_vacuum_tube2, rot2);
SubtractionSolid shell(rp_chamber_tee_inner, rp_chamber_tee_outer);
Volume rp_chamber_vol("rp_chamber_walls_vol", shell, aluminum);
Volume rp_vacuum_vol("rp_chamber_vacuum_vol", rp_chamber_tee_inner, vacuum);
pv = assembly.placeVolume(rp_chamber_vol);
auto vacuum_pv = assembly.placeVolume(rp_vacuum_vol);
vacuum_pv.addPhysVolID("element", 1);
rp_chamber_vol.setVisAttributes(vis0.c_str());
rp_vacuum_vol.setVisAttributes(dtor.invisible());
Box frame_box_full(0.5 * frame_x, 0.5 * frame_y, 0.5 * frame_z); // Make the enclosure box for the detector frame
Box det_cutoff(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), 0.5 * frame_z); // Make the cutout where whe detector will sit
SubtractionSolid frame_box(frame_box_full, det_cutoff, det_offset); // The frame with a hole for the detector
// --- d1 ---
Volume frame1_vol("xsensor_frame1", frame_box, frame_mat);
PlacedVolume frame1_pv = rp_vacuum_vol.placeVolume(frame1_vol, shift_pos1);
// Loop over layers
double layer_pos_offset = -0.5*layering.totalThickness();
int layer_num = 1;
for (xml_coll_t c(x_det, _U(layer)); c; ++c) {
xml_comp_t x_layer = c;
int repeat = x_layer.repeat(); // How many times does the layer repeat. Defined in the compact description.
for (int j = 0; j < repeat; j++) {
string layer_name = _toString(layer_num, "layer1_%d");
double layer_thickness = layering.layer(layer_num - 1)->thickness(); // Get the thickness of the current layer.
Position layer_pos = Position(0, 0, layer_pos_offset + 0.5*layer_thickness);
Box layer_box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), layer_thickness);
Volume layer_vol(layer_name, layer_box, vacuum);
DetElement layer(sdet, layer_name, det_id);
// Loop over slices within the layer
double slice_pos_offset = -0.5*layer_thickness;
int slice_num = 1;
for (xml_coll_t si(x_layer, _U(slice)); si; ++si) {
xml_comp_t x_slice = si;
string slice_name = _toString(slice_num, "slice1_%d");
double slice_thickness = x_slice.thickness();
Position slice_pos = Position(0, 0, slice_pos_offset + 0.5 * slice_thickness);
Box slice_box = x_slice.isSensitive() ? Box(0.5 * pixel_x, 0.5 * pixel_y, slice_thickness) : Box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), slice_thickness);
Volume slice_vol(slice_name, slice_box, dtor.material(x_slice.materialStr()));
string slice_vis = dd4hep::getAttrOrDefault(x_slice, _Unicode(vis), "BlueVis");
slice_vol.setVisAttributes(slice_vis.c_str());
DetElement slice(layer, slice_name, det_id);
if (x_slice.isSensitive()) {
slice_vol.setSensitiveDetector(sens);
}
// Place the slice.
PlacedVolume slice_pv = layer_vol.placeVolume(slice_vol, slice_pos);
slice_pv.addPhysVolID("slice", slice_num);
slice.setPlacement(slice_pv);
slice_pos_offset += slice_thickness; // Move the position offset for the next slice.
++slice_num;
}
// Place the layer.
PlacedVolume layer_pv = frame1_vol.placeVolume(layer_vol, layer_pos + det_offset);
layer_pv.addPhysVolID("layer", layer_num);
layer.setPlacement(layer_pv);
layer_pos_offset += layer_thickness;
++layer_num;
}
}
// --- d2 ---
Volume frame2_vol("xsensor_frame2", frame_box, frame_mat);
PlacedVolume frame2_pv = rp_vacuum_vol.placeVolume(frame2_vol, shift_pos2);
// Loop over layers
layer_pos_offset = -0.5*layering.totalThickness();
layer_num = 1;
for (xml_coll_t c(x_det, _U(layer)); c; ++c) {
xml_comp_t x_layer = c;
int repeat = x_layer.repeat(); // How many times does the layer repeat. Defined in the compact description.
for (int j = 0; j < repeat; j++) {
string layer_name = _toString(layer_num, "layer2_%d");
double layer_thickness = layering.layer(layer_num - 1)->thickness(); // Get the thickness of the current layer.
Position layer_pos = Position(0, 0, layer_pos_offset + 0.5*layer_thickness);
Box layer_box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), layer_thickness);
Volume layer_vol(layer_name, layer_box, vacuum);
DetElement layer(sdet, layer_name, det_id);
// Loop over slices within the layer
double slice_pos_offset = -0.5*layer_thickness;
int slice_num = 1;
for (xml_coll_t si(x_layer, _U(slice)); si; ++si) {
xml_comp_t x_slice = si;
string slice_name = _toString(slice_num, "slice2_%d");
double slice_thickness = x_slice.thickness();
Position slice_pos = Position(0, 0, slice_pos_offset + 0.5 * slice_thickness);
Box slice_box = x_slice.isSensitive() ? Box(0.5 * pixel_x, 0.5 * pixel_y, slice_thickness) : Box(0.5 * (pixel_x + delta), 0.5 * (pixel_y + delta), slice_thickness);
Volume slice_vol(slice_name, slice_box, dtor.material(x_slice.materialStr()));
string slice_vis = dd4hep::getAttrOrDefault(x_slice, _Unicode(vis), "BlueVis");
slice_vol.setVisAttributes(slice_vis.c_str());
DetElement slice(layer, slice_name, det_id);
if (x_slice.isSensitive()) {
slice_vol.setSensitiveDetector(sens);
}
// Place the slice.
PlacedVolume slice_pv = layer_vol.placeVolume(slice_vol, slice_pos);
slice_pv.addPhysVolID("slice", slice_num);
slice.setPlacement(slice_pv);
slice_pos_offset += slice_thickness; // Move the position offset for the next slice.
++slice_num;
}
// Place the layer.
PlacedVolume layer_pv = frame2_vol.placeVolume(layer_vol, layer_pos + det_offset);
layer_pv.addPhysVolID("layer", layer_num);
layer.setPlacement(layer_pv);
layer_pos_offset += layer_thickness;
++layer_num;
}
}
pv = dtor.pickMotherVolume(sdet).placeVolume(assembly, Position(0, 0, zoffset));
pv.addPhysVolID("system", det_id);
sdet.setPlacement(pv);
assembly->GetShape()->ComputeBBox();
return sdet;
}
DECLARE_DETELEMENT(RomanPot, build_detector)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment