From e3a4cc7fcaa71c1954394fc4bdee90e3bcc38502 Mon Sep 17 00:00:00 2001 From: taw27 Date: Sat, 8 Sep 2007 21:18:28 +0000 Subject: Fixed eFOM git-svn-id: svn://cook.msm.cam.ac.uk:745/diff-tomo/dtr@128 bf6ca9ba-c028-0410-8290-897cf20841d1 --- src/ipr.c | 30 ++++++++++++++++++++---------- 1 file changed, 20 insertions(+), 10 deletions(-) diff --git a/src/ipr.c b/src/ipr.c index 982d8c3..aa55bd5 100644 --- a/src/ipr.c +++ b/src/ipr.c @@ -219,23 +219,33 @@ static double ipr_efom(ReflectionContext *rctx, Basis *basis) { if ( cur->type == REFLECTION_NORMAL ) { /* Can this basis approximately account for this reflection? */ - double n; - double p, q, r, u, v, w, d, e, f; + double det; + double a11, a12, a13, a21, a22, a23, a31, a32, a33; double h, k, l; /* Set up the coordinate transform from hkl to xyz */ - p = basis->a.x; q = basis->a.y; r = basis->a.z; - u = basis->b.x; v = basis->b.y; w = basis->b.z; - d = basis->c.x; e = basis->c.y; f = basis->c.z; + a11 = basis->a.x; a12 = basis->a.y; a13 = basis->a.z; + a21 = basis->b.x; a22 = basis->b.y; a23 = basis->b.z; + a31 = basis->c.x; a32 = basis->c.y; a33 = basis->c.z; /* Invert the matrix to get hkl from xyz */ - n = p*(v*f - w*e) - q*(u*f - w*d) + r*(u*e - d*v); - h = (p*(v*f-w*e)*cur->x + -q*(u*f-w*d)*cur->y + r*(u*e-v*d)*cur->z) / n; - k = (-u*(q*f-r*e)*cur->x + v*(p*f-r*d)*cur->y + -w*(p*e-q*d)*cur->z) / n; - l = (d*(q*w-r*v)*cur->x + -e*(p*w-r*u)*cur->y + f*(p*v-q*u)*cur->z) / n; + det = a11*(a22*a33 - a23*a32) - a12*(a21*a33 - a23*a31) + a13*(a21*a32 - a22*a31); + h = ((a22*a33-a23*a32)*cur->x + (a23*a31-a21*a33)*cur->y + (a21*a32-a22*a31)*cur->z)/ det; + k = ((a13*a32-a12*a33)*cur->x + (a11*a33-a13*a31)*cur->y + (a12*a31-a11*a32)*cur->z) / det; + l = ((a12*a23-a13*a22)*cur->x + (a13*a21-a11*a23)*cur->y + (a11*a22-a12*a21)*cur->z) / det; + + //printf("%12.0f %12.0f %12.0f ----> %3.2f %3.2f %3.2f", cur->x, cur->y, cur->z, h, k, l); + //printf("%f %f %f %f %f %f\n", cur->x, h*basis->a.x + k*basis->b.x + l*basis->c.x, + // cur->y, h*basis->a.y + k*basis->b.y + l*basis->c.y, + // cur->z, h*basis->a.z + k*basis->b.z + l*basis->c.z); /* Calculate the deviations in terms of |a|, |b| and |c| */ + h = fabs(h); k = fabs(k); l = fabs(l); h -= floor(h); k -= floor(k); l -= floor(l); + if ( h == 1.0 ) h = 0.0; + if ( k == 1.0 ) k = 0.0; + if ( l == 1.0 ) l = 0.0; + //printf(" ----> %3.2f %3.2f %3.2f\n", h, k, l); if ( h*h + k*k + l*l <= 0.1*0.1*0.1 ) n_indexed++; n_counted++; @@ -274,7 +284,7 @@ static Basis *ipr_choose_initial_basis(ControlContext *ctx) { if ( best_basis ) free(best_basis); best_efom = efom; best_basis = basis; - printf("IP: %3i: eFOM=%f\n", i, efom); + printf("IP: %3i: eFOM = %f %%\n", i, efom*100); } -- cgit v1.2.3