indexamajig.c 18.1 KB
Newer Older
Thomas White's avatar
Thomas White committed
1
/*
2
 * indexamajig.c
Thomas White's avatar
Thomas White committed
3
4
5
 *
 * Find hits, index patterns, output hkl+intensity etc.
 *
Thomas White's avatar
Thomas White committed
6
 * (c) 2006-2010 Thomas White <taw@physics.org>
Thomas White's avatar
Thomas White committed
7
8
9
10
11
12
13
14
15
16
 *
 * Part of CrystFEL - crystallography with a FEL
 *
 */


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

17
#define _GNU_SOURCE 1
Thomas White's avatar
Thomas White committed
18
19
20
21
22
23
#include <stdarg.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <unistd.h>
#include <getopt.h>
Thomas White's avatar
Thomas White committed
24
#include <hdf5.h>
Thomas White's avatar
Thomas White committed
25
#include <gsl/gsl_errno.h>
26
#include <pthread.h>
Thomas White's avatar
Thomas White committed
27
#include <sys/time.h>
Thomas White's avatar
Thomas White committed
28
29
30

#include "utils.h"
#include "hdf5-file.h"
Thomas White's avatar
Thomas White committed
31
#include "index.h"
32
#include "peaks.h"
33
#include "diffraction.h"
Thomas White's avatar
Thomas White committed
34
#include "diffraction-gpu.h"
35
#include "detector.h"
36
#include "sfac.h"
Thomas White's avatar
Thomas White committed
37
#include "filters.h"
38
#include "reflections.h"
Thomas White's avatar
Thomas White committed
39
40


41
42
43
44
45
#define MAX_THREADS (96)

struct process_args
{
	char *filename;
Thomas White's avatar
Thomas White committed
46
	int id;
47
48
	pthread_mutex_t *output_mutex;  /* Protects stdout */
	pthread_mutex_t *gpu_mutex;     /* Protects "gctx" */
49
50
51
52
53
54
55
56
57
58
59
	UnitCell *cell;
	int config_cmfilter;
	int config_noisefilter;
	int config_writedrx;
	int config_dumpfound;
	int config_verbose;
	int config_alternate;
	int config_nearbragg;
	int config_gpu;
	int config_simulate;
	int config_nomatch;
60
	int config_unpolar;
61
62
63
64
65
66
67
68
69
	IndexingMethod indm;
	const double *intensities;
	const unsigned int *counts;
	struct gpu_context *gctx;
};

struct process_result
{
	int hit;
Thomas White's avatar
Thomas White committed
70
	int peaks_sane;
71
72
73
};


Thomas White's avatar
Thomas White committed
74
75
76
77
78
79
80
81
82
83
static void show_help(const char *s)
{
	printf("Syntax: %s [options]\n\n", s);
	printf(
"Process and index FEL diffraction images.\n"
"\n"
"  -h, --help              Display this help message.\n"
"\n"
"  -i, --input=<filename>  Specify file containing list of images to process.\n"
"                           '-' means stdin, which is the default.\n"
Thomas White's avatar
Thomas White committed
84
"\n"
Thomas White's avatar
Thomas White committed
85
86
87
"      --indexing=<method> Use 'method' for indexing.  Choose from:\n"
"                           none     : no indexing\n"
"                           dirax    : invoke DirAx\n"
Thomas White's avatar
Thomas White committed
88
89
"\n\nWith just the above options, this program does not do much of practical "
"use.\nYou should also enable some of the following:\n\n"
90
"      --near-bragg        Output a list of reflection intensities to stdout.\n"
Thomas White's avatar
Thomas White committed
91
92
93
94
95
96
97
98
99
100
"                           When pixels with fractional indices within 0.1 of\n"
"                           integer values (the Bragg condition) are found,\n"
"                           the integral of pixels within a ten pixel radius\n"
"                           of the nearest-to-Bragg pixel will be reported as\n"
"                           the intensity.  The centroid of the pixels will\n"
"                           be given as the coordinates, as well as the h,k,l\n"
"                           (integer) indices of the reflection.  If a peak\n"
"                           was located by the initial peak search close to\n"
"                           the \"near Bragg\" location, its coordinates will\n"
"                           be taken as the centre instead.\n"
101
"      --simulate          Simulate the diffraction pattern using the indexed\n"
Thomas White's avatar
Thomas White committed
102
103
104
105
"                           unit cell.  The simulated pattern will be saved\n"
"                           as \"simulated.h5\".  You can TRY to combine this\n"
"                           with \"-j <n>\" with n greater than 1, but it's\n"
"                           not a good idea.\n"
106
107
108
109
110
111
112
"      --filter-cm         Perform common-mode noise subtraction on images\n"
"                           before proceeding.  Intensities will be extracted\n"
"                           from the image as it is after this processing.\n"
"      --filter-noise      Apply an aggressive noise filter which sets all\n"
"                           pixels in each 3x3 region to zero if any of them\n"
"                           have negative values.  Intensity measurement will\n"
"                           be performed on the image as it was before this.\n"
113
114
115
116
117
"      --write-drx         Write 'xfel.drx' for visualisation of reciprocal\n"
"                           space.  Implied by any indexing method other than\n"
"                           'none'.  Beware: the units in this file are\n"
"                           reciprocal Angstroms.\n"
"      --dump-peaks        Write the results of the peak search to stdout.\n"
118
119
"                           The intensities in this list are from the\n"
"                           centroid/integration procedure.\n"
120
121
122
"      --no-match          Don't attempt to match the indexed cell to the\n"
"                           model, just proceed with the one generated by the\n"
"                           auto-indexing procedure.\n"
123
"      --unpolarized       Don't correct for the polarisation of the X-rays.\n"
Thomas White's avatar
Thomas White committed
124
125
126
127
128
"\n\nOptions for greater performance or verbosity:\n\n"
"      --verbose           Be verbose about indexing.\n"
"      --gpu               Use the GPU to speed up the simulation.\n"
"  -j <n>                  Run <n> analyses in parallel.  Default 1.\n"
"\n\nControl of model and data input:\n\n"
Thomas White's avatar
Thomas White committed
129
130
131
"     --intensities=<file> Specify file containing reflection intensities\n"
"                           to use when simulating.\n"
" -p, --pdb=<file>         PDB file from which to get the unit cell to match.\n"
132
" -x, --prefix=<p>         Prefix filenames from input file with 'p'.\n"
Thomas White's avatar
Thomas White committed
133
);
Thomas White's avatar
Thomas White committed
134
135
136
}


137
static struct image *get_simage(struct image *template, int alternate)
138
{
139
	struct image *image;
140
	struct panel panels[2];
141

142
143
	image = malloc(sizeof(*image));

144
	/* Simulate a diffraction pattern */
145
146
147
	image->twotheta = NULL;
	image->data = NULL;
	image->det = template->det;
148
149

	/* View head-on (unit cell is tilted) */
150
151
152
153
	image->orientation.w = 1.0;
	image->orientation.x = 0.0;
	image->orientation.y = 0.0;
	image->orientation.z = 0.0;
154
155
156

	/* Detector geometry for the simulation
	 * - not necessarily the same as the original. */
157
158
159
	image->width = 1024;
	image->height = 1024;
	image->det.n_panels = 2;
160

161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
	if ( alternate ) {

		/* Upper */
		panels[0].min_x = 0;
		panels[0].max_x = 1023;
		panels[0].min_y = 512;
		panels[0].max_y = 1023;
		panels[0].cx = 523.6;
		panels[0].cy = 502.5;
		panels[0].clen = 56.4e-2;  /* 56.4 cm */
		panels[0].res = 13333.3;   /* 75 microns/pixel */

		/* Lower */
		panels[1].min_x = 0;
		panels[1].max_x = 1023;
		panels[1].min_y = 0;
		panels[1].max_y = 511;
		panels[1].cx = 520.8;
179
		panels[1].cy = 525.0;
180
181
182
183
		panels[1].clen = 56.7e-2;  /* 56.7 cm */
		panels[1].res = 13333.3;   /* 75 microns/pixel */

		image->det.panels = panels;
184

185
186
187
188
189
190
	} else {

		/* Copy pointer to old geometry */
		image->det.panels = template->det.panels;

	}
191

192
	image->lambda = ph_en_to_lambda(eV_to_J(1.8e3));
193
	image->features = template->features;
194
	image->filename = template->filename;
195
	image->indexed_cell = template->indexed_cell;
Thomas White's avatar
Thomas White committed
196
	image->f0 = template->f0;
197

Thomas White's avatar
Thomas White committed
198
199
200
201
	/* Prevent muppetry */
	image->hits = NULL;
	image->n_hits = 0;

202
203
	return image;
}
204
205


206
static void simulate_and_write(struct image *simage, struct gpu_context **gctx,
207
208
                               const double *intensities,
                               const unsigned int *counts, UnitCell *cell)
209
{
210
	/* Set up GPU if necessary */
Thomas White's avatar
Thomas White committed
211
	if ( (gctx != NULL) && (*gctx == NULL) ) {
212
		*gctx = setup_gpu(0, simage, intensities, counts);
213
214
	}

Thomas White's avatar
Thomas White committed
215
	if ( (gctx != NULL) && (*gctx != NULL) ) {
216
		get_diffraction_gpu(*gctx, simage, 24, 24, 40, cell);
217
	} else {
218
		get_diffraction(simage, 24, 24, 40,
219
		                intensities, counts, cell, 0, GRADIENT_MOSAIC);
220
	}
221
	record_image(simage, 0);
222

223
	hdf5_write("simulated.h5", simage->data, simage->width, simage->height,
224
225
226
227
		   H5T_NATIVE_FLOAT);
}


228
static void *process_image(void *pargsv)
229
{
230
	struct process_args *pargs = pargsv;
231
232
233
234
235
	struct hdfile *hdfile;
	struct image image;
	struct image *simage;
	float *data_for_measurement;
	size_t data_size;
236
237
238
239
240
241
242
243
244
245
246
247
	const char *filename = pargs->filename;
	UnitCell *cell = pargs->cell;
	int config_cmfilter = pargs->config_cmfilter;
	int config_noisefilter = pargs->config_noisefilter;
	int config_writedrx = pargs->config_writedrx;
	int config_dumpfound = pargs->config_dumpfound;
	int config_verbose = pargs->config_verbose;
	int config_alternate  = pargs->config_alternate;
	int config_nearbragg = pargs->config_nearbragg;
	int config_gpu = pargs->config_gpu;
	int config_simulate = pargs->config_simulate;
	int config_nomatch = pargs->config_nomatch;
248
	int config_unpolar = pargs->config_unpolar;
249
250
251
252
253
	IndexingMethod indm = pargs->indm;
	const double *intensities = pargs->intensities;
	const unsigned int *counts = pargs->counts;
	struct gpu_context *gctx = pargs->gctx;
	struct process_result *result;
254
255
256
257

	image.features = NULL;
	image.data = NULL;
	image.indexed_cell = NULL;
Thomas White's avatar
Thomas White committed
258
	image.id = pargs->id;
259
	image.filename = filename;
Thomas White's avatar
Thomas White committed
260
261
	image.hits = NULL;
	image.n_hits = 0;
262

263
	STATUS("Processing '%s'\n", image.filename);
264

265
266
267
	result = malloc(sizeof(*result));
	if ( result == NULL ) return NULL;

268
269
	hdfile = hdfile_open(filename);
	if ( hdfile == NULL ) {
270
271
		result->hit = 0;
		return result;
272
273
	} else if ( hdfile_set_first_image(hdfile, "/") ) {
		ERROR("Couldn't select path\n");
274
275
		result->hit = 0;
		return result;
276
277
	}

278
279
	#include "geometry-lcls.tmp"

280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
	hdf5_read(hdfile, &image);

	if ( config_cmfilter ) {
		filter_cm(&image);
	}

	/* Take snapshot of image after CM subtraction but before
	 * the aggressive noise filter. */
	data_size = image.width*image.height*sizeof(float);
	data_for_measurement = malloc(data_size);

	if ( config_noisefilter ) {
		filter_noise(&image, data_for_measurement);
	} else {

		int x, y;

		for ( x=0; x<image.width; x++ ) {
		for ( y=0; y<image.height; y++ ) {
			float val;
			val = image.data[x+image.width*y];
			data_for_measurement[x+image.width*y] = val;
		}
		}

	}

	/* Perform 'fine' peak search */
	search_peaks(&image);
	if ( image_feature_count(image.features) < 5 ) goto done;

311
	if ( config_dumpfound ) dump_peaks(&image, pargs->output_mutex);
312
313
314
315
316
317
318
319
320
321
322
323
324
325

	/* Not indexing nor writing xfel.drx?
	 * Then there's nothing left to do. */
	if ( (!config_writedrx) && (indm == INDEXING_NONE) ) {
		goto done;
	}

	/* Calculate orientation matrix (by magic) */
	if ( config_writedrx || (indm != INDEXING_NONE) ) {
		index_pattern(&image, cell, indm, config_nomatch,
		              config_verbose);
	}

	/* No cell at this point?  Then we're done. */
Thomas White's avatar
Thomas White committed
326
	result->peaks_sane = 0;
327
328
	if ( image.indexed_cell == NULL ) goto done;

Thomas White's avatar
Thomas White committed
329
330
331
332
333
334
335
336
337
	/* Sanity check */
	if ( !peak_sanity_check(&image, image.indexed_cell) ) {
		STATUS("Failed peak sanity check.\n");
		result->peaks_sane = 0;
		goto done;
	} else {
		result->peaks_sane = 1;
	}

Thomas White's avatar
Thomas White committed
338
339
340
341
	/* Get rid of noise-filtered version at this point */
	free(image.data);
	image.data = data_for_measurement;

342
343
	/* Measure intensities if requested */
	if ( config_nearbragg ) {
344
		output_intensities(&image, image.indexed_cell,
345
		                   pargs->output_mutex, config_unpolar);
346
347
	}

348
349
	simage = get_simage(&image, config_alternate);

350
351
352
	/* Simulate if requested */
	if ( config_simulate ) {
		if ( config_gpu ) {
353
			pthread_mutex_lock(pargs->gpu_mutex);
354
			simulate_and_write(simage, &gctx, intensities,
355
			                   counts, image.indexed_cell);
356
			pthread_mutex_unlock(pargs->gpu_mutex);
357
358
		} else {
			simulate_and_write(simage, NULL, intensities,
359
			                   counts, image.indexed_cell);
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
		}
	}

	/* Finished with alternate image */
	if ( simage->twotheta != NULL ) free(simage->twotheta);
	if ( simage->data != NULL ) free(simage->data);
	free(simage);

	/* Only free cell if found */
	free(image.indexed_cell);

done:
	free(image.data);
	free(image.det.panels);
	image_feature_list_free(image.features);
Thomas White's avatar
Thomas White committed
375
	free(image.hits);
376
377
	hdfile_close(hdfile);

378
379
380
381
382
383
	if ( image.indexed_cell == NULL ) {
		result->hit = 0;
	} else {
		result->hit = 1;
	}
	return result;
384
385
386
}


Thomas White's avatar
Thomas White committed
387
388
389
int main(int argc, char *argv[])
{
	int c;
Thomas White's avatar
Thomas White committed
390
	struct gpu_context *gctx = NULL;
Thomas White's avatar
Thomas White committed
391
392
	char *filename = NULL;
	FILE *fh;
393
	char *rval = NULL;
Thomas White's avatar
Thomas White committed
394
395
	int n_images;
	int n_hits;
Thomas White's avatar
Thomas White committed
396
	int n_sane;
397
	int config_noindex = 0;
Thomas White's avatar
Thomas White committed
398
	int config_dumpfound = 0;
399
	int config_nearbragg = 0;
Thomas White's avatar
Thomas White committed
400
	int config_writedrx = 0;
401
	int config_simulate = 0;
402
403
	int config_cmfilter = 0;
	int config_noisefilter = 0;
404
	int config_nomatch = 0;
Thomas White's avatar
Thomas White committed
405
	int config_gpu = 0;
406
	int config_verbose = 0;
407
	int config_alternate = 0;
408
	int config_unpolar = 0;
Thomas White's avatar
Thomas White committed
409
410
	IndexingMethod indm;
	char *indm_str = NULL;
411
412
413
	UnitCell *cell;
	double *intensities = NULL;
	char *intfile = NULL;
414
	unsigned int *counts;
Thomas White's avatar
Thomas White committed
415
	char *pdb = NULL;
416
	char *prefix = NULL;
417
418
419
	int nthreads = 1;
	pthread_t workers[MAX_THREADS];
	struct process_args *worker_args[MAX_THREADS];
Thomas White's avatar
Thomas White committed
420
	int worker_active[MAX_THREADS];
421
	int i;
422
	pthread_mutex_t output_mutex = PTHREAD_MUTEX_INITIALIZER;
423
	pthread_mutex_t gpu_mutex = PTHREAD_MUTEX_INITIALIZER;
Thomas White's avatar
Thomas White committed
424
425
426
427
428

	/* Long options */
	const struct option longopts[] = {
		{"help",               0, NULL,               'h'},
		{"input",              1, NULL,               'i'},
Thomas White's avatar
Thomas White committed
429
		{"gpu",                0, &config_gpu,         1},
430
		{"no-index",           0, &config_noindex,     1},
Thomas White's avatar
Thomas White committed
431
		{"dump-peaks",         0, &config_dumpfound,   1},
432
		{"near-bragg",         0, &config_nearbragg,   1},
Thomas White's avatar
Thomas White committed
433
434
		{"write-drx",          0, &config_writedrx,    1},
		{"indexing",           1, NULL,               'z'},
435
		{"simulate",           0, &config_simulate,    1},
436
437
		{"filter-cm",          0, &config_cmfilter,    1},
		{"filter-noise",       0, &config_noisefilter, 1},
438
		{"no-match",           0, &config_nomatch,     1},
439
		{"verbose",            0, &config_verbose,     1},
440
		{"alternate",          0, &config_alternate,   1},
441
		{"intensities",        1, NULL,               'q'},
Thomas White's avatar
Thomas White committed
442
		{"pdb",                1, NULL,               'p'},
443
		{"prefix",             1, NULL,               'x'},
444
		{"unpolarized",        0, &config_unpolar,     1},
Thomas White's avatar
Thomas White committed
445
446
447
448
		{0, 0, NULL, 0}
	};

	/* Short options */
449
	while ((c = getopt_long(argc, argv, "hi:wp:j:x:", longopts, NULL)) != -1) {
Thomas White's avatar
Thomas White committed
450
451
452
453
454
455
456
457
458
459
460
461

		switch (c) {
		case 'h' : {
			show_help(argv[0]);
			return 0;
		}

		case 'i' : {
			filename = strdup(optarg);
			break;
		}

Thomas White's avatar
Thomas White committed
462
463
464
465
466
		case 'z' : {
			indm_str = strdup(optarg);
			break;
		}

467
468
469
470
471
		case 'q' : {
			intfile = strdup(optarg);
			break;
		}

Thomas White's avatar
Thomas White committed
472
473
474
475
476
		case 'p' : {
			pdb = strdup(optarg);
			break;
		}

477
478
479
480
481
		case 'x' : {
			prefix = strdup(optarg);
			break;
		}

482
483
484
485
486
		case 'j' : {
			nthreads = atoi(optarg);
			break;
		}

Thomas White's avatar
Thomas White committed
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
		case 0 : {
			break;
		}

		default : {
			return 1;
		}
		}

	}

	if ( filename == NULL ) {
		filename = strdup("-");
	}
	if ( strcmp(filename, "-") == 0 ) {
		fh = stdin;
	} else {
		fh = fopen(filename, "r");
	}
	if ( fh == NULL ) {
Thomas White's avatar
Thomas White committed
507
		ERROR("Failed to open input file '%s'\n", filename);
Thomas White's avatar
Thomas White committed
508
509
		return 1;
	}
Thomas White's avatar
Thomas White committed
510
	free(filename);
Thomas White's avatar
Thomas White committed
511

512
	if ( intfile != NULL ) {
513
		counts = new_list_count();
514
		intensities = read_reflections(intfile, counts);
515
516
	} else {
		intensities = NULL;
517
		counts = NULL;
518
519
	}

Thomas White's avatar
Thomas White committed
520
521
522
523
	if ( pdb == NULL ) {
		pdb = strdup("molecule.pdb");
	}

524
	if ( prefix == NULL ) {
Thomas White's avatar
Thomas White committed
525
		prefix = strdup("");
526
527
	}

528
	if ( (nthreads == 0) || (nthreads > MAX_THREADS) ) {
529
530
531
532
		ERROR("Invalid number of threads.\n");
		return 1;
	}

Thomas White's avatar
Thomas White committed
533
534
	if ( indm_str == NULL ) {
		STATUS("You didn't specify an indexing method, so I won't"
535
536
537
		       " try to index anything.\n"
		       "If that isn't what you wanted, re-run with"
		       " --indexing=<method>.\n");
Thomas White's avatar
Thomas White committed
538
539
540
541
542
543
544
545
546
		indm = INDEXING_NONE;
	} else if ( strcmp(indm_str, "none") == 0 ) {
		indm = INDEXING_NONE;
	} else if ( strcmp(indm_str, "dirax") == 0) {
		indm = INDEXING_DIRAX;
	} else {
		ERROR("Unrecognised indexing method '%s'\n", indm_str);
		return 1;
	}
547
	free(indm_str);
Thomas White's avatar
Thomas White committed
548

Thomas White's avatar
Thomas White committed
549
	cell = load_cell_from_pdb(pdb);
550
	if ( cell == NULL ) {
551
552
553
554
555
		if ( pdb == NULL ) {
			ERROR("Couldn't read unit cell (from molecule.pdb)\n");
		} else {
			ERROR("Couldn't read unit cell (from %s)\n", pdb);
		}
556
557
		return 1;
	}
558
	free(pdb);
559

Thomas White's avatar
Thomas White committed
560
	gsl_set_error_handler_off();
Thomas White's avatar
Thomas White committed
561
562
	n_images = 0;
	n_hits = 0;
Thomas White's avatar
Thomas White committed
563
	n_sane = 0;
564

565
566
567
568
569
570
	for ( i=0; i<nthreads; i++ ) {
		worker_args[i] = malloc(sizeof(struct process_args));
		worker_args[i]->filename = malloc(1024);
		worker_active[i] = 0;
	}

571
572
	/* Initially, fire off the full number of threads */
	for ( i=0; i<nthreads; i++ ) {
Thomas White's avatar
Thomas White committed
573
574

		char line[1024];
575
576
		struct process_args *pargs;
		int r;
Thomas White's avatar
Thomas White committed
577

578
		pargs = worker_args[i];
Thomas White's avatar
Thomas White committed
579

Thomas White's avatar
Thomas White committed
580
		rval = fgets(line, 1023, fh);
Thomas White's avatar
Thomas White committed
581
		if ( rval == NULL ) continue;
Thomas White's avatar
Thomas White committed
582
		chomp(line);
583
		snprintf(pargs->filename, 1023, "%s%s", prefix, line);
584

585
586
		n_images++;

587
		pargs->output_mutex = &output_mutex;
588
		pargs->gpu_mutex = &gpu_mutex;
589
590
591
592
593
594
595
596
597
598
		pargs->config_cmfilter = config_cmfilter;
		pargs->config_noisefilter = config_noisefilter;
		pargs->config_writedrx = config_writedrx;
		pargs->config_dumpfound = config_dumpfound;
		pargs->config_verbose = config_verbose;
		pargs->config_alternate = config_alternate;
		pargs->config_nearbragg = config_nearbragg;
		pargs->config_gpu = config_gpu;
		pargs->config_simulate = config_simulate;
		pargs->config_nomatch = config_nomatch;
599
		pargs->config_unpolar = config_unpolar;
600
601
602
603
604
		pargs->cell = cell;
		pargs->indm = indm;
		pargs->intensities = intensities;
		pargs->counts = counts;
		pargs->gctx = gctx;
Thomas White's avatar
Thomas White committed
605
		pargs->id = i;
606

Thomas White's avatar
Thomas White committed
607
		worker_active[i] = 1;
608
609
		r = pthread_create(&workers[i], NULL, process_image, pargs);
		if ( r != 0 ) {
Thomas White's avatar
Thomas White committed
610
			worker_active[i] = 0;
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
			ERROR("Couldn't start thread %i\n", i);
		}

	}

	/* Start new threads as old ones finish */
	do {

		int i;

		for ( i=0; i<nthreads; i++ ) {

			char line[1024];
			int r;
			struct process_result *result = NULL;
			struct timespec t;
Thomas White's avatar
Thomas White committed
627
			struct timeval tv;
628
629
			struct process_args *pargs;

Thomas White's avatar
Thomas White committed
630
631
			if ( !worker_active[i] ) continue;

632
633
			pargs = worker_args[i];

Thomas White's avatar
Thomas White committed
634
635
636
			gettimeofday(&tv, NULL);
			t.tv_sec = tv.tv_sec;
			t.tv_nsec = tv.tv_usec * 1000 + 20000;
637
638
639
640
641

			r = pthread_timedjoin_np(workers[i], (void *)&result,
			                         &t);
			if ( r != 0 ) continue; /* Not ready yet */

Thomas White's avatar
Thomas White committed
642
643
			worker_active[i] = 0;

644
645
			if ( result != NULL ) {
				n_hits += result->hit;
Thomas White's avatar
Thomas White committed
646
				n_sane += result->peaks_sane;
647
648
649
650
651
652
				free(result);
			}

			rval = fgets(line, 1023, fh);
			if ( rval == NULL ) break;
			chomp(line);
653
			snprintf(pargs->filename, 1023, "%s%s", prefix, line);
654

Thomas White's avatar
Thomas White committed
655
			worker_active[i] = 1;
656
657
658
			r = pthread_create(&workers[i], NULL, process_image,
			                   pargs);
			if ( r != 0 ) {
Thomas White's avatar
Thomas White committed
659
				worker_active[i] = 0;
660
661
662
663
664
				ERROR("Couldn't start thread %i\n", i);
			}

			n_images++;
		}
Thomas White's avatar
Thomas White committed
665

Thomas White's avatar
Thomas White committed
666
667
	} while ( rval != NULL );

Thomas White's avatar
Thomas White committed
668
	/* Catch all remaining threads */
669
	for ( i=0; i<nthreads; i++ ) {
Thomas White's avatar
Thomas White committed
670
671
672

		struct process_result *result = NULL;

Thomas White's avatar
Thomas White committed
673
		if ( !worker_active[i] ) goto free;
Thomas White's avatar
Thomas White committed
674
675
676
677
678
679
680
681
682
683

		pthread_join(workers[i], (void *)&result);

		worker_active[i] = 0;

		if ( result != NULL ) {
			n_hits += result->hit;
			free(result);
		}

Thomas White's avatar
Thomas White committed
684
	free:
685
		if ( worker_args[i]->filename != NULL ) {
686
687
			free(worker_args[i]->filename);
		}
688
		free(worker_args[i]);
Thomas White's avatar
Thomas White committed
689

690
691
692
693
	}

	free(prefix);
	free(cell);
Thomas White's avatar
Thomas White committed
694
695
696
	fclose(fh);

	STATUS("There were %i images.\n", n_images);
Thomas White's avatar
Thomas White committed
697
	STATUS("%i hits were found, of which %i were sane.\n", n_hits, n_sane);
Thomas White's avatar
Thomas White committed
698

Thomas White's avatar
Thomas White committed
699
700
701
702
	if ( gctx != NULL ) {
		cleanup_gpu(gctx);
	}

Thomas White's avatar
Thomas White committed
703
704
	return 0;
}