diff --git a/CaloMC/fcl/prolog.fcl b/CaloMC/fcl/prolog.fcl index 8020baa350..a84b0531d2 100644 --- a/CaloMC/fcl/prolog.fcl +++ b/CaloMC/fcl/prolog.fcl @@ -21,7 +21,7 @@ CaloMC : { digiSampling : @local::HitMakerDigiSampling nMaxFragment : 1000 minPeakADC : @local::HitMakerMinPeakADC - readoutPEPerMeV : @local::readoutPEPerMeV + readoutPEPerMeV : @local::readoutPEPerMeVCsI MeVToADC : @local::MeVToADC diagLevel : 0 } @@ -55,7 +55,8 @@ CaloMC : { @table::CaloMC digitizationEnd : @local::HitMakerDigitizationEnd digitizationBuffer : 10 crystalNonUniformity : 1.05 - readoutPEPerMeV : @local::readoutPEPerMeV + readoutPEPerMeVCsI : @local::readoutPEPerMeVCsI + readoutPEPerMeVLyso : @local::readoutPEPerMeVLyso LRUCorrection : true BirksCorrection : true PEStatCorrection : true @@ -81,7 +82,8 @@ CaloMC : { @table::CaloMC minPeakADC : @local::HitMakerMinPeakADC nBinsPeak : 2 bufferDigi : 16 - readoutPEPerMeV : @local::readoutPEPerMeV + readoutPEPerMeVCsI : @local::readoutPEPerMeVCsI + readoutPEPerMeVLyso : @local::readoutPEPerMeVLyso MeVToADC : @local::MeVToADC diagLevel : 0 } diff --git a/CaloMC/src/CaloDigiMaker_module.cc b/CaloMC/src/CaloDigiMaker_module.cc index cfec407e3c..6911e64c0f 100644 --- a/CaloMC/src/CaloDigiMaker_module.cc +++ b/CaloMC/src/CaloDigiMaker_module.cc @@ -25,6 +25,7 @@ #include "Offline/RecoDataProducts/inc/CaloDigi.hh" #include "Offline/SeedService/inc/SeedService.hh" #include "Offline/MCDataProducts/inc/ProtonBunchTimeMC.hh" +#include "Offline/DataProducts/inc/CaloSiPMId.hh" #include "CLHEP/Vector/ThreeVector.h" #include "CLHEP/Random/RandPoissonQ.h" @@ -60,13 +61,14 @@ namespace mu2e { fhicl::Atom pbtmcTag { Name("protonBunchTimeMC"), Comment("ProtonBunchTimeMC producer") }; fhicl::Atom pulseFileName { Name("pulseFileName"), Comment("Calo pulse file name") }; fhicl::Atom pulseHistName { Name("pulseHistName"), Comment("Calo pulse hist name") }; - fhicl::Atom digitizationStart { Name("digitizationStart"), Comment("Start of digitization window relative to nominal pb time") }; - fhicl::Atom digitizationEnd { Name("digitizationEnd"), Comment("End of digitization window relative to nominal pb time")}; + fhicl::Atom digitizationStart { Name("digitizationStart"), Comment("Start of digitization window relative to nominal pb time") }; + fhicl::Atom digitizationEnd { Name("digitizationEnd"), Comment("End of digitization window relative to nominal pb time")}; fhicl::Atom addNoise { Name("addNoise"), Comment("Add noise to waveform") }; fhicl::Atom addRandomNoise { Name("addRandomNoise"), Comment("Add random salt and pepper noise") }; - fhicl::Atom digiSampling { Name("digiSampling"), Comment("Digitization time sampling") }; - fhicl::Atom pePerMeV { Name("readoutPEPerMeV"), Comment("Number of pe / MeV for Readout") }; - fhicl::Atom MeVToADC { Name("MeVToADC"), Comment("MeV to ADC conversion factor") }; + fhicl::Atom digiSampling { Name("digiSampling"), Comment("Digitization time sampling") }; + fhicl::Atom pePerMeVCsI { Name("readoutPEPerMeVCsI"), Comment("Number of pe / MeV for Readout for CsI") }; + fhicl::Atom pePerMeVLyso { Name("readoutPEPerMeVLyso"), Comment("Number of pe / MeV for Readout for LYSO") }; + fhicl::Atom MeVToADC { Name("MeVToADC"), Comment("MeV to ADC conversion factor") }; fhicl::Atom nBits { Name("nBits"), Comment("ADC Number of bits") }; fhicl::Atom nBinsPeak { Name("nBinsPeak"), Comment("Window size for finding local maximum to digitize wf") }; fhicl::Atom minPeakADC { Name("minPeakADC"), Comment("Minimum ADC hits of local peak to digitize") }; @@ -85,7 +87,8 @@ namespace mu2e { bufferDigi_ (config().bufferDigi()), startTimeBuffer_ (config().digiSampling()*config().bufferDigi()), maxADCCounts_ ((1 << config().nBits()) -1), - pePerMeV_ (config().pePerMeV()), + pePerMeVCsI_ (config().pePerMeVCsI()), + pePerMeVLyso_ (config().pePerMeVLyso()), MeVToADC_ (config().MeVToADC()), pulseShape_ (CaloPulseShape(config().pulseFileName(),config().pulseHistName(),config().digiSampling())), wfExtractor_ (config().bufferDigi(),config().nBinsPeak(),config().minPeakADC(),config().bufferDigi()), @@ -122,15 +125,16 @@ namespace mu2e { const art::ProductToken caloShowerToken_; art::InputTag ewMarkerTag_; art::InputTag pbtmcTag_; - double digitizationStart_; - double digitizationEnd_; - double timeFromProtonsToDRMarker_; - double digiSampling_; + float digitizationStart_; + float digitizationEnd_; + float timeFromProtonsToDRMarker_; + float digiSampling_; unsigned bufferDigi_; - double startTimeBuffer_; + float startTimeBuffer_; int maxADCCounts_; - double pePerMeV_; - double MeVToADC_; + float pePerMeVCsI_; + float pePerMeVLyso_; + float MeVToADC_; CaloPulseShape pulseShape_; CaloWFExtractor wfExtractor_; CLHEP::HepRandomEngine& engine_; @@ -227,8 +231,10 @@ namespace mu2e { bool CaloDigiMaker::fillROHits(unsigned iRO, std::vector& waveform, const CaloShowerROCollection& CaloShowerROs, const ProtonBunchTimeMC& pbtmc) { - bool isEmpty(true); - float scaleFactor(MeVToADC_/pePerMeV_); + bool isEmpty = true; + bool isCaphri = CaloSiPMId(iRO).crystal().isCaphri(); + auto pePerMeV = isCaphri ? pePerMeVLyso_ : pePerMeVCsI_; + auto scaleFactor = MeVToADC_/pePerMeV; for (const auto& CaloShowerRO : CaloShowerROs) { diff --git a/CaloMC/src/CaloHitTruthMatch_module.cc b/CaloMC/src/CaloHitTruthMatch_module.cc index e10cc9e927..f8f34f97a7 100644 --- a/CaloMC/src/CaloHitTruthMatch_module.cc +++ b/CaloMC/src/CaloHitTruthMatch_module.cc @@ -150,8 +150,8 @@ namespace mu2e { makeTruthMatch(event, *caloHitMCs, *caloHitMCTruth, *caloShowerMCTruth, primaryParticles); - event.put(std::move(caloHitMCTruth)); event.put(std::move(caloHitMCs)); + event.put(std::move(caloHitMCTruth)); if (fillDetailedMC_) event.put(std::move(caloShowerMCTruth)); } @@ -242,7 +242,7 @@ namespace mu2e { //sort CaloEDepMC by decreasing energy, create caloHit and keep track of hit -> MChit association std::sort(edeps.begin(),edeps.end(),[](const auto& a, const auto& b){return a.energyDep() > b.energyDep();}); - caloHitMCs.emplace_back(CaloHitMC(std::move(edeps))); + caloHitMCs.emplace_back(CaloHitMC(std::move(edeps),hit.crystalID())); art::Ptr hitMCPtr = art::Ptr(hitMCProductID, caloHitMCs.size()-1, hitMCProductGetter); CaloHitTruthMatch.addSingle(hitPtr,hitMCPtr); diff --git a/CaloMC/src/CaloShowerROMaker_module.cc b/CaloMC/src/CaloShowerROMaker_module.cc index 11f719b69b..120afaddb6 100644 --- a/CaloMC/src/CaloShowerROMaker_module.cc +++ b/CaloMC/src/CaloShowerROMaker_module.cc @@ -96,7 +96,8 @@ namespace mu2e { fhicl::Atom digitizationEnd { Name("digitizationEnd"), Comment("Maximum time of hit to be digitized") }; fhicl::Atom digitizationBuffer { Name("digitizationBuffer"), Comment("Digi time buffer for photon propagation inside crystal") }; fhicl::Atom crystalNonUniformity { Name("crystalNonUniformity"), Comment("LRU parameter 0") }; - fhicl::Atom pePerMeV { Name("readoutPEPerMeV"), Comment("Number of pe / MeV for Readout") }; + fhicl::Atom pePerMeVCsI { Name("readoutPEPerMeVCsI"), Comment("Number of pe / MeV for Readout for CsI") }; + fhicl::Atom pePerMeVLyso { Name("readoutPEPerMeVLyso"), Comment("Number of pe / MeV for Readout for LYSO") }; fhicl::Atom LRUCorrection { Name("LRUCorrection"), Comment("Include LRU corrections") }; fhicl::Atom BirksCorrection { Name("BirksCorrection"), Comment("Include Birks corrections") }; fhicl::Atom PEStatCorrection { Name("PEStatCorrection"), Comment("Include PE Poisson fluctuations") }; @@ -112,7 +113,8 @@ namespace mu2e { digitizationEnd_ (config().digitizationEnd()), digitizationBuffer_ (config().digitizationBuffer()), crystalNonUniformity_(config().crystalNonUniformity()), - pePerMeV_ (config().pePerMeV()), + pePerMeVCsI_ (config().pePerMeVCsI()), + pePerMeVLyso_ (config().pePerMeVLyso()), LRUCorrection_ (config().LRUCorrection()), BirksCorrection_ (config().BirksCorrection()), PEStatCorrection_ (config().PEStatCorrection()), @@ -151,7 +153,8 @@ namespace mu2e { float digitizationEnd_; float digitizationBuffer_; float crystalNonUniformity_; - float pePerMeV_; + float pePerMeVCsI_; + float pePerMeVLyso_; bool LRUCorrection_; bool BirksCorrection_; bool PEStatCorrection_; @@ -274,11 +277,14 @@ namespace mu2e { if (BirksCorrection_) edep_corr = step.energyDepBirks(); if (LRUCorrection_) edep_corr = LRUCorrection(crystalID, posZ/cryhalflength, edep_corr); + bool isCaphri = CrystalId(crystalID).isCaphri(); + float pePerMeV = isCaphri ? pePerMeVLyso_ : pePerMeVCsI_; + // Generate individual PEs and their arrival times for (int i=0; i PETime(NPE,hitTime); @@ -293,7 +299,7 @@ namespace mu2e { if (diagLevel_ > 1) for (const auto& time : PETime) hTime_->Fill(2.0*cryhalflength-posZ,time-hitTime); diagSum.totNPE += NPE; - diagSum.totEdepNPE += double(NPE)/pePerMeV_/2.0; //average between the two RO + diagSum.totEdepNPE += double(NPE)/pePerMeV/2.0; //average between the two RO } diagSum.totEdep += step.energyDepG4(); diagSum.totEdepCorr += edep_corr; diff --git a/CaloReco/fcl/common.fcl b/CaloReco/fcl/common.fcl index 7e54a02f73..5efd0185b1 100644 --- a/CaloReco/fcl/common.fcl +++ b/CaloReco/fcl/common.fcl @@ -5,7 +5,8 @@ BEGIN_PROLOG HitMakerMinPeakADC : 16 pulseFileName : "OfflineData/ConditionsService/CsI-waveform-2020-02-12.root" pulseHistName : "h_waveform" - readoutPEPerMeV : 30.0 + readoutPEPerMeVCsI : 30.0 + readoutPEPerMeVLyso : 500.0 MeVToADC : 16.0 ADCToMeV : 0.0625 END_PROLOG diff --git a/MCDataProducts/inc/CaloHitMC.hh b/MCDataProducts/inc/CaloHitMC.hh index c4a8826b55..f7a860d072 100644 --- a/MCDataProducts/inc/CaloHitMC.hh +++ b/MCDataProducts/inc/CaloHitMC.hh @@ -18,8 +18,15 @@ namespace mu2e class CaloHitMC { public: - CaloHitMC() : edeps_() {}; - CaloHitMC(const std::vector& edeps) : edeps_(edeps) {}; + CaloHitMC() : edeps_(), crystalID_(-1) {}; + CaloHitMC(const std::vector& edeps, int crystalId) : + edeps_(edeps), + crystalID_(crystalId) + {}; + CaloHitMC(const std::vector&& edeps, int crystalId) : + edeps_(std::move(edeps)), + crystalID_(crystalId) + {}; void resetSim(SimParticleRemapping const& remap); @@ -30,9 +37,11 @@ namespace mu2e float time () const {return edeps_.empty() ? 0.0 : edeps_.at(0).time();} float totalEnergyDep () const; float totalEnergyDepG4() const; + int crystalID () const { return crystalID_;} private: std::vector edeps_; + int crystalID_; }; using CaloHitMCCollection = std::vector; diff --git a/MCDataProducts/inc/CaloMCTruthAssns.hh b/MCDataProducts/inc/CaloMCTruthAssns.hh index cff49e72aa..c27faa23a7 100644 --- a/MCDataProducts/inc/CaloMCTruthAssns.hh +++ b/MCDataProducts/inc/CaloMCTruthAssns.hh @@ -12,7 +12,7 @@ namespace mu2e { using CaloShowerMCTruthAssn = art::Assns>; - using CaloHitMCTruthAssn = art::Assns ; + using CaloHitMCTruthAssn = art::Assns ; using CaloClusterMCTruthAssn = art::Assns ; }