From aedd1501af8707b9f502a2fd765d743f53c3d078 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Thu, 7 Mar 2013 13:23:41 +0100 Subject: pr_gradient_check: Fix and improve --- tests/pr_gradient_check.c | 72 +++++++++++++++++++++++++++++------------------ 1 file changed, 44 insertions(+), 28 deletions(-) diff --git a/tests/pr_gradient_check.c b/tests/pr_gradient_check.c index c7eb0129..2631b1bb 100644 --- a/tests/pr_gradient_check.c +++ b/tests/pr_gradient_check.c @@ -155,6 +155,7 @@ static Crystal *new_shifted_crystal(Crystal *cr, int refine, double incr_val) case REF_R : cell = cell_new_from_cell(crystal_get_cell(cr)); + crystal_set_cell(cr_new, cell); crystal_set_profile_radius(cr_new, r + incr_val); break; @@ -218,8 +219,8 @@ static void calc_either_side(Crystal *cr, double incr_val, -static int test_gradients(Crystal *cr, double incr_val, int refine, - const char *str, PartialityModel pmodel) +static double test_gradients(Crystal *cr, double incr_val, int refine, + const char *str, PartialityModel pmodel) { Reflection *refl; RefListIterator *iter; @@ -227,9 +228,11 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, int i; int *valid; int nref; - int n_acc, n_valid; + int n_good, n_invalid, n_small, n_nan, n_bad; RefList *reflections; //FILE *fh; + int ntot = 0; + double total = 0.0; reflections = find_intersections(crystal_get_image(cr), cr); crystal_set_reflections(cr, reflections); @@ -261,7 +264,8 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, //fh = fopen("wrongness.dat", "a"); - n_valid = nref; n_acc = 0; + n_invalid = 0; n_good = 0; + n_nan = 0; n_small = 0; n_bad = 0; i = 0; for ( refl = first_refl(reflections, &iter); refl != NULL; @@ -275,7 +279,7 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, get_indices(refl, &h, &k, &l); if ( !valid[i] ) { - n_valid--; + n_invalid++; } else { double r1, r2, p; @@ -289,8 +293,20 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, get_partial(refl, &r1, &r2, &p, &cl, &ch); - if ( (fabs(cgrad) > 5e-8) && - !within_tolerance(grad, cgrad, 10.0) ) + if ( fabs(cgrad) < 5e-8 ) { + n_small++; + continue; + } + + if ( isnan(cgrad) ) { + n_nan++; + continue; + } + + total += fabs(cgrad - grad); + ntot++; + + if ( !within_tolerance(grad, cgrad, 5.0) ) { STATUS("!- %s %3i %3i %3i" @@ -298,6 +314,7 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, " %10.2e %10.2e\n", str, h, k, l, grad, cgrad, cgrad/grad, r1, r2); + n_bad++; } else { @@ -307,7 +324,7 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, // str, h, k, l, grad, cgrad, cgrad/grad, // r1, r2); - n_acc++; + n_good++; } @@ -323,13 +340,12 @@ static int test_gradients(Crystal *cr, double incr_val, int refine, } - STATUS("%s: %i out of %i valid gradients were accurate.\n", - str, n_acc, n_valid); + STATUS("%3s: %3i within 5%%, %3i outside, %3i nan, %3i invalid, " + "%3i small. ", str, n_good, n_bad, n_nan, n_invalid, n_small); //fclose(fh); - if ( n_acc != n_valid ) return 1; - - return 0; + STATUS("Fractional error = %f %%\n", 100.0*total/ntot); + return total / ntot; } @@ -345,8 +361,8 @@ int main(int argc, char *argv[]) Crystal *cr; struct quaternion orientation; int i; - int val; const PartialityModel pmodel = PMODEL_SPHERE; + double tot = 0.0; image.width = 1024; image.height = 1024; @@ -374,8 +390,6 @@ int main(int argc, char *argv[]) deg2rad(90.0), deg2rad(90.0)); - val = 0; - for ( i=0; i<1; i++ ) { UnitCell *rot; @@ -389,32 +403,34 @@ int main(int argc, char *argv[]) &bz, &cx, &cy, &cz); incr_val = incr_frac * image.div; - val += test_gradients(cr, incr_val, REF_DIV, "div", pmodel); + tot += test_gradients(cr, incr_val, REF_DIV, "div", pmodel); + + incr_val = incr_frac * crystal_get_profile_radius(cr); + tot += test_gradients(cr, incr_val, REF_R, "R", pmodel); incr_val = incr_frac * ax; - val += test_gradients(cr, incr_val, REF_ASX, "ax*", pmodel); + tot += test_gradients(cr, incr_val, REF_ASX, "ax*", pmodel); incr_val = incr_frac * ay; - val += test_gradients(cr, incr_val, REF_ASY, "ay*", pmodel); + tot += test_gradients(cr, incr_val, REF_ASY, "ay*", pmodel); incr_val = incr_frac * az; - val += test_gradients(cr, incr_val, REF_ASZ, "az*", pmodel); + tot += test_gradients(cr, incr_val, REF_ASZ, "az*", pmodel); incr_val = incr_frac * bx; - val += test_gradients(cr, incr_val, REF_BSX, "bx*", pmodel); + tot += test_gradients(cr, incr_val, REF_BSX, "bx*", pmodel); incr_val = incr_frac * by; - val += test_gradients(cr, incr_val, REF_BSY, "by*", pmodel); + tot += test_gradients(cr, incr_val, REF_BSY, "by*", pmodel); incr_val = incr_frac * bz; - val += test_gradients(cr, incr_val, REF_BSZ, "bz*", pmodel); + tot += test_gradients(cr, incr_val, REF_BSZ, "bz*", pmodel); incr_val = incr_frac * cx; - val += test_gradients(cr, incr_val, REF_CSX, "cx*", pmodel); + tot += test_gradients(cr, incr_val, REF_CSX, "cx*", pmodel); incr_val = incr_frac * cy; - val += test_gradients(cr, incr_val, REF_CSY, "cy*", pmodel); + tot += test_gradients(cr, incr_val, REF_CSY, "cy*", pmodel); incr_val = incr_frac * cz; - val += test_gradients(cr, incr_val, REF_CSZ, "cz*", pmodel); + tot += test_gradients(cr, incr_val, REF_CSZ, "cz*", pmodel); } - STATUS("Returning 0 by default: gradients only needed for experimental" - " features of CrystFEL.\n"); + STATUS("Mean fractional error = %f %%\n", 100.0*tot/10.0); return 0; } -- cgit v1.2.3