diff options
Diffstat (limited to 'src/reproject.c')
-rw-r--r-- | src/reproject.c | 46 |
1 files changed, 19 insertions, 27 deletions
diff --git a/src/reproject.c b/src/reproject.c index efd32c9..ffd49b9 100644 --- a/src/reproject.c +++ b/src/reproject.c @@ -25,7 +25,10 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect int i; double smax = 0.5e9; double tilt, omega; - + double xt, yt, zt; + double nx, ny, nz; + double kx, ky, kz; + tilt = deg2rad(image.tilt); omega = deg2rad(image.omega); @@ -34,11 +37,21 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect i = 0; reflection = rctx->reflections; + + /* Calculate the (normalised) incident electron wavevector */ + xt = 0; yt = sin(tilt); zt = cos(tilt); + nx = xt*cos(omega) + yt*-sin(omega); + ny = xt*sin(omega) + yt*cos(omega); + nz = zt; + kx = nx / image.lambda; + ky = ny / image.lambda; + kz = nz / image.lambda; /* This is the centre of the Ewald sphere */ + + reflection_add(ctx->reflectionctx, kx, ky, kz, 1, REFLECTION_VECTOR_MARKER_1); + do { double xl, yl, zl; - double nx, ny, nz; - double xt, yt, zt; double a, b, c; double A1, A2, s1, s2, s; @@ -47,14 +60,6 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect yl = reflection->y; zl = reflection->z; - /* Now calculate the (normalised) incident electron wavevector */ - xt = 0; - yt = sin(tilt); - zt = cos(tilt); - nx = xt*cos(omega) + yt*-sin(omega); - ny = xt*sin(omega) + yt*cos(omega); - nz = zt; - /* Next, solve the relrod equation to calculate the excitation error */ a = 1.0; b = -2.0*(xl*nx + yl*ny + zl*nz); @@ -71,7 +76,6 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect double xi, yi, zi; double gx, gy, gz; - double kx, ky, kz; double cx, cy, cz; double ux, uy, uz; double uyt, uzt; @@ -83,9 +87,6 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect xi = xl + s*nx; yi = yl + s*ny; zi = zl + s*nz; /* Calculate Bragg angle */ - kx = nx / image.lambda; - ky = ny / image.lambda; - kz = nz / image.lambda; /* This is the centre of the Ewald sphere */ gx = xi - kx; gy = yi - ky; gz = zi - kz; /* This is the vector from the centre of the sphere to the intersection */ @@ -95,6 +96,9 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect double psi, disc; + reflection_add(ctx->reflectionctx, xl, yl, zl, 1, REFLECTION_GENERATED); + reflection_add(ctx->reflectionctx, xi, yi, zi, 1, REFLECTION_MARKER); + /* Determine where "up" is */ uyt = cos(tilt); uzt = -sin(tilt); /* tilt it (uxt not needed) */ ux = uyt*-sin(omega); uy = uyt*cos(omega); uz = uzt; /* rotate it */ @@ -111,18 +115,6 @@ ImageReflection *reproject_get_reflections(ImageRecord image, size_t *n, Reflect else if ( (psi < M_PI_2) && (disc < M_PI_2) ) psi = 3*M_PI_2 - psi; /* Case 3 */ else if ( (psi < M_PI_2) && (disc >= M_PI_2) ) psi = 3*M_PI_2 + psi; /* Case 4 */ psi = M_PI_2 - psi; /* Anticlockwise from "+x" instead of clockwise from "up" */ - if ( (reflection->h == 2) && (reflection->k == 2) && (reflection->l == 0) ) { - printf("c = %f, %f, %f\n", cx, cy, cz); - reflection_add(ctx->reflectionctx, cx, cy, cz, 1, REFLECTION_VECTOR_MARKER_2); - printf("r = %f, %f, %f\n", rx, ry, rz); - reflection_add(ctx->reflectionctx, rx, ry, rz, 1, REFLECTION_VECTOR_MARKER_3); - reflection->type = REFLECTION_MARKER; - printf("n = %f, %f, %f\n", kx, ky, kz); - reflection_add(ctx->reflectionctx, kx, ky, kz, 1, REFLECTION_VECTOR_MARKER_1); - } else { - reflection_add(ctx->reflectionctx, xl, yl, zl, 1, REFLECTION_GENERATED); - reflection_add(ctx->reflectionctx, xi, yi, zi, 1, REFLECTION_MARKER); - } /* Calculate image coordinates from polar representation */ if ( image.fmode == FORMULATION_CLEN ) { |