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

Resolve "Detailed MRICH"

parent b9b25150
No related branches found
No related tags found
1 merge request!65Resolve "Detailed MRICH"
......@@ -4,11 +4,14 @@
<!-- Some information about detector -->
<info name="Athena Detector" title="Athena Detector"
author="Jihee Kim, Whitney Armstrong"
author="Athena Collaboration"
url="https://eicweb.phy.anl.gov/EIC/detectors/athena.git"
status="development"
version="v1 2021-03-16">
<comment>Athena Detector</comment>
<comment>Athena Detector
- https://eicweb.phy.anl.gov/EIC/detectors/athena.git
- https://eicweb.phy.anl.gov/EIC/detectors/ip6.git
</comment>
</info>
<define>
......@@ -39,11 +42,31 @@
4.0*eV 1.5
5.1*eV 1.5
"/>
<matrix name="ABSLENGTH__Pyrex" coldim="2" values="
1.0*eV 10.0*cm
4.0*eV 10.0*cm
5.1*eV 10.0*cm
"/>
<matrix name= "REFLECTIVITY_mirror" coldim="2" values="
1.0*eV 0.9
4.0*eV 0.9
5.1*eV 0.9
"/>
"/>
<matrix name="RINDEX__Aerogel" coldim="2" values="
1.0*eV 1.030
4.0*eV 1.030
5.1*eV 1.030
"/>
<matrix name="ABSLENGTH__Aerogel" coldim="2" values="
1.0*eV 4.0*cm
4.0*eV 4.0*cm
5.1*eV 4.0*cm
"/>
<matrix name="RINDEX__Acrylic" coldim="2" values="
1240*eV/1100 1.49
1240*eV/600 1.49
1240*eV/400 1.49
"/>
</properties>
<includes>
......@@ -57,37 +80,87 @@
<fraction n="0.754" ref="N"/>
<fraction n="0.234" ref="O"/>
<fraction n="0.012" ref="Ar"/>
<property name="RINDEX" ref="RINDEX__Vacuum"/>
<property name="RINDEX" ref="RINDEX__Air"/>
<property name="ABSLENGTH" coldim="2" values="1*eV 200*m 5*eV 200*m"/>
</material>
<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="PyrexGlass">
<material name="PyrexGlassOptical">
<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"/>
<property name="ABSLENGTH" ref="ABSLENGTH__Pyrex"/>
</material>
<material name="AerogelOptical">
<D value="0.2" unit="g / cm3"/>
<fraction n="0.625" ref="SiliconOxide"/>
<fraction n="0.374" ref="SiliconOxide"/>
<fraction n="0.1" ref="C"/>
<property name="RINDEX" ref="RINDEX__Aerogel"/>
<property name="ABSLENGTH" ref="ABSLENGTH__Aerogel"/>
</material>
<material name="AcrylicOptical">
<D type="density" value="1.18" unit="g/cm3"/>
<composite n="5" ref="C"/>
<composite n="2" ref="O"/>
<composite n="8" ref="H"/>
<property name="RINDEX" ref="RINDEX__Acrylic"/>
</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>
<!-- Define the dimensions of the world volume -->
-->
<opticalsurface finish="polished" model="unified" name="MRICH_MirrorOpticalSurface" type="dielectric_metal" value="0">
</opticalsurface>
<opticalsurface finish="polished" model="unified" name="MRICH_LensOpticalSurface" type="dielectric_dielectric" value="0">
<property name="REFLECTIVITY" coldim="2" values="1240*eV/1100 0.08 1240*eV/400 0.08"/>
<!--
<property name="RINDEX" coldim="2" values="2.034*eV 1.56 4.136*eV 1.56"/>
<property name="SPECULARLOBECONSTANT" coldim="2" values="2.034*eV 0.3 4.136*eV 0.3 "/>
<property name="SPECULARSPIKECONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
<property name="BACKSCATTERCONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
-->
</opticalsurface>
<opticalsurface finish="polished" model="unified" name="MRICH_PhotoSensorOpticalSurface" type="dielectric_dielectric" value="0">
<!--
<property name="RINDEX" coldim="2" values="2.034*eV 1.35 4.136*eV 1.40"/>
<property name="SPECULARLOBECONSTANT" coldim="2" values="2.034*eV 0.3 4.136*eV 0.3 "/>
<property name="SPECULARSPIKECONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
<property name="BACKSCATTERCONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
-->
</opticalsurface>
<opticalsurface finish="polished" model="unified" name="MRICH_AerogelOpticalSurface" type="dielectric_dielectric" value="0">
<!--
<property name="RINDEX" coldim="2" values="2.034*eV 1.010 4.136*eV 1.010"/>
<property name="SPECULARLOBECONSTANT" coldim="2" values="2.034*eV 0.3 4.136*eV 0.3 "/>
<property name="SPECULARSPIKECONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
<property name="BACKSCATTERCONSTANT" coldim="2" values="2.034*eV 0.2 4.136*eV 0.2 "/>
-->
</opticalsurface>
</surfaces>
<limits>
......@@ -189,7 +262,9 @@
<include ref="compact/tof_barrel.xml"/>
<!--include ref="compact/rwell_tracker_barrel.xml"/-->
<!--
<include ref="compact/rwell_tracker_barrel.xml"/>
-->
<include ref="compact/cb_DIRC.xml"/>
<!-- When changing magnet, also select dimensions in definitions.xml. -->
......@@ -206,7 +281,7 @@
<include ref="compact/hcal.xml"/>
<!--include ref="compact/ce_GEM.xml"/-->
<!--include ref="compact/gem_tracker_endcap.xml"/-->
<include ref="compact/ce_mrich.xml"/>
<include ref="compact/mrich.xml"/>
<include ref="compact/tof_endcap.xml"/>
<include ref="compact/forward_trd.xml"/>
<include ref="compact/gaseous_rich.xml"/>
......@@ -216,12 +291,12 @@
<include ref="ip6/far_forward_romanpots.xml"/>
<include ref="ip6/far_forward_detectors.xml"/>
<!--
<include ref="compact/mm_tracker_barrel.xml"/>
<include ref="compact/cb_VTX_Barrel.xml"/>
-->
<readouts>
</readouts>
......
#!/bin/bash
npsim --enableG4GPS \
--compact athena.xml \
--runType vis \
--macro macro/mrich_vis.mac \
-O derp.root
<?xml version="1.0" encoding="UTF-8"?>
<lccdd>
<define>
<constant name="ce_MRICHRMin" value="15*cm"/>
<constant name="ce_MRICHRMax" value="BarrelTracking_rmax"/>
<constant name="ce_MRICHLength" value="BackwardCherenkov_length"/>
<constant name="ce_MRICHZMin" value="-BarrelTracking_length/2.0 - BackwardTracking_length"/>
</define>
<detectors>
<detector id="ce_MRICH_ID" name="ce_MRICH" type="refdet_ce_MRICH" readout="MRICHHits" vis="BlueVis" material="Air">
<dimensions rmin="ce_MRICHRMin" rmax="ce_MRICHRMax" length="ce_MRICHLength" z="ce_MRICHZMin"/>
<modules material="Quartz" thickness="10*cm" width="10*cm" gap="1*cm" />
</detector>
</detectors>
<readouts>
<readout name="MRICHHits">
<segmentation type="CartesianGridXY" grid_size_x="3*mm" grid_size_y="3*mm" />
<id>system:8,sector:8,module:14,x:32:-16,y:-16</id>
</readout>
</readouts>
</lccdd>
......@@ -57,7 +57,7 @@
<constant name="TrackerEndcapInnerLayer_dz" value="TrackerEndcapInner_length/TrackerEndcapInner_nLayers"/>
<constant name="TrackerEndcapOuter_zmin" value="TrackerBarrelOuter_length/2.0"/>
<constant name="TrackerEndcapOuter_zmax" value="TrackerEndcapOuter_zmin + 72"/>
<constant name="TrackerEndcapOuter_zmax" value="TrackerEndcapOuter_zmin + 70"/>
<constant name="TrackerEndcapOuter_length" value="TrackerEndcapOuter_zmax- TrackerEndcapOuter_zmin"/>
<constant name="TrackerEndcapOuterLayer_dz" value="TrackerEndcapOuter_length/TrackerEndcapOuter_nLayers"/>
......
......@@ -5,8 +5,8 @@
<constant name="world_z" value="100*m"/>
<constant name="Pi" value="3.14159265359"/>
<constant name="mil" value="0.0254*mm"/>
<constant name="inch" value="2.54*cm"/>
<comment>
------------
......@@ -232,7 +232,7 @@
Unused IDs: 132-139
</comment>
<constant name="ce_MRICH_ID" value="130"/>
<constant name="MRICH_ID" value="130"/>
<constant name="ce_GEM_ID" value="131"/>
<comment>
......@@ -340,7 +340,7 @@
<constant name="ForwardTRD_length" value="10.0*cm"/>
<constant name="ForwardTOF_length" value="3.0*cm"/>
<constant name="BackwardCherenkov_length" value="15.0*cm"/>
<constant name="BackwardCherenkov_length" value="20.0*cm"/>
<constant name="BackwardTOF_length" value="3.0*cm"/>
<comment> Total length of PID detectors above</comment>
......
......@@ -6,9 +6,11 @@
<vis name="AnlOrange" alpha="1" r="255/256" g="121/256" b="0"/>
<vis name="AnlRed" alpha="1" r="205/256" g="32/256" b="44/256"/>
<vis name="AnlGold" alpha="1" r="248/256" g="178/256" b="0/256"/>
<vis name="AnlGold_1" alpha="1" r="248/256" g="178/256" b="0/256" visible="true" drawingStyle="wireframe"/>
<vis name="AnlBlue" alpha="1" r="0/256" g="96/256" b="156/256"/>
<vis name="AnlTeal" alpha="1" r="0/256" g="161/256" b="156/256"/>
<vis name="AnlGray" alpha="1" r="102/256" g="102/256" b="102/256"/>
<vis name="AnlGray_1" alpha="1" r="102/256" g="102/256" b="102/256" visible="true" drawingStyle="wireframe"/>
<vis name="AnlLight_Gray" alpha="1" r="209/256" g="209/256" b="209/256"/>
<vis name="AnlOff_White" alpha="1" r="242/256" g="242/256" b="242/256"/>
<vis name="AnlDelta_Red" alpha="1" r="161/256" g="43/256" b="47/256"/>
......
<lccdd>
<comment> MRICH </comment>
<define>
<constant name="MRICH_rmin" value="15*cm"/>
<constant name="MRICH_rmax" value="BarrelTracking_rmax"/>
<constant name="MRICH_length" value="BackwardCherenkov_length"/>
<constant name="MRICH_zmin" value="-BarrelTracking_length/2.0 - BackwardTracking_length"/>
<constant name="MRICHAerogel_thickness" value="30.0*mm"/>
<constant name="MRICHAerogel_width" value="126.5*mm"/>
<constant name="MRICHFoam_thickness" value="2*mm"/>
<constant name="MRICHFresnelLens_thickness" value="0.06*inch"/>
<constant name="MRICHAerogelLensGap_thickness" value="2*mm"/>
<constant name="MRICHLensMirrorGap_thickness" value="2*mm"/>
<constant name="MRICHPhotoDet_thickness" value="1.5*mm"/>
<constant name="MRICHPhotoDet_length" value="48.5*mm"/>
<constant name="MRICHGlassWindow_width" value="103.5*mm"/>
<constant name="MRICHGlassPhotoDet_thickness" value="2.0*mm"/>
<constant name="MRICHRearExtraSpace_thickness" value="10.0*mm"/>
<constant name="MRICHLensPhotoDet_length" value="136.4*mm"/>
<constant name="MRICHMirror_thickness" value="2.0*mm"/>
<constant name="MRICHMirror_length" value="MRICHLensPhotoDet_length - MRICHLensMirrorGap_thickness"/>
<constant name="MRICHMirror_width1" value="MRICHAerogel_width "/>
<constant name="MRICHMirror_width2" value="MRICHGlassWindow_width"/>
<constant name="MRICHFresnelLensEffectiveDiameter" value="6.0*inch"/>
<constant name="MRICHFresnelLensGroove_pitch" value="inch/125"/>
<constant name="MRICHCarbonFrame_thickness" value="1.0*mm"/>
<constant name="MRICHCarbonFrame_width" value="MRICHAerogel_width+2.0*MRICHFoam_thickness + 2.0*MRICHCarbonFrame_thickness"/>
<constant name="MRICHModules_nx" value="floor((MRICH_rmax-MRICH_rmin)/MRICHCarbonFrame_width)"/>
<constant name="MRICHModules_ny" value="floor((MRICH_rmax-MRICH_rmin)/MRICHCarbonFrame_width)"/>
<constant name="MRICHCarbonFrame_length"
value="MRICHAerogel_thickness
+ 2.0*MRICHCarbonFrame_thickness
+ 2.0*MRICHFoam_thickness
+ MRICHAerogelLensGap_thickness
+ MRICHFresnelLens_thickness
+ MRICHLensPhotoDet_length
+ MRICHGlassPhotoDet_thickness
+ MRICHRearExtraSpace_thickness "/>
</define>
<limits>
</limits>
<regions>
</regions>
<display>
</display>
<detectors>
<detector id="MRICH_ID" name="MRICH" type="athena_MRICH"
readout="MRICHHits"
projective="false"
vis="InvisibleWithDaughters" material="Air">
<dimensions rmin="MRICH_rmin" rmax="MRICH_rmax" length="MRICHCarbonFrame_length" zmin="-1.0*abs(MRICH_zmin)"/>
<module name="MRICH_module1" vis="InvisibleWithDaughters"
width="MRICHCarbonFrame_width"
height="MRICHCarbonFrame_width"
length="MRICHCarbonFrame_length">
<frame vis="AnlGray" thickness="MRICHCarbonFrame_thickness" material="CarbonFiber"/>
<aerogel vis="AnlTeal" length="MRICHAerogel_thickness" width="MRICHAerogel_width" material="AerogelOptical">
<frame vis="AnlGold_1" thickness="MRICHFoam_thickness" material="PolystyreneFoam" />
</aerogel>
<lens vis="AnlViolet" thickness="MRICHFresnelLens_thickness"
pitch="MRICHFresnelLensGroove_pitch" focal_length="6.0*inch"
effective_diameter="MRICHFresnelLensEffectiveDiameter"
width="MRICHAerogel_width"
material="AcrylicOptical"/>
<mirror vis="AnlGray" x1="MRICHMirror_width1" x2="MRICHMirror_width2" length="MRICHMirror_length"
surface="MRICH_MirrorOpticalSurface" thickness="MRICHMirror_thickness" material="AluminumOxide"/>
<photodet width="MRICHGlassWindow_width" thickness="MRICHGlassPhotoDet_thickness" material="PyrexGlassOptical">
<sensor nx="2" ny="2" thickness="MRICHPhotoDet_thickness" width="MRICHPhotoDet_length" material="SiliconOxide"/>
</photodet>
</module>
<!--
<layer id="1">
<array nx="4" ny="4" module="MRICH_module1">
<position x="0" y="0" z="0"/>
</array>
</layer>
-->
</detector>
</detectors>
<readouts>
<readout name="MRICHHits">
<segmentation type="CartesianGridXY" grid_size_x="3*mm" grid_size_y="3*mm" />
<id>system:8,module:14,sensor:8,x:32:-16,y:-16</id>
</readout>
</readouts>
</lccdd>
/run/initialize
#/vis/open OGL 800x800-0+0
/vis/open OGLIQt 800x800-0+0
#
/vis/drawVolume
/vis/viewer/set/viewpointThetaPhi 20 30
/vis/viewer/zoom 1.1
#/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 2
#/vis/filtering/trajectories/create/particleFilter
#/vis/filtering/trajectories/particleFilter-0/add gamma
#/vis/filtering/trajectories/particleFilter-0/invert true
#/vis/filtering/trajectories/particleFilter-0/verbose true
#/vis/filtering/trajectories/particleFilter-0/active true
/vis/ogl/set/displayListLimit 50000
/vis/viewer/flush
/gps/verbose 2
/gps/particle pi-
/gps/number 1
/gps/ene/type Gauss
/gps/ene/mono 8.0 GeV
/gps/ene/sigma 1.0 GeV
/gps/pos/type Volume
/gps/pos/shape Cylinder
/gps/pos/centre 0.0 0.0 0.0 cm
/gps/pos/radius 0.001 cm
/gps/pos/halfz 1 cm
/gps/position 0 0 0 cm
#/gps/direction 0 0.1 1.0
/gps/ang/type iso
/gps/ang/mintheta 10 degree
/gps/ang/maxtheta 20 degree
/run/beamOn 1
/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 2
#/vis/filtering/trajectories/create/particleFilter
#/vis/filtering/trajectories/particleFilter-0/add gamma
#/vis/filtering/trajectories/particleFilter-0/invert true
#/vis/filtering/trajectories/particleFilter-0/verbose true
#/vis/filtering/trajectories/particleFilter-0/active true
/vis/ogl/set/displayListLimit 500000
/vis/viewer/flush
/run/beamOn 1
#/control/execute macro/gun.mac
......@@ -5,9 +5,18 @@
// some utility functions that can be shared
namespace athena::geo {
typedef ROOT::Math::XYPoint Point;
using Point = ROOT::Math::XYPoint;
// fill rectangles in a ring
/** Fill rectangles in a ring (disk).
*
* @param ref 2D reference point.
* @param sx x side length
* @param sy y side length
* @param rmin inner radius of disk
* @param rmax outer radius of disk to fill
* @param phmin phi min
* @param phmax phi max
*/
std::vector<Point> fillRectangles(Point ref, double sx, double sy, double rmin, double rmax,
double phmin = -M_PI, double phmax = M_PI);
// fill squares in a ring
......
//
// Author : Whit Armstrong (warmstrong@anl.gov)
//
#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 "GeometryHelpers.h"
#include "Math/Vector3D.h"
#include "Math/AxisAngle.h"
#include "Math/VectorUtil.h"
using namespace std;
using namespace dd4hep;
using namespace dd4hep::rec;
using Placements = vector<PlacedVolume>;
static Ref_t createDetector(Detector& description, xml::Handle_t e, SensitiveDetector sens){
xml_det_t x_det = e;
Material air = description.material("AirOptical");
Material vacuum = description.vacuum();
string det_name = x_det.nameStr();
//xml::Component pos = x_det.position();
DetElement sdet(det_name, x_det.id());
Assembly assembly(det_name);
sens.setType("photoncounter");
OpticalSurfaceManager surfMgr = description.surfaceManager();
bool projective = getAttrOrDefault(x_det, _Unicode(projective), false);
PlacedVolume pv;
map<string, Volume> modules;
map<string, Placements> sensitives;
map<string, Volume> module_assemblies;
std::map<std::string,DetElement> module_assembly_delements;
int n_sensor = 1;
xml::Component dims = x_det.dimensions();
auto rmin = dims.rmin();
auto rmax = dims.rmax();
auto length = dims.length();
auto zmin = dims.zmin();
// expect only one module (for now)
xml_comp_t x_mod = x_det.child(_U(module));
string mod_name = x_mod.nameStr();
double mod_width = getAttrOrDefault(x_mod, _U(width), 130.0 * mm);
double mod_height = getAttrOrDefault(x_mod, _U(height), 130.0 * mm);
double mod_length = getAttrOrDefault(x_mod, _U(length), 130.0 * mm);
// various components
xml_comp_t x_frame = x_mod.child(_Unicode(frame));
xml_comp_t x_aerogel = x_mod.child(_Unicode(aerogel));
xml_comp_t x_lens = x_mod.child(_Unicode(lens));
xml_comp_t x_mirror = x_mod.child(_Unicode(mirror));
xml_comp_t x_photodet = x_mod.child(_Unicode(photodet));
// module
Box m_solid(mod_width / 2.0, mod_height / 2.0, mod_length / 2.0);
Volume m_volume(mod_name, m_solid, air);
m_volume.setVisAttributes(description.visAttributes(x_mod.visStr()));
DetElement mod_de( mod_name + std::string("_mod_") + std::to_string(1), 1);
// todo module frame
double frame_thickness = getAttrOrDefault(x_frame, _U(thickness), 2.0 * mm);
// aerogel box
xml_comp_t x_aerogel_frame = x_aerogel.child(_Unicode(frame));
double aerogel_width = getAttrOrDefault(x_aerogel, _U(width), 130.0 * mm);
double aerogel_length = getAttrOrDefault(x_aerogel, _U(length), 130.0 * mm);
double foam_thickness = getAttrOrDefault(x_aerogel_frame, _U(thickness), 2.0 * mm);
Material foam_mat = description.material(x_aerogel_frame.materialStr());
Material aerogel_mat = description.material(x_aerogel.materialStr());
auto aerogel_vis = getAttrOrDefault<std::string>(x_aerogel, _U(vis), std::string("InvisibleWithDaughters"));
auto foam_vis = getAttrOrDefault<std::string>(x_aerogel_frame, _U(vis), std::string("RedVis"));
// aerogel foam frame
Box foam_box(aerogel_width / 2.0 + foam_thickness, aerogel_width / 2.0 + foam_thickness, (aerogel_length + foam_thickness) / 2.0);
Box aerogel_sub_box(aerogel_width / 2.0, aerogel_width / 2.0, (aerogel_length + foam_thickness) / 2.0);
SubtractionSolid foam_frame_solid(foam_box, aerogel_sub_box, Position(0, 0, foam_thickness));
Volume foam_vol(mod_name+"_aerogel_frame", foam_frame_solid, foam_mat);
foam_vol.setVisAttributes(description.visAttributes(foam_vis));
double foam_frame_zpos = -mod_length / 2.0 + frame_thickness + (aerogel_length + foam_thickness) / 2.0;
m_volume.placeVolume(foam_vol,Position(0,0,foam_frame_zpos));
// aerogel
Box aerogel_box(aerogel_width / 2.0, aerogel_width / 2.0, (aerogel_length) / 2.0);
Volume aerogel_vol(mod_name+"_aerogel", aerogel_box, aerogel_mat);
aerogel_vol.setVisAttributes(description.visAttributes(aerogel_vis));
double aerogel_zpos = foam_frame_zpos + foam_thickness / 2.0;
pv = m_volume.placeVolume(aerogel_vol,Position(0,0,aerogel_zpos));
DetElement aerogel_de(mod_de, mod_name + std::string("_aerogel_de") + std::to_string(1), 1);
aerogel_de.setPlacement(pv);
auto aerogel_surf = surfMgr.opticalSurface(dd4hep::getAttrOrDefault(x_aerogel, _Unicode(surface), "MRICH_AerogelOpticalSurface"));
SkinSurface skin0(description, aerogel_de, Form("MRICH_aerogel_skin_surface_%d", 1), aerogel_surf, aerogel_vol);
skin0.isValid();
// Fresnel Lens
// - The lens has a constant groove pitch (delta r) as opposed to fixing the groove height.
// - The lens area outside of the effective diamtere is flat.
// - The grooves are not curved, rather they are polycone shaped, ie a flat approximating the curvature.
auto lens_vis = getAttrOrDefault<std::string>(x_lens, _U(vis), std::string("AnlBlue"));
double groove_pitch = getAttrOrDefault(x_lens, _Unicode(pitch), 0.5 * mm);
double lens_f = getAttrOrDefault(x_lens, _Unicode(focal_length), 6.0*2.54*cm);
double eff_diameter = getAttrOrDefault(x_lens, _Unicode(effective_diameter), 130.0 * mm);
double lens_width = getAttrOrDefault(x_lens, _Unicode(width), 6.7*2.54*cm);
double center_thickness = getAttrOrDefault(x_lens, _U(thickness), 2.0 * mm);
double n_acrylic = 1.49;
double lens_curvature = 1.0 / (lens_f*(n_acrylic - 1.0));
double full_ring_rmax = std::min(eff_diameter / 2.0, lens_width/2.0);
double N_grooves = std::ceil((full_ring_rmax) / groove_pitch);
double groove_last_rmin = (N_grooves - 1) * groove_pitch;
double groove_last_rmax = N_grooves * groove_pitch;
auto groove_sagitta = [&](double r) { return lens_curvature * std::pow(r, 2) / (1.0 + 1.0); };
double lens_thickness = groove_sagitta(groove_last_rmax) - groove_sagitta(groove_last_rmin) + center_thickness;
Material lens_mat = description.material(x_lens.materialStr());
Box lens_box(lens_width / 2.0, lens_width / 2.0, (center_thickness) / 2.0);
SubtractionSolid flat_lens(lens_box, Tube(0.0, full_ring_rmax, 2 * center_thickness));
Assembly lens_vol(mod_name + "_lens");
Volume flatpart_lens_vol( "flatpart_lens", flat_lens, lens_mat);
lens_vol.placeVolume(flatpart_lens_vol);//,Position(0,0,lens_zpos));
Solid fresnel_lens_solid;
int i_groove = 0;
double groove_rmax = groove_pitch;
double groove_rmin = 0;
while ( groove_rmax <= full_ring_rmax ) {
double dZ = groove_sagitta(groove_rmax) - groove_sagitta(groove_rmin);
//std::cout << " dZ = " << dZ << ", lens_thickness = " << lens_thickness << "\n";
//std::cout << "groove_rmin = " << groove_rmin << "\n";
//std::cout << "groove_rmax = " << groove_rmax << "\n";
Polycone groove_solid(0, 2.0 * M_PI,
{groove_rmin, groove_rmin, groove_rmin},
{groove_rmax, groove_rmax, groove_rmin},
{-lens_thickness/2.0, lens_thickness/2.0-dZ, lens_thickness/2.0});
Volume lens_groove_vol("lens_groove_" + std::to_string(i_groove), groove_solid, lens_mat);
lens_vol.placeVolume(lens_groove_vol); //,Position(0,0,lens_zpos));
//Volume groove_vol(groove_solid, lens_mat, par->name.c_str(), 0, 0, 0);
//new G4PVPlacement(0, par->pos, Groove_log[i], par->name.c_str(), motherLV, false, 0, OverlapCheck());
//phi1 = phi1 + halfpi; //g4 pre-defined: halfpi=pi/2
//Tube sub_cylinder(r0, r1, 3*eff_diameter);
//IntersectionSolid groove_solid(lens_box,lens_sphere, Position(0,0,-eff_diameter/2.0 + lens_thickness/2.0+(t-lens_t)/2.0 ));
//IntersectionSolid lens_ring(groove_solid, sub_cylinder);
//if (i_groove == 0) {
// fresnel_lens_solid = groove_solid;
//} else {
// fresnel_lens_solid = UnionSolid(fresnel_lens_solid, groove_solid);
//}
//r0 = r1;
//if(i_groove > 3) {
// SubtractionSolid flat_lens(lens_box,Tube(0.0, r0, 3*eff_diameter));
// fresnel_lens_solid = UnionSolid(fresnel_lens_solid, flat_lens);
// break; // temporary
//}
i_groove++;
groove_rmin = (i_groove )*groove_pitch;
groove_rmax = (i_groove+1)*groove_pitch;
}
//fresnel_lens_solid = UnionSolid(fresnel_lens_solid, flat_lens);
//Volume lens_vol(mod_name + "_lens", fresnel_lens_solid, lens_mat);
lens_vol.setVisAttributes(description.visAttributes(lens_vis));
double lens_zpos = aerogel_zpos +aerogel_length/ 2.0 + foam_thickness + lens_thickness/2.0;
pv = m_volume.placeVolume(lens_vol,Position(0,0,lens_zpos));
DetElement lens_de(mod_de, mod_name + std::string("_lens_de") + std::to_string(1), 1);
lens_de.setPlacement(pv);
auto surf = surfMgr.opticalSurface(dd4hep::getAttrOrDefault(x_lens, _Unicode(surface), "MRICH_LensOpticalSurface"));
SkinSurface skin(description, lens_de, Form("MRichFresnelLens_skin_surface_%d", 1), surf, lens_vol);
skin.isValid();
// mirror
auto mirror_vis = getAttrOrDefault<std::string>(x_mirror, _U(vis), std::string("AnlGray"));
double mirror_x1 = getAttrOrDefault(x_mirror, _U(x1), 100.0 * mm);
double mirror_x2 = getAttrOrDefault(x_mirror, _U(x2), 80.0 * mm);
double mirror_length = getAttrOrDefault(x_mirror, _U(length), 130.0 * mm);
double mirror_thickness = getAttrOrDefault(x_mirror, _U(thickness), 2.0 * mm);
double outer_x1 = (mirror_x1+mirror_thickness)/2.0;
double outer_x2 = (mirror_x2+mirror_thickness)/2.0;
Trd2 outer_mirror_trd(outer_x1, outer_x2, outer_x1, outer_x2, mirror_length/2.0);
Trd2 inner_mirror_trd(mirror_x1 / 2.0, mirror_x2 / 2.0, mirror_x1 / 2.0,mirror_x2 / 2.0, mirror_length/2.0+0.1*mm);
SubtractionSolid mirror_solid(outer_mirror_trd, inner_mirror_trd);
Material mirror_mat = description.material(x_mirror.materialStr());
Volume mirror_vol(mod_name+"_mirror", mirror_solid, mirror_mat);
double mirror_zpos = lens_zpos + lens_thickness/2.0 + foam_thickness + mirror_length/2.0;
pv = m_volume.placeVolume(mirror_vol,Position(0,0,mirror_zpos));
DetElement mirror_de(mod_de, mod_name + std::string("_mirror_de") + std::to_string(1), 1);
mirror_de.setPlacement(pv);
auto mirror_surf = surfMgr.opticalSurface(dd4hep::getAttrOrDefault(x_mirror, _Unicode(surface), "MRICH_MirrorOpticalSurface"));
SkinSurface skin1(description, mirror_de, Form("MRICH_mirror_skin_surface_%d", 1), mirror_surf, mirror_vol);
skin1.isValid();
// photon detector
xml_comp_t x_photodet_sensor = x_photodet.child(_Unicode(sensor));
auto photodet_vis = getAttrOrDefault<std::string>(x_photodet, _U(vis), std::string("AnlRed"));
double photodet_width = getAttrOrDefault(x_photodet, _U(width), 130.0 * mm);
double photodet_thickness = getAttrOrDefault(x_photodet, _U(thickness), 2.0 * mm);
double sensor_thickness = getAttrOrDefault(x_photodet_sensor, _U(thickness), 2.0 * mm);
Material photodet_mat = description.material(x_photodet.materialStr());
Material sensor_mat = description.material(x_photodet_sensor.materialStr());
int sensor_nx = getAttrOrDefault(x_photodet_sensor, _Unicode(nx), 2);
int sensor_ny = getAttrOrDefault(x_photodet_sensor, _Unicode(ny), 2);
Box window_box(photodet_width/2.0,photodet_width/2.0,photodet_thickness/2.0);
Volume window_vol(mod_name+"_window", window_box, photodet_mat);
window_vol.setSensitiveDetector(sens);
double window_zpos = mirror_zpos + mirror_length/2.0+photodet_thickness/2.0;
pv = m_volume.placeVolume(window_vol,Position(0,0,window_zpos));
DetElement comp_de(mod_de, std::string("mod_sensor_de_") + std::to_string(1) , 1);
comp_de.setPlacement(pv);
pv.addPhysVolID("sensor", n_sensor);
//for (size_t ic = 0; ic < sensVols.size(); ++ic) {
// PlacedVolume sens_pv = sensVols[ic];
// DetElement comp_de(mod_de, std::string("de_") + sens_pv.volume().name(), ic + 1);
// comp_de.setPlacement(sens_pv);
// // Acts::ActsExtension* sensorExtension = new Acts::ActsExtension();
// //// sensorExtension->addType("sensor", "detector");
// // comp_de.addExtension<Acts::ActsExtension>(sensorExtension);
// //// comp_de.setAttributes(description, sens_pv.volume(), x_layer.regionStr(),
// //// x_layer.limitsStr(),
// //// xml_det_t(xmleles[m_nam]).visStr());
//}
//DetElement window_de(sdet, mod_name + std::string("_window_de") + std::to_string(1), 1);
//window_de.setPlacement(pv);
window_vol.setSensitiveDetector(sens);
sensitives[mod_name].push_back(pv);
++n_sensor;
modules[mod_name] = m_volume;
module_assembly_delements[mod_name] = mod_de;
// end module
int i_mod = 1;
// detector envelope
Tube envShape(rmin, rmax, length / 2., 0., 2 * M_PI);
Volume envVol("MRICH_Envelope", envShape, air);
envVol.setVisAttributes(description.visAttributes(x_det.visStr()));
// place modules in the sectors (disk)
auto points = athena::geo::fillSquares({0., 0.}, mod_width, rmin, rmax);
// mod_name = ...
Placements& sensVols = sensitives[mod_name];
auto mod_v = modules[mod_name];
// determine module direction, always facing z = 0
double roty = dims.zmin() < 0. ? -M_PI : 0 ;
int imod = 1;
for (auto& p : points) {
ROOT::Math::XYZVector x_location(p.x(), p.y(), zmin+std::signbit(zmin)*mod_length/2.0);
ROOT::Math::XYZVector z_dir(0, 0, 1);
ROOT::Math::XYZVector x_dir(1, 0, 0);
ROOT::Math::XYZVector rot_axis = x_location.Cross(z_dir);
double rot_angle = ROOT::Math::VectorUtil::Angle(z_dir,x_location);
ROOT::Math::AxisAngle proj_rot(rot_axis,-1.0*rot_angle);
ROOT::Math::AxisAngle grid_fix_rot(x_location,0.0*rot_angle);
auto new_x_dir = grid_fix_rot*x_dir;
// operations are inversely ordered
Transform3D tr = Translation3D(p.x(), p.y(), 0.) // move to position
* RotationY(roty); // facing z = 0.
if(projective) {
tr = Translation3D(p.x(), p.y(), 0.) // move to position
* grid_fix_rot // keep the modules oriented vertially
* proj_rot // projective rotation
* RotationY(roty); // facing z = 0.
}
// mod placement
pv = envVol.placeVolume(mod_v, tr);
pv.addPhysVolID("module", i_mod);
auto mod_det_element = module_assembly_delements[mod_name].clone(mod_name + "__" + std::to_string(i_mod));
mod_det_element.setPlacement(pv);
sdet.add(mod_det_element);
i_mod++;
}
// place envelope
Volume motherVol = description.pickMotherVolume(sdet);
PlacedVolume envPV = motherVol.placeVolume(envVol, Position(0, 0, zmin));
envPV.addPhysVolID("system", x_det.id());
sdet.setPlacement(envPV);
return sdet;
}
//void addModules(Volume &mother, xml::DetElement &detElem, Detector &description, 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 = description.material(mods.materialStr());
// auto gasMat = description.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(description.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(athena_MRICH, 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 "GeometryHelpers.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 = athena::geo::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)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment