Remove all bandwidth and divergence stuff, fix pattern_sim

Bandwidth and divergence didn't work very well

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;
 
 
 	}