From 058468a0b315eecc553ca5dfe5fae1c3dc8639fd Mon Sep 17 00:00:00 2001 From: leonarski_f Date: Thu, 28 May 2026 16:02:52 +0200 Subject: [PATCH] jfjoch_process: Add orientation-only refinement --- common/IndexingSettings.h | 2 +- image_analysis/IndexAndRefine.cpp | 5 +++ .../geom_refinement/XtalOptimizer.cpp | 6 ++-- tools/jfjoch_process.cpp | 35 ++++++++++++------- 4 files changed, 31 insertions(+), 17 deletions(-) diff --git a/common/IndexingSettings.h b/common/IndexingSettings.h index 35eda787..7f603f93 100644 --- a/common/IndexingSettings.h +++ b/common/IndexingSettings.h @@ -7,7 +7,7 @@ #include enum class IndexingAlgorithmEnum {FFBIDX, FFT, FFTW, Auto, None}; -enum class GeomRefinementAlgorithmEnum {None, BeamCenter}; +enum class GeomRefinementAlgorithmEnum {None, OrientationOnly, BeamCenter}; class IndexingSettings { IndexingAlgorithmEnum algorithm; diff --git a/image_analysis/IndexAndRefine.cpp b/image_analysis/IndexAndRefine.cpp index 08a91678..f888528e 100644 --- a/image_analysis/IndexAndRefine.cpp +++ b/image_analysis/IndexAndRefine.cpp @@ -128,6 +128,11 @@ void IndexAndRefine::RefineGeometryIfNeeded(DataMessage &msg, IndexAndRefine::In switch (experiment.GetIndexingSettings().GetGeomRefinementAlgorithm()) { case GeomRefinementAlgorithmEnum::None: break; + case GeomRefinementAlgorithmEnum::OrientationOnly: + XtalOptimizerRotationOnly(data, msg.spots, 0.2); + XtalOptimizerRotationOnly(data, msg.spots, 0.1); + XtalOptimizerRotationOnly(data, msg.spots, 0.05); + break; case GeomRefinementAlgorithmEnum::BeamCenter: if (XtalOptimizer(data, msg.spots)) { outcome.experiment.BeamX_pxl(data.geom.GetBeamX_pxl()) diff --git a/image_analysis/geom_refinement/XtalOptimizer.cpp b/image_analysis/geom_refinement/XtalOptimizer.cpp index 3771216a..4df7a6ae 100644 --- a/image_analysis/geom_refinement/XtalOptimizer.cpp +++ b/image_analysis/geom_refinement/XtalOptimizer.cpp @@ -825,9 +825,9 @@ bool XtalOptimizerRotationOnly(XtalOptimizerData &data, ceres::AngleAxisToRotationMatrix(rot_aa, R_raw); // row-major 3x3 Eigen::Matrix3d R; - R << R_raw[0], R_raw[1], R_raw[2], - R_raw[3], R_raw[4], R_raw[5], - R_raw[6], R_raw[7], R_raw[8]; + R << R_raw[0], R_raw[3], R_raw[6], + R_raw[1], R_raw[4], R_raw[7], + R_raw[2], R_raw[5], R_raw[8]; const Eigen::Vector3d A(a0.x, a0.y, a0.z); const Eigen::Vector3d B(b0.x, b0.y, b0.z); diff --git a/tools/jfjoch_process.cpp b/tools/jfjoch_process.cpp index e7b5e325..a03c263b 100644 --- a/tools/jfjoch_process.cpp +++ b/tools/jfjoch_process.cpp @@ -55,7 +55,7 @@ void print_usage() { std::cout << " -X, --indexing-algorithm Indexing algorithm (FFBIDX|FFT|FFTW|Auto|None)" << std::endl; std::cout << " -S, --space-group Space group number - used for both indexing and scaling" << std::endl; std::cout << " -C, --unit-cell Fix reference unit cell: \"a,b,c,alpha,beta,gamma\"" << std::endl; - std::cout << " --no-beam-center-refine No least-square beam center refinement" << std::endl; + std::cout << " -r, --refine Geometry refinement algorithm (none|orientation|beam_and_lattice)" << std::endl; std::cout << std::endl; std::cout << " Scaling and merging" << std::endl; @@ -77,7 +77,6 @@ enum { OPT_SPOT_THRESHOLD, OPT_SPOT_RESOLUTION, OPT_MAX_SPOTS, - OPT_NO_BEAM_CENTER_REFINE, OPT_MIN_PARTIALITY, OPT_MIN_IMAGE_CC, OPT_SCALING_ITERATIONS, @@ -102,12 +101,12 @@ static option long_options[] = { {"refine-bfactor", no_argument, nullptr, 'B'}, {"wedge", optional_argument, nullptr, 'w'}, {"scale-merge", no_argument, nullptr, 'M'}, + {"refine", required_argument, nullptr, 'r'}, {"spot-sigma", required_argument, nullptr, OPT_SPOT_SIGMA}, {"spot-threshold", required_argument, nullptr, OPT_SPOT_THRESHOLD}, {"spot-high-resolution", required_argument, nullptr, OPT_SPOT_RESOLUTION}, {"max-spots", required_argument, nullptr, OPT_MAX_SPOTS}, - {"no-beam-center-refine", no_argument, nullptr, OPT_NO_BEAM_CENTER_REFINE}, {"min-partiality", required_argument, nullptr, OPT_MIN_PARTIALITY}, {"min-image-cc", required_argument, nullptr, OPT_MIN_IMAGE_CC}, {"scaling-iterations", required_argument, nullptr, OPT_SCALING_ITERATIONS}, @@ -199,7 +198,6 @@ int main(int argc, char **argv) { bool verbose = false; bool rotation_indexing = false; std::optional rotation_indexing_range; - bool refine_beam_center = true; bool run_scaling = false; bool anomalous_mode = false; std::optional space_group_number; @@ -216,6 +214,7 @@ int main(int argc, char **argv) { int64_t scaling_iter = 3; IndexingAlgorithmEnum indexing_algorithm = IndexingAlgorithmEnum::Auto; + GeomRefinementAlgorithmEnum refinement_algorithm = GeomRefinementAlgorithmEnum::BeamCenter; IntensityFormat intensity_format = IntensityFormat::MTZ; PartialityModel partiality_model = PartialityModel::Fixed; @@ -230,7 +229,7 @@ int main(int argc, char **argv) { int opt; int option_index = 0; - const char *short_opts = "vo:N:s:e:t:R::X:C:z:FABw::S:MP:"; + const char *short_opts = "vo:N:s:e:t:R::X:C:z:FABw::S:MP:r:"; while ((opt = getopt_long(argc, argv, short_opts, long_options, &option_index)) != -1) { switch (opt) { @@ -278,6 +277,23 @@ int main(int argc, char **argv) { } break; } + case 'r': { + std::string alg = optarg ? optarg : ""; + std::transform(alg.begin(), alg.end(), alg.begin(), + [](unsigned char c) { return static_cast(std::tolower(c)); }); + if (alg == "none") + refinement_algorithm = GeomRefinementAlgorithmEnum::None; + else if (alg == "beam_and_lattice") + refinement_algorithm = GeomRefinementAlgorithmEnum::BeamCenter; + else if (alg == "orientation") + refinement_algorithm = GeomRefinementAlgorithmEnum::OrientationOnly; + else { + logger.Error("Invalid geom refinement algorithm: {}", alg); + print_usage(); + exit(EXIT_FAILURE); + } + break; + } case 'C': { auto uc = parse_unit_cell_arg(optarg); if (!uc.has_value()) { @@ -345,9 +361,6 @@ int main(int argc, char **argv) { max_spot_count_override = atoll(optarg); logger.Info("Max spot count overridden to {}", max_spot_count_override.value()); break; - case OPT_NO_BEAM_CENTER_REFINE: - refine_beam_center = false; - break; case 'M': run_scaling = true; break; @@ -476,11 +489,7 @@ int main(int argc, char **argv) { indexing_settings.RotationIndexing(rotation_indexing); if (rotation_indexing_range.has_value()) indexing_settings.RotationIndexingMinAngularRange_deg(rotation_indexing_range.value()); - - if (refine_beam_center) - indexing_settings.GeomRefinementAlgorithm(GeomRefinementAlgorithmEnum::BeamCenter); - else - indexing_settings.GeomRefinementAlgorithm(GeomRefinementAlgorithmEnum::None); + indexing_settings.GeomRefinementAlgorithm(refinement_algorithm); experiment.ImportIndexingSettings(indexing_settings); ScalingSettings scaling_settings;