LArVoxelReadout.h

Go to the documentation of this file.

#ifndef LArG4_LArVoxelReadout_h
#define LArG4_LArVoxelReadout_h

#include <algorithm> // std::max()
#include <stddef.h>
#include <vector>

#include "Geant4/G4PVPlacement.hh"
#include "Geant4/G4VSensitiveDetector.hh"

#include "art/Framework/Services/Registry/ServiceHandle.h"
#include "larcore/Geometry/Geometry.h"
#include "larcoreobj/SimpleTypesAndConstants/geo_vectors.h"
#include "lardataobj/Simulation/SimChannel.h"
#include "larsim/Simulation/LArG4Parameters.h"
namespace detinfo {
  class DetectorClocksData;
  class DetectorPropertiesData;
}

// Forward declarations
class G4HCofThisEvent;
class G4TouchableHistory;
class G4Step;

namespace CLHEP {
  class HepRandomEngine;
}

namespace larg4 {

  struct TPCID_t {
    unsigned short int Cryostat, TPC;
    bool operator<(const TPCID_t& than) const
    {
      return (Cryostat < than.Cryostat) || ((Cryostat == than.Cryostat) && (TPC < than.TPC));
    } // operator< ()
  };  // TPCID_t

  template <class IDTYPE>
  class G4PVPlacementWithID : public G4PVPlacement {
  public:
    typedef IDTYPE ID_t;

    ID_t ID; 

    G4PVPlacementWithID(const G4Transform3D& Transform3D,
                        const G4String& pName,
                        G4LogicalVolume* pLogical,
                        G4VPhysicalVolume* pMother,
                        G4bool pMany,
                        G4int pCopyNo,
                        G4bool pSurfChk = false,
                        ID_t id = ID_t())
      : G4PVPlacement(Transform3D, pName, pLogical, pMother, pMany, pCopyNo, pSurfChk), ID(id)
    {}

    virtual ~G4PVPlacementWithID() {}
  }; // G4PVPlacementWithID<>

  typedef G4PVPlacementWithID<TPCID_t> G4PVPlacementInTPC;

  class LArVoxelReadout : public G4VSensitiveDetector {
  public:
    typedef std::map<unsigned int, sim::SimChannel> ChannelMap_t;

    struct Setup_t {

      CLHEP::HepRandomEngine* propGen = nullptr;

      double offPlaneMargin = 0.0;
    }; // struct Setup_t

    LArVoxelReadout(std::string const& name);

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

    void Setup(Setup_t const& setupData);

    void SetSingleTPC(unsigned int cryostat, unsigned int tpc);

    void SetDiscoverTPC();

    // Required for classes that inherit from G4VSensitiveDetector.
    //
    // Called at start and end of each event.
    virtual void Initialize(G4HCofThisEvent*);
    virtual void EndOfEvent(G4HCofThisEvent*);

    // Called to clear any accumulated information.
    virtual void clear();

    // The key method of this class.  It's called by Geant4 for each
    // step within the read-out geometry.  It accumulates the energy
    // in the G4Step in the LArVoxelList.
    virtual G4bool ProcessHits(G4Step*, G4TouchableHistory*);

    // Empty methods; they have to be defined, but they're rarely
    // used in Geant4 applications.
    virtual void DrawAll();
    virtual void PrintAll();

    // Independent method; clears the vector of SimChannels as well as the
    // channel number to SimChannel map.  Has to be separate from the
    // clear method above because that run is run for every G4 event, ie
    // each MCTruth in the art::Event, while we want to only run this at
    // the end of the G4 processing for each art::Event.
    void ClearSimChannels();

    std::vector<sim::SimChannel> GetSimChannels() const;

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

    const ChannelMap_t& GetSimChannelMap() const;
    ChannelMap_t& GetSimChannelMap();

    const ChannelMap_t& GetSimChannelMap(unsigned short cryo, unsigned short tpc) const;
    ChannelMap_t& GetSimChannelMap(unsigned short cryo, unsigned short tpc);

  private:
    // N.B. This code is not thread-safe, as it presupposes that there
    // is a "current" clock-data object.  Such a pattern should be
    // avoided for the users of larg4.
    friend class LArVoxelReadoutGeometry;

    void SetClockData(detinfo::DetectorClocksData const* const clockData) noexcept
    {
      fClockData = clockData;
    }

    void SetPropertiesData(detinfo::DetectorPropertiesData const* const detProp) noexcept
    {
      fDetProp = detProp;
    }

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

    void SetRandomEngines(CLHEP::HepRandomEngine* pPropGen);

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

    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
    {
      for (auto c : v)
        if (c == tpc) return true;
      return false;
    }

    // Used in electron-cluster calculations
    // External parameters for the electron-cluster calculation.
    // obtained from LArG4Parameters, LArProperties, and DetectorProperties services
    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;

    std::vector<std::vector<ChannelMap_t>> fChannelMaps; 
    art::ServiceHandle<geo::Geometry const> fGeoHandle;  
    art::ServiceHandle<sim::LArG4Parameters const>
      fLgpHandle;        
    unsigned int fTPC;   
    unsigned int fCstat; 
    bool bSingleTPC;     

    CLHEP::HepRandomEngine* fPropGen = nullptr; 

    detinfo::DetectorClocksData const* fClockData{nullptr};
    detinfo::DetectorPropertiesData const* fDetProp{nullptr};

    //these are the things for doing the separated EDeps
    void ProcessStep(G4Step*);

    G4ThreeVector fStepStart;
    G4ThreeVector fStepEnd;
    size_t fNSteps;
  };

}

#endif // LArG4_LArVoxelReadout_h