ImagePreprocessor: Separate azint from preprocessing (later put azint on GPU!)

This commit is contained in:
2026-04-22 13:34:07 +02:00
parent d0446b453b
commit 9338a876fd
12 changed files with 59 additions and 45 deletions
+1
View File
@@ -87,6 +87,7 @@ struct DataMessage {
std::optional<float> processing_time_s;
std::optional<float> spot_finding_time_s;
std::optional<float> azint_time_s;
std::optional<float> indexing_time_s;
std::optional<float> refinement_time_s;
std::optional<float> bragg_prediction_time_s;
+1
View File
@@ -163,6 +163,7 @@ See [DECTRIS documentation](https://github.com/dectris/documentation/tree/main/s
| b_factor | float | Estimated B-factor (Angstrom^2) | | |
| compression_time | float | Time spent on compression/decompressing image \[s\] | | |
| preprocessing_time | float | Time spent on preparing the image for analysis \[s\] | | |
| azint_time | float | Time spent on azimuthal integration \[s\] | | |
| spot_finding_time | float | Time spent on spot finding \[s\] | | |
| indexing_time | float | Time spent on indexing \[s\] | | |
| refinement_time | float | Time spent on refinement of indexing solution and experimental geometry \[s\] | | |
@@ -689,6 +689,8 @@ namespace {
message.processing_time_s = GetCBORFloat(value);
else if (key == "preprocessing_time")
message.preprocessing_time_s = GetCBORFloat(value);
else if (key == "azint_time")
message.azint_time_s = GetCBORFloat(value);
else if (key == "profile_radius")
message.profile_radius = GetCBORFloat(value);
else if (key == "mosaicity")
@@ -705,6 +705,7 @@ void CBORStream2Serializer::SerializeImageInternal(CborEncoder &mapEncoder, cons
CBOR_ENC(mapEncoder, "b_factor", message.b_factor);
CBOR_ENC(mapEncoder, "indexing_time", message.indexing_time_s);
CBOR_ENC(mapEncoder, "processing_time", message.processing_time_s);
CBOR_ENC(mapEncoder, "azint_time", message.azint_time_s);
CBOR_ENC(mapEncoder, "spot_finding_time", message.spot_finding_time_s);
CBOR_ENC(mapEncoder, "bragg_prediction_time", message.bragg_prediction_time_s);
CBOR_ENC(mapEncoder, "integration_time", message.integration_time_s);
+6 -3
View File
@@ -26,7 +26,7 @@ MXAnalysisWithoutFPGA::MXAnalysisWithoutFPGA(const DiffractionExperiment &in_exp
mask_resolution(experiment.GetPixelsNum(), false),
mask_high_res(-1),
mask_low_res(-1) {
preprocessor = std::make_unique<ImagePreprocessor>(in_experiment, in_integration, in_mask, spotFinder->GetInputBuffer());
preprocessor = std::make_unique<ImagePreprocessor>(in_experiment, in_mask, spotFinder->GetInputBuffer());
}
void MXAnalysisWithoutFPGA::Analyze(DataMessage &output,
@@ -48,6 +48,11 @@ void MXAnalysisWithoutFPGA::Analyze(DataMessage &output,
const auto preprocessing_end_time = std::chrono::steady_clock::now();
output.preprocessing_time_s = std::chrono::duration<float>(preprocessing_end_time - preprocessing_start_time).count();
const auto azint_start_time = std::chrono::steady_clock::now();
preprocessor->AzimIntegration(integration, profile);
const auto azint_end_time = std::chrono::steady_clock::now();
output.azint_time_s = std::chrono::duration<float>(azint_end_time - azint_start_time).count();
if (spot_finding_settings.enable) {
// Update resolution mask
if (mask_high_res != spot_finding_settings.high_resolution_limit
@@ -66,8 +71,6 @@ void MXAnalysisWithoutFPGA::Analyze(DataMessage &output,
*prediction);
}
preprocessor->Update(profile);
output.max_viable_pixel_value = ret.max_value;
output.min_viable_pixel_value = ret.min_value;
output.error_pixel_count = ret.error_pixel_count;
@@ -5,18 +5,12 @@
#include "ImagePreprocessor.h"
ImagePreprocessor::ImagePreprocessor(const DiffractionExperiment &experiment,
const AzimuthalIntegration &integration,
const PixelMask &mask,
std::vector<int32_t> &processed_image)
: npixels(experiment.GetPixelsNum()),
experiment(experiment),
integration(integration),
azint_sum(integration.GetBinNumber(), 0.0),
azint_sum2(integration.GetBinNumber(), 0.0),
azint_count(integration.GetBinNumber(), 0),
processed_image(processed_image),
mask_1bit(npixels, false),
azint_bins(integration.GetBinNumber()),
saturation_limit(experiment.GetSaturationLimit()) {
if (processed_image.size() != npixels)
@@ -26,11 +20,6 @@ ImagePreprocessor::ImagePreprocessor(const DiffractionExperiment &experiment,
mask_1bit[i] = (mask.GetMask().at(i) != 0);
}
void ImagePreprocessor::Update(AzimuthalIntegrationProfile &profile) const {
profile.Clear(integration);
profile.Add(azint_sum, azint_count);
}
ImageStatistics ImagePreprocessor::Analyze(const uint8_t *image_ptr, CompressedImageMode image_mode) {
switch (image_mode) {
case CompressedImageMode::Int8:
@@ -54,20 +43,11 @@ template<class T>
ImageStatistics ImagePreprocessor::Analyze(const uint8_t *input, T err_pixel_val, T sat_pixel_val) {
auto image = reinterpret_cast<const T *>(input);
for (int i = 0; i < azint_count.size(); i++) {
azint_sum[i] = 0.0f;
azint_sum2[i] = 0.0f;
azint_count[i] = 0;
}
ImageStatistics ret{};
if (sat_pixel_val > saturation_limit)
sat_pixel_val = static_cast<T>(saturation_limit);
auto &pixel_to_bin = integration.GetPixelToBin();
auto &corrections = integration.Corrections();
for (int i = 0; i < npixels; i++) {
if (mask_1bit[i] != 0) {
processed_image[i] = INT32_MIN;
@@ -86,14 +66,34 @@ ImageStatistics ImagePreprocessor::Analyze(const uint8_t *input, T err_pixel_val
ret.max_value = image[i];
if (image[i] < ret.min_value)
ret.min_value = image[i];
const uint16_t bin = pixel_to_bin[i];
if (bin < azint_bins) {
float val = image[i] * corrections[i];
azint_sum[bin] += val;
++azint_count[bin];
}
}
}
return ret;
}
}
void ImagePreprocessor::AzimIntegration(const AzimuthalIntegration &integration, AzimuthalIntegrationProfile &profile) {
const auto azint_bins = integration.GetBinNumber();
std::vector<float> azint_sum(azint_bins);
std::vector<float> azint_sum2(azint_bins);
std::vector<uint32_t> azint_count(azint_bins);
for (int i = 0; i < azint_count.size(); i++) {
azint_sum[i] = 0.0f;
azint_sum2[i] = 0.0f;
azint_count[i] = 0;
}
auto &pixel_to_bin = integration.GetPixelToBin();
auto &corrections = integration.Corrections();
for (int i = 0; i < npixels; i++) {
const uint16_t bin = pixel_to_bin[i];
if (bin < azint_bins) {
float val = static_cast<float>(processed_image[i]) * corrections[i];
azint_sum[bin] += val;
++azint_count[bin];
}
}
profile.Clear(integration);
profile.Add(azint_sum, azint_count);
}
@@ -23,27 +23,16 @@ class ImagePreprocessor {
protected:
const size_t npixels;
const DiffractionExperiment &experiment;
const AzimuthalIntegration &integration;
std::vector<float> azint_sum;
std::vector<float> azint_sum2;
std::vector<uint32_t> azint_count;
std::vector<int32_t> &processed_image;
std::vector<bool> mask_1bit;
uint16_t azint_bins;
const int64_t saturation_limit;
template <class T> ImageStatistics Analyze(const uint8_t *input, T err_value, T sat_value);
public:
ImagePreprocessor(const DiffractionExperiment &experiment,
const AzimuthalIntegration &integration,
const PixelMask &mask,
std::vector<int32_t> &processed_image);
ImageStatistics Analyze(const uint8_t *decompressed_image, CompressedImageMode image_mode);
void Update(AzimuthalIntegrationProfile &profile) const;
void AzimIntegration(const AzimuthalIntegration &integration, AzimuthalIntegrationProfile &profile);
};
+8
View File
@@ -114,6 +114,7 @@ void JFJochReceiverPlots::Setup(const DiffractionExperiment &experiment, const A
bragg_prediction_time.Clear(r);
preprocessing_time.Clear(r);
compression_time.Clear(r);
azint_time.Clear(r);
}
void JFJochReceiverPlots::Add(const DataMessage &msg, const AzimuthalIntegrationProfile &profile) {
@@ -144,6 +145,7 @@ void JFJochReceiverPlots::Add(const DataMessage &msg, const AzimuthalIntegration
bragg_prediction_time.AddElement(msg.number, msg.bragg_prediction_time_s);
preprocessing_time.AddElement(msg.number, msg.preprocessing_time_s);
compression_time.AddElement(msg.number, msg.compression_time_s);
azint_time.AddElement(msg.number, msg.azint_time_s);
if (msg.indexing_unit_cell) {
indexing_uc_a.AddElement(msg.number, msg.indexing_unit_cell->a);
@@ -379,6 +381,11 @@ MultiLinePlot JFJochReceiverPlots::GetPlots(const PlotRequest &request) {
if (!compression.x.empty())
ret.AddPlot(compression);
auto azint = azint_time.GetMeanPerBin(nbins, start, incr, request.fill_value);
azint.title = "azint";
if (!azint.x.empty())
ret.AddPlot(azint);
auto total = total_processing_time.GetMeanPerBin(nbins, start, incr, request.fill_value);
total.title = "total";
ret.AddPlot(total);
@@ -429,6 +436,7 @@ MeanProcessingTime JFJochReceiverPlots::GetMeanProcessingTime() const {
ret.bragg_prediction = bragg_prediction_time.Mean();
ret.processing = total_processing_time.Mean();
ret.preprocessing = preprocessing_time.Mean();
ret.azint = azint_time.Mean();
return ret;
}
+3 -1
View File
@@ -20,13 +20,14 @@
struct MeanProcessingTime {
float compression;
float preprocessing;
float azint;
float spot_finding;
float indexing;
float refinement;
float integration;
float bragg_prediction;
float processing;
float preprocessing;
};
class JFJochReceiverPlots {
@@ -95,6 +96,7 @@ class JFJochReceiverPlots {
StatusVector total_processing_time;
StatusVector preprocessing_time;
StatusVector compression_time;
StatusVector azint_time;
MultiLinePlot GetROIPlot(PlotType type, int64_t nbins, float start, float incr,
const std::optional<float> &fill_value) const;
+2 -1
View File
@@ -765,9 +765,10 @@ int main(int argc, char **argv) {
}
auto image_mean_time = plots.GetMeanProcessingTime();
logger.Info("Per-image time: (mean; microseconds): decompress {:.0f} preprocess {:.0f} spot finding {:.0f} indexing {:.0f} refinement {:.0f} prediction {:.0f} integration {:.0f} total {:.0f}",
logger.Info("Per-image time: (mean; microseconds): decompress {:.0f} preprocess {:.0f} azint {:.0f} spot finding {:.0f} indexing {:.0f} refinement {:.0f} prediction {:.0f} integration {:.0f} total {:.0f}",
image_mean_time.compression * 1e6,
image_mean_time.preprocessing * 1e6,
image_mean_time.azint * 1e6,
image_mean_time.spot_finding * 1e6,
image_mean_time.indexing * 1e6,
image_mean_time.refinement * 1e6,
+5
View File
@@ -102,6 +102,7 @@ void HDF5DataFilePluginMX::OpenFile(HDF5File &data_file, const DataMessage &msg,
bragg_prediction_time.reserve(images_per_file);
preprocessing_time.reserve(images_per_file);
compression_time.reserve(images_per_file);
azint_time.reserve(images_per_file);
}
void HDF5DataFilePluginMX::Write(const DataMessage &msg, uint64_t image_number) {
@@ -159,6 +160,7 @@ void HDF5DataFilePluginMX::Write(const DataMessage &msg, uint64_t image_number)
bragg_prediction_time[image_number] = msg.bragg_prediction_time_s.value_or(NAN);
preprocessing_time[image_number] = msg.preprocessing_time_s.value_or(NAN);
compression_time[image_number] = msg.compression_time_s.value_or(NAN);
azint_time[image_number] = msg.azint_time_s.value_or(NAN);
if (indexing) {
indexed[image_number] = msg.indexing_result.value_or(0);
@@ -226,6 +228,9 @@ void HDF5DataFilePluginMX::WriteFinal(HDF5File &data_file) {
data_file.SaveVector("/entry/MX/preprocessingTime", preprocessing_time.vec())->Units("s");
if (!compression_time.empty())
data_file.SaveVector("/entry/MX/compressionTime", compression_time.vec())->Units("s");
if (!azint_time.empty())
data_file.SaveVector("/entry/MX/azIntTime", azint_time.vec())->Units("s");
if (!strong_pixel_count.empty())
data_file.SaveVector("/entry/MX/strongPixels", strong_pixel_count.vec());
+1
View File
@@ -58,6 +58,7 @@ class HDF5DataFilePluginMX : public HDF5DataFilePlugin {
AutoIncrVector<float> bragg_prediction_time;
AutoIncrVector<float> preprocessing_time;
AutoIncrVector<float> compression_time;
AutoIncrVector<float> azint_time;
public:
explicit HDF5DataFilePluginMX(const StartMessage& msg);
void OpenFile(HDF5File &data_file, const DataMessage& msg, size_t images_per_file) override;