mx/Jungfraujoch
2
1
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases 67 Wiki Activity

XtalOptimizer: Include explicit rotation axis

Browse Source
This commit is contained in:
leonarski_f
2026-01-30 14:02:11 +01:00
parent 475ba92d01
commit a4ac885d1c
1 changed files with 26 additions and 14 deletions
Download Patch File Download Diff File Expand all files Collapse all files
image_analysis/geom_refinement/XtalOptimizer.cpp
+26 -14
View File
@@ -9,7 +9,8 @@ struct XtalResidual {
XtalResidual(double x, double y,
double lambda,
double pixel_size,
const Eigen::Matrix3d &gonio_back_rot,
Coord rot_axis,
double angle_rad,
double exp_h, double exp_k,
double exp_l,
gemmi::CrystalSystem symmetry)
@@ -19,7 +20,8 @@ struct XtalResidual {
exp_h(exp_h),
exp_k(exp_k),
exp_l(exp_l),
gonio_back_rot_(gonio_back_rot),
rot_axis(rot_axis),
angle_rad(angle_rad),
symmetry(symmetry) {
}
@@ -59,7 +61,15 @@ struct XtalResidual {
// Apply goniometer "back-to-start" rotation:
// brings observed reciprocal from image orientation into reference crystal frame
Eigen::Matrix<T, 3, 3> R_gonio_back = gonio_back_rot_.cast<T>();
T a[3];
a[0] = T(angle_rad*rot_axis.x);
a[1] = T(angle_rad*rot_axis.y);
a[2] = T(angle_rad*rot_axis.z);
T rot_arr[9];
ceres::AngleAxisToRotationMatrix(a, rot_arr);
Eigen::Matrix<T, 3, 3> R_gonio_back(rot_arr);
Eigen::Matrix<T, 3, 1> e_obs_recip = R_gonio_back * e_obs_recip_raw;
Eigen::Matrix<T, 3, 1> e_pred;
@@ -139,7 +149,8 @@ struct XtalResidual {
const double exp_h;
const double exp_k;
const double exp_l;
const Eigen::Matrix3d gonio_back_rot_;
const Coord rot_axis;
const double angle_rad;
gemmi::CrystalSystem symmetry;
};
@@ -419,9 +430,12 @@ bool XtalOptimizerInternal(XtalOptimizerData &data,
for (const auto &pt: spots) {
if (!data.index_ice_rings && pt.ice_ring)
continue;
Eigen::Matrix3d gonio_back_rot = Eigen::Matrix3d::Identity();
Coord axis = Coord(1,0,0);
float angle_rad = 0.0;
Coord recip = data.geom.DetectorToRecip(pt.x, pt.y);
if (data.axis) {
@@ -429,13 +443,8 @@ bool XtalOptimizerInternal(XtalOptimizerData &data,
auto rot = data.axis->GetTransformationAngle(angle_deg);
recip = rot * recip;
Eigen::Matrix3d Rg;
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < 3; ++j) {
Rg(i, j) = static_cast<double>(rot.arr()[i + j * 3]);
}
}
gonio_back_rot = Rg.transpose();
angle_rad = angle_deg * M_PI / 180.0;
axis = data.axis->GetAxis();
}
double h_fp = recip * vec0;
@@ -456,7 +465,8 @@ bool XtalOptimizerInternal(XtalOptimizerData &data,
new XtalResidual(pt.x, pt.y,
data.geom.GetWavelength_A(),
data.geom.GetPixelSize_mm(),
gonio_back_rot,
axis,
angle_rad,
h, k, l,
data.crystal_system)),
nullptr,
@@ -547,7 +557,9 @@ bool XtalOptimizerInternal(XtalOptimizerData &data,
data.geom.PoniRot1_rad(detector_rot[0]).PoniRot2_rad(detector_rot[1]);
if (!data.refine_unit_cell) {
Coord rot_vector{static_cast<float>(latt_vec0[0]), static_cast<float>(latt_vec0[1]), static_cast<float>(latt_vec0[2])};
Coord rot_vector{
static_cast<float>(latt_vec0[0]), static_cast<float>(latt_vec0[1]), static_cast<float>(latt_vec0[2])
};
Coord rot_vector_norm = rot_vector.Normalize();
std::cout << "X " << rot_vector.Length() * 180.0 / M_PI << " " << rot_vector_norm << std::endl;
}