From 72a68138d9bbec9da550c72ba9f17d28c4470311 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Mon, 6 Jul 2026 08:59:36 -0500 Subject: [PATCH 1/7] Modify light propagation correction to include particle propagation time --- .../LightPropagationCorrectionAna_module.cc | 10 +- .../LightPropagationCorrection_module.cc | 114 ++++++++++++++++-- .../LightPropagationCorrection_module.hh | 13 +- ...sbnd_lightpropagationcorrection_config.fcl | 4 + 4 files changed, 130 insertions(+), 11 deletions(-) diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrectionAna_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrectionAna_module.cc index 9b09202a2..2f29e3662 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrectionAna_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrectionAna_module.cc @@ -56,6 +56,8 @@ class LightPropagationCorrectionAna : public art::EDAnalyzer { double fNuToFLight; double fNuToFCharge; double fOpFlashT0Corrected; + double fParticlePropagationTime; + double fPhotonPropagationTime; unsigned int _eventID; unsigned int _runID; @@ -95,7 +97,6 @@ void LightPropagationCorrectionAna::analyze(art::Event const& e) return; } - std::vector> correctedopflash_v; art::fill_ptr_vector(correctedopflash_v, correctedopflash_h); @@ -121,11 +122,15 @@ void LightPropagationCorrectionAna::analyze(art::Event const& e) std::cout << "Corrected flash time light only " << correctedopflash->NuToFLight << std::endl; std::cout << "Corrected flash time tpc z corr " << correctedopflash->NuToFCharge << std::endl; std::cout << "Corrected flash time prop corr tpc z corr " << correctedopflash->OpFlashT0Corrected << std::endl; + std::cout << "Particle propagation time " << correctedopflash->ParticlePropagationTime << std::endl; + std::cout << "Photon propagation time " << correctedopflash->PhotonPropagationTime << std::endl; fOpFlashT0 = correctedopflash->OpFlashT0; fNuToFLight = correctedopflash->NuToFLight; fNuToFCharge = correctedopflash->NuToFCharge; fOpFlashT0Corrected = correctedopflash->OpFlashT0Corrected; + fParticlePropagationTime = correctedopflash->ParticlePropagationTime; + fPhotonPropagationTime = correctedopflash->PhotonPropagationTime; fTree->Fill(); } } @@ -144,7 +149,8 @@ void LightPropagationCorrectionAna::beginJob() fTree->Branch("fNuToFLight", &fNuToFLight, "NuToFLight/d"); fTree->Branch("fNuToFCharge", &fNuToFCharge, "NuToFCharge/d"); fTree->Branch("fOpFlashT0Corrected", &fOpFlashT0Corrected, "OpFlashT0Corrected/d"); - + fTree->Branch("fParticlePropagationTime", &fParticlePropagationTime, "ParticlePropagationTime/d"); + fTree->Branch("fPhotonPropagationTime", &fPhotonPropagationTime, "PhotonPropagationTime/d"); } DEFINE_ART_MODULE(LightPropagationCorrectionAna) \ No newline at end of file diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index 03931100c..22f2cc00d 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -21,14 +21,17 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet fVGroupVUV( p.get("VGroupVUV") ), fNuScoreThreshold( p.get("NuScoreThreshold") ), fFMScoreThreshold( p.get("FMScoreThreshold") ), + fMinHitPE ( p.get("MinHitPE") ), + fPreWindow ( p.get("PreWindow") ), + fPostWindow ( p.get("PostWindow") ), fDebug( p.get("Debug", false) ) // // More initializers here. { // 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 + fParticlePropagationTimePerChannel.resize(fNOpChannels, 0.0); // Initialize with zero or any default value + fPhotonPropagationTimePerChannel.resize(fNOpChannels, 0.0); // Initialize with zero or any default value for(unsigned int opch=0; opchTime(); double newFlashTime = 0.0; + double particlePropTime = 0.0; + double photonPropTime = 0.0; for(const auto& lflash : flash_v) { // Get Flash Barycenter double Ycenter, Zcenter, Ywidth, Zwidth; @@ -286,16 +291,21 @@ void sbnd::LightPropagationCorrection::produce(art::Event & e) double flasht0 = lflash.time; // Refine t0 calculation flasht0 = _flasht0calculator->GetFlashT0(lflash.time, GetAssociatedLiteHits(lflash, ophits)); + this->GetSelectedChannelsFlash(lflash.time, GetAssociatedLiteHits(lflash, ophits)); recob::OpFlash flash(flasht0, lflash.time_err, flasht0, ( flasht0) / 1600., lflash.channel_pe, 0, 0, 1, // this are just default values 100., -1., Ycenter, Ywidth, Zcenter, Zwidth); newFlashTime = flasht0; + particlePropTime = GetAverageParticlePropagationTime()/1000; + photonPropTime = GetAveragePhotonPropagationTime()/1000; sbn::CorrectedOpFlashTiming correctedOpFlashTiming; correctedOpFlashTiming.OpFlashT0 = originalFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; correctedOpFlashTiming.NuToFLight = (Zcenter/fSpeedOfLight)/1000; correctedOpFlashTiming.NuToFCharge = (fRecoVz/fSpeedOfLight)/1000; correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; + correctedOpFlashTiming.ParticlePropagationTime = particlePropTime; + correctedOpFlashTiming.PhotonPropagationTime = photonPropTime; correctedOpFlashTimes->emplace_back(std::move(correctedOpFlashTiming)); } @@ -418,7 +428,9 @@ 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 + fParticlePropagationTimePerChannel.assign(fNOpChannels, 0.0); // Reset the particle propagation time vector for each channel + fPhotonPropagationTimePerChannel.assign(fNOpChannels, 0.0); // Reset the photon propagation time vector for each channel fChargeBarycenterX.assign(2, 0.0); fChargeBarycenterY.assign(2, 0.0); fChargeBarycenterZ.assign(2, 0.0); @@ -431,13 +443,17 @@ void sbnd::LightPropagationCorrection::ResetSliceInfo() void sbnd::LightPropagationCorrection::GetPropagationTimeCorrectionPerChannel() { // Implementation - for(size_t opdet = 0; opdet < fOpDetID.size(); ++opdet) { +for(size_t opdet = 0; opdet < fOpDetID.size(); ++opdet) { double _opDetX = fOpDetX[opdet]; double _opDetY = fOpDetY[opdet]; double _opDetZ = fOpDetZ[opdet]; float minPropTime = 999999999.; + float minPartPropTime = 999999999.; + float minLightPropTime = 999999999.; for(size_t sp=0; sp 0; + if(!isInSameTPC) continue; // Skip points not in the same TPC double dx = fSpacePointX[sp] - _opDetX; double dy = fSpacePointY[sp] - _opDetY; double dz = fSpacePointZ[sp] - _opDetZ; @@ -454,9 +470,15 @@ void sbnd::LightPropagationCorrection::GetPropagationTimeCorrectionPerChannel() 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; + if(PropTime < minPropTime) { + minPropTime = PropTime; + minPartPropTime = partPropTime; + minLightPropTime = lightPropTime; + } } fTimeCorrectionPerChannel[opdet] = -minPropTime; + fParticlePropagationTimePerChannel[opdet] = minPartPropTime; + fPhotonPropagationTimePerChannel[opdet] = minLightPropTime; } } @@ -535,7 +557,7 @@ void sbnd::LightPropagationCorrection::FillCorrectionTree(double & newFlashTime, ::lightana::LiteOpHitArray_t sbnd::LightPropagationCorrection::GetAssociatedLiteHits(::lightana::LiteOpFlash_t lite_flash, ::lightana::LiteOpHitArray_t lite_hits_v) - { +{ ::lightana::LiteOpHitArray_t flash_hits_v; for(auto const& hitidx : lite_flash.asshit_idx) { @@ -543,4 +565,80 @@ ::lightana::LiteOpHitArray_t sbnd::LightPropagationCorrection::GetAssociatedLite } return flash_hits_v; - } \ No newline at end of file +} + +double sbnd::LightPropagationCorrection::GetAverageParticlePropagationTime() +{ + double sum = 0.0; + int n = 0; + for (size_t i=0; i > selected_hits; + double pe_sum = 0.0; + + // Limpiar vector de salida por si ya tenía contenido + fSelectedChannelList.clear(); + // Fill vector with selected hits in the specified window + for (auto const& hit : ophit_list) { + + if (hit.peak_time < flash_time + fPostWindow && + hit.peak_time > flash_time - fPreWindow && + hit.pe > fMinHitPE) { + + // (PE, peak_time, channel) + selected_hits.emplace_back(hit.pe, hit.peak_time, hit.channel); + pe_sum += hit.pe; + } + } + + // Sort by PE in descending order + std::sort( + selected_hits.begin(), + selected_hits.end(), + [](auto const& a, auto const& b) { + return std::get<0>(a) > std::get<0>(b); + } + ); + + double pe_count = 0.0; + + // Loop over selected ophits + for (size_t ix = 0; ix < selected_hits.size(); ix++) { + + double pe = std::get<0>(selected_hits[ix]); + size_t channel = std::get<2>(selected_hits[ix]); + + // Guardar channel ID + fSelectedChannelList.push_back(channel); + + pe_count += pe; + if (pe_count / pe_sum > .5) + break; + } +} diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index 26e6d47c7..f0c1ea94b 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -142,7 +142,9 @@ private: void FillLiteOpHit(std::vector const& , std::vector<::lightana::LiteOpHit_t>& ); void FillCorrectionTree(double & , recob::OpFlash const& , std::vector const& , std::vector const& ); ::lightana::LiteOpHitArray_t GetAssociatedLiteHits(::lightana::LiteOpFlash_t , ::lightana::LiteOpHitArray_t ); - + void GetSelectedChannelsFlash(double , ::lightana::LiteOpHitArray_t ); + double GetAverageParticlePropagationTime(); + double GetAveragePhotonPropagationTime(); geo::WireReadoutGeom const& fWireReadout = art::ServiceHandle()->Get(); opdet::sbndPDMapAlg fPDSMap; @@ -177,6 +179,9 @@ private: bool fSaveCorrectionTree; std::vector fTimeCorrectionPerChannel; + std::vector fParticlePropagationTimePerChannel; + std::vector fPhotonPropagationTimePerChannel; + std::vector fSelectedChannelList; double fRecoVx = 0.0; double fRecoVy = 0.0; double fRecoVz = 0.0; @@ -207,6 +212,12 @@ private: double fNuScoreThreshold; double fFMScoreThreshold; + // Parameters for the selected channel computation + int fMinHitPE; + double fPreWindow; + double fPostWindow; + + bool fDebug; art::ServiceHandle tfs; diff --git a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl index 3277a7bb8..2205072a7 100644 --- a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl +++ b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl @@ -35,6 +35,10 @@ LightPropagationCorrection: NuScoreThreshold: 0.1 FMScoreThreshold: 0.001 + MinHitPE: @local::FlashT0SelectedChannels.MinHitPE + PreWindow: @local::FlashT0SelectedChannels.PreWindow + PostWindow: @local::FlashT0SelectedChannels.PostWindow + FlashMatchingTool: "BarycenterFM" } From 6f3747649c389c72775a998aed6723a0068f494e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Mon, 6 Jul 2026 09:05:07 -0500 Subject: [PATCH 2/7] Add pdfraction as parameter --- .../LightPropagationCorrection_module.cc | 3 ++- .../LightPropagationCorrection_module.hh | 1 + .../job/sbnd_lightpropagationcorrection_config.fcl | 1 + 3 files changed, 4 insertions(+), 1 deletion(-) diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index 22f2cc00d..d1498c3c7 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -24,6 +24,7 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet fMinHitPE ( p.get("MinHitPE") ), fPreWindow ( p.get("PreWindow") ), fPostWindow ( p.get("PostWindow") ), + fPDFraction ( p.get("PDFraction") ), fDebug( p.get("Debug", false) ) // // More initializers here. @@ -638,7 +639,7 @@ void sbnd::LightPropagationCorrection::GetSelectedChannelsFlash( fSelectedChannelList.push_back(channel); pe_count += pe; - if (pe_count / pe_sum > .5) + if (pe_count / pe_sum > fPDFraction) break; } } diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index f0c1ea94b..356f7caa0 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -216,6 +216,7 @@ private: int fMinHitPE; double fPreWindow; double fPostWindow; + double fPDFraction; bool fDebug; diff --git a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl index 2205072a7..ecfd6e2e8 100644 --- a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl +++ b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl @@ -38,6 +38,7 @@ LightPropagationCorrection: MinHitPE: @local::FlashT0SelectedChannels.MinHitPE PreWindow: @local::FlashT0SelectedChannels.PreWindow PostWindow: @local::FlashT0SelectedChannels.PostWindow + PDFraction: @local::FlashT0SelectedChannels.PDFraction FlashMatchingTool: "BarycenterFM" } From 04271a37e2f1c539beef010aa888c41c44ebe6e5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Tue, 7 Jul 2026 04:11:57 -0500 Subject: [PATCH 3/7] Clean code adding correctopflash function --- .../LightPropagationCorrection_module.cc | 141 +++++++++--------- .../LightPropagationCorrection_module.hh | 5 +- 2 files changed, 71 insertions(+), 75 deletions(-) diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index d1498c3c7..f0cf6d3e3 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -117,8 +117,8 @@ void sbnd::LightPropagationCorrection::produce(art::Event & e) //PFP to space points art::FindManyP pfp_sp_assns(pfpHandle, e, fSpacePointLabel); //OpFlash to OpHit - art::FindManyP flashToOpHitAssns_tpc0(opflashListHandle_tpc0, e, fOpFlashLabel_tpc0); - art::FindManyP flashToOpHitAssns_tpc1(opflashListHandle_tpc1, e, fOpFlashLabel_tpc1); + flashToOpHitAssns_tpc0 = std::make_unique>( opflashListHandle_tpc0, e, fOpFlashLabel_tpc0); + flashToOpHitAssns_tpc1 = std::make_unique>(opflashListHandle_tpc1, e, fOpFlashLabel_tpc1); // PFP Metadata art::FindManyP pfp_to_metadata(pfpHandle, e, fReco2Label); @@ -254,70 +254,14 @@ void sbnd::LightPropagationCorrection::produce(art::Event & e) } else throw art::Exception(art::errors::LogicError) << " Flash matching tool " << fFlashMatchingTool << " not supported ." << std::endl; - - // Get the ophits associated to the flash - std::vector> ophitlist; - if(flashFM[0]->XCenter()<0) - { - ophitlist = flashToOpHitAssns_tpc0.at(flashFM[0].key()); - _mgr = _mgr_tpc0; // Use the TPC 0 flash finder manager - } - else - { - ophitlist = flashToOpHitAssns_tpc1.at(flashFM[0].key()); - _mgr = _mgr_tpc1; // Use the TPC 1 flash finder manager - } - - std::vector newOpHitList; - std::vector oldOpHitList; - for(const auto& ophit : ophitlist) { - oldOpHitList.push_back(*ophit); - } - // Get the list of the corrected OpHits - this->CorrectOpHitTime(ophitlist, newOpHitList); - // Create the list of ophit lite to be used in the flash finder - ::lightana::LiteOpHitArray_t ophits; - this->FillLiteOpHit(newOpHitList, ophits); - // Create the flash manager - auto const flash_v = _mgr.RecoFlash(ophits); - double originalFlashTime = flashFM[0]->Time(); - double newFlashTime = 0.0; - double particlePropTime = 0.0; - double photonPropTime = 0.0; - for(const auto& lflash : flash_v) { - // Get Flash Barycenter - double Ycenter, Zcenter, Ywidth, Zwidth; - _flashgeo->GetFlashLocation(lflash.channel_pe, Ycenter, Zcenter, Ywidth, Zwidth); - // Get flasht0 - double flasht0 = lflash.time; - // Refine t0 calculation - flasht0 = _flasht0calculator->GetFlashT0(lflash.time, GetAssociatedLiteHits(lflash, ophits)); - this->GetSelectedChannelsFlash(lflash.time, GetAssociatedLiteHits(lflash, ophits)); - recob::OpFlash flash(flasht0, lflash.time_err, flasht0, - ( flasht0) / 1600., lflash.channel_pe, - 0, 0, 1, // this are just default values - 100., -1., Ycenter, Ywidth, Zcenter, Zwidth); - newFlashTime = flasht0; - particlePropTime = GetAverageParticlePropagationTime()/1000; - photonPropTime = GetAveragePhotonPropagationTime()/1000; - sbn::CorrectedOpFlashTiming correctedOpFlashTiming; - correctedOpFlashTiming.OpFlashT0 = originalFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; - correctedOpFlashTiming.NuToFLight = (Zcenter/fSpeedOfLight)/1000; - correctedOpFlashTiming.NuToFCharge = (fRecoVz/fSpeedOfLight)/1000; - correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; - correctedOpFlashTiming.ParticlePropagationTime = particlePropTime; - correctedOpFlashTiming.PhotonPropagationTime = photonPropTime; - correctedOpFlashTimes->emplace_back(std::move(correctedOpFlashTiming)); - } - + sbn::CorrectedOpFlashTiming correctedOpFlashTiming; + CorrectOpFlash(flashFM[0], correctedOpFlashTiming); + correctedOpFlashTimes->emplace_back(std::move(correctedOpFlashTiming)); art::Ptr newCorrectedOpFlashTimingPtr = make_correctedopflashtime_ptr(correctedOpFlashTimes->size()-1); newCorrectedOpFlashTimingSliceAssn->addSingle(slice, newCorrectedOpFlashTimingPtr); newCorrectedOpFlashTimingOpFlashAssn->addSingle(flashFM[0], newCorrectedOpFlashTimingPtr); - if(fSaveCorrectionTree){ - this->FillCorrectionTree(newFlashTime, *flashFM[0], oldOpHitList, newOpHitList); - } } - if(fSaveCorrectionTree) fTree->Fill(); + ResetEventVars(); e.put(std::move(correctedOpFlashTimes)); e.put(std::move(newCorrectedOpFlashTimingSliceAssn)); @@ -599,11 +543,8 @@ void sbnd::LightPropagationCorrection::GetSelectedChannelsFlash( double flash_time, ::lightana::LiteOpHitArray_t ophit_list) { - // Ahora guardamos también el channel ID std::vector< std::tuple > selected_hits; double pe_sum = 0.0; - - // Limpiar vector de salida por si ya tenía contenido fSelectedChannelList.clear(); // Fill vector with selected hits in the specified window for (auto const& hit : ophit_list) { @@ -619,14 +560,7 @@ void sbnd::LightPropagationCorrection::GetSelectedChannelsFlash( } // Sort by PE in descending order - std::sort( - selected_hits.begin(), - selected_hits.end(), - [](auto const& a, auto const& b) { - return std::get<0>(a) > std::get<0>(b); - } - ); - + std::sort(selected_hits.begin(), selected_hits.end(), [](auto const& a, auto const& b) {return std::get<0>(a) > std::get<0>(b);}); double pe_count = 0.0; // Loop over selected ophits @@ -635,7 +569,7 @@ void sbnd::LightPropagationCorrection::GetSelectedChannelsFlash( double pe = std::get<0>(selected_hits[ix]); size_t channel = std::get<2>(selected_hits[ix]); - // Guardar channel ID + // Save channel ID fSelectedChannelList.push_back(channel); pe_count += pe; @@ -643,3 +577,62 @@ void sbnd::LightPropagationCorrection::GetSelectedChannelsFlash( break; } } + +void sbnd::LightPropagationCorrection::CorrectOpFlash(art::Ptr const& flash, sbn::CorrectedOpFlashTiming &correctedOpFlashTiming) +{ + // Get the ophits associated to the flash + std::vector> ophitlist; + if(flash->XCenter()<0) + { + ophitlist = flashToOpHitAssns_tpc0->at(flash.key()); + _mgr = _mgr_tpc0; // Use the TPC 0 flash finder manager + } + else + { + ophitlist = flashToOpHitAssns_tpc1->at(flash.key()); + _mgr = _mgr_tpc1; // Use the TPC 1 flash finder manager + } + std::vector newOpHitList; + std::vector oldOpHitList; + for(const auto& ophit : ophitlist) { + oldOpHitList.push_back(*ophit); + } + // Get the list of the corrected OpHits + this->CorrectOpHitTime(ophitlist, newOpHitList); + // Create the list of ophit lite to be used in the flash finder + ::lightana::LiteOpHitArray_t ophits; + this->FillLiteOpHit(newOpHitList, ophits); + // Create the flash manager + auto const flash_v = _mgr.RecoFlash(ophits); + double originalFlashTime = flash->Time(); + double newFlashTime = 0.0; + double particlePropTime = 0.0; + double photonPropTime = 0.0; + for(const auto& lflash : flash_v) { + // Get Flash Barycenter + double Ycenter, Zcenter, Ywidth, Zwidth; + _flashgeo->GetFlashLocation(lflash.channel_pe, Ycenter, Zcenter, Ywidth, Zwidth); + // Get flasht0 + double flasht0 = lflash.time; + // Refine t0 calculation + flasht0 = _flasht0calculator->GetFlashT0(lflash.time, GetAssociatedLiteHits(lflash, ophits)); + this->GetSelectedChannelsFlash(lflash.time, GetAssociatedLiteHits(lflash, ophits)); + recob::OpFlash flash(flasht0, lflash.time_err, flasht0, + ( flasht0) / 1600., lflash.channel_pe, + 0, 0, 1, // this are just default values + 100., -1., Ycenter, Ywidth, Zcenter, Zwidth); + newFlashTime = flasht0; + particlePropTime = GetAverageParticlePropagationTime()/1000; + photonPropTime = GetAveragePhotonPropagationTime()/1000; + correctedOpFlashTiming.OpFlashT0 = originalFlashTime; + correctedOpFlashTiming.NuToFLight = (Zcenter/fSpeedOfLight)/1000; + correctedOpFlashTiming.NuToFCharge = (fRecoVz/fSpeedOfLight)/1000; + correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime; + correctedOpFlashTiming.ParticlePropagationTime = particlePropTime; + correctedOpFlashTiming.PhotonPropagationTime = photonPropTime; + } + + if(fSaveCorrectionTree){ + this->FillCorrectionTree(newFlashTime, *flash, oldOpHitList, newOpHitList); + } +} diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index 356f7caa0..f8e1474e5 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -145,9 +145,13 @@ private: void GetSelectedChannelsFlash(double , ::lightana::LiteOpHitArray_t ); double GetAverageParticlePropagationTime(); double GetAveragePhotonPropagationTime(); + void CorrectOpFlash(art::Ptr const& , sbn::CorrectedOpFlashTiming &); geo::WireReadoutGeom const& fWireReadout = art::ServiceHandle()->Get(); opdet::sbndPDMapAlg fPDSMap; + // Flash ophit associations + std::unique_ptr> flashToOpHitAssns_tpc0; + std::unique_ptr> flashToOpHitAssns_tpc1; //Flash finder manager ::lightana::FlashFinderManager _mgr; ::lightana::FlashFinderManager _mgr_tpc0; @@ -218,7 +222,6 @@ private: double fPostWindow; double fPDFraction; - bool fDebug; art::ServiceHandle tfs; From b9efe8bcee7182d804bebb78d6ceb60c661f307f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Tue, 7 Jul 2026 07:27:37 -0500 Subject: [PATCH 4/7] Add lightprop correction in default MC reco2 workflow --- .../JobConfigurations/standard/reco/config/workflow_reco2.fcl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl b/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl index fd2083d46..1654bf746 100755 --- a/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl +++ b/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl @@ -21,6 +21,7 @@ #include "SBNDCVNMapper.fcl" #include "blipreco_configs.fcl" #include "lightcalo_sbnd.fcl" +#include "sbnd_lightpropagationcorrection_config.fcl" BEGIN_PROLOG @@ -82,6 +83,7 @@ sbnd_reco2_producers:{ opt0finderSCE: @local::sbnd_opt0_finder_one_to_many tpcpmtbarycentermatching: @local::TPCPMTBarycenterMatchProducerMC tpcpmtbarycentermatchingSCE: @local::TPCPMTBarycenterMatchProducerSCEMC + lightpropagationcorrectionSCE: @local::LightPropagationCorrectionSCE ### Uncalibrated calorimetry producer for calibration caloskimmer caloskimCalorimetry: @local::caloskim_calorimetry @@ -136,6 +138,7 @@ sbnd_reco2_producer_sequence: [ #, fmatchoparaSCE , opt0finderSCE , tpcpmtbarycentermatchingSCE + , lightpropagationcorrectionSCE , cvn , cvnSCE , blipreco From 990f5b3fa5e7b7c807db8caccd0ab1bd428c7584 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Tue, 21 Jul 2026 09:48:54 -0500 Subject: [PATCH 5/7] fcl refactoring --- sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd.fcl | 5 ++++- .../standard/caf/cafmakerjob_sbnd_data_sce.fcl | 1 - 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd.fcl b/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd.fcl index aa5aeafc7..ab521a814 100644 --- a/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd.fcl +++ b/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd.fcl @@ -156,9 +156,12 @@ physics.producers.cafmaker.FlashMatchSCECryoSuffixes: [""] physics.producers.cafmaker.TrackHitFillRRStartCut: 1000 physics.producers.cafmaker.TrackHitFillRREndCut: 1000 -# Include 3D barycenter flahs matching +# Include 3D barycenter flash matching physics.producers.cafmaker.TPCPMTBarycenterMatchLabel: "tpcpmtbarycentermatching" +# Include propagation correction on the flash timing +physics.producers.cafmaker.CorrectedOpFlashLabel: "lightpropagationcorrectionSCE" + # Include MC version of PMT Software Trigger physics.producers.cafmaker.SBNDSoftwareTriggerLabel: "pmtmetricmc" diff --git a/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd_data_sce.fcl b/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd_data_sce.fcl index 6aad2cdf6..552c9b465 100644 --- a/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd_data_sce.fcl +++ b/sbndcode/JobConfigurations/standard/caf/cafmakerjob_sbnd_data_sce.fcl @@ -14,7 +14,6 @@ physics.producers.cafmaker.SBNDCRTTrackMatchLabel: "crttrackmatchingSCE" physics.producers.cafmaker.OpT0Label: "opt0finderSCE" physics.producers.cafmaker.TPCPMTBarycenterMatchLabel: "tpcpmtbarycentermatchingSCE" physics.producers.cafmaker.CVNLabel: "cvnSCE" -physics.producers.cafmaker.CorrectedOpFlashLabel: "lightpropagationcorrectionSCE" physics.producers.cafmaker.LightCaloLabel: "lightcaloSCE" physics.producers.cnnid.WireLabel: "sptpc2d:dnnsp" From cc9820b5dffca36cf2f95bd7b74634682acc44a4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Wed, 29 Jul 2026 08:50:48 -0500 Subject: [PATCH 6/7] Change light prop module to not read RWM and ETT in MC --- .../standard/reco/config/workflow_reco2.fcl | 2 +- .../LightPropagationCorrection_module.cc | 42 ++++++++++--------- .../LightPropagationCorrection_module.hh | 2 + ...sbnd_lightpropagationcorrection_config.fcl | 4 ++ 4 files changed, 30 insertions(+), 20 deletions(-) diff --git a/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl b/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl index 1654bf746..1c2bf07d9 100755 --- a/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl +++ b/sbndcode/JobConfigurations/standard/reco/config/workflow_reco2.fcl @@ -83,7 +83,7 @@ sbnd_reco2_producers:{ opt0finderSCE: @local::sbnd_opt0_finder_one_to_many tpcpmtbarycentermatching: @local::TPCPMTBarycenterMatchProducerMC tpcpmtbarycentermatchingSCE: @local::TPCPMTBarycenterMatchProducerSCEMC - lightpropagationcorrectionSCE: @local::LightPropagationCorrectionSCE + lightpropagationcorrectionSCE: @local::LightPropagationCorrectionSCEMC ### Uncalibrated calorimetry producer for calibration caloskimmer caloskimCalorimetry: @local::caloskim_calorimetry diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index f0cf6d3e3..7a625ff5b 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -25,6 +25,7 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet fPreWindow ( p.get("PreWindow") ), fPostWindow ( p.get("PostWindow") ), fPDFraction ( p.get("PDFraction") ), + fIsMC ( p.get("IsMC") ), fDebug( p.get("Debug", false) ) // // More initializers here. @@ -201,21 +202,24 @@ void sbnd::LightPropagationCorrection::produce(art::Event & e) // 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> 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 c newFlashTime = flasht0; particlePropTime = GetAverageParticlePropagationTime()/1000; photonPropTime = GetAveragePhotonPropagationTime()/1000; - correctedOpFlashTiming.OpFlashT0 = originalFlashTime; + correctedOpFlashTiming.OpFlashT0 = originalFlashTime + fEventTriggerTime/1000 - fRWMTime/1000; correctedOpFlashTiming.NuToFLight = (Zcenter/fSpeedOfLight)/1000; correctedOpFlashTiming.NuToFCharge = (fRecoVz/fSpeedOfLight)/1000; - correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime; + correctedOpFlashTiming.OpFlashT0Corrected = newFlashTime + fEventTriggerTime/1000 - fRWMTime/1000;; correctedOpFlashTiming.ParticlePropagationTime = particlePropTime; correctedOpFlashTiming.PhotonPropagationTime = photonPropTime; } diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index f8e1474e5..9e3062c91 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -222,6 +222,8 @@ private: double fPostWindow; double fPDFraction; + bool fIsMC; // Flag to indicate if the data is from MC or real data + bool fDebug; art::ServiceHandle tfs; diff --git a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl index ecfd6e2e8..c462c5aea 100644 --- a/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl +++ b/sbndcode/LightPropagationCorrection/job/sbnd_lightpropagationcorrection_config.fcl @@ -41,6 +41,8 @@ LightPropagationCorrection: PDFraction: @local::FlashT0SelectedChannels.PDFraction FlashMatchingTool: "BarycenterFM" + + IsMC: false } LightPropagationCorrection.AlgoConfig_tpc0.TPC: 0 @@ -52,5 +54,7 @@ LightPropagationCorrectionSCE.SpacePointLabel: "pandoraSCE" LightPropagationCorrectionSCE.OpT0FinderModuleLabel: "opt0finderSCE" LightPropagationCorrectionSCE.TPCPMTBarycenterFMModuleLabel: "tpcpmtbarycentermatchingSCE" +LightPropagationCorrectionSCEMC: @local::LightPropagationCorrectionSCE +LightPropagationCorrectionSCEMC.IsMC: true END_PROLOG \ No newline at end of file From 8e335eba87bcb51afe5d00301209df3d984f7050 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Alejandro=20S=C3=A1nchez=20Castillo?= Date: Wed, 29 Jul 2026 10:29:15 -0500 Subject: [PATCH 7/7] Fix issue with nopchannels --- .../LightPropagationCorrection_module.cc | 1 + .../LightPropagationCorrection_module.hh | 1 + 2 files changed, 2 insertions(+) diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc index 7a625ff5b..8be1d4cd9 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.cc @@ -30,6 +30,7 @@ sbnd::LightPropagationCorrection::LightPropagationCorrection(fhicl::ParameterSet // // More initializers here. { + fNOpChannels = fWireReadout.NOpChannels(); // Initialize the TimeCorrectionVector PerChannel fTimeCorrectionPerChannel.resize(fNOpChannels, 0.0); // Initialize with zero or any default value fParticlePropagationTimePerChannel.resize(fNOpChannels, 0.0); // Initialize with zero or any default value diff --git a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh index 9e3062c91..51e0c8d38 100644 --- a/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh +++ b/sbndcode/LightPropagationCorrection/LightPropagationCorrection_module.hh @@ -232,6 +232,7 @@ private: int fEvent; int fRun; int fSubrun; + size_t fNOpChannels; double _fNuScore; double _fFMScore; double fEventTriggerTime=-999999.;