pattern_sim.c 16.5 KB
Newer Older
1
/*
Thomas White's avatar
Thomas White committed
2
3
4
 * pattern_sim.c
 *
 * Simulate diffraction patterns from small crystals
5
 *
Thomas White's avatar
Thomas White committed
6
 * (c) 2006-2010 Thomas White <taw@physics.org>
7
 *
Thomas White's avatar
Thomas White committed
8
 * Part of CrystFEL - crystallography with a FEL
9
10
11
12
13
14
15
16
17
18
19
20
21
 *
 */


#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <stdarg.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <unistd.h>
Thomas White's avatar
Thomas White committed
22
#include <getopt.h>
23

24
#include "image.h"
Thomas White's avatar
Thomas White committed
25
#include "diffraction.h"
26
#include "diffraction-gpu.h"
27
#include "cell.h"
28
29
#include "utils.h"
#include "hdf5-file.h"
Thomas White's avatar
Thomas White committed
30
#include "detector.h"
31
#include "peaks.h"
Thomas White's avatar
Thomas White committed
32
#include "sfac.h"
33
#include "reflections.h"
34
#include "beam-parameters.h"
Thomas White's avatar
Thomas White committed
35
#include "symmetry.h"
36
37


Thomas White's avatar
Thomas White committed
38
static void show_help(const char *s)
39
{
40
	printf("Syntax: %s [options]\n\n", s);
Thomas White's avatar
Thomas White committed
41
	printf(
Thomas White's avatar
Thomas White committed
42
"Simulate diffraction patterns from small crystals probed with femtosecond\n"
Thomas White's avatar
Thomas White committed
43
44
"pulses of X-rays from a free electron laser.\n"
"\n"
45
" -h, --help                Display this help message.\n"
Thomas White's avatar
Thomas White committed
46
"\n"
47
48
49
" -p, --pdb=<file>          PDB file from which to get the unit cell.\n"
"                            (The actual Bragg intensities come from the\n"
"                            intensities file)\n"
50
"     --simulation-details  Show technical details of the simulation.\n"
51
"     --gpu                 Use the GPU to speed up the calculation.\n"
Thomas White's avatar
Thomas White committed
52
53
"     --gpu-dev=<n>         Use GPU device <n>.  Omit this option to see the\n"
"                            available devices.\n"
54
55
" -g, --geometry=<file>     Get detector geometry from file.\n"
" -b, --beam=<file>         Get beam parameters from file.\n"
56
57
58
59
"\n"
"     --near-bragg          Output h,k,l,I near Bragg conditions.\n"
" -n, --number=<N>          Generate N images.  Default 1.\n"
"     --no-images           Do not output any HDF5 files.\n"
60
61
62
63
64
65
" -o, --output=<filename>   Output HDF5 filename.  Default: sim.h5.\n"
"                            If you choose to simulate more than one pattern,\n"
"                            the filename given will be postfixed with a\n"
"                            hyphen, the image number and '.h5'.  In this\n"
"                            case, the default value is 'sim', such that the\n"
"                            files produced are sim-1.h5, sim-2.h5 and so on.\n"
Thomas White's avatar
Thomas White committed
66
" -r, --random-orientation  Use a randomly generated orientation\n"
67
"                            (a new orientation will be used for each image).\n"
68
69
"     --powder=<file>       Write a summed pattern of all images simulated by\n"
"                            this invocation as the given filename.\n"
70
" -i, --intensities=<file>  Specify file containing reflection intensities\n"
Thomas White's avatar
Thomas White committed
71
"                            (and phases) to use.\n"
72
" -y, --symmetry=<sym>      The symmetry of the intensities file.\n"
Thomas White's avatar
Thomas White committed
73
" -t, --gradients=<method>  Use <method> for the calculation of shape\n"
74
75
"                            transform intensities.  Choose from:\n"
"                             mosaic      : Take the intensity of the nearest\n"
Thomas White's avatar
Thomas White committed
76
77
78
"                                           Bragg position.  This is the\n"
"                                           fastest method and the only one\n"
"                                           supported on the GPU.\n"
79
80
"                             interpolate : Interpolate trilinearly between\n"
"                                           six adjacent Bragg intensities.\n"
Thomas White's avatar
Thomas White committed
81
82
"                                           This method has intermediate\n"
"                                           accuracy.\n"
83
"                             phased      : As 'interpolate', but take phase\n"
Thomas White's avatar
Thomas White committed
84
85
86
"                                           values into account.  This is the\n"
"                                           most accurate method, but the\n"
"                                           slowest.\n"
87
88
89
"     --really-random       Use really random numbers for the orientation and\n"
"                            crystal size.  By default, the same sequence\n"
"                            will be used for each new run.\n"
90
91
92
93
94
"     --min-size=<s>        Generate random crystal sizes using <s> as the\n"
"                            minimum crystal size in nm. --max-size is also\n"
"                            required.\n"
"     --max-size=<s>        Use <s> as the maximum crystal size in nm.\n"
"                            --min-size is also required.\n"
Thomas White's avatar
Thomas White committed
95
96
97
98
99
"\n"
"By default, the simulation aims to be as accurate as possible.  For greater\n"
"speed, or for testing, you can choose to disable certain things using the\n"
"following options.\n"
"\n"
Thomas White's avatar
Thomas White committed
100
"     --no-noise            Do not calculate Poisson noise.\n"
101
);
Thomas White's avatar
Thomas White committed
102
103
104
}


105
static void show_details()
Thomas White's avatar
Thomas White committed
106
{
Thomas White's avatar
Thomas White committed
107
108
109
110
111
	printf(
"This program simulates diffraction patterns from small crystals illuminated\n"
"with femtosecond X-ray pulses from a free electron laser.\n"
"\n"
"The lattice transform from the specified number of unit cells is calculated\n"
112
113
"using the closed-form solution for a truncated lattice faceted on the\n"
"(001), (010) and (100) planes:\n"
Thomas White's avatar
Thomas White committed
114
"\n"
115
116
117
"I_latt(q) =  sin^2(pi*na*q.a)/sin^2(pi*q.a)\n"
"           * sin^2(pi*nb*q.b)/sin^2(pi*q.b)\n"
"           * sin^2(pi*nc*q.c)/sin^2(pi*q.c)\n"
Thomas White's avatar
Thomas White committed
118
119
120
121
"\n"
"na = number of unit cells in 'a' direction (likewise nb, nc)\n"
" q = reciprocal vector (1/d convention, not 2pi/d)\n"
"\n"
122
123
124
"This is multiplied by a model of the underlying molecular transform, I_mol(q).\n"
"This can be approximated to varying levels of accuracy by the methods given by\n"
"the '--gradients' option.\n"
Thomas White's avatar
Thomas White committed
125
126
127
"\n"
"Expected intensities at the CCD are then calculated using:\n"
"\n"
128
"I(q) = I0 * r^2 * I_latt(q) * I_mol(q) * S\n"
Thomas White's avatar
Thomas White committed
129
130
131
132
"\n"
"I0 = number of photons per unit area in the incident beam\n"
" r = Thomson radius\n"
" S = solid angle of corresponding pixel\n"
133
134
135
"\n"
"Polarisation is not currently included in pattern_sim, although it is included\n"
"in the analysis of Bragg peaks by 'indexamajig'.\n"
Thomas White's avatar
Thomas White committed
136
137
"\n"
"Poisson counts are generated from the expected intensities using Knuth's\n"
Thomas White's avatar
Thomas White committed
138
139
140
"algorithm.  When the intensity is sufficiently high that Knuth's algorithm\n"
"would result in machine precision problems, a normal distribution with\n"
"standard deviation sqrt(I) is used instead.\n"
141
);
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
}


static struct quaternion read_quaternion()
{
	do {

		int r;
		float w, x, y, z;
		char line[1024];
		char *rval;

		printf("Please input quaternion: w x y z\n");
		rval = fgets(line, 1023, stdin);
		if ( rval == NULL ) return invalid_quaternion();
		chomp(line);

		r = sscanf(line, "%f %f %f %f", &w, &x, &y, &z);
		if ( r == 4 ) {

			struct quaternion quat;

			quat.w = w;
			quat.x = x;
			quat.y = y;
			quat.z = z;

			return quat;

		} else {
172
			ERROR("Invalid rotation '%s'\n", line);
173
174
175
		}

	} while ( 1 );
176
177
178
}


179
180
181
182
183
184
185
186
187
188
189
190
static int random_ncells(double len, double min_nm, double max_nm)
{
	int min_cells, max_cells, wid;

	min_cells = min_nm*1e-9 / len;
	max_cells = max_nm*1e-9 / len;
	wid = max_cells - min_cells;

	return wid*random()/RAND_MAX + min_cells;
}


191
192
int main(int argc, char *argv[])
{
193
	int c;
194
	struct image image;
195
	struct gpu_context *gctx = NULL;
196
	double *powder;
197
198
	char *intfile = NULL;
	double *intensities;
199
	char *rval;
200
	double *phases;
201
	unsigned char *flags;
Thomas White's avatar
Thomas White committed
202
	int config_simdetails = 0;
203
204
	int config_nearbragg = 0;
	int config_randomquat = 0;
205
	int config_noimages = 0;
Thomas White's avatar
Thomas White committed
206
	int config_nonoise = 0;
207
	int config_nosfac = 0;
208
	int config_gpu = 0;
209
	int config_random = 0;
210
	char *powder_fn = NULL;
211
	char *filename = NULL;
212
	char *grad_str = NULL;
213
	char *outfile = NULL;
Thomas White's avatar
Thomas White committed
214
	char *geometry = NULL;
215
	char *beamfile = NULL;
216
	GradientMethod grad;
Thomas White's avatar
Thomas White committed
217
218
219
	int ndone = 0;    /* Number of simulations done (images or not) */
	int number = 1;   /* Number used for filename of image */
	int n_images = 1; /* Generate one image by default */
220
	int done = 0;
Thomas White's avatar
Thomas White committed
221
222
	UnitCell *input_cell;
	struct quaternion orientation;
Thomas White's avatar
Thomas White committed
223
	int gpu_dev = -1;
224
225
226
	int random_size = 0;
	double min_size = 0.0;
	double max_size = 0.0;
227
	char *sym = NULL;
Thomas White's avatar
Thomas White committed
228

229
	/* Long options */
Thomas White's avatar
Thomas White committed
230
	const struct option longopts[] = {
231
232
		{"help",               0, NULL,               'h'},
		{"simulation-details", 0, &config_simdetails,  1},
233
		{"gpu",                0, &config_gpu,         1},
234
235
236
		{"near-bragg",         0, &config_nearbragg,   1},
		{"random-orientation", 0, NULL,               'r'},
		{"number",             1, NULL,               'n'},
Thomas White's avatar
Thomas White committed
237
		{"no-images",          0, &config_noimages,    1},
Thomas White's avatar
Thomas White committed
238
		{"no-noise",           0, &config_nonoise,     1},
239
		{"intensities",        1, NULL,               'i'},
240
		{"symmetry",           1, NULL,               'y'},
241
		{"powder",             1, NULL,               'w'},
Thomas White's avatar
Thomas White committed
242
		{"gradients",          1, NULL,               't'},
243
		{"pdb",                1, NULL,               'p'},
244
		{"output",             1, NULL,               'o'},
Thomas White's avatar
Thomas White committed
245
		{"geometry",           1, NULL,               'g'},
246
		{"beam",               1, NULL,               'b'},
247
		{"really-random",      0, &config_random,      1},
Thomas White's avatar
Thomas White committed
248
		{"gpu-dev",            1, NULL,                2},
249
250
		{"min-size",           1, NULL,                3},
		{"max-size",           1, NULL,                4},
251
		{0, 0, NULL, 0}
Thomas White's avatar
Thomas White committed
252
	};
253

254
	/* Short options */
255
	while ((c = getopt_long(argc, argv, "hrn:i:t:p:o:g:b:y:",
Thomas White's avatar
Thomas White committed
256
	                        longopts, NULL)) != -1) {
257

Thomas White's avatar
Thomas White committed
258
		switch (c) {
Thomas White's avatar
Thomas White committed
259
		case 'h' :
Thomas White's avatar
Thomas White committed
260
261
			show_help(argv[0]);
			return 0;
262

Thomas White's avatar
Thomas White committed
263
		case 'r' :
264
265
266
			config_randomquat = 1;
			break;

Thomas White's avatar
Thomas White committed
267
		case 'n' :
268
269
270
271
272
			n_images = strtol(optarg, &rval, 10);
			if ( *rval != '\0' ) {
				ERROR("Invalid number of images.\n");
				return 1;
			}
273
274
			break;

Thomas White's avatar
Thomas White committed
275
		case 'i' :
276
277
278
			intfile = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
279
		case 't' :
280
281
282
			grad_str = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
283
		case 'p' :
284
285
286
			filename = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
287
		case 'o' :
288
289
290
			outfile = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
291
		case 'w' :
292
293
294
			powder_fn = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
295
		case 'g' :
Thomas White's avatar
Thomas White committed
296
297
298
			geometry = strdup(optarg);
			break;

299
300
301
302
		case 'b' :
			beamfile = strdup(optarg);
			break;

303
304
305
306
		case 'y' :
			sym = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
307
308
309
310
		case 2 :
			gpu_dev = atoi(optarg);
			break;

311
		case 3 :
312
313
314
315
316
			min_size = strtod(optarg, &rval);
			if ( *rval != '\0' ) {
				ERROR("Invalid minimum size.\n");
				return 1;
			}
317
318
319
320
			random_size++;
			break;

		case 4 :
321
322
323
324
325
			max_size = strtod(optarg, &rval);
			if ( *rval != '\0' ) {
				ERROR("Invalid maximum size.\n");
				return 1;
			}
326
327
328
			random_size++;
			break;

Thomas White's avatar
Thomas White committed
329
		case 0 :
Thomas White's avatar
Thomas White committed
330
			break;
331

Thomas White's avatar
Thomas White committed
332
		default :
Thomas White's avatar
Thomas White committed
333
334
			return 1;
		}
335
336
337

	}

338
	if ( config_random ) {
339
340
341
342
343
344
345
346
		FILE *fh;
		unsigned int seed;
		fh = fopen("/dev/urandom", "r");
		fread(&seed, sizeof(seed), 1, fh);
		fclose(fh);
		srand(seed);
	}

347
348
349
350
351
	if ( random_size == 1 ) {
		ERROR("You must specify both --min-size and --max-size.\n");
		return 1;
	}

352
353
354
355
	if ( filename == NULL ) {
		filename = strdup("molecule.pdb");
	}

356
357
358
359
360
361
362
363
	if ( outfile == NULL ) {
		if ( n_images == 1 ) {
			outfile = strdup("sim.h5");
		} else {
			outfile = strdup("sim");
		}
	}

364
365
	if ( sym == NULL ) sym = strdup("1");

Thomas White's avatar
Thomas White committed
366
	if ( config_simdetails ) {
367
		show_details();
Thomas White's avatar
Thomas White committed
368
369
370
		return 0;
	}

371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
	if ( grad_str == NULL ) {
		STATUS("You didn't specify a gradient calculation method, so"
		       " I'm using the 'mosaic' method, which is fastest.\n");
		grad = GRADIENT_MOSAIC;
	} else if ( strcmp(grad_str, "mosaic") == 0 ) {
		grad = GRADIENT_MOSAIC;
	} else if ( strcmp(grad_str, "interpolate") == 0) {
		grad = GRADIENT_INTERPOLATE;
	} else if ( strcmp(grad_str, "phased") == 0) {
		grad = GRADIENT_PHASED;
	} else {
		ERROR("Unrecognised gradient method '%s'\n", grad_str);
		return 1;
	}
	free(grad_str);

	if ( config_gpu && (grad != GRADIENT_MOSAIC) ) {
		ERROR("Only the mosaic method can be used for gradients when"
		      "calculating on the GPU.\n");
		return 1;
	}

Thomas White's avatar
Thomas White committed
393
394
395
396
397
	if ( geometry == NULL ) {
		ERROR("You need to specify a geometry file with --geometry\n");
		return 1;
	}

398
399
400
401
402
403
	if ( beamfile == NULL ) {
		ERROR("You need to specify a beam parameter file"
		      " with --beam\n");
		return 1;
	}

404
	if ( intfile == NULL ) {
405

406
407
408
409
410
411
		/* Gentle reminder */
		STATUS("You didn't specify the file containing the ");
		STATUS("reflection intensities (with --intensities).\n");
		STATUS("I'll simulate a flat intensity distribution.\n");
		intensities = NULL;
		phases = NULL;
412
413
		flags = NULL;

414
	} else {
415
416

		int i;
Thomas White's avatar
Thomas White committed
417
		ReflItemList *items;
418

419
420
421
422
423
		if ( grad == GRADIENT_PHASED ) {
			phases = new_list_phase();
		} else {
			phases = NULL;
		}
Thomas White's avatar
Thomas White committed
424
425
		intensities = new_list_intensity();
		phases = new_list_phase();
426
		flags = new_list_flag();
427
428
		items = read_reflections(intfile, intensities, phases,
		                         NULL, NULL);
429
		free(intfile);
430
431
432
433
434
435

		for ( i=0; i<num_items(items); i++ ) {
			struct refl_item *it = get_item(items, i);
			set_flag(flags, it->h, it->k, it->l, 1);
		}

Thomas White's avatar
Thomas White committed
436
437
438
439
440
441
442
		/* Check that the intensities have the correct symmetry */
		if ( check_symmetry(items, sym) ) {
			ERROR("The input reflection list does not appear to"
			      " have symmetry %s\n", sym);
			return 1;
		}

Thomas White's avatar
Thomas White committed
443
		delete_items(items);
444
445
	}

446
447
448
449
450
451
452
	image.det = get_detector_geometry(geometry);
	if ( image.det == NULL ) {
		ERROR("Failed to read detector geometry from '%s'\n", geometry);
		return 1;
	}
	free(geometry);

453
454
455
456
457
458
459
	image.beam = get_beam_parameters(beamfile);
	if ( image.beam == NULL ) {
		ERROR("Failed to read beam parameters from '%s'\n", beamfile);
		return 1;
	}
	free(beamfile);

Thomas White's avatar
Thomas White committed
460
	/* Define image parameters */
Thomas White's avatar
Thomas White committed
461
462
	image.width = image.det->max_fs + 1;
	image.height = image.det->max_ss + 1;
463
	image.lambda = ph_en_to_lambda(eV_to_J(image.beam->photon_energy));
Thomas White's avatar
Thomas White committed
464
465
466
467

	/* Load unit cell */
	input_cell = load_cell_from_pdb(filename);
	if ( input_cell == NULL ) {
Thomas White's avatar
Thomas White committed
468
469
		exit(1);
	}
Thomas White's avatar
Thomas White committed
470
471

	/* Initialise stuff */
472
	image.filename = NULL;
Thomas White's avatar
Thomas White committed
473
474
	image.features = NULL;
	image.flags = NULL;
Thomas White's avatar
Thomas White committed
475
476
	image.i0 = 1.0;
	image.i0_available = 1;
Thomas White's avatar
Thomas White committed
477

478
479
	powder = calloc(image.width*image.height, sizeof(*powder));

Thomas White's avatar
Thomas White committed
480
	/* Splurge a few useful numbers */
481
	STATUS("Wavelength is %f nm\n", image.lambda/1.0e-9);
Thomas White's avatar
Thomas White committed
482

483
484
	do {

Thomas White's avatar
Thomas White committed
485
		int na, nb, nc;
486
		double a, b, c, d;
Thomas White's avatar
Thomas White committed
487
		UnitCell *cell;
Thomas White's avatar
Thomas White committed
488

489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
		if ( random_size ) {

			double alen, blen, clen, dis;

			cell_get_parameters(input_cell, &alen, &blen, &clen,
			                    &dis, &dis, &dis);

			na = random_ncells(alen, min_size, max_size);
			nb = random_ncells(blen, min_size, max_size);
			nc = random_ncells(clen, min_size, max_size);

		} else {

			na = 8;
			nb = 8;
			nc = 8;

		}
Thomas White's avatar
Thomas White committed
507

508
509
		/* Read quaternion from stdin */
		if ( config_randomquat ) {
Thomas White's avatar
Thomas White committed
510
			orientation = random_quaternion();
511
		} else {
Thomas White's avatar
Thomas White committed
512
			orientation = read_quaternion();
513
514
		}

Thomas White's avatar
Thomas White committed
515
		STATUS("Orientation is %5.3f %5.3f %5.3f %5.3f\n",
Thomas White's avatar
Thomas White committed
516
517
		       orientation.w, orientation.x,
		       orientation.y, orientation.z);
Thomas White's avatar
Thomas White committed
518

Thomas White's avatar
Thomas White committed
519
		if ( !quaternion_valid(orientation) ) {
520
			ERROR("Orientation modulus is not zero!\n");
521
522
523
			return 1;
		}

Thomas White's avatar
Thomas White committed
524
525
		cell = cell_rotate(input_cell, orientation);

526
527
528
		/* Ensure no residual information */
		image.data = NULL;
		image.twotheta = NULL;
529

530
		cell_get_parameters(cell, &a, &b, &c, &d, &d, &d);
Thomas White's avatar
Thomas White committed
531
532
		STATUS("Particle size = %i x %i x %i"
		       " ( = %5.2f x %5.2f x %5.2f nm)\n",
533
534
	               na, nb, nc, na*a/1.0e-9, nb*b/1.0e-9, nc*c/1.0e-9);

535
		if ( config_gpu ) {
536
			if ( gctx == NULL ) {
537
				gctx = setup_gpu(config_nosfac,
538
539
				                 intensities, flags, sym,
				                 gpu_dev);
540
			}
541
			get_diffraction_gpu(gctx, &image, na, nb, nc, cell);
542
		} else {
Thomas White's avatar
Thomas White committed
543
			get_diffraction(&image, na, nb, nc, intensities, phases,
Thomas White's avatar
Thomas White committed
544
			                flags, cell, grad, sym);
545
		}
546
		if ( image.data == NULL ) {
547
548
549
550
			ERROR("Diffraction calculation failed.\n");
			goto skip;
		}

551
		record_image(&image, !config_nonoise);
Thomas White's avatar
Thomas White committed
552

Thomas White's avatar
Thomas White committed
553
		if ( config_nearbragg ) {
554

Thomas White's avatar
Thomas White committed
555
			/* Do EITHER: */
556

Thomas White's avatar
Thomas White committed
557
558
559
560
561
			//image.div = beam->divergence;
			//image.bw = beam->bandwidth;
			//image.profile_radius = 0.0001e9;
			//image.reflections = find_intersections(&image,
			//                               image.indexed_cell, 0);
562

Thomas White's avatar
Thomas White committed
563
564
			image.reflections = find_projected_peaks(&image,
			                            image.indexed_cell, 0, 0.1);
565

Thomas White's avatar
Thomas White committed
566
567
568
569
570
571
			integrate_reflections(&image, 0, 0);

			/* OR */

			//image.reflections = integrate_pixels(&image, 0, 0.1,
			//                                     config_polar);
572

Thomas White's avatar
Thomas White committed
573
574
		}

575
		if ( powder_fn != NULL ) {
576
577
578
579
580
581
582

			int x, y, w;

			w = image.width;

			for ( x=0; x<image.width; x++ ) {
			for ( y=0; y<image.height; y++ ) {
583
				powder[x+w*y] += (double)image.data[x+w*y];
584
585
586
587
			}
			}

			if ( !(ndone % 10) ) {
588
				hdf5_write(powder_fn, powder,
589
				           image.width, image.height,
590
				           H5T_NATIVE_DOUBLE);
591
592
593
			}
		}

594
595
596
597
		if ( !config_noimages ) {

			char filename[1024];

598
599
600
601
602
603
604
			if ( n_images != 1 ) {
				snprintf(filename, 1023, "%s-%i.h5",
				         outfile, number);
			} else {
				strncpy(filename, outfile, 1023);
			}

605
			number++;
Thomas White's avatar
Thomas White committed
606

607
608
			/* Write the output file */
			hdf5_write(filename, image.data,
609
			           image.width, image.height, H5T_NATIVE_FLOAT);
610
611

		}
612

613
614
615
		/* Clean up */
		free(image.data);
		free(image.twotheta);
Thomas White's avatar
Thomas White committed
616
		cell_free(cell);
Thomas White's avatar
Thomas White committed
617

618
skip:
Thomas White's avatar
Thomas White committed
619
		ndone++;
620

621
		if ( n_images && (ndone >= n_images) ) done = 1;
Thomas White's avatar
Thomas White committed
622

623
	} while ( !done );
624

625
626
627
628
	if ( gctx != NULL ) {
		cleanup_gpu(gctx);
	}

Thomas White's avatar
Thomas White committed
629
630
	free(image.det->panels);
	free(image.det);
631
	free(image.beam);
632
	free(powder);
Thomas White's avatar
Thomas White committed
633
	cell_free(input_cell);
634
	free(intensities);
Thomas White's avatar
Thomas White committed
635
636
	free(outfile);
	free(filename);
637
	free(sym);
638

639
	return 0;
640
}