diff --git a/compact/ce_ecal.xml b/compact/ce_ecal.xml new file mode 100644 index 0000000000000000000000000000000000000000..35b7ba0e6f4db72a7e2540a139f49790859ea59b --- /dev/null +++ b/compact/ce_ecal.xml @@ -0,0 +1,109 @@ +<lccdd> + + <define> + <constant name="CrystalEndcap_rmin" value="VertexTrackerBarrel_rmin"/> + <constant name="CrystalEndcap_rmax" value="EcalEndcapN_rmin"/> + <constant name="CrystalBox_x_length" value="20.0*mm"/> + <constant name="CrystalBox_y_length" value="20.0*mm"/> + <constant name="CrystalBox_z_length" value="200.0*mm"/> + <constant name="CrystalBox_offset" value="0.000001*mm"/> + <constant name="CrystalEndcap_x_pos" value="0.0*m"/> + <constant name="CrystalEndcap_y_pos" value="0.0*m"/> + <constant name="CrystalEndcap_z_pos" value="-EcalEndcapN_zmin"/> + </define> + + + <limits> + </limits> + + <regions> + </regions> + + <!-- Common Generic visualization attributes --> + <comment>Common Generic visualization attributes</comment> + <display> + </display> + + <detectors> + <comment> + ------------------------------- + Backwards Endcap EM Calorimeter + ------------------------------- + A layered EM calorimeter with tungsten and silicon (or scintillator) strips + </comment> + <detector id="ECalEndcapN_ID" + name="EcalEndcapN" + type="refdet_PolyhedraEndcapCalorimeter2" + reflect="true" + readout="EcalEndcapHits" + vis="EcalEndcapVis" + calorimeterType="EM_ENDCAP"> + <position x="0" y="0" z="0"/> + <dimensions + numsides="CaloSides" + zmin="EcalEndcapN_zmin" + rmin="EcalEndcapN_rmin" + rmax="EcalBarrel_rmin " /> + <layer repeat="EcalEndcapNLayer1_NRepeat"> + <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> + <slice material="Copper" thickness="EcalCopperThickness"/> + <slice material="Kapton" thickness="EcalKaptonThickness"/> + <slice material="Air" thickness="EcalAir1Thickness"/> + </layer> + <layer repeat="EcalEndcapNLayer2_NRepeat"> + <slice material="TungstenDens24" thickness="EcalThinTungstenThickness"/> + <slice material="Air" thickness="EcalAir2Thickness"/> + <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> + <slice material="Copper" thickness="EcalCopperThickness"/> + <slice material="Kapton" thickness="EcalKaptonThickness"/> + <slice material="Air" thickness="EcalAir1Thickness"/> + </layer> + <layer repeat="EcalEndcapNLayer3_NRepeat"> + <slice material="TungstenDens24" thickness="EcalThickTungstenThickness"/> + <slice material="Air" thickness="EcalAir2Thickness"/> + <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> + <slice material="Copper" thickness="EcalCopperThickness"/> + <slice material="Kapton" thickness="EcalKaptonThickness"/> + <slice material="Air" thickness="EcalAir1Thickness"/> + </layer> + </detector> + + <comment> + ------------------- + Crystal Calorimeter + ------------------- + Backward (negative pseudorapidity) angle electron detector + </comment> + <detector + id="CrystalEndcap_ID" + name="ElectronECAL" + type="CrystalCalorimeterEndcap" + readout="CrystalEcalHits" + vis="GreenVis"> + <position x="CrystalEndcap_x_pos" y="CrystalEndcap_y_pos" z="CrystalEndcap_z_pos" /> + <dimensions rmin="CrystalEndcap_rmin" rmax="CrystalEndcap_rmax" x="CrystalBox_x_length" y="CrystalBox_y_length" z="CrystalBox_z_length" delta="CrystalBox_offset"/> + </detector> + </detectors> + + <!-- Definition of the readout segmentation/definition --> + <readouts> + <!-- + <readout name="PlaneTrackerHits"> + <segmentation type="CartesianGridXY" grid_size_x="20.0*mm" grid_size_y="20.0*mm" /> + <id>system:5,module:4,x:32:-16,y:-16</id> + </readout> + --> + <readout name="CrystalEcalHits"> + <segmentation type="CartesianGridXY" grid_size_x="CrystalBox_x_length" grid_size_y="CrystalBox_y_length" /> + <id>system:8,sector:4,module:20,x:48:-8,y:-8</id> + </readout> + <readout name="EcalEndcapNHits"> + <segmentation type="CartesianGridXY" grid_size_x="3.5 * mm" grid_size_y="3.5 * mm"/> + <id>system:8,barrel:3,module:4,layer:6,slice:5,x:32:-16,y:-16</id> + </readout> + </readouts> + + <plugins> + </plugins> + +</lccdd> diff --git a/compact/ce_ecal_crystal_glass.xml b/compact/ce_ecal_crystal_glass.xml new file mode 100644 index 0000000000000000000000000000000000000000..cc36383be8b946cbd200512db58a969baa5d4387 --- /dev/null +++ b/compact/ce_ecal_crystal_glass.xml @@ -0,0 +1,91 @@ +<lccdd> + <define> + <constant name="CrystalModule_sx" value="20.0*mm"/> + <constant name="CrystalModule_sy" value="20.0*mm"/> + <constant name="CrystalModule_sz" value="20.0*cm"/> + <constant name="CrystalModule_wrap" value="0.5*mm"/> + <constant name="CrystalModule_dz" value="-8*cm"/> + <constant name="GlassModule_sx" value="40.0*mm"/> + <constant name="GlassModule_sy" value="40.0*mm"/> + <constant name="GlassModule_sz" value="40.0*cm"/> + <constant name="GlassModule_dz" value="0.0*cm"/> + <constant name="GlassModule_wrap" value="1.0*mm"/> + <constant name="EcalEndcapN_z0" value="-EcalEndcapN_zmin-max(CrystalModule_sz,GlassModule_sz)/2."/> + </define> + <detectors> + <comment> Backwards Endcap EM Calorimeter, placements generated by script </comment> + <detector id="ECalEndcapN_ID" name="EcalEndcapN" type="HomogeneousCalorimeter" readout="EcalEndcapNHits"> + <position x="0" y="0" z="EcalEndcapN_z0"/> + <rotation x="0" y="0" z="0"/> + <placements> + <lines sector="1" mirrorx="true" mirrory="true"> + <module sizex="CrystalModule_sx" sizey="CrystalModule_sy" sizez="CrystalModule_sz" material="PbWO4" vis="AnlTeal"/> + <wrapper thickness="CrystalModule_wrap" material="Epoxy" vis="WhiteVis"/> + <position z="CrystalModule_dz"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*1/2." begin="7" nmods="17"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*3/2." begin="7" nmods="17"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*5/2." begin="7" nmods="17"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*7/2." begin="6" nmods="18"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*9/2." begin="6" nmods="18"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*11/2." begin="5" nmods="19"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*13/2." begin="3" nmods="19"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*15/2." begin="0" nmods="22"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*17/2." begin="0" nmods="22"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*19/2." begin="0" nmods="22"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*21/2." begin="0" nmods="22"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*23/2." begin="0" nmods="22"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*25/2." begin="0" nmods="20"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*27/2." begin="0" nmods="20"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*29/2." begin="0" nmods="18"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*31/2." begin="0" nmods="18"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*33/2." begin="0" nmods="16"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*35/2." begin="0" nmods="16"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*37/2." begin="0" nmods="14"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*39/2." begin="0" nmods="14"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*41/2." begin="0" nmods="12"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*43/2." begin="0" nmods="12"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*45/2." begin="0" nmods="6"/> + <line x="(CrystalModule_sx + CrystalModule_wrap)/2." y="(CrystalModule_sy + CrystalModule_wrap)*47/2." begin="0" nmods="6"/> + </lines> + <lines sector="2" mirrorx="true" mirrory="true"> + <module sizex="GlassModule_sx" sizey="GlassModule_sy" sizez="GlassModule_sz" material="PbGlass" vis="AnlBlue"/> + <wrapper thickness="GlassModule_wrap" material="Epoxy" vis="WhiteVis"/> + <position z="GlassModule_dz"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*1/2." begin="12" nmods="11"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*3/2." begin="12" nmods="11"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*5/2." begin="12" nmods="11"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*7/2." begin="11" nmods="12"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*9/2." begin="11" nmods="12"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*11/2." begin="11" nmods="12"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*13/2." begin="10" nmods="12"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*15/2." begin="9" nmods="13"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*17/2." begin="8" nmods="13"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*19/2." begin="7" nmods="14"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*21/2." begin="6" nmods="14"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*23/2." begin="3" nmods="16"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*25/2." begin="0" nmods="19"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*27/2." begin="0" nmods="18"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*29/2." begin="0" nmods="18"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*31/2." begin="0" nmods="17"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*33/2." begin="0" nmods="17"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*35/2." begin="0" nmods="15"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*37/2." begin="0" nmods="13"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*39/2." begin="0" nmods="11"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*41/2." begin="0" nmods="10"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*43/2." begin="0" nmods="8"/> + <line x="(GlassModule_sx + GlassModule_wrap)/2." y="(GlassModule_sy + GlassModule_wrap)*45/2." begin="0" nmods="6"/> + </lines> + </placements> + </detector> + </detectors> + <readouts> + <comment>Effectively no segmentation, the segmentation is used to provide cell dimension info</comment> + <readout name="EcalEndcapNHits"> + <segmentation type="MultiSegmentation" key="sector"> + <segmentation name="CrystalSeg" key_value="1" type="CartesianGridXY" grid_size_x="CrystalModule_sx + CrystalModule_wrap" grid_size_y="CrystalModule_sy + CrystalModule_wrap"/> + <segmentation name="GlassSeg" key_value="2" type="CartesianGridXY" grid_size_x="GlassModule_sx + GlassModule_wrap" grid_size_y="GlassModule_sy + GlassModule_wrap"/> + </segmentation> + <id>system:8,sector:4,module:20,x:32:-16,y:-16</id> + </readout> + </readouts> +</lccdd> diff --git a/compact/definitions.xml b/compact/definitions.xml index efdbdb09e008beb4efd9f06ee6fa1fc56d0f764f..ca239b3c8a8d15c32e915c3d5475bff1b2dd1dbf 100644 --- a/compact/definitions.xml +++ b/compact/definitions.xml @@ -503,7 +503,7 @@ <constant name="EcalEndcapP_length" value="EndcapPTotalCal_length * EndcapP_CalDivide"/> <!--constant name="EcalEndcapN_length" value="EndcapNTotalCal_length * EndcapN_CalDivide"/--> - <constant name="EcalEndcapN_length" value="40*cm"/> + <constant name="EcalEndcapN_length" value="41*cm"/> <comment> These need to be set in sync with the forward and backward detectors </comment> diff --git a/compact/ecal.xml b/compact/ecal.xml index 9357337739bd00f52714255f35ebdbed1953a3ec..9761fd9ba7e3e6f685b54ff6aa1cf0c90e9c0034 100644 --- a/compact/ecal.xml +++ b/compact/ecal.xml @@ -1,18 +1,6 @@ <lccdd> <define> - <constant name="CrystalEndcap_rmin" value="VertexTrackerBarrel_rmin"/> - <constant name="CrystalEndcap_rmax" value="EcalEndcapN_rmin"/> - <constant name="CrystalBox_x_length" value="20.0*mm"/> - <constant name="CrystalBox_y_length" value="20.0*mm"/> - <constant name="CrystalBox_z_length" value="200.0*mm"/> - <constant name="CrystalBox_offset" value="0.000001*mm"/> - <constant name="CrystalEndcap_x_pos" value="0.0*m"/> - <constant name="CrystalEndcap_y_pos" value="0.0*m"/> - <constant name="CrystalEndcap_z_pos" value="-EcalEndcapN_zmin"/> - - - <constant name="EcalEndcapP_rmax" value="Solenoid_rmax "/> </define> @@ -29,6 +17,8 @@ </display> <include ref="ecal_barrel.xml"/> + <!--<include ref="ce_ecal.xml"/>--> + <include ref="ce_ecal_crystal_glass.xml"/> <detectors> <comment> @@ -41,7 +31,7 @@ name="EcalEndcapP" reflect="false" type="refdet_PolyhedraEndcapCalorimeter2" - readout="EcalEndcapHits" + readout="EcalEndcapPHits" vis="EcalEndcapVis" calorimeterType="EM_ENDCAP" > <position x="0" y="0" z="-0"/> @@ -73,65 +63,6 @@ <slice material="Air" thickness="EcalAir1Thickness"/> </layer> </detector> - - <comment> - ------------------------------- - Backwards Endcap EM Calorimeter - ------------------------------- - A layered EM calorimeter with tungsten and silicon (or scintillator) strips - </comment> - <detector id="ECalEndcapN_ID" - name="EcalEndcapN" - type="refdet_PolyhedraEndcapCalorimeter2" - reflect="true" - readout="EcalEndcapHits" - vis="EcalEndcapVis" - calorimeterType="EM_ENDCAP"> - <position x="0" y="0" z="0"/> - <dimensions - numsides="CaloSides" - zmin="EcalEndcapN_zmin" - rmin="EcalEndcapN_rmin" - rmax="EcalBarrel_rmin " /> - <layer repeat="EcalEndcapNLayer1_NRepeat"> - <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> - <slice material="Copper" thickness="EcalCopperThickness"/> - <slice material="Kapton" thickness="EcalKaptonThickness"/> - <slice material="Air" thickness="EcalAir1Thickness"/> - </layer> - <layer repeat="EcalEndcapNLayer2_NRepeat"> - <slice material="TungstenDens24" thickness="EcalThinTungstenThickness"/> - <slice material="Air" thickness="EcalAir2Thickness"/> - <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> - <slice material="Copper" thickness="EcalCopperThickness"/> - <slice material="Kapton" thickness="EcalKaptonThickness"/> - <slice material="Air" thickness="EcalAir1Thickness"/> - </layer> - <layer repeat="EcalEndcapNLayer3_NRepeat"> - <slice material="TungstenDens24" thickness="EcalThickTungstenThickness"/> - <slice material="Air" thickness="EcalAir2Thickness"/> - <slice material="Silicon" thickness="EcalSiliconThickness" sensitive="yes" limits="cal_limits"/> - <slice material="Copper" thickness="EcalCopperThickness"/> - <slice material="Kapton" thickness="EcalKaptonThickness"/> - <slice material="Air" thickness="EcalAir1Thickness"/> - </layer> - </detector> - - <comment> - ------------------- - Crystal Calorimeter - ------------------- - Backward (negative pseudorapidity) angle electron detector - </comment> - <detector - id="CrystalEndcap_ID" - name="ElectronECAL" - type="CrystalCalorimeterEndcap" - readout="CrystalEcalHits" - vis="GreenVis"> - <position x="CrystalEndcap_x_pos" y="CrystalEndcap_y_pos" z="CrystalEndcap_z_pos" /> - <dimensions rmin="CrystalEndcap_rmin" rmax="CrystalEndcap_rmax" x="CrystalBox_x_length" y="CrystalBox_y_length" z="CrystalBox_z_length" delta="CrystalBox_offset"/> - </detector> </detectors> <!-- Definition of the readout segmentation/definition --> @@ -142,11 +73,7 @@ <id>system:5,module:4,x:32:-16,y:-16</id> </readout> --> - <readout name="CrystalEcalHits"> - <segmentation type="CartesianGridXY" grid_size_x="CrystalBox_x_length" grid_size_y="CrystalBox_y_length" /> - <id>system:8,sector:4,module:20,x:48:-8,y:-8</id> - </readout> - <readout name="EcalEndcapHits"> + <readout name="EcalEndcapPHits"> <segmentation type="CartesianGridXY" grid_size_x="3.5 * mm" grid_size_y="3.5 * mm"/> <id>system:8,barrel:3,module:4,layer:6,slice:5,x:32:-16,y:-16</id> </readout> diff --git a/scripts/ce_ecal_placement.py b/scripts/ce_ecal_placement.py new file mode 100644 index 0000000000000000000000000000000000000000..138bc7c43ce70cf1ea018e24359d764dd1b642b4 --- /dev/null +++ b/scripts/ce_ecal_placement.py @@ -0,0 +1,237 @@ +''' + A script to calcualte placement of ecal endcap modules + lxml is not included in container, get it by simply typing 'pip install lxml' + Author: Chao Peng (ANL) + Date: 06/17/2021 +''' + +import numpy as np +import argparse +from lxml import etree as ET + + +# start pos, total number of blocks +CRYSTAL_ALIGNMENT = [ + (8, 17), + (8, 17), + (8, 17), + (7, 18), + (7, 18), + (6, 19), + (4, 19), + (1, 22), + (1, 22), + (1, 22), + (1, 22), + (1, 22), + (1, 20), + (1, 20), + (1, 18), + (1, 18), + (1, 16), + (1, 16), + (1, 14), + (1, 14), + (1, 12), + (1, 12), + (1, 6), + (1, 6), +] + +GLASS_ALIGNMENT = [ + (13, 11), + (13, 11), + (13, 11), + (12, 12), + (12, 12), + (12, 12), + (11, 12), + (10, 13), + (9, 13), + (8, 14), + (7, 14), + (4, 16), + (1, 19), + (1, 18), + (1, 18), + (1, 17), + (1, 17), + (1, 15), + (1, 13), + (1, 11), + (1, 10), + (1, 8), + (1, 6), +] + +# calculate positions of modules with a quad-alignment and module size +def individual_placement(alignment, module_x=20.5, module_y=20.5): + placements = [] + for row, (start, num) in enumerate(alignment): + for col in np.arange(start - 1, start + num - 1): + placements.append(((col + 0.5)*module_y, (row + 0.5)*module_x)) + placements = np.asarray(placements) + return np.vstack((placements, + np.vstack((placements.T[0]*-1., placements.T[1])).T, + np.vstack((placements.T[0], placements.T[1]*-1.)).T, + np.vstack((placements.T[0]*-1., placements.T[1]*-1.)).T)) + + + +if __name__ == '__main__': + parser = argparse.ArgumentParser() + parser.add_argument('-s', '--save', default='compact/ce_ecal_crystal_glass.xml', + help='path to save compact file.') + args = parser.parse_args() + + data = ET.Element('lccdd') + defines = ET.SubElement(data, 'define') + + # constants: name, value + consts = [ + ('CrystalModule_sx', '20.0*mm'), + ('CrystalModule_sy', '20.0*mm'), + ('CrystalModule_sz', '20.0*cm'), + ('CrystalModule_wrap', '0.5*mm'), + ('CrystalModule_dz', '-8*cm'), + ('GlassModule_sx', '40.0*mm'), + ('GlassModule_sy', '40.0*mm'), + ('GlassModule_sz', '40.0*cm'), + ('GlassModule_dz', '0.0*cm'), + ('GlassModule_wrap', '1.0*mm'), + ('EcalEndcapN_z0', '-EcalEndcapN_zmin-max(CrystalModule_sz,GlassModule_sz)/2.'), + ] + + for name, value in consts: + constant = ET.SubElement(defines, 'constant') + constant.set('name', name) + constant.set('value', value) + + # this value will be used multiple times, so define it here + readout_name = 'EcalEndcapNHits' + + # detector and its dimension/position/rotation + dets = ET.SubElement(data, 'detectors') + cmt = ET.SubElement(dets, 'comment') + cmt.text = ' Backwards Endcap EM Calorimeter, placements generated by script ' + + det = ET.SubElement(dets, 'detector') + det.set('id', 'ECalEndcapN_ID') + det.set('name', 'EcalEndcapN') + det.set('type', 'HomogeneousCalorimeter') + det.set('readout', readout_name) + + pos = ET.SubElement(det, 'position') + pos.set('x', '0') + pos.set('y', '0') + pos.set('z', 'EcalEndcapN_z0') + + rot = ET.SubElement(det, 'rotation') + rot.set('x', '0') + rot.set('y', '0') + rot.set('z', '0') + + # placements of modules + plm = ET.SubElement(det, 'placements') + + # crystals + # crystal = ET.SubElement(plm, 'individuals') + crystal = ET.SubElement(plm, 'lines') + crystal.set('sector', '1') + crystal_mod = ET.SubElement(crystal, 'module') + crystal_mod.set('sizex', 'CrystalModule_sx') + crystal_mod.set('sizey', 'CrystalModule_sy') + crystal_mod.set('sizez', 'CrystalModule_sz') + crystal_mod.set('material', 'PbWO4') + crystal_mod.set('vis', 'AnlTeal') + crystal_wrap = ET.SubElement(crystal, 'wrapper') + crystal_wrap.set('thickness', 'CrystalModule_wrap') + crystal_wrap.set('material', 'Epoxy') + crystal_wrap.set('vis', 'WhiteVis') + crystal_pos = ET.SubElement(crystal, 'position') + crystal_pos.set('z', 'CrystalModule_dz') + # crystal placements (for individuals) + # for m, (x, y) in enumerate(individual_placement(CRYSTAL_ALIGNMENT)): + # module = ET.SubElement(crystal, 'placement') + # module.set('x', '{:.3f}*mm'.format(x)) + # module.set('y', '{:.3f}*mm'.format(y)) + # module.set('z', '0') + # module.set('id', '{:d}'.format(m)) + + # crystal placements (for lines) + crystal.set('mirrorx', 'true') + crystal.set('mirrory', 'true') + for row, (begin, nmods) in enumerate(CRYSTAL_ALIGNMENT): + line = ET.SubElement(crystal, 'line') + line.set('x', '(CrystalModule_sx + CrystalModule_wrap)/2.') + line.set('y', '(CrystalModule_sy + CrystalModule_wrap)*{:d}/2.'.format(row*2 + 1)) + line.set('begin', '{:d}'.format(begin - 1)) + line.set('nmods', '{:d}'.format(nmods)) + + + # glass + # glass = ET.SubElement(plm, 'individuals') + glass = ET.SubElement(plm, 'lines') + glass.set('sector', '2') + glass_mod = ET.SubElement(glass, 'module') + glass_mod.set('sizex', 'GlassModule_sx') + glass_mod.set('sizey', 'GlassModule_sy') + glass_mod.set('sizez', 'GlassModule_sz') + # TODO: change glass material + glass_mod.set('material', 'PbGlass') + glass_mod.set('vis', 'AnlBlue') + glass_wrap = ET.SubElement(glass, 'wrapper') + glass_wrap.set('thickness', 'GlassModule_wrap') + glass_wrap.set('material', 'Epoxy') + glass_wrap.set('vis', 'WhiteVis') + glass_pos = ET.SubElement(glass, 'position') + glass_pos.set('z', 'GlassModule_dz') + # crystal placements (for individuals) + # for m, (x, y) in enumerate(individual_placement(GLASS_ALIGNMENT, 41.0, 41.0)): + # module = ET.SubElement(glass, 'placement') + # module.set('x', '{:.3f}*mm'.format(x)) + # module.set('y', '{:.3f}*mm'.format(y)) + # module.set('z', '0') + # module.set('id', '{:d}'.format(m)) + + # crystal placements (for lines) + glass.set('mirrorx', 'true') + glass.set('mirrory', 'true') + for row, (begin, nmods) in enumerate(GLASS_ALIGNMENT): + line = ET.SubElement(glass, 'line') + line.set('x', '(GlassModule_sx + GlassModule_wrap)/2.') + line.set('y', '(GlassModule_sy + GlassModule_wrap)*{:d}/2.'.format(row*2 + 1)) + line.set('begin', '{:d}'.format(begin - 1)) + line.set('nmods', '{:d}'.format(nmods)) + + # readout + readouts = ET.SubElement(data, 'readouts') + cmt = ET.SubElement(readouts, 'comment') + cmt.text = 'Effectively no segmentation, the segmentation is used to provide cell dimension info' + readout = ET.SubElement(readouts, 'readout') + readout.set('name', readout_name) + seg = ET.SubElement(readout, 'segmentation') + # need segmentation to provide cell dimension info + # seg.set('type', 'NoSegmentation') + seg.set('type', 'MultiSegmentation') + seg.set('key', 'sector') + crystal_seg = ET.SubElement(seg, 'segmentation') + crystal_seg.set('name', 'CrystalSeg') + crystal_seg.set('key_value', '1') + crystal_seg.set('type', 'CartesianGridXY') + crystal_seg.set('grid_size_x', 'CrystalModule_sx + CrystalModule_wrap') + crystal_seg.set('grid_size_y', 'CrystalModule_sy + CrystalModule_wrap') + glass_seg = ET.SubElement(seg, 'segmentation') + glass_seg.set('name', 'GlassSeg') + glass_seg.set('key_value', '2') + glass_seg.set('type', 'CartesianGridXY') + glass_seg.set('grid_size_x', 'GlassModule_sx + GlassModule_wrap') + glass_seg.set('grid_size_y', 'GlassModule_sy + GlassModule_wrap') + rid = ET.SubElement(readout, 'id') + rid.text = 'system:8,sector:4,module:20,x:32:-16,y:-16' + + + text = ET.tostring(data, pretty_print=True) + with open(args.save, 'wb') as f: + f.write(text) + diff --git a/src/HomogeneousCalorimeter_geo.cpp b/src/HomogeneousCalorimeter_geo.cpp index f3d74c7f418d8a3121048ba3348b6119c98272fa..5ed339dc3dddce25f725a1981b45f243eb44ca42 100644 --- a/src/HomogeneousCalorimeter_geo.cpp +++ b/src/HomogeneousCalorimeter_geo.cpp @@ -1,13 +1,16 @@ //========================================================================== // A general implementation for homogeneous calorimeter // it supports three types of placements -// 1. Module placement with module dimensions and positions +// 1. Individual module placement with module dimensions and positions // 2. Array placement with module dimensions and numbers of row and columns // 3. Disk placement with module dimensions and (Rmin, Rmax), and (Phimin, Phimax) +// 4. Lines placement with module dimensions and (mirrorx, mirrory) +// (NOTE: anchor point is the 0th block of the line instead of line center) //-------------------------------------------------------------------------- // Author: Chao Peng (ANL) // Date: 06/09/2021 //========================================================================== + #include "GeometryHelpers.h" #include "DD4hep/DetFactoryHelper.h" #include <XML/Helper.h> @@ -28,8 +31,7 @@ using namespace dd4hep::detail; * * * \code - * <detector id="1" name="HyCal" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis"> - * <dimensions shape="box" sizex="120*cm" sizey="120*cm" sizez="46*cm"/> + * <detector id="1" name="HyCal" type="HomogeneousCalorimeter" readout="EcalHits"> * <position x="0" y="0" z="0"/> * <rotation x="0" y="0" z="0"/> * <placements> @@ -65,24 +67,22 @@ using namespace dd4hep::detail; * </placements> * </detector> * - * <detector id="2" name="SomeBlocks" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis"> - * <dimensions shape="box" sizex="100*cm" sizey="100*cm" sizez="20.5*cm"/> + * <detector id="2" name="SomeBlocks" type="HomogeneousCalorimeter" readout="EcalHits"> * <position x="0" y="0" z="30*cm"/> * <rotation x="0" y="0" z="0"/> * <placements> - * <blocks sector="1"/> + * <individuals sector="1"/> * <module sizex="2.05*cm" sizey="2.05*cm" sizez="20*cm" vis="GreenVis" material="PbWO4"/> * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/> - * <placement x="1*cm" y="1*cm" z="0"/> - * <placement x="-1*cm" y="1*cm" z="0"/> - * <placement x="1*cm" y="-1*cm" z="0"/> - * <placement x="-1*cm" y="-1*cm" z="0"/> - * </blocks> + * <placement x="1*cm" y="1*cm" z="0" id="1"/> + * <placement x="-1*cm" y="1*cm" z="0" id="2"/> + * <placement x="1*cm" y="-1*cm" z="0" id="3"/> + * <placement x="-1*cm" y="-1*cm" z="0" id="4"/> + * </individuals> * </placements> * </detector> * - * <detector id="2" name="DiskShapeCalorimeter" type="HomogeneousCalorimeter" readout="EcalHits" vis="GreenVis"> - * <dimensions shape="disk" rmin="25*cm" rmax="125*cm" length="20.5*cm" phimin="0" phimax="360*degree"/> + * <detector id="2" name="DiskShapeCalorimeter" type="HomogeneousCalorimeter" readout="EcalHits"> * <position x="0" y="0" z="-30*cm"/> * <rotation x="0" y="0" z="0"/> * <placements> @@ -91,15 +91,30 @@ using namespace dd4hep::detail; * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/> * </placements> * </detector> + * + * <detector id="3" name="SomeLines" type="HomogeneousCalorimeter" readout="EcalHits"> + * <position x="0" y="0" z="60*cm"/> + * <rotation x="0" y="0" z="0"/> + * <placements> + * <lines sector="1" mirrorx="true" mirrory="true"/> + * <module sizex="2.05*cm" sizey="2.05*cm" sizez="20*cm" vis="GreenVis" material="PbWO4"/> + * <wrapper thickness="0.015*cm" material="Epoxy" vis="WhiteVis"/> + * <line x="10.25*mm" y="10.25*mm" begin="8" nmods="16"/> + * <line x="10.25*mm" y="30.75*mm" begin="8" nmods="16"/> + * <line x="10.25*mm" y="51.25*mm" begin="8" nmods="16"/> + * </individuals> + * </placements> + * </detector> * \endcode * * @{ */ // headers -static void add_blocks(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); -static void add_array(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); -static void add_disk(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); +static void add_individuals(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); +static void add_array(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); +static void add_disk(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); +static void add_lines(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int id); // helper function to get x, y, z if defined in a xml component template<class XmlComp> @@ -124,61 +139,23 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete int detID = detElem.id(); DetElement det(detName, detID); sens.setType("calorimeter"); - - // top-level children - xml::Component dims = detElem.dimensions(); - - // build envelop from dimensions - std::string shape = dd4hep::getAttrOrDefault(dims, _Unicode(shape), ""); - // no shape input, try to determine shape from dimension variables - if (shape.empty()) { - if (dims.hasAttr(_Unicode(rmin)) && dims.hasAttr(_Unicode(rmax)) && dims.hasAttr(_Unicode(length))) { - shape = "disk"; - } else if (dims.hasAttr(_Unicode(sizex)) && dims.hasAttr(_Unicode(sizey)) && dims.hasAttr(_Unicode(sizez))) { - shape = "box"; - } else { - std::cerr << func << " Error: Cannot determine shape of the calorimeter. " - "Add shape (box, or disk) into dimensions\n"; - return det; - } - } - - Volume envVol(detName + "_envelope"); - std::string filler = dd4hep::getAttrOrDefault(detElem, _Unicode(filler), "Air"); - envVol.setMaterial(desc.material(filler)); - envVol.setVisAttributes(desc.visAttributes(detElem.visStr())); - - // convert to lower case - std::transform(shape.begin(), shape.end(), shape.begin(), [](unsigned char c){ return std::tolower(c); }); - - if (shape == "box") { - double sx = dims.attr<double>(_Unicode(sizex)); - double sy = dims.attr<double>(_Unicode(sizey)); - double sz = dims.attr<double>(_Unicode(sizez)); - envVol.setSolid(Box(sx/2., sy/2., sz/2.)); - } else if (shape == "disk") { - double rmin = dims.rmin(); - double rmax = dims.rmax(); - double length = dims.length(); - double phimin = dd4hep::getAttrOrDefault<double>(dims, _Unicode(phimin), 0.); - double phimax = dd4hep::getAttrOrDefault<double>(dims, _Unicode(phimax), 2*M_PI); - envVol.setSolid(Tube(rmin, rmax, length/2., phimin, phimax)); - } else { - std::cerr << func << " Error: Unsupported shape " << shape << ", use box or disk\n"; - return det; - } + // envelope + Assembly assembly(detName); // module placement xml::Component plm = detElem.child(_Unicode(placements)); int sector = 1; - for (xml::Collection_t arr(plm, _Unicode(array)); arr; ++arr) { - add_array(desc, envVol, arr, sens, sector++); + for (xml::Collection_t mod(plm, _Unicode(individuals)); mod; ++mod) { + add_individuals(desc, assembly, mod, sens, sector++); } - for (xml::Collection_t mod(plm, _Unicode(blocks)); mod; ++mod) { - add_blocks(desc, envVol, mod, sens, sector++); + for (xml::Collection_t arr(plm, _Unicode(array)); arr; ++arr) { + add_array(desc, assembly, arr, sens, sector++); } for (xml::Collection_t disk(plm, _Unicode(disk)); disk; ++disk) { - add_disk(desc, envVol, disk, sens, sector++); + add_disk(desc, assembly, disk, sens, sector++); + } + for (xml::Collection_t lines(plm, _Unicode(lines)); lines; ++lines) { + add_lines(desc, assembly, lines, sens, sector++); } // detector position and rotation @@ -186,7 +163,7 @@ static Ref_t create_detector(Detector& desc, xml::Handle_t handle, SensitiveDete auto rot = get_xml_xyz(detElem, _Unicode(rotation)); Volume motherVol = desc.pickMotherVolume(det); Transform3D tr = Translation3D(pos.x(), pos.y(), pos.z()) * RotationZYX(rot.z(), rot.y(), rot.x()); - PlacedVolume envPV = motherVol.placeVolume(envVol, tr); + PlacedVolume envPV = motherVol.placeVolume(assembly, tr); envPV.addPhysVolID("system", detID); det.setPlacement(envPV); return det; @@ -213,6 +190,9 @@ Volume build_module(Detector &desc, xml::Collection_t &plm, SensitiveDetector &s } else { auto wrp = plm.child(_Unicode(wrapper)); auto thickness = wrp.attr<double>(_Unicode(thickness)); + if (thickness == 0.) { + return modVol; + } auto wrpMat = desc.material(wrp.attr<std::string>(_Unicode(material))); Box wrpShape((sx + thickness)/2., (sy + thickness)/2., sz/2.); Volume wrpVol("wrapper_vol", wrpShape, wrpMat); @@ -223,12 +203,35 @@ Volume build_module(Detector &desc, xml::Collection_t &plm, SensitiveDetector &s } } +// place modules, id must be provided +static void add_individuals(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) +{ + Position modSize; + auto modVol = build_module(desc, plm, sens, modSize); + int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid); + + for (xml::Collection_t pl(plm, _Unicode(placement)); pl; ++pl) { + Position pos(dd4hep::getAttrOrDefault<double>(pl, _Unicode(x), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(y), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(z), 0.)); + Position rot(dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotx), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(roty), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotz), 0.)); + auto mid = pl.attr<int>(_Unicode(id)); + Transform3D tr = Translation3D(pos.x(), pos.y(), pos.z()) + * RotationZYX(rot.z(), rot.y(), rot.x()); + auto modPV = env.placeVolume(modVol, tr); + modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", mid); + } +} + // place array of modules -static void add_array(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) +static void add_array(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) { Position modSize; auto modVol = build_module(desc, plm, sens, modSize); int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid); + int id_begin = dd4hep::getAttrOrDefault<int>(plm, _Unicode(id_begin), 1); int nrow = plm.attr<int>(_Unicode(nrow)); int ncol = plm.attr<int>(_Unicode(ncol)); @@ -252,42 +255,21 @@ static void add_array(Detector& desc, Volume &env, xml::Collection_t &plm, Sensi } double px = begx + modSize.x()*j; double py = begy - modSize.y()*i; - Transform3D tr = Translation3D(pos.x() + px, pos.y() + py, pos.z()) - * RotationZYX(rot.z(), rot.y(), rot.x()); + Transform3D tr = RotationZYX(rot.z(), rot.y(), rot.x()) + * Translation3D(pos.x() + px, pos.y() + py, pos.z()); auto modPV = env.placeVolume(modVol, tr); - modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", i*ncol + j); + modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", i*ncol + j + id_begin); } } } -// place modules -static void add_blocks(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) -{ - Position modSize; - auto modVol = build_module(desc, plm, sens, modSize); - int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid); - - int mid = 1; - for (xml::Collection_t pl(plm, _Unicode(placement)); pl; ++pl, ++mid) { - Position pos(dd4hep::getAttrOrDefault<double>(pl, _Unicode(x), 0.), - dd4hep::getAttrOrDefault<double>(pl, _Unicode(y), 0.), - dd4hep::getAttrOrDefault<double>(pl, _Unicode(z), 0.)); - Position rot(dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotx), 0.), - dd4hep::getAttrOrDefault<double>(pl, _Unicode(roty), 0.), - dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotz), 0.)); - Transform3D tr = Translation3D(pos.x(), pos.y(), pos.z()) - * RotationZYX(rot.z(), rot.y(), rot.x()); - auto modPV = env.placeVolume(modVol, tr); - modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", mid); - } -} - // place disk of modules -static void add_disk(Detector& desc, Volume &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) +static void add_disk(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) { Position modSize; auto modVol = build_module(desc, plm, sens, modSize); int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid); + int id_begin = dd4hep::getAttrOrDefault<int>(plm, _Unicode(id_begin), 1); double rmin = plm.attr<double>(_Unicode(rmin)); double rmax = plm.attr<double>(_Unicode(rmax)); double phimin = dd4hep::getAttrOrDefault<double>(plm, _Unicode(phimin), 0.); @@ -297,12 +279,57 @@ static void add_disk(Detector& desc, Volume &env, xml::Collection_t &plm, Sensit // placement to mother auto pos = get_xml_xyz(plm, _Unicode(position)); auto rot = get_xml_xyz(plm, _Unicode(rotation)); - int mid = 1; + int mid = 0; for (auto &p : points) { - Transform3D tr = Translation3D(pos.x() + p.x(), pos.y() + p.y(), pos.z()) - * RotationZYX(rot.z(), rot.y(), rot.x()); + Transform3D tr = RotationZYX(rot.z(), rot.y(), rot.x()) + * Translation3D(pos.x() + p.x(), pos.y() + p.y(), pos.z()); auto modPV = env.placeVolume(modVol, tr); - modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", mid++); + modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", id_begin + mid++); + } +} + +// place lines of modules (anchor point is the 0th module of this line) +static void add_lines(Detector& desc, Assembly &env, xml::Collection_t &plm, SensitiveDetector &sens, int sid) +{ + Position modSize; + auto modVol = build_module(desc, plm, sens, modSize); + int sector_id = dd4hep::getAttrOrDefault<int>(plm, _Unicode(sector), sid); + int id_begin = dd4hep::getAttrOrDefault<int>(plm, _Unicode(id_begin), 1); + bool mirrorx = dd4hep::getAttrOrDefault<bool>(plm, _Unicode(mirrorx), false); + bool mirrory = dd4hep::getAttrOrDefault<bool>(plm, _Unicode(mirrory), false); + + // line placement + int mid = 1; + for (xml::Collection_t pl(plm, _Unicode(line)); pl; ++pl) { + Position pos(dd4hep::getAttrOrDefault<double>(pl, _Unicode(x), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(y), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(z), 0.)); + Position rot(dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotx), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(roty), 0.), + dd4hep::getAttrOrDefault<double>(pl, _Unicode(rotz), 0.)); + int begin = dd4hep::getAttrOrDefault<int>(pl, _Unicode(begin), 0); + int nmods = pl.attr<int>(_Unicode(nmods)); + + std::vector<std::pair<double, double>> translations; + for (int i = 0; i < nmods; ++i) { + translations.push_back(std::pair<double, double>{pos.x() + (begin + i)*modSize.x(), pos.y()}); + if (mirrorx) { + translations.push_back(std::pair<double, double>{-pos.x() - (begin + i)*modSize.x(), pos.y()}); + } + if (mirrory) { + translations.push_back(std::pair<double, double>{pos.x() + (begin + i)*modSize.x(), -pos.y()}); + } + if (mirrorx && mirrory) { + translations.push_back(std::pair<double, double>{-pos.x() - (begin + i)*modSize.x(), -pos.y()}); + } + } + + for (auto &p : translations) { + Transform3D tr = RotationZYX(rot.z(), rot.y(), rot.x()) + * Translation3D(p.first, p.second, pos.z()); + auto modPV = env.placeVolume(modVol, tr); + modPV.addPhysVolID("sector", sector_id).addPhysVolID("module", id_begin + mid++); + } } } //@}