From f7be30c4e727ad517bd3bf1ee3dde914c01910ea Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 8 Oct 2010 17:15:32 +0200 Subject: facetron: Multithread the integration, as well --- src/facetron.c | 190 +++++++++++++++++++++++++++++++++++++++++++-------------- 1 file changed, 143 insertions(+), 47 deletions(-) diff --git a/src/facetron.c b/src/facetron.c index 486cde6d..55acb645 100644 --- a/src/facetron.c +++ b/src/facetron.c @@ -30,6 +30,7 @@ #include "reflections.h" #include "stream.h" #include "geometry.h" +#include "peaks.h" #define MAX_THREADS (256) @@ -44,10 +45,15 @@ struct process_args int finish; int done; + /* Analysis routine */ + void (*func)(struct process_args *); + /* Analysis parameters */ const char *sym; + pthread_mutex_t *list_lock; /* Protects 'obs', 'i_full' and 'cts' */ ReflItemList *obs; double *i_full; + unsigned int *cts; }; @@ -73,36 +79,93 @@ static void show_help(const char *s) } -static void refine_image(struct image *image, ReflItemList *obs, double *i_full) +static void refine_image(struct process_args *pargs) { - //struct hdfile *hdfile; + /* Do, er, something. */ +} + + +static double partiality(struct image *image, + signed int h, signed int k, signed int l) +{ + return 1.0; +} + + +static void integrate_image(struct process_args *pargs) +{ + struct reflhit *spots; + int j, n; + struct hdfile *hdfile; + struct image *image = pargs->image; - image->features = NULL; image->data = NULL; image->flags = NULL; - image->hits = NULL; - image->n_hits = 0; - - /* View head-on (unit cell is tilted) */ - image->orientation.w = 1.0; - image->orientation.x = 0.0; - image->orientation.y = 0.0; - image->orientation.z = 0.0; - - //hdfile = hdfile_open(pargs->filename); - //if ( hdfile == NULL ) { - // return; - //} else if ( hdfile_set_first_image(hdfile, "/") ) { - // ERROR("Couldn't select path\n"); - // hdfile_close(hdfile); - // return; - //} - - //hdf5_read(hdfile, &image, 1); + hdfile = hdfile_open(image->filename); + if ( hdfile == NULL ) { + ERROR("Couldn't open '%s'\n", image->filename); + return; + } else if ( hdfile_set_image(hdfile, "/data/data0") ) { + ERROR("Couldn't select path\n"); + hdfile_close(hdfile); + return; + } + + if ( hdf5_read(hdfile, pargs->image, 0) ) { + ERROR("Couldn't read '%s'\n", image->filename); + hdfile_close(hdfile); + return; + } + +goto skip; + + /* Figure out which spots should appear in this pattern, + * using a large divergence and bandwidth to avoid missing + * reflection tails. */ + spots = find_intersections(image, image->indexed_cell, + image->div, image->bw, &n, 0); + + /* For each reflection, estimate the partiality */ + for ( j=0; jlist_lock); + integrate_intensity(pargs->i_full, h, k, l, i_partial); + integrate_count(pargs->cts, h, k, l, 1); + if ( !find_item(pargs->obs, h, k, l) ) { + add_item(pargs->obs, h, k, l); + } + pthread_mutex_unlock(pargs->list_lock); + + } + +skip: free(image->data); if ( image->flags != NULL ) free(image->flags); - //hdfile_close(hdfile); + hdfile_close(hdfile); + //free(spots); } @@ -115,7 +178,7 @@ static void *worker_thread(void *pargsv) int wakeup; - refine_image(pargs->image, pargs->obs, pargs->i_full); + pargs->func(pargs); pthread_mutex_lock(&pargs->control_mutex); pargs->done = 1; @@ -139,12 +202,15 @@ static void *worker_thread(void *pargsv) } -static void refine_all(struct image *images, int n_total_patterns, - struct detector *det, const char *sym, - ReflItemList *obs, double *i_full, int nthreads) +static void munch_threads(struct image *images, int n_total_patterns, + struct detector *det, const char *sym, + ReflItemList *obs, double *i_full, unsigned int *cts, + int nthreads, void (*func)(struct process_args *), + const char *text) { pthread_t workers[MAX_THREADS]; struct process_args *worker_args[MAX_THREADS]; + pthread_mutex_t list_lock = PTHREAD_MUTEX_INITIALIZER; int worker_active[MAX_THREADS]; int i; int n_done = 0; @@ -157,6 +223,11 @@ static void refine_all(struct image *images, int n_total_patterns, worker_active[i] = 0; pthread_mutex_init(&worker_args[i]->control_mutex, NULL); worker_args[i]->sym = sym; + worker_args[i]->obs = obs; + worker_args[i]->i_full = i_full; + worker_args[i]->cts = cts; + worker_args[i]->list_lock = &list_lock; + worker_args[i]->func = func; } @@ -211,7 +282,7 @@ static void refine_all(struct image *images, int n_total_patterns, pargs->done = 0; n_done++; - progress_bar(n_done, n_total_patterns, "Refining"); + progress_bar(n_done, n_total_patterns, text); /* If there are no more patterns, "done" will remain * zero, so the last pattern will not be re-counted. */ @@ -243,43 +314,53 @@ static void refine_all(struct image *images, int n_total_patterns, if ( pargs->done ) { n_done++; - progress_bar(n_done, n_total_patterns, "Refining"); + progress_bar(n_done, n_total_patterns, text); } /* else this thread was not busy */ } + + for ( i=0; iindexed_cell, - image->div, image->bw, &n, 0); + /* Divide the totals to get the means */ + for ( i=0; ih, it->k, it->l); + total /= lookup_count(cts, it->h, it->k, it->l); + set_intensity(i_full, it->h, it->k, it->l, total); } + + free(cts); } @@ -420,7 +501,7 @@ int main(int argc, char *argv[]) UnitCell *cell; char *filename; - char fnamereal[1024]; + char *fnamereal; if ( find_chunk(fh, &cell, &filename) == 1 ) { ERROR("Couldn't get all of the filenames and cells" @@ -437,11 +518,17 @@ int main(int argc, char *argv[]) free(filename); filename = tmp; } + fnamereal = malloc(1024); snprintf(fnamereal, 1023, "%s%s", prefix, filename); images[i].filename = fnamereal; images[i].div = 0.5e-3; images[i].bw = 0.001; + images[i].orientation.w = 1.0; + images[i].orientation.x = 0.0; + images[i].orientation.y = 0.0; + images[i].orientation.z = 0.0; + images[i].det = det; free(filename); @@ -452,7 +539,8 @@ int main(int argc, char *argv[]) free(prefix); /* Make initial estimates */ - estimate_full(images, n_total_patterns, det, sym, obs, i_full); + estimate_full(images, n_total_patterns, det, sym, obs, i_full, + nthreads); /* Iterate */ for ( i=0; ipanels); + free(det); + for ( i=0; i