CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes | Static Private Attributes
DreamSD Class Reference

#include <DreamSD.h>

Inheritance diagram for DreamSD:
CaloSD SensitiveCaloDetector Observer< const BeginOfRun * > Observer< const BeginOfEvent * > Observer< const BeginOfTrack * > Observer< const EndOfTrack * > Observer< const EndOfEvent * > SensitiveDetector

Public Member Functions

 DreamSD (const std::string &, const edm::EventSetup &, const SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
 
uint32_t setDetUnitId (const G4Step *) override
 
 ~DreamSD () override
 
- Public Member Functions inherited from CaloSD
 CaloSD (const std::string &aSDname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p, const SimTrackManager *, float timeSlice=1., bool ignoreTkID=false)
 
void clear () override
 
void clearHits () override
 
void DrawAll () override
 
void EndOfEvent (G4HCofThisEvent *eventHC) override
 
void fillHits (edm::PCaloHitContainer &, const std::string &) override
 
void Initialize (G4HCofThisEvent *HCE) override
 
bool isItFineCalo (const G4VTouchable *touch)
 
void PrintAll () override
 
bool ProcessHits (G4GFlashSpot *aSpot, G4TouchableHistory *) override
 
G4bool ProcessHits (G4Step *step, G4TouchableHistory *) override
 
void reset () override
 
 ~CaloSD () override
 
- Public Member Functions inherited from SensitiveCaloDetector
 SensitiveCaloDetector (const std::string &iname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p)
 
- Public Member Functions inherited from SensitiveDetector
void EndOfEvent (G4HCofThisEvent *eventHC) override
 
const std::vector< std::string > & getNames () const
 
void Initialize (G4HCofThisEvent *eventHC) override
 
bool isCaloSD () const
 
 SensitiveDetector (const std::string &iname, const edm::EventSetup &es, const SensitiveDetectorCatalog &, edm::ParameterSet const &p, bool calo)
 
 ~SensitiveDetector () override
 
- Public Member Functions inherited from Observer< const BeginOfRun * >
 Observer ()
 
void slotForUpdate (const BeginOfRun * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const BeginOfEvent * >
 Observer ()
 
void slotForUpdate (const BeginOfEvent * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const BeginOfTrack * >
 Observer ()
 
void slotForUpdate (const BeginOfTrack * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const EndOfTrack * >
 Observer ()
 
void slotForUpdate (const EndOfTrack * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const EndOfEvent * >
 Observer ()
 
void slotForUpdate (const EndOfEvent * iT)
 
virtual ~Observer ()
 

Protected Member Functions

double getEnergyDeposit (const G4Step *) override
 
void initRun () override
 
- Protected Member Functions inherited from CaloSD
bool checkHit ()
 
CaloG4HitcreateNewHit (const G4Step *, const G4Track *)
 
virtual void endEvent ()
 
virtual double EnergyCorrected (const G4Step &step, const G4Track *)
 
virtual bool filterHit (CaloG4Hit *, double)
 
double getAttenuation (const G4Step *aStep, double birk1, double birk2, double birk3) const
 
virtual uint16_t getDepth (const G4Step *)
 
virtual bool getFromLibrary (const G4Step *step)
 
int getNumberOfHits ()
 
double getResponseWt (const G4Track *)
 
virtual int getTrackID (const G4Track *)
 
bool hitExists (const G4Step *)
 
void ignoreRejection ()
 
virtual void initEvent (const BeginOfEvent *)
 
void processHit (const G4Step *step)
 
void resetForNewPrimary (const G4Step *)
 
void setNumberCheckedHits (int val)
 
void setParameterized (bool val)
 
G4ThreeVector setToGlobal (const G4ThreeVector &, const G4VTouchable *) const
 
G4ThreeVector setToLocal (const G4ThreeVector &, const G4VTouchable *) const
 
virtual int setTrackID (const G4Step *)
 
void setUseMap (bool val)
 
void update (const ::EndOfEvent *) override
 
void update (const BeginOfEvent *) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const BeginOfRun *) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const BeginOfTrack *trk) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const EndOfTrack *trk) override
 This routine will be called when the appropriate signal arrives. More...
 
void updateHit (CaloG4Hit *)
 
- Protected Member Functions inherited from SensitiveDetector
TrackInformationcmsTrackInformation (const G4Track *aTrack)
 
Local3DPoint ConvertToLocal3DPoint (const G4ThreeVector &point) const
 
Local3DPoint FinalStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint InitialStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint LocalPostStepPosition (const G4Step *step) const
 
Local3DPoint LocalPreStepPosition (const G4Step *step) const
 
void NaNTrap (const G4Step *step) const
 
void setNames (const std::vector< std::string > &)
 
- Protected Member Functions inherited from Observer< const EndOfEvent * >
virtual void update (const EndOfEvent *)=0
 This routine will be called when the appropriate signal arrives. More...
 

Private Types

typedef std::map< G4LogicalVolume *, DoublesDimensionMap
 
typedef std::pair< double, double > Doubles
 

Private Member Functions

double cherenkovDeposit_ (const G4Step *aStep)
 Returns the total energy due to Cherenkov radiation. More...
 
double crystalLength (G4LogicalVolume *) const
 
double crystalWidth (G4LogicalVolume *) const
 
double curve_LY (const G4Step *, int)
 
void fillMap (const std::string &, double, double)
 
double getAverageNumberOfPhotons_ (const double charge, const double beta, const G4Material *aMaterial, const G4MaterialPropertyVector *rIndex)
 Returns average number of photons created by track. More...
 
double getPhotonEnergyDeposit_ (const G4ParticleMomentum &p, const G4ThreeVector &x, const G4Step *aStep)
 Returns energy deposit for a given photon. More...
 
void initMap (const std::string &, const edm::EventSetup &)
 
bool setPbWO2MaterialProperties_ (G4Material *aMaterial)
 Sets material properties at run-time... More...
 

Private Attributes

double birk1_
 
double birk2_
 
double birk3_
 
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
 Table of Cherenkov angle integrals vs photon momentum. More...
 
bool dd4hep_
 
bool doCherenkov_
 
G4MaterialPropertiesTable * materialPropertiesTable_
 
int nphotons_
 
bool readBothSide_
 
int side_
 
double slopeLY_
 
bool useBirk_
 
DimensionMap xtalLMap_
 

Static Private Attributes

static constexpr double k_ScaleFromDD4HepToG4 = 10.0
 
static constexpr double k_ScaleFromDDDToG4 = 1.0
 

Additional Inherited Members

- Protected Types inherited from SensitiveDetector
enum  coordinates { WorldCoordinates, LocalCoordinates }
 
- Protected Attributes inherited from CaloSD
CaloG4HitcurrentHit
 
CaloHitID currentID
 
float edepositEM
 
float edepositHAD
 
double eminHit
 
double energyCut
 
G4ThreeVector entranceLocal
 
G4ThreeVector entrancePoint
 
bool forceSave
 
float incidentEnergy
 
double kmaxIon
 
double kmaxNeutron
 
double kmaxProton
 
G4ThreeVector posGlobal
 
CaloHitID previousID
 
bool suppressHeavy
 
double tmaxHit
 

Detailed Description

Definition at line 15 of file DreamSD.h.

Member Typedef Documentation

◆ DimensionMap

typedef std::map<G4LogicalVolume *, Doubles> DreamSD::DimensionMap
private

Definition at line 32 of file DreamSD.h.

◆ Doubles

typedef std::pair<double, double> DreamSD::Doubles
private

Definition at line 31 of file DreamSD.h.

Constructor & Destructor Documentation

◆ DreamSD()

DreamSD::DreamSD ( const std::string &  name,
const edm::EventSetup es,
const SensitiveDetectorCatalog clg,
edm::ParameterSet const &  p,
const SimTrackManager manager 
)

Definition at line 37 of file DreamSD.cc.

42  : CaloSD(name, es, clg, p, manager) {
43  edm::ParameterSet m_EC = p.getParameter<edm::ParameterSet>("ECalSD");
44  useBirk_ = m_EC.getParameter<bool>("UseBirkLaw");
45  doCherenkov_ = m_EC.getParameter<bool>("doCherenkov");
46  birk1_ = m_EC.getParameter<double>("BirkC1") * (g / (MeV * cm2));
47  birk2_ = m_EC.getParameter<double>("BirkC2");
48  birk3_ = m_EC.getParameter<double>("BirkC3");
49  slopeLY_ = m_EC.getParameter<double>("SlopeLightYield");
50  readBothSide_ = m_EC.getUntrackedParameter<bool>("ReadBothSide", false);
51  dd4hep_ = m_EC.getUntrackedParameter<bool>("DD4Hep", false);
52 
53  chAngleIntegrals_.reset(nullptr);
54 
55  edm::LogVerbatim("EcalSim") << "Constructing a DreamSD with name " << GetName()
56  << "\nDreamSD:: Use of Birks law is set to " << useBirk_
57  << " with three constants kB = " << birk1_ << ", C1 = " << birk2_ << ", C2 = " << birk3_
58  << "\n Slope for Light yield is set to " << slopeLY_
59  << "\n Parameterization of Cherenkov is set to " << doCherenkov_
60  << ", readout both sides is " << readBothSide_ << " and dd4hep flag " << dd4hep_;
61 
62  initMap(name, es);
63 }

References birk1_, birk2_, birk3_, chAngleIntegrals_, dd4hep_, doCherenkov_, g, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), initMap(), MeV, Skims_PA_cff::name, AlCaHLTBitMon_ParallelJobs::p, readBothSide_, slopeLY_, and useBirk_.

◆ ~DreamSD()

DreamSD::~DreamSD ( )
inlineoverride

Definition at line 22 of file DreamSD.h.

22 {}

Member Function Documentation

◆ cherenkovDeposit_()

double DreamSD::cherenkovDeposit_ ( const G4Step *  aStep)
private

Returns the total energy due to Cherenkov radiation.

Definition at line 241 of file DreamSD.cc.

241  {
242  double cherenkovEnergy = 0;
244  return cherenkovEnergy;
245  G4Material *material = aStep->GetTrack()->GetMaterial();
246 
247  // Retrieve refractive index
248  G4MaterialPropertyVector *Rindex = materialPropertiesTable_->GetProperty("RINDEX");
249  if (Rindex == nullptr) {
250  edm::LogWarning("EcalSim") << "Couldn't retrieve refractive index";
251  return cherenkovEnergy;
252  }
253 
254  // V.Ivanchenko - temporary close log output for 9.5
255  // Material refraction properties
256  int Rlength = Rindex->GetVectorLength() - 1;
257  double Pmin = Rindex->Energy(0);
258  double Pmax = Rindex->Energy(Rlength);
259 #ifdef EDM_ML_DEBUG
260  edm::LogVerbatim("EcalSim") << "Material properties: \n Pmin = " << Pmin << " Pmax = " << Pmax;
261 #endif
262  // Get particle properties
263  const G4StepPoint *pPreStepPoint = aStep->GetPreStepPoint();
264  const G4StepPoint *pPostStepPoint = aStep->GetPostStepPoint();
265  const G4ThreeVector &x0 = pPreStepPoint->GetPosition();
266  G4ThreeVector p0 = aStep->GetDeltaPosition().unit();
267  const G4DynamicParticle *aParticle = aStep->GetTrack()->GetDynamicParticle();
268  const double charge = aParticle->GetDefinition()->GetPDGCharge();
269  // beta is averaged over step
270  double beta = 0.5 * (pPreStepPoint->GetBeta() + pPostStepPoint->GetBeta());
271  double BetaInverse = 1.0 / beta;
272 
273 #ifdef EDM_ML_DEBUG
274  edm::LogVerbatim("EcalSim") << "Particle properties: \n charge = " << charge << " beta = " << beta;
275 #endif
276 
277  // Now get number of photons generated in this step
278  double meanNumberOfPhotons = getAverageNumberOfPhotons_(charge, beta, material, Rindex);
279  if (meanNumberOfPhotons <= 0.0) { // Don't do anything
280 #ifdef EDM_ML_DEBUG
281  edm::LogVerbatim("EcalSim") << "Mean number of photons is zero: " << meanNumberOfPhotons << ", stopping here";
282 #endif
283  return cherenkovEnergy;
284  }
285 
286  // number of photons is in unit of Geant4...
287  meanNumberOfPhotons *= aStep->GetStepLength();
288 
289  // Now get a poisson distribution
290  int numPhotons = static_cast<int>(G4Poisson(meanNumberOfPhotons));
291  // edm::LogVerbatim("EcalSim") << "Number of photons = " << numPhotons;
292  if (numPhotons <= 0) {
293 #ifdef EDM_ML_DEBUG
294  edm::LogVerbatim("EcalSim") << "Poission number of photons is zero: " << numPhotons << ", stopping here";
295 #endif
296  return cherenkovEnergy;
297  }
298 
299  // Material refraction properties
300  double dp = Pmax - Pmin;
301  double maxCos = BetaInverse / (*Rindex)[Rlength];
302  double maxSin2 = (1.0 - maxCos) * (1.0 + maxCos);
303 
304  // Finally: get contribution of each photon
305  for (int iPhoton = 0; iPhoton < numPhotons; ++iPhoton) {
306  // Determine photon momentum
307  double randomNumber;
308  double sampledMomentum, sampledRI;
309  double cosTheta, sin2Theta;
310 
311  // sample a momentum (not sure why this is needed!)
312  do {
313  randomNumber = G4UniformRand();
314  sampledMomentum = Pmin + randomNumber * dp;
315  sampledRI = Rindex->Value(sampledMomentum);
316  cosTheta = BetaInverse / sampledRI;
317 
318  sin2Theta = (1.0 - cosTheta) * (1.0 + cosTheta);
319  randomNumber = G4UniformRand();
320 
321  } while (randomNumber * maxSin2 > sin2Theta);
322 
323  // Generate random position of photon on cone surface
324  // defined by Theta
325  randomNumber = G4UniformRand();
326 
327  double phi = twopi * randomNumber;
328  double sinPhi = sin(phi);
329  double cosPhi = cos(phi);
330 
331  // Create photon momentum direction vector
332  // The momentum direction is still w.r.t. the coordinate system where the
333  // primary particle direction is aligned with the z axis
334  double sinTheta = sqrt(sin2Theta);
335  double px = sinTheta * cosPhi;
336  double py = sinTheta * sinPhi;
337  double pz = cosTheta;
338  G4ThreeVector photonDirection(px, py, pz);
339 
340  // Rotate momentum direction back to global (crystal) reference system
341  photonDirection.rotateUz(p0);
342 
343  // Create photon position and momentum
344  randomNumber = G4UniformRand();
345  G4ThreeVector photonPosition = x0 + randomNumber * aStep->GetDeltaPosition();
346  G4ThreeVector photonMomentum = sampledMomentum * photonDirection;
347 
348  // Collect energy on APD
349  cherenkovEnergy += getPhotonEnergyDeposit_(photonMomentum, photonPosition, aStep);
350  }
351  return cherenkovEnergy;
352 }

References zMuMuMuonUserData::beta, ALCARECOTkAlJpsiMuMu_cff::charge, funct::cos(), Calorimetry_cff::dp, getAverageNumberOfPhotons_(), getPhotonEnergyDeposit_(), materialPropertiesTable_, phi, dataAnalyzerFineBiningParameters_cff::Pmax, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, funct::sin(), and mathSSE::sqrt().

Referenced by getEnergyDeposit().

◆ crystalLength()

double DreamSD::crystalLength ( G4LogicalVolume *  lv) const
private

Definition at line 221 of file DreamSD.cc.

221  {
222  double length = -1.;
223  DimensionMap::const_iterator ite = xtalLMap_.find(lv);
224  if (ite != xtalLMap_.end())
225  length = ite->second.first;
226  return length;
227 }

References xtalLMap_.

Referenced by curve_LY(), and getPhotonEnergyDeposit_().

◆ crystalWidth()

double DreamSD::crystalWidth ( G4LogicalVolume *  lv) const
private

Definition at line 230 of file DreamSD.cc.

230  {
231  double width = -1.;
232  DimensionMap::const_iterator ite = xtalLMap_.find(lv);
233  if (ite != xtalLMap_.end())
234  width = ite->second.second;
235  return width;
236 }

References ApeEstimator_cff::width, and xtalLMap_.

Referenced by getPhotonEnergyDeposit_().

◆ curve_LY()

double DreamSD::curve_LY ( const G4Step *  aStep,
int  flag 
)
private

Definition at line 194 of file DreamSD.cc.

194  {
195  auto const stepPoint = aStep->GetPreStepPoint();
196  auto const lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
197  G4String nameVolume = lv->GetName();
198 
199  double weight = 1.;
200  G4ThreeVector localPoint = setToLocal(stepPoint->GetPosition(), stepPoint->GetTouchable());
201  double crlength = crystalLength(lv);
202  double localz = localPoint.x();
203  double dapd = 0.5 * crlength - flag * localz; // Distance from closest APD
204  if (dapd >= -0.1 || dapd <= crlength + 0.1) {
205  if (dapd <= 100.)
206  weight = 1.0 + slopeLY_ - dapd * 0.01 * slopeLY_;
207  } else {
208  edm::LogWarning("EcalSim") << "DreamSD: light coll curve : wrong distance "
209  << "to APD " << dapd << " crlength = " << crlength << " crystal name = " << nameVolume
210  << " z of localPoint = " << localz << " take weight = " << weight;
211  }
212 #ifdef EDM_ML_DEBUG
213  edm::LogVerbatim("EcalSim") << "DreamSD, light coll curve : " << dapd << " crlength = " << crlength
214  << " crystal name = " << nameVolume << " z of localPoint = " << localz
215  << " take weight = " << weight;
216 #endif
217  return weight;
218 }

References crystalLength(), RemoveAddSevLevel::flag, CaloSD::setToLocal(), slopeLY_, and mps_merge::weight.

Referenced by getEnergyDeposit().

◆ fillMap()

void DreamSD::fillMap ( const std::string &  name,
double  length,
double  width 
)
private

Definition at line 176 of file DreamSD.cc.

176  {
177  const G4LogicalVolumeStore *lvs = G4LogicalVolumeStore::GetInstance();
178  std::vector<G4LogicalVolume *>::const_iterator lvcite;
179  G4LogicalVolume *lv = nullptr;
180  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++) {
181  if ((*lvcite)->GetName() == static_cast<G4String>(name)) {
182  lv = (*lvcite);
183  break;
184  }
185  }
186 #ifdef EDM_ML_DEBUG
187  edm::LogVerbatim("EcalSim") << "DreamSD::fillMap (for " << name << " Logical Volume " << lv << " Length " << length
188  << " Width " << width;
189 #endif
190  xtalLMap_.insert(std::pair<G4LogicalVolume *, Doubles>(lv, Doubles(length, width)));
191 }

References Skims_PA_cff::name, ApeEstimator_cff::width, and xtalLMap_.

Referenced by initMap().

◆ getAverageNumberOfPhotons_()

double DreamSD::getAverageNumberOfPhotons_ ( const double  charge,
const double  beta,
const G4Material *  aMaterial,
const G4MaterialPropertyVector *  rIndex 
)
private

Returns average number of photons created by track.

Definition at line 357 of file DreamSD.cc.

360  {
361  const double rFact = 369.81 / (eV * cm);
362 
363  if (beta <= 0.0)
364  return 0.0;
365 
366  double BetaInverse = 1. / beta;
367 
368  // Vectors used in computation of Cerenkov Angle Integral:
369  // - Refraction Indices for the current material
370  // - new G4PhysicsOrderedFreeVector allocated to hold CAI's
371 
372  // Min and Max photon momenta
373  int Rlength = Rindex->GetVectorLength() - 1;
374  double Pmin = Rindex->Energy(0);
375  double Pmax = Rindex->Energy(Rlength);
376 
377  // Min and Max Refraction Indices
378  double nMin = (*Rindex)[0];
379  double nMax = (*Rindex)[Rlength];
380 
381  // Max Cerenkov Angle Integral
382  double CAImax = chAngleIntegrals_.get()->GetMaxValue();
383 
384  double dp = 0., ge = 0., CAImin = 0.;
385 
386  // If n(Pmax) < 1/Beta -- no photons generated
387  if (nMax < BetaInverse) {
388  }
389 
390  // otherwise if n(Pmin) >= 1/Beta -- photons generated
391  else if (nMin > BetaInverse) {
392  dp = Pmax - Pmin;
393  ge = CAImax;
394  }
395  // If n(Pmin) < 1/Beta, and n(Pmax) >= 1/Beta, then
396  // we need to find a P such that the value of n(P) == 1/Beta.
397  // Interpolation is performed by the GetPhotonEnergy() and
398  // GetProperty() methods of the G4MaterialPropertiesTable and
399  // the GetValue() method of G4PhysicsVector.
400  else {
401  Pmin = Rindex->Value(BetaInverse);
402  dp = Pmax - Pmin;
403  // need boolean for current implementation of G4PhysicsVector
404  // ==> being phased out
405  double CAImin = chAngleIntegrals_->Value(Pmin);
406  ge = CAImax - CAImin;
407  }
408 
409  // Calculate number of photons
410  double numPhotons = rFact * charge / eplus * charge / eplus * (dp - ge * BetaInverse * BetaInverse);
411 
412 #ifdef EDM_ML_DEBUG
413  edm::LogVerbatim("EcalSim") << "@SUB=getAverageNumberOfPhotons\nCAImin = " << CAImin << "\nCAImax = " << CAImax
414  << "\ndp = " << dp << ", ge = " << ge << "\nnumPhotons = " << numPhotons;
415 #endif
416  return numPhotons;
417 }

References zMuMuMuonUserData::beta, chAngleIntegrals_, ALCARECOTkAlJpsiMuMu_cff::charge, Calorimetry_cff::dp, metDiagnosticParameterSet_cfi::nMax, metDiagnosticParameterSet_cfi::nMin, and dataAnalyzerFineBiningParameters_cff::Pmax.

Referenced by cherenkovDeposit_().

◆ getEnergyDeposit()

double DreamSD::getEnergyDeposit ( const G4Step *  aStep)
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 66 of file DreamSD.cc.

66  {
67  // take into account light collection curve for crystals
68  double weight = curve_LY(aStep, side_);
69  if (useBirk_)
71  double edep = aStep->GetTotalEnergyDeposit() * weight;
72 
73  // Get Cerenkov contribution
74  if (doCherenkov_) {
75  edep += cherenkovDeposit_(aStep);
76  }
77 #ifdef EDM_ML_DEBUG
78  edm::LogVerbatim("EcalSim") << "DreamSD:: " << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() << " Side "
79  << side_ << " Light Collection Efficiency " << weight << " Weighted Energy Deposit "
80  << edep / CLHEP::MeV << " MeV";
81 #endif
82  return edep;
83 }

References birk1_, birk2_, birk3_, cherenkovDeposit_(), curve_LY(), doCherenkov_, CaloSD::getAttenuation(), MeV, side_, useBirk_, and mps_merge::weight.

◆ getPhotonEnergyDeposit_()

double DreamSD::getPhotonEnergyDeposit_ ( const G4ParticleMomentum &  p,
const G4ThreeVector &  x,
const G4Step *  aStep 
)
private

Returns energy deposit for a given photon.

Definition at line 503 of file DreamSD.cc.

503  {
504  double energy = 0;
505 
506  // Crystal dimensions
507 
508  // edm::LogVerbatim("EcalSim") << p << x;
509 
510  // 1. Check if this photon goes straight to the APD:
511  // - assume that APD is at x=xtalLength/2.0
512  // - extrapolate from x=x0 to x=xtalLength/2.0 using momentum in x-y
513 
514  G4StepPoint *stepPoint = aStep->GetPreStepPoint();
515  G4LogicalVolume *lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
516  G4String nameVolume = lv->GetName();
517 
518  double crlength = crystalLength(lv);
519  double crwidth = crystalWidth(lv);
520  double dapd = 0.5 * crlength - x.x(); // Distance from closest APD
521  double y = p.y() / p.x() * dapd;
522 
523 #ifdef EDM_ML_DEBUG
524  edm::LogVerbatim("EcalSim") << "Distance to APD: " << dapd << " - y at APD: " << y;
525 #endif
526  // Not straight: compute probability
527  if (std::abs(y) > crwidth * 0.5) {
528  }
529 
530  // 2. Retrieve efficiency for this wavelength (in nm, from MeV)
531  double waveLength = p.mag() * 1.239e8;
532 
533  energy = p.mag() * PMTResponse::getEfficiency(waveLength);
534 #ifdef EDM_ML_DEBUG
535  edm::LogVerbatim("EcalSim") << "Wavelength: " << waveLength << " - Energy: " << energy;
536 #endif
537  return energy;
538 }

References funct::abs(), crystalLength(), crystalWidth(), HCALHighEnergyHPDFilter_cfi::energy, PMTResponse::getEfficiency(), AlCaHLTBitMon_ParallelJobs::p, x, and y.

Referenced by cherenkovDeposit_().

◆ initMap()

void DreamSD::initMap ( const std::string &  sd,
const edm::EventSetup es 
)
private

Definition at line 117 of file DreamSD.cc.

117  {
118  if (dd4hep_) {
120  es.get<IdealGeometryRecord>().get(cpv);
121  const cms::DDFilter filter("ReadOutName", sd);
122  cms::DDFilteredView fv((*cpv), filter);
123  while (fv.firstChild()) {
124  std::string name = static_cast<std::string>(dd4hep::dd::noNamespace(fv.name()));
125  std::vector<double> paras(fv.parameters());
126 #ifdef EDM_ML_DEBUG
127  edm::LogVerbatim("EcalSim") << "DreamSD::initMap (for " << sd << "): Solid " << name << " Shape "
128  << cms::dd::name(cms::DDSolidShapeMap, fv.shape()) << " Parameter 0 = " << paras[0];
129 #endif
130  // Set length to be the largest size, width the smallest
131  std::sort(paras.begin(), paras.end());
132  double length = 2.0 * k_ScaleFromDD4HepToG4 * paras.back();
133  double width = 2.0 * k_ScaleFromDD4HepToG4 * paras.front();
134  fillMap(name, length, width);
135  }
136  } else {
138  es.get<IdealGeometryRecord>().get(cpv);
139  DDSpecificsMatchesValueFilter filter{DDValue("ReadOutName", sd, 0)};
140  DDFilteredView fv((*cpv), filter);
141 
142  bool dodet = fv.firstChild();
143  while (dodet) {
144  const DDSolid &sol = fv.logicalPart().solid();
145  std::vector<double> paras(sol.parameters());
146  std::string name = sol.name().name();
147 #ifdef EDM_ML_DEBUG
148  edm::LogVerbatim("EcalSim") << "DreamSD::initMap (for " << sd << "): Solid " << name << " Shape " << sol.shape()
149  << " Parameter 0 = " << paras[0];
150 #endif
151  // Set length to be the largest size, width the smallest
152  std::sort(paras.begin(), paras.end());
153  double length = 2.0 * k_ScaleFromDDDToG4 * paras.back();
154  double width = 2.0 * k_ScaleFromDDDToG4 * paras.front();
155  fillMap(name, length, width);
156  dodet = fv.next();
157  }
158  }
159 #ifdef EDM_ML_DEBUG
160  edm::LogVerbatim("EcalSim") << "DreamSD: Length Table for ReadOutName = " << sd << ":";
161 #endif
162  DimensionMap::const_iterator ite = xtalLMap_.begin();
163  int i = 0;
164  for (; ite != xtalLMap_.end(); ite++, i++) {
165  G4String name = "Unknown";
166  if (ite->first != nullptr)
167  name = (ite->first)->GetName();
168 #ifdef EDM_ML_DEBUG
169  edm::LogVerbatim("EcalSim") << " " << i << " " << ite->first << " " << name << " L = " << ite->second.first
170  << " W = " << ite->second.second;
171 #endif
172  }
173 }

References dd4hep_, cms::DDSolidShapeMap, fillMap(), ALCARECOTkAlBeamHalo_cff::filter, DDFilteredView::firstChild(), cms::DDFilteredView::firstChild(), edm::EventSetup::get(), get, mps_fire::i, k_ScaleFromDD4HepToG4, k_ScaleFromDDDToG4, DDFilteredView::logicalPart(), Skims_PA_cff::name, cms::dd::name(), DDName::name(), DDBase< N, C >::name(), cms::DDFilteredView::name(), DDFilteredView::next(), DDSolid::parameters(), cms::DDFilteredView::parameters(), sd, DDSolid::shape(), cms::DDFilteredView::shape(), DDLogicalPart::solid(), AlCaHLTBitMon_QueryRunRegistry::string, ApeEstimator_cff::width, and xtalLMap_.

Referenced by DreamSD().

◆ initRun()

void DreamSD::initRun ( )
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 86 of file DreamSD.cc.

86  {
87  // Get the material and set properties if needed
88  DimensionMap::const_iterator ite = xtalLMap_.begin();
89  const G4LogicalVolume *lv = (ite->first);
90  G4Material *material = lv->GetMaterial();
91  edm::LogVerbatim("EcalSim") << "DreamSD::initRun: Initializes for material " << material->GetName() << " in "
92  << lv->GetName();
93  materialPropertiesTable_ = material->GetMaterialPropertiesTable();
95  if (!setPbWO2MaterialProperties_(material)) {
96  edm::LogWarning("EcalSim") << "Couldn't retrieve material properties table\n Material = " << material->GetName();
97  }
98  materialPropertiesTable_ = material->GetMaterialPropertiesTable();
99  }
100 }

References materialPropertiesTable_, setPbWO2MaterialProperties_(), and xtalLMap_.

◆ setDetUnitId()

uint32_t DreamSD::setDetUnitId ( const G4Step *  aStep)
overridevirtual

Implements CaloSD.

Definition at line 103 of file DreamSD.cc.

103  {
104  const G4VTouchable *touch = aStep->GetPreStepPoint()->GetTouchable();
105  uint32_t id = (touch->GetReplicaNumber(1)) * 10 + (touch->GetReplicaNumber(0));
106  side_ = readBothSide_ ? -1 : 1;
107  if (side_ < 0) {
108  ++id;
109  }
110 #ifdef EDM_ML_DEBUG
111  edm::LogVerbatim("EcalSim") << "DreamSD:: ID " << id;
112 #endif
113  return id;
114 }

References triggerObjects_cff::id, readBothSide_, and side_.

◆ setPbWO2MaterialProperties_()

bool DreamSD::setPbWO2MaterialProperties_ ( G4Material *  aMaterial)
private

Sets material properties at run-time...

Definition at line 422 of file DreamSD.cc.

422  {
423  std::string pbWO2Name("E_PbWO4");
424  if (pbWO2Name != aMaterial->GetName()) { // Wrong material!
425  edm::LogWarning("EcalSim") << "This is not the right material: "
426  << "expecting " << pbWO2Name << ", got " << aMaterial->GetName();
427  return false;
428  }
429 
430  G4MaterialPropertiesTable *table = new G4MaterialPropertiesTable();
431 
432  // Refractive index as a function of photon momentum
433  // FIXME: Should somehow put that in the configuration
434  const int nEntries = 14;
435  double PhotonEnergy[nEntries] = {1.7712 * eV,
436  1.8368 * eV,
437  1.90745 * eV,
438  1.98375 * eV,
439  2.0664 * eV,
440  2.15625 * eV,
441  2.25426 * eV,
442  2.3616 * eV,
443  2.47968 * eV,
444  2.61019 * eV,
445  2.75521 * eV,
446  2.91728 * eV,
447  3.09961 * eV,
448  3.30625 * eV};
449  double RefractiveIndex[nEntries] = {2.17728,
450  2.18025,
451  2.18357,
452  2.18753,
453  2.19285,
454  2.19813,
455  2.20441,
456  2.21337,
457  2.22328,
458  2.23619,
459  2.25203,
460  2.27381,
461  2.30282,
462  2.34666};
463 
464  table->AddProperty("RINDEX", PhotonEnergy, RefractiveIndex, nEntries);
465  aMaterial->SetMaterialPropertiesTable(table); // FIXME: could this leak? What does G4 do?
466 
467  // Calculate Cherenkov angle integrals:
468  // This is an ad-hoc solution (we hold it in the class, not in the material)
469  chAngleIntegrals_ = std::make_unique<G4PhysicsOrderedFreeVector>();
470 
471  int index = 0;
472  double currentRI = RefractiveIndex[index];
473  double currentPM = PhotonEnergy[index];
474  double currentCAI = 0.0;
475  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
476  double prevPM = currentPM;
477  double prevCAI = currentCAI;
478  double prevRI = currentRI;
479  while (++index < nEntries) {
480  currentRI = RefractiveIndex[index];
481  currentPM = PhotonEnergy[index];
482  currentCAI = 0.5 * (1.0 / (prevRI * prevRI) + 1.0 / (currentRI * currentRI));
483  currentCAI = prevCAI + (currentPM - prevPM) * currentCAI;
484 
485  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
486 
487  prevPM = currentPM;
488  prevCAI = currentCAI;
489  prevRI = currentRI;
490  }
491 
492 #ifdef EDM_ML_DEBUG
493  edm::LogVerbatim("EcalSim") << "Material properties set for " << aMaterial->GetName();
494 #endif
495  return true;
496 }

References chAngleIntegrals_, AlCaHLTBitMon_QueryRunRegistry::string, and TableParser::table.

Referenced by initRun().

Member Data Documentation

◆ birk1_

double DreamSD::birk1_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ birk2_

double DreamSD::birk2_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ birk3_

double DreamSD::birk3_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ chAngleIntegrals_

std::unique_ptr<G4PhysicsOrderedFreeVector> DreamSD::chAngleIntegrals_
private

Table of Cherenkov angle integrals vs photon momentum.

Definition at line 63 of file DreamSD.h.

Referenced by DreamSD(), getAverageNumberOfPhotons_(), and setPbWO2MaterialProperties_().

◆ dd4hep_

bool DreamSD::dd4hep_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and initMap().

◆ doCherenkov_

bool DreamSD::doCherenkov_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ k_ScaleFromDD4HepToG4

constexpr double DreamSD::k_ScaleFromDD4HepToG4 = 10.0
staticconstexprprivate

Definition at line 53 of file DreamSD.h.

Referenced by initMap().

◆ k_ScaleFromDDDToG4

constexpr double DreamSD::k_ScaleFromDDDToG4 = 1.0
staticconstexprprivate

Definition at line 52 of file DreamSD.h.

Referenced by initMap().

◆ materialPropertiesTable_

G4MaterialPropertiesTable* DreamSD::materialPropertiesTable_
private

Definition at line 64 of file DreamSD.h.

Referenced by cherenkovDeposit_(), and initRun().

◆ nphotons_

int DreamSD::nphotons_
private

Definition at line 66 of file DreamSD.h.

◆ readBothSide_

bool DreamSD::readBothSide_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and setDetUnitId().

◆ side_

int DreamSD::side_
private

Definition at line 60 of file DreamSD.h.

Referenced by getEnergyDeposit(), and setDetUnitId().

◆ slopeLY_

double DreamSD::slopeLY_
private

Definition at line 57 of file DreamSD.h.

Referenced by curve_LY(), and DreamSD().

◆ useBirk_

bool DreamSD::useBirk_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ xtalLMap_

DimensionMap DreamSD::xtalLMap_
private

Definition at line 58 of file DreamSD.h.

Referenced by crystalLength(), crystalWidth(), fillMap(), initMap(), and initRun().

DreamSD::chAngleIntegrals_
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
Table of Cherenkov angle integrals vs photon momentum.
Definition: DreamSD.h:63
DreamSD::readBothSide_
bool readBothSide_
Definition: DreamSD.h:55
DDAxes::y
ApeEstimator_cff.width
width
Definition: ApeEstimator_cff.py:24
DreamSD::Doubles
std::pair< double, double > Doubles
Definition: DreamSD.h:31
mps_fire.i
i
Definition: mps_fire.py:428
metDiagnosticParameterSet_cfi.nMin
nMin
Definition: metDiagnosticParameterSet_cfi.py:10
DreamSD::fillMap
void fillMap(const std::string &, double, double)
Definition: DreamSD.cc:176
multPhiCorr_741_25nsDY_cfi.py
py
Definition: multPhiCorr_741_25nsDY_cfi.py:12
zMuMuMuonUserData.beta
beta
Definition: zMuMuMuonUserData.py:10
mps_merge.weight
weight
Definition: mps_merge.py:88
DreamSD::k_ScaleFromDDDToG4
static constexpr double k_ScaleFromDDDToG4
Definition: DreamSD.h:52
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
DreamSD::cherenkovDeposit_
double cherenkovDeposit_(const G4Step *aStep)
Returns the total energy due to Cherenkov radiation.
Definition: DreamSD.cc:241
DreamSD::dd4hep_
bool dd4hep_
Definition: DreamSD.h:55
CaloSD::getAttenuation
double getAttenuation(const G4Step *aStep, double birk1, double birk2, double birk3) const
Definition: CaloSD.cc:517
DreamSD::getPhotonEnergyDeposit_
double getPhotonEnergyDeposit_(const G4ParticleMomentum &p, const G4ThreeVector &x, const G4Step *aStep)
Returns energy deposit for a given photon.
Definition: DreamSD.cc:503
MeV
const double MeV
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
DDAxes::x
DreamSD::crystalWidth
double crystalWidth(G4LogicalVolume *) const
Definition: DreamSD.cc:230
DreamSD::doCherenkov_
bool doCherenkov_
Definition: DreamSD.h:55
cms::DDFilteredView
Definition: DDFilteredView.h:70
DreamSD::materialPropertiesTable_
G4MaterialPropertiesTable * materialPropertiesTable_
Definition: DreamSD.h:64
edm::LogWarning
Log< level::Warning, false > LogWarning
Definition: MessageLogger.h:122
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
DDSolid::shape
DDSolidShape shape(void) const
The type of the solid.
Definition: DDSolid.cc:119
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
edm::EventSetup::get
T get() const
Definition: EventSetup.h:80
cms::DDFilter
Definition: DDFilteredView.h:59
DreamSD::initMap
void initMap(const std::string &, const edm::EventSetup &)
Definition: DreamSD.cc:117
cms::dd::name
std::string name(Mapping a, V value)
Definition: DDSolidShapes.h:31
CaloSD::setToLocal
G4ThreeVector setToLocal(const G4ThreeVector &, const G4VTouchable *) const
Definition: CaloSD.cc:373
Calorimetry_cff.dp
dp
Definition: Calorimetry_cff.py:157
DreamSD::useBirk_
bool useBirk_
Definition: DreamSD.h:55
DreamSD::getAverageNumberOfPhotons_
double getAverageNumberOfPhotons_(const double charge, const double beta, const G4Material *aMaterial, const G4MaterialPropertyVector *rIndex)
Returns average number of photons created by track.
Definition: DreamSD.cc:357
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
DDBase::name
const N & name() const
Definition: DDBase.h:59
DreamSD::xtalLMap_
DimensionMap xtalLMap_
Definition: DreamSD.h:58
HCALHighEnergyHPDFilter_cfi.energy
energy
Definition: HCALHighEnergyHPDFilter_cfi.py:5
ALCARECOTkAlBeamHalo_cff.filter
filter
Definition: ALCARECOTkAlBeamHalo_cff.py:27
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
DreamSD::setPbWO2MaterialProperties_
bool setPbWO2MaterialProperties_(G4Material *aMaterial)
Sets material properties at run-time...
Definition: DreamSD.cc:422
ALCARECOTkAlJpsiMuMu_cff.charge
charge
Definition: ALCARECOTkAlJpsiMuMu_cff.py:47
DreamSD::birk1_
double birk1_
Definition: DreamSD.h:56
edm::ParameterSet
Definition: ParameterSet.h:47
cms::DDSolidShapeMap
const std::array< const cms::dd::NameValuePair< DDSolidShape >, 19 > DDSolidShapeMap
Definition: DDSolidShapes.h:97
PMTResponse::getEfficiency
static double getEfficiency(const double &waveLengthNm)
Return efficiency for given photon wavelength (in nm)
Definition: PMTResponse.cc:6
CaloSD::CaloSD
CaloSD(const std::string &aSDname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p, const SimTrackManager *, float timeSlice=1., bool ignoreTkID=false)
Definition: CaloSD.cc:33
DDName::name
const std::string & name() const
Returns the name.
Definition: DDName.cc:41
get
#define get
DreamSD::birk2_
double birk2_
Definition: DreamSD.h:56
edm::ESTransientHandle
Definition: ESTransientHandle.h:41
DreamSD::k_ScaleFromDD4HepToG4
static constexpr double k_ScaleFromDD4HepToG4
Definition: DreamSD.h:53
DDAxes::phi
multPhiCorr_741_25nsDY_cfi.px
px
Definition: multPhiCorr_741_25nsDY_cfi.py:10
DreamSD::crystalLength
double crystalLength(G4LogicalVolume *) const
Definition: DreamSD.cc:221
dataAnalyzerFineBiningParameters_cff.Pmax
Pmax
Definition: dataAnalyzerFineBiningParameters_cff.py:10
DDValue
Definition: DDValue.h:21
edm::LogVerbatim
Log< level::Info, true > LogVerbatim
Definition: MessageLogger.h:128
triggerObjects_cff.id
id
Definition: triggerObjects_cff.py:31
DDSolid::parameters
const std::vector< double > & parameters(void) const
Give the parameters of the solid.
Definition: DDSolid.cc:121
Skims_PA_cff.name
name
Definition: Skims_PA_cff.py:17
DDSolid
A DDSolid represents the shape of a part.
Definition: DDSolid.h:39
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
DDSpecificsMatchesValueFilter
Definition: DDFilter.h:70
AlignmentPI::index
index
Definition: AlignmentPayloadInspectorHelper.h:46
DreamSD::curve_LY
double curve_LY(const G4Step *, int)
Definition: DreamSD.cc:194
DreamSD::side_
int side_
Definition: DreamSD.h:60
sd
double sd
Definition: CascadeWrapper.h:113
DDFilteredView
Definition: DDFilteredView.h:20
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
DreamSD::slopeLY_
double slopeLY_
Definition: DreamSD.h:57
TableParser.table
table
Definition: TableParser.py:111
metDiagnosticParameterSet_cfi.nMax
nMax
Definition: metDiagnosticParameterSet_cfi.py:11
weight
Definition: weight.py:1
IdealGeometryRecord
Definition: IdealGeometryRecord.h:25
g
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
RemoveAddSevLevel.flag
flag
Definition: RemoveAddSevLevel.py:116
DreamSD::birk3_
double birk3_
Definition: DreamSD.h:56