Skip to content

Commit

Permalink
Update FiberTubeCalorimeter example
Browse files Browse the repository at this point in the history
  • Loading branch information
MarkusFrankATcernch committed Jan 18, 2024
1 parent c0120cd commit 7994dd9
Show file tree
Hide file tree
Showing 3 changed files with 137 additions and 29 deletions.
15 changes: 8 additions & 7 deletions examples/ClientTests/compact/FiberTubeCalorimeter.xml
Original file line number Diff line number Diff line change
Expand Up @@ -70,11 +70,12 @@
</define>

<display>
<vis name="AbsVis" alpha="1.0" r="0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="AbsVis" alpha="0.5" r="0.0" g="1.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="CerenVis" alpha="0.5" r="0.0" g="0.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="ScintVis" alpha="0.5" r="1.0" g="0.0" b="0.0" showDaughters="true" visible="true"/>
<vis name="phdetVis" alpha="1.0" r="1.0" g="1.0" b="1.0" showDaughters="true" visible="true"/>
<vis name="towerVis" alpha="0.1" r="0.3" g="0.3" b="0.3" showDaughters="true" visible="true"/>
<vis name="layerVis" alpha="0.5" r="0.1" g="0.1" b="0.9" showDaughters="true" visible="true"/>
</display>

<readouts>
Expand All @@ -92,18 +93,18 @@
zmin is where the dtector goes
-->
<dimensions numsides="DRFiberNsizeSCE"
thickness="DRFiberAbswidthSCE"
thickness="DRFiberAbswidthSCE+holeoverSCE"
z_length="DRFiberlengthSCE"
gap="gapSCE"
zmin="-world_side/2.+2*killthick+edgeoffset+DRcrystallength+EcalHcalgap"
z1="killthick"/>
<structure>
<core1 name="scintillator" rmax="DRFiberFibwidthSCE" rmin="0.0" material="Polystyrene" vis="ScintVis" sensitive="yes"/>
<core2 name="quartz" rmax="DRFiberFibwidthSCE" rmin="0.0" material="Quartz" vis="CerenVis" sensitive="yes"/>
<hole name="hole" rmax="(DRFiberFibwidthSCE+holeoverSCE)" rmin="0.0" material="Air" vis="holeVis" sensitive="yes"/>
<absorb material="Brass" vis="AbsVis" sensitive="yes"/>
<phdet1 name="phdet1" rmax="DRFiberFibwidthSCE" rmin="0.0" material="killMedia2" vis="phdetVis" sensitive="yes"/>
<phdet2 name="phdet2" rmax="DRFiberFibwidthSCE" rmin="0.0" material="killMedia3" vis="phdetVis" sensitive="yes"/>
<core2 name="quartz" rmax="DRFiberFibwidthSCE" rmin="0.0" material="Quartz" vis="CerenVis" sensitive="yes"/>
<hole name="hole" rmax="DRFiberFibwidthSCE+holeoverSCE" rmin="DRFiberFibwidthSCE" material="Air" vis="holeVis" sensitive="yes"/>
<absorb material="Brass" vis="AbsVis" sensitive="yes"/>
<phdet1 name="phdet1" rmax="DRFiberFibwidthSCE" rmin="0.0" material="killMedia2" vis="phdetVis" sensitive="yes"/>
<phdet2 name="phdet2" rmax="DRFiberFibwidthSCE" rmin="0.0" material="killMedia3" vis="phdetVis" sensitive="yes"/>
</structure>
</detector>
</detectors>
Expand Down
95 changes: 95 additions & 0 deletions examples/ClientTests/scripts/FiberTubeCalorimeter.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,95 @@
# ==========================================================================
# AIDA Detector description implementation
# --------------------------------------------------------------------------
# Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
# All rights reserved.
#
# For the licensing terms see $DD4hepINSTALL/LICENSE.
# For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
#
# ==========================================================================
#
#
from __future__ import absolute_import, unicode_literals
import os
import time
import DDG4
from DDG4 import OutputLevel as Output
from g4units import GeV, MeV, m
#
#
"""
dd4hep simulation example setup using the python configuration
@author M.Frank
@version 1.0
"""


def run():
args = DDG4.CommandLine()
kernel = DDG4.Kernel()
install_dir = os.environ['DD4hepExamplesINSTALL']
kernel.loadGeometry(str("file:" + install_dir + "/examples/ClientTests/compact/FiberTubeCalorimeter.xml"))

DDG4.importConstants(kernel.detectorDescription(), debug=False)
geant4 = DDG4.Geant4(kernel, tracker='Geant4TrackerCombineAction')
geant4.printDetectors()
# Configure UI
if args.macro:
ui = geant4.setupCshUI(macro=args.macro)
else:
ui = geant4.setupCshUI()
if args.batch:
ui.Commands = ['/run/beamOn ' + str(args.events), '/ddg4/UI/terminate']

# Configure field
geant4.setupTrackingField(prt=True)
# Configure Event actions
prt = DDG4.EventAction(kernel, 'Geant4ParticlePrint/ParticlePrint')
prt.OutputLevel = Output.DEBUG
prt.OutputType = 3 # Print both: table and tree
kernel.eventAction().adopt(prt)

generator_output_level = Output.INFO

# Configure G4 geometry setup
seq, act = geant4.addDetectorConstruction("Geant4DetectorGeometryConstruction/ConstructGeo")
act.DebugMaterials = True
act.DebugElements = False
act.DebugVolumes = True
act.DebugShapes = True
seq, act = geant4.addDetectorConstruction("Geant4DetectorSensitivesConstruction/ConstructSD")

# Configure I/O
geant4.setupROOTOutput('RootOutput', 'FiberTubeCalorimeter_' + time.strftime('%Y-%m-%d_%H-%M'))

# Setup particle gun
gun = geant4.setupGun("Gun", particle='e+', energy=20 * GeV, multiplicity=1, position=(0.0, 0.0, -369.0))
gun.OutputLevel = generator_output_level

# And handle the simulation particles.
part = DDG4.GeneratorAction(kernel, "Geant4ParticleHandler/ParticleHandler")
kernel.generatorAction().adopt(part)
part.SaveProcesses = ['Decay']
part.MinimalKineticEnergy = 100 * MeV
part.OutputLevel = Output.INFO # generator_output_level
part.enableUI()
user = DDG4.Action(kernel, "Geant4TCUserParticleHandler/UserParticleHandler")
user.TrackingVolume_Zmax = 3.0 * m
user.TrackingVolume_Rmax = 3.0 * m
user.enableUI()
part.adopt(user)

geant4.setupCalorimeter('FiberTubeCalorimeter')

# Now build the physics list:
phys = geant4.setupPhysics('QGSP_BERT')
phys.dump()
geant4.execute()


if __name__ == "__main__":
run()
56 changes: 34 additions & 22 deletions examples/ClientTests/src/FiberTubeCalorimeter_geo.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -68,27 +68,17 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector s
// these refer to different fields in the xml file for this detector
xml_comp_t fX_struct( x_det.child( _Unicode(structure) ) );
xml_comp_t fX_absorb( fX_struct.child( _Unicode(absorb) ) );
xml_comp_t fX_core1( fX_struct.child( _Unicode(core1) ) );
xml_comp_t fX_core2( fX_struct.child( _Unicode(core2) ) );
xml_comp_t fX_hole( fX_struct.child( _Unicode(hole) ) );
xml_comp_t fX_core1( fX_struct.child( _Unicode(core1) ) );
xml_comp_t fX_core2( fX_struct.child( _Unicode(core2) ) );
xml_comp_t fX_hole( fX_struct.child( _Unicode(hole) ) );
xml_comp_t fX_phdet1( fX_struct.child( _Unicode(phdet1) ) );
xml_comp_t fX_phdet2( fX_struct.child( _Unicode(phdet2) ) );

// detector element for entire detector.
DetElement sdet (det_name, det_id);
Volume motherVol = description.pickMotherVolume(sdet);
Box env_box ((2*Ncount+1)*(hthick+agap+tol),(2*Ncount+1)*(hthick+agap+tol), (hzlength+hzph+tol));
Volume envelopeVol (det_name, env_box, air);
envelopeVol.setAttributes(description,x_det.regionStr(),x_det.limitsStr(),x_det.visStr());

Material mat;
Transform3D trafo;
PlacedVolume pv;
Solid sol;

pv = motherVol.placeVolume(envelopeVol, Position(0.,0.,azmin+hzlength+hzph+tol));
pv.addPhysVolID("system", det_id);
sdet.setPlacement(pv); // associate the placed volume to the detector element

sens.setType("calorimeter");

// scint fiber
Expand Down Expand Up @@ -182,24 +172,40 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector s
}

// setup the volumes with the shapes and properties in one horixontal layer
Volume tube_row_vol("layer", Box(hthick,hthick,hzlength+hzph), air);
double dx = 2*(Ncount + Ncount+1)/2e0 * (hthick+agap) + tol;
double dy = hthick + tol;
double dz = hzlength+hzph + tol;
Box tube_row_box(dx, dy, dz);
Volume tube_row_vol("layer", tube_row_box, air);
tube_row_vol.setVisAttributes(description, x_det.visStr());
tube_row_vol.setSensitiveDetector(sens);

cout << tube_row_vol.name()
<< " dx: " << tube_row_box.x()
<< " dy: " << tube_row_box.y()
<< " dz: " << tube_row_box.z() << endl;
tube_row_vol.setVisAttributes(description, "layerVis");

for (int ijk=-Ncount; ijk<Ncount+1; ijk++) {
double mod_x_off = (ijk)*2*(hthick+agap);
Transform3D tr(RotationZYX(0.,0.,0.), Position(mod_x_off,0.,0.));
double mod_x_off = (ijk) * 2 * (hthick + agap);
int towernum = Ncount + ijk + 1;
pv = tube_row_vol.placeVolume((towernum%2 == 0) ? quartz_abs_vol : scint_abs_vol, tr);
pv = tube_row_vol.placeVolume((towernum%2 == 0) ? quartz_abs_vol : scint_abs_vol, Position(mod_x_off,0.,0.));
pv.addPhysVolID("tube", towernum);
cout << "Placing row " << setw(5) << right << ijk
<< " x-offset: " << setw(7) << right<< mod_x_off
//Box bounding_box = pv.volume().solid().GetBoundingBox();
cout << "Placing row " << setw(5) << right << ijk
<< " x-offset: " << setw(7) << right << mod_x_off
<< " volume of type " << pv.volume().name()
<< endl;
}

dy = 2*(Ncount + Ncount+1)/2e0 * (hthick+agap) + tol;
DetElement sdet (det_name, det_id);
Box env_box (dx+tol, dy+tol, dz+tol);
Volume envelopeVol (det_name, env_box, air);
envelopeVol.setAttributes(description, x_det.regionStr(), x_det.limitsStr(), x_det.visStr());

// Now stack multiple horizontal layers to form the final box
for (int ijk=-Ncount; ijk<Ncount+1; ijk++) {
double mod_y_off = (ijk)*2*(hthick+agap);
double mod_y_off = (ijk) * 2 * (tube_row_box.y() + agap);
Transform3D tr(RotationZYX(0.,0.,0.), Position(0.,mod_y_off,0.));
pv = envelopeVol.placeVolume(tube_row_vol, tr);
pv.addPhysVolID("layer", Ncount+ijk+1);
Expand All @@ -212,6 +218,12 @@ static Ref_t create_detector(Detector& description, xml_h e, SensitiveDetector s
<< " volume of type " << pv.volume().name()
<< endl;
}

// detector element for entire detector.
Volume motherVol = description.pickMotherVolume(sdet);
pv = motherVol.placeVolume(envelopeVol, Position(0.,0.,azmin+hzlength+hzph+tol));
pv.addPhysVolID("system", det_id);
sdet.setPlacement(pv); // associate the placed volume to the detector element
return sdet;
}

Expand Down

0 comments on commit 7994dd9

Please sign in to comment.