aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/diffraction-gpu.c216
-rw-r--r--src/diffraction.c1
-rw-r--r--src/indexamajig.c3
-rw-r--r--src/pattern_sim.c5
4 files changed, 137 insertions, 88 deletions
diff --git a/src/diffraction-gpu.c b/src/diffraction-gpu.c
index 2b8e97e5..15ff39cc 100644
--- a/src/diffraction-gpu.c
+++ b/src/diffraction-gpu.c
@@ -167,57 +167,13 @@ static int set_arg_mem(struct gpu_context *gctx, int idx, cl_mem val)
}
-void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
- int na, int nb, int nc, UnitCell *ucell)
+static void do_panels(struct gpu_context *gctx, struct image *image,
+ int *n_inf, int *n_neg, int *n_nan, double k, double frac,
+ double xo, double yo)
{
- cl_int err;
- double ax, ay, az;
- double bx, by, bz;
- double cx, cy, cz;
int i;
- cl_float16 cell;
- cl_int4 ncells;
- int n_inf = 0;
- int n_neg = 0;
- int n_nan = 0;
- if ( gctx == NULL ) {
- ERROR("GPU setup failed.\n");
- return;
- }
-
- cell_get_cartesian(ucell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
- cell.s[0] = ax; cell.s[1] = ay; cell.s[2] = az;
- cell.s[3] = bx; cell.s[4] = by; cell.s[5] = bz;
- cell.s[6] = cx; cell.s[7] = cy; cell.s[8] = cz;
-
- ncells.s[0] = na;
- ncells.s[1] = nb;
- ncells.s[2] = nc;
- ncells.s[3] = 0; /* unused */
-
- /* Ensure all required LUTs are available */
- check_sinc_lut(gctx, na);
- check_sinc_lut(gctx, nb);
- check_sinc_lut(gctx, nc);
-
- if ( set_arg_float(gctx, 2, 1.0/image->lambda) ) return;
- if ( set_arg_mem(gctx, 9, gctx->intensities) ) return;
- if ( set_arg_mem(gctx, 10, gctx->sinc_luts[na-1]) ) return;
- if ( set_arg_mem(gctx, 11, gctx->sinc_luts[nb-1]) ) return;
- if ( set_arg_mem(gctx, 12, gctx->sinc_luts[nc-1]) ) return;
- if ( set_arg_mem(gctx, 13, gctx->flags) ) return;
-
- /* Unit cell */
- err = clSetKernelArg(gctx->kern, 8, sizeof(cl_float16), &cell);
- if ( err != CL_SUCCESS ) {
- ERROR("Couldn't set unit cell: %s\n", clError(err));
- return;
- }
-
- /* Allocate memory for the result */
- image->data = calloc(image->width * image->height, sizeof(float));
- image->twotheta = calloc(image->width * image->height, sizeof(double));
+ if ( set_arg_float(gctx, 1, k) ) return;
/* Iterate over panels */
for ( i=0; i<image->det->n_panels; i++ ) {
@@ -225,14 +181,12 @@ void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
size_t dims[2];
size_t ldims[2] = {1, 1};
struct panel *p;
- cl_mem tt;
- size_t tt_size;
cl_mem diff;
size_t diff_size;
float *diff_ptr;
- float *tt_ptr;
int pan_width, pan_height;
int fs, ss;
+ cl_int err;
p = &image->det->panels[i];
@@ -247,25 +201,19 @@ void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
ERROR("Couldn't allocate diffraction memory\n");
return;
}
- tt_size = pan_width * pan_height * sizeof(cl_float);
- tt = clCreateBuffer(gctx->ctx, CL_MEM_WRITE_ONLY, tt_size,
- NULL, &err);
- if ( err != CL_SUCCESS ) {
- ERROR("Couldn't allocate twotheta memory\n");
- return;
- }
if ( set_arg_mem(gctx, 0, diff) ) return;
- if ( set_arg_mem(gctx, 1, tt) ) return;
- if ( set_arg_int(gctx, 3, pan_width) ) return;
- if ( set_arg_float(gctx, 4, p->cnx) ) return;
- if ( set_arg_float(gctx, 5, p->cny) ) return;
- if ( set_arg_float(gctx, 6, p->res) ) return;
- if ( set_arg_float(gctx, 7, p->clen) ) return;
- if ( set_arg_float(gctx, 14, p->fsx) ) return;
- if ( set_arg_float(gctx, 15, p->fsy) ) return;
- if ( set_arg_float(gctx, 16, p->ssx) ) return;
- if ( set_arg_float(gctx, 17, p->ssy) ) return;
+ if ( set_arg_int(gctx, 2, pan_width) ) return;
+ if ( set_arg_float(gctx, 3, p->cnx) ) return;
+ if ( set_arg_float(gctx, 4, p->cny) ) return;
+ if ( set_arg_float(gctx, 5, p->res) ) return;
+ if ( set_arg_float(gctx, 6, p->clen) ) return;
+ if ( set_arg_float(gctx, 13, p->fsx) ) return;
+ if ( set_arg_float(gctx, 14, p->fsy) ) return;
+ if ( set_arg_float(gctx, 15, p->ssx) ) return;
+ if ( set_arg_float(gctx, 16, p->ssy) ) return;
+ if ( set_arg_float(gctx, 17, xo) ) return;
+ if ( set_arg_float(gctx, 18, yo) ) return;
dims[0] = pan_width;
dims[1] = pan_height;
@@ -288,43 +236,124 @@ void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
clError(err));
return;
}
- tt_ptr = clEnqueueMapBuffer(gctx->cq, tt, CL_TRUE, CL_MAP_READ,
- 0, tt_size, 0, NULL, NULL, &err);
- if ( err != CL_SUCCESS ) {
- ERROR("Couldn't map tt buffer\n");
- return;
- }
for ( fs=0; fs<pan_width; fs++ ) {
for ( ss=0; ss<pan_height; ss++ ) {
- float val, tt;
+ float val;
int tfs, tss;
val = diff_ptr[fs + pan_width*ss];
- if ( isinf(val) ) n_inf++;
- if ( val < 0.0 ) n_neg++;
- if ( isnan(val) ) n_nan++;
- tt = tt_ptr[fs + pan_width*ss];
+ if ( isinf(val) ) (*n_inf)++;
+ if ( val < 0.0 ) (*n_neg)++;
+ if ( isnan(val) ) (*n_nan)++;
tfs = p->min_fs + fs;
tss = p->min_ss + ss;
- image->data[tfs + image->width*tss] = val;
- image->twotheta[tfs + image->width*tss] = tt;
+ image->data[tfs + image->width*tss] += frac*val;
}
}
clEnqueueUnmapMemObject(gctx->cq, diff, diff_ptr,
0, NULL, NULL);
- clEnqueueUnmapMemObject(gctx->cq, tt, tt_ptr,
- 0, NULL, NULL);
clReleaseMemObject(diff);
- clReleaseMemObject(tt);
}
+}
+
+
+void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
+ int na, int nb, int nc, UnitCell *ucell)
+{
+ double ax, ay, az;
+ double bx, by, bz;
+ double cx, cy, cz;
+ cl_float16 cell;
+ cl_int4 ncells;
+ cl_int err;
+ int n_inf = 0;
+ int n_neg = 0;
+ int n_nan = 0;
+ int fs, ss;
+ const int nkstep = 10;
+ const int nxs = 4;
+ const int nys = 4;
+ int kstep;
+ double klow, khigh, kinc, w;
+ double frac;
+ int xs, ys;
+ int n, ntot;
+
+ if ( gctx == NULL ) {
+ ERROR("GPU setup failed.\n");
+ return;
+ }
+ cell_get_cartesian(ucell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
+ cell.s[0] = ax; cell.s[1] = ay; cell.s[2] = az;
+ cell.s[3] = bx; cell.s[4] = by; cell.s[5] = bz;
+ cell.s[6] = cx; cell.s[7] = cy; cell.s[8] = cz;
+
+ ncells.s[0] = na;
+ ncells.s[1] = nb;
+ ncells.s[2] = nc;
+ ncells.s[3] = 0; /* unused */
+
+ /* Ensure all required LUTs are available */
+ check_sinc_lut(gctx, na);
+ check_sinc_lut(gctx, nb);
+ check_sinc_lut(gctx, nc);
+
+ if ( set_arg_mem(gctx, 8, gctx->intensities) ) return;
+ if ( set_arg_mem(gctx, 9, gctx->sinc_luts[na-1]) ) return;
+ if ( set_arg_mem(gctx, 10, gctx->sinc_luts[nb-1]) ) return;
+ if ( set_arg_mem(gctx, 11, gctx->sinc_luts[nc-1]) ) return;
+ if ( set_arg_mem(gctx, 12, gctx->flags) ) return;
+
+ /* Unit cell */
+ err = clSetKernelArg(gctx->kern, 7, sizeof(cl_float16), &cell);
+ if ( err != CL_SUCCESS ) {
+ ERROR("Couldn't set unit cell: %s\n", clError(err));
+ return;
+ }
+
+ /* Allocate memory for the result */
+ image->data = calloc(image->width * image->height, sizeof(float));
+
+ w = image->lambda * image->bw;
+ klow = 1.0 / (image->lambda + (w/2.0)); /* Smallest k */
+ khigh = 1.0 / (image->lambda - (w/2.0)); /* Largest k */
+ kinc = (khigh-klow) / (nkstep+1);
+
+ ntot = nkstep * nxs * nys;
+ frac = 1.0 / ntot;
+
+ n = 0;
+ for ( xs=0; xs<nxs; xs++ ) {
+ for ( ys=0; ys<nys; ys++ ) {
+
+ double xo, yo;
+
+ xo = (1.0/nxs) * xs;
+ yo = (1.0/nys) * ys;
+
+ for ( kstep=0; kstep<nkstep; kstep++ ) {
+
+ double k;
+
+ k = klow + kstep*kinc;
+
+ do_panels(gctx, image, &n_inf, &n_neg, &n_nan, k, frac,
+ xo, yo);
+
+ n++;
+ progress_bar(n, ntot, "Simulating");
+
+ }
+ }
+ }
if ( n_neg + n_inf + n_nan ) {
ERROR("WARNING: The GPU calculation produced %i negative"
@@ -332,6 +361,25 @@ void get_diffraction_gpu(struct gpu_context *gctx, struct image *image,
n_neg, n_inf, n_nan);
}
+ /* Calculate "2theta" values for detector geometry */
+ image->twotheta = calloc(image->width * image->height, sizeof(double));
+ for ( fs=0; fs<image->width; fs++ ) {
+ for ( ss=0; ss<image->height; ss++ ) {
+
+ struct rvec q;
+ double twotheta, k;
+ int idx;
+
+ /* Calculate k this time round */
+ k = 1.0/image->lambda;
+
+ q = get_q(image, fs, ss, &twotheta, k);
+
+ idx = fs + image->width*ss;
+ image->twotheta[idx] = twotheta;
+
+ }
+ }
}
diff --git a/src/diffraction.c b/src/diffraction.c
index 870569b7..b76154c2 100644
--- a/src/diffraction.c
+++ b/src/diffraction.c
@@ -420,7 +420,6 @@ void get_diffraction(struct image *image, int na, int nb, int nc,
image->data[idx] = I_lattice * I_molecule;
image->twotheta[idx] = twotheta;
-
}
progress_bar(fs, image->width-1, "Calculating diffraction");
}
diff --git a/src/indexamajig.c b/src/indexamajig.c
index 8061b1ae..2c1af284 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -291,9 +291,10 @@ static void process_image(const struct index_args *iargs,
image.indexed_cell = NULL;
image.det = copy_geom(iargs->det);
image.copyme = iargs->copyme;
- image.beam = beam;
+ image.reflections = NULL;
image.id = cookie;
image.filename = pargs->filename;
+ image.beam = beam;
hdfile = hdfile_open(image.filename);
if ( hdfile == NULL ) return;
diff --git a/src/pattern_sim.c b/src/pattern_sim.c
index 03b009b4..e05dcdec 100644
--- a/src/pattern_sim.c
+++ b/src/pattern_sim.c
@@ -497,6 +497,8 @@ int main(int argc, char *argv[])
image.width = image.det->max_fs + 1;
image.height = image.det->max_ss + 1;
image.lambda = ph_en_to_lambda(eV_to_J(image.beam->photon_energy));
+ image.bw = image.beam->bandwidth;
+ image.div = image.beam->divergence;
/* Load unit cell */
input_cell = load_cell_from_pdb(filename);
@@ -617,8 +619,7 @@ int main(int argc, char *argv[])
number++;
/* Write the output file */
- hdf5_write(filename, image.data,
- image.width, image.height, H5T_NATIVE_FLOAT);
+ hdf5_write_image(filename, &image);
}
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-15 10:15:50 +0000