diff options
| -rw-r--r-- | src/detector.c | 6 | ||||
| -rw-r--r-- | src/diffraction.c | 46 | ||||
| -rw-r--r-- | src/diffraction.h | 2 | ||||
| -rw-r--r-- | src/intensities.c | 2 |
4 files changed, 39 insertions, 17 deletions
diff --git a/src/detector.c b/src/detector.c index ae7576b2..1bcddef1 100644 --- a/src/detector.c +++ b/src/detector.c @@ -181,8 +181,12 @@ void record_image(struct image *image, int do_water, int do_poisson, if ( do_water ) { + struct rvec q; + + q = get_q(image, x, y, 1, NULL, 1.0/image->lambda); + /* Add intensity contribution from water */ - water = water_intensity(get_q(image, x, y, 1, NULL), + water = water_intensity(q, ph_lambda_to_en(image->lambda), BEAM_RADIUS, WATER_RADIUS); intensity += water; diff --git a/src/diffraction.c b/src/diffraction.c index f262c9d1..bcd0e100 100644 --- a/src/diffraction.c +++ b/src/diffraction.c @@ -23,7 +23,9 @@ #include "sfac.h" -#define SAMPLING (1) +#define SAMPLING (2) +#define BWSAMPLING (10) +#define BANDWIDTH (0.015) static double lattice_factor(struct rvec q, double ax, double ay, double az, @@ -136,7 +138,7 @@ double water_intensity(struct rvec q, double en, struct rvec get_q(struct image *image, unsigned int xs, unsigned int ys, - unsigned int sampling, float *ttp) + unsigned int sampling, float *ttp, float k) { struct rvec q; float twothetax, twothetay, twotheta, r; @@ -144,7 +146,6 @@ struct rvec get_q(struct image *image, unsigned int xs, unsigned int ys, float ry = 0.0; int p; - const float k = 1.0/image->lambda; const unsigned int x = xs / sampling; const unsigned int y = ys / sampling; /* Integer part only */ @@ -184,6 +185,7 @@ void get_diffraction(struct image *image, int na, int nb, int nc, int no_sfac) double bx, by, bz; double cx, cy, cz; double a, b, c, d; + float kc; if ( image->molecule == NULL ) return; @@ -210,6 +212,8 @@ void get_diffraction(struct image *image, int na, int nb, int nc, int no_sfac) /* Needed later for Lorentz calculation */ image->twotheta = malloc(image->width * image->height * sizeof(double)); + kc = 1.0/image->lambda; /* Centre value */ + for ( xs=0; xs<image->width*SAMPLING; xs++ ) { for ( ys=0; ys<image->height*SAMPLING; ys++ ) { @@ -222,19 +226,33 @@ void get_diffraction(struct image *image, int na, int nb, int nc, int no_sfac) const unsigned int x = xs / SAMPLING; const unsigned int y = ys / SAMPLING; /* Integer part only */ - q = get_q(image, xs, ys, SAMPLING, &twotheta); - image->twotheta[x + image->width*y] = twotheta; + int kstep; - f_lattice = lattice_factor(q, ax,ay,az,bx,by,bz,cx,cy,cz, - na, nb, nc); - if ( no_sfac ) { - f_molecule = 10000.0; - } else { - f_molecule = molecule_factor(image->molecule, q, - ax,ay,az,bx,by,bz,cx,cy,cz); - } + for ( kstep=0; kstep<BWSAMPLING; kstep++ ) { + + float k; + + /* Calculate k this time round */ + k = kc + (kstep-(BWSAMPLING/2)) * + kc*(BANDWIDTH/BWSAMPLING); - image->sfacs[x + image->width*y] += sw * f_molecule * f_lattice; + q = get_q(image, xs, ys, SAMPLING, &twotheta, k); + image->twotheta[x + image->width*y] = twotheta; + + f_lattice = lattice_factor(q, ax, ay, az, + bx, by, bz, + cx, cy, cz, + na, nb, nc); + if ( no_sfac ) { + f_molecule = 10000.0; + } else { + f_molecule = molecule_factor(image->molecule, q, + ax,ay,az,bx,by,bz,cx,cy,cz); + } + + image->sfacs[x + image->width*y] += sw * f_molecule * f_lattice; + + } } progress_bar(xs, SAMPLING*image->width-1, "Calculating lattice factors"); diff --git a/src/diffraction.h b/src/diffraction.h index 36d1c21f..f26d9a4e 100644 --- a/src/diffraction.h +++ b/src/diffraction.h @@ -24,5 +24,5 @@ extern void get_diffraction(struct image *image, int na, int nb, int nc, extern double water_intensity(struct rvec q, double en, double beam_r, double water_r); extern struct rvec get_q(struct image *image, unsigned int xs, unsigned int ys, - unsigned int sampling, float *ttp); + unsigned int sampling, float *ttp, float k); #endif /* DIFFRACTION_H */ diff --git a/src/intensities.c b/src/intensities.c index 8b5f51c5..4ec496c3 100644 --- a/src/intensities.c +++ b/src/intensities.c @@ -73,7 +73,7 @@ void output_intensities(struct image *image, UnitCell *cell) int found = 0; int j; - q = get_q(image, x, y, 1, NULL); + q = get_q(image, x, y, 1, NULL, 1.0/image->lambda); hd = q.u * ax + q.v * ay + q.w * az; kd = q.u * bx + q.v * by + q.w * bz; |