aboutsummaryrefslogtreecommitdiff
path: root/tests
diff options
context:
space:
mode:
Diffstat (limited to 'tests')
-rw-r--r--tests/pr_gradient_check.c103
1 files changed, 64 insertions, 39 deletions
diff --git a/tests/pr_gradient_check.c b/tests/pr_gradient_check.c
index 2631b1bb..cf11d18c 100644
--- a/tests/pr_gradient_check.c
+++ b/tests/pr_gradient_check.c
@@ -29,6 +29,7 @@
#include <stdlib.h>
#include <stdio.h>
+#include <gsl/gsl_statistics.h>
#include <image.h>
#include <cell.h>
@@ -230,9 +231,14 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
int nref;
int n_good, n_invalid, n_small, n_nan, n_bad;
RefList *reflections;
- //FILE *fh;
+ FILE *fh;
int ntot = 0;
double total = 0.0;
+ char tmp[32];
+ double *vec1;
+ double *vec2;
+ int n_line;
+ double cc;
reflections = find_intersections(crystal_get_image(cr), cr);
crystal_set_reflections(cr, reflections);
@@ -240,7 +246,7 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
nref = num_reflections(reflections);
if ( nref < 10 ) {
ERROR("Too few reflections found. Failing test by default.\n");
- return -1;
+ return 0.0;
}
vals[0] = malloc(nref*sizeof(long double));
@@ -248,13 +254,13 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
vals[2] = malloc(nref*sizeof(long double));
if ( (vals[0] == NULL) || (vals[1] == NULL) || (vals[2] == NULL) ) {
ERROR("Couldn't allocate memory.\n");
- return -1;
+ return 0.0;
}
valid = malloc(nref*sizeof(int));
if ( valid == NULL ) {
ERROR("Couldn't allocate memory.\n");
- return -1;
+ return 0.0;
}
for ( i=0; i<nref; i++ ) valid[i] = 1;
@@ -262,10 +268,18 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
calc_either_side(cr, incr_val, valid, vals, refine);
- //fh = fopen("wrongness.dat", "a");
+ snprintf(tmp, 32, "gradient-test-%s.dat", str);
+ fh = fopen(tmp, "w");
+
+ vec1 = malloc(nref*sizeof(double));
+ vec2 = malloc(nref*sizeof(double));
+ if ( (vec1 == NULL) || (vec2 == NULL) ) {
+ ERROR("Couldn't allocate memory.\n");
+ return 0.0;
+ }
n_invalid = 0; n_good = 0;
- n_nan = 0; n_small = 0; n_bad = 0;
+ n_nan = 0; n_small = 0; n_bad = 0; n_line = 0;
i = 0;
for ( refl = first_refl(reflections, &iter);
refl != NULL;
@@ -293,27 +307,32 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
get_partial(refl, &r1, &r2, &p, &cl, &ch);
- if ( fabs(cgrad) < 5e-8 ) {
- n_small++;
+ if ( isnan(cgrad) ) {
+ n_nan++;
continue;
}
- if ( isnan(cgrad) ) {
- n_nan++;
+ fprintf(fh, "%e %Le\n", cgrad, grad);
+ vec1[n_line] = cgrad;
+ vec2[n_line] = grad;
+ n_line++;
+
+ if ( fabs(cgrad) < 5e-8 ) {
+ n_small++;
continue;
}
total += fabs(cgrad - grad);
ntot++;
- if ( !within_tolerance(grad, cgrad, 5.0) )
+ if ( !within_tolerance(grad, cgrad, 10.0) )
{
- STATUS("!- %s %3i %3i %3i"
- " %10.2Le %10.2e ratio = %5.2Lf"
- " %10.2e %10.2e\n",
- str, h, k, l, grad, cgrad, cgrad/grad,
- r1, r2);
+ //STATUS("!- %s %3i %3i %3i"
+ // " %10.2Le %10.2e ratio = %5.2Lf"
+ // " %10.2e %10.2e\n",
+ // str, h, k, l, grad, cgrad, cgrad/grad,
+ // r1, r2);
n_bad++;
} else {
@@ -328,12 +347,6 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
}
- //fprintf(fh, "%e %f\n",
- //resolution(image->indexed_cell, h, k, l),
- //rad2deg(tt),
- // cgrad,
- // fabs((grad-cgrad)/grad));
-
}
i++;
@@ -342,10 +355,11 @@ static double test_gradients(Crystal *cr, double incr_val, int refine,
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);
+ fclose(fh);
- STATUS("Fractional error = %f %%\n", 100.0*total/ntot);
- return total / ntot;
+ cc = gsl_stats_correlation(vec1, 1, vec2, 1, n_line);
+ STATUS("CC = %+f\n", cc);
+ return cc;
}
@@ -362,7 +376,7 @@ int main(int argc, char *argv[])
struct quaternion orientation;
int i;
const PartialityModel pmodel = PMODEL_SPHERE;
- double tot = 0.0;
+ int fail = 0;
image.width = 1024;
image.height = 1024;
@@ -393,6 +407,7 @@ int main(int argc, char *argv[])
for ( i=0; i<1; i++ ) {
UnitCell *rot;
+ double val;
orientation = random_quaternion();
rot = cell_rotate(cell, orientation);
@@ -403,34 +418,44 @@ int main(int argc, char *argv[])
&bz, &cx, &cy, &cz);
incr_val = incr_frac * image.div;
- tot += test_gradients(cr, incr_val, REF_DIV, "div", pmodel);
+ val = test_gradients(cr, incr_val, REF_DIV, "div", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * crystal_get_profile_radius(cr);
- tot += test_gradients(cr, incr_val, REF_R, "R", pmodel);
+ val = test_gradients(cr, incr_val, REF_R, "R", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * ax;
- tot += test_gradients(cr, incr_val, REF_ASX, "ax*", pmodel);
+ val = test_gradients(cr, incr_val, REF_ASX, "ax*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * ay;
- tot += test_gradients(cr, incr_val, REF_ASY, "ay*", pmodel);
+ val = test_gradients(cr, incr_val, REF_ASY, "ay*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * az;
- tot += test_gradients(cr, incr_val, REF_ASZ, "az*", pmodel);
+ val = test_gradients(cr, incr_val, REF_ASZ, "az*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * bx;
- tot += test_gradients(cr, incr_val, REF_BSX, "bx*", pmodel);
+ val = test_gradients(cr, incr_val, REF_BSX, "bx*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * by;
- tot += test_gradients(cr, incr_val, REF_BSY, "by*", pmodel);
+ val = test_gradients(cr, incr_val, REF_BSY, "by*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * bz;
- tot += test_gradients(cr, incr_val, REF_BSZ, "bz*", pmodel);
+ val = test_gradients(cr, incr_val, REF_BSZ, "bz*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * cx;
- tot += test_gradients(cr, incr_val, REF_CSX, "cx*", pmodel);
+ val = test_gradients(cr, incr_val, REF_CSX, "cx*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * cy;
- tot += test_gradients(cr, incr_val, REF_CSY, "cy*", pmodel);
+ val = test_gradients(cr, incr_val, REF_CSY, "cy*", pmodel);
+ if ( val > 0.1 ) fail = 1;
incr_val = incr_frac * cz;
- tot += test_gradients(cr, incr_val, REF_CSZ, "cz*", pmodel);
+ val = test_gradients(cr, incr_val, REF_CSZ, "cz*", pmodel);
+ if ( val > 0.1 ) fail = 1;
}
- STATUS("Mean fractional error = %f %%\n", 100.0*tot/10.0);
- return 0;
+ return fail;
}
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-11 08:28:15 +0000