From c82e90663fa1b0ced70f00996f9c21aae1904e9b Mon Sep 17 00:00:00 2001 From: arvindkhuntia Date: Tue, 22 Sep 2020 15:07:34 +0200 Subject: [PATCH 1/3] [FV0] pile-up simulation --- .../include/FV0Simulation/Digitizer.h | 63 ++++- .../include/FV0Simulation/FV0DigParam.h | 23 +- .../FIT/FV0/simulation/src/Digitizer.cxx | 237 +++++++++++++----- .../src/FV0DigitizerSpec.cxx | 8 +- 4 files changed, 246 insertions(+), 85 deletions(-) diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h index e4f9fa53d4a38..8e9d7fc7a9f98 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h @@ -35,7 +35,7 @@ class Digitizer public: Digitizer() - : mTimeStamp(0), mIntRecord(), mEventId(-1), mSrcId(-1), mMCLabels(), mPmtChargeVsTime(), mNBins(), mPmtResponseGlobal(), mPmtResponseTemp() + : mTimeStamp(0), mIntRecord(), mEventId(-1), mSrcId(-1), mMCLabels(), mCache(), mPmtChargeVsTime(), mNBins(), NTimeBinsPerBC(), mPmtResponseGlobal(), mPmtResponseTemp() { } @@ -53,25 +53,68 @@ class Digitizer void setSrcId(Int_t id) { mSrcId = id; } void setInteractionRecord(const InteractionTimeRecord& ir) { mIntRecord = ir; } - void process(const std::vector& hits); - void analyseWaveformsAndStore(std::vector& digitsBC, - std::vector& digitsCh, - dataformats::MCTruthContainer& labels); + void process(const std::vector& hits, std::vector& digitsBC, + std::vector& digitsCh, o2::dataformats::MCTruthContainer& labels); + + void flush(std::vector& digitsBC, + std::vector& digitsCh, + o2::dataformats::MCTruthContainer& labels); const InteractionRecord& getInteractionRecord() const { return mIntRecord; } InteractionRecord& getInteractionRecord(InteractionRecord& src) { return mIntRecord; } uint32_t getOrbit() const { return mIntRecord.orbit; } uint16_t getBC() const { return mIntRecord.bc; } + using ChannelBCDataF = std::vector; + + struct BCCache : public o2::InteractionRecord { + std::vector labels; + //std::array Cfd_times; + //std::array Charges; + std::array mPmtChargeVsTime = {}; + + void clear() + { + for (auto& channel : mPmtChargeVsTime) { + std::fill(channel.begin(), channel.end(), 0.); + } + labels.clear(); + } + + //bool IsCounted = false; + //bool isWritten = 0; + //int EvID; + // void setEvID (const int ev) {EvID = ev;} + //BCCache(); + BCCache& operator=(const o2::InteractionRecord& ir) + { + o2::InteractionRecord::operator=(ir); + return *this; + } + void print() const; + }; + private: + static constexpr int BCCacheMin = 0, BCCacheMax = 7, NBC2Cache = 1 + BCCacheMax - BCCacheMin; + void createPulse(float mipFraction, int parID, double hitTime, std::array const& cachedIR, + int nCachedIR, const int& detID); + long mTimeStamp; // TF (run) timestamp - InteractionRecord mIntRecord; // Interaction record (orbit, bc) -> InteractionTimeRecord + InteractionTimeRecord mIntRecord; // Interaction record (orbit, bc) -> InteractionTimeRecord Int_t mEventId; // ID of the current event Int_t mSrcId; // signal, background or QED - std::vector mMCLabels; + std::deque mMCLabels; + std::deque mCache; + + BCCache& setBCCache(const o2::InteractionRecord& ir); + BCCache* getBCCache(const o2::InteractionRecord& ir); + + void storeBC(const BCCache& bc, std::vector& digitsBC, std::vector& digitsCh, + o2::dataformats::MCTruthContainer& labels); std::array, Constants::nFv0Channels> mPmtChargeVsTime; // Charge time series aka analogue signal pulse from PM - UInt_t mNBins; // Number of bins in pulse series + UInt_t mNBins; // + UInt_t NTimeBinsPerBC; Float_t mBinSize; // Time width of the pulse bin - HPTDC resolution /// vectors to store the PMT signal from cosmic muons @@ -80,7 +123,9 @@ class Digitizer /// Internal helper methods related to conversion of energy-deposition into el. signal Int_t SimulateLightYield(Int_t pmt, Int_t nPhot) const; - Float_t SimulateTimeCfd(Int_t channel) const; + Float_t SimulateTimeCfd(const ChannelBCDataF& pulse) const; + Float_t IntegrateCharge(const ChannelBCDataF& pulse) const; + // Float_t SimulateTimeCfd(Int_t channel, Int_t iCache) const; /// Functions related to splitting ring-5 cell signal to two readout channels static float getDistFromCellCenter(UInt_t cellId, double hitx, double hity); diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h index c70c4464ccf6a..c03e650c2f609 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h @@ -13,13 +13,13 @@ #include "CommonUtils/ConfigurableParam.h" #include "CommonUtils/ConfigurableParamHelper.h" +#include "CommonConstants/PhysicsConstants.h" namespace o2 { namespace fv0 { // parameters of FV0 digitization / transport simulation - struct FV0DigParam : public o2::conf::ConfigurableParamHelper { float intrinsicTimeRes = 0.91; // time resolution float photoCathodeEfficiency = 0.23; // quantum efficiency = nOfPhotoE_emitted_by_photocathode / nIncidentPhotons @@ -28,23 +28,22 @@ struct FV0DigParam : public o2::conf::ConfigurableParamHelper { float pmtTransitTime = 9.5; // PMT response time (corresponds to 1.9 ns rise time) float pmtTransparency = 0.25; // Transparency of the first dynode of the PMT float shapeConst = 1.18059e-14; // Crystal ball const parameter - float shapeMean = 1.70518e+01; // Crystal ball mean parameter + float shapeMean = 9.5; // Crystal ball mean parameter float shapeAlpha = -6.56586e-01; // Crystal ball alpha parameter float shapeN = 2.36408e+00; // Crystal ball N parameter float shapeSigma = 3.55445; // Crystal ball sigma parameter - float timeShiftCfd = 3.3; // TODO: adjust after PM design for FV0 is fixed + //float timeShiftCfd = 3.3; // From the cosmic measurements of FV0 [keep it for reference] + float timeShiftCfd = 5.3; // TODO: adjust after FV0 with FEE measurements are done float singleMipThreshold = 3.0; // in [MeV] of deposited energy float waveformNbins = 10000; // number of bins for the analog pulse waveform float waveformBinWidth = 0.01302; // number of bins for the analog - float timeCompensate = 23.25; // in ns - float chargeIntBinMin = (timeCompensate - 6.0) / waveformBinWidth; //Charge integration offset (cfd mean time - 6 ns) - float chargeIntBinMax = (timeCompensate + 14.0) / waveformBinWidth; //Charge integration offset (cfd mean time + 14 ns) - bool isIntegrateFull = false; - float cfdCheckWindow = 2.5; // in ns - int avgPhElectron = 201; - - //Optimization-related, derived constants - float oneOverPmtTransitTime2 = 1.0 / (pmtTransitTime * pmtTransitTime); + float avgCfdTimeForMip = 8.63; // in ns to shift the CFD time to zero TODO do ring wise + int chargeIntBinMin = (avgCfdTimeForMip - 6.0) / waveformBinWidth; //Charge integration offset (cfd mean time - 6 ns) + int chargeIntBinMax = (avgCfdTimeForMip + 14.0) / waveformBinWidth; //Charge integration offset (cfd mean time + 14 ns) + bool isIntegrateFull = false; // Full charge integration widow in 25 ns + float cfdCheckWindow = 2.5; // time window for the cfd in ns to trigger the charge integration + int avgNumberPhElectronPerMip = 201; // avg number of photo-electrons per MIP + float globalTimeOfFlight = 315.0 / o2::constants::physics::LightSpeedCm2NS; //TODO check the correct value for distance of FV0 to IP O2ParamDef(FV0DigParam, "FV0DigParam"); }; diff --git a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx index 2cc2bafbf07d2..17add7ae15a9b 100644 --- a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx +++ b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx @@ -37,6 +37,8 @@ void Digitizer::init() mNBins = FV0DigParam::Instance().waveformNbins; //Will be computed using detector set-up from CDB mBinSize = FV0DigParam::Instance().waveformBinWidth; //Will be set-up from CDB + NTimeBinsPerBC = int(o2::constants::lhc::LHCBunchSpacingNS / mBinSize); + for (Int_t detID = 0; detID < Constants::nFv0Channels; detID++) { mPmtChargeVsTime[detID].resize(mNBins); } @@ -50,7 +52,7 @@ void Digitizer::init() FV0DigParam::Instance().shapeN); // PMT response per hit [Global] - float x = mBinSize; + float x = mBinSize / 2.0; /// Calculate at BinCenter mPmtResponseGlobal.resize(mNBins); for (Int_t j = 0; j < mPmtResponseGlobal.size(); ++j) { mPmtResponseGlobal[j] = signalShapeFn.Eval(x); @@ -61,26 +63,26 @@ void Digitizer::init() LOG(INFO) << "V0Digitizer::init -> finished"; } -void Digitizer::process(const std::vector& hits) +void Digitizer::process(const std::vector& hits, std::vector& digitsBC, + std::vector& digitsCh, o2::dataformats::MCTruthContainer& labels) { LOG(INFO) << "[FV0] Digitizer::process(): begin with " << hits.size() << " hits"; + flush(digitsBC, digitsCh, labels); // flush cached signal which cannot be affect by new event std::vector hitIdx(hits.size()); std::iota(std::begin(hitIdx), std::end(hitIdx), 0); - std::sort(std::begin(hitIdx), std::end(hitIdx), [&hits](int a, int b) { return hits[a].GetTrackID() < hits[b].GetTrackID(); }); - + std::sort(std::begin(hitIdx), std::end(hitIdx), + [&hits](int a, int b) { return hits[a].GetTrackID() < hits[b].GetTrackID(); }); Int_t parentIdPrev = -10; // use ordered hits for (auto ids : hitIdx) { const auto& hit = hits[ids]; Int_t detId = hit.GetDetectorID(); Double_t hitEdep = hit.GetHitValue() * 1e3; //convert to MeV - // TODO: check how big is inaccuracy if more than 1 'below-threshold' particles hit the same detector cell if (hitEdep < FV0DigParam::Instance().singleMipThreshold) { continue; } - float distanceFromXc = 0; if (Geometry::instance()->isRing5(detId)) { distanceFromXc = getDistFromCellCenter(detId, hit.GetX(), hit.GetY()); @@ -90,14 +92,13 @@ void Digitizer::process(const std::vector& hits) while (iChannelPerCell < 2) { // loop over 2 channels, into which signal from each cell in ring 5 is split if (Geometry::instance()->isRing5(detId)) { // The first channel number is located counter-clockwise from the cell center - // and remains identical to the detector number, the second one is clockwise and incremented by 8 + // and remains identical to the detector number, the second one is clockwise and incremented by 8 if (iChannelPerCell == 1) { detId += 8; } - // Split signal magnitude to fractions depending on the distance of the hit from the cell center hitEdep = (hit.GetHitValue() * 1e3) * getSignalFraction(distanceFromXc, iChannelPerCell == 0); - //LOG(INFO) << " detId: " << detId << "-" << iChannelPerCell << " hitEdep: " << hitEdep << " distanceFromXc: " << distanceFromXc; + // LOG(INFO) << " detId: " << detId << "-" << iChannelPerCell << " hitEdep: " << hitEdep << " distanceFromXc: " << distanceFromXc; ++iChannelPerCell; } else { iChannelPerCell = 2; // not a ring 5 cell -> don't repeat the loop @@ -105,64 +106,118 @@ void Digitizer::process(const std::vector& hits) Double_t const nPhotons = hitEdep * DP::N_PHOTONS_PER_MEV; float const nPhE = SimulateLightYield(detId, nPhotons); - float avgMip = float(nPhE / FV0DigParam::Instance().avgPhElectron); + float mipFraction = float(nPhE / FV0DigParam::Instance().avgNumberPhElectronPerMip); Float_t const hitTime = hit.GetTime() * 1e9; - //------------------------------------------------------------------------------ - mPmtResponseTemp = mPmtResponseGlobal; - Size_t NBinShift = std::lround(hitTime / FV0DigParam::Instance().waveformBinWidth); - for (int m = 0; m < NBinShift; m++) - mPmtResponseTemp.push_back(0); - // rotate the vector element to shift all the elements by hit time - std::rotate(mPmtResponseTemp.rbegin(), mPmtResponseTemp.rbegin() + NBinShift, mPmtResponseTemp.rend()); - mPmtResponseTemp.resize(FV0DigParam::Instance().waveformNbins); - - for (int mm = 0; mm < mPmtResponseTemp.size(); mm++) { - mPmtChargeVsTime[detId][mm] += (mPmtResponseTemp[mm] * avgMip); - } - // Charged particles in MCLabel - Int_t const parentId = hit.GetTrackID(); - if (parentId != parentIdPrev) { - mMCLabels.emplace_back(parentId, mEventId, mSrcId, detId); - parentIdPrev = parentId; - } + Float_t timeHit = hitTime; + timeHit += mIntRecord.getTimeNS(); + o2::InteractionTimeRecord irHit(timeHit); + + std::array cachedIR; + int nCachedIR = 0; + + for (int i = BCCacheMin; i < BCCacheMax + 1; i++) { + double tNS = timeHit + o2::constants::lhc::LHCBunchSpacingNS * i; + cachedIR[nCachedIR].setFromNS(tNS); + if (tNS < 0 && cachedIR[nCachedIR] > irHit) { + continue; // don't go to negative BC/orbit (it will wrap) + } + setBCCache(cachedIR[nCachedIR++]); // ensure existence of cached container + } //BCCache loop + createPulse(mipFraction, hit.GetTrackID(), hitTime, cachedIR, nCachedIR, detId); + } //while loop + } //hitloop +} + +void Digitizer::createPulse(float mipFraction, int parID, double hitTime, + std::array const& cachedIR, int nCachedIR, const int& detId) +{ + + bool added[nCachedIR]; + for (int ir = 0; ir < nCachedIR; ir++) { + added[ir] = false; + } + + for (int ir = 0; ir < NBC2Cache; ir++) { + auto bcCache = getBCCache(cachedIR[ir]); + for (int ich = 0; ich < Constants::nFv0Channels; ich++) { + (*bcCache).mPmtChargeVsTime[ich].resize(NTimeBinsPerBC); + } + } + mPmtResponseTemp = mPmtResponseGlobal; + ///Time of flight subtracted from Hit time //TODO have different TOF according to thr ring number + Size_t NBinShift = std::lround((hitTime - FV0DigParam::Instance().globalTimeOfFlight) / FV0DigParam::Instance().waveformBinWidth); + + for (int m = 0; m < NBinShift; m++) { + mPmtResponseTemp.push_back(0); + } + /// rotate the vector element to shift all the elements by hit time + std::rotate(mPmtResponseTemp.rbegin(), mPmtResponseTemp.rbegin() + NBinShift, mPmtResponseTemp.rend()); + mPmtResponseTemp.resize(FV0DigParam::Instance().waveformNbins); + + for (int ir = 0; ir < int(mPmtResponseTemp.size() / NTimeBinsPerBC); ir++) { + auto bcCache = getBCCache(cachedIR[ir]); + + for (int iBin = 0; iBin < NTimeBinsPerBC; iBin++) { + (*bcCache).mPmtChargeVsTime[detId][iBin] += (mPmtResponseTemp[ir * NTimeBinsPerBC + iBin] * mipFraction); + } + added[ir] = true; + } + ///Add MC labels to BCs for those contributed to the PMT signal + for (int ir = 0; ir < nCachedIR; ir++) { + if (added[ir]) { + auto bcCache = getBCCache(cachedIR[ir]); + (*bcCache).labels.emplace_back(parID, mEventId, mSrcId, detId); + } + } +} + +void Digitizer::flush(std::vector& digitsBC, std::vector& digitsCh, + o2::dataformats::MCTruthContainer& labels) +{ + int nCached = mCache.size(); + if (nCached < 1) { + return; + } + for (auto bc : mCache) { + if (mIntRecord.differenceInBC(bc) > NBC2Cache) { // Build events those are separated by NBC2Cache BCs from current BC + storeBC(bc, digitsBC, digitsCh, labels); + mCache.pop_front(); + } else { + return; } - } //hit loop + } } -void Digitizer::analyseWaveformsAndStore(std::vector& digitsBC, - std::vector& digitsCh, - dataformats::MCTruthContainer& labels) +void Digitizer::storeBC(const BCCache& bc, std::vector& digitsBC, std::vector& digitsCh, + o2::dataformats::MCTruthContainer& labels) + { - // Sum charge of all time bins to get total charge collected for a given channel - size_t const first = digitsCh.size(); + int first = digitsCh.size(); size_t nStored = 0; - for (Int_t ipmt = 0; ipmt < Constants::nFv0Channels; ++ipmt) { - float time = SimulateTimeCfd(ipmt) - FV0DigParam::Instance().timeCompensate; - if (time < -FV0DigParam::Instance().cfdCheckWindow || time > FV0DigParam::Instance().cfdCheckWindow) - continue; - Float_t totalCharge = 0.0f; - auto const& analogSignal = mPmtChargeVsTime[ipmt]; + double totalCharge = 0; + for (int iPmt = 0; iPmt < Constants::nFv0Channels; iPmt++) { + double cfdWithOffset = SimulateTimeCfd(bc.mPmtChargeVsTime[iPmt]); + double cfdZero = cfdWithOffset - FV0DigParam::Instance().avgCfdTimeForMip; - int chargeIntMin = FV0DigParam::Instance().isIntegrateFull ? 0 : FV0DigParam::Instance().chargeIntBinMin; - int chargeIntMax = FV0DigParam::Instance().isIntegrateFull ? FV0DigParam::Instance().waveformNbins : FV0DigParam::Instance().chargeIntBinMax; + if (cfdZero < -FV0DigParam::Instance().cfdCheckWindow || cfdZero > FV0DigParam::Instance().cfdCheckWindow) + continue; - for (Int_t iTimeBin = chargeIntMin; iTimeBin < chargeIntMax; ++iTimeBin) { - Float_t const timeBinCharge = mPmtChargeVsTime[ipmt][iTimeBin]; - totalCharge += timeBinCharge; - } + //LOG(INFO) << "time inside analyse and store =========> " << cfdZero <<" detid "<(std::round(time)), static_cast(std::round(totalCharge))); + cfdZero *= DP::INV_TIME_PER_TDCCHANNEL; + + digitsCh.emplace_back(iPmt, static_cast(std::round(cfdZero)), + static_cast(std::round(totalCharge))); + totalCharge = 0; ++nStored; } - - // Send MClabels and digitsBC to storage - size_t const nBC = digitsBC.size(); - digitsBC.emplace_back(first, nStored, mIntRecord); - for (auto const& lbl : mMCLabels) { + int nBC = digitsBC.size(); + digitsBC.emplace_back(first, nStored, bc); + for (const auto& lbl : bc.labels) { labels.addElement(nBC, lbl); } - mMCLabels.clear(); } // ------------------------------------------------------------------------------- @@ -181,16 +236,35 @@ Int_t Digitizer::SimulateLightYield(Int_t pmt, Int_t nPhot) const : gRandom->Gaus((p * nPhot) + 0.5, TMath::Sqrt(p * (1. - p) * nPhot))); return n; } +//--------------------------------------------------------------------------- +Float_t Digitizer::IntegrateCharge(const ChannelBCDataF& pulse) const +{ + int chargeIntMin = FV0DigParam::Instance().isIntegrateFull ? 0 : FV0DigParam::Instance().chargeIntBinMin; + int chargeIntMax = FV0DigParam::Instance().isIntegrateFull ? NTimeBinsPerBC : FV0DigParam::Instance().chargeIntBinMax; -Float_t Digitizer::SimulateTimeCfd(Int_t channel) const + Float_t totalCharge = 0.0f; + for (int iTimeBin = chargeIntMin; iTimeBin < chargeIntMax; iTimeBin++) { + Float_t const timeBinCharge = pulse[iTimeBin]; + //LOG(INFO)<= binShift - ? 5.0f * mPmtChargeVsTime[channel][iTimeBin - binShift] - mPmtChargeVsTime[channel][iTimeBin] - : -mPmtChargeVsTime[channel][iTimeBin]); + Float_t sigPrev = -pulse[0]; //[0]; + for (Int_t iTimeBin = 1; iTimeBin < NTimeBinsPerBC; ++iTimeBin) { + Float_t const sigCurrent = (iTimeBin >= binShift ? 5.0f * pulse[iTimeBin - binShift] - pulse[iTimeBin] : -pulse[iTimeBin]); if (sigPrev < 0.0f && sigCurrent >= 0.0f) { timeCfd = Float_t(iTimeBin) * mBinSize; break; @@ -200,7 +274,6 @@ Float_t Digitizer::SimulateTimeCfd(Int_t channel) const return timeCfd; } -// The Distance is positive for top half-sectors (when the hit position is above the cell center (has higher y)) float Digitizer::getDistFromCellCenter(UInt_t cellId, double hitx, double hity) { Geometry* geo = Geometry::instance(); @@ -212,7 +285,7 @@ float Digitizer::getDistFromCellCenter(UInt_t cellId, double hitx, double hity) double a = -(y0 - pCell->y) / (x0 - pCell->x); double b = 1; double c = -(y0 - a * x0); - // Return the distance from hit to this line + //Return the distance from hit to this line return (a * hitx + b * hity + c) / TMath::Sqrt(a * a + b * b); } @@ -225,3 +298,43 @@ float Digitizer::getSignalFraction(float distanceFromXc, bool isFirstChannel) return isFirstChannel ? (1. - fraction) : fraction; } } + +//_____________________________________________________________________________ +o2::fv0::Digitizer::BCCache& Digitizer::setBCCache(const o2::InteractionRecord& ir) +{ + if (mCache.empty() || mCache.back() < ir) { + mCache.emplace_back(); + auto& cb = mCache.back(); + cb = ir; + return cb; + } + if (mCache.front() > ir) { + mCache.emplace_front(); + auto& cb = mCache.front(); + cb = ir; + return cb; + } + + for (auto cb = mCache.begin(); cb != mCache.end(); cb++) { + if ((*cb) == ir) { + return *cb; + } + if (ir < (*cb)) { + auto cbnew = mCache.emplace(cb); // insert new element before cb + (*cbnew) = ir; + return (*cbnew); + } + } + return mCache.front(); +} +//_____________________________________________________________________________ +o2::fv0::Digitizer::BCCache* Digitizer::getBCCache(const o2::InteractionRecord& ir) +{ + // get pointer on existing cache + for (auto cb = mCache.begin(); cb != mCache.end(); cb++) { + if ((*cb) == ir) { + return &(*cb); + } + } + return nullptr; +} \ No newline at end of file diff --git a/Steer/DigitizerWorkflow/src/FV0DigitizerSpec.cxx b/Steer/DigitizerWorkflow/src/FV0DigitizerSpec.cxx index 4de8f3c51d9e5..7e501e6f2ed81 100644 --- a/Steer/DigitizerWorkflow/src/FV0DigitizerSpec.cxx +++ b/Steer/DigitizerWorkflow/src/FV0DigitizerSpec.cxx @@ -85,12 +85,16 @@ class FV0DPLDigitizerTask : public o2::base::BaseDPLDigitizer // call actual digitization procedure mDigitizer.setEventId(part.entryID); mDigitizer.setSrcId(part.sourceID); - mDigitizer.process(hits); + mDigitizer.process(hits, mDigitsBC, mDigitsCh, mLabels); } - mDigitizer.analyseWaveformsAndStore(mDigitsBC, mDigitsCh, mLabels); LOG(INFO) << "[FV0] Has " << mDigitsBC.size() << " BC elements, " << mDigitsCh.size() << " mDigitsCh elements"; } + o2::InteractionTimeRecord terminateIR; + terminateIR.orbit = 0xffffffff; // supply IR in the infinite future to flush all cached BC + mDigitizer.setInteractionRecord(terminateIR); + mDigitizer.flush(mDigitsBC, mDigitsCh, mLabels); + // here we have all digits and we can send them to consumer (aka snapshot it onto output) LOG(INFO) << "FV0: Sending " << mDigitsBC.size() << " digitsBC and " << mDigitsCh.size() << " digitsCh."; From b901683e1029e9cce4f354be1965b5f8930dec78 Mon Sep 17 00:00:00 2001 From: arvindkhuntia Date: Mon, 5 Oct 2020 13:40:02 +0200 Subject: [PATCH 2/3] Ruben's comment implemented --- .../include/FV0Simulation/Digitizer.h | 15 +++----- .../include/FV0Simulation/FV0DigParam.h | 1 + .../FIT/FV0/simulation/src/Digitizer.cxx | 35 ++++++++----------- 3 files changed, 20 insertions(+), 31 deletions(-) diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h index 8e9d7fc7a9f98..c03614b0d5ecd 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Digitizer.h @@ -35,7 +35,7 @@ class Digitizer public: Digitizer() - : mTimeStamp(0), mIntRecord(), mEventId(-1), mSrcId(-1), mMCLabels(), mCache(), mPmtChargeVsTime(), mNBins(), NTimeBinsPerBC(), mPmtResponseGlobal(), mPmtResponseTemp() + : mTimeStamp(0), mIntRecord(), mEventId(-1), mSrcId(-1), mMCLabels(), mCache(), mPmtChargeVsTime(), mNBins(), mNTimeBinsPerBC(), mPmtResponseGlobal(), mPmtResponseTemp() { } @@ -69,8 +69,6 @@ class Digitizer struct BCCache : public o2::InteractionRecord { std::vector labels; - //std::array Cfd_times; - //std::array Charges; std::array mPmtChargeVsTime = {}; void clear() @@ -81,11 +79,6 @@ class Digitizer labels.clear(); } - //bool IsCounted = false; - //bool isWritten = 0; - //int EvID; - // void setEvID (const int ev) {EvID = ev;} - //BCCache(); BCCache& operator=(const o2::InteractionRecord& ir) { o2::InteractionRecord::operator=(ir); @@ -96,8 +89,8 @@ class Digitizer private: static constexpr int BCCacheMin = 0, BCCacheMax = 7, NBC2Cache = 1 + BCCacheMax - BCCacheMin; - void createPulse(float mipFraction, int parID, double hitTime, std::array const& cachedIR, - int nCachedIR, const int& detID); + void createPulse(float mipFraction, int parID, const double hitTime, std::array const& cachedIR, + int nCachedIR, const int detID); long mTimeStamp; // TF (run) timestamp InteractionTimeRecord mIntRecord; // Interaction record (orbit, bc) -> InteractionTimeRecord @@ -114,7 +107,7 @@ class Digitizer std::array, Constants::nFv0Channels> mPmtChargeVsTime; // Charge time series aka analogue signal pulse from PM UInt_t mNBins; // - UInt_t NTimeBinsPerBC; + UInt_t mNTimeBinsPerBC; Float_t mBinSize; // Time width of the pulse bin - HPTDC resolution /// vectors to store the PMT signal from cosmic muons diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h index c03e650c2f609..4697f4a804d46 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h @@ -35,6 +35,7 @@ struct FV0DigParam : public o2::conf::ConfigurableParamHelper { //float timeShiftCfd = 3.3; // From the cosmic measurements of FV0 [keep it for reference] float timeShiftCfd = 5.3; // TODO: adjust after FV0 with FEE measurements are done float singleMipThreshold = 3.0; // in [MeV] of deposited energy + float singleHitTimeThreshold = 120.0;// in [ns] to skip very slow particles float waveformNbins = 10000; // number of bins for the analog pulse waveform float waveformBinWidth = 0.01302; // number of bins for the analog float avgCfdTimeForMip = 8.63; // in ns to shift the CFD time to zero TODO do ring wise diff --git a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx index 17add7ae15a9b..3490c19ed5d0e 100644 --- a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx +++ b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx @@ -37,7 +37,7 @@ void Digitizer::init() mNBins = FV0DigParam::Instance().waveformNbins; //Will be computed using detector set-up from CDB mBinSize = FV0DigParam::Instance().waveformBinWidth; //Will be set-up from CDB - NTimeBinsPerBC = int(o2::constants::lhc::LHCBunchSpacingNS / mBinSize); + mNTimeBinsPerBC = int(o2::constants::lhc::LHCBunchSpacingNS / mBinSize); for (Int_t detID = 0; detID < Constants::nFv0Channels; detID++) { mPmtChargeVsTime[detID].resize(mNBins); @@ -79,8 +79,9 @@ void Digitizer::process(const std::vector& hits, std::vector FV0DigParam::Instance().singleHitTimeThreshold) { continue; } float distanceFromXc = 0; @@ -107,11 +108,9 @@ void Digitizer::process(const std::vector& hits, std::vector cachedIR; int nCachedIR = 0; @@ -124,12 +123,13 @@ void Digitizer::process(const std::vector& hits, std::vector const& cachedIR, int nCachedIR, const int& detId) +void Digitizer::createPulse(float mipFraction, int parID,const double hitTime, + std::array const& cachedIR, int nCachedIR, const int detId) { bool added[nCachedIR]; @@ -140,25 +140,20 @@ void Digitizer::createPulse(float mipFraction, int parID, double hitTime, for (int ir = 0; ir < NBC2Cache; ir++) { auto bcCache = getBCCache(cachedIR[ir]); for (int ich = 0; ich < Constants::nFv0Channels; ich++) { - (*bcCache).mPmtChargeVsTime[ich].resize(NTimeBinsPerBC); + (*bcCache).mPmtChargeVsTime[ich].resize(mNTimeBinsPerBC); } } - mPmtResponseTemp = mPmtResponseGlobal; + ///Time of flight subtracted from Hit time //TODO have different TOF according to thr ring number Size_t NBinShift = std::lround((hitTime - FV0DigParam::Instance().globalTimeOfFlight) / FV0DigParam::Instance().waveformBinWidth); + mPmtResponseTemp.resize(FV0DigParam::Instance().waveformNbins, 0.); + std::memcpy(&mPmtResponseTemp[NBinShift], &mPmtResponseGlobal[0], sizeof(double)*(FV0DigParam::Instance().waveformNbins-NBinShift)); - for (int m = 0; m < NBinShift; m++) { - mPmtResponseTemp.push_back(0); - } - /// rotate the vector element to shift all the elements by hit time - std::rotate(mPmtResponseTemp.rbegin(), mPmtResponseTemp.rbegin() + NBinShift, mPmtResponseTemp.rend()); - mPmtResponseTemp.resize(FV0DigParam::Instance().waveformNbins); - - for (int ir = 0; ir < int(mPmtResponseTemp.size() / NTimeBinsPerBC); ir++) { + for (int ir = 0; ir < int(mPmtResponseTemp.size() / mNTimeBinsPerBC); ir++) { auto bcCache = getBCCache(cachedIR[ir]); - for (int iBin = 0; iBin < NTimeBinsPerBC; iBin++) { - (*bcCache).mPmtChargeVsTime[detId][iBin] += (mPmtResponseTemp[ir * NTimeBinsPerBC + iBin] * mipFraction); + for (int iBin = 0; iBin < mNTimeBinsPerBC; iBin++) { + (*bcCache).mPmtChargeVsTime[detId][iBin] += (mPmtResponseTemp[ir * mNTimeBinsPerBC + iBin] * mipFraction); } added[ir] = true; } @@ -240,7 +235,7 @@ Int_t Digitizer::SimulateLightYield(Int_t pmt, Int_t nPhot) const Float_t Digitizer::IntegrateCharge(const ChannelBCDataF& pulse) const { int chargeIntMin = FV0DigParam::Instance().isIntegrateFull ? 0 : FV0DigParam::Instance().chargeIntBinMin; - int chargeIntMax = FV0DigParam::Instance().isIntegrateFull ? NTimeBinsPerBC : FV0DigParam::Instance().chargeIntBinMax; + int chargeIntMax = FV0DigParam::Instance().isIntegrateFull ? mNTimeBinsPerBC : FV0DigParam::Instance().chargeIntBinMax; Float_t totalCharge = 0.0f; for (int iTimeBin = chargeIntMin; iTimeBin < chargeIntMax; iTimeBin++) { @@ -263,7 +258,7 @@ Float_t Digitizer::SimulateTimeCfd(/*Int_t channel, */ const ChannelBCDataF& pul Int_t const binShift = TMath::Nint(FV0DigParam::Instance().timeShiftCfd / mBinSize); Float_t sigPrev = -pulse[0]; //[0]; - for (Int_t iTimeBin = 1; iTimeBin < NTimeBinsPerBC; ++iTimeBin) { + for (Int_t iTimeBin = 1; iTimeBin < mNTimeBinsPerBC; ++iTimeBin) { Float_t const sigCurrent = (iTimeBin >= binShift ? 5.0f * pulse[iTimeBin - binShift] - pulse[iTimeBin] : -pulse[iTimeBin]); if (sigPrev < 0.0f && sigCurrent >= 0.0f) { timeCfd = Float_t(iTimeBin) * mBinSize; From eb795b1dcab0d2dca4e62e1577aa831d70cce601 Mon Sep 17 00:00:00 2001 From: arvindkhuntia Date: Mon, 5 Oct 2020 13:43:29 +0200 Subject: [PATCH 3/3] clang format --- .../FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h | 2 +- Detectors/FIT/FV0/simulation/src/Digitizer.cxx | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h index 4697f4a804d46..0901ff3c1b0ee 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/FV0DigParam.h @@ -35,7 +35,7 @@ struct FV0DigParam : public o2::conf::ConfigurableParamHelper { //float timeShiftCfd = 3.3; // From the cosmic measurements of FV0 [keep it for reference] float timeShiftCfd = 5.3; // TODO: adjust after FV0 with FEE measurements are done float singleMipThreshold = 3.0; // in [MeV] of deposited energy - float singleHitTimeThreshold = 120.0;// in [ns] to skip very slow particles + float singleHitTimeThreshold = 120.0; // in [ns] to skip very slow particles float waveformNbins = 10000; // number of bins for the analog pulse waveform float waveformBinWidth = 0.01302; // number of bins for the analog float avgCfdTimeForMip = 8.63; // in ns to shift the CFD time to zero TODO do ring wise diff --git a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx index 3490c19ed5d0e..3ce0a80e69022 100644 --- a/Detectors/FIT/FV0/simulation/src/Digitizer.cxx +++ b/Detectors/FIT/FV0/simulation/src/Digitizer.cxx @@ -128,7 +128,7 @@ void Digitizer::process(const std::vector& hits, std::vector const& cachedIR, int nCachedIR, const int detId) { @@ -147,7 +147,7 @@ void Digitizer::createPulse(float mipFraction, int parID,const double hitTime, ///Time of flight subtracted from Hit time //TODO have different TOF according to thr ring number Size_t NBinShift = std::lround((hitTime - FV0DigParam::Instance().globalTimeOfFlight) / FV0DigParam::Instance().waveformBinWidth); mPmtResponseTemp.resize(FV0DigParam::Instance().waveformNbins, 0.); - std::memcpy(&mPmtResponseTemp[NBinShift], &mPmtResponseGlobal[0], sizeof(double)*(FV0DigParam::Instance().waveformNbins-NBinShift)); + std::memcpy(&mPmtResponseTemp[NBinShift], &mPmtResponseGlobal[0], sizeof(double) * (FV0DigParam::Instance().waveformNbins - NBinShift)); for (int ir = 0; ir < int(mPmtResponseTemp.size() / mNTimeBinsPerBC); ir++) { auto bcCache = getBCCache(cachedIR[ir]);