aboutsummaryrefslogtreecommitdiff
path: root/src/diffraction.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/diffraction.c')
-rw-r--r--src/diffraction.c97
1 files changed, 22 insertions, 75 deletions
diff --git a/src/diffraction.c b/src/diffraction.c
index de994133..a37c2fcb 100644
--- a/src/diffraction.c
+++ b/src/diffraction.c
@@ -27,9 +27,6 @@
#include "pattern_sim.h"
-#define SAMPLING (4)
-#define BWSAMPLING (10)
-#define DIVSAMPLING (1)
#define SINC_LUT_ELEMENTS (4096)
@@ -355,11 +352,9 @@ void get_diffraction(struct image *image, int na, int nb, int nc,
double ax, ay, az;
double bx, by, bz;
double cx, cy, cz;
- float klow, khigh, bwstep;
double *lut_a;
double *lut_b;
double *lut_c;
- double divxlow, divylow, divxstep, divystep;
cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
@@ -369,15 +364,6 @@ void get_diffraction(struct image *image, int na, int nb, int nc,
/* Needed later for Lorentz calculation */
image->twotheta = malloc(image->width * image->height * sizeof(double));
- klow = 1.0/(image->lambda*(1.0 + image->beam->bandwidth/2.0));
- khigh = 1.0/(image->lambda*(1.0 - image->beam->bandwidth/2.0));
- bwstep = (khigh-klow) / BWSAMPLING;
-
- divxlow = -image->beam->divergence/2.0;
- divylow = -image->beam->divergence/2.0;
- divxstep = image->beam->divergence / DIVSAMPLING;
- divystep = image->beam->divergence / DIVSAMPLING;
-
lut_a = get_sinc_lut(na);
lut_b = get_sinc_lut(nb);
lut_c = get_sinc_lut(nc);
@@ -385,73 +371,34 @@ void get_diffraction(struct image *image, int na, int nb, int nc,
for ( fs=0; fs<image->width; fs++ ) {
for ( ss=0; ss<image->height; ss++ ) {
- int fs_step, ss_step, kstep;
- int divxval, divyval;
- int idx = fs + image->width*ss;
-
- for ( fs_step=0; fs_step<SAMPLING; fs_step++ ) {
- for ( ss_step=0; ss_step<SAMPLING; ss_step++ ) {
- for ( kstep=0; kstep<BWSAMPLING; kstep++ ) {
- for ( divxval=0; divxval<DIVSAMPLING; divxval++ ) {
- for ( divyval=0; divyval<DIVSAMPLING; divyval++ ) {
-
- double k;
- double intensity;
- double f_lattice, I_lattice;
- double I_molecule;
- struct rvec q, qn;
- double twotheta;
- const double dfs = (double)fs
- + ((double)fs_step / SAMPLING);
- const double dss = (double)ss
- + ((double)ss_step / SAMPLING);
-
- double xdiv = divxlow + divxstep*(double)divxval;
- double ydiv = divylow + divystep*(double)divyval;
-
- /* Calculate k this time round */
- k = klow + (double)kstep * bwstep;
-
- qn = get_q(image, dfs, dss, &twotheta, k);
-
- /* x divergence */
- q.u = qn.u*cos(xdiv) +qn.w*sin(xdiv);
- q.v = qn.v;
- q.w = -qn.u*sin(xdiv) +qn.w*cos(xdiv);
+ int idx;
+ double k;
+ double intensity;
+ double f_lattice, I_lattice;
+ double I_molecule;
+ struct rvec q;
+ double twotheta;
- qn = q;
+ /* Calculate k this time round */
+ k = 1.0/image->lambda;
- /* y divergence */
- q.v = qn.v*cos(ydiv) +qn.w*sin(ydiv);
- q.w = -qn.v*sin(ydiv) +qn.w*cos(ydiv);
+ q = get_q(image, fs, ss, &twotheta, k);
- f_lattice = lattice_factor(q, ax, ay, az,
- bx, by, bz,
- cx, cy, cz,
- lut_a, lut_b, lut_c);
+ f_lattice = lattice_factor(q, ax, ay, az,
+ bx, by, bz,
+ cx, cy, cz,
+ lut_a, lut_b, lut_c);
- I_molecule = molecule_factor(intensities,
- phases, flags, q,
- ax,ay,az,bx,by,bz,cx,cy,cz,
- m, sym);
+ I_molecule = molecule_factor(intensities,
+ phases, flags, q,
+ ax,ay,az,bx,by,bz,cx,cy,cz,
+ m, sym);
- I_lattice = pow(f_lattice, 2.0);
- intensity = I_lattice * I_molecule;
-
- image->data[idx] += intensity;
-
- if ( fs_step + ss_step + kstep == 0 ) {
- image->twotheta[idx] = twotheta;
- }
-
- }
- }
- }
- }
- }
+ I_lattice = pow(f_lattice, 2.0);
- image->data[idx] /= (SAMPLING*SAMPLING*BWSAMPLING
- *DIVSAMPLING*DIVSAMPLING);
+ idx = fs + image->width*ss;
+ image->data[idx] = I_lattice * I_molecule;
+ image->twotheta[idx] = twotheta;
}