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, geo::GeometryCore const *geom, geo::WireReadoutGeom const *wireReadoutGeom, sim::LArG4Parameters const *lgp)
	Constructor. Can detect which TPC to cover by the name. More...
	LArVoxelReadout (std::string const &name, geo::GeometryCore const *geom, geo::WireReadoutGeom const *wireReadoutGeom, sim::LArG4Parameters const *lgp, 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 ()
	Returns the accumulated channel -> SimChannel map for the single TPC. More...
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)
	Returns the accumulated channel -> SimChannel map for the specified TPC. More...

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...
geo::GeometryCore const *	fGeo
geo::WireReadoutGeom const *	fWireReadoutGeom
sim::LArG4Parameters const *	fLgp
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 159 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 162 of file LArVoxelReadout.h.

Constructor & Destructor Documentation

larg4::LArVoxelReadout::LArVoxelReadout	(	std::string const &	name,
		geo::GeometryCore const *	geom,
		geo::WireReadoutGeom const *	wireReadoutGeom,
		sim::LArG4Parameters const *	lgp
	)

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

Definition at line 51 of file LArVoxelReadout.cxx.

References ClearSimChannels(), fLgp, fWireReadoutGeom, SetDiscoverTPC(), and SetSingleTPC().

    : G4VSensitiveDetector(name), fGeo{geom}, fWireReadoutGeom{wireReadoutGeom}, fLgp{lgp}
  {
    // 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();
  }

larg4::LArVoxelReadout::LArVoxelReadout	(	std::string const &	name,
		geo::GeometryCore const *	geom,
		geo::WireReadoutGeom const *	wireReadoutGeom,
		sim::LArG4Parameters const *	lgp,
		unsigned int	cryostat,
		unsigned int	tpc
	)

Constructor. Sets which TPC to work on.

Definition at line 72 of file LArVoxelReadout.cxx.

References SetSingleTPC().

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

Member Function Documentation

void larg4::LArVoxelReadout::clear ( )

virtual

Definition at line 150 of file LArVoxelReadout.cxx.

{}

void larg4::LArVoxelReadout::ClearSimChannels

(

)

Definition at line 153 of file LArVoxelReadout.cxx.

References fChannelMaps, fGeo, geo::GeometryCore::Ncryostats(), and geo::GeometryCore::NTPC().

Referenced by LArVoxelReadout().

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

void larg4::LArVoxelReadout::DrawAll ( )

virtual

Definition at line 571 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 324 of file LArVoxelReadout.cxx.

References util::abs(), sim::SimChannel::AddIonizationElectrons(), geo::TPCGeo::DriftAxisWithSign(), e, energy, larg4::IonizationAndScintillation::EnergyDeposit(), fChannelMaps, fDriftVelocity, fElectronClusterSize, fElectronLifetime, fGeo, geo::WireReadoutGeom::FirstPlane(), fLongitudinalDiffusion, fMinNumberOfElCluster, fOffPlaneMargin, for(), fPropGen, fTransverseDiffusion, fWireReadoutGeom, detinfo::DetectorClocksData::G4ToElecTime(), geo::PlaneGeo::GetCenter(), larg4::IonizationAndScintillation::Instance(), geo::Iterable< IterationPolicy, Transform >::Iterate(), MF_LOG_DEBUG, geo::WireReadoutGeom::NearestChannel(), geo::Negative, geo::WireReadoutGeom::Nplanes(), larg4::IonizationAndScintillation::NumberIonizationElectrons(), larsim::Utils::SCE::out_of_bounds(), geo::WireReadoutGeom::Plane0Pitch(), geo::WireReadoutGeom::PlanePitch(), geo::Positive, RecoverOffPlaneDeposit(), geo::GeometryCore::TPC(), detinfo::DetectorClocksData::TPCClock(), and xx.

Referenced by ProcessHits().

  {
    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;

    geo::Point_t xyz1;

    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 {
      geo::TPCID const tpcid{cryostat, tpc};
      geo::TPCGeo const& tpcg = fGeo->TPC(tpcid);

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

      auto const& plane_center = fWireReadoutGeom->FirstPlane(tpcid).GetCenter();
      double XDrift = std::abs(stepMidPoint.x() / CLHEP::cm - plane_center.X());
      auto const [_, driftSign] = tpcg.DriftAxisWithSign();
      if (driftSign == geo::DriftSign::Negative)
        XDrift = stepMidPoint.x() / CLHEP::cm - plane_center.X();
      else if (driftSign == geo::DriftSign::Positive)
        XDrift = plane_center.X() - 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];
      unsigned int const nplanes = fWireReadoutGeom->Nplanes(tpcid);
      if (nplanes == 2) { // special case for ArgoNeuT (plane 0 is the second wire plane)
        double const pitch = fWireReadoutGeom->PlanePitch(tpcid);
        TDrift = ((XDrift - pitch) * RecipDriftVel[0] + pitch * 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 averageYtransversePos = (stepMidPoint.y() / CLHEP::cm) + posOffsets.Y();
      double const averageZtransversePos = (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], averageYtransversePos, TDiffSig);
        PropRand.fireArray(nClus, &ZDiff[0], averageZtransversePos, TDiffSig);
      }
      else {
        YDiff.assign(nClus, averageYtransversePos);
        ZDiff.assign(nClus, averageZtransversePos);
      }

      // make a collection of electrons for each plane
      double const plane_0_center_x = fWireReadoutGeom->FirstPlane(tpcid).GetCenter().X();
      for (auto const& plane : fWireReadoutGeom->Iterate<geo::PlaneGeo>(tpcid)) {

        double Plane0Pitch = fWireReadoutGeom->Plane0Pitch(tpcid, plane.ID().Plane);

        // "-" sign is because Plane0Pitch output is positive. Andrzej
        xyz1.SetX(plane_0_center_x - 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 (unsigned int ip = 0; ip < plane.ID().Plane; ++ip) {
            TDiff += fWireReadoutGeom->PlanePitch(tpcid, ip, ip + 1) *
                     RecipDriftVel[nplanes == 3 ? ip + 1 : ip + 2];
          }
          xyz1.SetY(YDiff[k]);
          xyz1.SetZ(ZDiff[k]);


          // 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.
              xyz1 = RecoverOffPlaneDeposit(xyz1, plane);
            } // if charge lands off plane
            uint32_t channel = fWireReadoutGeom->NearestChannel(xyz1, plane.ID());

            // 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 144 of file LArVoxelReadout.cxx.

References fNSteps, and MF_LOG_DEBUG.

  {
    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 164 of file LArVoxelReadout.cxx.

References bSingleTPC, fCstat, and fTPC.

Referenced by GetSimChannelMap().

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

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

(

)

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

Definition at line 170 of file LArVoxelReadout.cxx.

References bSingleTPC, fCstat, fTPC, and GetSimChannelMap().

  {
    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 176 of file LArVoxelReadout.cxx.

References fChannelMaps.

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

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

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

Definition at line 182 of file LArVoxelReadout.cxx.

References fChannelMaps.

  {
    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 188 of file LArVoxelReadout.cxx.

References bSingleTPC, fCstat, and fTPC.

  {
    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 194 of file LArVoxelReadout.cxx.

References fChannelMaps.

  {
    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 319 of file LArVoxelReadout.h.

Referenced by ProcessHits().

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

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

virtual

Definition at line 109 of file LArVoxelReadout.cxx.

References sim::LArG4Parameters::DisableWireplanes(), detinfo::DetectorPropertiesData::DriftVelocity(), detinfo::DetectorPropertiesData::Efield(), sim::LArG4Parameters::ElectronClusterSize(), detinfo::DetectorPropertiesData::ElectronLifetime(), fClockData, fDetProp, fDontDriftThem, fDriftVelocity, fElectronClusterSize, fElectronLifetime, fLgp, fLongitudinalDiffusion, fMinNumberOfElCluster, fNSteps, fSkipWireSignalInTPCs, fTransverseDiffusion, sim::LArG4Parameters::LongitudinalDiffusion(), MF_LOG_DEBUG, sim::LArG4Parameters::MinNumberOfElCluster(), sim::LArG4Parameters::NoElectronPropagation(), sim::LArG4Parameters::SkipWireSignalInTPCs(), detinfo::DetectorPropertiesData::Temperature(), and sim::LArG4Parameters::TransverseDiffusion().

  {
    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 = fLgp->ElectronClusterSize();
    fMinNumberOfElCluster = fLgp->MinNumberOfElCluster();
    fLongitudinalDiffusion = fLgp->LongitudinalDiffusion();
    fTransverseDiffusion = fLgp->TransverseDiffusion();
    fDontDriftThem = fLgp->DisableWireplanes();
    fSkipWireSignalInTPCs = fLgp->SkipWireSignalInTPCs();

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

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

    fNSteps = 0;
  }

void larg4::LArVoxelReadout::PrintAll ( )

virtual

Definition at line 572 of file LArVoxelReadout.cxx.

{}

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

virtual

Definition at line 207 of file LArVoxelReadout.cxx.

References bSingleTPC, DriftIonizationElectrons(), fClockData, fCstat, fDontDriftThem, fNSteps, fSkipWireSignalInTPCs, fTPC, larg4::ParticleListAction::GetCurrentOrigTrackID(), larg4::ParticleListAction::GetCurrentTrackID(), Has(), larg4::G4PVPlacementWithID< IDTYPE >::ID, larg4::IonizationAndScintillation::Instance(), MF_LOG_DEBUG, and larg4::IonizationAndScintillation::Reset().

  {
    // 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 296 of file LArVoxelReadout.cxx.

References util::abs(), geo::PlaneGeo::ComposePoint(), geo::PlaneGeo::DeltaFromActivePlane(), geo::PlaneGeo::DistanceFromPlane(), fOffPlaneMargin, and geo::PlaneGeo::PointWidthDepthProjection().

Referenced by DriftIonizationElectrons().

  {
    // 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(distance, landingPos + offPlane);
  }

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

inlineprivatenoexcept

Definition at line 247 of file LArVoxelReadout.h.

    {
      fClockData = clockData;
    }

void larg4::LArVoxelReadout::SetDiscoverTPC

(

)

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

Definition at line 99 of file LArVoxelReadout.cxx.

References bSingleTPC, fCstat, fTPC, and MF_LOG_DEBUG.

Referenced by LArVoxelReadout().

  {
    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 266 of file LArVoxelReadout.h.

Referenced by Setup().

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

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

inlineprivatenoexcept

Definition at line 252 of file LArVoxelReadout.h.

    {
      fDetProp = detProp;
    }

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

private

Sets the random generators to be used.

Definition at line 289 of file LArVoxelReadout.cxx.

References fPropGen.

Referenced by Setup().

  {
    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 91 of file LArVoxelReadout.cxx.

References bSingleTPC, fCstat, fTPC, and MF_LOG_DEBUG.

Referenced by larg4::LArVoxelReadoutGeometry::Construct(), and LArVoxelReadout().

  {
    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 84 of file LArVoxelReadout.cxx.

References larg4::LArVoxelReadout::Setup_t::offPlaneMargin, larg4::LArVoxelReadout::Setup_t::propGen, SetOffPlaneChargeRecoveryMargin(), and SetRandomEngines().

Referenced by larg4::LArVoxelReadoutGeometry::Construct().

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

Friends And Related Function Documentation

friend class LArVoxelReadoutGeometry

friend

Definition at line 245 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 346 of file LArVoxelReadout.h.

Referenced by GetSimChannelMap(), GetSimChannels(), ProcessHits(), SetDiscoverTPC(), and SetSingleTPC().

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

private

Maps of cryostat, tpc to channel data.

Definition at line 340 of file LArVoxelReadout.h.

Referenced by ClearSimChannels(), DriftIonizationElectrons(), GetSimChannelMap(), and GetSimChannels().

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

private

Definition at line 350 of file LArVoxelReadout.h.

Referenced by Initialize(), and ProcessHits().

unsigned int larg4::LArVoxelReadout::fCstat

private

and in which cryostat (if bSingleTPC is true)

Definition at line 345 of file LArVoxelReadout.h.

Referenced by GetSimChannelMap(), GetSimChannels(), ProcessHits(), SetDiscoverTPC(), and SetSingleTPC().

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

private

Definition at line 351 of file LArVoxelReadout.h.

Referenced by Initialize().

bool larg4::LArVoxelReadout::fDontDriftThem

private

Definition at line 335 of file LArVoxelReadout.h.

Referenced by Initialize(), and ProcessHits().

double larg4::LArVoxelReadout::fDriftVelocity[3]

private

Definition at line 329 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

double larg4::LArVoxelReadout::fElectronClusterSize

private

Definition at line 333 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

double larg4::LArVoxelReadout::fElectronLifetime

private

Definition at line 332 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

geo::GeometryCore const* larg4::LArVoxelReadout::fGeo

private

Definition at line 341 of file LArVoxelReadout.h.

Referenced by ClearSimChannels(), and DriftIonizationElectrons().

sim::LArG4Parameters const* larg4::LArVoxelReadout::fLgp

private

Definition at line 343 of file LArVoxelReadout.h.

Referenced by Initialize(), and LArVoxelReadout().

double larg4::LArVoxelReadout::fLongitudinalDiffusion

private

Definition at line 330 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

int larg4::LArVoxelReadout::fMinNumberOfElCluster

private

Definition at line 334 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

size_t larg4::LArVoxelReadout::fNSteps

private

Definition at line 358 of file LArVoxelReadout.h.

Referenced by EndOfEvent(), Initialize(), and ProcessHits().

double larg4::LArVoxelReadout::fOffPlaneMargin = 0.0

private

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

Definition at line 338 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and RecoverOffPlaneDeposit().

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

private

random engine for charge propagation

Definition at line 348 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and SetRandomEngines().

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

private

Definition at line 336 of file LArVoxelReadout.h.

Referenced by Initialize(), and ProcessHits().

G4ThreeVector larg4::LArVoxelReadout::fStepEnd

private

Definition at line 357 of file LArVoxelReadout.h.

G4ThreeVector larg4::LArVoxelReadout::fStepStart

private

Definition at line 356 of file LArVoxelReadout.h.

unsigned int larg4::LArVoxelReadout::fTPC

private

which TPC this LArVoxelReadout corresponds to

Definition at line 344 of file LArVoxelReadout.h.

Referenced by GetSimChannelMap(), GetSimChannels(), ProcessHits(), SetDiscoverTPC(), and SetSingleTPC().

double larg4::LArVoxelReadout::fTransverseDiffusion

private

Definition at line 331 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and Initialize().

geo::WireReadoutGeom const* larg4::LArVoxelReadout::fWireReadoutGeom

private

Definition at line 342 of file LArVoxelReadout.h.

Referenced by DriftIonizationElectrons(), and LArVoxelReadout().

---

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

• larsim/v10_02_02/source/larsim/LegacyLArG4/LArVoxelReadout.h
• larsim/v10_02_02/source/larsim/LegacyLArG4/LArVoxelReadout.cxx