From 3638ffdee226e248c351362fed7fdc213d4d39b0 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 8 Jul 2015 14:00:16 +0200 Subject: Remove unnecessary calculations from prediction refinement --- libcrystfel/src/predict-refine.c | 62 +++++++++++++++------------------------ tests/prediction_gradient_check.c | 38 +++++++++++++----------- 2 files changed, 45 insertions(+), 55 deletions(-) diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c index 6ef31805..6784aae2 100644 --- a/libcrystfel/src/predict-refine.c +++ b/libcrystfel/src/predict-refine.c @@ -349,23 +349,16 @@ void refine_radius(Crystal *cr, struct image *image) } -/* Returns d(xh-xpk)/dP, where P = any parameter */ -static double x_gradient(int param, struct reflpeak *rp, struct detector *det, - double lambda, UnitCell *cell) +/* Returns dx_h/dP, where P = any parameter */ +static double x_gradient(int param, Reflection *refl, UnitCell *cell, + struct panel *p, double lambda) { signed int h, k, l; - double xpk, ypk, xh, yh; - double fsh, ssh; double x, z, wn; double ax, ay, az, bx, by, bz, cx, cy, cz; - twod_mapping(rp->peak->fs, rp->peak->ss, &xpk, &ypk, rp->panel); - get_detector_pos(rp->refl, &fsh, &ssh); - twod_mapping(fsh, ssh, &xh, &yh, rp->panel); - get_indices(rp->refl, &h, &k, &l); - + get_indices(refl, &h, &k, &l); wn = 1.0 / lambda; - cell_get_reciprocal(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); x = h*ax + k*bx + l*cx; z = h*az + k*bz + l*cz; @@ -373,13 +366,13 @@ static double x_gradient(int param, struct reflpeak *rp, struct detector *det, switch ( param ) { case GPARAM_ASX : - return h * rp->panel->clen / (wn+z); + return h * p->clen / (wn+z); case GPARAM_BSX : - return k * rp->panel->clen / (wn+z); + return k * p->clen / (wn+z); case GPARAM_CSX : - return l * rp->panel->clen / (wn+z); + return l * p->clen / (wn+z); case GPARAM_ASY : return 0.0; @@ -391,13 +384,13 @@ static double x_gradient(int param, struct reflpeak *rp, struct detector *det, return 0.0; case GPARAM_ASZ : - return -h * x * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -h * x * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_BSZ : - return -k * x * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -k * x * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_CSZ : - return -l * x * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -l * x * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_DETX : return -1; @@ -415,23 +408,16 @@ static double x_gradient(int param, struct reflpeak *rp, struct detector *det, } -/* Returns d(yh-ypk)/dP, where P = any parameter */ -static double y_gradient(int param, struct reflpeak *rp, struct detector *det, - double lambda, UnitCell *cell) +/* Returns dy_h/dP, where P = any parameter */ +static double y_gradient(int param, Reflection *refl, UnitCell *cell, + struct panel *p, double lambda) { signed int h, k, l; - double xpk, ypk, xh, yh; - double fsh, ssh; double y, z, wn; double ax, ay, az, bx, by, bz, cx, cy, cz; - twod_mapping(rp->peak->fs, rp->peak->ss, &xpk, &ypk, rp->panel); - get_detector_pos(rp->refl, &fsh, &ssh); - twod_mapping(fsh, ssh, &xh, &yh, rp->panel); - get_indices(rp->refl, &h, &k, &l); - + get_indices(refl, &h, &k, &l); wn = 1.0 / lambda; - cell_get_reciprocal(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); y = h*ay + k*by + l*cy; z = h*az + k*bz + l*cz; @@ -448,22 +434,22 @@ static double y_gradient(int param, struct reflpeak *rp, struct detector *det, return 0.0; case GPARAM_ASY : - return h * rp->panel->clen / (wn+z); + return h * p->clen / (wn+z); case GPARAM_BSY : - return k * rp->panel->clen / (wn+z); + return k * p->clen / (wn+z); case GPARAM_CSY : - return l * rp->panel->clen / (wn+z); + return l * p->clen / (wn+z); case GPARAM_ASZ : - return -h * y * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -h * y * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_BSZ : - return -k * y * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -k * y * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_CSZ : - return -l * y * rp->panel->clen / (wn*wn + 2*wn*z + z*z); + return -l * y * p->clen / (wn*wn + 2*wn*z + z*z); case GPARAM_DETX : return 0; @@ -553,8 +539,8 @@ static int iterate(struct reflpeak *rps, int n, UnitCell *cell, /* Positional x terms */ for ( k=0; kdet, - image->lambda, cell); + gradients[k] = x_gradient(rv[k], rps[i].refl, cell, + rps[i].panel, image->lambda); } for ( k=0; kdet, - image->lambda, cell); + gradients[k] = y_gradient(rv[k], rps[i].refl, cell, + rps[i].panel, image->lambda); } for ( k=0; k #include #include -#include "../src/predict-refine.c" int checkrxy; +static void twod_mapping(double fs, double ss, double *px, double *py, + struct panel *p) +{ + double xs, ys; + + xs = fs*p->fsx + ss*p->ssx; + ys = fs*p->fsy + ss*p->ssy; + + *px = (xs + p->cnx) / p->res; + *py = (ys + p->cny) / p->res; +} + + static void scan(RefList *reflections, RefList *compare, int *valid, long double *vals[3], int idx) { @@ -260,7 +272,6 @@ static double test_gradients(Crystal *cr, double incr_val, int refine, long double grad1, grad2, grad; double cgrad; signed int h, k, l; - struct reflpeak rp; get_indices(refl, &h, &k, &l); @@ -283,27 +294,20 @@ static double test_gradients(Crystal *cr, double incr_val, int refine, } else { - struct imagefeature pk; struct image *image; - pk.fs = 0.0; - pk.ss = 0.0; - image = crystal_get_image(cr); - rp.refl = refl; - rp.peak = &pk; - rp.panel = &image->det->panels[0]; if ( checkrxy == 1 ) { - cgrad = x_gradient(refine, &rp, - crystal_get_image(cr)->det, - crystal_get_image(cr)->lambda, - crystal_get_cell(cr)); + cgrad = x_gradient(refine, refl, + crystal_get_cell(cr), + &image->det->panels[0], + crystal_get_image(cr)->lambda); } else { - cgrad = y_gradient(refine, &rp, - crystal_get_image(cr)->det, - crystal_get_image(cr)->lambda, - crystal_get_cell(cr)); + cgrad = y_gradient(refine, refl, + crystal_get_cell(cr), + &image->det->panels[0], + crystal_get_image(cr)->lambda); } } -- cgit v1.2.3