DetectorPropertiesStandard.h

Go to the documentation of this file.

// \file DetectorPropertiesStandard.h
//
// \brief service to contain information about detector electronics, etc
//
// \author brebel@fnal.gov
//
// Separation of service from Detector info class:
// jpaley@fnal.gov
#ifndef DETINFO_DETECTORPROPERTIESSTD_H
#define DETINFO_DETECTORPROPERTIESSTD_H

// LArSoft libraries
#include "larcorealg/CoreUtils/ProviderPack.h"
#include "larcorealg/Geometry/GeometryCore.h"
#include "larcorealg/Geometry/WireReadoutGeom.h"
#include "larcoreobj/SimpleTypesAndConstants/PhysicalConstants.h"
#include "lardataalg/DetectorInfo/DetectorClocks.h"
#include "lardataalg/DetectorInfo/DetectorProperties.h"
#include "lardataalg/DetectorInfo/DetectorPropertiesData.h"
#include "lardataalg/DetectorInfo/LArProperties.h"

// framework libraries
#include "fhiclcpp/fwd.h"
#include "fhiclcpp/types/Atom.h"
#include "fhiclcpp/types/OptionalAtom.h"
#include "fhiclcpp/types/Sequence.h"

// C/C++ standard libraries
#include <set>

namespace detinfo {

  class DetectorPropertiesStandard final : public DetectorProperties {
  public:
    using providers_type =
      lar::ProviderPack<geo::GeometryCore, geo::WireReadoutGeom, detinfo::LArProperties>;

    struct Configuration_t {
      using Name = fhicl::Name;
      using Comment = fhicl::Comment;

      fhicl::Sequence<double> Efield{
        Name("Efield"),
        Comment("electric field in front of each wire plane (the last one is "
                "the big one!) [kV/cm]")};

      fhicl::Atom<double> Electronlifetime{Name("Electronlifetime"),
                                           Comment("electron lifetime in liquid argon [us]")};
      fhicl::Atom<double> Temperature{Name("Temperature"), Comment("argon temperature [K]")};
      fhicl::Atom<double> ElectronsToADC{
        Name("ElectronsToADC"),
        Comment("conversion factor: (ADC counts)/(ionization electrons)")};
      fhicl::Atom<unsigned int> NumberTimeSamples{
        Name("NumberTimeSamples"),
        Comment("number of TPC readout TDC clock ticks per event")};
      fhicl::Atom<unsigned int> ReadOutWindowSize{
        Name("ReadOutWindowSize"),
        Comment("number of TPC readout TDC clock ticks per readout window")};

      // The following are not really "optional": the ones for the views which
      // are present are mandatory.
      fhicl::OptionalAtom<double> TimeOffsetU{
        Name("TimeOffsetU"),
        Comment("tick offset subtracted to to convert spacepoint coordinates "
                "to hit times on view U")};
      fhicl::OptionalAtom<double> TimeOffsetV{
        Name("TimeOffsetV"),
        Comment("tick offset subtracted to to convert spacepoint coordinates "
                "to hit times on view V")};
      fhicl::OptionalAtom<double> TimeOffsetZ{
        Name("TimeOffsetZ"),
        Comment("tick offset subtracted to to convert spacepoint coordinates "
                "to hit times on view Z")};
      fhicl::OptionalAtom<double> TimeOffsetY{
        Name("TimeOffsetY"),
        Comment("tick offset subtracted to to convert spacepoint coordinates "
                "to hit times on view Y")};
      fhicl::OptionalAtom<double> TimeOffsetX{
        Name("TimeOffsetX"),
        Comment("tick offset subtracted to to convert spacepoint coordinates "
                "to hit times on view X")};

      fhicl::Atom<double> SternheimerA{
        Name("SternheimerA"),
        Comment("parameter a of Sternheimer correction delta = 2log(10) x - "
                "cbar + { a (x1-x)^k } theta(x1-x), x = log10(p/m)")};
      fhicl::Atom<double> SternheimerK{
        Name("SternheimerK"),
        Comment("parameter k of Sternheimer correction delta = 2log(10) x - "
                "cbar + { a (x_1-x)^k } theta(x1-x), x = log10(p/m)")};
      fhicl::Atom<double> SternheimerX0{
        Name("SternheimerX0"),
        Comment("minimum x = log10(p/m) for the application of Sternheimer "
                "correction")};
      fhicl::Atom<double> SternheimerX1{
        Name("SternheimerX1"),
        Comment("parameter x_1 of Sternheimer correction delta = 2log(10) x - "
                "cbar + { a (x_1-x)^k } theta(x1-x), x = log10(p/m)")};
      fhicl::Atom<double> SternheimerCbar{
        Name("SternheimerCbar"),
        Comment("parameter cbar of Sternheimer correction delta = 2log(10) x - "
                "cbar + { a (x_1-x)^k } theta(x1-x), x = log10(p/m)")};
      fhicl::Atom<double> DriftVelFudgeFactor{
        Name("DriftVelFudgeFactor"),
        Comment("Allows a scaling factor to fudge the drift velocity "
                "calculation (as suggested by DriftVel Stancari")};

      fhicl::Atom<bool> UseIcarusMicrobooneDriftModel{
        Name("UseIcarusMicrobooneDriftModel"),
        Comment("Allows user to decide to use the ICARUS+MicroBooNE drift "
                "model for velocity calculation as in arXiv:2008.09765"),
        false};

      fhicl::Atom<bool> IncludeInterPlanePitchInXTickOffsets{
        Name("IncludeInterPlanePitchInXTickOffsets"),
        Comment("Historically, ConvertTicksToX has allowed for the drift time "
                "between the wire planes. This is appropriate for "
                "recob::RawDigits, and recob::Wires from the 1D unfolding, "
                "but is not appropriate for recob::Wires from WireCell. "
                "When disabled, the conversion propagates to the last plane, "
                "rather than the first plane and assumes a standard drift velocity. "
                "The default value is 'true', retaining the 'classic' behaviour"),
        true};

      fhicl::Atom<bool> SimpleBoundary{Name("SimpleBoundaryProcess"), Comment("")};

      fhicl::Atom<double> ModBoxAlpha{
        Name("ModBoxAlpha"),
        Comment("alpha parameter in the Modified Box recombination model."),
        util::kModBoxA};

      fhicl::Atom<double> ModBoxBeta{
        Name("ModBoxBeta"),
        Comment("beta parameter in the Modified Box recombination model."),
        util::kModBoxB};

    }; // Configuration_t

    DetectorPropertiesStandard(fhicl::ParameterSet const& pset,
                               geo::GeometryCore const* geo,
                               geo::WireReadoutGeom const* wireReadoutGeom,
                               detinfo::LArProperties const* lp,
                               std::set<std::string> const& ignore_params = {});

    DetectorPropertiesStandard(DetectorPropertiesStandard const&) = delete;
    virtual ~DetectorPropertiesStandard() = default;

    void SetNumberTimeSamples(unsigned int nsamp) { fNumberTimeSamples = nsamp; }

    // Accessors.

    double Efield(unsigned int planegap = 0) const override; 

    double DriftVelocity(double efield = 0.,
                         double temperature = 0.) const override; 

    double BirksCorrection(double dQdX) const override;
    double BirksCorrection(double dQdX, double EField) const override;
    double ModBoxCorrection(double dQdX) const override;
    double ModBoxCorrection(double dQdX, double EField) const override;

    double ElectronLifetime() const override
    {
      return fElectronlifetime; //< microseconds
    }

    double Density(double temperature = 0.) const override; 

    double Temperature() const override { return fTemperature; }

    double Eloss(double mom, double mass, double tcut) const override;

    double ElossVar(double mom, double mass) const override;

    double ElectronsToADC() const override { return fElectronsToADC; }
    unsigned int NumberTimeSamples() const override { return fNumberTimeSamples; }
    unsigned int ReadOutWindowSize() const override { return fReadOutWindowSize; }
    double TimeOffsetU() const override { return fTimeOffsetU; };
    double TimeOffsetV() const override { return fTimeOffsetV; };
    double TimeOffsetZ() const override { return fTimeOffsetZ; };
    double TimeOffsetY() const override { return fTimeOffsetY; };

    bool SimpleBoundary() const override { return fSimpleBoundary; }

    DetectorPropertiesData DataFor(detinfo::DetectorClocksData const& clock_data) const override;

  private:
    void ValidateAndConfigure(fhicl::ParameterSet const& p,
                              std::set<std::string> const& ignore_params);

    struct SternheimerParameters_t {
      double a;    
      double k;    
      double x0;   
      double x1;   
      double cbar; 
    };

    // service providers we depend on;
    // in principle could be replaced by a single providerpack_type.
    detinfo::LArProperties const* fLP;
    geo::GeometryCore const* fGeo;
    geo::WireReadoutGeom const* fChannelMap;

    std::vector<double> fEfield;     
    double fElectronlifetime;        
    double fTemperature;             
    double fElectronsToADC;          
    unsigned int fNumberTimeSamples; 
    unsigned int fReadOutWindowSize; 
    double fTimeOffsetU;             
    double fTimeOffsetV;             
    double fTimeOffsetZ;             
    double fTimeOffsetY;             
    double fTimeOffsetX;             
    double fDriftVelFudgeFactor;     

    bool fUseIcarusMicrobooneDriftModel; 

    bool fIncludeInterPlanePitchInXTickOffsets;

    SternheimerParameters_t fSternheimerParameters; 

    std::vector<std::vector<double>> fDriftDirection;

    bool fSimpleBoundary;

    double fModBoxA;
    double fModBoxB;

  }; // class DetectorPropertiesStandard
} // namespace detinfo

#endif // DETINFO_DETECTOR_PROPERTIES_H