From 49efa904b7ac7807584beaa5abc2b1605b1a0423 Mon Sep 17 00:00:00 2001 From: Marc Schefer Date: Tue, 11 Nov 2025 17:21:22 +0100 Subject: [PATCH 1/4] Improve debug log --- .../SEFramework/Pipeline/SourceGrouping.h | 1 + .../src/lib/Pipeline/SourceGrouping.cpp | 13 ++++++++ .../src/lib/Source/EntangledSource.cpp | 7 +++-- .../Source/SourceWithOnDemandProperties.cpp | 1 + .../SEImplementation/Grouping/AssocGrouping.h | 1 + .../Grouping/MoffatGrouping.h | 2 ++ .../Grouping/SplitSourcesGrouping.h | 2 +- .../src/lib/Grouping/AssocGrouping.cpp | 17 +++++++++- .../src/lib/Grouping/MoffatGrouping.cpp | 13 ++++++++ .../src/lib/Grouping/SplitSourcesGrouping.cpp | 31 ++++++++++++++----- .../Measurement/MultithreadedMeasurement.cpp | 14 ++++++++- .../FlexibleModelFittingIterativeTask.cpp | 7 +++++ .../src/lib/Prefetcher/Prefetcher.cpp | 5 +++ 13 files changed, 101 insertions(+), 13 deletions(-) diff --git a/SEFramework/SEFramework/Pipeline/SourceGrouping.h b/SEFramework/SEFramework/Pipeline/SourceGrouping.h index fc744f9f2..07d9c6419 100644 --- a/SEFramework/SEFramework/Pipeline/SourceGrouping.h +++ b/SEFramework/SEFramework/Pipeline/SourceGrouping.h @@ -121,6 +121,7 @@ class SourceGrouping : public SourceGroupingInterface { std::shared_ptr m_group_factory; std::list> m_source_groups; unsigned int m_hard_limit; + unsigned int m_total_sources_waiting = 0; }; /* End of SourceGrouping class */ diff --git a/SEFramework/src/lib/Pipeline/SourceGrouping.cpp b/SEFramework/src/lib/Pipeline/SourceGrouping.cpp index 8269f262a..417bdfd17 100644 --- a/SEFramework/src/lib/Pipeline/SourceGrouping.cpp +++ b/SEFramework/src/lib/Pipeline/SourceGrouping.cpp @@ -21,11 +21,15 @@ */ #include "SEFramework/Pipeline/SourceGrouping.h" +#include #include namespace SourceXtractor { +static Elements::Logging logger = Elements::Logging::getLogger("SourceGrouping"); + + SourceGrouping::SourceGrouping(std::shared_ptr grouping_criteria, std::shared_ptr group_factory, unsigned int hard_limit) @@ -33,6 +37,9 @@ SourceGrouping::SourceGrouping(std::shared_ptr grouping_criter } void SourceGrouping::receiveSource(std::unique_ptr source) { + m_total_sources_waiting++; + logger.debug() << "Receiving source, sources in grouping: " << m_total_sources_waiting; + // Pointer which points to the group of the source SourceGroupInterface* matched_group = nullptr; @@ -87,6 +94,8 @@ void SourceGrouping::receiveSource(std::unique_ptr source) { void SourceGrouping::receiveProcessSignal(const ProcessSourcesEvent& process_event) { std::vector>::iterator> groups_to_process; + logger.debug() << "Received processing signal, total sources waiting in grouping: " << m_total_sources_waiting; + // We iterate through all the SourceGroups we have for (auto group_it = m_source_groups.begin(); group_it != m_source_groups.end(); ++group_it) { // We look at its Sources and if we find at least one that needs to be processed we put it in groups_to_process @@ -101,9 +110,13 @@ void SourceGrouping::receiveProcessSignal(const ProcessSourcesEvent& process_eve // For each SourceGroup that we put in groups_to_process, for (auto& group : groups_to_process) { // we remove it from our list of stored SourceGroups and notify our observers + m_total_sources_waiting -= (*group)->size(); + logger.debug() << "Sending group size " << (*group)->size() << ", sources remaining in grouping: " << m_total_sources_waiting; sendSource(std::move(*group)); m_source_groups.erase(group); } + + logger.debug() << "Processing signal handled, total sources remaining in grouping: " << m_total_sources_waiting; } std::set SourceGrouping::requiredProperties() const { diff --git a/SEFramework/src/lib/Source/EntangledSource.cpp b/SEFramework/src/lib/Source/EntangledSource.cpp index 32db767e0..9e0a8a60e 100644 --- a/SEFramework/src/lib/Source/EntangledSource.cpp +++ b/SEFramework/src/lib/Source/EntangledSource.cpp @@ -21,9 +21,12 @@ #include "SEFramework/Source/SourceGroupWithOnDemandProperties.h" #include "SEFramework/Task/GroupTask.h" +#include "ElementsKernel/Logging.h" namespace SourceXtractor { +static Elements::Logging logger = Elements::Logging::getLogger("EntangledSource"); + SourceGroupWithOnDemandProperties::EntangledSource::EntangledSource(std::shared_ptr source, SourceGroupWithOnDemandProperties& group) : m_source(source), m_group(group) { // Normally, it should not be possible that the given source is of type @@ -38,7 +41,6 @@ SourceGroupWithOnDemandProperties::EntangledSource::EntangledSource(std::shared_ } const Property& SourceGroupWithOnDemandProperties::EntangledSource::getProperty(const PropertyId& property_id) const { - // If we already have the property stored in this object, returns it if (m_property_holder.isPropertySet(property_id)) { return m_property_holder.getProperty(property_id); @@ -63,7 +65,8 @@ const Property& SourceGroupWithOnDemandProperties::EntangledSource::getProperty( throw PropertyNotFoundException(property_id); } - // Use the task to make the property + //logger.debug() << "Computing property " << property_id.getString(); + // Use the task to make the property group_task->computeProperties(m_group); // The property should now be available either in this object or in the group object diff --git a/SEFramework/src/lib/Source/SourceWithOnDemandProperties.cpp b/SEFramework/src/lib/Source/SourceWithOnDemandProperties.cpp index 3ca7cc349..8dec32dce 100644 --- a/SEFramework/src/lib/Source/SourceWithOnDemandProperties.cpp +++ b/SEFramework/src/lib/Source/SourceWithOnDemandProperties.cpp @@ -46,6 +46,7 @@ const Property& SourceWithOnDemandProperties::getProperty(const PropertyId& prop // if not, get the task that makes it and execute, we should have it then auto task = m_task_provider->getTask(property_id); if (task) { + //logger.debug() << "Computing property " << property_id.getString(); task->computeProperties(const_cast(*this)); return m_property_holder.getProperty(property_id); } diff --git a/SEImplementation/SEImplementation/Grouping/AssocGrouping.h b/SEImplementation/SEImplementation/Grouping/AssocGrouping.h index 7e5e67a79..585c3dcc3 100644 --- a/SEImplementation/SEImplementation/Grouping/AssocGrouping.h +++ b/SEImplementation/SEImplementation/Grouping/AssocGrouping.h @@ -42,6 +42,7 @@ class AssocGrouping : public SourceGroupingInterface { std::shared_ptr m_group_factory; std::map> m_source_groups; unsigned int m_hard_limit; + unsigned int m_total_sources_waiting = 0; }; diff --git a/SEImplementation/SEImplementation/Grouping/MoffatGrouping.h b/SEImplementation/SEImplementation/Grouping/MoffatGrouping.h index 6220b157a..80db1cd50 100644 --- a/SEImplementation/SEImplementation/Grouping/MoffatGrouping.h +++ b/SEImplementation/SEImplementation/Grouping/MoffatGrouping.h @@ -60,6 +60,8 @@ class MoffatGrouping : public SourceGroupingInterface { size_t m_group_counter; std::map> m_groups; QuadTree> m_tree; + + unsigned int m_total_sources_waiting = 0; }; } diff --git a/SEImplementation/SEImplementation/Grouping/SplitSourcesGrouping.h b/SEImplementation/SEImplementation/Grouping/SplitSourcesGrouping.h index 81db92714..cba705852 100644 --- a/SEImplementation/SEImplementation/Grouping/SplitSourcesGrouping.h +++ b/SEImplementation/SEImplementation/Grouping/SplitSourcesGrouping.h @@ -42,7 +42,7 @@ class SplitSourcesGrouping : public SourceGroupingInterface { std::shared_ptr m_group_factory; std::map> m_source_groups; unsigned int m_hard_limit; - + unsigned int m_total_sources_waiting = 0; }; } diff --git a/SEImplementation/src/lib/Grouping/AssocGrouping.cpp b/SEImplementation/src/lib/Grouping/AssocGrouping.cpp index ce80b812c..283b087ac 100644 --- a/SEImplementation/src/lib/Grouping/AssocGrouping.cpp +++ b/SEImplementation/src/lib/Grouping/AssocGrouping.cpp @@ -17,9 +17,11 @@ #include "SEImplementation/Grouping/AssocGrouping.h" #include "SEImplementation/Plugin/AssocMode/AssocMode.h" - +#include namespace SourceXtractor { +static Elements::Logging logger = Elements::Logging::getLogger("AssocGrouping"); + AssocGrouping::AssocGrouping(std::shared_ptr group_factory, unsigned int hard_limit) : m_group_factory(group_factory), m_hard_limit(hard_limit) { @@ -31,6 +33,9 @@ std::set AssocGrouping::requiredProperties() const { /// Handles a new Source void AssocGrouping::receiveSource(std::unique_ptr source) { + m_total_sources_waiting++; + logger.debug() << "Receiving source, sources in grouping: " << m_total_sources_waiting; + auto source_id = source->getProperty().getGroupId(); if (m_source_groups.find(source_id) == m_source_groups.end()) { @@ -41,6 +46,9 @@ void AssocGrouping::receiveSource(std::unique_ptr source) { if (m_hard_limit > 0 && m_source_groups.at(source_id)->size() >= m_hard_limit) { // the stored group has reached the hard limit // send the current group to processing + m_total_sources_waiting -= m_source_groups.at(source_id)->size(); + logger.debug() << "Group " << source_id << " reached hard limit of " << m_hard_limit + << " sources, sending group to processing, sources remaining in grouping: " << m_total_sources_waiting; sendSource(std::move(m_source_groups.at(source_id))); // and replace it with a new empty one @@ -55,6 +63,8 @@ void AssocGrouping::receiveSource(std::unique_ptr source) { void AssocGrouping::receiveProcessSignal(const ProcessSourcesEvent& event) { std::vector groups_to_process; + logger.debug() << "Received processing signal, total sources waiting in grouping: " << m_total_sources_waiting; + // We iterate through all the SourceGroups we have for (auto const& it : m_source_groups) { // We look at its Sources and if we find at least one that needs to be processed @@ -69,9 +79,14 @@ void AssocGrouping::receiveProcessSignal(const ProcessSourcesEvent& event) { // For each SourceGroup that we put in groups_to_process, for (auto group_id : groups_to_process) { // we remove it from our list of stored SourceGroups and notify our observers + m_total_sources_waiting -= m_source_groups[group_id]->size(); + logger.debug() << "Sending group size " << m_source_groups[group_id]->size() << ", sources remaining in grouping: " << m_total_sources_waiting; + sendSource(std::move(m_source_groups[group_id])); m_source_groups.erase(group_id); } + + logger.debug() << "Processing signal handled, total sources remaining in grouping: " << m_total_sources_waiting; } } // SourceXtractor namespace diff --git a/SEImplementation/src/lib/Grouping/MoffatGrouping.cpp b/SEImplementation/src/lib/Grouping/MoffatGrouping.cpp index 79e6aad70..2e946bf38 100644 --- a/SEImplementation/src/lib/Grouping/MoffatGrouping.cpp +++ b/SEImplementation/src/lib/Grouping/MoffatGrouping.cpp @@ -27,8 +27,12 @@ #include "SEImplementation/Plugin/PixelCentroid/PixelCentroid.h" #include "SEImplementation/Plugin/PeakValue/PeakValue.h" +#include + namespace SourceXtractor { +static Elements::Logging logger = Elements::Logging::getLogger("MoffatGrouping"); + template <> struct QuadTreeTraits> { static double getCoord(std::shared_ptr t, size_t index) { @@ -58,6 +62,8 @@ std::set MoffatGrouping::requiredProperties() const { /// Handles a new Source void MoffatGrouping::receiveSource(std::unique_ptr source) { + m_total_sources_waiting++; + logger.debug() << "Receiving source, sources in grouping: " << m_total_sources_waiting; // Encapsulates the source unique_ptr auto& centroid = source->getProperty(); @@ -110,6 +116,8 @@ void MoffatGrouping::receiveSource(std::unique_ptr source) { void MoffatGrouping::receiveProcessSignal(const ProcessSourcesEvent& event) { std::vector groups_to_process; + logger.debug() << "Received processing signal, total sources waiting in grouping: " << m_total_sources_waiting; + // We iterate through all the SourceGroups we have for (auto const& it : m_groups) { // We look at its Sources and if we find at least one that needs to be processed @@ -125,6 +133,8 @@ void MoffatGrouping::receiveProcessSignal(const ProcessSourcesEvent& event) { for (auto group_id : groups_to_process) { processGroup(group_id); } + + logger.debug() << "Processing signal handled, total sources remaining in grouping: " << m_total_sources_waiting; } void MoffatGrouping::processGroup(unsigned int group_id) { @@ -136,6 +146,9 @@ void MoffatGrouping::processGroup(unsigned int group_id) { m_tree.remove(source_info); } + m_total_sources_waiting -= new_group->size(); + logger.debug() << "Sending group size " << new_group->size() << ", sources remaining in grouping: " << m_total_sources_waiting; + sendSource(std::move(new_group)); m_groups.erase(group_id); } diff --git a/SEImplementation/src/lib/Grouping/SplitSourcesGrouping.cpp b/SEImplementation/src/lib/Grouping/SplitSourcesGrouping.cpp index 5cd807296..a14d8c768 100644 --- a/SEImplementation/src/lib/Grouping/SplitSourcesGrouping.cpp +++ b/SEImplementation/src/lib/Grouping/SplitSourcesGrouping.cpp @@ -17,9 +17,11 @@ #include "SEImplementation/Grouping/SplitSourcesGrouping.h" #include "SEImplementation/Property/SourceId.h" - +#include namespace SourceXtractor { +static Elements::Logging logger = Elements::Logging::getLogger("SplitSourcesGrouping"); + SplitSourcesGrouping::SplitSourcesGrouping(std::shared_ptr group_factory, unsigned int hard_limit) : m_group_factory(group_factory), m_hard_limit(hard_limit) { @@ -31,6 +33,9 @@ std::set SplitSourcesGrouping::requiredProperties() const { /// Handles a new Source void SplitSourcesGrouping::receiveSource(std::unique_ptr source) { + m_total_sources_waiting++; + logger.debug() << "Receiving source, sources in grouping: " << m_total_sources_waiting; + auto source_id = source->getProperty().getDetectionId(); if (m_source_groups.find(source_id) == m_source_groups.end()) { @@ -39,14 +44,17 @@ void SplitSourcesGrouping::receiveSource(std::unique_ptr source } if (m_hard_limit > 0 && m_source_groups.at(source_id)->size() >= m_hard_limit) { - // the stored group has reached the hard limit - // send the current group to processing - sendSource(std::move(m_source_groups.at(source_id))); + // the stored group has reached the hard limit + // send the current group to processing + m_total_sources_waiting -= m_source_groups.at(source_id)->size(); + logger.debug() << "Group " << source_id << " reached hard limit of " << m_hard_limit + << " sources, sending group to processing, sources remaining in grouping: " << m_total_sources_waiting; + sendSource(std::move(m_source_groups.at(source_id))); - // and replace it with a new empty one - auto new_group = m_group_factory->createSourceGroup(); - m_source_groups[source_id] = std::move(new_group); - } + // and replace it with a new empty one + auto new_group = m_group_factory->createSourceGroup(); + m_source_groups[source_id] = std::move(new_group); + } m_source_groups.at(source_id)->addSource(std::move(source)); } @@ -55,6 +63,8 @@ void SplitSourcesGrouping::receiveSource(std::unique_ptr source void SplitSourcesGrouping::receiveProcessSignal(const ProcessSourcesEvent& event) { std::vector groups_to_process; + logger.debug() << "Received processing signal, total sources waiting in grouping: " << m_total_sources_waiting; + // We iterate through all the SourceGroups we have for (auto const& it : m_source_groups) { // We look at its Sources and if we find at least one that needs to be processed @@ -69,9 +79,14 @@ void SplitSourcesGrouping::receiveProcessSignal(const ProcessSourcesEvent& event // For each SourceGroup that we put in groups_to_process, for (auto group_id : groups_to_process) { // we remove it from our list of stored SourceGroups and notify our observers + m_total_sources_waiting -= m_source_groups[group_id]->size(); + logger.debug() << "Sending group size " << m_source_groups[group_id]->size() << ", sources remaining in grouping: " << m_total_sources_waiting; + sendSource(std::move(m_source_groups[group_id])); m_source_groups.erase(group_id); } + + logger.debug() << "Processing signal handled, total sources remaining in grouping: " << m_total_sources_waiting; } } // SourceXtractor namespace diff --git a/SEImplementation/src/lib/Measurement/MultithreadedMeasurement.cpp b/SEImplementation/src/lib/Measurement/MultithreadedMeasurement.cpp index b29fdecff..64fc5b3af 100644 --- a/SEImplementation/src/lib/Measurement/MultithreadedMeasurement.cpp +++ b/SEImplementation/src/lib/Measurement/MultithreadedMeasurement.cpp @@ -30,7 +30,7 @@ using namespace SourceXtractor; -static Elements::Logging logger = Elements::Logging::getLogger("Multithreading"); +static Elements::Logging logger = Elements::Logging::getLogger("MultithreadedMeasurement"); MultithreadedMeasurement::~MultithreadedMeasurement() { @@ -44,6 +44,7 @@ void MultithreadedMeasurement::startThreads() { } void MultithreadedMeasurement::stopThreads() { + logger.debug() << "Stopping worker threads"; m_input_done = true; m_thread_pool->block(); m_output_thread->join(); @@ -51,6 +52,8 @@ void MultithreadedMeasurement::stopThreads() { } void MultithreadedMeasurement::synchronizeThreads() { + logger.debug() << "Synchronizing worker threads"; + // Wait until all worker threads are done m_thread_pool->block(); @@ -74,6 +77,9 @@ void MultithreadedMeasurement::synchronizeThreads() { } void MultithreadedMeasurement::receiveSource(std::unique_ptr source_group) { + logger.debug() << "Receiving source group with " << source_group->size() << " sources, job queue size: " + << m_thread_pool->queued() << " output queue size: " << m_output_queue.size(); + // Force computation of SourceID here, where the order is still deterministic for (auto& source : *source_group) { source.getProperty(); @@ -127,8 +133,14 @@ void MultithreadedMeasurement::outputThreadLoop() { // Process the output queue while (!m_output_queue.empty()) { + logger.debug() << "Output queue: " << m_output_queue.size() << " Sending group to output, group has " + << m_output_queue.front().second->size() << " sources"; + sendSource(std::move(m_output_queue.front().second)); m_output_queue.pop_front(); + + logger.debug() << "Thread pool: active threads: " << m_thread_pool->activeThreads() << ", running jobs: " << m_thread_pool->running() + << " queued jobs: " << m_thread_pool->queued(); } if (m_input_done && m_thread_pool->running() + m_thread_pool->queued() == 0 && diff --git a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp index 1e99d9eae..08ecaf19f 100644 --- a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp +++ b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp @@ -376,6 +376,9 @@ std::shared_ptr> FlexibleModelFittingIterativeTask::createW } void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& group) const { + + auto start_time = std::chrono::high_resolution_clock::now(); + FittingState fitting_state; for (auto& source : group) { @@ -478,6 +481,10 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& } updateCheckImages(group, 1.0, fitting_state); + + auto end_time = std::chrono::high_resolution_clock::now(); + auto duration = std::chrono::duration_cast(end_time - start_time).count(); + logger.debug() << "Flexible model fitting completed in " << duration << " ms, group size: " << group.size(); } diff --git a/SEImplementation/src/lib/Prefetcher/Prefetcher.cpp b/SEImplementation/src/lib/Prefetcher/Prefetcher.cpp index 42f6325ff..0eaf239f1 100644 --- a/SEImplementation/src/lib/Prefetcher/Prefetcher.cpp +++ b/SEImplementation/src/lib/Prefetcher/Prefetcher.cpp @@ -144,16 +144,20 @@ void Prefetcher::receiveProcessSignal(const ProcessSourcesEvent& message) { } void Prefetcher::wait() { + logger.debug() << "Waiting for prefetcher to finish"; m_stop = true; m_output_thread->join(); } void Prefetcher::synchronize() { + logger.debug() << "Waiting for prefetcher queue to be empty"; + // Wait until the output queue is empty while (true) { { std::unique_lock output_lock(m_queue_mutex); if (m_received.empty()) { + logger.debug() << "Prefetcher queue is empty"; break; } else if (m_thread_pool->checkForException(false)) { @@ -164,6 +168,7 @@ void Prefetcher::synchronize() { throw Elements::Exception() << "No active threads and the queue is not empty! Please, report this as a bug"; } } + logger.debug() << "Queue: " << m_received.size() << " activeThreads: " << m_thread_pool->activeThreads(); std::this_thread::sleep_for(std::chrono::milliseconds(100)); } } From e9ae7d234dc485f37c911a89280eff09d1617572 Mon Sep 17 00:00:00 2001 From: Marc Schefer Date: Tue, 25 Nov 2025 16:18:58 +0100 Subject: [PATCH 2/4] more detailed profiling in model fitting and possible optimizations --- .../FlexibleModelFittingIterativeTask.h | 8 +- .../FlexibleModelFittingIterativeTask.cpp | 184 +++++++++++++++--- 2 files changed, 165 insertions(+), 27 deletions(-) diff --git a/SEImplementation/SEImplementation/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.h b/SEImplementation/SEImplementation/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.h index 9c0ecc3cc..2d56a17e8 100644 --- a/SEImplementation/SEImplementation/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.h +++ b/SEImplementation/SEImplementation/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.h @@ -85,6 +85,10 @@ class FlexibleModelFittingIterativeTask : public GroupTask { std::vector iterations_per_meta; std::vector fitting_areas_x; std::vector fitting_areas_y; + + // Cached images + std::vector>> cached_image_copies; + std::vector>> cached_weight_images; }; struct FittingState { @@ -111,8 +115,8 @@ class FlexibleModelFittingIterativeTask : public GroupTask { void updateCheckImages(SourceGroupInterface& group, double pixel_scale, FittingState& state) const; SeFloat computeChiSquared(SourceGroupInterface& group, SourceInterface& source, int index, double pixel_scale, FlexibleModelFittingParameterManager& manager, int& total_data_points, FittingState& state) const; - SeFloat computeChiSquaredForFrame(std::shared_ptr> image, - std::shared_ptr> model, std::shared_ptr> weights, int& data_points) const; + SeFloat computeChiSquaredForFrame(std::shared_ptr> image, + std::shared_ptr> model, std::shared_ptr> weights, int& data_points) const; int fitSourcePrepareParameters(FlexibleModelFittingParameterManager& parameter_manager, ModelFitting::EngineParameterManager& engine_parameter_manager, SourceInterface& source, int index, FittingState& state) const; diff --git a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp index 08ecaf19f..c79dd58bf 100644 --- a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp +++ b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp @@ -416,13 +416,29 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& } } + // pre compute images copies and weight images for each frame + for (auto frame : m_frames) { + int frame_index = frame->getFrameNb(); + if (isFrameValid(source, frame_index)) { + initial_state.cached_image_copies.push_back(createImageCopy(source, frame_index)); + initial_state.cached_weight_images.push_back(createWeightImage(source, frame_index)); + } else { + initial_state.cached_image_copies.push_back(nullptr); + initial_state.cached_weight_images.push_back(nullptr); + } + } + fitting_state.source_states.emplace_back(std::move(initial_state)); } + + // TODO Sort sources by flux to fit brightest sources first? // iterate over the whole group, fitting sources one at a time + auto time_before_iterations = std::chrono::high_resolution_clock::now(); + double prev_chi_squared = 999999.9; for (int iteration = 0; iteration < m_meta_iterations; iteration++) { int index = 0; @@ -440,12 +456,14 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& chi_squared /= fitting_state.source_states.size(); if (fabs(chi_squared - prev_chi_squared) / chi_squared < m_meta_iteration_stop) { - break; + break; } prev_chi_squared = chi_squared; } + auto time_after_iterations = std::chrono::high_resolution_clock::now(); + // Remove parameters that couldn't be fit from the output @@ -462,6 +480,8 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& } } + auto time_finalization_start = std::chrono::high_resolution_clock::now(); + // output a property for every source size_t index = 0; for (auto& source : group) { @@ -480,11 +500,22 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& index++; } - updateCheckImages(group, 1.0, fitting_state); + auto time_check_images = std::chrono::high_resolution_clock::now(); + //updateCheckImages(group, 1.0, fitting_state); auto end_time = std::chrono::high_resolution_clock::now(); - auto duration = std::chrono::duration_cast(end_time - start_time).count(); - logger.debug() << "Flexible model fitting completed in " << duration << " ms, group size: " << group.size(); + + auto total_duration = std::chrono::duration_cast(end_time - start_time).count(); + auto before_iterations_duration = std::chrono::duration_cast(time_before_iterations - start_time).count(); + auto iteration_duration = std::chrono::duration_cast(time_after_iterations - time_before_iterations).count(); + auto after_iterations_duration = std::chrono::duration_cast(time_finalization_start - time_after_iterations).count(); + auto finalization_duration = std::chrono::duration_cast(time_check_images - time_finalization_start).count(); + auto check_images_duration = std::chrono::duration_cast(end_time - time_check_images).count(); + + logger.debug() << "Flexible model fitting completed in " << total_duration << " ms, group size: " << group.size() + << " (preparation: " << before_iterations_duration << " ms, iterations: " << iteration_duration << " ms, post-processing: " + << after_iterations_duration << " ms, finalization: " << finalization_duration << " ms, check images: " + << check_images_duration << " ms)"; } @@ -574,8 +605,16 @@ int FlexibleModelFittingIterativeTask::fitSourcePrepareModels(FlexibleModelFitti ResidualEstimator& res_estimator, int& good_pixels, SourceGroupInterface& group, SourceInterface& source, int index, FittingState& state, double down_scaling) const { + double pixel_scale = 1.0; + double frame_duration = 0.0; + double image_duration = 0.0; + double deblend_duration = 0.0; + double weight_duration = 0.0; + double residuals_duration = 0.0; + + int valid_frames = 0; for (auto frame : m_frames) { int frame_index = frame->getFrameNb(); @@ -583,11 +622,19 @@ int FlexibleModelFittingIterativeTask::fitSourcePrepareModels(FlexibleModelFitti if (isFrameValid(source, frame_index)) { valid_frames++; + auto start_time = std::chrono::high_resolution_clock::now(); + auto stamp_rect = getFittingRect(source, frame_index); auto frame_model = createFrameModel(source, pixel_scale, parameter_manager, frame, stamp_rect, down_scaling); - auto image = createImageCopy(source, frame_index); + auto frame_time = std::chrono::high_resolution_clock::now(); + frame_duration += std::chrono::duration_cast(frame_time - start_time).count(); + + auto image = VectorImage::create(*state.source_states[index].cached_image_copies[frame_index]); + auto image_time = std::chrono::high_resolution_clock::now(); + image_duration += std::chrono::duration_cast(image_time - frame_time).count(); + auto deblend_image = createDeblendImage(group, source, index, frame, state); for (int y = 0; y < image->getHeight(); ++y) { for (int x = 0; x < image->getWidth(); ++x) { @@ -595,23 +642,39 @@ int FlexibleModelFittingIterativeTask::fitSourcePrepareModels(FlexibleModelFitti } } - auto weight = createWeightImage(source, frame_index); + auto deblend_time = std::chrono::high_resolution_clock::now(); + deblend_duration += std::chrono::duration_cast(deblend_time - image_time).count(); + auto weight_image = state.source_states[index].cached_weight_images[frame_index]; + + // count number of pixels that can be used for fitting - for (int y = 0; y < weight->getHeight(); ++y) { - for (int x = 0; x < weight->getWidth(); ++x) { - good_pixels += (weight->at(x, y) != 0.); + for (int y = 0; y < weight_image->getHeight(); ++y) { + for (int x = 0; x < weight_image->getWidth(); ++x) { + good_pixels += (weight_image->at(x, y) != 0.); } } - + + auto weight_time = std::chrono::high_resolution_clock::now(); + weight_duration += std::chrono::duration_cast(weight_time - deblend_time).count(); // Setup residuals + auto data_vs_model = - createDataVsModelResiduals(image, std::move(frame_model), weight, - //LogChiSquareComparator(m_modified_chi_squared_scale)); + createDataVsModelResiduals(image, std::move(frame_model), weight_image, AsinhChiSquareComparator(m_modified_chi_squared_scale)); res_estimator.registerBlockProvider(std::move(data_vs_model)); + + auto end_time = std::chrono::high_resolution_clock::now(); + residuals_duration += std::chrono::duration_cast(end_time - weight_time).count(); + } } + logger.debug() << "Frame model preparation times (ms): " + << " model: " << frame_duration + << ", image copy: " << image_duration + << ", deblend: " << deblend_duration + << ", weight: " << weight_duration + << ", setup residuals: " << residuals_duration; return valid_frames; } @@ -686,6 +749,8 @@ void FlexibleModelFittingIterativeTask::fitSourceUpdateState( void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, SourceInterface& source, int index, FittingState& state) const { + auto start_time = std::chrono::high_resolution_clock::now(); + ////////////////////////////////////////////// // Determine size of fitted area and if needed downsize factor @@ -708,6 +773,8 @@ void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, S << " scaling factor: " << down_scaling; } + auto init_time = std::chrono::high_resolution_clock::now(); + ////////////////////////////////////////////// // Prepare parameters @@ -716,6 +783,8 @@ void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, S int n_free_parameters = fitSourcePrepareParameters( parameter_manager, engine_parameter_manager, source, index, state); + auto prep_param_time = std::chrono::high_resolution_clock::now(); + /////////////////////////////////////////////////////////////////////////////////// // Add models for all frames ResidualEstimator res_estimator {}; @@ -723,6 +792,8 @@ void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, S int valid_frames = fitSourcePrepareModels( parameter_manager, res_estimator, n_good_pixels, group, source, index, state, down_scaling); + auto prep_model_time = std::chrono::high_resolution_clock::now(); + /////////////////////////////////////////////////////////////////////////////// // Check that we had enough data for the fit @@ -751,6 +822,8 @@ void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, S prior->setupPrior(parameter_manager, source, res_estimator); } + auto prep_prior_time = std::chrono::high_resolution_clock::now(); + ///////////////////////////////////////////////////////////////////////////////// // Model fitting @@ -764,15 +837,32 @@ void FlexibleModelFittingIterativeTask::fitSource(SourceGroupInterface& group, S } auto duration = solution.duration; + auto fit_time = std::chrono::high_resolution_clock::now(); + //////////////////////////////////////////////////////////////////////////////////// // compute chi squared - SeFloat avg_reduced_chi_squared = fitSourceComputeChiSquared(parameter_manager, group, source, index, state); + SeFloat avg_reduced_chi_squared = fitSourceComputeChiSquared(parameter_manager, group, source, index, state); + + auto chi_squared_time = std::chrono::high_resolution_clock::now(); + //////////////////////////////////////////////////////////////////////////////////// // update state with results fitSourceUpdateState(parameter_manager, source, avg_reduced_chi_squared, duration, iterations, stop_reason, flags, solution, index, state); + + auto end_time = std::chrono::high_resolution_clock::now(); + + logger.debug() << " Flexible model fitting for source " + << " completed in " << std::chrono::duration_cast(end_time - start_time).count() << " ms (init: " + << std::chrono::duration_cast(init_time - start_time).count() << " ms, param prep: " + << std::chrono::duration_cast(prep_param_time - init_time).count() << " ms, model prep: " + << std::chrono::duration_cast(prep_model_time - prep_param_time).count() << " ms, prior prep: " + << std::chrono::duration_cast(prep_prior_time - prep_model_time).count() << " ms, fit: " + << std::chrono::duration_cast(fit_time - prep_prior_time).count() << " ms, chi_squared: " + << std::chrono::duration_cast(chi_squared_time - fit_time).count() << " ms, update state: " + << std::chrono::duration_cast(end_time - chi_squared_time).count() << " ms)"; } void FlexibleModelFittingIterativeTask::updateCheckImages(SourceGroupInterface& group, @@ -802,6 +892,7 @@ void FlexibleModelFittingIterativeTask::updateCheckImages(SourceGroupInterface& index++; } + index = 0; for (auto& src : group) { for (auto frame : m_frames) { int frame_index = frame->getFrameNb(); @@ -812,7 +903,7 @@ void FlexibleModelFittingIterativeTask::updateCheckImages(SourceGroupInterface& auto frame_model = createFrameModel(src, pixel_scale, parameter_manager, frame, stamp_rect); auto final_stamp = frame_model.getImage(); - auto weight_image = createWeightImage(src, frame_index); + auto weight_image = state.source_states[index].cached_weight_images[frame_index]; { auto debug_image = CheckImages::getInstance().getModelFittingImage(frame_index); @@ -849,22 +940,20 @@ void FlexibleModelFittingIterativeTask::updateCheckImages(SourceGroupInterface& } } + index++; } } -SeFloat FlexibleModelFittingIterativeTask::computeChiSquaredForFrame(std::shared_ptr> image, - std::shared_ptr> model, std::shared_ptr> weights, int& data_points) const { +SeFloat FlexibleModelFittingIterativeTask::computeChiSquaredForFrame(std::shared_ptr> image, + std::shared_ptr> model, std::shared_ptr> weights, int& data_points) const { double reduced_chi_squared = 0.0; data_points = 0; - ImageAccessor imageAccessor(image), modelAccessor(model); - ImageAccessor weightAccessor(weights); - for (int y=0; y < image->getHeight(); y++) { for (int x=0; x < image->getWidth(); x++) { - double tmp = imageAccessor.getValue(x, y) - modelAccessor.getValue(x, y); - reduced_chi_squared += tmp * tmp * weightAccessor.getValue(x, y) * weightAccessor.getValue(x, y); - if (weightAccessor.getValue(x, y) > 0) { + double tmp = image->at(x, y) - model->at(x, y); + reduced_chi_squared += tmp * tmp * weights->at(x, y) * weights->at(x, y); + if (weights->at(x, y) > 0) { data_points++; } } @@ -874,6 +963,15 @@ SeFloat FlexibleModelFittingIterativeTask::computeChiSquaredForFrame(std::shared SeFloat FlexibleModelFittingIterativeTask::computeChiSquared(SourceGroupInterface& group, SourceInterface& source, int index, double pixel_scale, FlexibleModelFittingParameterManager& manager, int& total_data_points, FittingState& state) const { + + double frame_model_time = 0.0; + double image_time = 0.0; + double image_copy_time = 0.0; + double deblend_time = 0.0; + double deblend_image_time = 0.0; + double weight_time = 0.0; + double chi_squared_time = 0.0; + SeFloat total_chi_squared = 0; total_data_points = 0; int valid_frames = 0; @@ -881,29 +979,65 @@ SeFloat FlexibleModelFittingIterativeTask::computeChiSquared(SourceGroupInterfac int frame_index = frame->getFrameNb(); // Validate that each frame covers the model fitting region if (isFrameValid(source, frame_index)) { + auto start_time = std::chrono::high_resolution_clock::now(); + valid_frames++; auto stamp_rect = getFittingRect(source, frame_index); auto frame_model = createFrameModel(source, pixel_scale, manager, frame, stamp_rect); + + auto create_time = std::chrono::high_resolution_clock::now(); + frame_model_time += std::chrono::duration_cast(create_time - start_time).count(); + auto final_stamp = frame_model.getImage(); - auto image = createImageCopy(source, frame_index); + auto image_time_point = std::chrono::high_resolution_clock::now(); + image_time += std::chrono::duration_cast(image_time_point - create_time).count(); + + auto image = VectorImage::create(*state.source_states[index].cached_image_copies[frame_index]); + + auto image_copy_time_point = std::chrono::high_resolution_clock::now(); + image_copy_time += std::chrono::duration_cast(image_copy_time_point - image_time_point).count(); + + auto deblend_image = createDeblendImage(group, source, index, frame, state); + + auto deblend_time_point = std::chrono::high_resolution_clock::now(); + deblend_time += std::chrono::duration_cast(deblend_time_point - image_copy_time_point).count(); + for (int y = 0; y < image->getHeight(); ++y) { for (int x = 0; x < image->getWidth(); ++x) { image->at(x, y) -= deblend_image->at(x, y); } } - auto weight = createWeightImage(source, frame_index); + auto deblend_image_time_point = std::chrono::high_resolution_clock::now(); + deblend_image_time += std::chrono::duration_cast(deblend_image_time_point - deblend_time_point).count(); + + auto weight_image = state.source_states[index].cached_weight_images[frame_index]; + + auto weight_time_point = std::chrono::high_resolution_clock::now(); + weight_time += std::chrono::duration_cast(weight_time_point - deblend_image_time_point).count(); int data_points = 0; - SeFloat chi_squared = computeChiSquaredForFrame(image, final_stamp, weight, data_points); + SeFloat chi_squared = computeChiSquaredForFrame(image, final_stamp, weight_image, data_points); + + auto chi_squared_time_point = std::chrono::high_resolution_clock::now(); + chi_squared_time += std::chrono::duration_cast(chi_squared_time_point - weight_time_point).count(); total_data_points += data_points; total_chi_squared += chi_squared; } } + logger.debug() << "Frame model chi squared computation times (ms): " + << " model: " << frame_model_time + << ", image copy: " << image_copy_time + << ", image: " << image_time + << ", deblend: " << deblend_time + << ", deblend image: " << deblend_image_time + << ", weight: " << weight_time + << ", chi_squared: " << chi_squared_time; + return total_chi_squared; } From fe030b651702c1548672557dbe114a91d2e4ea6c Mon Sep 17 00:00:00 2001 From: Marc Schefer Date: Wed, 3 Dec 2025 13:03:59 +0100 Subject: [PATCH 3/4] re-enables OpenCV FFT --- SEBenchmarks/src/program/BenchConvolution.cpp | 2 +- SEFramework/CMakeLists.txt | 5 ++-- .../Convolution/OpenCVConvolution.h | 29 ++++++++++++++---- SEFramework/src/lib/FFT/FFT.cpp | 12 ++++++++ SEImplementation/CMakeLists.txt | 3 +- .../Image/DownSampledImagePsf.h | 11 +++++-- .../SEImplementation/Image/ImagePsf.h | 28 ++++++++++++++++- .../src/lib/Image/DownSampledImagePsf.cpp | 3 ++ .../FlexibleModelFittingIterativeTask.cpp | 30 +++++++++++++++++++ 9 files changed, 111 insertions(+), 12 deletions(-) diff --git a/SEBenchmarks/src/program/BenchConvolution.cpp b/SEBenchmarks/src/program/BenchConvolution.cpp index 60562272c..ead3ece7b 100644 --- a/SEBenchmarks/src/program/BenchConvolution.cpp +++ b/SEBenchmarks/src/program/BenchConvolution.cpp @@ -106,7 +106,7 @@ class BenchConvolution : public Elements::Program { #ifdef WITH_OPENCV logger.info() << "Timing OpenCV implementation"; - auto opencv_result = benchmark(image, kernel, repeat, measures); + auto opencv_result = benchmark>(image, kernel, repeat, measures); #endif if (krn_size <= 10 || img_size <= 20) { diff --git a/SEFramework/CMakeLists.txt b/SEFramework/CMakeLists.txt index 775552542..ccaa175c4 100644 --- a/SEFramework/CMakeLists.txt +++ b/SEFramework/CMakeLists.txt @@ -33,6 +33,7 @@ find_package(CCfits) find_package(BoostDLL) find_package(FFTW COMPONENTS single double) find_package(WCSLIB REQUIRED) +find_package(OpenCV REQUIRED) #=============================================================================== @@ -59,8 +60,8 @@ elements_add_library(SEFramework src/lib/CoordinateSystem/*.cpp LINK_LIBRARIES ElementsKernel SEUtils Table Configuration MathUtils FilePool - WCSLIB ${CCFITS_LIBRARIES} ${FFTW_LIBRARIES} - INCLUDE_DIRS SEUtils ${CCFITS_INCLUDE_DIRS} ${BoostDLL_INCLUDE_DIRS} ${FFTW_INCLUDE_DIRS} + WCSLIB ${CCFITS_LIBRARIES} ${FFTW_LIBRARIES} ${OpenCV_LIBRARIES} + INCLUDE_DIRS SEUtils ${CCFITS_INCLUDE_DIRS} ${BoostDLL_INCLUDE_DIRS} ${FFTW_INCLUDE_DIRS} ${OpenCV_INCLUDE_DIRS} PUBLIC_HEADERS SEFramework) #=============================================================================== diff --git a/SEFramework/SEFramework/Convolution/OpenCVConvolution.h b/SEFramework/SEFramework/Convolution/OpenCVConvolution.h index 16f9fc943..e5b858656 100644 --- a/SEFramework/SEFramework/Convolution/OpenCVConvolution.h +++ b/SEFramework/SEFramework/Convolution/OpenCVConvolution.h @@ -28,10 +28,11 @@ namespace SourceXtractor { +template class OpenCVConvolution { public: - explicit OpenCVConvolution(std::shared_ptr> img) - : m_kernel(img->getWidth(), img->getHeight(), CV_32F) { + explicit OpenCVConvolution(std::shared_ptr> img) + : m_kernel(img->getWidth(), img->getHeight(), CV_32F) , m_original_kernel{VectorImage::create(*MirrorImage::create(img))} { cv::Mat aux(img->getWidth(), img->getHeight(), CV_32F); std::copy(img->getData().begin(), img->getData().end(), aux.begin()); cv::flip(aux, m_kernel, -1); @@ -39,13 +40,26 @@ class OpenCVConvolution { virtual ~OpenCVConvolution() = default; - void convolve(std::shared_ptr> image) const { + void convolve(std::shared_ptr> image) const { + auto vector_image = VectorImage::create(image); + cv::Mat image_cv (image->getHeight(), image->getWidth(), CV_32F); - std::copy(image->getData().begin(), image->getData().end(), image_cv.begin()); + // for (int y = 0; y < image->getHeight(); ++y) { + // for (int x = 0; x < image->getWidth(); ++x) { + // image_cv.at(y, x) = image->getValue(x, y); + // } + // } + std::copy(vector_image->getData().begin(), vector_image->getData().end(), image_cv.begin()); cv::filter2D(image_cv, image_cv, -1, m_kernel); - std::copy(image_cv.begin(), image_cv.end(), image->getData().begin()); + for (int y = 0; y < image->getHeight(); ++y) { + for (int x = 0; x < image->getWidth(); ++x) { + image->setValue(x, y, image_cv.at(y, x)); + } + } + + //std::copy(image_cv.begin(), image_cv.end(), image->getData().begin()); } std::size_t getWidth() const { @@ -56,8 +70,13 @@ class OpenCVConvolution { return m_kernel.size[1]; } + std::shared_ptr> getKernel() const { + return m_original_kernel; + } + private: cv::Mat m_kernel;// (size, size, CV_32F); + std::shared_ptr> m_original_kernel; }; } // end of SourceXtractor diff --git a/SEFramework/src/lib/FFT/FFT.cpp b/SEFramework/src/lib/FFT/FFT.cpp index da53fda05..79c78ca4d 100644 --- a/SEFramework/src/lib/FFT/FFT.cpp +++ b/SEFramework/src/lib/FFT/FFT.cpp @@ -26,6 +26,8 @@ #include #include +#include + namespace SourceXtractor { /** @@ -130,6 +132,9 @@ auto FFT::createForwardPlan(int width, int height, std::vector& inout) -> } // No available plan yet, so get one from FFTW + + std::cout << "!!!!!!!!!!Creating new FFTW forward plan for size " << width << "x" << height << std::endl; + boost::upgrade_to_unique_lock write_lock{read_lock}; boost::lock_guard lock_planner{fftw_global_plan_mutex}; @@ -162,7 +167,11 @@ auto FFT::createInversePlan(int width, int height, std::vector& inout) -> static boost::shared_mutex mutex; static std::map, plan_ptr_t> plan_cache; + + auto before_time = std::chrono::high_resolution_clock::now(); boost::upgrade_lock read_lock{mutex}; + auto after_time = std::chrono::high_resolution_clock::now(); + std::cout << "!!!!Time to get read lock: " << std::chrono::duration_cast(after_time - before_time).count() << " us" << std::endl; auto pi = plan_cache.find(std::make_tuple(width, height)); if (pi != plan_cache.end()) { @@ -170,6 +179,9 @@ auto FFT::createInversePlan(int width, int height, std::vector& inout) -> } // No available plan yet, so get one from FFTW + + std::cout << "!!!!!!!!!!Creating new FFTW inverse plan for size " << width << "x" << height << std::endl; + boost::upgrade_to_unique_lock write_lock{read_lock}; boost::lock_guard lock_planner{fftw_global_plan_mutex}; diff --git a/SEImplementation/CMakeLists.txt b/SEImplementation/CMakeLists.txt index 1b52a93d3..4dd316c54 100644 --- a/SEImplementation/CMakeLists.txt +++ b/SEImplementation/CMakeLists.txt @@ -37,6 +37,7 @@ find_package(BoostPython ${PYTHON_EXPLICIT_VERSION}) find_package(Boost 1.53 REQUIRED) find_package(OnnxRuntime) find_package(Log4CPP REQUIRED) +find_package(OpenCV REQUIRED) if (${Boost_VERSION} LESS "106700") message(WARNING " @@ -104,7 +105,7 @@ elements_add_library(SEImplementation BoostPython PythonLibs PythonInterp Pyston ${OPT_LIBRARIES} - INCLUDE_DIRS ${BoostPython_INCLUDE_DIRS} ${PYTHON_INCLUDE_DIRS} ${OPT_INCLUDES} + INCLUDE_DIRS ${BoostPython_INCLUDE_DIRS} ${PYTHON_INCLUDE_DIRS} ${OPT_INCLUDES} ${OpenCV_INCLUDE_DIRS} PUBLIC_HEADERS SEImplementation) #=============================================================================== diff --git a/SEImplementation/SEImplementation/Image/DownSampledImagePsf.h b/SEImplementation/SEImplementation/Image/DownSampledImagePsf.h index 79d082767..bea676337 100644 --- a/SEImplementation/SEImplementation/Image/DownSampledImagePsf.h +++ b/SEImplementation/SEImplementation/Image/DownSampledImagePsf.h @@ -52,9 +52,12 @@ class DownSampledImagePsf { std::shared_ptr> getScaledKernel(SeFloat scale) const; void convolve(std::shared_ptr> image) const; - // For ConvolutionContext optimization + +#ifdef USE_DFT_CONVOLUTION_FOR_PSF + // // For ConvolutionContext optimization std::unique_ptr::ConvolutionContext> prepare(const std::shared_ptr>& model_ptr) const; void convolve(std::shared_ptr> image, std::unique_ptr::ConvolutionContext>& context) const; +#endif private: double m_down_scaling; @@ -66,15 +69,19 @@ class DownSampledImagePsf { namespace ModelFitting { + /** * Specialization of PsfTraits, as DFTConvolution has the concept of context */ -template<> +#ifdef USE_DFT_CONVOLUTION_FOR_PSF + template<> struct PsfTraits { using context_t = typename std::unique_ptr; static constexpr bool has_context = true; }; +#endif + } // end of ModelFitting #endif /* _SEIMPLEMENTATION_IMAGE_DOWNSAMPLEDIMAGEPSF_H_ */ diff --git a/SEImplementation/SEImplementation/Image/ImagePsf.h b/SEImplementation/SEImplementation/Image/ImagePsf.h index 97902cdc5..8f79fbd80 100644 --- a/SEImplementation/SEImplementation/Image/ImagePsf.h +++ b/SEImplementation/SEImplementation/Image/ImagePsf.h @@ -31,12 +31,35 @@ #include "SEFramework/Image/ProcessedImage.h" #include "SEFramework/Convolution/Convolution.h" +//#define USE_DFT_CONVOLUTION_FOR_PSF +//#define USE_DIRECT_CONVOLUTION_FOR_PSF +#define USE_OPENCV_CONVOLUTION_FOR_PSF + +#ifdef USE_DIRECT_CONVOLUTION_FOR_PSF +#include "SEFramework/Convolution/DirectConvolution.h" +#endif + +#ifdef USE_OPENCV_CONVOLUTION_FOR_PSF +#include "SEFramework/Convolution/OpenCVConvolution.h" +#endif namespace SourceXtractor { +#ifdef USE_DIRECT_CONVOLUTION_FOR_PSF +class ImagePsf: public DirectConvolution> { +private: + typedef DirectConvolution> base_t; +#endif +#ifdef USE_DFT_CONVOLUTION_FOR_PSF class ImagePsf: public DFTConvolution> { private: - typedef DFTConvolution> base_t; +typedef DFTConvolution> base_t; +#endif +#ifdef USE_OPENCV_CONVOLUTION_FOR_PSF +class ImagePsf: public OpenCVConvolution { +private: + typedef OpenCVConvolution base_t; +#endif public: @@ -78,12 +101,15 @@ namespace ModelFitting { /** * Specialization of PsfTraits, as DFTConvolution has the concept of context */ +#ifdef USE_DFT_CONVOLUTION_FOR_PSF template<> struct PsfTraits { using context_t = typename std::unique_ptr; static constexpr bool has_context = true; }; +#endif + } // end of ModelFitting diff --git a/SEImplementation/src/lib/Image/DownSampledImagePsf.cpp b/SEImplementation/src/lib/Image/DownSampledImagePsf.cpp index 7656410b4..e375fb8af 100644 --- a/SEImplementation/src/lib/Image/DownSampledImagePsf.cpp +++ b/SEImplementation/src/lib/Image/DownSampledImagePsf.cpp @@ -105,6 +105,8 @@ void DownSampledImagePsf::convolve(std::shared_ptr> image) } } +#ifdef USE_DFT_CONVOLUTION_FOR_PSF + std::unique_ptr::ConvolutionContext> DownSampledImagePsf::prepare( const std::shared_ptr>& model_ptr) const { if (m_psf != nullptr) { @@ -121,6 +123,7 @@ void DownSampledImagePsf::convolve(std::shared_ptr> image, } } +#endif } diff --git a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp index c79dd58bf..30127461b 100644 --- a/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp +++ b/SEImplementation/src/lib/Plugin/FlexibleModelFitting/FlexibleModelFittingIterativeTask.cpp @@ -522,6 +522,9 @@ void FlexibleModelFittingIterativeTask::computeProperties(SourceGroupInterface& std::shared_ptr> FlexibleModelFittingIterativeTask::createDeblendImage( SourceGroupInterface& group, SourceInterface& source, int source_index, std::shared_ptr frame, FittingState& state) const { + + auto start_time = std::chrono::high_resolution_clock::now(); + int frame_index = frame->getFrameNb(); auto rect = getFittingRect(source, frame_index); @@ -553,22 +556,49 @@ std::shared_ptr> FlexibleModelFittingIterativeTask::createD index++; } + auto parameter_time = std::chrono::high_resolution_clock::now(); + + double frame_model_duration = 0.0; + double get_image_duration = 0.0; + double deblend_image_duration = 0.0; + auto deblend_image = VectorImage::create(rect.getWidth(), rect.getHeight()); index = 0; for (auto& src : group) { if (index != source_index) { + auto source_start_time = std::chrono::high_resolution_clock::now(); + auto frame_model = createFrameModel(src, pixel_scale, parameter_manager, frame, rect); + + auto frame_model_time = std::chrono::high_resolution_clock::now(); + frame_model_duration += std::chrono::duration_cast(frame_model_time - source_start_time).count(); + auto final_stamp = frame_model.getImage(); + auto get_image_time = std::chrono::high_resolution_clock::now(); + get_image_duration += std::chrono::duration_cast(get_image_time - frame_model_time).count(); + for (int y = 0; y < final_stamp->getHeight(); ++y) { for (int x = 0; x < final_stamp->getWidth(); ++x) { deblend_image->at(x, y) += final_stamp->at(x, y); } } + + auto deblend_image_time = std::chrono::high_resolution_clock::now(); + deblend_image_duration += std::chrono::duration_cast(deblend_image_time - get_image_time).count(); } index++; } + auto end_time = std::chrono::high_resolution_clock::now(); + + // logger.debug() << "Deblend image preparation times (ms): " + // << " parameter setup: " << std::chrono::duration_cast(parameter_time - start_time).count() + // << ", frame models: " << frame_model_duration + // << ", get images: " << get_image_duration + // << ", deblend image sum: " << deblend_image_duration + // << ", total: " << std::chrono::duration_cast(end_time - start_time).count(); + return deblend_image; } From 8eb42f1b7dca6d2f4c36eb72f8368b71e4b69843 Mon Sep 17 00:00:00 2001 From: Marc Schefer Date: Wed, 3 Dec 2025 13:07:28 +0100 Subject: [PATCH 4/4] add opencv-devel dependency for github CI --- .github/workflows/dependencies.txt | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/dependencies.txt b/.github/workflows/dependencies.txt index 3a1576bb4..7c92c8fef 100644 --- a/.github/workflows/dependencies.txt +++ b/.github/workflows/dependencies.txt @@ -17,3 +17,4 @@ levmar-devel log4cpp-devel onnxruntime-devel wcslib-devel +opencv-devel