LArSoft  v06_85_00
Liquid Argon Software toolkit - http://larsoft.org/
LArVoxelReadout.cxx
Go to the documentation of this file.
1 
8 // C/C++ standard library
9 #include <cstdio> // std::sscanf()
10 #include <cmath> // std::ceil()
11 #include <string>
12 #include <map>
13 #include <utility> // std::move()
14 #include <cassert>
15 
16 // GEANT
17 #include "Geant4/G4HCofThisEvent.hh"
18 #include "Geant4/G4TouchableHistory.hh"
19 #include "Geant4/G4TouchableHandle.hh"
20 #include "Geant4/G4Step.hh"
21 #include "Geant4/G4StepPoint.hh"
22 #include "Geant4/G4ThreeVector.hh"
23 #include "Geant4/G4VPhysicalVolume.hh"
24 
25 // framework libraries
26 #include "cetlib_except/exception.h"
29 
30 // LArSoft code
34 #include "larcoreobj/SimpleTypesAndConstants/RawTypes.h" // raw::ChannelID_t
35 
36 // CLHEP
37 #include "CLHEP/Random/RandGauss.h"
38 #include "CLHEP/Random/RandPoisson.h"
39 #include "CLHEP/Random/RandFlat.h"
40 
41 namespace larg4 {
42 
43 
44  //---------------------------------------------------------------------------------------
45  // Constructor. Note that we force the name of this sensitive
46  // detector to be the value expected by LArVoxelListAction.
47  LArVoxelReadout::LArVoxelReadout(std::string const& name)
48  : G4VSensitiveDetector(name)
49  {
50  // Initialize values for the electron-cluster calculation.
52 
53  const detinfo::DetectorClocks* ts = lar::providerFrom<detinfo::DetectorClocksService>();
54  fClock = ts->TPCClock();
55 
56  // the standard name contains cryostat and TPC;
57  // if we don't find it, we will detect the TPC at each Geant hit
58  unsigned int cryostat, tpc;
59  if (std::sscanf(name.c_str(),"%*19s%1u%*4s%u",&cryostat, &tpc) == 2)
60  SetSingleTPC(cryostat, tpc);
61  else
63 
64  } // LArVoxelReadout::LArVoxelReadout()
65 
66  //---------------------------------------------------------------------------------------
67  // Constructor. Note that we force the name of this sensitive
68  // detector to be the value expected by LArVoxelListAction.
70  (std::string const& name, unsigned int cryostat, unsigned int tpc)
71  : LArVoxelReadout(name)
72  { SetSingleTPC(cryostat, tpc); }
73 
74 
75  //---------------------------------------------------------------------------------------
76  void LArVoxelReadout::Setup(Setup_t const& setupData) {
77 
79  SetRandomEngines(setupData.propGen);
80 
81  } // LArVoxelReadout::Setup()
82 
83 
84  //---------------------------------------------------------------------------------------
85  void LArVoxelReadout::SetSingleTPC(unsigned int cryostat, unsigned int tpc) {
86  bSingleTPC = true;
87  fCstat = cryostat;
88  fTPC = tpc;
89  LOG_DEBUG("LArVoxelReadout")
90  << GetName() << "covers C=" << fCstat << " T=" << fTPC;
91  } // LArVoxelReadout::SetSingleTPC()
92 
94  bSingleTPC = false;
95  fCstat = 0;
96  fTPC = 0;
97  LOG_DEBUG("LArVoxelReadout") << GetName() << " autodetects TPC";
98  } // LArVoxelReadout::SetDiscoverTPC()
99 
100 
101  //---------------------------------------------------------------------------------------
103 
104  //---------------------------------------------------------------------------------------
105  // Called at the start of each event.
107  {
108  // for c2: larp is unused
109  //auto const * larp = lar::providerFrom<detinfo::LArPropertiesService>();
110  auto const * detprop = lar::providerFrom<detinfo::DetectorPropertiesService>();
111  fElectronLifetime = detprop->ElectronLifetime();
112  for (int i = 0; i<3; ++i)
113  fDriftVelocity[i] = detprop->DriftVelocity(detprop->Efield(i),
114  detprop->Temperature())/1000.;
115 
122 
123  LOG_DEBUG("LArVoxelReadout") << " e lifetime: " << fElectronLifetime
124  << "\n Temperature: " << detprop->Temperature()
125  << "\n Drift velocity: " << fDriftVelocity[0]
126  <<" "<<fDriftVelocity[1]<<" "<<fDriftVelocity[2];
127 
129 
130  fNSteps=0;
131  }
132 
133  //---------------------------------------------------------------------------------------
134  // Called at the end of each event.
136  {
137  mf::LogWarning("LArVoxelReadout") << "Total number of steps was " << fNSteps << std::endl;
138 
139  } // LArVoxelReadout::EndOfEvent()
140 
141  //---------------------------------------------------------------------------------------
143  {
144  }
145 
146  //--------------------------------------------------------------------------------------
149  size_t cryo = 0;
150  for (auto& cryoData: fChannelMaps) { // each, a vector of maps
151  cryoData.resize(fGeoHandle->NTPC(cryo++));
152  for (auto& channelsMap: cryoData) channelsMap.clear(); // each, a map
153  } // for cryostats
154  } // LArVoxelReadout::ClearSimChannels()
155 
156 
158  {
159  if (bSingleTPC) return GetSimChannelMap(fCstat, fTPC);
160  throw cet::exception("LArVoxelReadout") << "TPC not specified";
161  } // LArVoxelReadout::GetSimChannelMap() const
162 
164  if (bSingleTPC) return GetSimChannelMap(fCstat, fTPC);
165  throw cet::exception("LArVoxelReadout") << "TPC not specified";
166  } // LArVoxelReadout::GetSimChannelMap()
167 
168 
170  (unsigned short cryo, unsigned short tpc) const
171  { return fChannelMaps.at(cryo).at(tpc); }
172 
174  (unsigned short cryo, unsigned short tpc)
175  { return fChannelMaps.at(cryo).at(tpc); }
176 
177 
178  std::vector<sim::SimChannel> LArVoxelReadout::GetSimChannels() const {
179  if (bSingleTPC) return GetSimChannels(fCstat, fTPC);
180  throw cet::exception("LArVoxelReadout") << "TPC not specified";
181  } // LArVoxelReadout::GetSimChannels()
182 
183  std::vector<sim::SimChannel> LArVoxelReadout::GetSimChannels
184  (unsigned short cryo, unsigned short tpc) const
185  {
186  std::vector<sim::SimChannel> channels;
187  const ChannelMap_t& chmap = fChannelMaps.at(cryo).at(tpc);
188  channels.reserve(chmap.size());
189  for(const auto& chpair: chmap) channels.push_back(chpair.second);
190  return channels;
191  } // LArVoxelReadout::GetSimChannels(short, short)
192 
193 
194  //---------------------------------------------------------------------------------------
195  // Called for each step.
196  G4bool LArVoxelReadout::ProcessHits( G4Step* step, G4TouchableHistory* pHistory)
197  {
198  // All work done for the "parallel world" "box of voxels" in
199  // LArVoxelReadoutGeometry makes this a fairly simple routine.
200  // First, the usual check for non-zero energy:
201 
202  // Only process the hit if the step is inside the active volume and
203  // it has deposited energy. The hit being inside the active volume
204  // is virtually sure to happen because the LArVoxelReadoutGeometry
205  // that this class makes use of only has voxels for inside the TPC.
206 
207  // The step can be no bigger than the size of the voxel,
208  // because of the geometry set up in LArVoxelGeometry and the
209  // transportation set up in PhysicsList. Find the mid-point
210  // of the step.
211 
212  if ( step->GetTotalEnergyDeposit() > 0 ){
213 
214  // Make sure we have the IonizationAndScintillation singleton
215  // reset to this step
217  fNSteps++;
218  if( !fDontDriftThem ){
219 
220  G4ThreeVector midPoint = 0.5*( step->GetPreStepPoint()->GetPosition()
221  + step->GetPostStepPoint()->GetPosition() );
222  double g4time = step->GetPreStepPoint()->GetGlobalTime();
223 
224  // Find the Geant4 track ID for the particle responsible for depositing the
225  // energy. if we are only storing primary EM shower particles, and this energy
226  // is from a secondary etc EM shower particle, the ID returned is the primary
227  const int trackID = ParticleListAction::GetCurrentTrackID();
228 
229  // Find out which TPC we are in.
230  // If this readout object covers just one, we already know it.
231  // Otherwise, we have to ask Geant where we are.
232  unsigned short int cryostat = 0, tpc = 0;
233  if (bSingleTPC) {
234  cryostat = fCstat;
235  tpc = fTPC;
236  }
237  else {
238  // detect the TPC we are in
239  const G4VTouchable* pTouchable = step->GetPreStepPoint()->GetTouchable();
240  if (!pTouchable) {
241  throw cet::exception
242  ("LArG4") << "Untouchable step in LArVoxelReadout::ProcessHits()";
243  }
244 
245  // one of the ancestors of the touched volume is supposed to be
246  // actually a G4PVPlacementInTPC that knows which TPC it covers;
247  // currently, it's the 4th in the ladder:
248  // [0] voxel [1] voxel tower [2] voxel plane [3] the full box;
249  G4int depth = 0;
250  while (depth < pTouchable->GetHistoryDepth()) {
251  const G4PVPlacementInTPC* pPVinTPC =
252  dynamic_cast<const G4PVPlacementInTPC*>
253  (pTouchable->GetVolume(depth++));
254  if (!pPVinTPC) continue;
255  cryostat = pPVinTPC->ID.Cryostat;
256  tpc = pPVinTPC->ID.TPC;
257  if (Has(fSkipWireSignalInTPCs, tpc))
258  {
259  return true;
260  }
261  break;
262  } // while
263  if (depth < pTouchable->GetHistoryDepth()) {
264  // this is a fundamental error where the step does not happen in
265  // any TPC; this should not happen in the readout geometry!
266  throw cet::exception
267  ("LArG4") << "No TPC ID found in LArVoxelReadout::ProcessHits()";
268  } // if
269  LOG_DEBUG("LArVoxelReadoutHit") << " hit in C=" << cryostat << " T=" << tpc;
270  } // if more than one TPC
271 
272  // Note that if there is no particle ID for this energy deposit, the
273  // trackID will be sim::NoParticleId.
274 
275  DriftIonizationElectrons(midPoint, g4time, trackID, cryostat, tpc);
276  } // end we are drifting
277  } // end there is non-zero energy deposition
278 
279  return true;
280  }
281 
282 
283  //----------------------------------------------------------------------------
284  void LArVoxelReadout::SetRandomEngines(CLHEP::HepRandomEngine* pPropGen) {
285  assert(pPropGen); // random engine must be present
286  fPropGen = pPropGen;
287  } // LArVoxelReadout::SetRandomEngines()
288 
289 
290  //----------------------------------------------------------------------------
292  (geo::Point_t const& pos, geo::PlaneGeo const& plane) const
293  {
294  //
295  // translate the landing position on the two frame coordinates
296  // ("width" and "depth")
297  //
298  auto const landingPos = plane.PointWidthDepthProjection(pos);
299 
300  //
301  // compute the distance of the landing position on the two frame
302  // coordinates ("width" and "depth");
303  // keep the point within 10 micrometers (0.001 cm) from the border
304  //
305  auto const offPlane = plane.DeltaFromActivePlane(landingPos, 0.001);
306 
307  //
308  // if both the distances are below the margin, move the point to
309  // the border
310  //
311 
312  // nothing to recover: landing is inside
313  if ((offPlane.X() == 0.0) && (offPlane.Y() == 0.0)) return pos;
314 
315  // won't recover: too far
316  if ((std::abs(offPlane.X()) > fOffPlaneMargin)
317  || (std::abs(offPlane.Y()) > fOffPlaneMargin))
318  return pos;
319 
320  // we didn't fully decompose because it might be unnecessary;
321  // now we need the full thing
322  auto const distance = plane.DistanceFromPlane(pos);
323 
324  return plane.ComposePoint<geo::Point_t>(distance, landingPos + offPlane);
325 
326  } // LArVoxelReadout::RecoverOffPlaneDeposit()
327 
328 
329  //----------------------------------------------------------------------------
330  // energy is passed in with units of MeV, dx has units of cm
331  void LArVoxelReadout::DriftIonizationElectrons(G4ThreeVector stepMidPoint,
332  const double simTime,
333  int trackID,
334  unsigned short int cryostat, unsigned short int tpc)
335  {
336  auto const * ts = lar::providerFrom<detinfo::DetectorClocksService>();
337 
338  // this must be always true, unless caller has been sloppy
339  assert(fPropGen); // No propagation random generator provided?!
340 
341  CLHEP::RandGauss PropRand(*fPropGen);
342 
343  // This routine gets called frequently, once per every particle
344  // traveling through every voxel. Use whatever tricks we can to
345  // increase its execution speed.
346 
347  static double LifetimeCorr_const = -1000. * fElectronLifetime;
348  static double LDiff_const = std::sqrt(2.*fLongitudinalDiffusion);
349  static double TDiff_const = std::sqrt(2.*fTransverseDiffusion);
350  static double RecipDriftVel[3] = {1./fDriftVelocity[0],
351  1./fDriftVelocity[1],
352  1./fDriftVelocity[2]};
353 
354  struct Deposit_t {
355  double energy = 0.;
356  double electrons = 0.;
357 
358  void add(double more_energy, double more_electrons)
359  { energy += more_energy; electrons += more_electrons; }
360  }; // Deposit_t
361 
362  // Map of electrons to store - catalogued by map[channel][tdc]
363  std::map<raw::ChannelID_t, std::map<unsigned int, Deposit_t>> DepositsToStore;
364 
365  double xyz1[3] = {0.};
366 
367  double const xyz[3] = {stepMidPoint.x() / CLHEP::cm,
368  stepMidPoint.y() / CLHEP::cm,
369  stepMidPoint.z() / CLHEP::cm};
370 
371  // Already know which TPC we're in because we have been told
372 
373  try{
374  const geo::TPCGeo &tpcg = fGeoHandle->TPC(tpc, cryostat);
375 
376  // X drift distance - the drift direction can be either in
377  // the positive or negative direction, so use std::abs
378 
380  double XDrift = std::abs(stepMidPoint.x()/CLHEP::cm - tpcg.PlaneLocation(0)[0]);
381  //std::cout<<tpcg.DriftDirection()<<std::endl;
382  if (tpcg.DriftDirection() == geo::kNegX)
383  XDrift = stepMidPoint.x()/CLHEP::cm - tpcg.PlaneLocation(0)[0];
384  else if (tpcg.DriftDirection() == geo::kPosX)
385  XDrift = tpcg.PlaneLocation(0)[0] - stepMidPoint.x()/CLHEP::cm;
386 
387  if(XDrift < 0.) return;
388 
389  // Get SCE {x,y,z} offsets for particular location in TPC
390  geo::Vector_t posOffsets;
391  auto const* SCE = lar::providerFrom<spacecharge::SpaceChargeService>();
392  if (SCE->EnableSimSpatialSCE() == true)
393  {
394  posOffsets = SCE->GetPosOffsets({ xyz[0], xyz[1], xyz[2] });
395  }
396  posOffsets.SetX(-posOffsets.X());
397 
398  // Drift time (nano-sec)
399  double TDrift;
400  XDrift += posOffsets.X();
401 
402  // Space charge distortion could push the energy deposit beyond the wire
403  // plane (see issue #15131). Given that we don't have any subtlety in the
404  // simulation of this region, bringing the deposit exactly on the plane
405  // should be enough for the time being.
406  if (XDrift < 0.) XDrift = 0.;
407 
408  TDrift = XDrift * RecipDriftVel[0];
409  if (tpcg.Nplanes() == 2){// special case for ArgoNeuT (plane 0 is the second wire plane)
410  TDrift = ((XDrift - tpcg.PlanePitch(0,1)) * RecipDriftVel[0]
411  + tpcg.PlanePitch(0,1) * RecipDriftVel[1]);
412  }
413 
414  const double lifetimecorrection = TMath::Exp(TDrift / LifetimeCorr_const);
417 
418  // if we have no electrons (too small energy or too large recombination)
419  // we are done already here
420  if (nIonizedElectrons <= 0) {
421  LOG_DEBUG("LArVoxelReadout")
422  << "No electrons drifted to readout, " << energy << " MeV lost.";
423  return;
424  }
425  // includes the effect of lifetime
426  const double nElectrons = nIonizedElectrons * lifetimecorrection;
427 
428  // Longitudinal & transverse diffusion sigma (cm)
429  double SqrtT = std::sqrt(TDrift);
430  double LDiffSig = SqrtT * LDiff_const;
431  double TDiffSig = SqrtT * TDiff_const;
432  double electronclsize = fElectronClusterSize;
433 
434  int nClus = (int) std::ceil(nElectrons / electronclsize);
435  if (nClus < fMinNumberOfElCluster)
436  {
437  electronclsize = nElectrons / fMinNumberOfElCluster;
438  if (electronclsize < 1.0)
439  {
440  electronclsize = 1.0;
441  }
442  nClus = (int) std::ceil(nElectrons / electronclsize);
443  }
444 
445  // Compute arrays of values as quickly as possible.
446  std::vector< double > XDiff(nClus);
447  std::vector< double > YDiff(nClus);
448  std::vector< double > ZDiff(nClus);
449  std::vector< double > nElDiff(nClus, electronclsize);
450  std::vector< double > nEnDiff(nClus);
451 
452  // fix the number of electrons in the last cluster, that has smaller size
453  nElDiff.back() = nElectrons - (nClus-1)*electronclsize;
454 
455  for(size_t xx = 0; xx < nElDiff.size(); ++xx){
456  if(nElectrons > 0) nEnDiff[xx] = energy/nElectrons*nElDiff[xx];
457  else nEnDiff[xx] = 0.;
458  }
459 
460  double const avegageYtransversePos
461  = (stepMidPoint.y()/CLHEP::cm) + posOffsets.Y();
462  double const avegageZtransversePos
463  = (stepMidPoint.z()/CLHEP::cm) + posOffsets.Z();
464 
465  // Smear drift times by x position and drift time
466  if (LDiffSig > 0.0)
467  PropRand.fireArray( nClus, &XDiff[0], 0., LDiffSig);
468  else
469  XDiff.assign(nClus, 0.0);
470 
471  if (TDiffSig > 0.0) {
472  // Smear the Y,Z position by the transverse diffusion
473  PropRand.fireArray( nClus, &YDiff[0], avegageYtransversePos, TDiffSig);
474  PropRand.fireArray( nClus, &ZDiff[0], avegageZtransversePos, TDiffSig);
475  }
476  else {
477  YDiff.assign(nClus, avegageYtransversePos);
478  ZDiff.assign(nClus, avegageZtransversePos);
479  }
480 
481  // make a collection of electrons for each plane
482  for(size_t p = 0; p < tpcg.Nplanes(); ++p){
483 
484  geo::PlaneGeo const& plane = tpcg.Plane(p);
485 
486  double Plane0Pitch = tpcg.Plane0Pitch(p);
487 
488  // "-" sign is because Plane0Pitch output is positive. Andrzej
489  xyz1[0] = tpcg.PlaneLocation(0)[0] - Plane0Pitch;
490 
491  // Drift nClus electron clusters to the induction plane
492  for(int k = 0; k < nClus; ++k){
493  // Correct drift time for longitudinal diffusion and plane
494  double TDiff = TDrift + XDiff[k] * RecipDriftVel[0];
495  // Take into account different Efields between planes
496  // Also take into account special case for ArgoNeuT where Nplanes = 2.
497  for (size_t ip = 0; ip<p; ++ip){
498  TDiff += tpcg.PlanePitch(ip,ip+1) * RecipDriftVel[tpcg.Nplanes()==3?ip+1:ip+2];
499  }
500  xyz1[1] = YDiff[k];
501  xyz1[2] = ZDiff[k];
502 
504 
505  // grab the nearest channel to the xyz1 position
506  try{
507  if (fOffPlaneMargin != 0) {
508  // get the effective position where to consider the charge landed;
509  //
510  // Some optimisations are possible; in particular, this method
511  // could be extended to inform us if the point was too far.
512  // Currently, if that is the case the code will proceed, find the
513  // point is off plane, emit a warning and skip the deposition.
514  //
515  auto const landingPos
516  = RecoverOffPlaneDeposit({ xyz1[0], xyz1[1], xyz1[2] }, plane);
517  xyz1[0] = landingPos.X();
518  xyz1[1] = landingPos.Y();
519  xyz1[2] = landingPos.Z();
520 
521  } // if charge lands off plane
522  uint32_t channel = fGeoHandle->NearestChannel(xyz1, p, tpc, cryostat);
523 
526  // Add potential decay/capture/etc delay effect, simTime.
527  unsigned int tdc = fClock.Ticks(ts->G4ToElecTime(TDiff + simTime));
528 
529  // Add electrons produced by each cluster to the map
530  DepositsToStore[channel][tdc].add(nEnDiff[k], nElDiff[k]);
531  }
532  catch(cet::exception &e){
533  mf::LogWarning("LArVoxelReadout") << "unable to drift electrons from point ("
534  << xyz[0] << "," << xyz[1] << "," << xyz[2]
535  << ") with exception " << e;
536  }
537  } // end loop over clusters
538  } // end loop over planes
539 
540  // Now store them in SimChannels
541  ChannelMap_t& ChannelDataMap = fChannelMaps[cryostat][tpc];
542 
543  // browse deposited data on each channel: (channel; deposit data in time)
544  for(auto const& deposit_per_channel: DepositsToStore){
545 
546  raw::ChannelID_t channel = deposit_per_channel.first;
547 
548  // find whether we already have this channel
549  auto iChannelData = ChannelDataMap.find(channel);
550 
551  // channelData is the SimChannel these deposits are going to be added to
552  // If there is such a channel already, use it (first beanch).
553  // If it's a new channel, the inner assignment creates a new SimChannel
554  // in the map, and we save its reference in channelData
555  sim::SimChannel& channelData
556  = (iChannelData == ChannelDataMap.end())
557  ? (ChannelDataMap[channel] = sim::SimChannel(channel))
558  : iChannelData->second;
559 
560  // go through all deposits, one for each TDC: (TDC, deposit data)
561  for(auto const& deposit_per_tdc: deposit_per_channel.second) {
562  channelData.AddIonizationElectrons(trackID,
563  deposit_per_tdc.first,
564  deposit_per_tdc.second.electrons,
565  xyz,
566  deposit_per_tdc.second.energy);
567 
568  } // for deposit on TDCs
569  } // for deposit on channels
570 
571  } // end try intended to catch points where TPC can't be found
572  catch(cet::exception &e){
573  mf::LogWarning("LArVoxelReadout") << "step cannot be found in a TPC\n"
574  << e;
575  }
576 
577  return;
578  }
579 
580  //---------------------------------------------------------------------------------------
581  // Never used but still have to be defined for G4
584 
585 } // namespace larg4
CLHEP::HepRandomEngine * propGen
Random engine for charge propagation.
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:167
unsigned int fTPC
which TPC this LArVoxelReadout corresponds to
double PlanePitch(unsigned int p1=0, unsigned int p2=1) const
Returns the center of the TPC volume in world coordinates [cm].
Definition: TPCGeo.cxx:418
Double_t xx
Definition: macro.C:12
virtual const ::detinfo::ElecClock & TPCClock() const =0
Lends a constant TPC clock with time set to trigger time.
Collection of what it takes to set a LArVoxelReadout up.
Energy deposited on a readout channel by simulated tracks.
Definition: SimChannel.h:143
A G4PVPlacement with an additional identificator.
unsigned int Nplanes() const
Number of planes in this tpc.
Definition: TPCGeo.h:145
Geometry information for a single TPC.
Definition: TPCGeo.h:37
const ChannelMap_t & GetSimChannelMap() const
Returns the accumulated channel -> SimChannel map for the single TPC.
LArVoxelReadout(std::string const &name)
Constructor. Can detect which TPC to cover by the name.
art::ServiceHandle< sim::LArG4Parameters > fLgpHandle
Handle to the LArG4 parameters service.
Point ComposePoint(WireDecomposedVector_t const &decomp) const
Returns the 3D vector from composition of projection and distance.
Definition: PlaneGeo.h:818
Geant4 interface.
void Setup(Setup_t const &setupData)
Reads all the configuration elements from setupData
WidthDepthProjection_t PointWidthDepthProjection(geo::Point_t const &point) const
Returns the projection of the specified point on the plane.
Definition: PlaneGeo.h:897
void SetRandomEngines(CLHEP::HepRandomEngine *pPropGen)
Sets the random generators to be used.
int Ticks() const
Current clock tick (that is, the number of tick Time() falls in).
Definition: ElecClock.h:235
art::ServiceHandle< geo::Geometry > fGeoHandle
Handle to the Geometry service.
Drift towards negative X values.
Definition: geo_types.h:109
unsigned int Ncryostats() const
Returns the number of cryostats in the detector.
bool bSingleTPC
true if this readout is associated with a single TPC
ID_t ID
Physical Volume identificator.
geo::Point_t RecoverOffPlaneDeposit(geo::Point_t const &pos, geo::PlaneGeo const &plane) const
Returns the point on the specified plane closest to position.
double TransverseDiffusion() const
void SetSingleTPC(unsigned int cryostat, unsigned int tpc)
Associates this readout to one specific TPC.
std::map< unsigned int, sim::SimChannel > ChannelMap_t
Type of map channel -> sim::SimChannel.
Use Geant4&#39;s user "hooks" to maintain a list of particles generated by Geant4.
bool NoElectronPropagation() const
std::vector< unsigned short int > fSkipWireSignalInTPCs
double ElectronClusterSize() const
Geometry information for a single wire plane.The plane is represented in the geometry by a solid whic...
Definition: PlaneGeo.h:78
double energy
Definition: plottest35.C:25
static IonizationAndScintillation * Instance()
std::vector< sim::SimChannel > GetSimChannels() const
Creates a list with the accumulated information for the single TPC.
bool Has(std::vector< unsigned short int > v, unsigned short int tpc) const
virtual G4bool ProcessHits(G4Step *, G4TouchableHistory *)
double Plane0Pitch(unsigned int p) const
Returns the center of the TPC volume in world coordinates [cm].
Definition: TPCGeo.cxx:354
DriftDirection_t DriftDirection() const
Returns an enumerator value describing the drift direction.
Definition: TPCGeo.h:127
std::vector< std::vector< ChannelMap_t > > fChannelMaps
Maps of cryostat, tpc to channel data.
virtual void EndOfEvent(G4HCofThisEvent *)
double fOffPlaneMargin
Charge deposited within this many [cm] from the plane is lead onto it.
void SetDiscoverTPC()
Sets this readout to discover the TPC of each processed hit.
Conversion of times between different formats and references.
const std::vector< unsigned short int > SkipWireSignalInTPCs() const
unsigned int NTPC(unsigned int cstat=0) const
Returns the total number of TPCs in the specified cryostat.
void AddIonizationElectrons(TrackID_t trackID, TDC_t tdc, double numberElectrons, double const *xyz, double energy)
Add ionization electrons and energy to this channel.
Definition: SimChannel.cxx:52
void DriftIonizationElectrons(G4ThreeVector stepMidPoint, const double simTime, int trackID, unsigned short int cryostat, unsigned short int tpc)
A Geant4 sensitive detector that accumulates voxel information.
double offPlaneMargin
Margin for charge recovery (see LArVoxelReadout).
WidthDepthProjection_t DeltaFromActivePlane(WidthDepthProjection_t const &proj, double wMargin, double dMargin) const
Returns a projection vector that, added to the argument, gives a projection inside (or at the border ...
Definition: PlaneGeo.cxx:593
Drift towards positive X values.
Definition: geo_types.h:108
int MinNumberOfElCluster() const
double DistanceFromPlane(geo::Point_t const &point) const
Returns the distance of the specified point from the wire plane.
Definition: PlaneGeo.h:557
unsigned int fCstat
and in which cryostat (if bSingleTPC is true)
void SetOffPlaneChargeRecoveryMargin(double margin)
Sets the margin for recovery of charge drifted off-plane.
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
CLHEP::HepRandomEngine * fPropGen
random engine for charge propagation
virtual void Initialize(G4HCofThisEvent *)
#define LOG_DEBUG(id)
::detinfo::ElecClock fClock
TPC electronics clock.
TPCGeo const & TPC(unsigned int const tpc=0, unsigned int const cstat=0) const
Returns the specified TPC.
PlaneGeo const & Plane(geo::View_t view) const
Return the plane in the tpc with View_t view.
Definition: TPCGeo.cxx:298
double LongitudinalDiffusion() const
unsigned int ChannelID_t
Type representing the ID of a readout channel.
Definition: RawTypes.h:27
Float_t e
Definition: plot.C:34
raw::ChannelID_t NearestChannel(geo::Point_t const &worldLoc, geo::PlaneID const &planeid) const
Returns the ID of the channel nearest to the specified position.
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:187
const double * PlaneLocation(unsigned int p) const
Returns the coordinates of the center of the specified plane [cm].
Definition: TPCGeo.cxx:412
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
bool DisableWireplanes() const