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

BraggPrediction: Remove an assumption on S0 being (0,0,1/lambda) + check with the viewer predictions are OK

Browse Source
This commit is contained in:
leonarski_f
2025-10-04 18:31:35 +02:00
parent 9cad8dd9a5
commit d801b17b81
1 changed files with 22 additions and 25 deletions
Download Patch File Download Diff File Expand all files Collapse all files
image_analysis/bragg_integration/BraggPrediction.cpp
+22 -25
View File
@@ -8,7 +8,6 @@ std::multimap<float, Reflection> CalcBraggPredictions(const DiffractionExperimen
float high_res_A,
float ewald_dist_cutoff,
int max_hkl) {
// Assumption: Scattering vector is (0,0,1)
std::multimap<float, Reflection> reflections;
auto geom = experiment.GetDiffractionGeometry();
@@ -18,11 +17,12 @@ std::multimap<float, Reflection> CalcBraggPredictions(const DiffractionExperimen
float one_over_dmax = 1.0f / high_res_A;
float one_over_dmax_sq = one_over_dmax * one_over_dmax;
float one_over_wavelength = 1.0f/ geom.GetWavelength_A();
float one_over_wavelength = 1.0f / geom.GetWavelength_A();
Coord Astar = lattice.Astar();
Coord Bstar = lattice.Bstar();
Coord Cstar = lattice.Cstar();
Coord S0 = geom.GetScatteringVector();
for (int h = -max_hkl; h < max_hkl; h++) {
// Precompute A* h contribution
@@ -45,38 +45,35 @@ std::multimap<float, Reflection> CalcBraggPredictions(const DiffractionExperimen
float recip_y = AhBk_y + Cstar.y * l;
float recip_z = AhBk_z + Cstar.z * l;
// Check resolution (actually through 1/d^2)
float dot_xy = recip_y * recip_y + recip_z * recip_z;
float dot = dot_xy + recip_x * recip_x;
float dot = recip_x * recip_x + recip_y * recip_y + recip_z * recip_z;
if (dot > one_over_dmax_sq)
continue;
// Get Distance from Ewald sphere, to check if reflection is in diffracting condition
float S_z = recip_z + one_over_wavelength;
float S_len = sqrtf(dot_xy + S_z * S_z);
float S_x = recip_x + S0.x;
float S_y = recip_y + S0.y;
float S_z = recip_z + S0.z;
float S_len = sqrtf(S_x * S_x + S_y * S_y + S_z * S_z);
float dist_ewald_sphere = std::fabs(S_len - one_over_wavelength);
if (!std::isfinite(dist_ewald_sphere) || dist_ewald_sphere > ewald_dist_cutoff )
continue;
if (dist_ewald_sphere <= ewald_dist_cutoff ) {
Coord recip(recip_x, recip_y, recip_z);
Coord recip(recip_x, recip_y, recip_z);
auto [x,y] = geom.RecipToDector(recip);
auto [x,y] = geom.RecipToDector(recip);
if ((x < 0) || (x >= det_width_pxl) || (y < 0) || (y >= det_height_pxl))
continue;
if ((x < 0) || (x >= det_width_pxl) || (y < 0) || (y >= det_height_pxl))
continue;
float d = 1.0f / sqrtf(dot);
float d = 1.0f / sqrtf(dot);
reflections.insert(std::make_pair(dist_ewald_sphere, Reflection{
.h = h,
.k = k,
.l = l,
.predicted_x = x,
.predicted_y = y,
.d = d
}));
reflections.insert(std::make_pair(dist_ewald_sphere, Reflection{
.h = h,
.k = k,
.l = l,
.predicted_x = x,
.predicted_y = y,
.d = d
}));
}
}
}
}