SiLinearChargeDivider Class Reference

#include <SiLinearChargeDivider.h>

Inheritance diagram for SiLinearChargeDivider:

Public Member Functions
SiChargeDivider::ionization_type	divide (const PSimHit *, const LocalVector &, double, const StripGeomDetUnit &det, CLHEP::HepRandomEngine *) override
void	setParticleDataTable (const ParticleDataTable *pdt) override
	SiLinearChargeDivider (const edm::ParameterSet &conf)
	~SiLinearChargeDivider () override
Public Member Functions inherited from SiChargeDivider
virtual	~SiChargeDivider ()

Private Member Functions
float	driftXPos (const Local3DPoint &pos, const LocalVector &drift, double thickness)
void	fluctuateEloss (double const particleMass, float momentum, float eloss, float length, int NumberOfSegmentation, float elossVector[], CLHEP::HepRandomEngine *)
void	readPulseShape (const std::string &pulseShapeFileName)
float	TimeResponse (const PSimHit *hit, const StripGeomDetUnit &det)

Private Attributes
const int	chargedivisionsPerStrip
const double	cosmicShift
const double	deltaCut
std::unique_ptr< SiG4UniversalFluctuation >	fluctuate
const bool	fluctuateCharge
const bool	peakMode
double	pulseResolution
unsigned int	pulset0Idx
std::vector< double >	pulseValues
const ParticleDataTable *	theParticleDataTable

Additional Inherited Members
Public Types inherited from SiChargeDivider
typedef std::vector< EnergyDepositUnit >	ionization_type

Detailed Description

Concrete implementation of SiChargeDivider. It divides the charge on the line connecting entry and exit point of the SimTrack in the Silicon. Effects that are considered here are:

• fluctuations of the charge deposit fronm segment to segment
• pulse shape ( in peak of deconvolution mode)

Definition at line 27 of file SiLinearChargeDivider.h.

Constructor & Destructor Documentation

SiLinearChargeDivider()

SiLinearChargeDivider::SiLinearChargeDivider

(

const edm::ParameterSet &

conf

)

Definition at line 9 of file SiLinearChargeDivider.cc.

References edm::FileInPath::fullPath(), edm::ParameterSet::getParameter(), peakMode, and readPulseShape().

    :  // Run APV in peak instead of deconvolution mode, which degrades the time resolution.
      peakMode(conf.getParameter<bool>("APVpeakmode")),
      // Enable interstrip Landau fluctuations within a cluster.
      fluctuateCharge(conf.getParameter<bool>("LandauFluctuations")),
      // Number of segments per strip into which charge is divided during
      // simulation. If large, precision of simulation improves.
      chargedivisionsPerStrip(conf.getParameter<int>("chargeDivisionsPerStrip")),
      // delta cutoff in MeV, has to be same as in Geant (0.120425 MeV corresponding to 100um range for electrons)
      deltaCut(conf.getParameter<double>("DeltaProductionCut")),
      //Offset for digitization during the MTCC and in general for taking cosmic particle
      //The value to be used it must be evaluated and depend on the volume defnition used
      //for the cosimc generation (Considering only the tracker the value is 11 ns)
      cosmicShift(conf.getUntrackedParameter<double>("CosmicDelayShift")),
      theParticleDataTable(nullptr),
      // Geant4 engine used to fluctuate the charge from segment to segment
      fluctuate(new SiG4UniversalFluctuation())

{
  readPulseShape(conf.getParameter<edm::FileInPath>(peakMode ? "APVShapePeakFile" : "APVShapeDecoFile").fullPath());
}

~SiLinearChargeDivider()

SiLinearChargeDivider::~SiLinearChargeDivider ( )

override

Definition at line 31 of file SiLinearChargeDivider.cc.

{}

Member Function Documentation

divide()

SiChargeDivider::ionization_type SiLinearChargeDivider::divide ( const PSimHit *  hit, const LocalVector &  driftdir, double  moduleThickness, const StripGeomDetUnit &  det, CLHEP::HepRandomEngine *  engine  )

overridevirtual

Implements SiChargeDivider.

Definition at line 71 of file SiLinearChargeDivider.cc.

References cms::cuda::assert(), chargedivisionsPerStrip, driftXPos(), simKBmtfDigis_cfi::eLoss, nano_mu_digi_cff::float, fluctuateCharge, fluctuateEloss(), mps_fire::i, StripTopology::localPitch(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), StripGeomDetUnit::specificTopology(), theParticleDataTable, and TimeResponse().

                                                                                             {
  // signal after pulse shape correction
  float const decSignal = TimeResponse(hit, det);

  // if out of time go home!
  if (0 == decSignal)
    return ionization_type();

  // Get the nass if the particle, in MeV.
  // Protect from particles with Mass = 0, assuming then the pion mass
  assert(theParticleDataTable != nullptr);
  ParticleData const* particle = theParticleDataTable->particle(hit->particleType());
  double const particleMass = particle ? particle->mass() * 1000 : 139.57;
  double const particleCharge = particle ? particle->charge() : 1.;

  if (!particle) {
    LogDebug("SiLinearChargeDivider") << "Cannot find particle of type " << hit->particleType()
                                      << " in the PDT we assign to this particle the mass and charge of the Pion";
  }

  int NumberOfSegmentation =
      // if neutral: just one deposit....
      (fabs(particleMass) < 1.e-6 || particleCharge == 0)
          ? 1
          :
          // computes the number of segments from number of segments per strip times number of strips.
          (int)(1 + chargedivisionsPerStrip *
                        fabs(driftXPos(hit->exitPoint(), driftdir, moduleThickness) -
                             driftXPos(hit->entryPoint(), driftdir, moduleThickness)) /
                        det.specificTopology().localPitch(hit->localPosition()));

  // Eloss in GeV
  float eLoss = hit->energyLoss();

  // Prepare output
  ionization_type _ionization_points;
  _ionization_points.resize(NumberOfSegmentation);

  // Fluctuate charge in track subsegments
  LocalVector direction = hit->exitPoint() - hit->entryPoint();
  if (NumberOfSegmentation <= 1) {
    // here I need a random... not 0.5
    _ionization_points[0] = EnergyDepositUnit(eLoss * decSignal / eLoss, hit->entryPoint() + 0.5f * direction);
  } else {
    float eLossVector[NumberOfSegmentation];
    if (fluctuateCharge) {
      fluctuateEloss(particleMass, hit->pabs(), eLoss, direction.mag(), NumberOfSegmentation, eLossVector, engine);
      // Save the energy of each segment
      for (int i = 0; i != NumberOfSegmentation; i++) {
        // take energy value from vector eLossVector,
        _ionization_points[i] =
            EnergyDepositUnit(eLossVector[i] * decSignal / eLoss,
                              hit->entryPoint() + float((i + 0.5) / NumberOfSegmentation) * direction);
      }
    } else {
      // Save the energy of each segment
      for (int i = 0; i != NumberOfSegmentation; i++) {
        // take energy value from eLoss average over n.segments.
        _ionization_points[i] =
            EnergyDepositUnit(decSignal / float(NumberOfSegmentation),
                              hit->entryPoint() + float((i + 0.5) / NumberOfSegmentation) * direction);
      }
    }
  }
  return _ionization_points;
}

driftXPos()

float SiLinearChargeDivider::driftXPos ( const Local3DPoint &  pos, const LocalVector &  drift, double  thickness  )

inlineprivate

Definition at line 57 of file SiLinearChargeDivider.h.

References shallow::drift(), Calorimetry_cff::thickness, PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by divide().

                                                                                              {
    return pos.x() + (thickness / 2. - pos.z()) * drift.x() / drift.z();
  }

fluctuateEloss()

void SiLinearChargeDivider::fluctuateEloss ( double const  particleMass, float  momentum, float  eloss, float  length, int  NumberOfSegmentation, float  elossVector[], CLHEP::HepRandomEngine *  engine  )

private

Definition at line 142 of file SiLinearChargeDivider.cc.

References deltaCut, fluctuate, mps_fire::i, cuy::ii, and particleFlowDisplacedVertex_cfi::ratio.

Referenced by divide().

                                                                         {
  // Generate charge fluctuations.
  float sum = 0.;
  double deltaCutoff;
  double mom = particleMomentum * 1000.;
  double seglen = length / NumberOfSegs * 10.;
  double segeloss = (1000. * eloss) / NumberOfSegs;
  for (int i = 0; i < NumberOfSegs; i++) {
    // The G4 routine needs momentum in MeV, mass in MeV, delta-cut in MeV,
    // track segment length in mm, segment eloss in MeV
    // Returns fluctuated eloss in MeV
    // the cutoff is sometimes redefined inside, so fix it.
    deltaCutoff = deltaCut;
    sum += (elossVector[i] =
                fluctuate->SampleFluctuations(mom, particleMass, deltaCutoff, seglen, segeloss, engine) / 1000.);
  }

  if (sum > 0.) {  // If fluctuations give eloss>0.
    // Rescale to the same total eloss
    float ratio = eloss / sum;
    for (int ii = 0; ii < NumberOfSegs; ii++)
      elossVector[ii] = ratio * elossVector[ii];
  } else {  // If fluctuations gives 0 eloss
    float averageEloss = eloss / NumberOfSegs;
    for (int ii = 0; ii < NumberOfSegs; ii++)
      elossVector[ii] = averageEloss;
  }
  return;
}

readPulseShape()

void SiLinearChargeDivider::readPulseShape ( const std::string &  pulseShapeFileName )

private

Definition at line 33 of file SiLinearChargeDivider.cc.

References funct::abs(), HLT_2024v13_cff::distance, geometryDiff::epsilon, Exception, mps_splice::line, pulseResolution, pulset0Idx, pulseValues, AlCaHLTBitMon_QueryRunRegistry::string, and relativeConstraints::value.

Referenced by SiLinearChargeDivider().

                                                                              {
  // Pulse shape file format: empty lines and comments (lines starting with '#') are ignored
  // one line "resolution: value" is interpreted as the resolution
  // all other lines are read as consecutive values of the shape
  std::ifstream shapeFile(pulseShapeFileName.c_str());
  if (!shapeFile.good()) {
    throw cms::Exception("FileError") << "Problem opening APV Shape file: " << pulseShapeFileName;
  }
  pulseResolution = -1.;
  std::string line;
  const std::string resoPrefix{"resolution: "};
  while (std::getline(shapeFile, line)) {
    if ((!line.empty()) && (line.substr(1) != "#")) {
      std::istringstream lStr{line};
      if (line.substr(0, resoPrefix.size()) == resoPrefix) {
        lStr.seekg(resoPrefix.size());
        lStr >> pulseResolution;
      } else {
        double value;
        while (lStr >> value) {
          pulseValues.push_back(value);
        }
      }
    }
  }
  if (pulseValues.empty() || (pulseResolution == -1.)) {
    throw cms::Exception("WrongAPVPulseShape") << "Problem reading from APV pulse shape file " << pulseShapeFileName
                                               << ": " << (pulseValues.empty() ? "no values" : "no resolution");
  }
  const auto maxIt = std::max_element(pulseValues.begin(), pulseValues.end());
  if (std::abs((*maxIt) - 1.) > std::numeric_limits<double>::epsilon()) {
    throw cms::Exception("WrongAPVPulseShape")
        << "The max value of the APV pulse shape stored in the text file used in "
           "SimGeneral/MixingModule/python/SiStripSimParameters_cfi.py is not equal to 1. Need to be fixed.";
  }
  pulset0Idx = std::distance(pulseValues.begin(), maxIt);
}

setParticleDataTable()

void SiLinearChargeDivider::setParticleDataTable ( const ParticleDataTable *  pdt )

inlineoverridevirtual

Implements SiChargeDivider.

Definition at line 40 of file SiLinearChargeDivider.h.

References theParticleDataTable.

{ theParticleDataTable = pdt; }

TimeResponse()

float SiLinearChargeDivider::TimeResponse ( const PSimHit *  hit, const StripGeomDetUnit &  det  )

private

Definition at line 178 of file SiLinearChargeDivider.cc.

References cosmicShift, createfilelist::int, mag(), pulseResolution, pulset0Idx, pulseValues, GeomDet::surface(), Surface::toGlobal(), and x.

Referenced by divide().

                                                                                         {
  // x is difference between the tof and the tof for a photon (reference)
  // converted into a bin number
  const auto dTOF = det.surface().toGlobal(hit->localPosition()).mag() / 30. + cosmicShift - hit->tof();
  const int x = int(dTOF / pulseResolution) + pulset0Idx;
  if (x < 0 || x >= int(pulseValues.size()))
    return 0;
  return hit->energyLoss() * pulseValues[std::size_t(x)];
}

Member Data Documentation

chargedivisionsPerStrip

const int SiLinearChargeDivider::chargedivisionsPerStrip

private

Definition at line 46 of file SiLinearChargeDivider.h.

Referenced by divide().

cosmicShift

const double SiLinearChargeDivider::cosmicShift

private

Definition at line 48 of file SiLinearChargeDivider.h.

Referenced by TimeResponse().

deltaCut

const double SiLinearChargeDivider::deltaCut

private

Definition at line 47 of file SiLinearChargeDivider.h.

Referenced by fluctuateEloss().

fluctuate

std::unique_ptr<SiG4UniversalFluctuation> SiLinearChargeDivider::fluctuate

private

Definition at line 55 of file SiLinearChargeDivider.h.

Referenced by fluctuateEloss().

fluctuateCharge

const bool SiLinearChargeDivider::fluctuateCharge

private

Definition at line 45 of file SiLinearChargeDivider.h.

Referenced by divide().

peakMode

const bool SiLinearChargeDivider::peakMode

private

Definition at line 44 of file SiLinearChargeDivider.h.

Referenced by SiLinearChargeDivider().

pulseResolution

double SiLinearChargeDivider::pulseResolution

private

Definition at line 50 of file SiLinearChargeDivider.h.

Referenced by readPulseShape(), and TimeResponse().

pulset0Idx

unsigned int SiLinearChargeDivider::pulset0Idx

private

Definition at line 51 of file SiLinearChargeDivider.h.

Referenced by readPulseShape(), and TimeResponse().

pulseValues

std::vector<double> SiLinearChargeDivider::pulseValues

private

Definition at line 52 of file SiLinearChargeDivider.h.

Referenced by readPulseShape(), and TimeResponse().

theParticleDataTable

const ParticleDataTable* SiLinearChargeDivider::theParticleDataTable

private

Definition at line 49 of file SiLinearChargeDivider.h.

Referenced by divide(), and setParticleDataTable().