Commit a162a01f authored by Thomas White's avatar Thomas White
Browse files

Improved framework for detector geometry

parent c85a8a30
......@@ -16,13 +16,25 @@
#ifndef DETECTOR_H
#define DETECTOR_H
struct image;
#include "image.h"
/* Position of central beam for upper and lower CCDs */
#define UPPER_CX (491.9)
#define UPPER_CY (440.7)
#define LOWER_CX (492.0)
#define LOWER_CY (779.7)
struct panel
{
int min_x; /* Smallest x value considered to be in this panel */
int max_x; /* Largest x value considered to be in this panel */
int min_y; /* ... and so on */
int max_y;
float cx; /* Location of centre */
float cy;
};
struct detector
{
struct panel *panels;
int n_panels;
};
extern void record_image(struct image *image, int do_water, int do_poisson,
int do_bloom);
......
......@@ -68,20 +68,25 @@ void get_ewald(struct image *image)
for ( x=0; x<image->width; x++ ) {
for ( y=0; y<image->height; y++ ) {
double rx, ry, r;
double rx = 0.0;
double ry = 0.0;
double r;
double twothetax, twothetay, twotheta;
double qx, qy, qz;
struct rvec q;
int p;
/* Calculate q vectors for Ewald sphere */
if ( y >= 512 ) {
/* Top CCD */
rx = ((double)x - UPPER_CX) / image->resolution;
ry = ((double)y - UPPER_CY) / image->resolution;
} else {
/* Bottom CCD */
rx = ((double)x - LOWER_CX) / image->resolution;
ry = ((double)y - LOWER_CY) / image->resolution;
for ( p=0; p<image->det.n_panels; p++ ) {
if ( (x >= image->det.panels[p].min_x)
&& (x <= image->det.panels[p].max_x)
&& (y >= image->det.panels[p].min_y)
&& (y <= image->det.panels[p].max_y) ) {
rx = ((double)x - image->det.panels[p].cx)
/ image->resolution;
ry = ((double)y - image->det.panels[p].cy)
/ image->resolution;
}
}
r = sqrt(pow(rx, 2.0) + pow(ry, 2.0));
......
......@@ -23,6 +23,7 @@
#include "utils.h"
#include "cell.h"
#include "detector.h"
/* How is the scaling of the image described? */
......@@ -79,6 +80,7 @@ struct image {
double *twotheta;
struct molecule *molecule;
UnitCell *indexed_cell;
struct detector det;
struct quaternion orientation;
......
......@@ -111,19 +111,30 @@ void index_pattern(struct image *image, IndexingMethod indm)
for ( i=0; i<image_feature_count(image->features); i++ ) {
struct imagefeature *f;
double rx = 0.0;
double ry = 0.0;
int p;
int found = 0;
f = image_get_feature(image->features, i);
if ( f == NULL ) continue;
if ( f->y >=512 ) {
/* Top half of CCD */
map_position(image, f->x-UPPER_CX, f->y-UPPER_CY,
&f->rx, &f->ry, &f->rz);
} else {
/* Lower half of CCD */
map_position(image, f->x-LOWER_CX, f->y-LOWER_CY,
&f->rx, &f->ry, &f->rz);
for ( p=0; p<image->det.n_panels; p++ ) {
if ( (f->x >= image->det.panels[p].min_x)
&& (f->x <= image->det.panels[p].max_x)
&& (f->y >= image->det.panels[p].min_y)
&& (f->y <= image->det.panels[p].max_y) ) {
rx = ((double)f->x - image->det.panels[p].cx);
ry = ((double)f->y - image->det.panels[p].cy);
found = 1;
}
}
if ( !found ) {
ERROR("No mapping found for %f,%f\n", f->x, f->y);
continue;
}
map_position(image, rx, ry, &f->rx, &f->ry, &f->rz);
}
write_drx(image);
......
......@@ -166,6 +166,24 @@ int main(int argc, char *argv[])
image.data = NULL;
image.indexed_cell = NULL;
/* Set up detector configuration */
image.det.n_panels = 2;
image.det.panels = malloc(2*sizeof(struct panel));
/* Upper panel */
image.det.panels[0].min_x = 0;
image.det.panels[0].max_x = 1023;
image.det.panels[0].min_y = 512;
image.det.panels[0].max_y = 1023;
image.det.panels[0].cx = 491.9;
image.det.panels[0].cy = 440.7;
/* Lower panel */
image.det.panels[1].min_x = 0;
image.det.panels[1].max_x = 1023;
image.det.panels[1].min_y = 0;
image.det.panels[1].max_y = 511;
image.det.panels[1].cx = 492.0;
image.det.panels[1].cy = 779.7;
STATUS("Processing '%s'\n", line);
n_images++;
......
......@@ -220,6 +220,24 @@ int main(int argc, char *argv[])
image.lambda = ph_en_to_lambda(eV_to_J(2.0e3)); /* Wavelength */
image.molecule = load_molecule();
/* Set up detector configuration */
image.det.n_panels = 2;
image.det.panels = malloc(2*sizeof(struct panel));
/* Upper panel */
image.det.panels[0].min_x = 0;
image.det.panels[0].max_x = 1023;
image.det.panels[0].min_y = 512;
image.det.panels[0].max_y = 1023;
image.det.panels[0].cx = 491.9;
image.det.panels[0].cy = 440.7;
/* Lower panel */
image.det.panels[1].min_x = 0;
image.det.panels[1].max_x = 1023;
image.det.panels[1].min_y = 0;
image.det.panels[1].max_y = 511;
image.det.panels[1].cx = 492.0;
image.det.panels[1].cy = 779.7;
/* Splurge a few useful numbers */
STATUS("Wavelength is %f nm\n", image.lambda/1.0e-9);
......
......@@ -432,7 +432,11 @@ void dump_peaks(struct image *image)
double q, rx, ry, rz;
int x, y;
double rcx = 0.0;
double rcy = 0.0;
int found = 0;
struct imagefeature *f;
int p;
f = image_get_feature(image->features, i);
if ( f == NULL ) continue;
......@@ -440,16 +444,24 @@ void dump_peaks(struct image *image)
x = f->x;
y = f->y;
if ( f->y >=512 ) {
/* Top half of CCD */
map_position(image, f->x-UPPER_CX, f->y-UPPER_CY,
&rx, &ry, &rz);
} else {
/* Lower half of CCD */
map_position(image, f->x-LOWER_CX, f->y-LOWER_CY,
&rx, &ry, &rz);
for ( p=0; p<image->det.n_panels; p++ ) {
if ( (x >= image->det.panels[p].min_x)
&& (x <= image->det.panels[p].max_x)
&& (y >= image->det.panels[p].min_y)
&& (y <= image->det.panels[p].max_y) ) {
rcx = ((double)x - image->det.panels[p].cx)
/ image->resolution;
rcy = ((double)y - image->det.panels[p].cy)
/ image->resolution;
found = 1;
}
}
if ( !found ) {
ERROR("No mapping found for %i,%i\n", x, y);
continue;
}
map_position(image, rcx, rcy, &rx, &ry, &rz);
q = modulus(rx, ry, rz);
printf("%i\t%i\t%f\t%i\n", x, y, q/1.0e9,
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment