aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2021-05-28 15:05:35 +0200
committerThomas White <taw@physics.org>2023-07-28 13:22:05 +0200
commit5007a946fab23f0363020918dd5acf48294b9dde (patch)
tree10ec423d665da5827e844621fa2e57976ad72362 /libcrystfel
parent6efd03ddd0699eb4ddc5a1f1e1885cae373b49d7 (diff)
Initial re-parameterisation of prediction refinement with fs/ss coordinates
The gradients need to be updated.
Diffstat (limited to 'libcrystfel')
-rw-r--r--libcrystfel/src/geometry.c29
-rw-r--r--libcrystfel/src/geometry.h10
-rw-r--r--libcrystfel/src/predict-refine.c51
-rw-r--r--libcrystfel/src/predict-refine.h9
4 files changed, 39 insertions, 60 deletions
diff --git a/libcrystfel/src/geometry.c b/libcrystfel/src/geometry.c
index 8e8a445d..5289856b 100644
--- a/libcrystfel/src/geometry.c
+++ b/libcrystfel/src/geometry.c
@@ -1086,9 +1086,11 @@ void polarisation_correction(RefList *list, UnitCell *cell,
}
-/* Returns dx_h/dP, where P = any parameter */
-double x_gradient(int param, Reflection *refl, UnitCell *cell,
- struct detgeom_panel *p)
+/* Returns dfs_refl/dP, where P = any parameter *
+ * fs_refl is the fast scan position of refl in panel 'p',
+ * in metres (not pixels) */
+double fs_gradient(int param, Reflection *refl, UnitCell *cell,
+ struct detgeom_panel *p)
{
signed int h, k, l;
double xl, zl, kpred;
@@ -1140,16 +1142,20 @@ double x_gradient(int param, Reflection *refl, UnitCell *cell,
case GPARAM_CLEN :
return -xl / (kpred+zl);
+ default:
+ ERROR("Positional gradient requested for parameter %i?\n", param);
+ abort();
+
}
- ERROR("Positional gradient requested for parameter %i?\n", param);
- abort();
}
-/* Returns dy_h/dP, where P = any parameter */
-double y_gradient(int param, Reflection *refl, UnitCell *cell,
- struct detgeom_panel *p)
+/* Returns dss_refl/dP, where P = any parameter
+ * ss_refl is the slow scan position of refl in panel 'p',
+ * in metres (not pixels) */
+double ss_gradient(int param, Reflection *refl, UnitCell *cell,
+ struct detgeom_panel *p)
{
signed int h, k, l;
double yl, zl, kpred;
@@ -1201,8 +1207,9 @@ double y_gradient(int param, Reflection *refl, UnitCell *cell,
case GPARAM_CLEN :
return -yl / (kpred+zl);
- }
+ default :
+ ERROR("Positional gradient requested for parameter %i?\n", param);
+ abort();
- ERROR("Positional gradient requested for parameter %i?\n", param);
- abort();
+ }
}
diff --git a/libcrystfel/src/geometry.h b/libcrystfel/src/geometry.h
index 21bd4d39..a7057f12 100644
--- a/libcrystfel/src/geometry.h
+++ b/libcrystfel/src/geometry.h
@@ -8,7 +8,7 @@
* Copyright © 2012 Richard Kirian
*
* Authors:
- * 2010-2020 Thomas White <taw@physics.org>
+ * 2010-2021 Thomas White <taw@physics.org>
* 2012 Richard Kirian
*
* This file is part of CrystFEL.
@@ -120,10 +120,10 @@ extern void polarisation_correction(RefList *list, UnitCell *cell,
extern double sphere_fraction(double rlow, double rhigh, double pr);
extern double gaussian_fraction(double rlow, double rhigh, double pr);
-extern double x_gradient(int param, Reflection *refl, UnitCell *cell,
- struct detgeom_panel *p);
-extern double y_gradient(int param, Reflection *refl, UnitCell *cell,
- struct detgeom_panel *p);
+extern double fs_gradient(int param, Reflection *refl, UnitCell *cell,
+ struct detgeom_panel *p);
+extern double ss_gradient(int param, Reflection *refl, UnitCell *cell,
+ struct detgeom_panel *p);
#ifdef __cplusplus
}
diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c
index 5210846a..1ed3bfed 100644
--- a/libcrystfel/src/predict-refine.c
+++ b/libcrystfel/src/predict-refine.c
@@ -68,19 +68,6 @@ static const enum gparam rv[] =
};
-void twod_mapping(double fs, double ss, double *px, double *py,
- struct detgeom_panel *p, double dx, double dy)
-{
- double xs, ys;
-
- xs = fs*p->fsx + ss*p->ssx; /* pixels */
- ys = fs*p->fsy + ss*p->ssy; /* pixels */
-
- *px = (xs + p->cnx) * p->pixel_pitch + dx; /* metres */
- *py = (ys + p->cny) * p->pixel_pitch + dy; /* metres */
-}
-
-
static double r_dev(struct reflpeak *rp)
{
/* Excitation error term */
@@ -88,29 +75,19 @@ static double r_dev(struct reflpeak *rp)
}
-double x_dev(struct reflpeak *rp, struct detgeom *det,
- double dx, double dy)
+double fs_dev(struct reflpeak *rp, struct detgeom *det)
{
- /* Peak position term */
- double xpk, ypk, xh, yh;
double fsh, ssh;
- twod_mapping(rp->peak->fs, rp->peak->ss, &xpk, &ypk, rp->panel, dx, dy);
get_detector_pos(rp->refl, &fsh, &ssh);
- twod_mapping(fsh, ssh, &xh, &yh, rp->panel, dx, dy);
- return xh-xpk;
+ return fsh - rp->peak->fs;
}
-double y_dev(struct reflpeak *rp, struct detgeom *det,
- double dx, double dy)
+double ss_dev(struct reflpeak *rp, struct detgeom *det)
{
- /* Peak position term */
- double xpk, ypk, xh, yh;
double fsh, ssh;
- twod_mapping(rp->peak->fs, rp->peak->ss, &xpk, &ypk, rp->panel, dx, dy);
get_detector_pos(rp->refl, &fsh, &ssh);
- twod_mapping(fsh, ssh, &xh, &yh, rp->panel, dx, dy);
- return yh-ypk;
+ return ssh - rp->peak->ss;
}
@@ -397,10 +374,10 @@ static int iterate(struct reflpeak *rps, int n, UnitCell *cell,
}
- /* Positional x terms */
+ /* Positional fs terms */
for ( k=0; k<num_params; k++ ) {
- gradients[k] = x_gradient(rv[k], rps[i].refl, cell,
- rps[i].panel);
+ gradients[k] = fs_gradient(rv[k], rps[i].refl, cell,
+ rps[i].panel);
}
for ( k=0; k<num_params; k++ ) {
@@ -422,17 +399,17 @@ static int iterate(struct reflpeak *rps, int n, UnitCell *cell,
}
- v_c = x_dev(&rps[i], image->detgeom, *total_x, *total_y);
+ v_c = fs_dev(&rps[i], image->detgeom, *total_x, *total_y);
v_c *= -gradients[k];
v_curr = gsl_vector_get(v, k);
gsl_vector_set(v, k, v_curr + v_c);
}
- /* Positional y terms */
+ /* Positional ss terms */
for ( k=0; k<num_params; k++ ) {
- gradients[k] = y_gradient(rv[k], rps[i].refl, cell,
- rps[i].panel);
+ gradients[k] = ss_gradient(rv[k], rps[i].refl, cell,
+ rps[i].panel);
}
for ( k=0; k<num_params; k++ ) {
@@ -454,7 +431,7 @@ static int iterate(struct reflpeak *rps, int n, UnitCell *cell,
}
- v_c = y_dev(&rps[i], image->detgeom, *total_x, *total_y);
+ v_c = ss_dev(&rps[i], image->detgeom, *total_x, *total_y);
v_c *= -gradients[k];
v_curr = gsl_vector_get(v, k);
gsl_vector_set(v, k, v_curr + v_c);
@@ -535,13 +512,13 @@ static double pred_residual(struct reflpeak *rps, int n, struct detgeom *det,
r = 0.0;
for ( i=0; i<n; i++ ) {
- r += pow(x_dev(&rps[i], det, dx, dy), 2.0);
+ r += pow(fs_dev(&rps[i], det, dx, dy), 2.0);
}
res += r;
r = 0.0;
for ( i=0; i<n; i++ ) {
- r += pow(y_dev(&rps[i], det, dx, dy), 2.0);
+ r += pow(ss_dev(&rps[i], det, dx, dy), 2.0);
}
res += r;
diff --git a/libcrystfel/src/predict-refine.h b/libcrystfel/src/predict-refine.h
index d1e5fbe7..e74e00a8 100644
--- a/libcrystfel/src/predict-refine.h
+++ b/libcrystfel/src/predict-refine.h
@@ -46,17 +46,12 @@ struct reflpeak {
* Prediction refinement: refinement of indexing solutions before integration.
*/
-extern void twod_mapping(double fs, double ss, double *px, double *py,
- struct detgeom_panel *p, double dx, double dy);
-
extern int refine_prediction(struct image *image, Crystal *cr, Mille *mille);
extern int refine_radius(Crystal *cr, struct image *image);
-extern double x_dev(struct reflpeak *rp, struct detgeom *det,
- double dx, double dy);
+extern double fs_dev(struct reflpeak *rp, struct detgeom *det);
-extern double y_dev(struct reflpeak *rp, struct detgeom *det,
- double dx, double dy);
+extern double ss_dev(struct reflpeak *rp, struct detgeom *det);
#endif /* PREDICT_REFINE_H */