LDMX-Software
Updates to the hcal prototype geometry
The current ldmx-hcal-prototype-v1 geometry differs from what was actually used at the test beam so we need to add a couple of new GDML descriptions (corresponding to different configurations that were used) which also brings along some necessary updates to LDMX-sw (mostly for the Hcal submodule but also in DetDescr).
- The ordering of horizontal and vertical layers in the prototype differs from the regular Hcal
- Some measurements were carried out with the horizontal bars shifted which means that the center of the bars cannot be assumed to be at x/y = 0
@bryngemark I've pushed a version which has @PeterGy's latest work on https://github.com/LDMX-Software/ldmx-sw/tree/iss1063 in the Detectors/data/ldmx-hcal-prototype-v2.0 directory. There are some different versions of the TS available in the trigger_physvols file, but everything except plastic should be commented out by default.
from LDMX.Framework import ldmxcfg
from LDMX.SimCore import generators
from LDMX.SimCore import simulator
process = ldmxcfg.Process('test_TS')
process.maxEvents = 1000
process.run = 10
process.outputFiles = ['test_TS.root']
gun = generators.gun('particle_gun')
gun.particle = 'e-'
gun.direction = [0., 0., 1.]
gun.position = [0., 0., -1000.]
gun.energy = 4.
simulation = simulator.simulator('test_TS')
simulation.generators=[gun]
simulation.setDetector('ldmx-hcal-prototype-v2.0')
from LDMX.TrigScint.trigScint import TrigScintQIEDigiProducer
from LDMX.TrigScint.trigScint import TrigScintRecHitProducer
tsDigis = TrigScintQIEDigiProducer.up()
tsDigis.number_of_strips = 12
tsRecHits = TrigScintRecHitProducer.up()
process.sequence=[simulation,
tsDigis,
tsRecHits
]
I stole this from what @PeterGy has been using and it does ~something~ although you and Peter are probably better than me at knowing if it is sensible :)
I've uploaded the basic configuration (plastic TS, no shifted bars) of the April testbeam to https://github.com/LDMX-Software/ldmx-sw/tree/iss1063 together with an updated DetDescr readout geometry. See updated test config below.
However, if you try to run this by default you will crash when running the HcalRecProducer. This is because the ordering of the layers has changed. In the April testbeam, odd layers were vertical while in the standard geometries odd layers in the Hcal are horizontal. To fix this, there are two places (at least) that need to be able to take this into account.
In https://github.com/LDMX-Software/SimCore/blob/1fd615f297b6daad5a80ba41e29596f1c64cf633/src/SimCore/HcalSD.cxx#L136-L146
int stripID = -1;
// Odd layers have bars horizontal
// Even layers have bars vertical
// 5cm wide bars are HARD-CODED
if (section == ldmx::HcalID::BACK && layer % 2 == 1)
stripID = int((localPosition.y() + scint->GetYHalfLength()) / 50.0);
else if (section == ldmx::HcalID::BACK && layer % 2 == 0)
stripID = int((localPosition.x() + scint->GetXHalfLength()) / 50.0);
else
stripID = int((localPosition.z() + scint->GetZHalfLength()) / 50.0);
in https://github.com/LDMX-Software/Hcal/blob/49ee21d7634d353ee538470aaa4d89d1d2cebccf/src/Hcal/HcalRecProducer.cxx#L254-L259
// set position along the bar
if ((id_posend.layer() % 2) == 1) {
position.SetX(position_bar);
} else {
position.SetY(position_bar);
}
in https://github.com/LDMX-Software/ldmx-sw/blob/f49a7f8daae371425a538d629dc895944465de34/DetDescr/src/DetDescr/HcalGeometry.cxx#L110-L116
/**
Now compute, y(x) position for even(odd) layers, relative to the
center of detector. Strips enumeration starts from -y(-x)
stripcenter will be large for +y(+x) and the halfwidth of the strip
needs to be subtracted The halfwidth of the scintillator is given by
ZeroStrip_. The x(y) position is set to the center of the strip (0).
*/
if (layer % 2 == 1) {
y = stripcenter - getZeroStrip(section, layer);
x = 0;
} else {
x = stripcenter - getZeroStrip(section, layer);
y = 0;
}
and in https://github.com/LDMX-Software/Hcal/blob/49ee21d7634d353ee538470aaa4d89d1d2cebccf/src/Hcal/HcalDigiProducer.cxx#L171
distance_along_bar = (layer % 2) ? position[0] : position[1];
If you want to start testing running simulations with what's currently in this branch, you can just manually patch these three places and you should (probably) be fine to start testing things.
The easiest way to support these geometry differences that I can think of would be to rely on the HcalGeometry condition as the single source of truth. Then you could just have a member function in https://github.com/LDMX-Software/ldmx-sw/blob/trunk/DetDescr/include/DetDescr/HcalGeometry.h along the lines of (although maybe with some more elegant but equivalent expression that I couldn't think of right now)
bool layerIsHorizontal(const int layer) const {
if (someRuntimeConfig) {
return layer % 2 == 0;
} else {
return layer % 2 == 1;
}
}
then the places that currently has the layer ordering hard coded could rely on the geometry condition and do something along the lines of
// set position along the bar
if (hcalGeometry.layerIsHorizontal(id_posend.layer())) {
position.SetX(position_bar);
} else {
position.SetY(position_bar);
}
For the Hcal parts, this trivial to do. However, this would add a new dependency for SimCore on import LDMX.Hcal.HcalGeometry. Maybe that's ok? There's already a dependence on the corresponding the Ecal geometry in EcalSD.
from LDMX.Framework import ldmxcfg
from LDMX.SimCore import generators
from LDMX.SimCore import simulator
process = ldmxcfg.Process('test')
process.maxEvents = 1000
process.run = 10
process.outputFiles = ['test.root']
gun = generators.gun('particle_gun')
gun.particle = 'e-'
gun.direction = [0., 0., 1.]
gun.position = [0., 0., -1000.]
gun.energy = 4.
simulation = simulator.simulator('test_TS')
simulation.generators=[gun]
simulation.setDetector('ldmx-hcal-prototype-v2.0')
from LDMX.TrigScint.trigScint import TrigScintQIEDigiProducer
from LDMX.TrigScint.trigScint import TrigScintRecHitProducer
tsDigis = TrigScintQIEDigiProducer.up()
tsDigis.number_of_strips = 12
tsRecHits = TrigScintRecHitProducer.up()
import LDMX.Hcal.HcalGeometry
from LDMX.Hcal import digi
from LDMX.Hcal import hcal
from LDMX.Hcal import hcal_hardcoded_conditions
digis = digi.HcalDigiProducer()
reco = digi.HcalRecProducer()
process.sequence=[simulation,
digis,
reco,
tsDigis,
tsRecHits
]
I went ahead and tested doing things this way and it wasn't particularly complicated to do. So if it makes sense, I can go ahead and implement it.
Another option is to have the parity of the layers configurable, as in horizontalParity_=1, and then check,
if layer % 2 == horizontalParity_
...
Could still be set in the python script or the conditions but you wouldn't have to check for which configuration you're using.
Oh yeah that's much better. The question of putting it as part of the condition or as a part of the configuration is hard, i would lean towards putting it in the condition since, in principle at least, once you choose your geometry it would be a bug to select the wrong parameter. On the other hand though, having the partity change invisibly when you change detector geometry might be the kind of thing that makes someone pull their hair out if it isn't obvious
ldmx-hcal-prototype-v2 is the resolution of this issue, no? Let's close this issue
