Commit 8e22aed6 authored by Alexander Kiselev's avatar Alexander Kiselev
Browse files

SiPM sensor geometry for dRICH implemented

parent bde09ca5
......@@ -378,7 +378,7 @@ StatusCode Jug::PID::IRTAlgorithm::execute( void )
double eVenergy = 1E9*hit.energy();
//printf("%f\n", eVenergy);
if (!QE_pass(eVenergy, m_rngUni()) ||
m_rngUni() > m_GeometricEfficiency.value()*m_SafetyFactor.value()) {
m_rngUni() > /*m_GeometricEfficiency.value()**/m_SafetyFactor.value()) {
// Prefer to create the stepping history in a clean way, namely use only those
// photons which were rejected; FIXME: optionally include all; anyway, now
......@@ -396,20 +396,45 @@ StatusCode Jug::PID::IRTAlgorithm::execute( void )
//if (vs.z < 2500 || vs.z > 2700) continue;
//if (vs.z < 2200 || vs.z > 2400) continue;
useful_photons_found = true;
// Photon accepted; add it to the internal event structure;
auto photon = new OpticalPhoton();
// Hit location after pixelization;
TVector3 lxyPixel;
// FIXME: '0' - for now assume a single photon detector type for the whole ERICH (DRICH),
// which must be a reasonable assumption; eventually may want to generalize;
const auto pd = m_IrtDet->m_PhotonDetectors[0];
photon->SetPhotonDetector(pd);
uint64_t vcopy = hit.cellID() & m_ReadoutCellMask;
{
uint64_t vcopy = hit.cellID() & m_ReadoutCellMask;
const auto irt = pd->GetIRT(vcopy);
auto sensor = dynamic_cast<const FlatSurface*>(irt->tail()->GetSurface());
photon->SetVolumeCopy(vcopy);
double pitch = m_SensorPixelPitch.value(), sens = m_SensorPixelSize.value();//pitch - gap;
const auto &x = hit.position();
TVector3 lpt(x.x, x.y, x.z);
double lxy[2] = {sensor->GetLocalX(lpt), sensor->GetLocalY(lpt)}, lxyPixels[2] = {0.0};
//printf("%7.2f %7.2f\n", lxy[0], lxy[1]);
int ixy[2];
for(unsigned iq=0; iq<2; iq++) {
ixy[iq] = (int)floor(lxy[iq]/pitch);//m_SensorPixelPitch.value());
// "+0.5" assumes even number of pixels in each projection;
lxyPixels[iq] = pitch*(ixy[iq] + 0.5);
} //for iq
// Well, just ignore this small fraction of unlucky photons (do not even use them for track pinning);
if (fabs(lxy[0] - lxyPixels[0]) > sens/2 || fabs(lxy[1] - lxyPixels[1]) > sens/2) continue;
lxyPixel = sensor->GetSpacePoint(lxyPixels[0], lxyPixels[1]);
}
useful_photons_found = true;
// Photon accepted; add it to the internal event structure;
auto photon = new OpticalPhoton();
photon->SetVolumeCopy(vcopy);
photon->SetDetectionPosition(lxyPixel);
photon->SetPhotonDetector(pd);
{
// Start vertex and momentum;
photon->SetVertexPosition(vtx);
{
......@@ -427,15 +452,6 @@ StatusCode Jug::PID::IRTAlgorithm::execute( void )
photon->SetVertexRefractiveIndex(ri);
}
}
// Hit location;
{
// FIXME: take sensor plane orientation into account, and also use floor()
// rather than gaussian smearing;
double sigma = 3.4/sqrt(12.0);
const auto &x = hit.position();
photon->SetDetectionPosition(TVector3(x.x + sigma*m_rngUni(), x.y + sigma*m_rngUni(), x.z));
}
// At this point all of the CherenkovPhoton class internal variables are actually
// set the same way as in a standalone G4 code, except for the parent momentum at vertex;
......@@ -486,6 +502,7 @@ StatusCode Jug::PID::IRTAlgorithm::execute( void )
if (b1 && b2) {
const unsigned zbins = radiator->GetTrajectoryBinCount();
TVector3 nn = (to - from).Unit();
// FIXME: hardcoded;
from += (0.010)*nn;
to -= (0.010)*nn;
......
......@@ -103,9 +103,11 @@ namespace Jug::PID {
Gaudi::Property<std::string> m_Detector {this, "Detector", ""};
// The famous ~0.85 for the S13660-3050AE-08 SiPM and its QE table expected; well, as long
// as the fill factor is actually mimiced in E(D)Rich_geo.cpp, geometric efficiency should be 1.0 here;
Gaudi::Property<double> m_GeometricEfficiency {this, "GeometricEfficiency", 1.0};
//Gaudi::Property<double> m_GeometricEfficiency {this, "GeometricEfficiency", 1.0};
// Another famous factor (life is always below expectations) which people like to see ~0.7;
Gaudi::Property<double> m_SafetyFactor {this, "SafetyFactor", 1.0};
Gaudi::Property<double> m_SensorPixelSize {this, "SensorPixelSize", 3.0};
Gaudi::Property<double> m_SensorPixelPitch {this, "SensorPixelPitch", 3.2};
// Array of {e, qe} points and the desired number of bins in a fast lookup table;
Gaudi::Property<std::vector<std::pair<double, double>>> m_QE_input_data {this, "QEcurve"};
Gaudi::Property<unsigned> m_QEbins {this, "QEbins", 0};
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment