aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--src/hdf5-file.c3
-rw-r--r--src/image.h13
-rw-r--r--src/index.c25
-rw-r--r--src/pattern_sim.c1
-rw-r--r--src/relrod.c202
-rw-r--r--src/relrod.h24
6 files changed, 5 insertions, 263 deletions
diff --git a/src/hdf5-file.c b/src/hdf5-file.c
index 0a63a853..3fc5de8c 100644
--- a/src/hdf5-file.c
+++ b/src/hdf5-file.c
@@ -192,9 +192,6 @@ int hdf5_read(struct hdfile *f, struct image *image)
image->height = f->nx;
image->width = f->ny;
- /* Always camera length/lambda formulation for FEL */
- image->fmode = FORMULATION_CLEN;
-
/* Read wavelength from file */
image->lambda = get_wavelength(f);
if ( image->lambda < 0.0 ) {
diff --git a/src/image.h b/src/image.h
index 543d924c..326e70b7 100644
--- a/src/image.h
+++ b/src/image.h
@@ -26,13 +26,6 @@
#include "detector.h"
-/* How is the scaling of the image described? */
-typedef enum {
- FORMULATION_CLEN,
- FORMULATION_PIXELSIZE
-} FormulationMode;
-
-
/* Structure describing a feature in an image */
struct imagefeature {
@@ -83,12 +76,6 @@ struct image {
struct quaternion orientation;
- /* Image parameters can be given as camera length or pixel size.
- * If FORMULATION_CLEN, then camera_len and resolution must be given.
- * If FORMULATION_PIXELSIZE, then pixel_size only is needed.*/
- FormulationMode fmode;
- double pixel_size;
-
/* Wavelength must always be given */
double lambda; /* Wavelength in m */
diff --git a/src/index.c b/src/index.c
index e9f5cf44..4e8460f6 100644
--- a/src/index.c
+++ b/src/index.c
@@ -57,26 +57,11 @@ int map_position(struct image *image, double dx, double dy,
return 0;
}
- if ( image->fmode == FORMULATION_CLEN ) {
-
- /* Convert pixels to metres */
- x /= image->det.panels[p].res;
- y /= image->det.panels[p].res; /* Convert pixels to metres */
- d = sqrt((x*x) + (y*y));
- twotheta = atan2(d, image->det.panels[p].clen);
-
- } else if (image->fmode == FORMULATION_PIXELSIZE ) {
-
- /* Convert pixels to metres^-1 */
- x *= image->pixel_size;
- y *= image->pixel_size; /* Convert pixels to metres^-1 */
- d = sqrt((x*x) + (y*y));
- twotheta = atan2(d, k);
-
- } else {
- ERROR("Unrecognised formulation mode in mapping_scale.\n");
- return -1;
- }
+ /* Convert pixels to metres */
+ x /= image->det.panels[p].res;
+ y /= image->det.panels[p].res; /* Convert pixels to metres */
+ d = sqrt((x*x) + (y*y));
+ twotheta = atan2(d, image->det.panels[p].clen);
psi = atan2(y, x);
diff --git a/src/pattern_sim.c b/src/pattern_sim.c
index e389385d..340acd8b 100644
--- a/src/pattern_sim.c
+++ b/src/pattern_sim.c
@@ -223,7 +223,6 @@ int main(int argc, char *argv[])
/* Define image parameters */
image.width = 1024;
image.height = 1024;
- image.fmode = FORMULATION_CLEN;
image.lambda = ph_en_to_lambda(eV_to_J(2.0e3)); /* Wavelength */
image.molecule = load_molecule();
diff --git a/src/relrod.c b/src/relrod.c
deleted file mode 100644
index 4c850310..00000000
--- a/src/relrod.c
+++ /dev/null
@@ -1,202 +0,0 @@
-/*
- * relrod.c
- *
- * Calculate reflection positions via line-sphere intersection test
- *
- * (c) 2006-2010 Thomas White <taw@physics.org>
- *
- * Part of CrystFEL - crystallography with a FEL
- *
- */
-
-
-#include <stdlib.h>
-#include <math.h>
-#include <stdio.h>
-
-#include "image.h"
-#include "utils.h"
-#include "cell.h"
-
-
-static void mapping_rotate(double x, double y, double z,
- double *ddx, double *ddy, double *ddz,
- double omega, double tilt)
-{
- double nx, ny, nz;
- double x_temp, y_temp, z_temp;
-
- /* First: rotate image clockwise until tilt axis is aligned
- * horizontally. */
- nx = x*cos(omega) + y*sin(omega);
- ny = -x*sin(omega) + y*cos(omega);
- nz = z;
-
- /* Now, tilt about the x-axis ANTICLOCKWISE around +x, i.e. the
- * "wrong" way. This is because the crystal is rotated in the
- * experiment, not the Ewald sphere. */
- x_temp = nx; y_temp = ny; z_temp = nz;
- nx = x_temp;
- ny = cos(tilt)*y_temp + sin(tilt)*z_temp;
- nz = -sin(tilt)*y_temp + cos(tilt)*z_temp;
-
- /* Finally, reverse the omega rotation to restore the location of the
- * image in 3D space */
- x_temp = nx; y_temp = ny; z_temp = nz;
- nx = x_temp*cos(-omega) + y_temp*sin(-omega);
- ny = -x_temp*sin(-omega) + y_temp*cos(-omega);
- nz = z_temp;
-
- *ddx = nx;
- *ddy = ny;
- *ddz = nz;
-}
-
-
-void get_reflections(struct image *image, UnitCell *cell, double smax)
-{
- ImageFeatureList *flist;
- double tilt, omega, wavenumber;
- double nx, ny, nz; /* "normal" vector */
- double kx, ky, kz; /* Electron wavevector ("normal" times 1/lambda) */
- double ux, uy, uz; /* "up" vector */
- double rx, ry, rz; /* "right" vector */
- double asx, asy, asz; /* "a*" lattice parameter */
- double bsx, bsy, bsz; /* "b*" lattice parameter */
- double csx, csy, csz; /* "c*" lattice parameter */
- signed int h, k, l;
- double res_max;
-
- /* Get the reciprocal unit cell */
- cell_get_reciprocal(cell, &asx, &asy, &asz,
- &bsx, &bsy, &bsz,
- &csx, &csy, &csz);
-
- /* Prepare list and some parameters */
- flist = image_feature_list_new();
- tilt = image->tilt;
- omega = image->omega;
- wavenumber = 1.0/image->lambda;
-
- /* Calculate (roughly) the maximum resolution */
- if ( image->fmode == FORMULATION_CLEN ) {
- double w2, h2;
- w2 = image->width/2; h2 = image->height/2;
- w2 = pow(w2, 2.0); h2 = pow(h2, 2.0);
- res_max = sqrt(w2 + h2) / image->resolution;
- res_max *= (wavenumber / image->camera_len);
- } else {
- ERROR("Unrecognised formulation mode in get_reflections"
- " (resolution cutoff calculation)\n");
- return;
- }
- res_max = pow(res_max, 2.0);
-
- /* Calculate the (normalised) incident electron wavevector */
- mapping_rotate(0.0, 0.0, 1.0, &nx, &ny, &nz, omega, tilt);
- kx = nx / image->lambda;
- ky = ny / image->lambda;
- kz = nz / image->lambda; /* This is the centre of the Ewald sphere */
-
- /* Determine where "up" is */
- mapping_rotate(0.0, 1.0, 0.0, &ux, &uy, &uz, omega, tilt);
-
- /* Determine where "right" is */
- mapping_rotate(1.0, 0.0, 0.0, &rx, &ry, &rz, omega, tilt);
-
- for ( h=-50; h<50; h++ ) {
- for ( k=-50; k<50; k++ ) {
- for ( l=-50; l<50; l++ ) {
-
- double xl, yl, zl;
- double a, b, c;
- double s1, s2, s, t;
- double g_sq, gn;
-
- /* Get the coordinates of the reciprocal lattice point */
- xl = h*asx + k*bsx + l*csx;
- yl = h*asy + k*bsy + l*csy;
- zl = h*asz + k*bsz + l*csz;
- g_sq = modulus_squared(xl, yl, zl);
- gn = xl*nx + yl*ny + zl*nz;
-
- /* Early bailout if resolution is clearly too high */
- if ( g_sq > res_max ) continue;
-
- /* Next, solve the relrod equation to calculate
- * the excitation error */
- a = 1.0;
- b = 2.0*(wavenumber + gn);
- c = -2.0*gn*wavenumber + g_sq;
- t = -0.5*(b + sign(b)*sqrt(b*b - 4.0*a*c));
- s1 = t/a;
- s2 = c/t;
- if ( fabs(s1) < fabs(s2) ) s = s1; else s = s2;
-
- /* Skip this reflection if s is large */
- if ( fabs(s) <= smax ) {
-
- double xi, yi, zi;
- double gx, gy, gz;
- double theta;
- double x, y;
- double dx, dy, psi;
-
- /* Determine the intersection point */
- xi = xl + s*nx; yi = yl + s*ny; zi = zl + s*nz;
-
- /* Calculate Bragg angle */
- gx = xi - kx;
- gy = yi - ky;
- gz = zi - kz; /* This is the vector from the centre of
- * the sphere to the intersection */
- theta = angle_between(-kx, -ky, -kz, gx, gy, gz);
-
- /* Calculate azimuth of point in image
- * (anticlockwise from +x) */
- dx = xi*rx + yi*ry + zi*rz;
- dy = xi*ux + yi*uy + zi*uz;
- psi = atan2(dy, dx);
-
- /* Get image coordinates from polar
- * representation */
- if ( image->fmode == FORMULATION_CLEN ) {
- x = image->camera_len*tan(theta)*cos(psi);
- y = image->camera_len*tan(theta)*sin(psi);
- x *= image->resolution;
- y *= image->resolution;
- } else if ( image->fmode==FORMULATION_PIXELSIZE ) {
- x = tan(theta)*cos(psi) / image->lambda;
- y = tan(theta)*sin(psi) / image->lambda;
- x /= image->pixel_size;
- y /= image->pixel_size;
- } else {
- ERROR("Unrecognised formulation mode "
- "in get_reflections\n");
- return;
- }
-
- x += image->x_centre;
- y += image->y_centre;
-
- /* Sanity check */
- if ( (x>=0) && (x<image->width)
- && (y>=0) && (y<image->height) ) {
-
- /* Record the reflection.
- * Intensity should be multiplied by relrod
- * spike function, except reprojected
- * reflections aren't used quantitatively for
- * anything. */
- image_add_feature(flist, x, y, image, 1.0);
-
- } /* else it's outside the picture somewhere */
-
- } /* else reflection is not excited in this orientation */
-
- }
- }
- }
-
- image->rflist = flist;
-}
diff --git a/src/relrod.h b/src/relrod.h
deleted file mode 100644
index 9e664392..00000000
--- a/src/relrod.h
+++ /dev/null
@@ -1,24 +0,0 @@
-/*
- * relrod.h
- *
- * Calculate reflection positions via line-sphere intersection test
- *
- * (c) 2006-2010 Thomas White <taw@physics.org>
- *
- * Part of CrystFEL - crystallography with a FEL
- *
- */
-
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
-
-#ifndef RELROD_H
-#define RELROD_H
-
-#include "image.h"
-#include "cell.h"
-
-extern void get_reflections(struct image *image, UnitCell *cell, double smax);
-
-#endif /* RELROD_H */
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-15 16:45:57 +0000