Remove RotationScanAccumulator
This commit is contained in:
@@ -31,8 +31,6 @@ ADD_LIBRARY(JFJochImageAnalysis STATIC
|
||||
RotationIndexer.h
|
||||
RotationParameters.cpp
|
||||
RotationParameters.h
|
||||
RotationSpotAccumulator.cpp
|
||||
RotationSpotAccumulator.h
|
||||
WriteMmcif.cpp
|
||||
WriteMmcif.h)
|
||||
|
||||
|
||||
@@ -1,141 +0,0 @@
|
||||
//
|
||||
// Created by jungfrau on 2/4/26.
|
||||
//
|
||||
|
||||
#include "RotationSpotAccumulator.h"
|
||||
|
||||
RotationSpot3D::RotationSpot3D(SpotToSave s, float phi_deg, int64_t image) {
|
||||
x = s.x * s.intensity;
|
||||
y = s.y * s.intensity;
|
||||
phi = phi_deg * s.intensity;
|
||||
intensity = s.intensity;
|
||||
maxc = s.maxc;
|
||||
first_image = image;
|
||||
last_image = image;
|
||||
}
|
||||
|
||||
void RotationSpot3D::Add(const SpotToSave& s, float phi_deg, int64_t image) {
|
||||
x += s.x * s.intensity;
|
||||
y += s.y * s.intensity;
|
||||
phi += phi_deg * s.intensity;
|
||||
intensity += s.intensity;
|
||||
maxc = std::max(maxc, s.maxc);
|
||||
ice_ring = ice_ring || s.ice_ring;
|
||||
first_image = std::min(first_image, image);
|
||||
last_image = std::max(last_image, image);
|
||||
}
|
||||
|
||||
float RotationSpot3D::calcX() const {
|
||||
return x / intensity;
|
||||
}
|
||||
|
||||
float RotationSpot3D::calcY() const {
|
||||
return y / intensity;
|
||||
}
|
||||
|
||||
float RotationSpot3D::calcPhi() const {
|
||||
return phi / intensity;
|
||||
}
|
||||
|
||||
float RotationSpot3D::calcI() const {
|
||||
return intensity;
|
||||
}
|
||||
|
||||
int64_t RotationSpot3D::calcMAXC() const {
|
||||
return maxc;
|
||||
}
|
||||
|
||||
int64_t RotationSpot3D::getLastImage() const {
|
||||
return last_image;
|
||||
}
|
||||
|
||||
void RotationSpotAccumulator::AddImage(int64_t image, const std::vector<SpotToSave>& in_spots) {
|
||||
std::lock_guard<std::mutex> lock(mutex_);
|
||||
spots[image] = in_spots;
|
||||
}
|
||||
|
||||
|
||||
std::vector<RotationSpot3D> RotationSpotAccumulator::FinalizeAll() {
|
||||
std::lock_guard<std::mutex> lock(mutex_);
|
||||
|
||||
int64_t image0 = INT64_MIN;
|
||||
|
||||
std::vector<RotationSpot3D> completed_;
|
||||
std::vector<RotationSpot3D> pending_;
|
||||
|
||||
for (auto& [image, v]: spots) {
|
||||
float image_angle = axis_.GetIncrement_deg() * image + axis_.GetWedge_deg() / 2.0f;
|
||||
|
||||
std::vector<RotationSpot3D> tmp;
|
||||
tmp.reserve(pending_.size() + v.size());
|
||||
|
||||
if (image == image0 + 1) {
|
||||
for (const auto &p: pending_) {
|
||||
if (p.getLastImage() != image - 1) {
|
||||
completed_.push_back(p);
|
||||
} else {
|
||||
tmp.push_back(p);
|
||||
}
|
||||
}
|
||||
|
||||
// Spatial hash for tmp (pending)
|
||||
const float cell_size = config_.grid_cell_size;
|
||||
const float inv_cell = (cell_size > 0.0f) ? (1.0f / cell_size) : 1.0f;
|
||||
|
||||
std::unordered_map<int64_t, std::vector<int>> grid;
|
||||
grid.reserve(tmp.size() * 2);
|
||||
|
||||
auto cell_key = [](int32_t cx, int32_t cy) -> int64_t {
|
||||
return (static_cast<int64_t>(cx) << 32) ^ static_cast<uint32_t>(cy);
|
||||
};
|
||||
|
||||
for (int i = 0; i < static_cast<int>(tmp.size()); ++i) {
|
||||
const int32_t cx = static_cast<int32_t>(std::floor(tmp[i].calcX() * inv_cell));
|
||||
const int32_t cy = static_cast<int32_t>(std::floor(tmp[i].calcY() * inv_cell));
|
||||
grid[cell_key(cx, cy)].push_back(i);
|
||||
}
|
||||
|
||||
for (const auto &s: v) {
|
||||
bool matched = false;
|
||||
|
||||
const int32_t scx = static_cast<int32_t>(std::floor(s.x * inv_cell));
|
||||
const int32_t scy = static_cast<int32_t>(std::floor(s.y * inv_cell));
|
||||
|
||||
for (int dx = -1; dx <= 1 && !matched; ++dx) {
|
||||
for (int dy = -1; dy <= 1 && !matched; ++dy) {
|
||||
auto it = grid.find(cell_key(scx + dx, scy + dy));
|
||||
if (it == grid.end())
|
||||
continue;
|
||||
|
||||
for (int idx : it->second) {
|
||||
float dist_x = s.x - tmp[idx].calcX();
|
||||
float dist_y = s.y - tmp[idx].calcY();
|
||||
|
||||
if (dist_x * dist_x + dist_y * dist_y < config_.xy_tolerance_pxl_sq) {
|
||||
tmp[idx].Add(s, image_angle, image);
|
||||
matched = true;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (!matched)
|
||||
tmp.emplace_back(s, image_angle, image);
|
||||
}
|
||||
} else {
|
||||
for (const auto &p: pending_)
|
||||
completed_.push_back(p);
|
||||
|
||||
for (const auto &s: v)
|
||||
tmp.emplace_back(s, image_angle, image);
|
||||
}
|
||||
pending_ = std::move(tmp);
|
||||
image0 = image;
|
||||
}
|
||||
|
||||
for (const auto &p: pending_)
|
||||
completed_.push_back(p);
|
||||
|
||||
return completed_;
|
||||
}
|
||||
@@ -1,70 +0,0 @@
|
||||
// RotationSpotAccumulator.h
|
||||
// SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
|
||||
// SPDX-License-Identifier: GPL-3.0-only
|
||||
|
||||
#pragma once
|
||||
|
||||
#include <vector>
|
||||
#include <map>
|
||||
#include <unordered_map>
|
||||
#include <unordered_set>
|
||||
#include <cmath>
|
||||
#include <mutex>
|
||||
#include <optional>
|
||||
|
||||
#include "../common/SpotToSave.h"
|
||||
#include "../common/GoniometerAxis.h"
|
||||
#include "../common/DiffractionGeometry.h"
|
||||
|
||||
// Forward declaration
|
||||
class RotationSpot3D {
|
||||
float x = 0;
|
||||
float y = 0;
|
||||
float phi = 0;
|
||||
float intensity = 0;
|
||||
int64_t maxc = 0;
|
||||
bool ice_ring = false;
|
||||
|
||||
int64_t first_image = INT64_MAX;
|
||||
int64_t last_image = INT64_MIN;
|
||||
|
||||
public:
|
||||
RotationSpot3D(SpotToSave s, float phi_deg, int64_t image);
|
||||
void Add(const SpotToSave& s, float phi_deg, int64_t image);
|
||||
|
||||
float calcX() const;
|
||||
float calcY() const;
|
||||
float calcPhi() const;
|
||||
float calcI() const;
|
||||
int64_t calcMAXC() const;
|
||||
int64_t getLastImage() const;
|
||||
};
|
||||
|
||||
class RotationSpotAccumulator {
|
||||
public:
|
||||
struct Config {
|
||||
float xy_tolerance_pxl_sq = 2.0f * 2.0f;
|
||||
float grid_cell_size = 10.0f; // pixels
|
||||
};
|
||||
|
||||
private:
|
||||
mutable std::mutex mutex_;
|
||||
|
||||
const DiffractionGeometry& geom_;
|
||||
const GoniometerAxis& axis_;
|
||||
Config config_;
|
||||
|
||||
std::map<int64_t, std::vector<SpotToSave>> spots;
|
||||
|
||||
public:
|
||||
RotationSpotAccumulator(const DiffractionGeometry& geom,
|
||||
const GoniometerAxis& axis,
|
||||
Config config)
|
||||
: geom_(geom), axis_(axis), config_(config) {}
|
||||
|
||||
// Add spots from one image (thread-safe, handles out-of-order)
|
||||
void AddImage(int64_t image, const std::vector<SpotToSave>& in_spots);
|
||||
|
||||
// Finalize all remaining spots (call at end of dataset)
|
||||
std::vector<RotationSpot3D> FinalizeAll();
|
||||
};
|
||||
@@ -64,7 +64,6 @@ ADD_EXECUTABLE(jfjoch_test
|
||||
RotationIndexerTest.cpp
|
||||
TopPixelsTest.cpp
|
||||
HKLKeyTest.cpp
|
||||
RotationScanAccumulatorTest.cpp
|
||||
)
|
||||
|
||||
target_link_libraries(jfjoch_test Catch2WithMain JFJochBroker JFJochReceiver JFJochReader JFJochWriter JFJochImageAnalysis JFJochCommon JFJochHLSSimulation JFJochPreview)
|
||||
|
||||
@@ -1,70 +0,0 @@
|
||||
// SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
|
||||
// SPDX-License-Identifier: GPL-3.0-only
|
||||
|
||||
#include <iostream>
|
||||
|
||||
#include <catch2/catch_all.hpp>
|
||||
#include "../image_analysis//RotationSpotAccumulator.h"
|
||||
|
||||
TEST_CASE("RotationSpotAccumulator") {
|
||||
DiffractionGeometry geometry;
|
||||
geometry.BeamX_pxl(1000).BeamY_pxl(1000).DetectorDistance_mm(100);
|
||||
|
||||
GoniometerAxis axis("U", 0, 0.1, Coord{1, 0, 0}, std::nullopt);
|
||||
|
||||
RotationSpotAccumulator::Config config{};
|
||||
RotationSpotAccumulator accumulator(geometry, axis, config);
|
||||
|
||||
accumulator.AddImage(2, {
|
||||
SpotToSave{.x = 500, .y = 500, .intensity = 100},
|
||||
SpotToSave{.x = 101, .y = 101, .intensity = 50}
|
||||
});
|
||||
|
||||
accumulator.AddImage(0, {
|
||||
SpotToSave{.x = 100, .y = 100, .intensity = 50},
|
||||
SpotToSave{.x = 200, .y = 200, .intensity = 50}
|
||||
});
|
||||
|
||||
accumulator.AddImage(1, {
|
||||
SpotToSave{.x = 500, .y = 500, .intensity = 100},
|
||||
SpotToSave{.x = 101, .y = 100, .intensity = 100}
|
||||
});
|
||||
|
||||
accumulator.AddImage(4, {
|
||||
SpotToSave{.x = 500, .y = 500, .intensity = 100}
|
||||
});
|
||||
|
||||
auto v = accumulator.FinalizeAll();
|
||||
|
||||
struct Expect {
|
||||
float x, y, phi, I;
|
||||
};
|
||||
|
||||
std::vector<Expect> expected = {
|
||||
{200.0f, 200.0f, 0.05f, 50.0f},
|
||||
{500.0f, 500.0f, 0.2f, 200.0f},
|
||||
{100.75f, 100.25f, 0.15f, 200.0f},
|
||||
{500.0f, 500.0f, 0.45f, 100.0f}
|
||||
};
|
||||
|
||||
auto key = [](const Expect& e) { return std::pair<float,float>(e.x, e.y); };
|
||||
|
||||
std::vector<Expect> got;
|
||||
got.reserve(v.size());
|
||||
for (const auto& s : v) {
|
||||
got.push_back({s.calcX(), s.calcY(), s.calcPhi(), s.calcI()});
|
||||
}
|
||||
|
||||
std::sort(expected.begin(), expected.end(),
|
||||
[&](const Expect& a, const Expect& b){ return key(a) < key(b); });
|
||||
std::sort(got.begin(), got.end(),
|
||||
[&](const Expect& a, const Expect& b){ return key(a) < key(b); });
|
||||
|
||||
REQUIRE(got.size() == expected.size());
|
||||
for (size_t i = 0; i < got.size(); ++i) {
|
||||
CHECK(got[i].x == Catch::Approx(expected[i].x));
|
||||
CHECK(got[i].y == Catch::Approx(expected[i].y));
|
||||
CHECK(got[i].phi == Catch::Approx(expected[i].phi));
|
||||
CHECK(got[i].I == Catch::Approx(expected[i].I));
|
||||
}
|
||||
}
|
||||
Reference in New Issue
Block a user