diff --git a/sbndcode/CRT/CRTReco/CRTStripHitProducer_module.cc b/sbndcode/CRT/CRTReco/CRTStripHitProducer_module.cc index 1f212abd1..4a2be4d7f 100644 --- a/sbndcode/CRT/CRTReco/CRTStripHitProducer_module.cc +++ b/sbndcode/CRT/CRTReco/CRTStripHitProducer_module.cc @@ -23,10 +23,9 @@ #include "sbnobj/SBND/CRT/FEBData.hh" #include "sbnobj/SBND/CRT/CRTStripHit.hh" -#include "sbnobj/SBND/Timing/DAQTimestamp.hh" +#include "sbnobj/SBND/Timing/FrameShiftInfo.hh" #include "sbndcode/Geometry/GeometryWrappers/CRTGeoService.h" -#include "sbndcode/Decoders/PTB/sbndptb.h" #include "sbndcode/Timing/SBNDRawTimingObj.h" #include "sbndcode/ChannelMaps/CRT/CRTChannelMapService.h" @@ -50,20 +49,23 @@ class sbnd::crt::CRTStripHitProducer : public art::EDProducer { void produce(art::Event& e) override; std::vector CreateStripHits(art::Ptr &data, const uint32_t ref_time_s, - const uint32_t ref_time_ns); + const uint32_t ref_time_ns, const uint32_t etrig_time_s, + const uint32_t etrig_time_ns); std::set UnixSet(const std::vector> &datas); - bool SPECTDCReference(art::Event& e, const uint64_t &raw_ts, uint64_t &ref_time); - bool PTBHLTReference(art::Event& e, const uint64_t &raw_ts, uint64_t &ref_time, uint32_t &hlt_code); private: art::ServiceHandle fCRTGeoService; art::ServiceHandle fCRTChannelMapService; - std::string fFEBDataModuleLabel; + art::InputTag fFEBDataModuleLabel; + art::InputTag fFrameShiftModuleLabel; uint16_t fADCThreshold; std::vector fErrorCoeff; bool fAllowFlag1; + bool fApplyETrigWindow; + double fETrigMin; + double fETrigMax; bool fApplyTs0Window; double fTs0Min; double fTs0Max; @@ -72,24 +74,19 @@ class sbnd::crt::CRTStripHitProducer : public art::EDProducer { double fTs1Max; bool fCorrectForDifferentSecond; bool fReferenceTs0; - int fTimingType; - std::string fDAQHeaderModuleLabel; - std::string fDAQHeaderInstanceLabel; - uint32_t fRawTSCorrection; - uint32_t fMaxAllowedRefTimeDiff; - std::string fSPECTDCModuleLabel; - uint32_t fSPECTDCETrigChannel; - std::string fPTBModuleLabel; - std::vector fAllowedPTBHLTs; }; sbnd::crt::CRTStripHitProducer::CRTStripHitProducer(fhicl::ParameterSet const& p) : EDProducer{p} - , fFEBDataModuleLabel(p.get("FEBDataModuleLabel")) + , fFEBDataModuleLabel(p.get("FEBDataModuleLabel")) + , fFrameShiftModuleLabel(p.get("FrameShiftModuleLabel")) , fADCThreshold(p.get("ADCThreshold")) , fErrorCoeff(p.get>("ErrorCoeff")) , fAllowFlag1(p.get("AllowFlag1")) + , fApplyETrigWindow(p.get("ApplyETrigWindow")) + , fETrigMin(p.get("ETrigMin", 0)) + , fETrigMax(p.get("ETrigMax", std::numeric_limits::max())) , fApplyTs0Window(p.get("ApplyTs0Window")) , fTs0Min(p.get("Ts0Min", 0)) , fTs0Max(p.get("Ts0Max", std::numeric_limits::max())) @@ -98,15 +95,6 @@ sbnd::crt::CRTStripHitProducer::CRTStripHitProducer(fhicl::ParameterSet const& p , fTs1Max(p.get("Ts1Max", std::numeric_limits::max())) , fCorrectForDifferentSecond(p.get("CorrectForDifferentSecond")) , fReferenceTs0(p.get("ReferenceTs0")) - , fTimingType(p.get("TimingType", 0)) - , fDAQHeaderModuleLabel(p.get("DAQHeaderModuleLabel", "")) - , fDAQHeaderInstanceLabel(p.get("DAQHeaderInstanceLabel", "")) - , fRawTSCorrection(p.get("RawTSCorrection", 0)) - , fMaxAllowedRefTimeDiff(p.get("MaxAllowedRefTimeDiff", 0)) - , fSPECTDCModuleLabel(p.get("SPECTDCModuleLabel", "")) - , fSPECTDCETrigChannel(p.get("SPECTDCETrigChannel", 4)) - , fPTBModuleLabel(p.get("PTBModuleLabel", "")) - , fAllowedPTBHLTs(p.get>("AllowedPTBHLTs", {})) { produces>(); produces>(); @@ -127,66 +115,37 @@ void sbnd::crt::CRTStripHitProducer::produce(art::Event& e) std::vector> FEBDataVec; art::fill_ptr_vector(FEBDataVec, FEBDataHandle); - uint64_t raw_ts = 0, ref_time = 0; - uint32_t ref_time_s = 0, ref_time_ns = 0; + uint32_t ref_time_s = 0, ref_time_ns = 0, etrig_time_s = 0, etrig_time_ns = 0; if(fReferenceTs0) { - art::Handle DAQHeaderHandle; - e.getByLabel(fDAQHeaderModuleLabel, fDAQHeaderInstanceLabel, DAQHeaderHandle); + art::Handle frameShiftHandle; + if (!e.getByLabel(fFrameShiftModuleLabel, frameShiftHandle)) + throw std::runtime_error("Frame Shift Info object is invalid, check data quality"); - if(DAQHeaderHandle.isValid()) - { - artdaq::RawEvent rawHeaderEvent = artdaq::RawEvent(*DAQHeaderHandle); - raw_ts = rawHeaderEvent.timestamp() - fRawTSCorrection; - } + uint64_t ref_time = frameShiftHandle->FrameDefault(); - int timingType = fTimingType; - int timingCh = 0; + ref_time_s = ref_time / sbnd::timing::kSecondInNanoseconds; + ref_time_ns = ref_time % sbnd::timing::kSecondInNanoseconds; - if(timingType == 0) - { - uint64_t spec_tdc_ref_time = 0; - - if(SPECTDCReference(e, raw_ts, spec_tdc_ref_time)) - { - ref_time = spec_tdc_ref_time; - timingCh = fSPECTDCETrigChannel; - } - else - ++timingType; - } + timingReferenceInfo->timingType = frameShiftHandle->TimingTypeDefault(); + timingReferenceInfo->timingChannel = frameShiftHandle->TimingChannelDefault(); - if(timingType == 1) - { - uint64_t ptb_hlt_ref_time = 0; - uint32_t hlt_code = 0; - - if(PTBHLTReference(e, raw_ts, ptb_hlt_ref_time, hlt_code)) - { - ref_time = ptb_hlt_ref_time; - timingCh = hlt_code; - } - else - ++timingType; - } + if(timingReferenceInfo->timingType == sbnd::timing::kNoShiftType) + { + std::set unix_set = UnixSet(FEBDataVec); + ref_time_s = unix_set.size() ? *unix_set.rbegin(): 0; + } - if(timingType == 2) - { - std::set unix_set = UnixSet(FEBDataVec); - ref_time = unix_set.size() ? *unix_set.rbegin() * static_cast(1e9) : 0; - } - - ref_time_s = ref_time / static_cast(1e9); - ref_time_ns = ref_time % static_cast(1e9); + uint64_t etrig_time = frameShiftHandle->FrameEtrig(); - timingReferenceInfo->timingType = timingType; - timingReferenceInfo->timingChannel = timingCh; + etrig_time_s = etrig_time / sbnd::timing::kSecondInNanoseconds; + etrig_time_ns = etrig_time % sbnd::timing::kSecondInNanoseconds; } for(auto data : FEBDataVec) { - std::vector newStripHits = CreateStripHits(data, ref_time_s, ref_time_ns); + std::vector newStripHits = CreateStripHits(data, ref_time_s, ref_time_ns, etrig_time_s, etrig_time_ns); for(auto hit : newStripHits) { @@ -203,7 +162,8 @@ void sbnd::crt::CRTStripHitProducer::produce(art::Event& e) } std::vector sbnd::crt::CRTStripHitProducer::CreateStripHits(art::Ptr &data, const uint32_t ref_time_s, - const uint32_t ref_time_ns) + const uint32_t ref_time_ns, const uint32_t etrig_time_s, + const uint32_t etrig_time_ns) { std::vector stripHits; @@ -249,6 +209,27 @@ std::vector sbnd::crt::CRTStripHitProducer::CreateStripH if(fReferenceTs0) t0 -= ref_time_ns; + if(fApplyETrigWindow) + { + double t0_etrig = (int)data->Ts0() + module.t0DelayCorrection; + const int64_t unix_diff = static_cast(ref_time_s) - static_cast(etrig_time_s); + + if(unix_diff < -1 || unix_diff > 1) + { + throw std::runtime_error("Unix timestamps differ by more than 1 (" + std::to_string(unix_diff) + ")"); + } + + t0_etrig -= etrig_time_ns; + + if(unix_diff == 1) + t0_etrig -= 1e9; + else if(unix_diff == -1) + t0_etrig += 1e9; + + if(t0_etrig < fETrigMin || t0_etrig > fETrigMax) + return stripHits; + } + if(fApplyTs0Window && (t0 < fTs0Min || t0 > fTs0Max)) return stripHits; @@ -309,92 +290,4 @@ std::set sbnd::crt::CRTStripHitProducer::UnixSet(const std::vector> TDCHandle; - e.getByLabel(fSPECTDCModuleLabel, TDCHandle); - - if(!TDCHandle.isValid() || TDCHandle->size() == 0) - return found; - - std::vector> TDCVec; - art::fill_ptr_vector(TDCVec, TDCHandle); - - uint64_t min_diff = std::numeric_limits::max(); - uint64_t min_diff_ts = 0; - - for(auto ts : TDCVec) - { - if(ts->Channel() == fSPECTDCETrigChannel) - { - uint64_t diff = raw_ts > ts->Timestamp() ? raw_ts - ts->Timestamp() : ts->Timestamp() - raw_ts; - - if(diff < min_diff) - { - min_diff = diff; - min_diff_ts = ts->Timestamp(); - found = true; - } - } - } - - if(min_diff > fMaxAllowedRefTimeDiff) - return false; - - ref_time = min_diff_ts; - - return found; -} - -bool sbnd::crt::CRTStripHitProducer::PTBHLTReference(art::Event& e, const uint64_t &raw_ts, uint64_t &ref_time, uint32_t &hlt_code) -{ - bool found = false; - - art::Handle> PTBHandle; - e.getByLabel(fPTBModuleLabel, PTBHandle); - - if(!PTBHandle.isValid() || PTBHandle->size() == 0) - return found; - - std::vector> PTBVec; - art::fill_ptr_vector(PTBVec, PTBHandle); - - uint64_t min_diff = std::numeric_limits::max(); - uint64_t min_diff_ts = 0; - - for(auto ptb : PTBVec) - { - for(auto hlt : ptb->GetHLTriggers()) - { - uint64_t hlt_timestamp = (hlt.timestamp * 20); - std::bitset<32> hlt_word_bitset = std::bitset<32>(hlt.trigger_word); - - for(uint32_t allowed_hlt : fAllowedPTBHLTs) - { - if(hlt_word_bitset[allowed_hlt]) - { - uint64_t diff = raw_ts > hlt_timestamp ? raw_ts - hlt_timestamp : hlt_timestamp - raw_ts; - - if(diff < min_diff) - { - min_diff = diff; - min_diff_ts = hlt_timestamp; - hlt_code = allowed_hlt; - found = true; - } - } - } - } - } - - if(min_diff > fMaxAllowedRefTimeDiff) - return false; - - ref_time = min_diff_ts; - - return found; -} - DEFINE_ART_MODULE(sbnd::crt::CRTStripHitProducer) diff --git a/sbndcode/CRT/CRTReco/crtrecoproducers_sbnd.fcl b/sbndcode/CRT/CRTReco/crtrecoproducers_sbnd.fcl index 9aa8e2a26..0093419b9 100644 --- a/sbndcode/CRT/CRTReco/crtrecoproducers_sbnd.fcl +++ b/sbndcode/CRT/CRTReco/crtrecoproducers_sbnd.fcl @@ -5,6 +5,7 @@ BEGIN_PROLOG crtstriphitproducer_sbnd: { FEBDataModuleLabel: "crtsim" + FrameShiftModuleLabel: "" ADCThreshold: 60 ErrorCoeff: [ 0.26, -0.27, 0.025 ] AllowFlag1: false @@ -17,20 +18,13 @@ crtstriphitproducer_sbnd: crtstriphitproducer_data_sbnd: @local::crtstriphitproducer_sbnd crtstriphitproducer_data_sbnd.FEBDataModuleLabel: "crtdecoder" +crtstriphitproducer_data_sbnd.FrameShiftModuleLabel: "frameshift" crtstriphitproducer_data_sbnd.AllowFlag1: true crtstriphitproducer_data_sbnd.CorrectForDifferentSecond: true crtstriphitproducer_data_sbnd.ReferenceTs0: true -crtstriphitproducer_data_sbnd.TimingType: 0 -crtstriphitproducer_data_sbnd.DAQHeaderModuleLabel: "daq" -crtstriphitproducer_data_sbnd.DAQHeaderInstanceLabel: "RawEventHeader" -crtstriphitproducer_data_sbnd.RawTSCorrection: 367000 -crtstriphitproducer_data_sbnd.MaxAllowedRefTimeDiff: 3000000 -crtstriphitproducer_data_sbnd.SPECTDCModuleLabel: "tdcdecoder" -crtstriphitproducer_data_sbnd.PTBModuleLabel: "ptbdecoder" -crtstriphitproducer_data_sbnd.AllowedPTBHLTs: [ 1, 2, 3, 4, 5, 14, 15 ] -crtstriphitproducer_data_sbnd.ApplyTs0Window: true -crtstriphitproducer_data_sbnd.Ts0Min: -1.71e6 -crtstriphitproducer_data_sbnd.Ts0Max: 1.51e6 +crtstriphitproducer_data_sbnd.ApplyETrigWindow: true +crtstriphitproducer_data_sbnd.ETrigMin: -1.71e6 +crtstriphitproducer_data_sbnd.ETrigMax: 1.51e6 crtclusterproducer_sbnd: { diff --git a/sbndcode/CRT/CRTVeto/CRTVetoProducer_module.cc b/sbndcode/CRT/CRTVeto/CRTVetoProducer_module.cc index 1abe23046..a362e5e1b 100644 --- a/sbndcode/CRT/CRTVeto/CRTVetoProducer_module.cc +++ b/sbndcode/CRT/CRTVeto/CRTVetoProducer_module.cc @@ -59,8 +59,6 @@ class sbnd::crt::CRTVetoProducer : public art::EDProducer { void produce(art::Event& e) override; - bool SPECTDCReference(art::Event& e, const uint64_t &raw_ts, uint64_t &ref_time, const uint32_t channel); - // Functions to check if a point is in a specific Tagger bool inSouth(const CRTTagger t); bool inNorth(const CRTTagger t); @@ -76,16 +74,9 @@ class sbnd::crt::CRTVetoProducer : public art::EDProducer { std::string fCRTClusterModuleLabel; std::string fCRTSpacePointModuleLabel; std::string fCRTTimingReferenceInfoLabel; - std::string fSPECTDCModuleLabel; bool fIsData; bool fDebug; - // raw timestamp correction variables - std::string fDAQHeaderModuleLabel; - std::string fDAQHeaderInstanceLabel; - uint32_t fRawTSCorrection; - uint32_t fMaxAllowedRefTimeDiff; - unsigned int fEventID; unsigned int fRun; unsigned int fSubRun; @@ -100,13 +91,8 @@ sbnd::crt::CRTVetoProducer::CRTVetoProducer(fhicl::ParameterSet const& p) , fCRTClusterModuleLabel(p.get("CRTClusterModuleLabel")) , fCRTSpacePointModuleLabel(p.get("CRTSpacePointModuleLabel")) , fCRTTimingReferenceInfoLabel(p.get("CRTTimingReferenceInfoLabel")) - , fSPECTDCModuleLabel(p.get("SPECTDCModuleLabel", "")) , fIsData(p.get("IsData")) , fDebug(p.get("Debug", false)) - , fDAQHeaderModuleLabel(p.get("DAQHeaderModuleLabel", "")) - , fDAQHeaderInstanceLabel(p.get("DAQHeaderInstanceLabel", "")) - , fRawTSCorrection(p.get("RawTSCorrection", 0)) - , fMaxAllowedRefTimeDiff(p.get("MaxAllowedRefTimeDiff", 0)) { produces>(); produces>(); @@ -127,40 +113,20 @@ void sbnd::crt::CRTVetoProducer::produce(art::Event& e) std::cout << std::endl; } - uint64_t raw_ts = 0; - int64_t ref_time = 0; - int64_t ref_time_ns = 0; - - // Data needs to be handled differently than MC - + // Data needs to be handled differently than MC if (fIsData) { - - // First Need to calculate the raw_ts in the same way as the CRTStripProducer - art::Handle DAQHeaderHandle; - e.getByLabel(fDAQHeaderModuleLabel, fDAQHeaderInstanceLabel, DAQHeaderHandle); - - if(DAQHeaderHandle.isValid()) - { - if (fDebug) - std::cout << "DAQHeader is valid --> Calculate raw ts" << std::endl; - artdaq::RawEvent rawHeaderEvent = artdaq::RawEvent(*DAQHeaderHandle); - raw_ts = rawHeaderEvent.timestamp() - fRawTSCorrection; - } art::Handle TimingReferenceHandle; e.getByLabel(fCRTTimingReferenceInfoLabel, TimingReferenceHandle); if(TimingReferenceHandle.isValid()) { uint16_t TType = TimingReferenceHandle->timingType; - uint16_t TChannel = TimingReferenceHandle->timingChannel; - - // the CRTStripProducer should reference the Event Trigger - // If we don't find this --> return false - if (TType != 0 || TChannel != 4) { + // Check if the timing is valid + if (TType != 0 && TType != 1) { if (fDebug) { std::cout << std::endl; - std::cout << "ETrig was not referenced for this event --> Don't Flag Event!" << std::endl; + std::cout << "TType is not valid for this event --> Don't Flag Event!" << std::endl; std::cout << std::endl; } auto crtveto = std::make_unique(false, false, false, false, false); @@ -175,29 +141,8 @@ void sbnd::crt::CRTVetoProducer::produce(art::Event& e) } - // The Strip Hits were referenced properly --> Change to RWM - uint64_t spec_tdc_ref_time_etrig = 0; - uint64_t spec_tdc_ref_time_rwm = 0; - if(SPECTDCReference(e, raw_ts, spec_tdc_ref_time_etrig, 4) && SPECTDCReference(e, raw_ts, spec_tdc_ref_time_rwm, 2)) { - - ref_time += spec_tdc_ref_time_rwm; - ref_time -= spec_tdc_ref_time_etrig; - - } - else { - if (fDebug) { - std::cout << std::endl; - std::cout << std::endl; - std::cout << "Was NOT able to reference time to the RWM !!!" << std::endl; - std::cout << std::endl; - std::cout << std::endl; - } - } - ref_time_ns = ref_time % static_cast(1e9); - } // end of check for Data - // Get the Clusters for this event art::Handle> CRTClusterHandle; e.getByLabel(fCRTClusterModuleLabel, CRTClusterHandle); @@ -236,7 +181,6 @@ void sbnd::crt::CRTVetoProducer::produce(art::Event& e) // Get the time and tagger for this space point double t = sp->Ts0(); // Starts in ns - t -= ref_time_ns; t /= 1000; // convert to us @@ -328,46 +272,6 @@ void sbnd::crt::CRTVetoProducer::produce(art::Event& e) } // end of produce -bool sbnd::crt::CRTVetoProducer::SPECTDCReference(art::Event& e, const uint64_t &raw_ts, uint64_t &ref_time, const uint32_t channel) -{ - bool found = false; - - art::Handle> TDCHandle; - e.getByLabel(fSPECTDCModuleLabel, TDCHandle); - - if(!TDCHandle.isValid() || TDCHandle->size() == 0) - return found; - - std::vector> TDCVec; - art::fill_ptr_vector(TDCVec, TDCHandle); - - uint64_t min_diff = std::numeric_limits::max(); - uint64_t min_diff_ts = 0; - - for(auto ts : TDCVec) - { - // For the CRT Veto, we want to use the RWM --> Channel 2 - if(ts->Channel() == channel) - { - uint64_t diff = raw_ts > ts->Timestamp() ? raw_ts - ts->Timestamp() : ts->Timestamp() - raw_ts; - - if(diff < min_diff) - { - min_diff = diff; - min_diff_ts = ts->Timestamp(); - found = true; - } - } - } - - if(min_diff > fMaxAllowedRefTimeDiff) - return false; - - ref_time = min_diff_ts; - - return found; -} - bool sbnd::crt::CRTVetoProducer::inSouth(const CRTTagger t) { if (t == kSouthTagger) { diff --git a/sbndcode/CRT/CRTVeto/crtvetoproducer_sbnd.fcl b/sbndcode/CRT/CRTVeto/crtvetoproducer_sbnd.fcl index 7be12fb67..ef3d16a40 100644 --- a/sbndcode/CRT/CRTVeto/crtvetoproducer_sbnd.fcl +++ b/sbndcode/CRT/CRTVeto/crtvetoproducer_sbnd.fcl @@ -15,25 +15,28 @@ crtvetoproducer_data_sbnd: @local::crtvetoproducer_sbn crtvetoproducer_data_sbnd.IsData: true # Values for testing -crtvetoproducer_data_sbnd.WindowStart: -1.185 -crtvetoproducer_data_sbnd.WindowEnd: 0.515 -crtvetoproducer_data_sbnd.SPECTDCModuleLabel: "tdcdecoder" +#crtvetoproducer_data_sbnd.WindowStart: -1.185 +#crtvetoproducer_data_sbnd.WindowEnd: 0.515 +crtvetoproducer_data_sbnd.WindowStart: 0.605 +crtvetoproducer_data_sbnd.WindowEnd: 2.305 +#crtvetoproducer_data_sbnd.SPECTDCModuleLabel: "tdcdecoder" -crtvetoproducer_data_sbnd.DAQHeaderModuleLabel: "daq" -crtvetoproducer_data_sbnd.DAQHeaderInstanceLabel: "RawEventHeader" -crtvetoproducer_data_sbnd.RawTSCorrection: 367000 -crtvetoproducer_data_sbnd.MaxAllowedRefTimeDiff: 3000000 +#crtvetoproducer_data_sbnd.DAQHeaderModuleLabel: "daq" +#crtvetoproducer_data_sbnd.DAQHeaderInstanceLabel: "RawEventHeader" +crtvetoproducer_data_sbnd.FrameShiftModuleLabel: "frameshift" +#crtvetoproducer_data_sbnd.RawTSCorrection: 367000 +#crtvetoproducer_data_sbnd.MaxAllowedRefTimeDiff: 3000000 # Top Hat Configuration for Debugging crtvetoproducer_data_sbnd_tophat: @local::crtvetoproducer_sbnd crtvetoproducer_data_sbnd_tophat.IsData: true # Values for testing -crtvetoproducer_data_sbnd_tophat.SPECTDCModuleLabel: "tdcdecoder" -crtvetoproducer_data_sbnd_tophat.DAQHeaderModuleLabel: "daq" -crtvetoproducer_data_sbnd_tophat.DAQHeaderInstanceLabel: "RawEventHeader" -crtvetoproducer_data_sbnd_tophat.RawTSCorrection: 367000 -crtvetoproducer_data_sbnd_tophat.MaxAllowedRefTimeDiff: 3000000 +#crtvetoproducer_data_sbnd_tophat.SPECTDCModuleLabel: "tdcdecoder" +#crtvetoproducer_data_sbnd_tophat.DAQHeaderModuleLabel: "daq" +#crtvetoproducer_data_sbnd_tophat.DAQHeaderInstanceLabel: "RawEventHeader" +#crtvetoproducer_data_sbnd_tophat.RawTSCorrection: 367000 +#crtvetoproducer_data_sbnd_tophat.MaxAllowedRefTimeDiff: 3000000 crtvetoproducer_data_sbnd_tophat.WindowStart: -2 crtvetoproducer_data_sbnd_tophat.WindowEnd: 2 diff --git a/sbndcode/Decoders/CRT/run_crt_ptb_decoder_sbnd.fcl b/sbndcode/Decoders/CRT/run_crt_ptb_decoder_sbnd.fcl deleted file mode 100644 index 7cb1b72b6..000000000 --- a/sbndcode/Decoders/CRT/run_crt_ptb_decoder_sbnd.fcl +++ /dev/null @@ -1,42 +0,0 @@ -#include "crt_decoder_sbnd.fcl" -#include "services_sbnd.fcl" -#include "rootoutput_sbnd.fcl" -#include "crt_channel_map_service_sbnd.fcl" -#include "SBNDPTBDecoderDefaults.fcl" - -process_name: Decoders - -services: -{ - @table::sbnd_services - CRTChannelMapService: @local::crt_channel_map_sbnd -} - -source: -{ - module_type: RootInput -} - -outputs: -{ - out1: - { - @table::sbnd_rootoutput - dataTier: "reconstructed" - } -} - -physics: -{ - producers: - { - crtdecoder: @local::crt_decoder - ptbdecoder: @local::SBNDPTBDecoderDefaults - } - - decoders: [ crtdecoder, ptbdecoder ] - trigger_paths: [ decoders ] - - stream1: [ out1 ] - end_paths: [ stream1 ] -} diff --git a/sbndcode/Decoders/CRT/run_crt_ptb_tdc_decoder_sbnd.fcl b/sbndcode/Decoders/CRT/run_crt_ptb_tdc_decoder_sbnd.fcl index ca58e58b8..8c4f82ba9 100644 --- a/sbndcode/Decoders/CRT/run_crt_ptb_tdc_decoder_sbnd.fcl +++ b/sbndcode/Decoders/CRT/run_crt_ptb_tdc_decoder_sbnd.fcl @@ -4,6 +4,7 @@ #include "crt_channel_map_service_sbnd.fcl" #include "SBNDPTBDecoderDefaults.fcl" #include "spectdc_decoder_sbnd.fcl" +#include "frameshift_sbnd_data.fcl" process_name: Decoders @@ -34,9 +35,10 @@ physics: crtdecoder: @local::crt_decoder ptbdecoder: @local::SBNDPTBDecoderDefaults tdcdecoder: @local::spec_tdc_decoder_sbnd + frameshift: @local::frameshift_data } - decoders: [ crtdecoder, ptbdecoder, tdcdecoder ] + decoders: [ crtdecoder, ptbdecoder, tdcdecoder, frameshift ] trigger_paths: [ decoders ] stream1: [ out1 ] diff --git a/sbndcode/Decoders/CRT/run_crt_ptb_tdc_pmt_decoder_sbnd.fcl b/sbndcode/Decoders/CRT/run_crt_ptb_tdc_pmt_decoder_sbnd.fcl index 6bd07c8b2..1a9e747df 100644 --- a/sbndcode/Decoders/CRT/run_crt_ptb_tdc_pmt_decoder_sbnd.fcl +++ b/sbndcode/Decoders/CRT/run_crt_ptb_tdc_pmt_decoder_sbnd.fcl @@ -5,6 +5,7 @@ #include "SBNDPTBDecoderDefaults.fcl" #include "spectdc_decoder_sbnd.fcl" #include "pmtdecoder.fcl" +#include "frameshift_sbnd_data.fcl" process_name: Decoders @@ -37,9 +38,10 @@ physics: ptbdecoder: @local::SBNDPTBDecoderDefaults tdcdecoder: @local::spec_tdc_decoder_sbnd pmtdecoder: @local::pmtdecoder + frameshift: @local::frameshift_data } - decoders: [ crtdecoder, ptbdecoder, tdcdecoder, pmtdecoder ] + decoders: [ crtdecoder, ptbdecoder, tdcdecoder, frameshift, pmtdecoder ] trigger_paths: [ decoders ] stream1: [ out1 ] diff --git a/sbndcode/Decoders/PMT/CMakeLists.txt b/sbndcode/Decoders/PMT/CMakeLists.txt index 02ea063b9..0b5ec3709 100644 --- a/sbndcode/Decoders/PMT/CMakeLists.txt +++ b/sbndcode/Decoders/PMT/CMakeLists.txt @@ -3,7 +3,7 @@ cet_build_plugin( SBNDPMTDecoder art::module LIBRARIES sbndaq_artdaq_core::sbndaq-artdaq-core_Overlays_SBND sbnobj::SBND_Timing - sbndcode_OpDetSim + sbndcode_OpDetSim lardata::Utilities canvas::canvas lardataobj::RecoBase diff --git a/sbndcode/Decoders/PMT/SBNDPMTDecoder_module.cc b/sbndcode/Decoders/PMT/SBNDPMTDecoder_module.cc index 18a56d8d5..c31c20d7a 100644 --- a/sbndcode/Decoders/PMT/SBNDPMTDecoder_module.cc +++ b/sbndcode/Decoders/PMT/SBNDPMTDecoder_module.cc @@ -29,9 +29,13 @@ #include "sbndaq-artdaq-core/Overlays/Common/CAENV1730Fragment.hh" #include "sbndaq-artdaq-core/Overlays/SBND/PTBFragment.hh" #include "sbndaq-artdaq-core/Overlays/FragmentType.hh" -#include "sbnobj/SBND/Timing/DAQTimestamp.hh" +#include "sbnobj/SBND/Timing/FrameShiftInfo.hh" #include "sbndcode/Decoders/PTB/sbndptb.h" #include "sbndcode/Timing/SBNDRawTimingObj.h" +#include "sbndcode/Calibration/PDSDatabaseInterface/PMTCalibrationDatabase.h" +#include "sbndcode/Calibration/PDSDatabaseInterface/IPMTCalibrationDatabaseService.h" +#include "sbndcode/OpDetSim/TriggerEmulationService.h" +#include "sbndcode/OpDetSim/sbndPDMapAlg.hh" #include "art_root_io/TFileService.h" #include "TH1D.h" @@ -48,8 +52,7 @@ #include #include #include - -#include "sbndcode/OpDetSim/TriggerEmulationService.h" +#include namespace sbndaq { class SBNDPMTDecoder; @@ -78,27 +81,17 @@ class sbndaq::SBNDPMTDecoder : public art::EDProducer { uint32_t get_boardid(artdaq::Fragment & frag); void get_timing(artdaq::Fragment & frag, uint16_t & postpercent, uint32_t & ttt, uint32_t & len, int & tick); -private: - uint fdebug; - uint ftiming_type; + std::vector make_chmap(sbndDB::PMTCalibrationDatabase const* pmt_calib_db) const; - std::vector fcaen_fragment_name; - std::string fcaen_module_label; +private: + uint8_t fdebug; - std::vector fignore_fragid; + std::vector fcaen_fragment_name; + art::InputTag fframeshift_module_label; + + std::vector fignore_fragid; uint32_t fnominal_length; - - uint32_t fraw_ts_correction; - - std::string fspectdc_product_name; - uint32_t fspectdc_ftrig_ch; - uint32_t fspectdc_etrig_ch; - - std::string fptb_product_name; - std::vector fptb_allowed_hlts; - uint fptb_raw_diff_max; - - uint fallowed_time_diff; + uint32_t fallowed_time_diff; std::string fpmt_instance_name; std::string fflt_instance_name; @@ -113,51 +106,42 @@ class sbndaq::SBNDPMTDecoder : public art::EDProducer { std::vector fignore_timing_ch; int fn_maxflashes; - uint fn_caenboards; + uint8_t fn_caenboards; + uint8_t fn_caenpmtchannels; + uint16_t ftiming_caen_offset; uint16_t fthreshold_ftrig; uint16_t fdefault_postpercent; // should be a number between 0 and 100 - uint ffragid_offset; - uint fhist_evt; + uint16_t ffragid_offset; + uint32_t fhist_evt; - std::vector fset_fragid_map; + std::vector fset_fragid_map; bool fuse_set_map; - - std::vector fch_map; - int fmon_threshold; // histogram info std::stringstream histname; //raw waveform hist name art::ServiceHandle tfs; - uint evt_counter = 0; + uint32_t evt_counter = 0; + + sbndDB::PMTCalibrationDatabase const* fpmt_calib_db; + opdet::sbndPDMapAlg opdetmap; //map for photon detector types }; sbndaq::SBNDPMTDecoder::SBNDPMTDecoder(fhicl::ParameterSet const& p) : EDProducer{p} // , { - fdebug = p.get("debug",0); + fdebug = p.get("debug",0); - fcaen_fragment_name = p.get>("caen_fragment_name"); - fcaen_module_label = p.get("caen_module_label","daq"); + fcaen_fragment_name = p.get>("caen_fragment_name"); + fframeshift_module_label = p.get("frameshift_module_label"); - fignore_fragid = p.get>("ignore_fragid",{}); + fignore_fragid = p.get>("ignore_fragid",{}); fnominal_length = p.get("nominal_length",5000); + fallowed_time_diff = p.get("allowed_time_diff",3000); // us!!! - ftiming_type = p.get("timing_type",0); - - fraw_ts_correction = p.get("raw_ts_correction",367000); // ns - - fspectdc_product_name = p.get("spectdc_product_name","tdcdecoder"); - fspectdc_ftrig_ch = p.get("spectdc_ftrig_ch",3); - fspectdc_etrig_ch = p.get("spectdc_etrig_ch",4); - - fptb_product_name = p.get("ptb_product_name","ptbdecoder"); - fptb_allowed_hlts = p.get>("ptb_allowed_hlts",{}); - fptb_raw_diff_max = p.get("ptb_raw_diff_max",3000000); // ns - - fallowed_time_diff = p.get("allowed_time_diff",3000); // us!!! + // configure output labels fpmt_instance_name = p.get("pmtInstanceName","PMTChannels"); fflt_instance_name = p.get("ftrigInstanceName","FTrigChannels"); ftim_instance_name = p.get("timingInstanceName","TimingChannels"); @@ -170,17 +154,18 @@ sbndaq::SBNDPMTDecoder::SBNDPMTDecoder(fhicl::ParameterSet const& p) fignore_timing_ch = p.get>("ignore_timing_ch",{}); fn_maxflashes = p.get("n_maxflashes",30); - fn_caenboards = p.get("n_caenboards",8); + fn_caenboards = p.get("n_caenboards",8); + fn_caenpmtchannels = p.get("n_caenpmtchannels",15); + ftiming_caen_offset = p.get("timing_caen_offset",900); fthreshold_ftrig = p.get("threshold_ftrig",16350); fdefault_postpercent = p.get("default_postpercent",80); - ffragid_offset = p.get("fragid_offset",40960); - fhist_evt = p.get("hist_evt",1); + ffragid_offset = p.get("fragid_offset",40960); + fhist_evt = p.get("hist_evt",1); - fset_fragid_map = p.get>("set_fragid_map",{}); + fset_fragid_map = p.get>("set_fragid_map",{}); fuse_set_map = p.get("use_set_map",false); - fch_map = p.get>("ch_map",{}); - + fpmt_calib_db = lar::providerFrom(); fmon_threshold = p.get("mon_threshold", 15); produces< std::vector< raw::OpDetWaveform > >(fpmt_instance_name); @@ -223,13 +208,15 @@ void sbndaq::SBNDPMTDecoder::produce(art::Event& evt) evt_counter++; std::vector> board_frag_v(fn_caenboards); + std::vector ch_map = make_chmap(fpmt_calib_db); + uint ncont = 0; // counter for number of containers bool found_caen = false; // get CAEN fragments - for (const std::string &caen_name : fcaen_fragment_name){ + for (const art::InputTag &caen_name : fcaen_fragment_name){ art::Handle> fragmentHandle; - evt.getByLabel(fcaen_module_label, caen_name, fragmentHandle); + evt.getByLabel(caen_name, fragmentHandle); if (!fragmentHandle.isValid() || fragmentHandle->size() == 0){ if (fdebug>0) std::cout << "No CAEN V1730 fragments with label " << caen_name << " found." << std::endl; @@ -285,144 +272,20 @@ void sbndaq::SBNDPMTDecoder::produce(art::Event& evt) return; } - // create a timing type per event so the default doesn't get overwritten uint64_t event_trigger_time = 0; // in ns - auto timing_type = ftiming_type; - int timing_ch = -1; - - // get the raw event header timestamp - art::Handle header_handle; - uint64_t raw_timestamp = 0; - evt.getByLabel("daq", "RawEventHeader", header_handle); - if ((header_handle.isValid())){ - auto rawheader = artdaq::RawEvent(*header_handle); - raw_timestamp = rawheader.timestamp() - fraw_ts_correction; // includes sec + ns portion - if (fdebug>1) - std::cout << "Raw timestamp (w/ correction) -> " << "ts (ns): " << raw_timestamp % uint64_t(1e9) << ", sec (s): " << raw_timestamp / uint64_t(1e9) << std::endl; - } + + art::Handle frameShiftHandle; + if (!evt.getByLabel(fframeshift_module_label, frameShiftHandle)) + throw std::runtime_error("Frame Shift Info object is invalid, check data quality"); - // get spec tdc product - if (timing_type==0){ - art::Handle> tdcHandle; - evt.getByLabel(fspectdc_product_name,tdcHandle); - bool found_ett = false; - std::vector tdc_etrig_v; - uint64_t min_raw_tdc_diff = uint64_t(1e12); - - if (!tdcHandle.isValid() || tdcHandle->size() == 0){ - if (fdebug>0) std::cout << "No SPECTDC products found." << std::endl; - timing_type++; - } - else{ - if (fdebug>1) std::cout << "SPECTDC (decoded) products found: " << std::endl; - const std::vector tdc_v(*tdcHandle); - - for (size_t i=0; i1)){ - std::cout << " TDC CH " << ch << " -> " - << "name: " << name - << ", ts (ns): " << ts%uint64_t(1e9) - << ", sec (s): " << ts/uint64_t(1e9) - << ", offset: " << offset - << std::endl; - } - if (ch==fspectdc_etrig_ch){ - found_ett = true; - tdc_etrig_v.push_back(ts); - timing_ch = fspectdc_etrig_ch; - } - } - if (found_ett==false) - timing_type++; - else{ - if (tdc_etrig_v.size()==1){ - event_trigger_time = tdc_etrig_v.front()%uint64_t(1e9); - } - else{ // finding the closest ETRIG to the raw timestamp - for (size_t i=0; i < tdc_etrig_v.size(); i++){ - event_trigger_time = tdc_etrig_v.front()%uint64_t(1e9); - auto tdc_etrig = tdc_etrig_v[i]; - uint64_t diff = (tdc_etrig < (raw_timestamp)) ? (raw_timestamp - tdc_etrig) : (tdc_etrig - raw_timestamp); - if (diff < min_raw_tdc_diff){ - event_trigger_time = tdc_etrig%uint64_t(1e9); - min_raw_tdc_diff = diff; - } - } - } - } - } // end else statement for finding tdc - } - if (timing_type==1){ - // get ptb product - art::Handle> ptbHandle; - evt.getByLabel(fptb_product_name,ptbHandle); - bool found_ptb_hlt = false; - if ((!ptbHandle.isValid() || ptbHandle->size() == 0)){ - if (fdebug>0) std::cout << "No PTB products found." << std::endl; - timing_type++; - } - else{ - const std::vector ptb_v(*ptbHandle); - - if (fdebug>1) std::cout << "PTB (decoded) HLTs: " << std::endl; - for (size_t i=0; i ptb_hlt_bitset = std::bitset<32>(trig.trigger_word); - for (auto allowed_hlt : fptb_allowed_hlts){ - if (ptb_hlt_bitset[allowed_hlt]){ - auto diff = (ptb_timestamp < raw_timestamp) ? (raw_timestamp - ptb_timestamp) : (ptb_timestamp - raw_timestamp); - if (diff > fptb_raw_diff_max) continue; - - found_ptb_hlt = true; - event_trigger_time = ptb_timestamp%(uint(1e9)); - timing_ch = allowed_hlt; - if (fdebug>1){ - std::cout << " PTB HLT " << allowed_hlt << " -> " - << "ts (ns): " << ptb_timestamp%(uint(1e9)) - << std::endl; - } - } - } - } - } // end hlt check - } // end of loop over ptb products - if (found_ptb_hlt==false){ - timing_type++; - } - } // end handle validity check - } // end of timing type 1 - if (timing_type==2){ - // CAEN only configuration - // if neither PTB or SPEC products are found, the timestamp will be equal to - // the start of the waveform (according to CAEN TTT and wvfm length) - event_trigger_time = 0; - timing_ch=16; - } + event_trigger_time = frameShiftHandle->FrameDefault()%1'000'000'000; + evtTimingInfo->timingType = frameShiftHandle->TimingTypeDefault(); + evtTimingInfo->timingChannel = frameShiftHandle->TimingChannelDefault(); - if (fdebug>1){ - std::cout << "Using " ; - if (timing_type==0) - std::cout << "SPECTDC for event reference time, "; - else if (timing_type==1) - std::cout << "PTB for event reference time, "; - else if (timing_type==2) - std::cout << "CAEN-only for event reference time, "; - std::cout << "channel/word " << timing_ch << ", event trigger time: " << event_trigger_time << " ns" << std::endl; - } - evtTimingInfo->timingType = timing_type; - evtTimingInfo->timingChannel = timing_ch; + // if frameshift has no shift, the timestamp will be equal to + // the start of the waveform (according to CAEN TTT and wvfm length) + if(evtTimingInfo->timingType == sbnd::timing::kNoShiftType) + event_trigger_time = 0; bool extended_flag = false; @@ -561,7 +424,7 @@ void sbndaq::SBNDPMTDecoder::produce(art::Event& evt) wvfmHist->SetBinContent(n + 1, (double)combined_wvfm[n]); } double time_diff = (int(iwvfm_start) - int(event_trigger_time))*1e-3; // us - if ((std::abs(time_diff) > fallowed_time_diff) && (timing_type<2)){ + if ((std::abs(time_diff) > fallowed_time_diff) && (event_trigger_time!=0)){ // second rollover between reference time and waveform start if (std::abs(time_diff + 1e6) < fallowed_time_diff) time_diff += 1e6; // us @@ -571,26 +434,23 @@ void sbndaq::SBNDPMTDecoder::produce(art::Event& evt) else if (fdebug>1) std::cout << "WARNING: TIME DIFFERENCE IS GREATER THAN " << fallowed_time_diff << " us. Event timestamp: " << event_trigger_time << ", waveform timestamp: " << iwvfm_start << std::endl; } - uint ch; - if (i == 15) - ch = fragid; - else - ch = fch_map.at(fragid*15 + i); - raw::OpDetWaveform waveform(time_diff, ch, combined_wvfm); if ((i == 15) && (foutput_ftrig_wvfm)){ - fltwvfmVec->push_back(waveform); + // for the timing ch of each board, the opdetwaveform channel is equal to the board id + fltwvfmVec->emplace_back(time_diff, fragid, combined_wvfm); art::Ptr wvfmPtr = make_fltwvfm_ptr(fltwvfmVec->size()-1); fltTimingAssns->addSingle(brdTimingInfoPtr, wvfmPtr); } else if ((fragid == 8)){ // fyi: this hardcodes the timing caen board index + // for the timing caen, the opdetwaveform channel is offset (900) + chidx if ((foutput_timing_wvfm) && (std::find(fignore_timing_ch.begin(), fignore_timing_ch.end(), i) == fignore_timing_ch.end())){ - timwvfmVec->push_back(waveform); + timwvfmVec->emplace_back(time_diff, ftiming_caen_offset + i, combined_wvfm); art::Ptr wvfmPtr = make_timwvfm_ptr(timwvfmVec->size()-1); timTimingAssns->addSingle(brdTimingInfoPtr, wvfmPtr); } } else{ - pmtwvfmVec->push_back(waveform); + // for normal pmt chs, the opdetwaveform channel is from the ch_map + pmtwvfmVec->emplace_back(time_diff, ch_map.at(fragid*15+i), combined_wvfm); art::Ptr wvfmPtr = make_pmtwvfm_ptr(pmtwvfmVec->size()-1); pmtTimingAssns->addSingle(brdTimingInfoPtr, wvfmPtr); } @@ -615,14 +475,12 @@ void sbndaq::SBNDPMTDecoder::produce(art::Event& evt) } } } // end board loop + board_frag_v.clear(); // loop through flashes art::ServiceHandle tfs; art::ServiceHandle fTriggerService; int PMTPerBoard = fTriggerService->getPMTPerBoard(); - // int fTotalCAENBoards = fTriggerService->getTotalCAENBoards(); - //std::vector< std::vector > MonPulsesAll; - //MonPulseAll.clear(); std::vector MonPulsesFlat; std::vector pulseSizes; MonPulsesFlat.clear(); @@ -812,4 +670,20 @@ uint32_t sbndaq::SBNDPMTDecoder::get_boardid(artdaq::Fragment & frag){ return boardid; } +std::vector sbndaq::SBNDPMTDecoder::make_chmap(sbndDB::PMTCalibrationDatabase const* pmt_calib_db) const { + std::vector ch_map(fn_caenboards*fn_caenpmtchannels,9999); + if (pmt_calib_db==nullptr){ + throw std::runtime_error("PMT Calibration Database pointer is null."); + } + auto nopdets = opdetmap.size(); + for (size_t idet=0; idetgetCAENDigitizer(idet) - 1; + int channel = pmt_calib_db->getCAENDigitizerChannel(idet); + ch_map.at(board*15 + channel) = idet; + } + return ch_map; +} + DEFINE_ART_MODULE(sbndaq::SBNDPMTDecoder) diff --git a/sbndcode/Decoders/PMT/pmtdecoder.fcl b/sbndcode/Decoders/PMT/pmtdecoder.fcl index ea2f7cc1b..e7fd09d65 100644 --- a/sbndcode/Decoders/PMT/pmtdecoder.fcl +++ b/sbndcode/Decoders/PMT/pmtdecoder.fcl @@ -5,31 +5,16 @@ pmtdecoder: module_type: "SBNDPMTDecoder" debug: 0 # 0=no output, 1=standard, 2=verbose - # output timing system configuration - ## priority system - ## SPEC TDC ETT [0] -> PTB ETT [1] -> CAEN-only [2] - timing_type: 0 # default is SPEC TDC ETT - - # configure CAEN inputs - caen_fragment_name: [ "CAENV1730", "ContainerCAENV1730" ] + caen_fragment_name: [ "daq:CAENV1730", "daq:ContainerCAENV1730" ] + frameshift_module_label: "frameshift" ignore_fragid: [] # boards to ignore for testing purposes, fragid==8 is the timing board, this is 0-indexed! - nominal_length: 5000 # expected waveform length - raw_ts_correction: 367000 # ns, correction to the rawheader timestamp - - # configure SPEC TDC inputs - spectdc_product_name: "tdcdecoder" # used if the tdc decoder has been run - spectdc_ftrig_ch: 3 # channel number of the spec tdc that stores the flash trigger - spectdc_etrig_ch: 4 # channel number of the spec tdc that stores the event trigger - - # configure PTB inputs; used in case spec etrig is missing - ptb_product_name: "ptbdecoder" - ptb_allowed_hlts: [1,2,3,4,5,14,15] # beam, offbeam, and crossing muon hlts - ptb_trig_raw_diff_max: 3000000 # the maximum allowed time difference in ns + nominal_length: 5000 # ticks, expected waveform length # configure output options # use time difference between reference + CAEN used to determine rollover and/or weirdness allowed_time_diff: 3000 # us!! not nanoseconds + # configure the output data product instance names pmtInstanceName: "PMTChannels" ftrigInstanceName: "FTrigChannels" @@ -41,7 +26,9 @@ pmtdecoder: # more detailed configurables n_maxflashes: 30 # maximum number of flashes anticipated, this is used for memory allocation n_caenboards: 9 # number of caen boards; this is used for memory allocation - # 1-8 = PMT caens, 9th = timing caen + # 1-8 = PMT caens, 9th = timing caen + n_caenpmtchannels: 15 # number of PMT(!) channels per caen board; for memory allocation of chmap + timing_caen_offset: 900 # the channel offset for the timing caen board threshold_ftrig: 16350 # the threshold for the flash trigger channel fragid_offset: 40960 # the offset to subtract to get single digit fragids; if use_set_map is true, should set this to 0 @@ -54,19 +41,6 @@ pmtdecoder: ## to use this, must set `fragid_offset` to 0!!! use_set_map: false set_fragid_map: [3,7,5,9,11,13,15,17,19] # use if `use_set_map` is true: digitizers 0-8 are mapped to these fragids - - # ordered in sets of 15 by increasing fragID - # the PMT_to_Digitizers.txt file in the same directory is the reference for this - # text file is obtained from docdb 33668-v9 - ch_map: [17,16,71,70,41,40,15,14,69,68,13,12,67,66,65, # digitzer 0 - 39,38,11,10,9,8,63,62,37,36,7,6,61,60,64, # digitizer 1 - 95,94,149,148,119,118,93,92,147,146,91,90,145,144,143, # digitizer 2 - 117,116,89,88,87,86,141,140,115,114,85,84,139,138,142, # digitizer 3 - 173,172,227,226,197,196,171,170,225,224,169,168,223,222,167, # digitizer 4 - 195,194,221,220,165,164,219,218,193,192,163,162,217,216,166, # digitizer 5 - 251,250,305,304,275,274,249,248,303,302,247,246,301,300,245, # digitizer 6 - 273,272,299,298,243,242,297,296,271,270,241,240,295,294,244, # digitizer 7 - 900,901,902,903,904,905,906,907,908,909,910,911,912,913,914] # digitizer 8 (TIMING CAEN) } END_PROLOG diff --git a/sbndcode/Decoders/PMT/run_pmtdecoder.fcl b/sbndcode/Decoders/PMT/run_pmtdecoder.fcl index b0c8f5a32..de73c2809 100644 --- a/sbndcode/Decoders/PMT/run_pmtdecoder.fcl +++ b/sbndcode/Decoders/PMT/run_pmtdecoder.fcl @@ -1,10 +1,13 @@ +#include "services_sbnd.fcl" #include "pmtdecoder.fcl" #include "spectdc_decoder_sbnd.fcl" #include "SBNDPTBDecoderDefaults.fcl" +#include "frameshift_sbnd_data.fcl" process_name: PMTDecoder services: { + @table::sbnd_services # includes PMTCalibDB, TriggerEmulation TFileService : {fileName: "decoder_hist.root"} TriggerEmulationService: { service_type: TriggerEmulationService @@ -29,6 +32,7 @@ outputs: { module_type: RootOutput fileName: "pmtdecoder-art.root" + dataTier: "reconstructed" } } @@ -39,6 +43,7 @@ physics: tdcdecoder: @local::spec_tdc_decoder_sbnd pmtdecoder: @local::pmtdecoder ptbdecoder: @local::SBNDPTBDecoderDefaults + frameshift: @local::frameshift_data } filters:{} @@ -46,7 +51,7 @@ physics: analyzers:{} - reco: [tdcdecoder, ptbdecoder, pmtdecoder] + reco: [tdcdecoder, ptbdecoder, frameshift, pmtdecoder] # define the output stream, there could be more than one if using filters diff --git a/sbndcode/Decoders/SPECTDC/run_ptb_tdc_decoder_frameshift.fcl b/sbndcode/Decoders/SPECTDC/run_ptb_tdc_decoder_frameshift.fcl new file mode 100644 index 000000000..e9b6a2050 --- /dev/null +++ b/sbndcode/Decoders/SPECTDC/run_ptb_tdc_decoder_frameshift.fcl @@ -0,0 +1,63 @@ +#include "spectdc_decoder_sbnd.fcl" +#include "SBNDPTBDecoderDefaults.fcl" + +#include "frameshift_sbnd_data.fcl" + +process_name: Decoders + +services: { + TFileService : {fileName: "decoder_hist.root"} + TriggerEmulationService: { + service_type: TriggerEmulationService + MonWidth: 12 + TotalCAENBoards: 8 + PMTPerBoard: 15 + Baseline: 14250 + MC: false + } +} + +#source is a root file +source: +{ + module_type: RootInput + maxEvents: -1 # Number of events to create +} + +outputs: +{ + out1: + { + module_type: RootOutput + fileName: "decoder-frameshift-art.root" + } +} + +physics: +{ + + producers:{ + tdcdecoder: @local::spec_tdc_decoder_sbnd + ptbdecoder: @local::SBNDPTBDecoderDefaults + frameshift: @local::frameshift_data + } + + filters:{} + + analyzers:{} + + + reco: [tdcdecoder, ptbdecoder, frameshift] + + + # define the output stream, there could be more than one if using filters + stream1: [out1] + + # trigger_paths is a keyword and contains the paths that modify the art::event, + # ie filters and producers + trigger_paths: [reco] + + # end_paths is a keyword and contains the paths that do not modify the art::Event, + # ie analyzers and output streams. these all run simultaneously + end_paths: [stream1] +} diff --git a/sbndcode/Decoders/XARAPUCA/SBNDXARAPUCADecoder_module.cc b/sbndcode/Decoders/XARAPUCA/SBNDXARAPUCADecoder_module.cc index ac94c56e1..b3cf6cd01 100644 --- a/sbndcode/Decoders/XARAPUCA/SBNDXARAPUCADecoder_module.cc +++ b/sbndcode/Decoders/XARAPUCA/SBNDXARAPUCADecoder_module.cc @@ -10,9 +10,10 @@ * @brief Defines and implements the SBNDXARAPUCADecoder class which inherits from an art::EDProducer * as the decoder for V1740B digitizers, intended for the X-ARAPUCAs. * @details The current version of the SBND X-ARAPUCAs decoder implements the updates shown in - * the SBN Document 38475-v1 in the SBN Document Database. + * the SBN Document 43891-v1 in the SBN Document Database. * @note A Python version of the binary decoding is available for testing purposes. You can find * it [here: V1740 binary decoder](https://github.com/aliciavr/V1740_binary_decoder). + * @version 5.0 */ #include "art/Framework/Core/EDProducer.h" @@ -32,12 +33,11 @@ #include "sbndaq-artdaq-core/Overlays/FragmentType.hh" #include "sbndaq-artdaq-core/Overlays/SBND/CAENV1740Fragment.hh" -#include "sbndcode/Decoders/PTB/sbndptb.h" -#include "sbndaq-artdaq-core/Overlays/SBND/PTBFragment.hh" - #include "sbndcode/Timing/SBNDRawTimingObj.h" #include "sbnobj/SBND/Timing/DAQTimestamp.hh" +#include "sbnobj/SBND/Timing/FrameShiftInfo.hh" + #include "lardataobj/RawData/OpDetWaveform.h" #include "art_root_io/TFileService.h" @@ -63,6 +63,15 @@ namespace sbndaq { */ class sbndaq::SBNDXARAPUCADecoder : public art::EDProducer { public: + +/** + * @brief Constructor of the class SBNDXARAPUCADecoder. Given a fhicl::ParameterSet object it gets the configuration needed by this + * art::EDProducer from the FHiCL file. + * @param[in] p A set of parameters containing configuration values for the decoder from the FHiCL file. + * @details This constructor accesses the parameters in the FHiCL configuration file and set the configuration to be run for + * this module. It sets settings for CAEN fragments, board information, timing, debug and verbosity options, and creates + * instances of the products associated with this module: a vector of raw::OpDetWaveform objects and a raw::TimingReferenceInfo object. + */ explicit SBNDXARAPUCADecoder(fhicl::ParameterSet const& p); SBNDXARAPUCADecoder(SBNDXARAPUCADecoder const&) = delete; @@ -71,10 +80,20 @@ class sbndaq::SBNDXARAPUCADecoder : public art::EDProducer { SBNDXARAPUCADecoder& operator=(SBNDXARAPUCADecoder&&) = delete; // Main class method. + /** + * @brief Main function of the art::EDProducer: the produce function analyzes every event producing waveforms after the decoding. + * @param[in] e The event to be processed. + * @details It is the main function of the art::EDProducer module: it processes each event, searches for CAEN V1740 fragments, decodes + * them and creates the output products. + * 1. It gets the raw timestamp from the artdaq::RawEventHeader product. + * 2. It tries to get a valid timing reference from the selected time frame (FRAMESHIFT or CAEN-only). If the selected + * timing frame is not found, it tries to get the next one in the priority list. + * 3. It searches for CAEN V1740 fragments in the event, decodes them and creates the output products: a vector of raw::OpDetWaveform. + * 4. It dumps the products in the event. + */ void produce(art::Event& e) override; private: - constexpr static uint64_t NANOSEC_IN_SEC = 1'000'000'000; /**< Number of nanoseconds in one second. */ constexpr static uint64_t MICROSEC_IN_NANOSEC = 1'000; /**< Number of nanoseconds in one microsecond. */ constexpr static double NANOSEC_TO_MICROSEC = 1E-3; /**< Conversion factor from nanoseconds to microseconds. */ @@ -87,12 +106,10 @@ class sbndaq::SBNDXARAPUCADecoder : public art::EDProducer { constexpr static uint32_t NUM_CONSECUTIVE_SAMPLES = 3; /**< Number of consecutive samples for the same channel when decoding waveforms. */ constexpr static uint16_t NANOSEC_PER_TICK = 8; /**< Number of nanoseconds per Trigger Time Tag (TTT) tick. */ - constexpr static uint16_t SPEC_TDC_TIMING = 0; /**< Timing reference frame: ETRIG (Event Trigger) timestamp from the SPEC-TDC. */ - constexpr static uint16_t PTB_TIMING = 1; /**< Timing reference frame: HLT (High Level Trigger) timestamp from the PTB (Penn Trigger Board). */ - constexpr static uint16_t CAEN_ONLY_TIMING = 2; /**< Timing reference frame: CAEN-only. */ + constexpr static uint16_t FRAMESHIFT_TIMING = 0; /**< Timing reference frame: FrameShift timestamp. */ + constexpr static uint16_t CAEN_ONLY_TIMING = 1; /**< Timing reference frame: CAEN-only. */ - constexpr static uint16_t HLT_NOT_FOUND = 999; /**< Random value to indicate no HL trigger found. */ - constexpr static uint16_t HLT_TOO_FAR = 1000; /**< Random value to indicate HL trigger found but too far from the raw timestamp. */ + constexpr static uint16_t TTT_DEFAULT = 0; /**< Default integer value for the nominal TTT. */ unsigned int fevent_counter; /**< Event counter. */ @@ -100,64 +117,300 @@ class sbndaq::SBNDXARAPUCADecoder : public art::EDProducer { unsigned int ffragment_id_offset; /**< Offset value to get the slot of every CAEN board. */ std::vector fboard_id_list; /**< Possible CAEN board identifiers based on the slot they are installed. */ - std::string fcaen_module_label; /**< Label identifying the module where the CAEN fragments come from. */ - std::vector fcaen_fragment_names; /**< Fragments accepted by this decoder. */ + std::vector fcaen_fragment_names; /**< Fragments accepted by this decoder. */ std::string fwaveforms_instance_name; /**< Name assigned to the product instance containing the waveforms for each channel and board generated by this art::EDProducer. */ std::string ftiming_ref_instance_name; /**< Name assigned to the product instance containing the timing reference information for each event generated by this art::EDProducer. */ - uint16_t ftiming_priority; /**< Timing priority configured: 0 SPEC-TDC RWM, 1 SPEC-TDC ETT, 2 CAEN-only. */ + bool fframeshift_mode; /**< FrameShift mode: 0 use timing priority, 1 use FrameShift timestamp as reference. */ uint16_t factive_timing_frame; /**< Active timing frame while processing each event. */ unsigned int fns_per_sample; /**< Number of nanoseconds per sample. */ - std::string fraw_event_header_product_name; /**< Label identifying the module that created the artdaq::detail::RawEventHeader data product. */ - uint16_t fraw_trig_max_diff; /**< Maximum difference (in ns) allowed between trigger timestamps and raw timestamp. */ - uint32_t fraw_timestamp_correction; /**< Correction (in ns) to be applied to the raw timestamp. */ - - std::string fspectdc_product_name; /**< Label identifying the module containing the SPEC-TDC decoder products. */ - uint32_t fspectdc_ftrig_ch; /**< Channel assigned by the SPEC-TDC to the flash triggers. */ - uint32_t fspectdc_etrig_ch; /**< Channel assigned by the SPEC-TDC to the event triggers. */ - - std::string fptb_product_name; /**< Label identifying the module containing the PTB decoder products. */ - std::vector fallowed_hl_triggers; /**< List of allowed high-level triggers from the PTB. */ + art::InputTag fframeshift_product_name; /**< Label identifying the module containing the FrameShift module products. */ art::ServiceHandle tfs; /**< ServiceHandle object to store the histograms in the decoder_hist.root output file. */ int fstore_debug_waveforms; /**< Number of waveforms to store in the ServiceHandle object for debugging purposes (0: none, -1: all, n: first n waveforms each event). */ - bool fdebug_tdc_handle; /**< If `true` SPEC-TDC information is printed. */ - bool fdebug_ptb_handle; /**< If `true` PTB information is printed. */ + bool fcombine_ext_frag; /**< If `true` combines extended fragments into a single raw::OpDetWaveform object. */ + int32_t fallowed_jittering; /**< Allowed jittering (in ns) between fragments to be combined. */ + + bool fdebug_frameshift_handle; /**< If `true` FrameShift information is printed. */ bool fdebug_fragments_handle; /**< If `true` V1740B CAEN fragments art::Handle information is printed. */ bool fdebug_timing; /**< If `true` timing data is printed. */ bool fdebug_buffer; /**< If `true` the buffer status is printed. */ bool fdebug_waveforms; /**< If `true` waveforms decoding data is printed. */ + bool fdebug_extended_fragments; /**< If `true` extended fragments information is printed. */ + bool fdebug_jittering; /**< If `true` trigger jittering information is printed. */ bool fverbose; /**< If `true` it increases verbosity of console output for detailed processing steps. */ - - // Class methods. - void decode_fragment(uint64_t timestamp, std::vector & fragment_indices, const artdaq::Fragment& fragment, std::vector & prod_wvfms); - bool get_ptb_hlt_timestamp(art::Event& e, uint64_t corr_raw_timestamp, uint64_t & timestamp, uint16_t & hlt_code); - bool get_spec_tdc_etrig_timestamp(art::Event& e, uint64_t corr_raw_timestamp, uint64_t & timestamp); + // Auxiliary methods. + +/** + * @brief Returns the signed difference between two timestamps (t1 - t2). + * @param[in] t1 The first timestamp in nanoseconds. + * @param[in] t2 The second timestamp in nanoseconds. + * @return The signed difference between the two timestamps (t1 - t2) in nanoseconds. + * @pre The timestamp is given in nanoseconds and the maximum value is less than 2^63 ns. + */ + int64_t signed_difference_UTC_timestamp(uint64_t t1, uint64_t t2); + +/** + * @brief Returns the absolute difference between two timestamps. + * @param[in] t1 The first timestamp in nanoseconds. + * @param[in] t2 The second timestamp in nanoseconds. + * @return The absolute difference between the two timestamps in nanoseconds. + * @pre The timestamp is given in nanoseconds and the maximum value is less than 2^63 ns. + */ + uint64_t abs_difference(uint64_t t1, uint64_t t2); + +/** + * @brief Formats a timestamp in nanoseconds into a easily readable string format. + * @param[in] timestamp The timestamp in nanoseconds to be formatted. + * @return A string representation of the timestamp in the format "(seconds)nanoseconds ns". + * @details The function divides the input timestamp by 1,000,000,000 to obtain the seconds component and uses the modulus operator to get the remaining nanoseconds. + */ + std::string print_timestamp(uint64_t timestamp); - void save_prod_wvfm(size_t board_idx, size_t ch, double ini_wvfm_timestamp, const std::vector > & wvfms, std::vector & prod_wvfms); - void save_debug_wvfm(size_t board_idx, size_t fragment_idx, int ch, double ini_wvfm_timestamp, double end_wvfm_timestamp, const std::vector > & wvfms); + // Timing. + +/** + * @brief This function shifts the initial and end timestamps of a waveform based on the provided timing information. + * + * @param[in] TTT_ticks The trigger time tag in ticks from the CAEN V1740 header. + * @param[in] TTT_end_ns The end time of the Trigger Time Tag (TTT) in nanoseconds. + * @param[in] frag_timestamp The timestamp of the fragment in nanoseconds. + * @param[in] timestamp The reference timestamp (from FRAMESHIFT or CAEN-only) in nanoseconds. + * @param[in] num_samples_per_wvfm The number of samples per waveform. + * @param[out] ini_wvfm_timestamp The initial timestamp of the waveform in microseconds (output). + * @param[out] end_wvfm_timestamp The end timestamp of the waveform in microseconds (output). + * + * @details + * This function calculates the initial and end timestamps for a waveform based on the provided timing information. It takes into account + * potential rollovers in the trigger time tag and adjusts the timestamps accordingly. The function supports different timing frames, including + * FRAMESHIFT and CAEN-only timing. The calculated timestamps are returned in microseconds. + * + */ + void shift_time(uint64_t TTT_ticks, int64_t TTT_end_ns, uint64_t frag_timestamp, uint64_t timestamp, uint32_t num_samples_per_wvfm, double& ini_wvfm_timestamp, double& end_wvfm_timestamp, size_t& appended_samples); +/** + * @brief Searches for the frameshift timestamp in the event if any frameshift product is found and returns it as a TimingReferenceInfo object. + * @param[in] e The event to be processed. + * @param[in,out] timestamp The frameshift timestamp (if found). + * @param[in,out] timing_type The type of the frameshift timestamp (if found). + * @param[in,out] timing_channel The channel of the frameshift timestamp (if found). + * @return A unique pointer to a TimingReferenceInfo object containing the frameshift timestamp information (if found), or nullptr if no valid frameshift product is found. + * @details It searches for the frameshift products in the event and looks for the default frameshift timestamp, type and channel. If any valid frameshift product is found, it returns its information as a TimingReferenceInfo object. + */ + std::unique_ptr get_frameshift_timestamp(art::Event& e, uint64_t& timestamp); + + // Waveforms decoding. + +/** + * @brief Decodes the waveforms from a CAEN V1740 fragment (binary decoding stage). + * + * @param[in] fragment The input CAEN V1740 fragment containing raw data to be decoded. + * @param[out] wvfms A 2D vector where the decoded waveforms will be stored. Each inner vector corresponds to a channel's waveform. + * @param[in] header_size The size of the event header in bytes. + * @param[in] num_channels The number of channels in the fragment. + * @param[in] num_samples_per_wvfm The number of samples per waveform for each channel. + * @param[in] num_words_per_wvfms The total number of 32-bit words containing waveform data in the fragment. + * @param[in] num_samples_per_group The number of samples per group of channels (8 channels per group). + * + * @details This function reads raw 32-bit words from the fragment, storing them in a buffer, and extracts 12-bit samples. + * It assigns each sample to the appropriate channel and sample index based on its position in the sequence. The decoding process includes: + * - Initializing a buffer to hold incoming data and tracking the number of bits currently stored. + * - Iterating over each 32-bit word in the waveform data section of the fragment, adding it to the buffer. + * - Extracting 12-bit samples from the buffer as long as there are enough bits available, and assigning them to their respective channels and sample indices using calculated formulas. + * The function ensures that all samples are correctly mapped to their channels, taking into account the grouping of channels and the interleaving of samples. + * @note The function assumes that the input fragment is valid and contains the expected structure for a CAEN V1740 device. + * @see SBN Document 38475-v1 for more details on the binary decoding. + */ + void decode_waveforms(const artdaq::Fragment& fragment, std::vector>& wvfms, size_t header_size, uint32_t num_channels, uint32_t num_samples_per_wvfm, uint32_t num_words_per_wvfms, uint32_t num_samples_per_group); + +/** + * @brief Extract a sample from a 64-bit buffer using the specified bit positions. + * + * @param[in] buffer An unsigned 64-bit integer which represents a temporal buffer for the read words and where the samples are extracted from. + * @param[in] msb An unsigned 32-bit integer representing the most significative bit (MSB) where the readout from the buffer paramter. + * @param[in] lsb An unsigned 32-bit integer representing the less significative bit (LSB) from we end read + * + * @details The function shifts the buffer to the right by the number of positions specified by `lsb` so that the least significant bit of the + * sample aligns with bit 0. It then applies a mask to isolate the bits between `lsb` and `msb`, inclusive. + * + * @return The extracted sample as a 16-bit unsigned integer. + */ uint16_t get_sample(uint64_t buffer, uint32_t msb, uint32_t lsb); + +/** + * @brief Read a 32-bit word from the data pointer and advances the pointer. + * + * @param[in, out] data_ptr A reference to a pointer pointing to the current position in the data. + * + * @details This function retrieves a 32-bit word from the memory location pointed to by `data_ptr`. After reading, it advances `data_ptr` to + * the next 32-bit word location. + * + * @return The 32-bit word read from the location pointed to by `data_ptr`. + */ uint32_t read_word(const uint32_t* & data_ptr); + +/** + * @brief Generates a unique global channel identifier using the board slot and the channel number of that board. + * + * @param[in] board Index of the board in `fboard_id_list` from which to derive the board slot. + * @param[in] board_channel The specific channel number on the given board. + * + * @details This function computes a `channel_id` by combining the board slot and the specific + * channel number on that board. + * The unique identifier `channel_id` (\f$ CH_{ID} $\f) is computed as follows: + * \f[ + * CH\_{ID} = B\_{ID} \times 100 + CH\_{B} + * \f] + * + * Where: + * - \f$ B\_{ID} \f$ is the fragment ID retrieved from `fboard_id_list` based on the slot. + * - \f$ CH\_B \f$ is the channel number on that board. + * + * @return A unique identifier for the specified channel as an unsigned integer. + */ unsigned int get_channel_id(unsigned int board, unsigned int board_channel); - std::string print_timestamp(uint64_t timestamp); - int64_t signed_difference(uint64_t t1, uint64_t t2); - uint64_t abs_difference(uint64_t t1, uint64_t t2); -}; + // Combines decoded waveforms. /** - * @brief Constructor of the class SBNDXARAPUCADecoder. Given a fhicl::ParameterSet object it gets the configuration needed by this - * art::EDProducer from the FHiCL file. - * @param[in] p A set of parameters containing configuration values for the decoder from the FHiCL file. - * @details This constructor accesses the parameters in the FHiCL configuration file and set the configuration to be run for - * this module. It sets settings for CAEN fragments, board information, timing, debug and verbosity options, and creates - * instances of the products associated with this module: a vector of raw::OpDetWaveform objects and a raw::TimingReferenceInfo object. + * @brief Combines the waveforms from the current fragment with the previously stored waveforms. + * @param[in,out] wvfms A 2D vector containing the (combined if it was needed before) waveforms. + * @param[in] fragment_wvfms A 2D vector containing the waveforms from the current fragment to be combined. + * @param[in] num_channels The number of channels per board. + * + * @details + * The function performs the following steps: + * 1. Checks if the `wvfms` vector is empty. If it is, it resizes it to accommodate `num_channels` channels. + * 2. Iterates over each channel from 0 to `num_channels - 1`. + * 3. For each channel, it appends the samples from `fragment_wvfms` to the corresponding channel in `wvfms`. + * + * This approach ensures that waveforms from multiple fragments are concatenated correctly, maintaining the order of samples for each channel. + * + * @pre The `fragment_wvfms` vector should contain waveforms for all channels of the board being processed. + * @pre The `wvfms` vector should be either empty or already contain waveforms for all channels of the board being processed. + */ + void append_waveforms(std::vector>& wvfms, const std::vector>& fragment_wvfms, uint32_t num_channels); + + // Dumps and saves waveforms. +/** + * @brief Dump products into a `raw::OpDetWaveform` object and into a debug histogram file. + * @param[in,out] prod_wvfms A reference to the vector where the produced `raw::OpDetWaveform` objects are dumped into products. + * @param[in,out] wvfms A 2D vector containing the (combined if needed) waveforms. + * @param[in,out] fragment_indices A reference to a vector keeping track of the number of fragments decoded per board. + * @param[in] board_idx The board index (position in the list of boards). + * @param[in] num_channels The number of channels per board. + * @param[in] ini_wvfm_timestamp The initial timestamp of the waveform in microseconds. + * @param[in] end_wvfm_timestamp The final timestamp of the waveform in microseconds. + * + * @details + * The function performs the following steps: + * 1. Determines the number of debug waveforms to be stored based on the configuration parameter `fstore_debug_waveforms`. + * - If `fstore_debug_waveforms` is set to -1, all channels are considered for debug storage. + * - Otherwise, it takes the minimum between `num_channels` and `fstore_debug_waveforms`. + * 2. Iterates over each channel up to the determined number of debug waveforms: + * - Calls `save_prod_wvfm` to convert and store the waveform in the products. + * - Calls `save_debug_wvfm` to save the waveform as a histogram for debugging purposes. + * 3. For any remaining channels beyond the debug limit, it only calls `save_prod_wvfm` to store the waveform in the products. + * 4. Clears the `wvfms` vector to free up memory after processing. + * + * @pre ini_wvfm_timestamp and end_wvfm_timestamp is assumed to be given in microseconds when the timing frame is not CAEN_ONLY_TIMING. + * @pre The `wvfms` vector should contain waveforms for all channels of the board specified by `board_idx`. + * @pre The `fragment_indices` vector should have been initialized with a size equal to the number of boards being processed. + * @pre The `prod_wvfms` vector should be ready to accept new `raw::OpDetWaveform` objects. + * + * @see raw::OpDetWaveform + * @see save_prod_wvfm + * @see save_debug_wvfm + */ + void dump_waveforms(std::vector & prod_wvfms, std::vector>& wvfms, size_t first_frag_idx, size_t board_index, uint32_t num_channels, double ini_wvfm_timestamp, double end_wvfm_timestamp); + +/** + * @brief Converts a production waveform from raw ADC data into a `raw::OpDetWaveform` object, + * assigning it a global channel ID and timestamp, and appends it to the output collection. + * + * @param[in] b The board index (position in the list of boards). + * @param[in] ch The channel index (channel number from which the waveform is extracted). + * @param[in] ini_wvfm_timestamp The initial timestamp of the waveform in microseconds. + * @param[in] wvfms A 2D vector containing the waveforms. Each inner vector corresponds to a channel, + * and each element of the inner vector represents a sample (ADC count). + * @param[out] prod_wvfms A reference to the vector where the produced `raw::OpDetWaveform` objects + * will be stored. + * + * @details + * The function performs the following steps: + * 1. Determines the global channel ID for the specified board and channel using `get_channel_id`. + * 2. Creates a `raw::OpDetWaveform` object with the initial timestamp, global channel ID, and waveform data. + * 3. Appends the `raw::OpDetWaveform` object to the provided `prod_wvfms` collection. + * + * @pre ini_wvfm_timestamp is assumed to be given in microseconds when the timing frame is not CAEN_ONLY_TIMING. + * + * @see raw::OpDetWaveform */ + void save_prod_wvfm(size_t board_idx, size_t ch, double ini_wvfm_timestamp, const std::vector > & wvfms, std::vector & prod_wvfms); + +/** + * @brief Saves debug waveform as a histogram for debugging purposes, based on the waveform data provided. + * The histogram is labeled with event number, fragment, board slot ID and channel information. + * + * @param[in] board_idx The board index (position in the list of boards). + * @param[in] fragment_idx The fragment index (order in decoding). + * @param[in] ch The channel index (channel number from which the waveform is extracted). + * @param[in] ini_wvfm_timestamp The initial timestamp of the waveform. + * @param[in] end_wvfm_timestamp The final timestamp of the waveform. + * @param[in] wvfms A 2D vector containing the waveforms. Each inner vector corresponds to a channel, + * and each element of the inner vector represents a sample (ADC count). + * + * @details + * The function generates a histogram with: + * - X-axis representing time in microseconds if the timing frame is from the FRAMESHIFT reference, otherwise in samples. + * - Y-axis representing ADC counts. + * The histogram is stored using ROOT's `TFileService`. + * + * @pre ini_wvfm_timestamp and end_wvfm_timestamp are assumed to be given in microseconds when the timing frame is not CAEN_ONLY_TIMING. + * + * @see TH1I, TFileService + */ + void save_debug_wvfm(size_t board_idx, size_t fragment_idx, int ch, double ini_wvfm_timestamp, double end_wvfm_timestamp, const std::vector > & wvfms); + + // Main processing method. + +/** + * @brief Decodes a CAEN V1740 fragment, extracts the waveforms and combines them into the output product. + * @param[in] timestamp The valid timestamp used as reference for the event. + * @param[in,out] first_frag_timestamp The nominal timestamp calculated for the fragment being processed. + * @param[in, out] first_frag_idx The first index of the fragments being combined. + * @param[in,out] first_TTT The nominal Trigger Time Tag (TTT) calculated for the fragment being processed. + * @param[in,out] fragment_indices A vector containing the indices of the fragments already processed for each board. + * @param[in] fragment The artdaq::Fragment object to be processed. + * @param[in,out] prod_wvfms The vector of raw::OpDetWaveform objects to be filled with the decoded waveforms. + * @param[in,out] wvfms A 2D vector containing the waveforms for all channels and boards decoded so far. + * @param[in] last_one A boolean flag indicating if the fragment being processed is the last one in the event. + + * @details This method decodes a CAEN V1740 fragment, extracts the waveforms for each channel, shifts them in time according to the + * timing reference and combines them into the output product. + * 1. It checks if the fragment ID corresponds to a valid board. + * 2. It accesses the metadata and header of the fragment to get information about the number of channels, samples, words, etc. + * 3. It calculates the number of samples per waveform and checks if it is a nominal or extended fragment. + * 4. It decodes the waveforms for each channel in the fragment. + * 5. It shifts the waveforms in time according to the timing reference. + * 6. It combines the decoded waveforms into the output product. + * 7. It dumps the waveforms if required. + * + * @pre The art::Event object containing the fragment has been processed to get a valid timing reference. + * @post The vector of raw::OpDetWaveform objects is filled with the decoded waveforms from the fragment. + * @see shift_time + * @see decode_waveforms + * @see append_waveforms + * @see dump_waveforms + */ + void decode_fragment(uint64_t timestamp, uint64_t& first_frag_timestamp, size_t& first_frag_idx, int32_t& first_TTT, int32_t& prev_TTT, std::vector & fragment_indices, const artdaq::Fragment& fragment, std::vector & prod_wvfms, std::vector>& wvfms, bool last_one, size_t& appended_samples); + +}; + sbndaq::SBNDXARAPUCADecoder::SBNDXARAPUCADecoder(fhicl::ParameterSet const& p) : EDProducer{p} { @@ -165,8 +418,7 @@ sbndaq::SBNDXARAPUCADecoder::SBNDXARAPUCADecoder(fhicl::ParameterSet const& p) fevent_counter = 0; // Gets the CAEN fragments information. - fcaen_module_label = p.get ("caen_module_label", "daq"); - fcaen_fragment_names = p.get > ("caen_fragment_names", {"CAENV1740", "ContainerCAENV1740"}); + fcaen_fragment_names = p.get > ("caen_fragment_names", {"CAENV1740", "ContainerCAENV1740"}); // Gets the CAEN boards information. ffragment_id_offset = p.get ("fragment_id_offset", 41216); @@ -178,33 +430,26 @@ sbndaq::SBNDXARAPUCADecoder::SBNDXARAPUCADecoder(fhicl::ParameterSet const& p) ftiming_ref_instance_name = p.get ("timing_ref_instance_name", ""); // Gets timing configuration. - ftiming_priority = p.get ("timing_priority", SPEC_TDC_TIMING); factive_timing_frame = CAEN_ONLY_TIMING; fns_per_sample = p.get ("ns_per_sample", 16); - // Gets the raw timestamp configuration. - fraw_event_header_product_name = p.get ("raw_timestamp_product_name", "daq"); - fraw_trig_max_diff = p.get ("raw_trig_max_diff", 3000000); - fraw_timestamp_correction = p.get ("raw_timestamp_correction", 367000); - - // SPEC-TDC access configuration. - fspectdc_product_name = p.get("spectdc_product_name", "tdcdecoder"); - fspectdc_ftrig_ch = p.get("spectdc_ftrig_ch", 3); - fspectdc_etrig_ch = p.get("spectdc_etrig_ch", 4); - - // PTB access configuration. - fptb_product_name = p.get("ptb_product_name", "ptbdecoder"); - fallowed_hl_triggers = p.get>("allowed_hl_triggers", {1, 2, 3, 4, 5, 14, 15}); + // FrameShift access configuration. + fframeshift_product_name = p.get("frameshift_product_name", "frameshift"); // Gets the number of waveforms to store in the debug output file. fstore_debug_waveforms = p.get ("store_debug_waveforms", 0); + // Gets the combination of extended fragments option. + fcombine_ext_frag = p.get ("combine_ext_frag", true); + fallowed_jittering = p.get ("allowed_jittering", 64); + // Gets the debug and verbose options. - fdebug_ptb_handle = p.get ("debug_ptb_handle", false); - fdebug_tdc_handle = p.get ("debug_tdc_handle", false); + fdebug_frameshift_handle = p.get ("debug_frameshift_handle", false); fdebug_fragments_handle = p.get ("debug_fragments_handle", false); fdebug_timing = p.get ("debug_timing", false); fdebug_waveforms = p.get ("debug_waveforms", false); + fdebug_extended_fragments = p.get ("debug_extended_fragments", false); + fdebug_jittering = p.get ("debug_jittering", false); fdebug_buffer = p.get ("debug_buffer", false); fverbose = p.get ("verbose", false); @@ -213,19 +458,7 @@ sbndaq::SBNDXARAPUCADecoder::SBNDXARAPUCADecoder(fhicl::ParameterSet const& p) produces (ftiming_ref_instance_name); } - /** - * @brief Main function of the art::EDProducer: the produce function analyzes every event producing waveforms after the decoding. - * @param[in] e The event to be processed. - * @details It is the main function of the art::EDProducer module: it processes each event, searches for CAEN V1740 fragments, decodes - * them and creates the output products. - * 1. It gets the raw timestamp from the artdaq::RawEventHeader product. - * 2. It tries to get a valid timing reference from the selected time frame (SPEC-TDC ETRIG, PTB HLT or CAEN-only). If the selected - * timing frame is not found, it tries to get the next one in the priority list. - * 3. It searches for CAEN V1740 fragments in the event, decodes them and creates the output products: a vector of raw::OpDetWaveform. - * 4. It dumps the products in the event. - */ -void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) -{ +void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) { if (fverbose) std::cout << "\n > SBNDXARAPUCADecoder::produce: entering the produce function." << std::endl; // Advances the event counter. @@ -233,80 +466,48 @@ void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) // Initializes the output instance products. auto prod_wvfms = std::make_unique > (); - auto prod_event_timing_ref_info = std::make_unique (); - prod_event_timing_ref_info->timingType = CAEN_ONLY_TIMING; // CAEN-only timing frame as default. - prod_event_timing_ref_info->timingChannel = 0; + // CAEN-only timing frame as default. + auto prod_event_timing_ref_info = std::make_unique(CAEN_ONLY_TIMING, 0); // CAEN-only timing frame as default. if (fverbose) std::cout << "\n > SBNDXARAPUCADecoder::produce: products initialized." << std::endl; - // Gets the raw timestamp from the artdaq::RawEventHeader product. - uint64_t corr_raw_timestamp = 0; - art::Handle header_handle; - e.getByLabel(fraw_event_header_product_name, "RawEventHeader", header_handle); - if (header_handle.isValid()) { - auto raw_header = artdaq::RawEvent(*header_handle); - corr_raw_timestamp = raw_header.timestamp() - fraw_timestamp_correction; - - if (fverbose | fdebug_tdc_handle | fdebug_ptb_handle | fdebug_timing) std::cout << "\n > SBNDXARAPUCADecoder::produce: RAW timestamp from the event header: " << print_timestamp(corr_raw_timestamp) << "." << std::endl; - if (fdebug_tdc_handle | fdebug_ptb_handle | fdebug_timing) { - std::cout << "\t RAW timestamp: BEFORE correction: " << print_timestamp(raw_header.timestamp()) << " - AFTER correction: " << print_timestamp(corr_raw_timestamp) << "." << std::endl; - } - } - - // [0] ETRIG timestamp, [1] HLT timestamp, [2] CAEN-only timestamp. - uint64_t timestamp = 0; - - // Gets the SPEC-TDC product (if any). - if (ftiming_priority == SPEC_TDC_TIMING) { - if (get_spec_tdc_etrig_timestamp(e, corr_raw_timestamp, timestamp)) { - prod_event_timing_ref_info->timingType = SPEC_TDC_TIMING; - prod_event_timing_ref_info->timingChannel = fspectdc_etrig_ch; - factive_timing_frame = SPEC_TDC_TIMING; - if (fverbose | fdebug_tdc_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: SPEC-TDC ETRIG reference, " << print_timestamp(timestamp) << ", found close enough to the raw timestamp. Using SPEC-TDC timing frame as reference." << std::endl; - } else { - ftiming_priority = PTB_TIMING; - if (fverbose | fdebug_tdc_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: SPEC-TDC ETRIG reference not found. Trying to find PTB reference..." << std::endl; - } - } - - // Gets the PTB product (if any). - if (ftiming_priority == PTB_TIMING) { - uint16_t hlt_code = HLT_NOT_FOUND; - if (get_ptb_hlt_timestamp(e, corr_raw_timestamp, timestamp, hlt_code)) { - prod_event_timing_ref_info->timingType = PTB_TIMING; - prod_event_timing_ref_info->timingChannel = hlt_code; - factive_timing_frame = PTB_TIMING; - if (fverbose | fdebug_ptb_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: PTB HLT reference, " << print_timestamp(timestamp) << ", found close enough to the raw timestamp. Using PTB timing frame as reference." << std::endl; - } else { - ftiming_priority = CAEN_ONLY_TIMING; - if (fverbose | fdebug_ptb_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: PTB reference not found. Searching for CAEN-only reference..." << std::endl; + uint64_t timestamp = 0; + + // ============= FRAMESHIFT TIMING ============ // + if (auto const fs_timing_ref_info = get_frameshift_timestamp(e, timestamp)) { + *prod_event_timing_ref_info = *fs_timing_ref_info; + factive_timing_frame = FRAMESHIFT_TIMING; + if (fverbose | fdebug_frameshift_handle) { + std::cout << "\n > SBNDXARAPUCADecoder::produce: FRAMESHIFT reference, " << print_timestamp(timestamp) << ", found. Using FRAMESHIFT timing frame as reference." << std::endl; + std::cout << "\t FRAMESHIFT reference details: timing type: " << prod_event_timing_ref_info->timingType << ", timing channel: " << prod_event_timing_ref_info->timingChannel << "." << std::endl; } + } else { + if (fverbose | fdebug_frameshift_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: FrameShift reference not found. Using CAEN-only timing frame as reference." << std::endl; } + // ============================================== // - if (fverbose | fdebug_tdc_handle | fdebug_ptb_handle | fdebug_timing) { - std::cout << "\n > SBNDXARAPUCADecoder::produce: selected timing frame: "; - - if (factive_timing_frame == SPEC_TDC_TIMING) { - std::cout << "SPEC-TDC."; - } else if (factive_timing_frame == PTB_TIMING) { - std::cout << "PTB."; - } else if (ftiming_priority == CAEN_ONLY_TIMING) { - std::cout << "CAEN-only."; - } - - std::cout << " Timestamp: " << print_timestamp(timestamp) <<", timing type: " << prod_event_timing_ref_info->timingType << ", channel: " << prod_event_timing_ref_info->timingChannel << "." << std::endl; - } // Flag to track if valid CAEN fragments are found. bool found_caen = false; - std::vector fragment_indices(fnum_caen_boards, 0); + std::vector fragment_indices(fnum_caen_boards, -1); + std::vector> wvfms; + if (fdebug_extended_fragments) { + std::cout << " Waveforms size: " << wvfms.size() << std::endl; + } + + size_t appended_samples = 0; + size_t first_frag_idx = 0; + int32_t first_TTT = TTT_DEFAULT; + int32_t prev_TTT = 0; + uint64_t first_frag_timestamp = TTT_DEFAULT; + bool last_one = false; if (fverbose | fdebug_fragments_handle) std::cout << "\n > SBNDXARAPUCADecoder::produce: searching for V1740 fragments..." << std::endl; - for (const std::string &caen_name: fcaen_fragment_names) { + for (const art::InputTag &caen_name: fcaen_fragment_names) { art::Handle > fragment_handle; - e.getByLabel(fcaen_module_label, caen_name, fragment_handle); + e.getByLabel(caen_name, fragment_handle); // The art::Handle object is not valid. if (!fragment_handle.isValid()) { @@ -338,7 +539,8 @@ void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) for (size_t f = 0; f < num_caen_fragments; f++) { const artdaq::Fragment fragment = *container_fragment[f].get(); - decode_fragment(timestamp, fragment_indices, fragment, *prod_wvfms); + last_one = (f == (num_caen_fragments - 1)); + decode_fragment(timestamp, first_frag_timestamp, first_frag_idx, first_TTT, prev_TTT, fragment_indices, fragment, *prod_wvfms, wvfms, last_one, appended_samples); } // End CAEN V1740 fragments loop. } } // End Container fragments loop. @@ -350,7 +552,8 @@ void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) // It searches for all CAEN V1740 fragments. for (size_t f = 0; f < frag_handle_size; f++) { const artdaq::Fragment fragment = fragment_handle->at(f); - decode_fragment(timestamp, fragment_indices, fragment, *prod_wvfms); + last_one = f == (frag_handle_size - 1); + decode_fragment(timestamp, first_frag_timestamp, first_frag_idx, first_TTT, prev_TTT, fragment_indices, fragment, *prod_wvfms, wvfms, last_one, appended_samples); } // End CAEN V1740 fragments loop. } } // End extracting CAEN V1740 fragments. @@ -371,232 +574,9 @@ void sbndaq::SBNDXARAPUCADecoder::produce(art::Event& e) } } -/** -* @brief Searches for the SPEC-TDC ETRIG timestamp closest to the raw timestamp if any SPEC-TDC ETRIG product is found in the event. -* @param[in] e The event to be processed. -* @param[in] corr_raw_timestamp The corrected raw timestamp from the artdaq::RawEventHeader product. -* @param[in,out] timestamp The closest ETRIG timestamp to the raw timestamp (if found). -* @return A boolean indicating if a valid ETRIG timestamp was found close enough to the raw timestamp. -* @details It searches for the SPEC-TDC products in the event and looks for the ETRIG timestamps. If any ETRIG -* timestamp is found, it checks which one is the closest to the raw timestamp and if it is close enough (i.e. -* within fraw_trig_max_diff) it returns it as output. -*/ -bool sbndaq::SBNDXARAPUCADecoder::get_spec_tdc_etrig_timestamp(art::Event& e, uint64_t corr_raw_timestamp, uint64_t& timestamp) { - bool ett_found = false; - - art::Handle > tdc_handle; - e.getByLabel(fspectdc_product_name, tdc_handle); - - if (fverbose | fdebug_tdc_handle) std::cout << "\n > SBNDXARAPUCADecoder::get_spec_tdc_etrig_timestamp: searching for SPEC-TDC products..." << std::endl; - - // The art::Handle object is not valid. - if (!tdc_handle.isValid()) { - if (fdebug_tdc_handle) std::cout << "\n\tTDC-SPEC handle not valid for " << fspectdc_product_name << "." << std::endl; - - // The art::Handle object is empty. - } else if (tdc_handle->empty()) { - if (fdebug_tdc_handle) std::cout << "\n\tTDC-SPEC handle is empty." << std::endl; - - // The art::Handle object is valid and not empty. - } else { - if (fdebug_tdc_handle) std::cout << " \n\tDecoded TDC-SPEC products found: " << tdc_handle->size() << " products." << std::endl; - - unsigned int num_event_triggers = 0; - unsigned int num_flash_triggers = 0; - - uint64_t min_diff = std::numeric_limits::max(); - uint64_t closest_tdc_timestamp = 0; - - if (fdebug_tdc_handle) std::cout << "\t\tTDC Channel \t TDC Name \t\t TDC Timestamp [ns] \t\t TDC Offset [ns]"; - - // Loops over the SPEC-TDC products in the event. - for (size_t t = 0; t < tdc_handle->size(); t++) { - const uint64_t tdc_timestamp = tdc_handle->at(t).Timestamp(); // Timestamp of the signal [ns]. - const uint32_t tdc_channel = tdc_handle->at(t).Channel(); // Hardware channel. - - if (fdebug_tdc_handle) { - const std::string tdc_name = tdc_handle->at(t).Name(); // Name of channel input. - const uint64_t tdc_offset = tdc_handle->at(t).Offset(); // Channel specific offset [ns]. - std::cout << "\n\t\t\t " << tdc_channel << " \t\t\t" << tdc_name.c_str() - << "\t\t" << print_timestamp(tdc_timestamp) - << "\t\t " << tdc_offset; - } - - // Gets the closest ETRIG timestamp to the raw timestamp and counts the number of ETRIG and FTRIG. - if (tdc_channel == fspectdc_etrig_ch) { - num_event_triggers++; - - //uint64_t diff = (tdc_timestamp < corr_raw_timestamp) ? (corr_raw_timestamp - tdc_timestamp) : (tdc_timestamp - corr_raw_timestamp); - uint64_t diff = abs_difference(tdc_timestamp, corr_raw_timestamp); - if (diff < min_diff) { - min_diff = diff; - closest_tdc_timestamp = tdc_timestamp; - } - - if (fdebug_tdc_handle) std::cout << "\t\t" << diff << " ns away from the raw timestamp."; - - } else if (tdc_channel == fspectdc_ftrig_ch) { - num_flash_triggers++; - } - - } // End SPEC-TDC products extraction loop. - - if (fdebug_tdc_handle) { - std::cout << "\n\tNumber of event triggers (ETRIG): " << num_event_triggers << "." << std::endl; - std::cout << "\tNumber of flash triggers (FTRIG): " << num_flash_triggers << "." << std::endl; - if (closest_tdc_timestamp) { - std::cout << "\n\tClosest ETRIG timestamp to the raw timestamp: " << print_timestamp(closest_tdc_timestamp) << " with a difference of " << min_diff << " ns." << std::endl; - } else { - std::cout << "\n\tNo ETRIG timestamps found." << std::endl; - } - } - - // Checks if the closest ETRIG timestamp is close enough to the raw timestamp. - if (min_diff > fraw_trig_max_diff) { - closest_tdc_timestamp = 0; - } else { - ett_found = true; - timestamp = closest_tdc_timestamp; - } - } - - tdc_handle.removeProduct(); - return ett_found; -} - -/** - * @brief Searches for the PTB HLT timestamp closest to the raw timestamp if any PTB HLT product is found in the event. - * @param[in] e The event to be processed. - * @param[in] corr_raw_timestamp The corrected raw timestamp from the artdaq::RawEventHeader product. - * @param[in,out] timestamp The closest HLT timestamp to the raw timestamp (if found). - * @param[in,out] hlt_code The code of the closest HL trigger to the raw timestamp (if found). - * @return A boolean indicating if a valid HLT timestamp was found close enough to the raw timestamp. - * @details It searches for the PTB products in the event and looks for the HLT timestamps allowed in the fallowed_hl_triggers - * list. If any HLT timestamp is found, it checks which one is the closest to the raw timestamp and if it is close enough - * (i.e. within fraw_trig_max_diff) it returns it as output along with its trigger code. - */ -bool sbndaq::SBNDXARAPUCADecoder::get_ptb_hlt_timestamp(art::Event& e, uint64_t corr_raw_timestamp, uint64_t& timestamp, uint16_t& hlt_code) { - bool hlt_found = false; - - art::Handle > ptb_handle; - e.getByLabel(fptb_product_name, ptb_handle); - - if (fverbose | fdebug_ptb_handle) std::cout << "\n > SBNDXARAPUCADecoder::get_ptb_hlt_timestamp: PTB timing frame selected. Searching for PTB products..." << std::endl; - - // The art::Handle object is not valid. - if (!ptb_handle.isValid()) { - if (fdebug_ptb_handle) std::cout << "\n\tPTB handle not valid for " << fptb_product_name << "." << std::endl; - - // The art::Handle object is empty. - } else if (ptb_handle->empty()) { - if (fdebug_ptb_handle) std::cout << "\n\tPTB handle is empty." << std::endl; - - // The art::Handle object is valid and not empty. - } else { - if (fdebug_ptb_handle) std::cout << " \n\tDecoded PTB products found: " << ptb_handle->size() << " products."; - - uint64_t min_diff = std::numeric_limits::max(); - uint64_t closest_ptb_timestamp = 0; - uint16_t hlt_type = HLT_NOT_FOUND; - - // Loops over the PTB products in the event. - for (size_t i = 0; i < ptb_handle->size(); i++) { - std::vector hlt_triggers = ptb_handle->at(i).GetHLTriggers(); - - if (!hlt_triggers.empty()) { - if (fdebug_ptb_handle) std::cout << "\n\t\t Number of HL triggers: " << hlt_triggers.size() << "." << std::endl; - - // Loops over the HL triggers in the PTB product. - for (size_t j = 0; j < hlt_triggers.size(); j++) { - raw::ptb::Trigger hlt_trig = hlt_triggers.at(j); - uint64_t raw_ptb_timestamp = hlt_trig.timestamp; - uint64_t ptb_timestamp = raw_ptb_timestamp * 20; - std::bitset<32> hlt_word_bitset = std::bitset<32>(hlt_trig.trigger_word); - - if (fdebug_ptb_handle) { - std::cout << "\t\t\t HLT " << j << ": "; - std::cout << "Timestamp: " << print_timestamp(ptb_timestamp); - std::cout << "\t\tTrigger word: [" << hlt_word_bitset << "]"; - } - - uint32_t allowed_hlt = 0; // Random trigger code. In this case it is used as flag to check if an allowed HLT trigger is found. - // Loops over the allowed HL triggers. - for (size_t k = 0; k < fallowed_hl_triggers.size(); k++) { - if (hlt_word_bitset[fallowed_hl_triggers[k]]) { - allowed_hlt = fallowed_hl_triggers[k]; - uint64_t diff = abs_difference(ptb_timestamp, corr_raw_timestamp); - if (diff < min_diff) { - min_diff = diff; - closest_ptb_timestamp = ptb_timestamp; - hlt_type = allowed_hlt; - } - if (fdebug_ptb_handle) std::cout << " - Allowed HLT: " << allowed_hlt << " - " << diff << " ns away from the raw timestamp." << std::endl; - } - } - if ((fdebug_ptb_handle) && !allowed_hlt) std::cout << std::endl; - } - - } else { - if (fdebug_ptb_handle) std::cout << "\n\t\t No allowed HL triggers found." << std::endl; - } - } // End PTB products extraction loop. - - if (fdebug_ptb_handle) { - if (closest_ptb_timestamp) { - std::cout << "\n\t Closest HLT trigger to the raw timestamp: " << print_timestamp(closest_ptb_timestamp) << " - HLT type: " << hlt_type << ". With a difference of " << min_diff << " ns." << std::endl; - } else { - std::cout << "\n\t No allowed HL triggers found." << std::endl; - } - } - - // Checks if the closest HLT timestamp is close enough to the raw timestamp. - if (min_diff > fraw_trig_max_diff) { - closest_ptb_timestamp = 0; - hlt_code = HLT_TOO_FAR; - } else { - hlt_found = true; - timestamp = closest_ptb_timestamp; - hlt_code = hlt_type; - } - - } - - ptb_handle.removeProduct(); - return hlt_found; -} +// =============== Main fragment processing method =============== // -/** - * @brief This function processes a single fragment from a CAEN V1740 and stores the decoded waveforms in the provided product - * container. - * - * @param[in,out] fragment_indices A 1D vector tracking the number of fragments processed for each board. - * @param[in] fragment The input CAEN V1740 fragment containing raw data to be decoded. - * @param[out] prod_wvfms Vector where the decoded waveforms will be stored as raw::OpDetWaveform objects. - * - * @details - * - Identifies the board index corresponding to the fragment ID. - * - Verifies the fragment ID against known boards and ensures it is within a valid range. - * - - * - Decodes the fragment reading raw 32-bit words from the fragment, storing them in a buffer, - * and extracting 12-bit samples. The decoding includes: - * - Accessing metadata for the number of channels and samples. - * - Getting the initial timestamp with respect to the selected timing frame. - * - Binary decoding of the raw waveform. To assign efficiently each sequential sample extracted. These indices formulas are - * applied: - * - The board channel index: - * \f[ - * c = \left( \frac{S}{3} \mod 8 \right) + g \times 8 - * \f] - * - The channel sample index: - * \f[ - * s = (S \mod 3) + \left( \frac{S}{24} \times 3 \right) \mod s_{w}} - * \f] - * Where the group index is computed as \f$ \frac{S}{s_{g}} \f$. - * - Mapping samples to corresponding channels. - * - Populates the output vector (`prod_wvfms`) with decoded waveforms and optionally generates debug waveforms output. - * - */ -void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, std::vector & fragment_indices, const artdaq::Fragment& fragment, std::vector & prod_wvfms) { +void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, uint64_t& first_frag_timestamp, size_t& first_frag_idx, int32_t& first_TTT, int32_t& prev_TTT, std::vector & fragment_indices, const artdaq::Fragment& fragment, std::vector & prod_wvfms, std::vector>& wvfms, bool last_one, size_t& appended_samples) { auto fragment_id = fragment.fragmentID() - ffragment_id_offset; auto it = std::find(fboard_id_list.begin(), fboard_id_list.end(), fragment_id); size_t board_idx; @@ -614,6 +594,7 @@ void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, std::vecto } if (valid_fragment) { + fragment_indices[board_idx]++; if (fverbose) std::cout << "\n > SBNDXARAPUCADecoder::decode_fragment: decoding V1740 CAEN fragment " << fragment_indices[board_idx] << " from the board " << board_idx << " (slot " << fboard_id_list[board_idx] << "):" << std::endl; // =============== Accesses Event metadata and Event header for this fragment =============== // @@ -629,7 +610,8 @@ void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, std::vecto // Gets the number of words of the header and the waveforms. uint32_t num_words_per_event = header.eventSize; - uint32_t num_words_per_header = sizeof(CAENV1740EventHeader) / sizeof(uint32_t); + size_t header_size = sizeof(CAENV1740EventHeader); + uint32_t num_words_per_header = header_size / sizeof(uint32_t); uint32_t num_words_per_wvfms = (num_words_per_event - num_words_per_header); uint32_t num_bits_per_all_wvfms = num_words_per_wvfms * BITS_PER_WORD; @@ -637,14 +619,15 @@ void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, std::vecto uint32_t num_remaining_bits = num_bits_per_all_wvfms % BITS_PER_SAMPLE; uint32_t num_samples_per_wvfm = num_samples_per_all_wvfms / num_channels; uint32_t num_samples_per_group = num_samples_per_wvfm * NUM_CHANNELS_PER_GROUP; + uint32_t num_nominal_samples_per_wvfm = metadata->nSamples; if (fverbose | fdebug_waveforms) { - if (metadata->nSamples == num_samples_per_wvfm) { + if (num_nominal_samples_per_wvfm == num_samples_per_wvfm) { std::cout << " > SBNDXARAPUCADecoder::decode_fragment: [NOMINAL FRAGMENT] " << num_samples_per_wvfm << " samples/waveform." << " (" << num_samples_per_group << " samples per group - 8 channels per group -)." << std::endl; } else { std::cout << " > SBNDXARAPUCADecoder::decode_fragment: [EXTENDED FRAGMENT] " << num_samples_per_wvfm << " samples/waveform." << " (" << num_samples_per_group << " samples per group - 8 channels per group -)." << std::endl; } - std::cout << " > SBNDXARAPUCADecoder::decode_fragment: nominal number of samples per waveform: " << metadata->nSamples << "." << std::endl; + std::cout << " > SBNDXARAPUCADecoder::decode_fragment: nominal number of samples per waveform: " << num_nominal_samples_per_wvfm << "." << std::endl; std::cout << " > SBNDXARAPUCADecoder::decode_fragment: number of words for this fragment: " << num_words_per_event << " (Header: " << num_words_per_header << ", Waveform: " << num_words_per_wvfms << ") words." << std::endl; } @@ -658,169 +641,308 @@ void sbndaq::SBNDXARAPUCADecoder::decode_fragment(uint64_t timestamp, std::vecto // =============== Extracts timing information for this fragment =============== // - // Gets the timing information of the CAEN fragment. - int64_t pulse_duration_ns = num_samples_per_wvfm * fns_per_sample; // ns. uint64_t frag_timestamp = fragment.timestamp(); // ns. - int64_t frag_timestamp_s = frag_timestamp / NANOSEC_IN_SEC; // s. - int64_t frag_timestamp_ns = frag_timestamp % NANOSEC_IN_SEC; // ns. - uint32_t TTT_ticks = header.triggerTime(); int64_t TTT_end_ns = TTT_ticks * NANOSEC_PER_TICK; // ns. - // Gets the full TTT timestamp. - uint64_t full_TTT = 0; - // If the fragment timestamp is greater than the TTT end timestamp, it means that rollover occurred. - if (frag_timestamp_ns > TTT_end_ns) { - if (fverbose | fdebug_timing) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: CAEN TTT rollover occurred w.r.t. the fragment timestamp (FTS)." << std::endl; - full_TTT = (frag_timestamp_s + 1) * NANOSEC_IN_SEC + TTT_end_ns; - } else { - full_TTT = frag_timestamp_s * NANOSEC_IN_SEC + TTT_end_ns; + if (fdebug_timing) { + std::cout << "\t\t TTT header.extendedTriggerTime() [TTT_ticks] = " << header.extendedTriggerTime() << " ticks. \t TTT_end_ns = " << print_timestamp(header.extendedTriggerTime() * NANOSEC_PER_TICK) << "." << std::endl; + std::cout << "\t\t TTT header.triggerTimeRollOver(): " << header.triggerTimeRollOver() << std::endl; } - int64_t ref_timestamp = 0; - + //// =============== Start decoding the waveforms =============== // + std::vector > fragment_wvfms(num_channels, std::vector(num_samples_per_wvfm, 0)); + decode_waveforms(fragment, fragment_wvfms, header_size, num_channels, num_samples_per_wvfm, num_words_per_wvfms, num_samples_per_group); + double ini_wvfm_timestamp = 0; double end_wvfm_timestamp = 0; + + if (fverbose) { + std::cout << " > SBNDXARAPUCADecoder::decode_fragment: processing the decoded fragment and combines the extended ones to their nominal ones if needed." << std::endl; + } - // If an ETRIG or HLT timestamp was found it restarts the time from it. Otherwise the CAEN time frame is assigned. - if (factive_timing_frame != CAEN_ONLY_TIMING) { - ref_timestamp = signed_difference(full_TTT, timestamp); // ns. + if (fcombine_ext_frag) { - ini_wvfm_timestamp = (ref_timestamp - pulse_duration_ns) * NANOSEC_TO_MICROSEC; // us. - end_wvfm_timestamp = ref_timestamp * NANOSEC_TO_MICROSEC; // us. - } else { - ref_timestamp = full_TTT; // ns. + // Combination variables and checks. - ini_wvfm_timestamp = ((ref_timestamp - pulse_duration_ns) % NANOSEC_IN_SEC) * NANOSEC_TO_MICROSEC; // ns. - end_wvfm_timestamp = (ref_timestamp % NANOSEC_IN_SEC) * NANOSEC_TO_MICROSEC; // us. - } + bool is_nominal = (num_nominal_samples_per_wvfm == num_samples_per_wvfm); + bool is_first = (fragment_indices[board_idx] == 0); - if (fdebug_timing) { - std::cout << std::fixed << std::setprecision(0); - std::cout << "\t\t ns/tick = " << NANOSEC_PER_TICK << ", ns/sample = " << fns_per_sample << std::endl; - std::cout << "\t\t TTT header.TriggerTime() [TTT_ticks] = " << TTT_ticks << " ticks. \t TTT_end_ns = " << print_timestamp(TTT_end_ns) << "." << std::endl; - std::cout << "\t\t TTT header.extendedTriggerTime() [TTT_ticks] = " << header.extendedTriggerTime() << " ticks. \t TTT_end_ns = " << print_timestamp(header.extendedTriggerTime() * NANOSEC_PER_TICK) << "." << std::endl; - std::cout << "\t\t TTT header.triggerTimeRollOver(): " << header.triggerTimeRollOver() << std::endl; - std::cout << "\t\t Full Fragment timestamp: " << print_timestamp(frag_timestamp) << " = " << frag_timestamp_s << " s " << frag_timestamp_ns << " ns." << std::endl; - std::cout << "\t\t Full TTT - fragment timestamp = "<< abs_difference(full_TTT, frag_timestamp) << " ns." << " Post-percent: " << (double(abs_difference(full_TTT, frag_timestamp)) / double(pulse_duration_ns)) * 100 << "%." << std::endl; - if (factive_timing_frame == SPEC_TDC_TIMING) { - std::cout << "\t ETRIG (SPEC-TDC) timestamp of the fragment: " << std::endl; - std::cout << "\t\t Full UTC ETRIG timestamp: " << print_timestamp(timestamp) << "." << std::endl; - std::cout << "\t\t ETRIG SPEC-TDC difference applied to the CAEN frame (full timestamps): " << print_timestamp(full_TTT) << " - " << print_timestamp(timestamp) << " = " << ref_timestamp << " ns." << std::endl; - } else if (factive_timing_frame == PTB_TIMING) { - std::cout << "\t HLT ETRIG (PTB) timestamp of the fragment: " << std::endl; - std::cout << "\t\t Full UTC HLT ETRIG timestamp: " << print_timestamp(timestamp) << "." << std::endl; - std::cout << "\t\t HLT ETRIG (PTB) difference applied to the CAEN frame (full timestamps): " << print_timestamp(full_TTT) << " - " << print_timestamp(timestamp) << " = " << ref_timestamp << " ns." << std::endl; - } else if (factive_timing_frame == CAEN_ONLY_TIMING) { - std::cout << "\t CAEN trigger timestamp (TTT) of the fragment: " << std::endl; - std::cout << "\t\t Full UTC TTT timestamp: " << print_timestamp(full_TTT) << " = " << full_TTT / NANOSEC_IN_SEC << " s " << TTT_end_ns << " ns." << std::endl; + int32_t frag_length = num_samples_per_wvfm * fns_per_sample; // ns. + + // - Computation of the difference w.r.t. previous TTT difference. + int32_t TTT_dif = 0; + if (prev_TTT > TTT_end_ns) { + if (fverbose | fdebug_timing) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: CAEN TTT rollover occurred." << std::endl; + TTT_dif = prev_TTT - (NANOSEC_IN_SEC - TTT_end_ns); + } else { + TTT_dif = TTT_end_ns - prev_TTT; } - } - if (fverbose | fdebug_timing) { - std::cout << std::fixed << std::setprecision(3) << std::endl; - if (factive_timing_frame == SPEC_TDC_TIMING) { - std::cout << " > SBNDXARAPUCADecoder::decode_fragment: SPEC-TDC time window of " << end_wvfm_timestamp - ini_wvfm_timestamp << " us: [" << ini_wvfm_timestamp << ", " << end_wvfm_timestamp << "] us." << std::endl; - } else if (factive_timing_frame == PTB_TIMING) { - std::cout << " > SBNDXARAPUCADecoder::decode_fragment: PTB time window of " << end_wvfm_timestamp - ini_wvfm_timestamp << " us: [" << ini_wvfm_timestamp << ", " << end_wvfm_timestamp << "] us." << std::endl; - } else { // CAEN_ONLY_TIMING - std::cout << " > SBNDXARAPUCADecoder::decode_fragment: CAEN time window of " << end_wvfm_timestamp - ini_wvfm_timestamp << " us: [" << ini_wvfm_timestamp << ", " << end_wvfm_timestamp << "] us." << std::endl; + // - Computation of the difference between TTTs difference and the length of the fragment in ns (to detect jittering between fragments). + int32_t J = TTT_dif - frag_length; + + bool dump_comb_wvfms = !is_first && (is_nominal || (!is_nominal && (TTT_dif > (frag_length + fallowed_jittering)))); + + if (is_first) { + first_frag_idx = fragment_indices[board_idx]; + first_TTT = TTT_end_ns; + first_frag_timestamp = frag_timestamp; } - } - // =============== Start decoding the waveforms =============== // - if (fverbose) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: binary decoding of the waveforms starting... " << std::endl; - - std::vector > wvfms(num_channels, std::vector(num_samples_per_wvfm, 0)); - - // Absolute sample number [0, TOTAL_NUM_SAMPLES] where TOTAL_NUM_SAMPLES is the total number of samples stored for an event. - uint32_t S = 0; - // Buffer variables. - uint64_t buffer = 0; - uint32_t bits_in_buffer = 0; - - // Data pointer to the beggining of the waveforms stores in the event. - const uint32_t* data_ptr = reinterpret_cast(fragment.dataBeginBytes() + sizeof(CAENV1740EventHeader)); - // Accesses each word, stores it in the buffer and then the samples are extracted from the buffer. - for (size_t j = 0; j < num_words_per_wvfms; j++) { - uint64_t word = read_word(data_ptr); - - // Adds the new word to the buffer and increments the number of bits stored in it. - if (fdebug_buffer) std::cout << buffer << "[word: " << word << "]" << std::endl; - buffer |= word << bits_in_buffer; - bits_in_buffer += BITS_PER_WORD; // bytes * 8 bits/byte - if (fdebug_buffer) std::cout << " +" << buffer << " [bits in buffer: "<< bits_in_buffer << "]" << std::endl; - - // Obtains 12-bit sequences from the buffer and assigns each sample to the channel and channel sample it belongs to. - while (bits_in_buffer >= BITS_PER_SAMPLE) { - // Computes board channel, channel sample and group channel and assigns the sample to those indices. - uint32_t g = (S / num_samples_per_group); // Group index. - uint32_t c = ((S / NUM_CONSECUTIVE_SAMPLES) % NUM_CHANNELS_PER_GROUP) + g * NUM_GROUPS; // Channel index. - uint32_t s = (S % NUM_CONSECUTIVE_SAMPLES) + ((S / NUM_SAMPLES_PER_ROUND) * NUM_CONSECUTIVE_SAMPLES) % num_samples_per_wvfm; // Sample/channel index. - uint16_t sample = get_sample(buffer, BITS_PER_SAMPLE - 1, 0); - wvfms[c][s] = sample; - if (fdebug_waveforms) std::cout << "\tSample: " << sample << "\tg: " << g << "\tch: " << c << "\ts:" << s << "\tS: " << S << std::endl; - - // Updates the buffer status removing the read bits and decreasing the number of bits stored in it. - buffer >>= BITS_PER_SAMPLE; - bits_in_buffer -= BITS_PER_SAMPLE; - if (fdebug_buffer) std::cout << " -" << buffer << " [bits in buffer: "<< bits_in_buffer << "]" << std::endl; - - // Increments the absolute sample step. - S++; + if (dump_comb_wvfms) { + shift_time(TTT_ticks, first_TTT, first_frag_timestamp, timestamp, num_nominal_samples_per_wvfm, ini_wvfm_timestamp, end_wvfm_timestamp, appended_samples); + dump_waveforms(prod_wvfms, wvfms, first_frag_idx, board_idx, num_channels, ini_wvfm_timestamp, end_wvfm_timestamp); + first_frag_idx = fragment_indices[board_idx]; + first_TTT = TTT_end_ns; + first_frag_timestamp = frag_timestamp; + if (fdebug_jittering) std::cout << "NC - "; + } else { + if (fdebug_jittering) std::cout << "C - "; } - } + append_waveforms(wvfms, fragment_wvfms, num_channels); + appended_samples = wvfms[0].size() - num_nominal_samples_per_wvfm; + + if (last_one) { + if (fdebug_extended_fragments) std::cout << "\t\t LAST fragment " << std::endl; + shift_time(TTT_ticks, first_TTT, first_frag_timestamp, timestamp, num_nominal_samples_per_wvfm, ini_wvfm_timestamp, end_wvfm_timestamp, appended_samples); + dump_waveforms(prod_wvfms, wvfms, first_frag_idx, board_idx, num_channels, ini_wvfm_timestamp, end_wvfm_timestamp); + } + + if (fdebug_jittering) { + if (is_first) { + std::cout << "\t\t FIRST "; + } + if (is_nominal) { + int32_t TTT_ini_ns = TTT_end_ns - static_cast(num_nominal_samples_per_wvfm * fns_per_sample); + std::cout << "NOMINAL fragment " << fboard_id_list[board_idx] << " - " << fragment_indices[board_idx] << " (" << num_samples_per_wvfm << "): [" << TTT_ini_ns << ", " << TTT_end_ns << "] diff with prev: " << TTT_dif << " ns, length: " << frag_length << " ns, diff_with_length (jitt): " << J << std::endl; + } else { + int32_t TTT_ini_ns = TTT_end_ns - static_cast(num_nominal_samples_per_wvfm * fns_per_sample); + std::cout << "EXTENDED fragment " << fboard_id_list[board_idx] << " - " << fragment_indices[board_idx] << " (" << num_samples_per_wvfm << "): [" << TTT_ini_ns << ", " << TTT_end_ns << "] diff with prev: " << TTT_dif << " ns, length: " << frag_length << " ns, diff_with_length (jitt): " << J << std::endl; + } + } + + // WARNINGS + if (TTT_dif < frag_length - fallowed_jittering) { + if (!is_first) { + if (!is_nominal) std::cout << "\t\t WARNING: EXTENDED fragment too close! This fragment is being combined." << std::endl; + if (is_nominal) std::cout << "\t\t WARNING: NOMINAL fragment too close! This fragment is not being combined." << std::endl; + } + } + + prev_TTT = TTT_end_ns; - // The decoded waveforms are dumped into two products: - // - A xarapucadecoder-art.root file with the OpDetWaveforms as the product of this producer for further analysis. - // - A decoder_hist.root file gathering a waveform histograms. - if (fverbose) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: binary decoding complete, dumping products..." << std::endl; + } else { // Combination of extended fragments disabled. + first_frag_idx = fragment_indices[board_idx]; + shift_time(TTT_ticks, TTT_end_ns, frag_timestamp, timestamp, num_samples_per_wvfm, ini_wvfm_timestamp, end_wvfm_timestamp, appended_samples); + dump_waveforms(prod_wvfms, fragment_wvfms, first_frag_idx, board_idx, num_channels, ini_wvfm_timestamp, end_wvfm_timestamp); + } + } + +} + +// =============== Timing functions =============== // + +void sbndaq::SBNDXARAPUCADecoder::shift_time(uint64_t TTT_ticks, int64_t TTT_end_ns, uint64_t frag_timestamp, uint64_t timestamp, uint32_t num_samples_per_wvfm, double& ini_wvfm_timestamp, double& end_wvfm_timestamp, size_t& appended_samples) { + + int64_t pulse_duration_ns = num_samples_per_wvfm * fns_per_sample; // ns. + int64_t frag_timestamp_s = frag_timestamp / NANOSEC_IN_SEC; // s. + int64_t frag_timestamp_ns = frag_timestamp % NANOSEC_IN_SEC; // ns. + + // Gets the UNIX TTT timestamp. + uint64_t TTT_UNIX_timestamp = 0; + // If the fragment timestamp is greater than the TTT end timestamp, it means that rollover occurred. + if (frag_timestamp_ns > TTT_end_ns) { + if (fverbose | fdebug_timing) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: CAEN TTT rollover occurred w.r.t. the fragment timestamp (FTS)." << std::endl; + TTT_UNIX_timestamp = (frag_timestamp_s + 1) * NANOSEC_IN_SEC + TTT_end_ns; + } else { + TTT_UNIX_timestamp = frag_timestamp_s * NANOSEC_IN_SEC + TTT_end_ns; + } + + int64_t ref_timestamp = 0; + int64_t appended_ns = static_cast(appended_samples) * static_cast(fns_per_sample); // this is just for the debug output. + + // If an frameshift timestamp was found it restarts the time from it. Otherwise the CAEN time frame is assigned. + if (factive_timing_frame != CAEN_ONLY_TIMING) { + ref_timestamp = signed_difference_UTC_timestamp(TTT_UNIX_timestamp, timestamp); // ns. + ini_wvfm_timestamp = (ref_timestamp - pulse_duration_ns) * NANOSEC_TO_MICROSEC; // us. + end_wvfm_timestamp = (ref_timestamp + appended_ns) * NANOSEC_TO_MICROSEC; // us. + } else { + ref_timestamp = TTT_UNIX_timestamp; // ns. + + ini_wvfm_timestamp = ((ref_timestamp - pulse_duration_ns) % NANOSEC_IN_SEC) * NANOSEC_TO_MICROSEC; // ns. + end_wvfm_timestamp = ((ref_timestamp + appended_ns)% NANOSEC_IN_SEC) * NANOSEC_TO_MICROSEC; // us. + } + appended_samples = 0; + if (fdebug_timing) { + std::cout << std::fixed << std::setprecision(0); + std::cout << "\t\t ns/tick = " << NANOSEC_PER_TICK << ", ns/sample = " << fns_per_sample << std::endl; + std::cout << "\t\t TTT header.TriggerTime() [TTT_ticks] = " << TTT_ticks << " ticks. \t TTT_end_ns = " << print_timestamp(TTT_end_ns) << "." << std::endl; + std::cout << "\t\t Full Fragment timestamp: " << print_timestamp(frag_timestamp) << " = " << frag_timestamp_s << " s " << frag_timestamp_ns << " ns." << std::endl; + std::cout << "\t\t Full TTT - fragment timestamp = "<< abs_difference(TTT_UNIX_timestamp, frag_timestamp) << " ns." << " Post-percent: " << (double(abs_difference(TTT_UNIX_timestamp, frag_timestamp)) / double(pulse_duration_ns)) * 100 << "%." << std::endl; + if (factive_timing_frame == FRAMESHIFT_TIMING) { + std::cout << "\t FRAMESHIFT timestamp of the fragment: " << std::endl; + std::cout << "\t\t Full UTC FRAMESHIFT timestamp: " << print_timestamp(timestamp) << "." << std::endl; + std::cout << "\t\t FRAMESHIFT difference applied to the CAEN frame (UNIX timestamps): " << print_timestamp(TTT_UNIX_timestamp) << " - " << print_timestamp(timestamp) << " = " << ref_timestamp << " ns." << std::endl; + } else { + std::cout << "\t CAEN trigger timestamp (TTT) of the fragment: " << std::endl; + std::cout << "\t\t Full UTC TTT timestamp: " << print_timestamp(TTT_UNIX_timestamp) << " = " << TTT_UNIX_timestamp / NANOSEC_IN_SEC << " s " << TTT_end_ns << " ns." << std::endl; + } + } + + if (fverbose | fdebug_timing) { + std::cout << std::fixed << std::setprecision(3); + if (factive_timing_frame == FRAMESHIFT_TIMING) { + std::cout << " > SBNDXARAPUCADecoder::decode_fragment: FRAMESHIFT time window of " << end_wvfm_timestamp - ini_wvfm_timestamp << " us: [" << ini_wvfm_timestamp << ", " << end_wvfm_timestamp << "] us." << std::endl; + } else { // CAEN_ONLY_TIMING + std::cout << " > SBNDXARAPUCADecoder::decode_fragment: CAEN time window of " << end_wvfm_timestamp - ini_wvfm_timestamp << " us: [" << ini_wvfm_timestamp << ", " << end_wvfm_timestamp << "] us." << std::endl; + } + std::cout << " > SBNDXARAPUCADecoder::decode_fragment: TTT_end_ticks = " << TTT_ticks << " ticks. \t TTT_end_ns = " << print_timestamp(TTT_end_ns) << "." << std::endl; + } +} + +std::unique_ptr sbndaq::SBNDXARAPUCADecoder::get_frameshift_timestamp(art::Event& e, uint64_t& timestamp) { + std::unique_ptr fs_timing_ref_info = nullptr; + + art::Handle frameshiftHandle; + + if (fverbose || fdebug_frameshift_handle) std::cout << "\n > SBNDXARAPUCADecoder::get_frameshift_timestamp: Frameshift timing frame selected. Searching for Frameshift products..." << std::endl; + + if (e.getByLabel(fframeshift_product_name, frameshiftHandle)) { - uint32_t num_debug_wvfms; + uint16_t timing_type = frameshiftHandle->TimingTypeDefault(); + uint16_t timing_channel = frameshiftHandle->TimingChannelDefault(); - if (fstore_debug_waveforms == -1) { - num_debug_wvfms = num_channels; + if (timing_type != sbnd::timing::kNoShiftType) { // Only if timing type is valid, the timestamp is considered valid. + timestamp = frameshiftHandle->FrameDefault(); + fs_timing_ref_info = std::make_unique(timing_type, timing_channel); } else { - num_debug_wvfms = std::min(num_channels, fstore_debug_waveforms); + timestamp = 0; } + + } else { + if (fdebug_frameshift_handle) std::cout << "\n\tFrameshift handle not valid for " << fframeshift_product_name << "." << std::endl; + throw std::runtime_error("XA Decoder: Frame Shift Info object is invalid, check data quality"); + } - uint32_t ch; + return fs_timing_ref_info; +} + +// =============== Decodes the waveforms =============== // - for (ch = 0; ch < num_debug_wvfms; ch++) { - save_prod_wvfm(board_idx, ch, ini_wvfm_timestamp, wvfms, prod_wvfms); - save_debug_wvfm(board_idx, fragment_indices[board_idx], ch, ini_wvfm_timestamp, end_wvfm_timestamp, wvfms); +void sbndaq::SBNDXARAPUCADecoder::decode_waveforms(const artdaq::Fragment& fragment, std::vector>& wvfms, size_t header_size, uint32_t num_channels, uint32_t num_samples_per_wvfm, uint32_t num_words_per_wvfms, uint32_t num_samples_per_group) { + // =============== Start decoding the waveforms =============== // + if (fverbose) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: binary decoding of the waveforms starting... " << std::endl; + + // Absolute sample number [0, TOTAL_NUM_SAMPLES] where TOTAL_NUM_SAMPLES is the total number of samples stored for an event. + uint32_t S = 0; + // Buffer variables. + uint64_t buffer = 0; + uint32_t bits_in_buffer = 0; + + // Data pointer to the beggining of the waveforms stores in the event. + const uint32_t* data_ptr = reinterpret_cast(fragment.dataBeginBytes() + header_size); + // Accesses each word, stores it in the buffer and then the samples are extracted from the buffer. + for (size_t j = 0; j < num_words_per_wvfms; j++) { + uint64_t word = read_word(data_ptr); + + // Adds the new word to the buffer and increments the number of bits stored in it. + if (fdebug_buffer) std::cout << buffer << "[word: " << word << "]" << std::endl; + buffer |= word << bits_in_buffer; + bits_in_buffer += BITS_PER_WORD; // bytes * 8 bits/byte + if (fdebug_buffer) std::cout << " +" << buffer << " [bits in buffer: "<< bits_in_buffer << "]" << std::endl; + + // Obtains 12-bit sequences from the buffer and assigns each sample to the channel and channel sample it belongs to. + while (bits_in_buffer >= BITS_PER_SAMPLE) { + // Computes board channel, channel sample and group channel and assigns the sample to those indices. + uint32_t g = (S / num_samples_per_group); // Group index. + uint32_t c = ((S / NUM_CONSECUTIVE_SAMPLES) % NUM_CHANNELS_PER_GROUP) + g * NUM_GROUPS; // Channel index. + uint32_t s = (S % NUM_CONSECUTIVE_SAMPLES) + ((S / NUM_SAMPLES_PER_ROUND) * NUM_CONSECUTIVE_SAMPLES) % num_samples_per_wvfm; // Sample/channel index. + uint16_t sample = get_sample(buffer, BITS_PER_SAMPLE - 1, 0); + wvfms[c][s] = sample; + if (fdebug_waveforms) std::cout << "\tSample: " << sample << "\tg: " << g << "\tch: " << c << "\ts:" << s << "\tS: " << S << std::endl; + + // Updates the buffer status removing the read bits and decreasing the number of bits stored in it. + buffer >>= BITS_PER_SAMPLE; + bits_in_buffer -= BITS_PER_SAMPLE; + if (fdebug_buffer) std::cout << " -" << buffer << " [bits in buffer: "<< bits_in_buffer << "]" << std::endl; + + // Increments the absolute sample step. + S++; } + } +} + +uint16_t sbndaq::SBNDXARAPUCADecoder::get_sample(uint64_t buffer, uint32_t msb, uint32_t lsb) { + uint64_t mask = (1U << (msb - lsb + 1)) - 1; + uint64_t sample = buffer >> lsb; + return sample & mask; +} + +uint32_t sbndaq::SBNDXARAPUCADecoder::read_word(const uint32_t* & data_ptr) { + uint32_t word = *data_ptr; + data_ptr += 1; + return word; +} - for (;ch < num_channels; ch++) { - save_prod_wvfm(board_idx, ch, ini_wvfm_timestamp, wvfms, prod_wvfms); +unsigned int sbndaq::SBNDXARAPUCADecoder::get_channel_id(unsigned int board, unsigned int board_channel) { + unsigned int channel_id = fboard_id_list[board] * 100 + board_channel; + return channel_id; +} + +// =============== Combines the waveforms from extended fragments =============== // + +void sbndaq::SBNDXARAPUCADecoder::append_waveforms(std::vector>& wvfms, const std::vector>& fragment_wvfms, uint32_t num_channels) { + if (fdebug_extended_fragments) std::cout << " > SBNDXARAPUCADecoder::append_waveforms: combining waveforms from extended fragments..." << std::endl; + if (wvfms.empty()) { + if (fdebug_extended_fragments) { + std::cout << "\t\t Empty waveforms, resizing to " << num_channels << " channels." << std::endl; } + wvfms.resize(num_channels); + if (fverbose | fdebug_extended_fragments) { + std::cout << "\t\t Waveforms size BEFORE combining = " << wvfms.size() << " x " << wvfms[0].size() << std::endl; + } + } + for (uint32_t ch = 0; ch < num_channels; ch++) { + wvfms[ch].insert(wvfms[ch].end(), fragment_wvfms[ch].begin(), fragment_wvfms[ch].end()); + } + if (fverbose | fdebug_extended_fragments) { + std::cout << "\t\t Waveforms size AFTER combining = " << wvfms.size() << " x " << wvfms[0].size() << std::endl; + } +} + +// =============== Dumps the decoded waveforms =============== // + +void sbndaq::SBNDXARAPUCADecoder::dump_waveforms(std::vector & prod_wvfms, std::vector>& wvfms, size_t first_frag_idx, size_t board_idx, uint32_t num_channels, double ini_wvfm_timestamp, double end_wvfm_timestamp) { + + // The decoded waveforms are dumped into two products: + // - A xarapucadecoder-art.root file with the OpDetWaveforms as the product of this producer for further analysis. + // - A decoder_hist.root file gathering a waveform histograms. + if (fverbose) std::cout << " > SBNDXARAPUCADecoder::decode_fragment: binary decoding complete, dumping products..." << std::endl; + + uint32_t num_debug_wvfms; - fragment_indices[board_idx]++; + if (fstore_debug_waveforms == -1) { + num_debug_wvfms = num_channels; + } else { + num_debug_wvfms = std::min(num_channels, fstore_debug_waveforms); + } + + uint32_t ch; + + for (ch = 0; ch < num_debug_wvfms; ch++) { + save_prod_wvfm(board_idx, ch, ini_wvfm_timestamp, wvfms, prod_wvfms); + save_debug_wvfm(board_idx, first_frag_idx, ch, ini_wvfm_timestamp, end_wvfm_timestamp, wvfms); } + + for (;ch < num_channels; ch++) { + save_prod_wvfm(board_idx, ch, ini_wvfm_timestamp, wvfms, prod_wvfms); + } + + wvfms.clear(); + } -/** - * @brief Converts a production waveform from raw ADC data into a `raw::OpDetWaveform` object, - * assigning it a global channel ID and timestamp, and appends it to the output collection. - * - * @param[in] b The board index (position in the list of boards). - * @param[in] ch The channel index (channel number from which the waveform is extracted). - * @param[in] ini_wvfm_timestamp The initial timestamp of the waveform in microseconds. - * @param[in] wvfms A 2D vector containing the waveforms. Each inner vector corresponds to a channel, - * and each element of the inner vector represents a sample (ADC count). - * @param[out] prod_wvfms A reference to the vector where the produced `raw::OpDetWaveform` objects - * will be stored. - * - * @details - * The function performs the following steps: - * 1. Determines the global channel ID for the specified board and channel using `get_channel_id`. - * 2. Creates a `raw::OpDetWaveform` object with the initial timestamp, global channel ID, and waveform data. - * 3. Appends the `raw::OpDetWaveform` object to the provided `prod_wvfms` collection. - * - * @pre ini_wvfm_timestamp is assumed to be given in microseconds when the timing frame is not CAEN_ONLY_TIMING. - * - * @see raw::OpDetWaveform - */ void sbndaq::SBNDXARAPUCADecoder::save_prod_wvfm(size_t board_idx, size_t ch, double ini_wvfm_timestamp, const std::vector > & wvfms, std::vector & prod_wvfms) { unsigned int channel_id = get_channel_id(board_idx, ch); raw::OpDetWaveform waveform(ini_wvfm_timestamp, channel_id, wvfms[ch]); @@ -828,29 +950,6 @@ void sbndaq::SBNDXARAPUCADecoder::save_prod_wvfm(size_t board_idx, size_t ch, do prod_wvfms.push_back(waveform); } -/** - * @brief Saves debug waveform as a histogram for debugging purposes, based on the waveform data provided. - * The histogram is labeled with event number, fragment, board slot ID and channel information. - * - * @param[in] board_idx The board index (position in the list of boards). - * @param[in] fragment_idx The fragment index (order in decoding). - * @param[in] ch The channel index (channel number from which the waveform is extracted). - * @param[in] ini_wvfm_timestamp The initial timestamp of the waveform. - * @param[in] end_wvfm_timestamp The final timestamp of the waveform. - * @param[in] wvfms A 2D vector containing the waveforms. Each inner vector corresponds to a channel, - * and each element of the inner vector represents a sample (ADC count). - * - * @details - * The function generates a histogram with: - * - X-axis representing time in microseconds if the timing frame is SPEC_TDC_TIMING or PTB_TIMING, otherwise in samples. - * - Y-axis representing ADC counts. - * The histogram is stored using ROOT's `TFileService`. - * - * @pre ini_wvfm_timestamp and end_wvfm_timestamp are assumed to be given in microseconds when the timing frame is not CAEN_ONLY_TIMING. - * - * @see TH1I, TFileService - */ - void sbndaq::SBNDXARAPUCADecoder::save_debug_wvfm(size_t board_idx, size_t fragment_idx, int ch, double ini_wvfm_timestamp, double end_wvfm_timestamp, const std::vector > & wvfms) { std::stringstream hist_name(""); hist_name << "Event " << fevent_counter << " CH " << ch << " [frag " << fragment_idx << ", board " << board_idx << " (slot " << fboard_id_list[board_idx] << ") ] waveform"; @@ -858,11 +957,9 @@ void sbndaq::SBNDXARAPUCADecoder::save_debug_wvfm(size_t board_idx, size_t fragm if (factive_timing_frame != CAEN_ONLY_TIMING) { TH1I* hist = tfs->make(hist_name.str().c_str(), hist_name.str().c_str(), wvfms[ch].size(), ini_wvfm_timestamp, end_wvfm_timestamp); hist->GetYaxis()->SetTitle("ADCs"); - if (factive_timing_frame == SPEC_TDC_TIMING) { - hist->GetXaxis()->SetTitle("Time w.r.t. ETRIG (SPEC-TDC) [us]"); - } else if (factive_timing_frame == PTB_TIMING) { - hist->GetXaxis()->SetTitle("Time w.r.t. HLT ETRIG (PTB) [us]"); - } + if (factive_timing_frame == FRAMESHIFT_TIMING) { + hist->GetXaxis()->SetTitle("Time w.r.t. Frameshift reference [us]"); + } for (size_t i = 0; i < wvfms[ch].size(); i++) { hist->SetBinContent(i+1, wvfms[ch][i]); } @@ -877,72 +974,9 @@ void sbndaq::SBNDXARAPUCADecoder::save_debug_wvfm(size_t board_idx, size_t fragm } -/** - * @brief Extract a sample from a 64-bit buffer using the specified bit positions. - * - * @param[in] buffer An unsigned 64-bit integer which represents a temporal buffer for the read words and where the samples are extracted from. - * @param[in] msb An unsigned 32-bit integer representing the most significative bit (MSB) where the readout from the buffer paramter. - * @param[in] lsb An unsigned 32-bit integer representing the less significative bit (LSB) from we end read - * - * @details The function shifts the buffer to the right by the number of positions specified by `lsb` so that the least significant bit of the - * sample aligns with bit 0. It then applies a mask to isolate the bits between `lsb` and `msb`, inclusive. - * - * @return The extracted sample as a 16-bit unsigned integer. - */ -uint16_t sbndaq::SBNDXARAPUCADecoder::get_sample(uint64_t buffer, uint32_t msb, uint32_t lsb) { - uint64_t mask = (1U << (msb - lsb + 1)) - 1; - uint64_t sample = buffer >> lsb; - return sample & mask; -} - -/** - * @brief Read a 32-bit word from the data pointer and advances the pointer. - * - * @param[in, out] data_ptr A reference to a pointer pointing to the current position in the data. - * - * @details This function retrieves a 32-bit word from the memory location pointed to by `data_ptr`. After reading, it advances `data_ptr` to - * the next 32-bit word location. - * - * @return The 32-bit word read from the location pointed to by `data_ptr`. - */ -uint32_t sbndaq::SBNDXARAPUCADecoder::read_word(const uint32_t* & data_ptr) { - uint32_t word = *data_ptr; - data_ptr += 1; - return word; -} - -/** - * @brief Generates a unique global channel identifier using the board slot and the channel number of that board. - * - * @param[in] board Index of the board in `fboard_id_list` from which to derive the board slot. - * @param[in] board_channel The specific channel number on the given board. - * - * @details This function computes a `channel_id` by combining the board slot and the specific - * channel number on that board. - * The unique identifier `channel_id` (\f$ CH_{ID} $\f) is computed as follows: - * \f[ - * CH\_{ID} = B\_{ID} \times 100 + CH\_{B} - * \f] - * - * Where: - * - \f$ B\_{ID} \f$ is the fragment ID retrieved from `fboard_id_list` based on the slot. - * - \f$ CH\_B \f$ is the channel number on that board. - * - * @return A unique identifier for the specified channel as an unsigned integer. - */ -unsigned int sbndaq::SBNDXARAPUCADecoder::get_channel_id(unsigned int board, unsigned int board_channel) { - unsigned int channel_id = fboard_id_list[board] * 100 + board_channel; - return channel_id; -} +// ==================== Auxiliary functions ==================== // -/** - * @brief Returns the signed difference between two timestamps (t1 - t2). - * @param[in] t1 The first timestamp in nanoseconds. - * @param[in] t2 The second timestamp in nanoseconds. - * @return The signed difference between the two timestamps (t1 - t2) in nanoseconds. - * @pre The timestamp is given in nanoseconds and the maximum value is less than 2^63 ns. - */ -int64_t sbndaq::SBNDXARAPUCADecoder::signed_difference(uint64_t t1, uint64_t t2) { +int64_t sbndaq::SBNDXARAPUCADecoder::signed_difference_UTC_timestamp(uint64_t t1, uint64_t t2) { int64_t signed_diff = 0; if (t2 < t1) { @@ -954,13 +988,6 @@ int64_t sbndaq::SBNDXARAPUCADecoder::signed_difference(uint64_t t1, uint64_t t2) return signed_diff; } -/** - * @brief Returns the absolute difference between two timestamps. - * @param[in] t1 The first timestamp in nanoseconds. - * @param[in] t2 The second timestamp in nanoseconds. - * @return The absolute difference between the two timestamps in nanoseconds. - * @pre The timestamp is given in nanoseconds and the maximum value is less than 2^63 ns. - */ uint64_t sbndaq::SBNDXARAPUCADecoder::abs_difference(uint64_t t1, uint64_t t2) { uint64_t abs_diff = 0; @@ -973,12 +1000,6 @@ uint64_t sbndaq::SBNDXARAPUCADecoder::abs_difference(uint64_t t1, uint64_t t2) { return abs_diff; } -/** -* @brief Formats a timestamp in nanoseconds into a easily readable string format. -* @param[in] timestamp The timestamp in nanoseconds to be formatted. -* @return A string representation of the timestamp in the format "(seconds)nanoseconds ns". -* @details The function divides the input timestamp by 1,000,000,000 to obtain the seconds component and uses the modulus operator to get the remaining nanoseconds. -*/ std::string sbndaq::SBNDXARAPUCADecoder::print_timestamp(uint64_t timestamp) { return "(" + std::to_string(timestamp / NANOSEC_IN_SEC) + ")" + std::to_string(timestamp % NANOSEC_IN_SEC) + " ns"; } diff --git a/sbndcode/Decoders/XARAPUCA/run_xarapucadecoder.fcl b/sbndcode/Decoders/XARAPUCA/run_xarapucadecoder.fcl index 710a35c71..939d8e28e 100644 --- a/sbndcode/Decoders/XARAPUCA/run_xarapucadecoder.fcl +++ b/sbndcode/Decoders/XARAPUCA/run_xarapucadecoder.fcl @@ -1,6 +1,7 @@ #include "xarapucadecoder.fcl" #include "spectdc_decoder_sbnd.fcl" #include "SBNDPTBDecoderDefaults.fcl" +#include "frameshift_sbnd_data.fcl" process_name: XARAPUCADecoder @@ -17,7 +18,7 @@ outputs: { out1: { module_type: RootOutput fileName: "xarapucadecoder-art.root" # Output file containing the module products. - dataTier: "decoded-raw" + dataTier: "reconstructed" } } @@ -26,6 +27,7 @@ physics: { producers: { tdcdecoder: @local::spec_tdc_decoder_sbnd ptbdecoder: @local::SBNDPTBDecoderDefaults + frameshift: @local::frameshift_data xarapucadecoder: @local::xarapucadecoder } @@ -34,7 +36,7 @@ physics: { analyzers:{} - reco: [tdcdecoder, ptbdecoder, xarapucadecoder] + reco: [tdcdecoder, ptbdecoder, frameshift, xarapucadecoder] # Defines the output stream, there could be more than one if using filters. @@ -47,4 +49,4 @@ physics: { # end_paths is a keyword and contains the paths that do not modify the art::Event, # ie analyzers and output streams. These all run simultaneously. end_paths: [stream1] -} \ No newline at end of file +} diff --git a/sbndcode/Decoders/XARAPUCA/xarapucadecoder.fcl b/sbndcode/Decoders/XARAPUCA/xarapucadecoder.fcl index 3456450f9..2825c9316 100644 --- a/sbndcode/Decoders/XARAPUCA/xarapucadecoder.fcl +++ b/sbndcode/Decoders/XARAPUCA/xarapucadecoder.fcl @@ -7,39 +7,32 @@ xarapucadecoder: # INPUT configuration: # - CAEN Fragments. - caen_module_label: "daq" # Label identifying the CAEN fragments module. - caen_fragment_names: [ "CAENV1740", "ContainerCAENV1740" ] # Valid CAEN fragments for this decoder. + caen_fragment_names: [ "daq:CAENV1740", "daq:ContainerCAENV1740" ] # Valid CAEN fragments for this decoder. # - Boards configuration. num_caen_boards: 4 # Maximum number of boards to be considered. fragment_id_offset: 41216 # Offset to be applied to the fragment IDs (0xA100). board_id_list: [7, 13, 16, 19] # Slots where the boards are installed. # - Timing configuration. - timing_priority: 0 # 0 SPEC-TDC ETRIG, 1 PTB HLT (1-5, 14,15), 2 CAEN-only. ns_per_sample: 16 # Nanoseconds per sample. - # - RAW timestamp configuration: - raw_event_header_product_name: "daq" # Label identifying the module that created the artdaq::detail::RawEventHeader data product. - raw_trig_max_diff: 3000000 # Maximum difference (in ns) allowed between trigger timestamps and raw timestamp. - raw_timestamp_correction: 367000 # Correction (in ns) to be applied to the raw timestamp. - # - SPEC-TDC access configuration - spectdc_product_name: "tdcdecoder" # Name for getting SPEC-TDC Decoder products (if any). - spectdc_ftrig_ch: 3 # Channel assigned to flash triggers. - spectdc_etrig_ch: 4 # Channel assigned to ETT (Event Trigger Timestamp) triggers. - # - PTB access configuration - ptb_product_name: "ptbdecoder" # Name for getting PTB Decoder products (if any). - allowed_hl_triggers: [1, 2, 3, 4, 5, 14, 15] # List of allowed HLT trigger types (1-5, 14, 15). + # - FrameShift access configuration + frameshift_product_name: "frameshift" # Name for getting FrameShift module products (if any). # OUTPUT configuration: # - Output data product instance name. waveforms_instance_name: "XARAPUCAChannels" # Name for the instance product containing the raw decoded waveforms. timing_ref_instance_name: "" # Name for the instance product containing the timing reference information. store_debug_waveforms: 0 # Number of waveforms to store (0: none, -1: all, n: first n waveforms). + # Combination of extended fragments. + combine_ext_frag: true # (De)activates the combination of extended fragments into a single raw::OpDetWaveform object. + allowed_jittering: 64 # Allowed jittering (in ns) between fragments to be combined. # - Debug options. - debug_tdc_handle: false # (De)activates SPEC-TDC art::Handle information printing. - debug_ptb_handle: false # (De)activates PTB art::Handle information printing. + debug_frameshift_handle: false # (De)activates FrameShift art::Handle information printing. debug_fragments_handle: false # (De)activates V1740B CAEN fragments art::Handle information printing. debug_timing: false # (De)activates timing data printing. debug_buffer: false # (De)activates buffer status printing. debug_waveforms: false # (De)activates waveforms decoding printing. + debug_extended_fragments: false # (De)activates extended fragments information printing. + debug_jittering: false # (De)activates trigger jittering information printing. # - Verbose option. verbose: false # (De)activates verbosity. } diff --git a/sbndcode/Decoders/run_decoders_job.fcl b/sbndcode/Decoders/run_decoders_job.fcl index adb9f661d..a4fa855f5 100644 --- a/sbndcode/Decoders/run_decoders_job.fcl +++ b/sbndcode/Decoders/run_decoders_job.fcl @@ -37,6 +37,7 @@ #include "crt_decoder_sbnd.fcl" #include "crt_channel_map_service_sbnd.fcl" +#include "frameshift_sbnd_data.fcl" process_name: DECODE @@ -61,6 +62,7 @@ physics: xarapucadecoder: @local::xarapucadecoder ptbdecoder: @local::SBNDPTBDecoderDefaults crtdecoder: @local::crt_decoder + frameshift: @local::frameshift_data } filters: @@ -72,8 +74,8 @@ physics: analyzers: {} - dropchoppy: [ daq, ntbdecoder, crtdecoder, ptbdecoder, tdcdecoder, pmtdecoder, xarapucadecoder, filterchoppy ] - savechoppy: [ daq, ntbdecoder, crtdecoder, ptbdecoder, tdcdecoder, pmtdecoder, xarapucadecoder, "!filterchoppy" ] + dropchoppy: [ daq, ntbdecoder, crtdecoder, ptbdecoder, tdcdecoder, frameshift, pmtdecoder, xarapucadecoder, filterchoppy ] + savechoppy: [ daq, ntbdecoder, crtdecoder, ptbdecoder, tdcdecoder, frameshift, pmtdecoder, xarapucadecoder, "!filterchoppy" ] trigger_paths: [ dropchoppy, savechoppy ] stream1: [ out1 ] diff --git a/sbndcode/JobConfigurations/standard/reco/reco2_data.fcl b/sbndcode/JobConfigurations/standard/reco/reco2_data.fcl index 13452e2ff..32c9d7be8 100755 --- a/sbndcode/JobConfigurations/standard/reco/reco2_data.fcl +++ b/sbndcode/JobConfigurations/standard/reco/reco2_data.fcl @@ -1,7 +1,6 @@ #include "opt0finder_sbnd_data.fcl" #include "sbnd_tpcpmt3dbarycentermatching_config.fcl" #include "sbnd_lightpropagationcorrection_config.fcl" -#include "frameshift_sbnd_data.fcl" #include "standard_reco2_sbnd.fcl" services: @@ -33,15 +32,13 @@ physics.producers: ### SCE-aware shower reco for data pandoraSCEShower: @local::sbnd_sce_incremental_pandoraModularShowerCreationData pandoraSCEShowerSBN: @local::sbnd_sce_sbn_pandoraModularShowerCreationData - - frameshift: @local::frameshift_data } physics.reco2: [ pandora, pandoraTrack, blipPandoraCopy, blipPandoraTrackCopy, pandoraShower, pandoraShowerSBN, pandoraCaloData, pandoraPidData, pandoraLikePidData, cvn, opt0finder, crtveto, crtspacepointmatching, crttrackmatching, tpcpmtbarycentermatching, pandoraSCE, pandoraSCETrack, pandoraSCEShower, pandoraSCEShowerSBN, pandoraSCECaloData, pandoraSCEPidData, pandoraSCELikePidData, cvnSCE, opt0finderSCE, tpcpmtbarycentermatchingSCE, crtspacepointmatchingSCE, crttrackmatchingSCE, - caloskimCalorimetry, blipreco, lightpropagationcorrectionSCE, frameshift] + caloskimCalorimetry, blipreco, lightpropagationcorrectionSCE] physics.analyzers.caloskim.G4producer: "" diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index 4905589f3..cad824889 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -13,7 +13,6 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet fOpFlashLabel_tpc1 ( p.get("OpFlashLabel_tpc1") ), fSpacePointLabel( p.get("SpacePointLabel") ), fOpHitsModuleLabel( p.get("OpHitsModuleLabel") ), - fSPECTDCLabel( p.get("SPECTDCLabel") ), fFlashMatchingTool( p.get("FlashMatchingTool") ), fSaveCorrectionTree( p.get("SaveCorrectionTree") ), fSpeedOfLight( p.get("SpeedOfLight") ), @@ -25,12 +24,11 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet // // More initializers here. { + fNOpChannels = fWireReadout.NOpChannels(); // Initialize the TimeCorrectionVector PerChannel - for(size_t i = 0; i < fWireReadout.NOpChannels(); ++i) { - fTimeCorrectionPerChannel.push_back(0.0); // Initialize with zero or any default value - } + fTimeCorrectionPerChannel.resize(fNOpChannels, 0.0); // Initialize with zero or any default value - for(unsigned int opch=0; opchGetPropagationTimeCorrectionPerChannel(); - // Get the SPECTDC product required to go to the RWM reference frame - - art::Handle> tdcHandle; - e.getByLabel(fSPECTDCLabel, tdcHandle); - if (!tdcHandle.isValid() || tdcHandle->size() == 0){ - std::cout << "No SPECTDC products found. Skip this event." << std::endl; - ResetSliceInfo(); - continue; - } - else{ - const std::vector tdc_v(*tdcHandle); - for (size_t i=0; i> flashFM; if(fFlashMatchingTool == "OpT0Finder" ){ @@ -292,10 +271,10 @@ void sbnd::LightPropagationCorrection::produce(art::Event & e) 100., -1., Ycenter, Ywidth, Zcenter, Zwidth); newFlashTime = flasht0; sbn::CorrectedOpFlashTiming correctedOpFlashTiming; - correctedOpFlashTiming.OpFlashT0 = originalFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; + correctedOpFlashTiming.OpFlashT0 = originalFlashTime; correctedOpFlashTiming.NuToFLight = (Zcenter/fSpeedOfLight)/1000; correctedOpFlashTiming.NuToFCharge = (fRecoVz/fSpeedOfLight)/1000; - correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; + correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime; correctedOpFlashTimes->emplace_back(std::move(correctedOpFlashTiming)); } @@ -322,8 +301,6 @@ void sbnd::LightPropagationCorrection::beginJob() fTree->Branch("eventID", &fEvent, "eventID/i"); fTree->Branch("runID", &fRun, "runID/i"); fTree->Branch("subrunID", &fSubrun, "subrunID/i"); - fTree->Branch("RWMTime", &fRWMTime); - fTree->Branch("EventTriggerTime", &fEventTriggerTime); fTree->Branch("NuScore", &fNuScore); fTree->Branch("FMScore", &fFMScore); fTree->Branch("OpFlashTimeOld", &fOpFlashTimeOld); @@ -357,8 +334,6 @@ void sbnd::LightPropagationCorrection::ResetEventVars() fRun = 0; fSubrun = 0; _fNuScore = 0.0; - fRWMTime=-99999.; - fEventTriggerTime=-99999.; fNuScore.clear(); fFMScore.clear(); fOpFlashTimeOld.clear(); @@ -418,7 +393,7 @@ void sbnd::LightPropagationCorrection::ResetSliceInfo() fSpacePointY.clear(); fSpacePointZ.clear(); fSpacePointIntegral.clear(); - fTimeCorrectionPerChannel.resize(312, 0.0); // Reset the time correction vector for each channel + fTimeCorrectionPerChannel.assign(fNOpChannels, 0.0); // Reset the time correction vector for each channel fChargeBarycenterX.assign(2, 0.0); fChargeBarycenterY.assign(2, 0.0); fChargeBarycenterZ.assign(2, 0.0); @@ -442,16 +417,16 @@ void sbnd::LightPropagationCorrection::GetPropagationTimeCorrectionPerChannel() double dy = fSpacePointY[sp] - _opDetY; double dz = fSpacePointZ[sp] - _opDetZ; double distanceToOpDet = std::sqrt(dx*dx + dy*dy + dz*dz); - //double spToCathode = abs(fSpacePointX[sp]); // Distance from space point to cathode in mm - //double cathodeToOpDet = std::sqrt(_opDetX*_opDetX + dy*dy + dz*dz); // Distance from cathode to OpDet in mm - //float lightPropTimeVIS = spToCathode/fVGroupVUV + cathodeToOpDet/fVGroupVIS; // Speed - double cathodeToOpDet = std::sqrt(_opDetX*_opDetX + (dy/2)*(dy/2) + (dz/2)*(dz/2)); // Distance from cathode to OpDet in mm double spToCathode = std::sqrt( fSpacePointX[sp]*fSpacePointX[sp] + (dy/2)*(dy/2) + (dz/2)*(dz/2)); // Distance from space point to cathode in mm float lightPropTimeVIS = spToCathode/fVGroupVUV + cathodeToOpDet/fVGroupVIS; // Speed float lightPropTimeVUV = distanceToOpDet / fVGroupVUV; // Speed of light in mm/ns for VUV - float lightPropTime = std::min(lightPropTimeVIS, lightPropTimeVUV); + float lightPropTime = 0; + if(fPDSMap.pdType(opdet)=="pmt_coated" || fPDSMap.pdType(opdet)=="xarapuca_vuv") + lightPropTime = std::min(lightPropTimeVIS, lightPropTimeVUV); + else if(fPDSMap.pdType(opdet)=="pmt_uncoated" || fPDSMap.pdType(opdet)=="xarapuca_vis") + lightPropTime = lightPropTimeVIS; float partPropTime = std::sqrt((fSpacePointX[sp]-fRecoVx)*(fSpacePointX[sp]-fRecoVx) + (fSpacePointY[sp]-fRecoVy)*(fSpacePointY[sp]-fRecoVy) + (fSpacePointZ[sp]-fRecoVz)*(fSpacePointZ[sp]-fRecoVz))/fSpeedOfLight; float PropTime = lightPropTime + partPropTime; if(PropTime < minPropTime) minPropTime = PropTime; @@ -543,4 +518,4 @@ ::lightana::LiteOpHitArray_t sbnd::LightPropagationCorrection::GetAssociatedLite } return flash_hits_v; - } \ No newline at end of file + } diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index e5427282a..d0ae118fd 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -145,7 +145,7 @@ private: geo::WireReadoutGeom const& fWireReadout = art::ServiceHandle()->Get(); - + opdet::sbndPDMapAlg fPDSMap; //Flash finder manager ::lightana::FlashFinderManager _mgr; ::lightana::FlashFinderManager _mgr_tpc0; @@ -171,7 +171,6 @@ private: std::string fSpacePointLabel; std::string fOpHitsModuleLabel; std::string fOpFlashNewLabel; - std::string fSPECTDCLabel; std::string fFlashMatchingTool; bool fSaveCorrectionTree; @@ -215,10 +214,9 @@ private: int fEvent; int fRun; int fSubrun; + size_t fNOpChannels; double _fNuScore; double _fFMScore; - double fEventTriggerTime=-999999.; - double fRWMTime=-999999.; std::vector fNuScore; std::vector fFMScore; std::vector fOpFlashTimeOld; diff --git a/sbndcode/Timing/CMakeLists.txt b/sbndcode/Timing/CMakeLists.txt index 1ae0bbd44..b7023b7ca 100644 --- a/sbndcode/Timing/CMakeLists.txt +++ b/sbndcode/Timing/CMakeLists.txt @@ -1,9 +1,13 @@ -art_dictionary(DICTIONARY_LIBRARIES -lardataobj::RawData) +art_make_library( + LIBRARIES art::Persistency_Common + messagefacility::MF_MessageLogger + ) + +art_dictionary(DICTIONARY_LIBRARIES lardataobj::RawData) install_headers() -# install_fhicl() -# install_source() +#install_fhicl() +install_source() add_subdirectory(FrameShift) add_subdirectory(WaveformAlignment) diff --git a/sbndcode/Timing/FrameShift/FrameShift_module.cc b/sbndcode/Timing/FrameShift/FrameShift_module.cc index c675eb5a8..6643c1001 100644 --- a/sbndcode/Timing/FrameShift/FrameShift_module.cc +++ b/sbndcode/Timing/FrameShift/FrameShift_module.cc @@ -36,13 +36,21 @@ #include -namespace sbnd { +namespace sbnd::timing { class FrameShift; } -class sbnd::FrameShift : public art::EDProducer { +class sbnd::timing::FrameShift : public art::EDProducer { public: + /** + * @brief Constructor for the FrameShift art::EDProducer module. + * @param[in] p The FHiCL parameter set containing all configuration values for this module. + * @details Reads all FHiCL configuration parameters, including DAQ header module/instance + * labels, PTB and TDC decode labels, all HLT codes for each trigger stream (Beam, Offbeam, + * Xmuon), calibration shifts between TDC and PTB timing references, and debug/output options. + * Declares FrameShiftInfo and TimingInfo as output data products. + */ explicit FrameShift(fhicl::ParameterSet const& p); // The compiler-generated destructor is fine for non-base // classes without bare pointers or other resource use. @@ -54,713 +62,785 @@ class sbnd::FrameShift : public art::EDProducer { FrameShift& operator=(FrameShift&&) = delete; // Required functions. + /** + * @brief Main per-event processing function. + * @param[in,out] e The art::Event being processed. + * @details For each event this function: + * 1. Resets all per-event variables. + * 2. Retrieves the raw DAQ header timestamp. + * 3. Retrieves and unpacks PTB HLT trigger words and timestamps. + * 4. Retrieves SPECTDC timestamps for all relevant channels. + * 5. Identifies the ETRIG timestamps from TDC channel 4 and the PTB HLTs. + * 6. Selects a global ETRIG reference (TDC > PTB > raw header). + * 7. Picks the TDC and PTB timestamps closest to the global ETRIG reference. + * 8. Computes frame shift values for the CRT T1, beam gate, and event trigger frames. + * 9. Selects the default frame for the trigger stream (gate for Beam/Offbeam, ETRIG for Xmuon). + * 10. Produces a FrameShiftInfo and TimingInfo object and puts them in the event. + */ void produce(art::Event& e) override; // Selected optional functions. + /** + * @brief Called once before event processing begins; creates the output TTree and registers + * all branches if fMakeTree is true. + * @details Branches cover run/subrun/event identifiers, raw per-channel TDC timestamp + * vectors, the selected per-event TDC and PTB timestamps, and the three computed frame + * shift values (CRT T1, gate, ETRIG) together with their timing type and channel metadata. + */ void beginJob() override; + /** + * @brief Resets all per-event member variables to their default or invalid sentinel values. + * @details Sets run/subrun/event to -1, clears all TDC and PTB timestamp vectors, resets + * all selected timestamps to kInvalidTimestamp, resets the stream flags (_isBeam, + * _isOffbeam, _isXmuon) to false, and initialises all frame shift values together with + * their timing type and channel metadata to their respective kInvalid* sentinels. + * Called at the beginning of every produce() call. + */ void ResetEventVars(); - double SubtractUTCTimestmap(const uint64_t& ts1, const uint64_t& ts2); -private: - - // Declare member data here. - //DAQ Header - std::string fDAQHeaderModuleLabel; - std::string fDAQHeaderInstanceLabel; - uint64_t fRawTSCorrection = 0; //Correction for the Event timestamp to account for NTB - uint64_t fMaxAllowedRefTimeDiff; - uint64_t _raw_ts = 0; //ns - - //Timing Reference - art::InputTag fTimingRefPmtLabel; - art::InputTag fTimingRefCrtLabel; + /** + * @brief Returns the timestamp from a vector that is closest to a reference timestamp. + * @param[in] timestamps Vector of timestamps to search through (in nanoseconds). + * @param[in] reference The reference timestamp to compare against (in nanoseconds). + * @return The timestamp with the smallest absolute difference to @p reference, + * or kInvalidTimestamp if the vector is empty. + */ + uint64_t FindClosest(const std::vector ×tamps, const uint64_t &reference); + /** + * @brief Formats a 64-bit nanosecond timestamp into a human-readable string. + * @param[in] timestamp The timestamp in nanoseconds. + * @return A string of the form "(s) = X, (ns) = Y" where X is the whole-second component + * and Y is the sub-second nanosecond remainder. + */ + std::string PrintFormatTimestamp(const uint64_t ×tamp); + + /** + * @brief Reads the artdaq RawEventHeader from the event and stores the corrected timestamp. + * @param[in] e The art::Event to read from. + * @details Retrieves the artdaq::detail::RawEventHeader using fDAQHeaderModuleLabel and + * fDAQHeaderInstanceLabel, then subtracts fRawTSCorrection to account for the NTB offset, + * storing the result in _raw_ts. Throws a cet::exception if the handle is not valid. + */ + void GetRawTimestamp(const art::Event &e); + /** + * @brief Reads decoded SPECTDC DAQTimestamp objects from the event and sorts them by channel. + * @param[in] e The art::Event to read from. + * @details Retrieves std::vector using fTdcDecodeLabel and appends each + * timestamp to the appropriate per-channel vector: + * - Channel 0: CRT T1 reset (_tdc_ch0) + * - Channel 1: BES (_tdc_ch1) + * - Channel 2: RWM (_tdc_ch2) + * - Channel 4: Event trigger (_tdc_ch4) + * Logs a non-fatal info message if no timestamps are found. + */ + void GetTDCTimestamps(const art::Event &e); + /** + * @brief Reads decoded PTB HLT words from the event and unpacks per-bit trigger timestamps. + * @param[in] e The art::Event to read from. + * @details Retrieves std::vector using fPtbDecodeLabel. For each HLT + * trigger word the bitmask is unpacked so that every set bit becomes a separate entry in + * _ptb_hlt_trunmask and _ptb_hlt_unmask_timestamp, enabling per-HLT-code lookup later. + * HLT timestamps are converted from PTB clock ticks to nanoseconds (multiplied by 20). + * Throws a cet::exception if the handle is not valid or empty. + */ + void GetPTBTimestamps(const art::Event &e); + /** + * @brief Identifies the ETRIG timestamps from TDC channel 4 and PTB HLT words closest to + * the raw DAQ header timestamp. + * @details Uses the already-populated _tdc_ch4 and _ptb_hlt_unmask_timestamp collections. + * - Selects the TDC ETRIG timestamp closest to _raw_ts via FindClosest. + * - Collects all PTB HLT entries whose code appears in fPtbEtrigHlts, then selects the one + * closest to _raw_ts, storing the result in _hlt_etrig and _hlt_etrig_ts. + * Throws a cet::exception if either TDC or PTB ETRIG timestamps are missing. + */ + void FindETRIGs(); + + /** + * @brief Selects the global ETRIG reference timestamp with priority: + * TDC ETRIG > PTB HLT ETRIG > raw DAQ header timestamp. + * @return The chosen global ETRIG reference timestamp in nanoseconds. + * @details This timestamp is used as the anchor for selecting the closest individual TDC + * channel and PTB HLT timestamps in DecideRelevantTDCTimestamps and + * DecideRelevantPTBTimestamps. + */ + uint64_t DecideGlobalEtrigTimestamp(); + /** + * @brief For each active TDC channel, selects the timestamp closest to the global ETRIG reference. + * @param[in] global_etrig_ts The global ETRIG reference timestamp in nanoseconds. + * @details Searches _tdc_ch0 (CRT T1), _tdc_ch1 (BES), and _tdc_ch2 (RWM) for the closest + * timestamp to @p global_etrig_ts via FindClosest, storing results in _tdc_crtt1_ts, + * _tdc_bes_ts, and _tdc_rwm_ts respectively. Channels with no recorded timestamps are left + * at kInvalidTimestamp. + */ + void DecideRelevantTDCTimestamps(uint64_t global_etrig_ts); + /** + * @brief Identifies the trigger stream and retrieves the corresponding gate and CRT T1 + * HLT timestamps from the PTB record. + * @param[in] global_etrig_ts The global ETRIG reference timestamp in nanoseconds + * (used for API symmetry with DecideRelevantTDCTimestamps). + * @details Matches _hlt_etrig against fBeamEtrigHlt, fOffbeamEtrigHlt, and fXmuonEtrigHlt + * in that order to set _isBeam, _isOffbeam, or _isXmuon and assign the appropriate + * _hlt_gate and _hlt_crtt1 codes. Then scans _ptb_hlt_trunmask to find matching gate and + * CRT T1 entries and fills _hlt_gate_ts and _hlt_crtt1_ts. + * Throws a cet::exception if the ETRIG HLT code does not match any known stream. + */ + void DecideRelevantPTBTimestamps(uint64_t global_etrig_ts); - uint16_t _pmt_timing_type; // e.g. SPECTDC = 0; PTB HLT = 1; Do nothing = 2 - uint16_t _pmt_timing_ch; - uint16_t _crt_timing_type; - uint16_t _crt_timing_ch; +private: - //PTB - art::InputTag fPtbDecodeLabel; - - std::vector fPtbEtrigHlts; - - std::vector fBeamEtrigHlt; - std::vector fOffbeamEtrigHlt; - std::vector fXmuonEtrigHlt; - - int fBeamCrtT1Hlt; - int fOffbeamCrtT1Hlt; - - int fBeamGateHlt; - int fOffbeamGateHlt; - - std::vector _ptb_hlt_trigger; - std::vector _ptb_hlt_timestamp; - std::vector _ptb_hlt_unmask_timestamp; - std::vector _ptb_hlt_trunmask; - - std::vector _ptb_hlt_etrig; - std::vector _ptb_hlt_etrig_ts; - - bool _isBeam; - bool _isOffbeam; - bool _isXmuon; - int _hlt_etrig; - uint64_t _hlt_etrig_ts; - int _hlt_gate; - uint64_t _hlt_gate_ts; - int _hlt_crtt1; - uint64_t _hlt_crtt1_ts; + // FCL Parameters + art::InputTag fDAQHeaderLabel; ///< Instance label for DAQ header product. + uint64_t fRawTSCorrection; ///< Correction to the DAQ header timestamp to account for NTB offset. + art::InputTag fPtbDecodeLabel; ///< Module label for decoded PTB data. + std::vector fPtbEtrigHlts; ///< All allowed PTB HLT codes for event triggers. + std::vector fBeamEtrigHlt; ///< PTB HLT codes associated with beam trigger. + std::vector fOffbeamEtrigHlt; ///< PTB HLT codes associated with offbeam trigger. + std::vector fXmuonEtrigHlt; ///< PTB HLT codes associated with crossing muon trigger. + int fBeamCrtT1Hlt; ///< PTB HLT code for beam CRT T1 reset. + int fOffbeamCrtT1Hlt; ///< PTB HLT code for offbeam CRT T1 reset. + int fBeamGateHlt; ///< PTB HLT code for beam gate opening. + int fOffbeamGateHlt; ///< PTB HLT code for offbeam gate opening. + art::InputTag fTdcDecodeLabel; ///< Module label for decoded TDC data. + uint64_t fShiftData2MC; ///< Value to shift Data to MC so that data agree with MC [ns]. + ///< TODO: Derive this value and verify if it is consistent across PMT/CRT. + ///< TODO: Get this value from database instead of FHiCL parameter. + uint64_t fShiftTdcRwm2PtbGate; ///< Value to move TDC RWM (ch2) frame to agree with HLT gate frame. + ///< TODO: Get this value from database instead of FHiCL parameter. + uint64_t fShiftTdcEtrig2PtbEtrig; ///< Value to move TDC ETRIG (ch4) to agree with HLT ETRIG. + ///< TODO: Get from database. + uint64_t fShiftTdcCrtt12PtbCrtt1; ///< Value to move TDC CRT T1 reset (ch0) to agree with HLT CRT T1. + ///< TODO: Get from database. + bool fMakeTree; ///< Whether to produce a TTree in the hist file. + bool fDebugDAQHeader; ///< Whether to print debug statements relevant to DAQ header. + bool fDebugPtb; ///< Whether to print debug statements relevant to PTB. + bool fDebugTdc; ///< Whether to print debug statements relevant to TDC. + bool fDebugFrame; ///< Whether to print debug statements relevant to frame shifts. + + // Global Variables, set in processing + int _run, _subrun, _event; ///< Stores the unique run/subrun/event number combination for this event. + + uint64_t _raw_ts; ///< Stores DAQ header timestamp. + std::vector _ptb_hlt_trigger; ///< Stores full trigger word from PTB HLT (indexed per HLT). + std::vector _ptb_hlt_timestamp; ///< Stores timestamp associated with PTB HLT (indexed per HLT). + std::vector _ptb_hlt_unmask_timestamp; ///< Stores timestamp associated with PTB HLT for each element (indexed per HLT element once unmasked). + std::vector _ptb_hlt_trunmask; ///< Stores PTB HLT code for each element (indexed per HLT element once unmasked). + std::vector _ptb_hlt_etrig; ///< Stores all PTB HLTs that are an allowed type for this stream. + std::vector _ptb_hlt_etrig_ts; ///< Stores the associated timestamps for the above HLTs. + + bool _isBeam; ///< Event is a beam trigger (BNB or BNBLight). + bool _isOffbeam; ///< Event is an offbeam trigger (OffBeam or OffBeamLight). + bool _isXmuon; ///< Event is a crossing muon trigger. + + int _hlt_etrig; ///< Stores the HLT event trigger code once selected as closest to DAQ header timestamp. + uint64_t _hlt_etrig_ts; ///< Stores the HLT event trigger timestamp once selected as closest to DAQ header timestamp. + int _hlt_gate; ///< Stores the HLT gate opening code once selected as closest to DAQ header timestamp. + uint64_t _hlt_gate_ts; ///< Stores the HLT gate opening timestamp once selected as closest to DAQ header timestamp. + int _hlt_crtt1; ///< Stores the HLT CRT T1 code once selected as closest to DAQ header timestamp. + uint64_t _hlt_crtt1_ts; ///< Stores the HLT CRT T1 timestamp once selected as closest to DAQ header timestamp. - //TDC - art::InputTag fTdcDecodeLabel; - - std::vector _tdc_ch0; - std::vector _tdc_ch1; - std::vector _tdc_ch2; - std::vector _tdc_ch4; - - uint64_t _tdc_crtt1_ts; - uint64_t _tdc_bes_ts; - uint64_t _tdc_rwm_ts; - uint64_t _tdc_etrig_ts; - - - //Frame Shift - uint16_t _global_timing_type; // e.g. SPECTDC = 0; PTB HLT = 1; - uint64_t _global_frame; - - double _frame_tdc_crtt1; //ns - double _frame_tdc_bes; //ns - double _frame_tdc_rwm; //ns - double _frame_hlt_crtt1; //ns - double _frame_hlt_gate; //ns - - //Value to apply at CAF depending on which stream - double _frame_apply_at_caf; - - //Value to shift Data to MC -- so that data agree with MC [ns] - //TODO: Derive this value and verify if it is consistent across pmt/crt - //TODO: Get this value from database instead of fhicl parameter - double _frame_data2mc; //ns - - //Value to move RWM frame to agree with HLT Gate Frame - //This is derived by subtracting: TDC RWM - HLT Gate. - //Using the MC2025B dataset, this distribution has a mean of 1738 ns and std of 9 ns. - //TODO: Get this value from database instead of fhicl parameter - double fShiftRWM2Gate; //ns - - //Debug - bool fDebugDAQHeader; - bool fDebugTimingRef; - bool fDebugPtb; - bool fDebugTdc; - bool fDebugFrame; - - //---TREE PARAMETERS - TTree *fTree; - art::ServiceHandle tfs; - int _run, _subrun, _event; + std::vector _tdc_ch0; ///< Stores all timestamps recorded in TDC channel 0 (CRT T1). + std::vector _tdc_ch1; ///< Stores all timestamps recorded in TDC channel 1 (BES). + std::vector _tdc_ch2; ///< Stores all timestamps recorded in TDC channel 2 (RWM). + //std::vector _tdc_ch3; ///< Stores all timestamps recorded in TDC channel 3 (flash trigger). + std::vector _tdc_ch4; ///< Stores all timestamps recorded in TDC channel 4 (event trigger). + + uint64_t _tdc_crtt1_ts; ///< Stores the TDC CRT T1 timestamp selected as closest to DAQ header timestamp. + uint64_t _tdc_bes_ts; ///< Stores the TDC BES timestamp selected as closest to DAQ header timestamp. + uint64_t _tdc_rwm_ts; ///< Stores the TDC RWM timestamp selected as closest to DAQ header timestamp. + uint64_t _tdc_etrig_ts; ///< Stores the TDC event trigger timestamp selected as closest to DAQ header timestamp. + + uint64_t _frame_crtt1; ///< Stores the shift required to move from the PPS frame to the CRT T1 frame. + uint16_t _timing_type_crtt1; ///< Stores the type (TDC/PTB) used to produce the CRT T1 shift. + uint16_t _timing_channel_crtt1; ///< Stores the channel (TDC) or code (PTB) used to produce the CRT T1 shift. + uint64_t _frame_gate; ///< Stores the shift required to move from the PPS frame to the gate opening frame. + uint16_t _timing_type_gate; ///< Stores the type (TDC/PTB) used to produce the gate shift. + uint16_t _timing_channel_gate; ///< Stores the channel (TDC) or code (PTB) used to produce the gate shift. + uint64_t _frame_etrig; ///< Stores the shift required to move from the PPS frame to the event trigger frame. + uint16_t _timing_type_etrig; ///< Stores the type (TDC/PTB) used to produce the ETRIG shift. + uint16_t _timing_channel_etrig; ///< Stores the channel (TDC) or code (PTB) used to produce the ETRIG shift. + + uint64_t _frame_default; ///< Stores the shift required to move from the PPS frame to the default frame for this stream. + uint16_t _timing_type_default; ///< Stores the type (TDC/PTB) used to produce the default shift. + uint16_t _timing_channel_default; ///< Stores the channel (TDC) or code (PTB) used to produce the default shift. + + // Tree production + TTree *fTree; ///< Pointer to the output tree produced when fMakeTree is enabled. + art::ServiceHandle tfs; ///< Service handle used to create ROOT output objects. }; -sbnd::FrameShift::FrameShift(fhicl::ParameterSet const& p) - : EDProducer{p} // , - // More initializers here. +sbnd::timing::FrameShift::FrameShift(fhicl::ParameterSet const& p) + : EDProducer{p} { - fDAQHeaderInstanceLabel = p.get("DAQHeaderInstanceLabel"); - fDAQHeaderModuleLabel = p.get("DAQHeaderModuleLabel"); + fDAQHeaderLabel = p.get("DAQHeaderLabel"); fRawTSCorrection = p.get("RawTSCorrection"); - fMaxAllowedRefTimeDiff = p.get("MaxAllowedRefTimeDiff"); - - fTimingRefPmtLabel = p.get("TimingRefPmtLabel"); - fTimingRefCrtLabel = p.get("TimingRefCrtLabel"); - - fTdcDecodeLabel = p.get("TdcDecodeLabel"); fPtbDecodeLabel = p.get("PtbDecodeLabel"); - fPtbEtrigHlts = p.get>("PtbEtrigHlts"); - fBeamEtrigHlt = p.get>("BeamEtrigHlt"); fOffbeamEtrigHlt = p.get>("OffbeamEtrigHlt"); fXmuonEtrigHlt = p.get>("XmuonEtrigHlt"); - fBeamCrtT1Hlt = p.get("BeamCrtT1Hlt"); fOffbeamCrtT1Hlt = p.get("OffbeamCrtT1Hlt"); - fBeamGateHlt = p.get("BeamGateHlt"); fOffbeamGateHlt = p.get("OffbeamGateHlt"); - + fTdcDecodeLabel = p.get("TdcDecodeLabel"); + fShiftData2MC = p.get("ShiftData2MC"); + fShiftTdcRwm2PtbGate = p.get("ShiftTdcRwm2PtbGate"); + fShiftTdcEtrig2PtbEtrig = p.get("ShiftTdcEtrig2PtbEtrig"); + fShiftTdcCrtt12PtbCrtt1 = p.get("ShiftTdcCrtt12PtbCrtt1"); + fMakeTree = p.get("MakeTree", false); fDebugDAQHeader = p.get("DebugDAQHeader", false); - fDebugTimingRef = p.get("DebugTimingRef", false); fDebugPtb = p.get("DebugPtb", false); fDebugTdc = p.get("DebugTdc", false); fDebugFrame = p.get("DebugFrame", false); - //TODO: Get from database instead of fhicl parameters - _frame_data2mc = p.get("ShiftData2MC"); - - //TODO: Get from database instead of fhicl parameters - fShiftRWM2Gate = p.get("ShiftRWM2Gate"); - - produces< sbnd::timing::FrameShiftInfo >(); - produces< sbnd::timing::TimingInfo >(); + produces(); + produces(); } -void sbnd::FrameShift::produce(art::Event& e) +void sbnd::timing::FrameShift::produce(art::Event& e) { - ResetEventVars(); _run = e.id().run(); _subrun = e.id().subRun(); _event = e.id().event(); - if (fDebugTdc | fDebugPtb | fDebugTimingRef | fDebugFrame) - std::cout <<"#----------RUN " << _run << " SUBRUN " << _subrun << " EVENT " << _event <<"----------#\n"; - - //---------------------------TimingRef-----------------------------// - art::Handle timingRefPmtHandle; - e.getByLabel(fTimingRefPmtLabel, timingRefPmtHandle); - - if (!timingRefPmtHandle.isValid()){ - throw cet::exception("FrameShift") << "No raw::TimingReferenceInfo found w/ tag " << fTimingRefPmtLabel << ". Check data quality!"; - } - else{ - raw::TimingReferenceInfo const& pmt_timing(*timingRefPmtHandle); - _pmt_timing_type = pmt_timing.timingType; - _pmt_timing_ch = pmt_timing.timingChannel; - } - timingRefPmtHandle.removeProduct(); - - art::Handle timingRefCrtHandle; - e.getByLabel(fTimingRefCrtLabel, timingRefCrtHandle); - - if (!timingRefCrtHandle.isValid()){ - throw cet::exception("FrameShift") << "No raw::TimingReferenceInfo found w/ tag " << fTimingRefCrtLabel << ". Check data quality!"; - } - else{ - raw::TimingReferenceInfo const& crt_timing(*timingRefCrtHandle); - _crt_timing_type = crt_timing.timingType; - _crt_timing_ch = crt_timing.timingChannel; - } - timingRefCrtHandle.removeProduct(); - - //TODO: Add XA check! - - //Check what frame it is decoded to? - if (_pmt_timing_type != _crt_timing_type){ - throw cet::exception("FrameShift") << "Timing Reference for PMT and CRT are not the same! PMT type = " - << _pmt_timing_type << ", CRT type = " << _crt_timing_type; - }else{ - _global_timing_type = _pmt_timing_type; - } - - if (fDebugTimingRef){ - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "Timing Reference For Decoding Pmt" << std::endl; - std::cout << " Type = " << _pmt_timing_type << " (SPECTDC = 0; PTB HLT = 1; Nothing = 2)." << std::endl; - std::cout << " Channel = " << _pmt_timing_ch << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "Timing Reference For Decoding Crt" << std::endl; - std::cout << " Type = " << _crt_timing_type << " (SPECTDC = 0; PTB HLT = 1; Nothing = 2)." << std::endl; - std::cout << " Channel = " << _crt_timing_ch << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "Global Timing Type = " << _global_timing_type << " (SPECTDC = 0; PTB HLT = 1; Nothing = 2)." << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - } - - //----------------------------DAQ Header-----------------------------// + if(fDebugTdc | fDebugPtb | fDebugFrame) + std::cout <<"#----------RUN " << _run << " SUBRUN " << _subrun << " EVENT " << _event <<"----------#\n"; + + GetRawTimestamp(e); + GetPTBTimestamps(e); + GetTDCTimestamps(e); + FindETRIGs(); + + uint64_t global_etrig_ts = DecideGlobalEtrigTimestamp(); + DecideRelevantTDCTimestamps(global_etrig_ts); + DecideRelevantPTBTimestamps(global_etrig_ts); + + if(_isBeam || _isOffbeam) + { + //Frame CRT T1 + if(_tdc_crtt1_ts != kInvalidTimestamp) + { + _frame_crtt1 = _tdc_crtt1_ts - fShiftTdcCrtt12PtbCrtt1; + _timing_type_crtt1 = kSPECTDCType; + _timing_channel_crtt1 = 0; + } + else if(_hlt_crtt1_ts != kInvalidTimestamp) + { + _frame_crtt1 = _hlt_crtt1_ts; + _timing_type_crtt1 = kPTBHLTType; + _timing_channel_crtt1 = _hlt_crtt1; + } + else + { + _frame_crtt1 = 0; + _timing_type_crtt1 = kNoShiftType; + _timing_channel_crtt1 = kInvalidChannel; + } + + //Frame Beam Gate + if(_isBeam && _tdc_rwm_ts != kInvalidTimestamp) + { + _frame_gate = _tdc_rwm_ts - fShiftTdcRwm2PtbGate; //TODO: Add shift from Data to MD: + fShiftData2MC + _timing_type_gate = kSPECTDCType; + _timing_channel_gate = 2; + } + else if(_hlt_gate_ts != kInvalidTimestamp) + { + _frame_gate = _hlt_gate_ts; + _timing_type_gate = kPTBHLTType; + _timing_channel_gate = _hlt_gate; + } + else + { + _frame_gate = 0; + _timing_type_gate = kNoShiftType; + _timing_channel_gate = kInvalidChannel; + } + } + + //Frame ETRIG + if(_tdc_etrig_ts != kInvalidTimestamp) + { + _frame_etrig = _tdc_etrig_ts - fShiftTdcEtrig2PtbEtrig; + _timing_type_etrig = kSPECTDCType; + _timing_channel_etrig = 4; + } + else if(_hlt_etrig_ts != kInvalidTimestamp) + { + _frame_etrig = _hlt_etrig_ts; + _timing_type_etrig = kPTBHLTType; + _timing_channel_etrig = _hlt_etrig; + } + else + { + _frame_etrig = 0; + _timing_type_etrig = kNoShiftType; + _timing_channel_etrig = kInvalidChannel; + } + + //-----------------------Pick default frame-----------------------// + // The follow picks which frame to apply at downstream stage and store it as frame_default, based on the stream + // + // 1. Beam Stream: recontruct the beam spill + porch relative to RWM record in TDC since it has better resolution. + // Then, apply a constant to shift the RWM fram + // i.e. frame_tdc_rwm + shift rwm to beam gate HLT 26 + // + // 2. Offbeam Stream: reconstruct the porch relative to the beam gate opening frame + // i.e. frame_hlt_gate equiv + // + // 3. Xmuon Stream: shift to etrig + + if(_isBeam || _isOffbeam) + { + //Pick default stream -- beam gate + _frame_default = _frame_gate; + _timing_type_default = _timing_type_gate; + _timing_channel_default = _timing_channel_gate; + } + else if(_isXmuon) + { + //Pick default stream -- ETRIG + _frame_default = _frame_etrig; + _timing_type_default = _timing_type_etrig; + _timing_channel_default = _timing_channel_etrig; + } + + if(fDebugFrame) + { + std::cout << "--------------------------------------" << std::endl; + std::cout << "Frame Shift Results:" << std::endl; + std::cout << "Frame CRT T1 type " << _timing_type_crtt1 <<", channel " << _timing_channel_crtt1 << ": " + << PrintFormatTimestamp(_frame_crtt1) << std::endl; + std::cout << "Frame Beam Gate type " << _timing_type_gate <<", channel " << _timing_channel_gate << ": " + << PrintFormatTimestamp(_frame_gate) << std::endl; + std::cout << "Frame ETRIG type " << _timing_type_etrig <<", channel " << _timing_channel_etrig << ": " + << PrintFormatTimestamp(_frame_etrig) << std::endl; + std::cout << "Default Frame type " << _timing_type_default <<", channel " << _timing_channel_default << ": " + << PrintFormatTimestamp(_frame_default) << std::endl; + std::cout << "--------------------------------------" << std::endl; + } + + //Put product in event + std::unique_ptr newFrameShiftInfo(new FrameShiftInfo(_frame_crtt1, _timing_type_crtt1, _timing_channel_crtt1, + _frame_gate, _timing_type_gate, _timing_channel_gate, + _frame_etrig, _timing_type_etrig,_timing_channel_etrig, + _frame_default, _timing_type_default, _timing_channel_default)); + + std::unique_ptr newTimingInfo(new TimingInfo(_raw_ts, _tdc_crtt1_ts, _tdc_bes_ts, _tdc_rwm_ts, _tdc_etrig_ts, + _hlt_crtt1_ts, _hlt_etrig_ts, _hlt_gate_ts)); + + e.put(std::move(newTimingInfo)); + e.put(std::move(newFrameShiftInfo)); + + //Fill the tree + if(fMakeTree) + fTree->Fill(); +} + +uint64_t sbnd::timing::FrameShift::FindClosest(const std::vector ×tamps, const uint64_t &reference) +{ + uint64_t min_diff = kInvalidTimestamp; + uint64_t closest = kInvalidTimestamp; + + for(const uint64_t ×tamp : timestamps) + { + uint64_t diff = reference > timestamp ? reference - timestamp : timestamp - reference; + + if(diff < min_diff) + { + min_diff = diff; + closest = timestamp; + } + } + + return closest; +} + +std::string sbnd::timing::FrameShift::PrintFormatTimestamp(const uint64_t ×tamp) +{ + std::stringstream ss; + ss << "(s) = " << timestamp / kSecondInNanoseconds << ", (ns) = " << timestamp % kSecondInNanoseconds; + + return ss.str(); +} + +void sbnd::timing::FrameShift::GetRawTimestamp(const art::Event &e) +{ art::Handle DAQHeaderHandle; - e.getByLabel(fDAQHeaderModuleLabel, fDAQHeaderInstanceLabel, DAQHeaderHandle); + e.getByLabel(fDAQHeaderLabel, DAQHeaderHandle); - if (!DAQHeaderHandle.isValid()){ - throw cet::exception("FrameShift") << "No artdaq::detail::RawEventHeader found w/ tag " << fDAQHeaderModuleLabel << ". Check data quality!"; - } - else{ + if(!DAQHeaderHandle.isValid()) + throw cet::exception("FrameShift") << "No artdaq::detail::RawEventHeader found w/ tag " << fDAQHeaderLabel << ". Check data quality!"; + else + { artdaq::RawEvent rawHeaderEvent = artdaq::RawEvent(*DAQHeaderHandle); _raw_ts = rawHeaderEvent.timestamp() - fRawTSCorrection; - } - if (fDebugDAQHeader){ - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "DAQ Header Timestamp: " - << " (s) = " << _raw_ts/uint64_t(1e9) - << ", (ns) = " << _raw_ts%uint64_t(1e9) - << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - } - - //---------------------------TDC-----------------------------// - art::Handle> tdcHandle; - e.getByLabel(fTdcDecodeLabel, tdcHandle); + } + + if(fDebugDAQHeader) + { + std::cout << "----------------------------------------------------" << std::endl; + std::cout << "DAQ Header Timestamp: " << PrintFormatTimestamp(_raw_ts) << std::endl; + } +} - if (!tdcHandle.isValid() || tdcHandle->size() == 0){ - mf::LogInfo("FrameShift") << "No sbnd::timing::DAQTimestamp found w/ tag " << fTdcDecodeLabel << ". Check data quality!\n"; - //throw cet::exception("FrameShift") << "No sbnd::timing::DAQTimestamp found w/ tag " << fTdcDecodeLabel << ". Check data quality!"; - } - else{ - std::vector tdc_v(*tdcHandle); - for (size_t i=0; i> tdcHandle; + e.getByLabel(fTdcDecodeLabel, tdcHandle); + + if(!tdcHandle.isValid() || tdcHandle->size() == 0) + mf::LogInfo("FrameShift") << "No DAQTimestamp found w/ tag " << fTdcDecodeLabel << ". Check data quality!\n"; + else + { + for(auto const& tdc : *tdcHandle) + { + const uint32_t ch = tdc.Channel(); + const uint64_t ts = tdc.Timestamp(); + //Also TODO: Make use of picosecond timestamps + //const uint64_t ts_ps = tdc.TimestampPs(); + + switch(ch) + { + case 0: + _tdc_ch0.push_back(ts); + break; + case 1: + _tdc_ch1.push_back(ts); + break; + case 2: + _tdc_ch2.push_back(ts); + break; + //case 3: + // _tdc_ch3.push_back(ts); + // break; + case 4: + _tdc_ch4.push_back(ts); + break; + } + } } - } - tdcHandle.removeProduct(); +} - //---------------------------PTB-----------------------------// +void sbnd::timing::FrameShift::GetPTBTimestamps(const art::Event &e) +{ art::Handle> ptbHandle; - std::vector> ptb_v; e.getByLabel(fPtbDecodeLabel, ptbHandle); - if ((!ptbHandle.isValid() || ptbHandle->size() == 0)){ - throw cet::exception("FrameShift") << "No raw::ptb::sbndptb found w/ tag " << fPtbDecodeLabel << ". Check data quality!"; - } - else{ - art::fill_ptr_vector(ptb_v, ptbHandle); - // HLT Words - unsigned nHLTs = 0; - - for(auto const& ptb : ptb_v) - nHLTs += ptb->GetNHLTriggers(); + if((!ptbHandle.isValid() || ptbHandle->size() == 0)) + throw cet::exception("FrameShift") << "No raw::ptb::sbndptb found w/ tag " << fPtbDecodeLabel << ". Check data quality!"; + else + { + // HLT Words + unsigned nHLTs = 0; + + for(auto const& ptb : *ptbHandle) + nHLTs += ptb.GetNHLTriggers(); - _ptb_hlt_trigger.resize(nHLTs); - _ptb_hlt_timestamp.resize(nHLTs); - _ptb_hlt_trunmask.resize(nHLTs); - _ptb_hlt_unmask_timestamp.resize(nHLTs); - - unsigned hlt_i = 0; //For multiple upbits in trigger words for unmasking - unsigned h_i = 0; //For trigger with bitmask - - for(auto const& ptb : ptb_v){ - for(unsigned i = 0; i < ptb->GetNHLTriggers(); ++i){ - _ptb_hlt_trigger[h_i] = ptb->GetHLTrigger(i).trigger_word; - _ptb_hlt_timestamp[h_i] = ptb->GetHLTrigger(i).timestamp * uint64_t(20); //Units can be found in the Decoder Module - h_i++; + _ptb_hlt_trigger.resize(nHLTs); + _ptb_hlt_timestamp.resize(nHLTs); + _ptb_hlt_trunmask.resize(nHLTs); + _ptb_hlt_unmask_timestamp.resize(nHLTs); + + unsigned hlt_i = 0; //For multiple upbits in trigger words for unmasking + unsigned h_i = 0; //For trigger with bitmask + + for(auto const& ptb : *ptbHandle) + { + for(unsigned i = 0; i < ptb.GetNHLTriggers(); ++i) + { + _ptb_hlt_trigger[h_i] = ptb.GetHLTrigger(i).trigger_word; + _ptb_hlt_timestamp[h_i] = ptb.GetHLTrigger(i).timestamp * uint64_t(20); //Units can be found in the Decoder Module + h_i++; - int val = ptb->GetHLTrigger(i).trigger_word; - int upBit[32]; + int val = ptb.GetHLTrigger(i).trigger_word; + int upBit[32]; - for (int u=0; u<32; u++){ //setting default values for maximum of 32 bits - upBit[u]=-1; - } + for(int u=0; u<32; u++) + upBit[u] = -1; //setting default values for maximum of 32 bits - int numOfTrig =0; - for(int b=0; b<32;b++){ - if ((val & 0x01) ==1){ - upBit[numOfTrig] = b; - numOfTrig++; - } - val = val >> 1; - } + int numOfTrig = 0; + for(int b=0; b<32; b++) + { + if((val & 0x01) == 1) + { + upBit[numOfTrig] = b; + numOfTrig++; + } + + val = val >> 1; + } - if (numOfTrig ==1){ - _ptb_hlt_unmask_timestamp[hlt_i] = _ptb_hlt_timestamp[h_i-1]; - _ptb_hlt_trunmask[hlt_i] = upBit[0]; - hlt_i++; - }//End of if statement for single upbit - else if (numOfTrig > 1){ - nHLTs += (numOfTrig -1); - _ptb_hlt_unmask_timestamp.resize(nHLTs); - _ptb_hlt_trunmask.resize(nHLTs); + if(numOfTrig == 1) + { + _ptb_hlt_unmask_timestamp[hlt_i] = _ptb_hlt_timestamp[h_i-1]; + _ptb_hlt_trunmask[hlt_i] = upBit[0]; + hlt_i++; + } + else if(numOfTrig > 1) + { + nHLTs += (numOfTrig - 1); + _ptb_hlt_unmask_timestamp.resize(nHLTs); + _ptb_hlt_trunmask.resize(nHLTs); - for (int mult =0; mult < numOfTrig; mult++){ - _ptb_hlt_trunmask[hlt_i] = upBit[mult]; - _ptb_hlt_unmask_timestamp[hlt_i] = _ptb_hlt_timestamp[h_i-1]; - hlt_i++; - } //End of loop over multiple upbits - } //End of else statement for multiple triggers - } //End of loop over nHLTriggers - } //End of loop over ptb in ptb_v - } - ptbHandle.removeProduct(); - - if (fDebugPtb){ - for (size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++){ - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "HLT " << _ptb_hlt_trunmask[i] - << " sec (s) = " << _ptb_hlt_unmask_timestamp[i]/uint64_t(1e9) - << ", ts (ns) = " << _ptb_hlt_unmask_timestamp[i]%uint64_t(1e9) - < 1){ - //Get the one closest to the global frame - uint64_t min_diff = std::numeric_limits::max(); - for(auto ts : _tdc_ch0){ - uint64_t diff = _global_frame > ts ? _global_frame - ts : ts - _global_frame; - if(diff < min_diff) - { - min_diff = diff; - _tdc_crtt1_ts = ts; - } - } - } - - // ch1: BES - if (_tdc_ch1.size() == 1){ - _tdc_bes_ts = _tdc_ch1[0]; - } - else if(_tdc_ch1.size() > 1){ - //Get the one closest to the global frame - uint64_t min_diff = std::numeric_limits::max(); - for(auto ts : _tdc_ch1){ - uint64_t diff = _global_frame > ts ? _global_frame - ts : ts - _global_frame; - if(diff < min_diff) - { - min_diff = diff; - _tdc_bes_ts = ts; - } + for(int mult = 0; mult < numOfTrig; mult++) + { + _ptb_hlt_unmask_timestamp[hlt_i] = _ptb_hlt_timestamp[h_i-1]; + _ptb_hlt_trunmask[hlt_i] = upBit[mult]; + hlt_i++; + } + } + } //End of loop over nHLTriggers + } //End of loop over ptb in ptb_v } - } - // ch2: RWM - if (_tdc_ch2.size() == 1){ - _tdc_rwm_ts = _tdc_ch2[0]; - } - else if(_tdc_ch2.size() > 1){ - //Get the one closest to the global frame - uint64_t min_diff = std::numeric_limits::max(); - for(auto ts : _tdc_ch2){ - uint64_t diff = _global_frame > ts ? _global_frame - ts : ts - _global_frame; - if(diff < min_diff) - { - min_diff = diff; - _tdc_rwm_ts = ts; - } + if(fDebugPtb) + { + for(size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++) + { + std::cout << "----------------------------------------------------" << std::endl; + std::cout << "HLT " << _ptb_hlt_trunmask[i] << ": " << PrintFormatTimestamp(_ptb_hlt_unmask_timestamp[i]) << std::endl; + } } - } +} - // ch4: ETRIG - if (_tdc_ch4.size() == 0){ +void sbnd::timing::FrameShift::FindETRIGs() +{ + if(_tdc_ch4.size() == 0) throw cet::exception("FrameShift") << "No TDC ETRIG timestamps found! Check data quality!"; - } - else if (_tdc_ch4.size() == 1){ - _tdc_etrig_ts = _tdc_ch4[0]; - } - else if(_tdc_ch4.size() > 1){ - uint64_t min_diff = std::numeric_limits::max(); - for(auto ts : _tdc_ch4){ - uint64_t diff = _raw_ts > ts ? _raw_ts - ts : ts - _raw_ts; - if(diff < min_diff) - { - min_diff = diff; - _tdc_etrig_ts = ts; - } - } - } - - if (fDebugTdc){ - std::cout << "----------------------------------------------------" << std::endl; - std::cout << "TDC Channel 0 (CRTT1) Timestamp: " - << " (s) = " << _tdc_crtt1_ts/uint64_t(1e9) - << ", (ns) = " << _tdc_crtt1_ts%uint64_t(1e9) - << std::endl; - std::cout << "TDC Channel 1 (BES) Timestamp: " - << " (s) = " << _tdc_bes_ts/uint64_t(1e9) - << ", (ns) = " << _tdc_bes_ts%uint64_t(1e9) - << std::endl; - std::cout << "TDC Channel 2 (RWM) Timestamp: " - << " (s) = " << _tdc_rwm_ts/uint64_t(1e9) - << ", (ns) = " << _tdc_rwm_ts%uint64_t(1e9) - << std::endl; - std::cout << "TDC Channel 4 (ETRIG) Timestamp: " - << " (s) = " << _tdc_etrig_ts/uint64_t(1e9) - << ", (ns) = " << _tdc_etrig_ts%uint64_t(1e9) - << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - } - - //---------------------------PTB Frame-----------------------------// + else + _tdc_etrig_ts = FindClosest(_tdc_ch4, _raw_ts); - //HLT ETRIG //Grab all the ETRIG HLTS -- there might be more than 1 - for (size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++){ - for (size_t j = 0; j < fPtbEtrigHlts.size(); j++){ - if (_ptb_hlt_trunmask[i] == fPtbEtrigHlts[j]){ - _ptb_hlt_etrig.push_back(_ptb_hlt_trunmask[i]); - _ptb_hlt_etrig_ts.push_back(_ptb_hlt_unmask_timestamp[i]); - } + for(size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++) + { + for(size_t j = 0; j < fPtbEtrigHlts.size(); j++) + { + if(_ptb_hlt_trunmask[i] == fPtbEtrigHlts[j]) + { + _ptb_hlt_etrig.push_back(_ptb_hlt_trunmask[i]); + _ptb_hlt_etrig_ts.push_back(_ptb_hlt_unmask_timestamp[i]); + } + } } - } - if (_ptb_hlt_etrig.size() == 0){ + if(_ptb_hlt_etrig.size() == 0) throw cet::exception("FrameShift") << "No HLT ETRIG timestamps found! Check data quality!"; - } - else if (_ptb_hlt_etrig.size() == 1){ - _hlt_etrig = _ptb_hlt_etrig[0]; - _hlt_etrig_ts = _ptb_hlt_etrig_ts[0]; - } - else if(_ptb_hlt_etrig.size() > 1){ - uint64_t min_diff = std::numeric_limits::max(); - for(size_t i = 0; i < _ptb_hlt_etrig.size(); i++){ - uint64_t diff = _raw_ts > _ptb_hlt_etrig_ts[i] ? _raw_ts - _ptb_hlt_etrig_ts[i] : _ptb_hlt_etrig_ts[i] - _raw_ts; - if(diff < min_diff) - { - min_diff = diff; - _hlt_etrig = _ptb_hlt_etrig[i]; - _hlt_etrig_ts = _ptb_hlt_etrig_ts[i]; - } + else + { + uint64_t min_diff = std::numeric_limits::max(); + for(size_t i = 0; i < _ptb_hlt_etrig.size(); i++) + { + uint64_t diff = _raw_ts > _ptb_hlt_etrig_ts[i] ? _raw_ts - _ptb_hlt_etrig_ts[i] : _ptb_hlt_etrig_ts[i] - _raw_ts; + + if(diff < min_diff) + { + min_diff = diff; + _hlt_etrig = _ptb_hlt_etrig[i]; + _hlt_etrig_ts = _ptb_hlt_etrig_ts[i]; + } + } } - } +} - //Check which Gate/CRT T1 HLT to use based on ETRIG HLT - //Order to check: Beam -> Offbeam -> Xmuon - for (size_t i = 0; i < fBeamEtrigHlt.size(); i++){ - if (_hlt_etrig == fBeamEtrigHlt[i]){ - _hlt_gate = fBeamGateHlt; - _hlt_crtt1 = fBeamCrtT1Hlt; - _isBeam = true; - break; - } - } - if (!_isBeam){ - for (size_t i = 0; i < fOffbeamEtrigHlt.size(); i++){ - if (_hlt_etrig == fOffbeamEtrigHlt[i]){ - _hlt_gate = fOffbeamGateHlt; - _hlt_crtt1 = fOffbeamCrtT1Hlt; - _isOffbeam = true; - break; - } - } - } - //TODO: Grab Xmuon HLTs - - //Get Gate and CRT T1 HLT timestamps - //TODO: What if there is no Gate or CRT T1 HLT? - for (size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++){ - if(_ptb_hlt_trunmask[i] == _hlt_gate){ - _hlt_gate_ts = _ptb_hlt_unmask_timestamp[i]; +uint64_t sbnd::timing::FrameShift::DecideGlobalEtrigTimestamp() +{ + //Decide which global frame to use as reference + // Prioritise TDC ETRIG then PTB ETRIG then the raw header. + uint64_t global_etrig_ts = kInvalidTimestamp; + + if(_tdc_etrig_ts != kInvalidTimestamp) + global_etrig_ts = _tdc_etrig_ts; + else if(_hlt_etrig_ts != kInvalidTimestamp) + global_etrig_ts = _hlt_etrig_ts; + else + global_etrig_ts = _raw_ts; + + if(fDebugFrame) + { + std::cout << "----------------------------------------------------" << std::endl; + if(_tdc_etrig_ts != kInvalidTimestamp) + std::cout << "Using TDC ETRIG as Global Frame Reference" << std::endl; + else if(_hlt_etrig_ts != kInvalidTimestamp) + std::cout << "Using PTB HLT ETRIG as Global Frame Reference" << std::endl; + else + std::cout << "Using DAQ Header Timestamp as Global Frame Reference" << std::endl; + + std::cout << "Global Frame Timestamp: " << PrintFormatTimestamp(global_etrig_ts) << std::endl; } - if(_ptb_hlt_trunmask[i] == _hlt_crtt1){ - _hlt_crtt1_ts = _ptb_hlt_unmask_timestamp[i]; + + return global_etrig_ts; +} + +void sbnd::timing::FrameShift::DecideRelevantTDCTimestamps(uint64_t global_etrig_ts) +{ + // ch0: CRT T1 + if(_tdc_ch0.size() != 0) + _tdc_crtt1_ts = FindClosest(_tdc_ch0, global_etrig_ts); + + // ch1: BES + if(_tdc_ch1.size() != 0) + _tdc_bes_ts = FindClosest(_tdc_ch1, global_etrig_ts); + + // ch2: RWM + if(_tdc_ch2.size() != 0) + _tdc_rwm_ts = FindClosest(_tdc_ch2, global_etrig_ts); + + if(fDebugTdc) + { + std::cout << "----------------------------------------------------" << std::endl; + std::cout << "TDC Channel 0 (CRTT1) Timestamp: " << PrintFormatTimestamp(_tdc_crtt1_ts) << std::endl; + std::cout << "TDC Channel 1 (BES) Timestamp: " << PrintFormatTimestamp(_tdc_bes_ts) << std::endl; + std::cout << "TDC Channel 2 (RWM) Timestamp: "<< PrintFormatTimestamp(_tdc_rwm_ts) << std::endl; + std::cout << "TDC Channel 4 (ETRIG) Timestamp: "<< PrintFormatTimestamp(_tdc_etrig_ts) << std::endl; + std::cout << "----------------------------------------------------" << std::endl; + //for (size_t i = 0; i < _tdc_ch3.size(); i++) + // { + // std::cout << "TDC Channel 3 (FTRIG) Timestamp " << i << ": " << PrintFormatTimestamp(_tdc_ch3[i]) << std::endl; + // } + //std::cout << "----------------------------------------------------" << std::endl; } - } - - if (fDebugPtb){ - - std::cout << "----------------------------------------------------" << std::endl; - if (_isBeam) std::cout << "This is Beam Stream!" << std::endl; - if (_isOffbeam) std::cout << "This is Offbeam Stream!" << std::endl; - std::cout << "HLT ETRIG = " << _hlt_etrig - << ", Timestamp: " - << " (s) = " << _hlt_etrig_ts/uint64_t(1e9) - << ", (ns) = " << _hlt_etrig_ts%uint64_t(1e9) - << std::endl; - std::cout << "HLT Gate = " << _hlt_gate - << ", Timestamp: " - << " (s) = " << _hlt_gate_ts/uint64_t(1e9) - << ", (ns) = " << _hlt_gate_ts%uint64_t(1e9) - << std::endl; - std::cout << "HLT CRT T1 = " << _hlt_crtt1 - << ", Timestamp: " - << " (s) = " << _hlt_crtt1_ts/uint64_t(1e9) - << ", (ns) = " << _hlt_crtt1_ts%uint64_t(1e9) - << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - } - - //--------------------------Get the Global Frame--------------------------// - - if (_global_timing_type == 0){ - _global_frame = _tdc_etrig_ts; - }else if(_global_timing_type ==1){ - _global_frame = _hlt_etrig_ts; - } - //TODO: Add for global timing type 2 (nothing) - - uint64_t diff = _raw_ts > _global_frame ? _raw_ts - _global_frame : _global_frame - _raw_ts; - if (diff > fMaxAllowedRefTimeDiff){ - throw cet::exception("FrameShift") << "Global timestamp is too far from the raw timestamp! " - << "Raw TS = (s) " << _raw_ts/uint64_t(1e9) << " (ns) " << _raw_ts%uint64_t(1e9) << std::endl - << ", Gobal Frame TS = (s) " << _global_frame/uint64_t(1e9) << " (ns) " << _global_frame%uint64_t(1e9) - << ", Diff = " << diff - << ". Check data quality!"; - } - if (fDebugFrame){ - std::cout << "----------------------------------------------------" << std::endl; - if(_global_timing_type == 0){ - std::cout << "Global Frame Type = SPECTDC (0)" << std::endl; - }else if(_global_timing_type == 1){ - std::cout << "Global Frame Type = PTB HLT (1)" << std::endl; - }else{ - std::cout << "Global Frame Type = Nothing (2)" << std::endl; +} + +void sbnd::timing::FrameShift::DecideRelevantPTBTimestamps(uint64_t global_etrig_ts) +{ + //Check which Gate/CRT T1 HLT to use based on ETRIG HLT + //Order to check: Beam -> Offbeam -> Xmuon + + for(const int &beam_etrig_hlt : fBeamEtrigHlt) + { + if(_hlt_etrig == beam_etrig_hlt) + { + _hlt_gate = fBeamGateHlt; + _hlt_crtt1 = fBeamCrtT1Hlt; + _isBeam = true; + break; + } } - std::cout << "Global Frame Timestamp: " - << " (s) = " << _global_frame/uint64_t(1e9) - << ", (ns) = " << _global_frame%uint64_t(1e9) - << std::endl; - std::cout << "----------------------------------------------------" << std::endl; - } - - //-----------------------Compute Frame Shift-----------------------// - - //for beam stream - if (_isBeam){ - if(_tdc_crtt1_ts != 0) _frame_tdc_crtt1 = SubtractUTCTimestmap(_global_frame, _tdc_crtt1_ts); - if(_tdc_bes_ts != 0) _frame_tdc_bes = SubtractUTCTimestmap(_global_frame, _tdc_bes_ts); - if(_tdc_rwm_ts != 0) _frame_tdc_rwm = SubtractUTCTimestmap(_global_frame, _tdc_rwm_ts); - if(_hlt_crtt1_ts != 0) _frame_hlt_crtt1 = SubtractUTCTimestmap(_global_frame, _hlt_crtt1_ts); - if(_hlt_gate_ts != 0) _frame_hlt_gate = SubtractUTCTimestmap(_global_frame, _hlt_gate_ts); - } - //for offbeam stream - else if (_isOffbeam){ - if(_tdc_crtt1_ts != 0) _frame_tdc_crtt1 = SubtractUTCTimestmap(_global_frame, _tdc_crtt1_ts); - if(_hlt_crtt1_ts != 0) _frame_hlt_crtt1 = SubtractUTCTimestmap(_global_frame, _hlt_crtt1_ts); - if(_hlt_gate_ts != 0) _frame_hlt_gate = SubtractUTCTimestmap(_global_frame, _hlt_gate_ts); - } - //TODO: for cosmic stream - - if (fDebugFrame){ - std::cout << "--------------------------------------" << std::endl; - std::cout << "Frame TDC CRT T1 = " << std::setprecision(9) << _frame_tdc_crtt1 << " (ns)" << std::endl; - std::cout << "Frame TDC BES = " << std::setprecision(9) << _frame_tdc_bes << " (ns)" << std::endl; - std::cout << "Frame TDC RWM = " << std::setprecision(9) << _frame_tdc_rwm << " (ns)" << std::endl; - std::cout << "--------------------------------------" << std::endl; - std::cout << "Frame HLT CRT T1 = " << std::setprecision(9) << _frame_hlt_crtt1 << " (ns)" << std::endl; - std::cout << "Frame HLT Gate = " << std::setprecision(9) << _frame_hlt_gate << " (ns)" << std::endl; - std::cout << "--------------------------------------" << std::endl; - } - - //-----------------------Pick which one to use in CAF-----------------------// - // The follow picks which frame to apply at the CAF stage and store it as frameApplyAtCaf, based on the stream - // This is so that no decision making is applied CAF. - // - // 1. Beam Stream: recontruct the beam spill + porch relative to RWM record in TDC since it has better resolution. - // Then, apply a constant to shift the RWM fram - // i.e. frame_tdc_rwm + shift rwm to beam gate HLT 26 - // - // 2. Offbeam Stream: reconstruct the porch relative to the beam gate opening frame - // i.e. frame_hlt_gate equiv - // - // 3. Xmuon Stream: do nothing for now - // - // TODO: Align Data Beam and Offbeam with MC with frame_data2mc - // + MC: t = 0 = first proton in spill - // + Data: t = 0 = abitrary. All subsystem electronics time is reference to the last PPS - - if(_isBeam){ - if (_frame_tdc_rwm != 0){ - _frame_tdc_rwm += fShiftRWM2Gate; //align RWM signal to gate opening frame - _frame_apply_at_caf = _frame_tdc_rwm; // +frame_data2mc - }else if(_frame_hlt_gate != 0){ - _frame_apply_at_caf = _frame_hlt_gate; //+frame_data2mc - }else{ - _frame_apply_at_caf = 0; + + if(!_isBeam) + { + for(const int &offbeam_etrig_hlt : fOffbeamEtrigHlt) + { + if(_hlt_etrig == offbeam_etrig_hlt) + { + _hlt_gate = fOffbeamGateHlt; + _hlt_crtt1 = fOffbeamCrtT1Hlt; + _isOffbeam = true; + break; + } + } } - } - else if(_isOffbeam){ - if(_frame_hlt_gate != 0){ - _frame_apply_at_caf = _frame_hlt_gate; //+frame_data2mc - }else{ - _frame_apply_at_caf = 0; + + if(!_isBeam & !_isOffbeam) + { + for(const int &xmuon_etrig_hlt : fXmuonEtrigHlt) + { + if(_hlt_etrig == xmuon_etrig_hlt) + { + _isXmuon = true; + break; + } + } } - } - else if(_isXmuon){ - _frame_apply_at_caf = 0; - } - //Put product in event - std::unique_ptr< sbnd::timing::FrameShiftInfo > newFrameShiftInfo(new sbnd::timing::FrameShiftInfo(_global_timing_type, _frame_tdc_crtt1, _frame_tdc_bes, _frame_tdc_rwm, _frame_hlt_crtt1, _frame_hlt_gate, _frame_apply_at_caf)); - std::unique_ptr< sbnd::timing::TimingInfo > newTimingInfo(new sbnd::timing::TimingInfo(_raw_ts, _tdc_crtt1_ts, _tdc_bes_ts, _tdc_rwm_ts, _tdc_etrig_ts, _hlt_crtt1_ts, _hlt_etrig_ts, _hlt_gate_ts)); + if( !_isBeam & !_isOffbeam & !_isXmuon) + throw cet::exception("FrameShift") << "ETRIG HLT " << _hlt_etrig << " does not match any known Beam/Offbeam/Xmuon ETRIG HLT! Check data quality!"; - e.put(std::move(newTimingInfo)); - e.put(std::move(newFrameShiftInfo)); + //Get Gate and CRT T1 HLT timestamps + for(size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++) + { + if(_ptb_hlt_trunmask[i] == _hlt_gate) + _hlt_gate_ts = _ptb_hlt_unmask_timestamp[i]; + if(_ptb_hlt_trunmask[i] == _hlt_crtt1) + _hlt_crtt1_ts = _ptb_hlt_unmask_timestamp[i]; + } - //Fill the tree - fTree->Fill(); + if(fDebugPtb) + { + std::cout << "----------------------------------------------------" << std::endl; + if(_isBeam) std::cout << "This is Beam Stream!" << std::endl; + if(_isOffbeam) std::cout << "This is Offbeam Stream!" << std::endl; + if(_isXmuon) std::cout << "This is Crossing Muon Stream!" << std::endl; + std::cout << "HLT ETRIG = " << _hlt_etrig << ", Timestamp: " << PrintFormatTimestamp(_hlt_etrig_ts) << std::endl; + std::cout << "HLT Gate = " << _hlt_gate << ", Timestamp: " << PrintFormatTimestamp(_hlt_gate_ts) << std::endl; + std::cout << "HLT CRT T1 = " << _hlt_crtt1 << ", Timestamp: " << PrintFormatTimestamp(_hlt_crtt1_ts) << std::endl; + std::cout << "----------------------------------------------------" << std::endl; + } } -void sbnd::FrameShift::beginJob() +void sbnd::timing::FrameShift::beginJob() { - // Implementation of optional member function here. - //Event Tree - fTree = tfs->make("events", ""); + if(fMakeTree) + { + fTree = tfs->make("events", ""); - fTree->Branch("run", &_run); - fTree->Branch("subrun", &_subrun); - fTree->Branch("event", &_event); + fTree->Branch("run", &_run); + fTree->Branch("subrun", &_subrun); + fTree->Branch("event", &_event); - //TDC - fTree->Branch("tdc_ch0", &_tdc_ch0); - fTree->Branch("tdc_ch1", &_tdc_ch1); - fTree->Branch("tdc_ch2", &_tdc_ch2); - fTree->Branch("tdc_ch4", &_tdc_ch4); - - fTree->Branch("tdc_crtt1_ts", &_tdc_crtt1_ts); - fTree->Branch("tdc_bes_ts", &_tdc_bes_ts); - fTree->Branch("tdc_rwm_ts", &_tdc_rwm_ts); - fTree->Branch("tdc_etrig_ts", &_tdc_etrig_ts); - - //PTB - //fTree->Branch("ptb_hlt_trigger", &_ptb_hlt_trigger); - //fTree->Branch("ptb_hlt_timestamp", &_ptb_hlt_timestamp); - fTree->Branch("ptb_hlt_unmask_timestamp", &_ptb_hlt_unmask_timestamp); - fTree->Branch("ptb_hlt_trunmask", &_ptb_hlt_trunmask); - - fTree->Branch("hlt_etrig", &_hlt_etrig); - fTree->Branch("hlt_etrig_ts", &_hlt_etrig_ts); - fTree->Branch("hlt_gate", &_hlt_gate); - fTree->Branch("hlt_gate_ts", &_hlt_gate_ts); - fTree->Branch("hlt_crtt1", &_hlt_crtt1); - fTree->Branch("hlt_crtt1_ts", &_hlt_crtt1_ts); - - //Timing Reference - fTree->Branch("crt_timing_ch", &_crt_timing_ch); - fTree->Branch("crt_timing_type", &_crt_timing_type); - fTree->Branch("pmt_timing_ch", &_pmt_timing_ch); - fTree->Branch("pmt_timing_type", &_pmt_timing_type); - - //Frame Shift - fTree->Branch("global_timing_type", &_global_timing_type); - fTree->Branch("global_frame", &_global_frame); - fTree->Branch("frame_tdc_crtt1", &_frame_tdc_crtt1); - fTree->Branch("frame_tdc_bes", &_frame_tdc_bes); - fTree->Branch("frame_tdc_rwm", &_frame_tdc_rwm); - fTree->Branch("frame_hlt_crtt1", &_frame_hlt_crtt1); - fTree->Branch("frame_hlt_gate", &_frame_hlt_gate); - fTree->Branch("frame_data2mc", &_frame_data2mc); - fTree->Branch("frame_apply_at_caf", &_frame_apply_at_caf); + //TDC + fTree->Branch("tdc_ch0", &_tdc_ch0); + fTree->Branch("tdc_ch1", &_tdc_ch1); + fTree->Branch("tdc_ch2", &_tdc_ch2); + fTree->Branch("tdc_ch4", &_tdc_ch4); + + fTree->Branch("tdc_crtt1_ts", &_tdc_crtt1_ts); + fTree->Branch("tdc_bes_ts", &_tdc_bes_ts); + fTree->Branch("tdc_rwm_ts", &_tdc_rwm_ts); + fTree->Branch("tdc_etrig_ts", &_tdc_etrig_ts); + + //PTB + fTree->Branch("ptb_hlt_unmask_timestamp", &_ptb_hlt_unmask_timestamp); + fTree->Branch("ptb_hlt_trunmask", &_ptb_hlt_trunmask); + + fTree->Branch("hlt_etrig", &_hlt_etrig); + fTree->Branch("hlt_etrig_ts", &_hlt_etrig_ts); + fTree->Branch("hlt_gate", &_hlt_gate); + fTree->Branch("hlt_gate_ts", &_hlt_gate_ts); + fTree->Branch("hlt_crtt1", &_hlt_crtt1); + fTree->Branch("hlt_crtt1_ts", &_hlt_crtt1_ts); + + //Frame Shift + fTree->Branch("frame_crtt1", &_frame_crtt1); + fTree->Branch("timing_type_crtt1", &_timing_type_crtt1); + fTree->Branch("timing_channel_crtt1", &_timing_channel_crtt1); + + fTree->Branch("frame_gate", &_frame_gate); + fTree->Branch("timing_type_gate", &_timing_type_gate); + fTree->Branch("timing_channel_gate", &_timing_channel_gate); + + fTree->Branch("frame_etrig", &_frame_etrig); + fTree->Branch("timing_type_etrig", &_timing_type_etrig); + fTree->Branch("timing_channel_etrig", &_timing_channel_etrig); + + fTree->Branch("frame_default", &_frame_default); + fTree->Branch("timing_type_default", &_timing_type_default); + fTree->Branch("timing_channel_default", &_timing_channel_default); + } } -void sbnd::FrameShift::ResetEventVars() +void sbnd::timing::FrameShift::ResetEventVars() { _run = -1; _subrun = -1; _event = -1; - _raw_ts = 0; + _raw_ts = kInvalidTimestamp; _tdc_ch0.clear(); _tdc_ch1.clear(); _tdc_ch2.clear(); + //_tdc_ch3.clear(); _tdc_ch4.clear(); - _tdc_crtt1_ts = 0; - _tdc_bes_ts = 0; - _tdc_rwm_ts = 0; - _tdc_etrig_ts = 0; + _tdc_crtt1_ts = kInvalidTimestamp; + _tdc_bes_ts = kInvalidTimestamp; + _tdc_rwm_ts = kInvalidTimestamp; + _tdc_etrig_ts = kInvalidTimestamp; _ptb_hlt_trigger.clear(); _ptb_hlt_timestamp.clear(); @@ -773,57 +853,29 @@ void sbnd::FrameShift::ResetEventVars() _isBeam = false; _isOffbeam = false; _isXmuon = false; - _hlt_etrig = -1; - _hlt_etrig_ts = 0; - _hlt_gate = -1; - _hlt_gate_ts = 0; - _hlt_crtt1= -1; - _hlt_crtt1_ts = 0; - - _crt_timing_ch = 99; - _crt_timing_type = 99; - _pmt_timing_ch = 99; - _pmt_timing_type = 99; - - _global_timing_type = 99; - _global_frame = 0; - - _frame_tdc_crtt1 = 0; - _frame_tdc_bes = 0; - _frame_tdc_rwm = 0; - _frame_hlt_crtt1 = 0; - _frame_hlt_gate = 0; - _frame_apply_at_caf = 0; -} -double sbnd::FrameShift::SubtractUTCTimestmap(const uint64_t& ts1, const uint64_t& ts2) -{ - //Subtract two timestamps in UTC format - //ts1 and ts2 are in nanoseconds in uint64_t format - //Return the difference in nanoseconds as double - - double ts1_s = ts1 / uint64_t(1e9); - double ts1_ns = ts1 % uint64_t(1e9); - double ts2_s = ts2 / uint64_t(1e9); - double ts2_ns = ts2 % uint64_t(1e9); - - //std::cout << std::setprecision(15) << "ts1_s = " << ts1_s << std::setprecision(9)<< ", ts1_ns = " << ts1_ns << std::endl; - //std::cout << std::setprecision(15) << "ts2_s = " << ts2_s << std::setprecision(9)<< ", ts2_ns = " << ts2_ns << std::endl; + _hlt_etrig = std::numeric_limits::max(); + _hlt_etrig_ts = kInvalidTimestamp; + _hlt_gate = std::numeric_limits::max(); + _hlt_gate_ts = kInvalidTimestamp; + _hlt_crtt1= std::numeric_limits::max(); + _hlt_crtt1_ts = kInvalidTimestamp; - double diff_s = 0; - double diff_ns = 0; + _frame_crtt1 = kInvalidTimestamp; + _timing_type_crtt1 = kInvalidType; + _timing_channel_crtt1 = kInvalidChannel; - //If the same PPS, just subtract the nanoseconds - if(ts1_s == ts2_s){ - diff_ns = ts1_ns - ts2_ns; - } - //If ts1 is later than ts2, then subtract the seconds and add the nanoseconds - else{ - diff_s = ts1_s - ts2_s; - diff_ns = diff_s * (double)1e9 + ts1_ns - ts2_ns; - } - - return diff_ns; + _frame_gate = kInvalidTimestamp; + _timing_type_gate = kInvalidType; + _timing_channel_gate = kInvalidChannel; + + _frame_etrig = kInvalidTimestamp; + _timing_type_etrig = kInvalidType; + _timing_channel_etrig = kInvalidChannel; + + _frame_default = kInvalidTimestamp; + _timing_type_default = kInvalidType; + _timing_channel_default = kInvalidChannel; } -DEFINE_ART_MODULE(sbnd::FrameShift) +DEFINE_ART_MODULE(sbnd::timing::FrameShift) diff --git a/sbndcode/Timing/FrameShift/frameshift_sbnd_data.fcl b/sbndcode/Timing/FrameShift/frameshift_sbnd_data.fcl index 739fbc792..9a1e7dad1 100644 --- a/sbndcode/Timing/FrameShift/frameshift_sbnd_data.fcl +++ b/sbndcode/Timing/FrameShift/frameshift_sbnd_data.fcl @@ -1,40 +1,39 @@ BEGIN_PROLOG + frameshift_data: { - module_type: "FrameShift" - - DAQHeaderInstanceLabel: "RawEventHeader" - DAQHeaderModuleLabel: "daq" - RawTSCorrection: 367000 //367 us - MaxAllowedRefTimeDiff: 3000000 //3 ms + module_type: "FrameShift" - TimingRefPmtLabel: "pmtdecoder" - TimingRefCrtLabel: "crtstrips" - - TdcDecodeLabel: "tdcdecoder" - PtbDecodeLabel: "ptbdecoder" + DAQHeaderLabel: "daq:RawEventHeader" + RawTSCorrection: 367000 #367 us + + TdcDecodeLabel: "tdcdecoder" + PtbDecodeLabel: "ptbdecoder" - PtbEtrigHlts: [1, 2, 3, 4, 5, 14, 15] - - BeamEtrigHlt: [1,2] - OffbeamEtrigHlt: [3,4] - XmuonEtrigHlt: [5, 14, 15] + PtbEtrigHlts: [1, 2, 3, 4, 5, 14, 15, 16, 17] + + BeamEtrigHlt: [1,2,16] + OffbeamEtrigHlt: [3,4,17] + XmuonEtrigHlt: [5, 14, 15] - BeamCrtT1Hlt: 20 - OffbeamCrtT1Hlt: 21 + BeamCrtT1Hlt: 20 + OffbeamCrtT1Hlt: 21 - BeamGateHlt: 26 - OffbeamGateHlt: 27 + BeamGateHlt: 26 + OffbeamGateHlt: 27 - DebugDAQHeader: false - DebugTimingRef: false - DebugPtb: false - DebugTdc: false - DebugFrame: false - - ShiftData2MC: 0 - ShiftRWM2Gate: 1738 -} - + DebugDAQHeader: false + DebugPtb: false + DebugTdc: false + DebugFrame: false + ShiftData2MC: 0 + ShiftTdcRwm2PtbGate: 1738 #ns + #TDC RWM arrives 1738 ns later than HLT Gate on average, so shift TDC RWM earlier + ShiftTdcEtrig2PtbEtrig: 257 #ns + #TDC ETRIG arrives 257 ns later than HLT ETRIG on average, so shift TDC ETRIG earlier + ShiftTdcCrtt12PtbCrtt1: 295 #ns + #TDC CRTT1 arrives 295 ns later than HLT CRTT1 on average, so shift TDC CRTT1 earlier +} + END_PROLOG diff --git a/sbndcode/Timing/TimingUtils.cc b/sbndcode/Timing/TimingUtils.cc new file mode 100644 index 000000000..1156f11ee --- /dev/null +++ b/sbndcode/Timing/TimingUtils.cc @@ -0,0 +1,27 @@ +#include "sbndcode/Timing/TimingUtils.h" + +namespace sbnd { + double TimingUtils::SubtractUTCTimestmap(uint64_t ts1, uint64_t ts2) + { + + double ts1_s = ts1 / uint64_t(1e9); + double ts1_ns = ts1 % uint64_t(1e9); + double ts2_s = ts2 / uint64_t(1e9); + double ts2_ns = ts2 % uint64_t(1e9); + + double diff_s = 0; + double diff_ns = 0; + + //If the same PPS, just subtract the nanoseconds + if(ts1_s == ts2_s){ + diff_ns = ts1_ns - ts2_ns; + } + //If ts1 is later than ts2, then subtract the seconds and add the nanoseconds + else{ + diff_s = ts1_s - ts2_s; + diff_ns = diff_s * (double)1e9 + ts1_ns - ts2_ns; + } + + return diff_ns; + } +} \ No newline at end of file diff --git a/sbndcode/Timing/TimingUtils.h b/sbndcode/Timing/TimingUtils.h new file mode 100644 index 000000000..271bc5e6f --- /dev/null +++ b/sbndcode/Timing/TimingUtils.h @@ -0,0 +1,30 @@ +#ifndef TIMINGUTILS_H_SEEN +#define TIMINGUTILS_H_SEEN + +/** + * @file TimingUtils.h + * @brief Common utility functions for timing reconstruction. + */ + +#include "canvas/Persistency/Common/Ptr.h" + +namespace sbnd { + + /** + * @namespace TimingUtils + * @brief Collection of helper utilities for timing calculations. + */ + namespace TimingUtils{ + + /** + * @brief Subtracts two UTC timestamps. + * @param[in] ts1 First timestamp in nanoseconds. + * @param[in] ts2 Second timestamp in nanoseconds. + * @return Difference between `ts1 - ts2` in nanoseconds as double. + */ + double SubtractUTCTimestmap(uint64_t ts1, uint64_t ts2); + + } +} + +#endif \ No newline at end of file diff --git a/sbndcode/Timing/WaveformAlignment/WaveformAlignment_module.cc b/sbndcode/Timing/WaveformAlignment/WaveformAlignment_module.cc index d041b580f..423b7a0eb 100644 --- a/sbndcode/Timing/WaveformAlignment/WaveformAlignment_module.cc +++ b/sbndcode/Timing/WaveformAlignment/WaveformAlignment_module.cc @@ -1,10 +1,11 @@ -//////////////////////////////////////////////////////////////////////// -// Class: WaveformAlignment -// Plugin Type: Producer -// File: WaveformAlignment_module.cc -// -// Author: Lan Nguyen (vclnguyen@ucsb.edu) -//////////////////////////////////////////////////////////////////////// +/** + * @file WaveformAlignment_module.cc + * @brief art::EDProducer that aligns PMT waveforms in time + * by computing per-CAEN-digitiser jitter shifts + * from FTRIG waveforms digitised in channel 15 of each CAEN + * + * @author Lan Nguyen (vclnguyen@ucsb.edu) + */ #include "art/Framework/Core/EDProducer.h" #include "art/Framework/Core/ModuleMacros.h" @@ -39,6 +40,7 @@ #include "lardataobj/RawData/OpDetWaveform.h" #include "sbnobj/SBND/Timing/DAQTimestamp.hh" +#include "sbnobj/SBND/Timing/FrameShiftInfo.hh" #include "sbndcode/Timing/SBNDRawTimingObj.h" #include "sbndcode/Decoders/PTB/sbndptb.h" #include "sbndcode/Calibration/PDSDatabaseInterface/PMTCalibrationDatabase.h" @@ -53,9 +55,8 @@ namespace sbnd { class sbnd::WaveformAlignment : public art::EDProducer { public: + /** @brief Constructor. Reads FHiCL configuration and declares all output data products. */ explicit WaveformAlignment(fhicl::ParameterSet const& p); - // The compiler-generated destructor is fine for non-base - // classes without bare pointers or other resource use. // Plugins should not be copied or assigned. WaveformAlignment(WaveformAlignment const&) = delete; @@ -64,168 +65,241 @@ class sbnd::WaveformAlignment : public art::EDProducer { WaveformAlignment& operator=(WaveformAlignment&&) = delete; // Required functions. + /** + * @brief Per-event processing: reads FrameShift, TDC, PTB, and FTRIG data; + * computes per-board digitiser jitter shifts; applies timing corrections + * to PMT and timing waveforms; and writes aligned products to the event. + * @param[in,out] e The art event being processed. + */ void produce(art::Event& e) override; // Selected optional functions. + /** @brief Called once at the start of the job. Creates the output TTree and registers all branches. */ void beginJob() override; + /** @brief Called once at the end of the job. */ void endJob() override; + /** @brief Resets all per-event member variables to their invalid/default sentinel values. */ void ResetEventVars(); + /** + * @brief Fits the rising edge of an FTRIG waveform with a sigmoid function to find the transition midpoint. + * @param[in] wf Pointer to the FTRIG OpDetWaveform. + * @param[in] flashId Index of the current flash (used for plot naming). + * @param[in] boardId CAEN digitiser board ID (8 = timing board, 0–7 = PMT boards). + * @return Pair of (midpoint tick position, fit convergence flag). + */ std::pair FitFtrig(art::Ptr wf, const int flashId, const int boardId); + /** + * @brief Finds the timing frame timestamp closest to a given FTRIG UTC timestamp. + * @param[in] tsFtrig UTC timestamp of the FTRIG rising edge (ns). + * @param[in] frameVec Vector of reference frame UTC timestamps (ns). + * @param[in] frameOffset Known cable/path offset between the frame and the FTRIG (ns). + * @return UTC timestamp of the nearest frame in @p frameVec. + */ uint64_t FindNearestFrame(const uint64_t tsFtrig, const std::vector frameVec, const double frameOffset); + /** + * @brief Checks whether a UTC second-boundary rollover occurred between two timestamps. + * @param[in] tsFtrig UTC timestamp of the FTRIG (ns). + * @param[in] tsFrame UTC timestamp of the reference frame (ns). + * @return @c true if the two timestamps belong to different UTC seconds. + */ bool IsSecondRollover(const uint64_t tsFtrig, const uint64_t tsFrame); + /** + * @brief Computes the digitiser jitter shift using an integer-nanosecond reference frame. + * @param[in] decodeFrame UTC frame timestamp from the decoded waveform header (ns). + * @param[in] decodeTs Decoded waveform start timestamp (µs). + * @param[in] tickFtrig FTRIG rising-edge tick from the sigmoid fit. + * @param[in] frameVec Vector of reference frame UTC timestamps (ns). + * @param[in] frameOffset Known cable/path offset between frame and FTRIG (ns). + * @param[in] isFrameEarly @c true if the reference frame arrives before the FTRIG. + * @return Jitter shift in nanoseconds. + */ double ComputeShift(const uint64_t decodeFrame, const double decodeTs, const double tickFtrig , const std::vector frameVec , const double frameOffset , const bool isFrameEarly ); + /** + * @brief Computes the digitiser jitter shift using a picosecond-resolution reference frame. + * @param[in] decodeFrame UTC frame timestamp from the decoded waveform header (ns). + * @param[in] decodeTs Decoded waveform start timestamp (µs). + * @param[in] tickFtrig FTRIG rising-edge tick from the sigmoid fit. + * @param[in] frameVec Vector of reference frame UTC timestamps (ns). + * @param[in] frameVecPs Picosecond-precision counterpart of @p frameVec (ps). + * @param[in] frameOffset Known cable/path offset between frame and FTRIG (ns). + * @param[in] isFrameEarly @c true if the reference frame arrives before the FTRIG. + * @return Jitter shift in nanoseconds. + */ double ComputeShiftPs(const uint64_t decodeFrame, const double decodeTs, const double tickFtrig , const std::vector frameVec , const std::vector frameVecPs , const double frameOffset , const bool isFrameEarly ); + /** + * @brief Validates that the computed jitter shift is within the allowed acceptance window. + * @param[in] shift Computed shift value (ns). + * @param[in] boardId CAEN board ID (8 = timing board uses tighter bounds than PMT boards). + * @return @c true if the shift lies within the configured bounds. + */ bool CheckShift(const double shift, const int boardId); + /** + * @brief Saves a two-panel before/after comparison plot of FTRIG waveforms across all boards. + * @param[in] flashId Index of the flash to plot. + */ void PlotFtrigCompare(const int flashId); + /** + * @brief Saves a ROOT canvas showing the FTRIG waveform and its sigmoid rising-edge fit. + * @param[in] g TGraph of the FTRIG waveform ADC values. + * @param[in] fitf Fitted sigmoid TF1. + * @param[in] converged Whether the fit converged. + * @param[in] boardId CAEN board ID (for plot file naming). + * @param[in] flashId Flash index (for plot file naming). + */ void PlotFtrigFit(TGraph *g, TF1 *fitf, const bool converged, const int boardId, const int flashId); + /** @brief Clears all per-flash plotting containers (tick positions and waveform coordinate vectors). */ void ResetPlotVars(); private: //---GLOBAL PARAMETERS enum TimingStatus { - kUndefined = -99, //Something went wrong - kFitFailure = -1, //Fit does not converge - default 0 - kMissingTDC = -11, //Missing hardware information from TDC - kOutOfBoundTDC = -12, //Cannot find correct frame from TDC - kMissingPTB = -21, //Missing hardware information from PTB - kOutOfBoundPTB = -22, //Cannot find correct frame from PTB - kOutOfBoundCAEN = -3, //Cannot align with other CAEN digitisers since - kGoodTDC = 1, //Align using TDC O(2ns) - kGoodPTB = 2, //Align using PTB O(2ns) - kGoodCAEN = 3 //Align with CAEN earlier in the chain O(8ns) + kUndefined = -99, ///< Something went wrong + kFitFailure = -1, ///< Fit does not converge — default shift is 0 + kMissingTDC = -11, ///< Missing hardware information from TDC + kOutOfBoundTDC = -12, ///< Cannot find correct frame from TDC + kMissingPTB = -21, ///< Missing hardware information from PTB + kOutOfBoundPTB = -22, ///< Cannot find correct frame from PTB + kOutOfBoundCAEN = -3, ///< Cannot align with other CAEN digitisers + kGoodTDC = 1, ///< Aligned using TDC (~2 ns precision) + kGoodPTB = 2, ///< Aligned using PTB (~2 ns precision) + kGoodCAEN = 3 ///< Aligned with CAEN earlier in the chain (~8 ns precision) }; // Plotting - std::stringstream _plotName; //raw waveform hist name - std::map tickVec; - std::map> xVec; - std::map> yVec; - std::map boardMidX; - std::map boardMidY; + std::stringstream _plotName; ///< String stream used to build histogram/canvas names + std::map tickVec; ///< Per-board FTRIG rising-edge tick position from the fit + std::map> xVec; ///< Per-board waveform time axis values (ns) + std::map> yVec; ///< Per-board waveform ADC values + std::map boardMidX; ///< Per-board x-coordinate of the fitted midpoint + std::map boardMidY; ///< Per-board y-coordinate of the fitted midpoint // Shifting - std::vector boardJitter[9]; - std::vector boardStatus[9]; - std::vector boardFrame; + std::vector boardJitter[9]; ///< Per-board jitter shift values for each flash (ns); index = board ID + std::vector boardStatus[9]; ///< Per-board alignment status codes (TimingStatus) for each flash + std::vector boardFrame; ///< UTC timestamp(s) of the earliest PMT board FTRIG rising edge (ns) // Board counter - int nTotalBoard; - int nTimingBoard; - int nPmtBoard; - std::vector boardId_v; + int nTotalBoard; ///< Total number of CAEN boards seen in the event + int nTimingBoard; ///< Number of timing CAEN boards (board ID = fPmtBoard.back()) + int nPmtBoard; ///< Number of PMT CAEN boards + std::vector boardId_v; ///< Ordered list of board IDs found in the FTRIG data // Flash counter - int nFtrigFlash; + int nFtrigFlash; ///< Number of FTRIG flashes per board in the event //---TREE PARAMETERS - TTree *fTree; - art::ServiceHandle tfs; + TTree *fTree; ///< Output ROOT TTree + art::ServiceHandle tfs; ///< art TFileService for ROOT output - int _run, _subrun, _event; + int _run; ///< Run number + int _subrun; ///< Subrun number + int _event; ///< Event number // PMT Timing - uint16_t _pmt_timing_type; - uint16_t _pmt_timing_ch; + uint16_t _pmt_timing_type; ///< Timing reference type from FrameShiftInfo (0=TDC, 1=PTB, 2=CAEN) + uint16_t _pmt_timing_ch; ///< Timing reference channel from FrameShiftInfo // TDC stuff - std::vector _tdc_ch3; - std::vector _tdc_ch3_ps; - std::vector _tdc_ch4; + std::vector _tdc_ch3; ///< TDC channel 3 timestamps (ns) — FTRIG from CAEN timing board + std::vector _tdc_ch3_ps; ///< TDC channel 3 timestamps (ps) — picosecond precision if available //PTB stuff - std::vector _ptb_hlt_trigger; - std::vector _ptb_hlt_timestamp; - std::vector _ptb_hlt_unmask_timestamp; - std::vector _ptb_hlt_trunmask; + std::vector _ptb_hlt_trigger; ///< PTB HLT trigger bitmask words + std::vector _ptb_hlt_timestamp; ///< PTB HLT timestamps before bitmask expansion (ns) + std::vector _ptb_hlt_unmask_timestamp; ///< PTB HLT timestamps after per-bit unmasking (ns) + std::vector _ptb_hlt_trunmask; ///< PTB HLT trigger bit indices after unmasking //---SERVICE - sbndDB::PMTCalibrationDatabase const* fPMTCalibrationDatabaseService; + sbndDB::PMTCalibrationDatabase const* fPMTCalibrationDatabaseService; ///< PMT calibration database service //---FHICL CONFIG PARAMETERS // Product label - art::InputTag fTimingRefLabel; - art::InputTag fTdcDecodeLabel; - art::InputTag fPtbDecodeLabel; - art::InputTag fPmtDecodeLabel; - art::InputTag fPmtBoardLabel; - art::InputTag fTimingDecodeLabel; - art::InputTag fTimingBoardLabel; - art::InputTag fFtrigDecodeLabel; - art::InputTag fFtrigBoardLabel; + art::InputTag fFrameShiftLabel; ///< Input tag for the FrameShiftInfo product + art::InputTag fTdcDecodeLabel; ///< Input tag for the TDC DAQTimestamp product + art::InputTag fPtbDecodeLabel; ///< Input tag for the PTB sbndptb product + art::InputTag fPmtDecodeLabel; ///< Input tag for the PMT OpDetWaveform product + art::InputTag fPmtBoardLabel; ///< Input tag for the PMT BoardTimingInfo associations + art::InputTag fTimingDecodeLabel; ///< Input tag for the timing channel OpDetWaveform product + art::InputTag fTimingBoardLabel; ///< Input tag for the timing BoardTimingInfo associations + art::InputTag fFtrigDecodeLabel; ///< Input tag for the FTRIG OpDetWaveform product + art::InputTag fFtrigBoardLabel; ///< Input tag for the FTRIG BoardTimingInfo associations //PTB - std::vector fPtbAllowedHlts; + std::vector fPtbAllowedHlts; ///< List of PTB HLT bit indices accepted as valid timing references + + //FrameShift + uint64_t _frameDefault; ///< Default UTC frame timestamp from FrameShiftInfo (ns) // DAQ - double fWfLength; - int fnChperBoard; - std::vector fPmtBoard; - std::vector fPmtCol; + double fWfLength; ///< Expected waveform length in ticks + int fnChperBoard; ///< Number of channels per CAEN board + std::vector fPmtBoard; ///< CAEN board IDs for PMT boards (0–7) and timing board (8) + std::vector fPmtCol; ///< ROOT colour index per board for overlay plots // Shift configuration - int fTickLbPmt; - int fTickUbPmt; + int fTickLbPmt; ///< Lower search bound (tick) for the PMT FTRIG rising edge + int fTickUbPmt; ///< Upper search bound (tick) for the PMT FTRIG rising edge - int fTickLbTiming; - int fTickUbTiming; + int fTickLbTiming; ///< Lower search bound (tick) for the timing FTRIG rising edge + int fTickUbTiming; ///< Upper search bound (tick) for the timing FTRIG rising edge - int fFitBound; - int fFitTries; - std::vector fGradientInitialGuess; - - double fPmtJitterLowerBound; - double fPmtJitterUpperBound; - double fTimingJitterLowerBound; - double fTimingJitterUpperBound; + int fFitBound; ///< Half-width of the sigmoid fit range in ticks + int fFitTries; ///< Maximum number of re-fit attempts on non-convergence + std::vector fGradientInitialGuess; ///< Initial guess for the sigmoid gradient parameter per board + + double fPmtJitterLowerBound; ///< Minimum acceptable jitter shift for PMT boards (ns) + double fPmtJitterUpperBound; ///< Maximum acceptable jitter shift for PMT boards (ns) + double fTimingJitterLowerBound; ///< Minimum acceptable jitter shift for the timing board (ns) + double fTimingJitterUpperBound; ///< Maximum acceptable jitter shift for the timing board (ns) //Cable length - std::vector fTdc3CaenOffset; - std::vector fPtbCaenOffset; - double fPPSPath; + std::vector fTdc3CaenOffset; ///< Per-board cable-length offset between TDC ch3 and each CAEN board (ns) + std::vector fPtbCaenOffset; ///< Per-board cable-length offset between PTB and each CAEN board (ns) + double fPPSPath; ///< PPS cable propagation delay to be subtracted from all timestamps (ns) // Which Correction To Apply - bool fCorrectPPS; - bool fCorrectPMT2DigitiserCable; - bool fCorrectDigitiserJitter; - bool fCorrectPMTResponse; + bool fCorrectPPS; ///< Apply PPS cable path correction + bool fCorrectPMT2DigitiserCable; ///< Apply PMT-to-digitiser cable length (transit time) correction + bool fCorrectDigitiserJitter; ///< Apply per-board digitiser jitter shift correction + bool fCorrectPMTResponse; ///< Apply per-channel PMT response time correction // Debug - bool fDebugTimeRef; - bool fDebugFtrig; - bool fDebugPmt; - bool fDebugTiming; + bool fDebugTimeRef; ///< Enable verbose printout for timing reference selection + bool fDebugFtrig; ///< Enable verbose printout for FTRIG fitting and shift computation + bool fDebugPmt; ///< Enable verbose printout for PMT waveform alignment + bool fDebugTiming; ///< Enable verbose printout for timing channel waveform alignment // New product labels - std::string fFtrigNewLabel; - std::string fFtrigBoardNewLabel; + std::string fFtrigNewLabel; ///< Instance label for the output aligned FTRIG waveforms + std::string fFtrigBoardNewLabel; ///< Instance label for the output FTRIG–board associations - std::string fPmtNewLabel; - std::string fPmtBoardNewLabel; - std::string fPmtAlignNewLabel; + std::string fPmtNewLabel; ///< Instance label for the output aligned PMT waveforms + std::string fPmtBoardNewLabel; ///< Instance label for the output PMT–board associations + std::string fPmtAlignNewLabel; ///< Instance label for the output PMT–BoardAlignment associations - std::string fTimingNewLabel; - std::string fTimingBoardNewLabel; - std::string fTimingAlignNewLabel; + std::string fTimingNewLabel; ///< Instance label for the output aligned timing waveforms + std::string fTimingBoardNewLabel; ///< Instance label for the output timing–board associations + std::string fTimingAlignNewLabel; ///< Instance label for the output timing–BoardAlignment associations // Plotting - bool fSaveGoodFit; - bool fSaveBadFit; - bool fSaveCompare; - std::string fSavePath; + bool fSaveGoodFit; ///< Save PNG plots for waveforms where the sigmoid fit converged + bool fSaveBadFit; ///< Save PNG plots for waveforms where the sigmoid fit failed + bool fSaveCompare; ///< Save before/after comparison plots of FTRIG waveforms + std::string fSavePath; ///< Directory path for saved plot files }; @@ -233,7 +307,7 @@ sbnd::WaveformAlignment::WaveformAlignment(fhicl::ParameterSet const& p) : EDProducer{p} // // More initializers here. { - fTimingRefLabel = p.get("TimingRefLabel", "pmtdecoder"); + fFrameShiftLabel = p.get("FrameShiftLabel", "frameshift"); fTdcDecodeLabel = p.get("TdcDecodeLabel", "tdcdecoder"); fPtbDecodeLabel = p.get("PtbDecodeLabel", "ptbdecoder"); @@ -368,27 +442,16 @@ void sbnd::WaveformAlignment::produce(art::Event& e) if (fDebugTimeRef | fDebugFtrig | fDebugPmt | fDebugTiming) std::cout <<"#----------RUN " << _run << " SUBRUN " << _subrun << " EVENT " << _event <<"----------#\n"; - //------------------------PMT Timing--------------------------// - art::Handle timingRefHandle; - e.getByLabel(fTimingRefLabel, timingRefHandle); - - if (!timingRefHandle.isValid()){ - throw cet::exception("WaveformAlignment") << "No raw::TimingReferenceInfo found w/ tag " << fTimingRefLabel << ". Check data quality!"; + //------------------------Frame Shift--------------------------// + try { + auto const& frameShift = e.getProduct(fFrameShiftLabel); + _frameDefault = frameShift.FrameDefault(); + _pmt_timing_type = frameShift.TimingTypeDefault(); } - else{ - raw::TimingReferenceInfo const& pmt_timing(*timingRefHandle); - - _pmt_timing_type = pmt_timing.timingType; - _pmt_timing_ch = pmt_timing.timingChannel; - - if (fDebugTimeRef){ - std::cout << "Timing Reference For Decoding PMT" << std::endl; - std::cout << " Type = " << _pmt_timing_type << " (SPECTDC = 0; PTB HLT = 1; CAEN-only = 3)." << std::endl; - std::cout << " Channel = " << _pmt_timing_ch << " (TDC ETRIG = 4; PTB BNB Beam+Light = 2)." << std::endl; - } + catch(art::Exception const&) { + throw cet::exception("WaveformAlignment") << "No sbnd::timing::FrameShiftInfo found w/ tag " << fFrameShiftLabel << ". Check data quality!"; } - timingRefHandle.removeProduct(); - + //---------------------------TDC-----------------------------// art::Handle> tdcHandle; e.getByLabel(fTdcDecodeLabel, tdcHandle); @@ -408,12 +471,9 @@ void sbnd::WaveformAlignment::produce(art::Event& e) const uint32_t ch = tdc.Channel(); const uint64_t ts = tdc.Timestamp(); - if ((ch == 3) || (ch ==4)){ - if (ch == 3) { - _tdc_ch3.push_back(ts); - _tdc_ch3_ps.push_back(tdc.TimestampPs()); - } - if (ch == 4) _tdc_ch4.push_back(ts); + if (ch == 3){ + _tdc_ch3.push_back(ts); + _tdc_ch3_ps.push_back(tdc.TimestampPs()); if (fDebugTimeRef) std::cout << "Ch " << name << " sec (s) = " << ts/uint64_t(1e9) @@ -560,21 +620,21 @@ void sbnd::WaveformAlignment::produce(art::Event& e) }else{ //TDC as reference frame if timing type is 0 if (_pmt_timing_type == 0){ - if (_tdc_ch4.size() == 0 || _tdc_ch3.size() == 0 || _tdc_ch3_ps.size() == 0){ + if ( _tdc_ch3.size() == 0 || _tdc_ch3_ps.size() == 0){ status = kMissingTDC; }else{ double tempShift = 0; if (_tdc_ch3_ps.front() == 0){ //old decode version does not have picosecond - tempShift = ComputeShift( _tdc_ch4.front() - , wf->TimeStamp() + tempShift = ComputeShift( _frameDefault //decodeFrame + , wf->TimeStamp() //decodeTs , tickFtrig , _tdc_ch3 , fTdc3CaenOffset[wf->ChannelNumber()] , false ); }else{ - tempShift = ComputeShiftPs( _tdc_ch4.front() - , wf->TimeStamp() + tempShift = ComputeShiftPs( _frameDefault //decodeFrame + , wf->TimeStamp() //decodeTs , tickFtrig , _tdc_ch3 , _tdc_ch3_ps @@ -597,20 +657,10 @@ void sbnd::WaveformAlignment::produce(art::Event& e) //OR timing type is 0 but cannot find the TDC frame //TDC status error will be replaced with PTB status error if((_pmt_timing_type == 1) | (status == kMissingTDC) | (status == kOutOfBoundTDC)){ - //Find the HLT of ETRIG - uint64_t tsFrame = std::numeric_limits::min(); - for (size_t i = 0; i < _ptb_hlt_unmask_timestamp.size(); i++){ - for (size_t j = 0; j < fPtbAllowedHlts.size(); j++){ - if (_ptb_hlt_trunmask[i] == fPtbAllowedHlts[j]){ - tsFrame = _ptb_hlt_unmask_timestamp[i]; - break; - } - } - } - if(tsFrame == std::numeric_limits::min() || _ptb_hlt_unmask_timestamp.size() == 0){ + if( _ptb_hlt_unmask_timestamp.size() == 0){ status = kMissingPTB; }else{ - double tempShift = ComputeShift( tsFrame + double tempShift = ComputeShift( _frameDefault , wf->TimeStamp() , tickFtrig , _ptb_hlt_unmask_timestamp @@ -877,7 +927,6 @@ void sbnd::WaveformAlignment::beginJob() fTree->Branch("event", &_event); fTree->Branch("tdc_ch3", &_tdc_ch3); fTree->Branch("tdc_ch3_ps", &_tdc_ch3_ps); - fTree->Branch("tdc_ch4", &_tdc_ch4); fTree->Branch("pmt_timing_type", &_pmt_timing_type); fTree->Branch("pmt_timing_ch", &_pmt_timing_ch); @@ -911,12 +960,13 @@ void sbnd::WaveformAlignment::ResetEventVars() { _run = -1; _subrun = -1; _event = -1; + _frameDefault = sbnd::timing::kInvalidTimestamp; + _pmt_timing_type = std::numeric_limits::max(); _pmt_timing_ch = std::numeric_limits::max(); _tdc_ch3.clear(); _tdc_ch3_ps.clear(); - _tdc_ch4.clear(); _ptb_hlt_trigger.clear(); _ptb_hlt_timestamp.clear(); @@ -1168,7 +1218,7 @@ bool sbnd::WaveformAlignment::CheckShift(const double shift, const int boardId) if (boardId == 8){ if ((shift > fTimingJitterUpperBound) | (shift < fTimingJitterLowerBound)){ if (fDebugFtrig) std::cout << " board id = " << boardId - << " has jittering " << shift + << " has jittering " << std::setprecision(9) << shift << " out of bound [" << fTimingJitterLowerBound << ", " << fTimingJitterUpperBound << "]" << std::endl; return false; @@ -1176,7 +1226,7 @@ bool sbnd::WaveformAlignment::CheckShift(const double shift, const int boardId) }else{ if ((shift > fPmtJitterUpperBound) | (shift < fPmtJitterLowerBound)){ if (fDebugFtrig) std::cout << " board id = " << boardId - << " has jittering " << shift + << " has jittering " << std::setprecision(9)<< shift << " out of bound [" << fPmtJitterLowerBound << ", " << fPmtJitterUpperBound << "]" << std::endl; return false; diff --git a/sbndcode/Timing/WaveformAlignment/wfalign_sbnd_data.fcl b/sbndcode/Timing/WaveformAlignment/wfalign_sbnd_data.fcl index 1cf6d4d96..87bfd1f5d 100644 --- a/sbndcode/Timing/WaveformAlignment/wfalign_sbnd_data.fcl +++ b/sbndcode/Timing/WaveformAlignment/wfalign_sbnd_data.fcl @@ -20,7 +20,7 @@ wfalign_data: #Configure input projects #Types of High Level Trigger for data stream - PtbAllowedHlts: [1, 2, 3, 4, 5, 14, 15] + PtbAllowedHlts: [1, 2, 3, 4, 5, 14, 15, 16, 17] #Time interval between TDC ch3 (FTRIG) and ch15 CAEN Tdc3CaenOffset: [101, 101, 101, 101, 101, 101, 101, 101, 101] #ns #Time interval between HLT issued in PTB and ch15 CAEN