larg4::LArVoxelReadout Class Reference

Transports energy depositions from GEANT4 to TPC channels. More...

#include <LArVoxelReadout.h>

Inheritance diagram for larg4::LArVoxelReadout:

Classes
struct	Setup_t
	Collection of what it takes to set a LArVoxelReadout up. More...

Public Types
typedef std::map< unsigned int, sim::SimChannel >	ChannelMap_t
	Type of map channel -> sim::SimChannel. More...

Public Member Functions
	LArVoxelReadout (std::string const &name)
	Constructor. Can detect which TPC to cover by the name. More...
	LArVoxelReadout (std::string const &name, unsigned int cryostat, unsigned int tpc)
	Constructor. Sets which TPC to work on. More...
void	Setup (Setup_t const &setupData)
	Reads all the configuration elements from setupData More...
void	SetSingleTPC (unsigned int cryostat, unsigned int tpc)
	Associates this readout to one specific TPC. More...
void	SetDiscoverTPC ()
	Sets this readout to discover the TPC of each processed hit. More...
virtual void	Initialize (G4HCofThisEvent *)
virtual void	EndOfEvent (G4HCofThisEvent *)
virtual void	clear ()
virtual G4bool	ProcessHits (G4Step *, G4TouchableHistory *)
virtual void	DrawAll ()
virtual void	PrintAll ()
void	ClearSimChannels ()
std::vector< sim::SimChannel >	GetSimChannels () const
	Creates a list with the accumulated information for the single TPC. More...
std::vector< sim::SimChannel >	GetSimChannels (unsigned short cryo, unsigned short tpc) const
	Creates a list with the accumulated information for specified TPC. More...
const ChannelMap_t &	GetSimChannelMap () const
	Returns the accumulated channel -> SimChannel map for the single TPC. More...
ChannelMap_t &	GetSimChannelMap ()
const ChannelMap_t &	GetSimChannelMap (unsigned short cryo, unsigned short tpc) const
	Returns the accumulated channel -> SimChannel map for the specified TPC. More...
ChannelMap_t &	GetSimChannelMap (unsigned short cryo, unsigned short tpc)

Private Member Functions
void	SetClockData (detinfo::DetectorClocksData const *const clockData) noexcept
void	SetPropertiesData (detinfo::DetectorPropertiesData const *const detProp) noexcept
void	SetOffPlaneChargeRecoveryMargin (double margin)
	Sets the margin for recovery of charge drifted off-plane. More...
void	SetRandomEngines (CLHEP::HepRandomEngine *pPropGen)
	Sets the random generators to be used. More...
geo::Point_t	RecoverOffPlaneDeposit (geo::Point_t const &pos, geo::PlaneGeo const &plane) const
	Returns the point on the specified plane closest to position. More...
void	DriftIonizationElectrons (detinfo::DetectorClocksData const &clockData, G4ThreeVector stepMidPoint, const double simTime, int trackID, unsigned short int cryostat, unsigned short int tpc, int origTrackID)
bool	Has (std::vector< unsigned short int > v, unsigned short int tpc) const
void	ProcessStep (G4Step *)

Private Attributes
double	fDriftVelocity [3]
double	fLongitudinalDiffusion
double	fTransverseDiffusion
double	fElectronLifetime
double	fElectronClusterSize
int	fMinNumberOfElCluster
bool	fDontDriftThem
std::vector< unsigned short int >	fSkipWireSignalInTPCs
double	fOffPlaneMargin = 0.0
	Charge deposited within this many [cm] from the plane is lead onto it. More...
std::vector< std::vector < ChannelMap_t > >	fChannelMaps
	Maps of cryostat, tpc to channel data. More...
art::ServiceHandle < geo::Geometry const >	fGeoHandle
	Handle to the Geometry service. More...
art::ServiceHandle < sim::LArG4Parameters const >	fLgpHandle
	Handle to the LArG4 parameters service. More...
unsigned int	fTPC
	which TPC this LArVoxelReadout corresponds to More...
unsigned int	fCstat
	and in which cryostat (if bSingleTPC is true) More...
bool	bSingleTPC
	true if this readout is associated with a single TPC More...
CLHEP::HepRandomEngine *	fPropGen = nullptr
	random engine for charge propagation More...
detinfo::DetectorClocksData const *	fClockData {nullptr}
detinfo::DetectorPropertiesData const *	fDetProp {nullptr}
G4ThreeVector	fStepStart
G4ThreeVector	fStepEnd
size_t	fNSteps

Friends
class	LArVoxelReadoutGeometry

Detailed Description

Transports energy depositions from GEANT4 to TPC channels.

This class acts on single energy depositions from GEANT4, simulating the transportation of the ensuing ionisation electrons to the readout channels:

1. the number of ionisation electrons is read from the current larg4::IonizationAndScintillation instance
2. space charge displacement is optionally applied
3. lifetime correction is applied
4. charge is split in small electron clusters
5. each cluster is subject to longitudinal and transverse diffusion
6. each cluster is assigned to one TPC channel for each wire plane
7. optionally, charge is forced to stay on the planes; otherwise charge drifting outside the plane is lost

For each energy deposition, entries on the appropriate sim::SimChannel are added, with the information of the position where the energy deposit happened (in global coordinates, centimeters), the ID of the Geant4 track which produced the deposition, and the quantized time of arrival to the channel (in global TDC tick units). At most one entry is added for each electron cluster, but entries from the same energy deposit can be compacted if falling on the same TDC tick.

The main entry point of this class is the method ProcessHits().

Options

A few optional behaviours are supported:

• lead off-plane charge to the planes: regulated by RecoverOffPlaneDeposit(), if charge which reaches a wire plane is actually off it by less than the chosen margin, it's accounted for by that plane; by default the margin is 0 and all the charge off the plane is lost (with a warning)

Definition at line 170 of file LArVoxelReadout.h.

Member Typedef Documentation

typedef std::map<unsigned int, sim::SimChannel> larg4::LArVoxelReadout::ChannelMap_t

Type of map channel -> sim::SimChannel.

Definition at line 173 of file LArVoxelReadout.h.

Constructor & Destructor Documentation

larg4::LArVoxelReadout::LArVoxelReadout

(

std::string const &

name

)

Constructor. Can detect which TPC to cover by the name.

Definition at line 49 of file LArVoxelReadout.cxx.

                                                        : G4VSensitiveDetector(name)
  {
    // Initialize values for the electron-cluster calculation.
    ClearSimChannels();

    // the standard name contains cryostat and TPC;
    // if we don't find it, we will detect the TPC at each Geant hit
    unsigned int cryostat, tpc;
    if (std::sscanf(name.c_str(), "%*19s%1u%*4s%u", &cryostat, &tpc) == 2)
      SetSingleTPC(cryostat, tpc);
    else
      SetDiscoverTPC();

  } // LArVoxelReadout::LArVoxelReadout()

larg4::LArVoxelReadout::LArVoxelReadout	(	std::string const &	name,
		unsigned int	cryostat,
		unsigned int	tpc
	)

Constructor. Sets which TPC to work on.

Definition at line 67 of file LArVoxelReadout.cxx.

    : LArVoxelReadout(name)
  {
    SetSingleTPC(cryostat, tpc);
  }

Member Function Documentation

void larg4::LArVoxelReadout::clear ( )

virtual

Definition at line 146 of file LArVoxelReadout.cxx.

  {}

void larg4::LArVoxelReadout::ClearSimChannels

(

)

Definition at line 151 of file LArVoxelReadout.cxx.

  {
    fChannelMaps.resize(fGeoHandle->Ncryostats());
    size_t cryo = 0;
    for (auto& cryoData : fChannelMaps) { // each, a vector of maps
      cryoData.resize(fGeoHandle->NTPC(cryo++));
      for (auto& channelsMap : cryoData)
        channelsMap.clear(); // each, a map
    }                        // for cryostats
  }                          // LArVoxelReadout::ClearSimChannels()

void larg4::LArVoxelReadout::DrawAll ( )

virtual

Definition at line 595 of file LArVoxelReadout.cxx.

  {}

void larg4::LArVoxelReadout::DriftIonizationElectrons ( detinfo::DetectorClocksData const &  clockData, G4ThreeVector  stepMidPoint, const double  simTime, int  trackID, unsigned short int  cryostat, unsigned short int  tpc, int  origTrackID  )

private

Todo:

think about effects of drift between planes

Todo:

think about effects of drift between planes

Todo:

check on what happens if we allow the tdc value to be

Todo:

beyond the end of the expected number of ticks

Definition at line 340 of file LArVoxelReadout.cxx.

  {
    auto const tpcClock = clockData.TPCClock();

    // this must be always true, unless caller has been sloppy
    assert(fPropGen); // No propagation random generator provided?!

    CLHEP::RandGauss PropRand(*fPropGen);

    // This routine gets called frequently, once per every particle
    // traveling through every voxel. Use whatever tricks we can to
    // increase its execution speed.

    static double LifetimeCorr_const = -1000. * fElectronLifetime;
    static double LDiff_const = std::sqrt(2. * fLongitudinalDiffusion);
    static double TDiff_const = std::sqrt(2. * fTransverseDiffusion);
    static double RecipDriftVel[3] = {
      1. / fDriftVelocity[0], 1. / fDriftVelocity[1], 1. / fDriftVelocity[2]};

    struct Deposit_t {
      double energy = 0.;
      double electrons = 0.;

      void
      add(double more_energy, double more_electrons)
      {
        energy += more_energy;
        electrons += more_electrons;
      }
    }; // Deposit_t

    // Map of electrons to store - catalogued by map[channel][tdc]
    std::map<raw::ChannelID_t, std::map<unsigned int, Deposit_t>> DepositsToStore;

    double xyz1[3] = {0.};

    double const xyz[3] = {
      stepMidPoint.x() / CLHEP::cm, stepMidPoint.y() / CLHEP::cm, stepMidPoint.z() / CLHEP::cm};

    // Already know which TPC we're in because we have been told

    try {
      const geo::TPCGeo& tpcg = fGeoHandle->TPC(tpc, cryostat);

      // X drift distance - the drift direction can be either in
      // the positive or negative direction, so use std::abs

      /// \todo think about effects of drift between planes
      double XDrift = std::abs(stepMidPoint.x() / CLHEP::cm - tpcg.PlaneLocation(0)[0]);
      //std::cout<<tpcg.DriftDirection()<<std::endl;
      if (tpcg.DriftDirection() == geo::kNegX)
        XDrift = stepMidPoint.x() / CLHEP::cm - tpcg.PlaneLocation(0)[0];
      else if (tpcg.DriftDirection() == geo::kPosX)
        XDrift = tpcg.PlaneLocation(0)[0] - stepMidPoint.x() / CLHEP::cm;

      if (XDrift < 0.) return;

      // Get SCE {x,y,z} offsets for particular location in TPC
      geo::Vector_t posOffsets;
      auto const* SCE = lar::providerFrom<spacecharge::SpaceChargeService>();
      if (SCE->EnableSimSpatialSCE() == true) {
        posOffsets = SCE->GetPosOffsets({xyz[0], xyz[1], xyz[2]});
        if (larsim::Utils::SCE::out_of_bounds(posOffsets)) {
          return;
        }
      }
      posOffsets.SetX(-posOffsets.X());

      // Drift time (nano-sec)
      double TDrift;
      XDrift += posOffsets.X();

      // Space charge distortion could push the energy deposit beyond the wire
      // plane (see issue #15131). Given that we don't have any subtlety in the
      // simulation of this region, bringing the deposit exactly on the plane
      // should be enough for the time being.
      if (XDrift < 0.) XDrift = 0.;

      TDrift = XDrift * RecipDriftVel[0];
      if (tpcg.Nplanes() == 2) { // special case for ArgoNeuT (plane 0 is the second wire plane)
        TDrift = ((XDrift - tpcg.PlanePitch(0, 1)) * RecipDriftVel[0] +
                  tpcg.PlanePitch(0, 1) * RecipDriftVel[1]);
      }

      const double lifetimecorrection = TMath::Exp(TDrift / LifetimeCorr_const);
      const int nIonizedElectrons =
        larg4::IonizationAndScintillation::Instance()->NumberIonizationElectrons();
      const double energy = larg4::IonizationAndScintillation::Instance()->EnergyDeposit();

      // if we have no electrons (too small energy or too large recombination)
      // we are done already here
      if (nIonizedElectrons <= 0) {
        MF_LOG_DEBUG("LArVoxelReadout")
          << "No electrons drifted to readout, " << energy << " MeV lost.";
        return;
      }
      // includes the effect of lifetime
      const double nElectrons = nIonizedElectrons * lifetimecorrection;

      // Longitudinal & transverse diffusion sigma (cm)
      double SqrtT = std::sqrt(TDrift);
      double LDiffSig = SqrtT * LDiff_const;
      double TDiffSig = SqrtT * TDiff_const;
      double electronclsize = fElectronClusterSize;

      int nClus = (int)std::ceil(nElectrons / electronclsize);
      if (nClus < fMinNumberOfElCluster) {
        electronclsize = nElectrons / fMinNumberOfElCluster;
        if (electronclsize < 1.0) { electronclsize = 1.0; }
        nClus = (int)std::ceil(nElectrons / electronclsize);
      }

      // Compute arrays of values as quickly as possible.
      std::vector<double> XDiff(nClus);
      std::vector<double> YDiff(nClus);
      std::vector<double> ZDiff(nClus);
      std::vector<double> nElDiff(nClus, electronclsize);
      std::vector<double> nEnDiff(nClus);

      // fix the number of electrons in the last cluster, that has smaller size
      nElDiff.back() = nElectrons - (nClus - 1) * electronclsize;

      for (size_t xx = 0; xx < nElDiff.size(); ++xx) {
        if (nElectrons > 0)
          nEnDiff[xx] = energy / nElectrons * nElDiff[xx];
        else
          nEnDiff[xx] = 0.;
      }

      double const avegageYtransversePos = (stepMidPoint.y() / CLHEP::cm) + posOffsets.Y();
      double const avegageZtransversePos = (stepMidPoint.z() / CLHEP::cm) + posOffsets.Z();

      // Smear drift times by x position and drift time
      if (LDiffSig > 0.0)
        PropRand.fireArray(nClus, &XDiff[0], 0., LDiffSig);
      else
        XDiff.assign(nClus, 0.0);

      if (TDiffSig > 0.0) {
        // Smear the Y,Z position by the transverse diffusion
        PropRand.fireArray(nClus, &YDiff[0], avegageYtransversePos, TDiffSig);
        PropRand.fireArray(nClus, &ZDiff[0], avegageZtransversePos, TDiffSig);
      }
      else {
        YDiff.assign(nClus, avegageYtransversePos);
        ZDiff.assign(nClus, avegageZtransversePos);
      }

      // make a collection of electrons for each plane
      for (size_t p = 0; p < tpcg.Nplanes(); ++p) {

        geo::PlaneGeo const& plane = tpcg.Plane(p);

        double Plane0Pitch = tpcg.Plane0Pitch(p);

        // "-" sign is because Plane0Pitch output is positive. Andrzej
        xyz1[0] = tpcg.PlaneLocation(0)[0] - Plane0Pitch;

        // Drift nClus electron clusters to the induction plane
        for (int k = 0; k < nClus; ++k) {
          // Correct drift time for longitudinal diffusion and plane
          double TDiff = TDrift + XDiff[k] * RecipDriftVel[0];
          // Take into account different Efields between planes
          // Also take into account special case for ArgoNeuT where Nplanes = 2.
          for (size_t ip = 0; ip < p; ++ip) {
            TDiff +=
              tpcg.PlanePitch(ip, ip + 1) * RecipDriftVel[tpcg.Nplanes() == 3 ? ip + 1 : ip + 2];
          }
          xyz1[1] = YDiff[k];
          xyz1[2] = ZDiff[k];

          /// \todo think about effects of drift between planes

          // grab the nearest channel to the xyz1 position
          try {
            if (fOffPlaneMargin != 0) {
              // get the effective position where to consider the charge landed;
              //
              // Some optimisations are possible; in particular, this method
              // could be extended to inform us if the point was too far.
              // Currently, if that is the case the code will proceed, find the
              // point is off plane, emit a warning and skip the deposition.
              //
              auto const landingPos = RecoverOffPlaneDeposit({xyz1[0], xyz1[1], xyz1[2]}, plane);
              xyz1[0] = landingPos.X();
              xyz1[1] = landingPos.Y();
              xyz1[2] = landingPos.Z();

            } // if charge lands off plane
            uint32_t channel = fGeoHandle->NearestChannel(xyz1, p, tpc, cryostat);

            /// \todo check on what happens if we allow the tdc value to be
            /// \todo beyond the end of the expected number of ticks
            // Add potential decay/capture/etc delay effect, simTime.
            unsigned int tdc = tpcClock.Ticks(clockData.G4ToElecTime(TDiff + simTime));

            // Add electrons produced by each cluster to the map
            DepositsToStore[channel][tdc].add(nEnDiff[k], nElDiff[k]);
          }
          catch (cet::exception& e) {
            MF_LOG_DEBUG("LArVoxelReadout")
              << "unable to drift electrons from point (" << xyz[0] << "," << xyz[1] << ","
              << xyz[2] << ") with exception " << e;
          }
        } // end loop over clusters
      }   // end loop over planes

      // Now store them in SimChannels
      ChannelMap_t& ChannelDataMap = fChannelMaps[cryostat][tpc];

      // browse deposited data on each channel: (channel; deposit data in time)
      for (auto const& deposit_per_channel : DepositsToStore) {

        raw::ChannelID_t channel = deposit_per_channel.first;

        // find whether we already have this channel
        auto iChannelData = ChannelDataMap.find(channel);

        // channelData is the SimChannel these deposits are going to be added to
        // If there is such a channel already, use it (first beanch).
        // If it's a new channel, the inner assignment creates a new SimChannel
        // in the map, and we save its reference in channelData
        sim::SimChannel& channelData = (iChannelData == ChannelDataMap.end()) ?
                                         (ChannelDataMap[channel] = sim::SimChannel(channel)) :
                                         iChannelData->second;

        // go through all deposits, one for each TDC: (TDC, deposit data)
        for (auto const& deposit_per_tdc : deposit_per_channel.second) {
          channelData.AddIonizationElectrons(trackID,
                                             deposit_per_tdc.first,
                                             deposit_per_tdc.second.electrons,
                                             xyz,
                                             deposit_per_tdc.second.energy,
                                             origTrackID);

        } // for deposit on TDCs
      }   // for deposit on channels

    } // end try intended to catch points where TPC can't be found
    catch (cet::exception& e) {
      MF_LOG_DEBUG("LArVoxelReadout") << "step cannot be found in a TPC\n" << e;
    }

    return;
  }

void larg4::LArVoxelReadout::EndOfEvent ( G4HCofThisEvent *  )

virtual

Definition at line 139 of file LArVoxelReadout.cxx.

  {
    MF_LOG_DEBUG("LArVoxelReadout") << "Total number of steps was " << fNSteps << std::endl;
  }

const LArVoxelReadout::ChannelMap_t & larg4::LArVoxelReadout::GetSimChannelMap

(

)

const

Returns the accumulated channel -> SimChannel map for the single TPC.

Definition at line 163 of file LArVoxelReadout.cxx.

  {
    if (bSingleTPC) return GetSimChannelMap(fCstat, fTPC);
    throw cet::exception("LArVoxelReadout") << "TPC not specified";
  }

LArVoxelReadout::ChannelMap_t & larg4::LArVoxelReadout::GetSimChannelMap

(

)

Definition at line 170 of file LArVoxelReadout.cxx.

  {
    if (bSingleTPC) return GetSimChannelMap(fCstat, fTPC);
    throw cet::exception("LArVoxelReadout") << "TPC not specified";
  }

const LArVoxelReadout::ChannelMap_t & larg4::LArVoxelReadout::GetSimChannelMap	(	unsigned short	cryo,
		unsigned short	tpc
	)		const

Returns the accumulated channel -> SimChannel map for the specified TPC.

Definition at line 177 of file LArVoxelReadout.cxx.

  {
    return fChannelMaps.at(cryo).at(tpc);
  }

LArVoxelReadout::ChannelMap_t & larg4::LArVoxelReadout::GetSimChannelMap	(	unsigned short	cryo,
		unsigned short	tpc
	)

Definition at line 183 of file LArVoxelReadout.cxx.

  {
    return fChannelMaps.at(cryo).at(tpc);
  }

std::vector< sim::SimChannel > larg4::LArVoxelReadout::GetSimChannels

(

)

const

Creates a list with the accumulated information for the single TPC.

Definition at line 189 of file LArVoxelReadout.cxx.

  {
    if (bSingleTPC) return GetSimChannels(fCstat, fTPC);
    throw cet::exception("LArVoxelReadout") << "TPC not specified";
  }

std::vector< sim::SimChannel > larg4::LArVoxelReadout::GetSimChannels	(	unsigned short	cryo,
		unsigned short	tpc
	)		const

Creates a list with the accumulated information for specified TPC.

Definition at line 196 of file LArVoxelReadout.cxx.

  {
    std::vector<sim::SimChannel> channels;
    const ChannelMap_t& chmap = fChannelMaps.at(cryo).at(tpc);
    channels.reserve(chmap.size());
    for (const auto& chpair : chmap)
      channels.push_back(chpair.second);
    return channels;
  }

bool larg4::LArVoxelReadout::Has ( std::vector< unsigned short int >  v, unsigned short int  tpc  ) const

inlineprivate

Definition at line 329 of file LArVoxelReadout.h.

    {
      for (auto c : v)
        if (c == tpc) return true;
      return false;
    }

void larg4::LArVoxelReadout::Initialize ( G4HCofThisEvent *  )

virtual

Definition at line 103 of file LArVoxelReadout.cxx.

  {
    assert(fClockData != nullptr &&
           "You must use set the clock data pointer at the beginning "
           "of each module's event-level call.  This might be done through the "
           "LArVoxelReadoutGeometry::Sentry class.");
    assert(fDetProp != nullptr &&
           "You must use set the detector-properties data pointer at the beginning "
           "of each module's event-level call.  This might be done through the "
           "LArVoxelReadoutGeometry::Sentry class.");

    fElectronLifetime = fDetProp->ElectronLifetime();
    for (int i = 0; i < 3; ++i)
      fDriftVelocity[i] =
        fDetProp->DriftVelocity(fDetProp->Efield(i), fDetProp->Temperature()) / 1000.;

    fElectronClusterSize = fLgpHandle->ElectronClusterSize();
    fMinNumberOfElCluster = fLgpHandle->MinNumberOfElCluster();
    fLongitudinalDiffusion = fLgpHandle->LongitudinalDiffusion();
    fTransverseDiffusion = fLgpHandle->TransverseDiffusion();
    fDontDriftThem = fLgpHandle->DisableWireplanes();
    fSkipWireSignalInTPCs = fLgpHandle->SkipWireSignalInTPCs();

    MF_LOG_DEBUG("LArVoxelReadout")
      << " e lifetime: " << fElectronLifetime << "\n Temperature: " << fDetProp->Temperature()
      << "\n Drift velocity: " << fDriftVelocity[0] << " " << fDriftVelocity[1] << " "
      << fDriftVelocity[2];

    fDontDriftThem = (fDontDriftThem || fLgpHandle->NoElectronPropagation());

    fNSteps = 0;
  }

void larg4::LArVoxelReadout::PrintAll ( )

virtual

Definition at line 598 of file LArVoxelReadout.cxx.

  {}

G4bool larg4::LArVoxelReadout::ProcessHits ( G4Step *  step, G4TouchableHistory *  pHistory  )

virtual

Definition at line 209 of file LArVoxelReadout.cxx.

  {
    // All work done for the "parallel world" "box of voxels" in
    // LArVoxelReadoutGeometry makes this a fairly simple routine.
    // First, the usual check for non-zero energy:

    // Only process the hit if the step is inside the active volume and
    // it has deposited energy.  The hit being inside the active volume
    // is virtually sure to happen because the LArVoxelReadoutGeometry
    // that this class makes use of only has voxels for inside the TPC.

    // The step can be no bigger than the size of the voxel,
    // because of the geometry set up in LArVoxelGeometry and the
    // transportation set up in PhysicsList.  Find the mid-point
    // of the step.

    if (step->GetTotalEnergyDeposit() > 0) {

      // Make sure we have the IonizationAndScintillation singleton
      // reset to this step
      larg4::IonizationAndScintillation::Instance()->Reset(step);
      fNSteps++;
      if (!fDontDriftThem) {

        G4ThreeVector midPoint =
          0.5 * (step->GetPreStepPoint()->GetPosition() + step->GetPostStepPoint()->GetPosition());
        double g4time = step->GetPreStepPoint()->GetGlobalTime();

        // Find the Geant4 track ID for the particle responsible for depositing the
        // energy.  if we are only storing primary EM shower particles, and this energy
        // is from a secondary etc EM shower particle, the ID returned is the primary
        const int trackID = ParticleListAction::GetCurrentTrackID();
        // For all particles but shower daughters, origTrackID is the same as trackID.
        // For shower daughters it contains their original trackID instead of the
        // shower primary's trackID.
        const int origTrackID = ParticleListAction::GetCurrentOrigTrackID();

        // Find out which TPC we are in.
        // If this readout object covers just one, we already know it.
        // Otherwise, we have to ask Geant where we are.
        unsigned short int cryostat = 0, tpc = 0;
        if (bSingleTPC) {
          cryostat = fCstat;
          tpc = fTPC;
        }
        else {
          // detect the TPC we are in
          const G4VTouchable* pTouchable = step->GetPreStepPoint()->GetTouchable();
          if (!pTouchable) {
            throw cet::exception("LArG4") << "Untouchable step in LArVoxelReadout::ProcessHits()";
          }

          // one of the ancestors of the touched volume is supposed to be
          // actually a G4PVPlacementInTPC that knows which TPC it covers;
          // currently, it's the 4th in the ladder:
          // [0] voxel [1] voxel tower [2] voxel plane [3] the full box;
          G4int depth = 0;
          while (depth < pTouchable->GetHistoryDepth()) {
            const G4PVPlacementInTPC* pPVinTPC =
              dynamic_cast<const G4PVPlacementInTPC*>(pTouchable->GetVolume(depth++));
            if (!pPVinTPC) continue;
            cryostat = pPVinTPC->ID.Cryostat;
            tpc = pPVinTPC->ID.TPC;
            if (Has(fSkipWireSignalInTPCs, tpc)) { return true; }
            break;
          } // while
          if (depth < pTouchable->GetHistoryDepth()) {
            // this is a fundamental error where the step does not happen in
            // any TPC; this should not happen in the readout geometry!
            throw cet::exception("LArG4") << "No TPC ID found in LArVoxelReadout::ProcessHits()";
          } // if
          MF_LOG_DEBUG("LArVoxelReadoutHit") << " hit in C=" << cryostat << " T=" << tpc;
        } // if more than one TPC

        // Note that if there is no particle ID for this energy deposit, the
        // trackID will be sim::NoParticleId.

        DriftIonizationElectrons(*fClockData, midPoint, g4time, trackID, cryostat, tpc, origTrackID);
      } // end we are drifting
    }   // end there is non-zero energy deposition

    return true;
  }

void larg4::LArVoxelReadout::ProcessStep ( G4Step *  )

private

geo::Point_t larg4::LArVoxelReadout::RecoverOffPlaneDeposit ( geo::Point_t const &  pos, geo::PlaneGeo const &  plane  ) const

private

Returns the point on the specified plane closest to position.

Parameters

pos	the position to be tested (global coordinates, centimeters)
plane	the plane to test the position against

Returns

a position on plane, unless pos is too far from it

This method considers the distance of the position pos from the active part of the plane (see geo::Plane::DeltaFromActivePlane()). If the position is less than a configurable margin far from the plane, the closest point on the plane to that position is returned. Otherwise, the position itself is returned.

Ionization charge may be drifted so that when it arrives to the plane, it actually does not hit the area covered by wires. This can happen for many reasons:

• space charge distortion led the point outside the fiducial volume (this may be prevented by specific code)
• diffusion pushes the charge outside the instrumented region
• the geometry of the wire planes is such that planes have different coverage and what one plane can cover, the next can't

The "recovery" consists in forcing the charge to the instrumented area covered by the plane wires. The distance of the drifted charge from each plane border is computed and compared to the margin. If that distance is smaller than the margin, it is neglected and the charge is assigned a new position on that border.

This method provides the position that should be used for the charge deposition.

This is a simplistic approach to the simulation of border effects, assuming that in fact the electric field, which is continuous and pointing to the collection wires, will drive the charge to the wires even when they are "off track". No correction is applied for the additional time that such deviation would take.

Definition at line 303 of file LArVoxelReadout.cxx.

  {
    //
    // translate the landing position on the two frame coordinates
    // ("width" and "depth")
    //
    auto const landingPos = plane.PointWidthDepthProjection(pos);

    //
    // compute the distance of the landing position on the two frame
    // coordinates ("width" and "depth");
    // keep the point within 10 micrometers (0.001 cm) from the border
    //
    auto const offPlane = plane.DeltaFromActivePlane(landingPos, 0.001);

    //
    // if both the distances are below the margin, move the point to
    // the border
    //

    // nothing to recover: landing is inside
    if ((offPlane.X() == 0.0) && (offPlane.Y() == 0.0)) return pos;

    // won't recover: too far
    if ((std::abs(offPlane.X()) > fOffPlaneMargin) || (std::abs(offPlane.Y()) > fOffPlaneMargin))
      return pos;

    // we didn't fully decompose because it might be unnecessary;
    // now we need the full thing
    auto const distance = plane.DistanceFromPlane(pos);

    return plane.ComposePoint<geo::Point_t>(distance, landingPos + offPlane);
  } // LArVoxelReadout::RecoverOffPlaneDeposit()

void larg4::LArVoxelReadout::SetClockData ( detinfo::DetectorClocksData const *const  clockData )

inlineprivatenoexcept

Definition at line 251 of file LArVoxelReadout.h.

    {
      fClockData = clockData;
    }

void larg4::LArVoxelReadout::SetDiscoverTPC

(

)

Sets this readout to discover the TPC of each processed hit.

Definition at line 92 of file LArVoxelReadout.cxx.

  {
    bSingleTPC = false;
    fCstat = 0;
    fTPC = 0;
    MF_LOG_DEBUG("LArVoxelReadout") << GetName() << " autodetects TPC";
  }

void larg4::LArVoxelReadout::SetOffPlaneChargeRecoveryMargin ( double  margin )

inlineprivate

Sets the margin for recovery of charge drifted off-plane.

Parameters

margin

the extent of the margin on each frame coordinate [cm]

This method sets the margin for the recovery of off-plane ionization charge. See RecoverOffPlaneDeposit() for a description of that feature.

This method is used by LArVoxelReadout::Setup().

Definition at line 272 of file LArVoxelReadout.h.

    {
      fOffPlaneMargin = std::max(margin, 0.0);
    }

void larg4::LArVoxelReadout::SetPropertiesData ( detinfo::DetectorPropertiesData const *const  detProp )

inlineprivatenoexcept

Definition at line 257 of file LArVoxelReadout.h.

    {
      fDetProp = detProp;
    }

void larg4::LArVoxelReadout::SetRandomEngines ( CLHEP::HepRandomEngine *  pPropGen )

private

Sets the random generators to be used.

Definition at line 295 of file LArVoxelReadout.cxx.

  {
    assert(pPropGen); // random engine must be present
    fPropGen = pPropGen;
  }

void larg4::LArVoxelReadout::SetSingleTPC	(	unsigned int	cryostat,
		unsigned int	tpc
	)

Associates this readout to one specific TPC.

Definition at line 83 of file LArVoxelReadout.cxx.

  {
    bSingleTPC = true;
    fCstat = cryostat;
    fTPC = tpc;
    MF_LOG_DEBUG("LArVoxelReadout") << GetName() << "covers C=" << fCstat << " T=" << fTPC;
  }

void larg4::LArVoxelReadout::Setup

(

Setup_t const &

setupData

)

Reads all the configuration elements from setupData

Definition at line 75 of file LArVoxelReadout.cxx.

  {
    SetOffPlaneChargeRecoveryMargin(setupData.offPlaneMargin);
    SetRandomEngines(setupData.propGen);
  }

Friends And Related Function Documentation

friend class LArVoxelReadoutGeometry

friend

Definition at line 248 of file LArVoxelReadout.h.

Member Data Documentation

bool larg4::LArVoxelReadout::bSingleTPC

private

true if this readout is associated with a single TPC

Definition at line 356 of file LArVoxelReadout.h.

std::vector<std::vector<ChannelMap_t> > larg4::LArVoxelReadout::fChannelMaps

private

Maps of cryostat, tpc to channel data.

Definition at line 350 of file LArVoxelReadout.h.

detinfo::DetectorClocksData const* larg4::LArVoxelReadout::fClockData {nullptr}

private

Definition at line 360 of file LArVoxelReadout.h.

unsigned int larg4::LArVoxelReadout::fCstat

private

and in which cryostat (if bSingleTPC is true)

Definition at line 355 of file LArVoxelReadout.h.

detinfo::DetectorPropertiesData const* larg4::LArVoxelReadout::fDetProp {nullptr}

private

Definition at line 361 of file LArVoxelReadout.h.

bool larg4::LArVoxelReadout::fDontDriftThem

private

Definition at line 345 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fDriftVelocity[3]

private

Definition at line 339 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fElectronClusterSize

private

Definition at line 343 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fElectronLifetime

private

Definition at line 342 of file LArVoxelReadout.h.

art::ServiceHandle<geo::Geometry const> larg4::LArVoxelReadout::fGeoHandle

private

Handle to the Geometry service.

Definition at line 351 of file LArVoxelReadout.h.

art::ServiceHandle<sim::LArG4Parameters const> larg4::LArVoxelReadout::fLgpHandle

private

Handle to the LArG4 parameters service.

Definition at line 353 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fLongitudinalDiffusion

private

Definition at line 340 of file LArVoxelReadout.h.

int larg4::LArVoxelReadout::fMinNumberOfElCluster

private

Definition at line 344 of file LArVoxelReadout.h.

size_t larg4::LArVoxelReadout::fNSteps

private

Definition at line 368 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fOffPlaneMargin = 0.0

private

Charge deposited within this many [cm] from the plane is lead onto it.

Definition at line 348 of file LArVoxelReadout.h.

CLHEP::HepRandomEngine* larg4::LArVoxelReadout::fPropGen = nullptr

private

random engine for charge propagation

Definition at line 358 of file LArVoxelReadout.h.

std::vector<unsigned short int> larg4::LArVoxelReadout::fSkipWireSignalInTPCs

private

Definition at line 346 of file LArVoxelReadout.h.

G4ThreeVector larg4::LArVoxelReadout::fStepEnd

private

Definition at line 367 of file LArVoxelReadout.h.

G4ThreeVector larg4::LArVoxelReadout::fStepStart

private

Definition at line 366 of file LArVoxelReadout.h.

unsigned int larg4::LArVoxelReadout::fTPC

private

which TPC this LArVoxelReadout corresponds to

Definition at line 354 of file LArVoxelReadout.h.

double larg4::LArVoxelReadout::fTransverseDiffusion

private

Definition at line 341 of file LArVoxelReadout.h.

---

The documentation for this class was generated from the following files:

• LArVoxelReadout.h
• LArVoxelReadout.cxx