diff options
Diffstat (limited to 'src/facetron.c')
| -rw-r--r-- | src/facetron.c | 45 |
1 files changed, 33 insertions, 12 deletions
diff --git a/src/facetron.c b/src/facetron.c index 8d66d549..806ce54f 100644 --- a/src/facetron.c +++ b/src/facetron.c @@ -31,6 +31,7 @@ #include "geometry.h" #include "peaks.h" #include "thread-pool.h" +#include "beam-parameters.h" static void show_help(const char *s) @@ -45,6 +46,9 @@ static void show_help(const char *s) " (must be a file, not e.g. stdin)\n" " -o, --output=<filename> Output filename. Default: facetron.hkl.\n" " -g. --geometry=<file> Get detector geometry from file.\n" +" -b, --beam=<file> Get beam parameters from file (provides initial\n" +" values for parameters, and nominal wavelengths\n" +" if no per-shot value is found in an HDF5 file.\n" " -x, --prefix=<p> Prefix filenames from input file with <p>.\n" " --basename Remove the directory parts of the filenames.\n" " --no-check-prefix Don't attempt to correct the --prefix.\n" @@ -91,6 +95,7 @@ static void integrate_image(int mytask, void *tasks) int j, n; struct hdfile *hdfile; struct image *image = pargs->image; + double nominal_photon_energy = pargs->image->beam->photon_energy; hdfile = hdfile_open(image->filename); if ( hdfile == NULL ) { @@ -102,34 +107,31 @@ static void integrate_image(int mytask, void *tasks) return; } - if ( hdf5_read(hdfile, pargs->image, 0) ) { + /* FIXME: Nominal photon energy */ + if ( hdf5_read(hdfile, pargs->image, 0, nominal_photon_energy) ) { ERROR("Couldn't read '%s'\n", image->filename); hdfile_close(hdfile); return; } /* Figure out which spots should appear in this pattern */ - spots = find_intersections(image, image->indexed_cell, - image->div, image->bw, &n, 0); + spots = find_intersections(image, image->indexed_cell, &n, 1); /* For each reflection, estimate the partiality */ for ( j=0; j<n; j++ ) { signed int h, k, l; float i_partial; - double p; float xc, yc; + float i_full_est; h = spots[j].h; k = spots[j].k; l = spots[j].l; - /* Calculated partiality of this spot in this pattern */ - p = partiality(image, h, k, l); - /* Don't attempt to use spots with very small * partialities, since it won't be accurate. */ - if ( p < 0.1 ) continue; + if ( spots[j].p < 0.1 ) continue; /* Actual measurement of this reflection from this * pattern? */ @@ -139,8 +141,10 @@ static void integrate_image(int mytask, void *tasks) continue; } + i_full_est = i_partial * spots[j].p; + pthread_mutex_lock(pargs->list_lock); - integrate_intensity(pargs->i_full, h, k, l, i_partial); + integrate_intensity(pargs->i_full, h, k, l, i_full_est); integrate_count(pargs->cts, h, k, l, 1); if ( !find_item(pargs->obs, h, k, l) ) { add_item(pargs->obs, h, k, l); @@ -249,6 +253,7 @@ int main(int argc, char *argv[]) int n_total_patterns; struct image *images; int n_iter = 10; + struct beam_params *beam = NULL; /* Long options */ const struct option longopts[] = { @@ -256,6 +261,7 @@ int main(int argc, char *argv[]) {"input", 1, NULL, 'i'}, {"output", 1, NULL, 'o'}, {"geometry", 1, NULL, 'g'}, + {"beam", 1, NULL, 'b'}, {"prefix", 1, NULL, 'x'}, {"basename", 0, &config_basename, 1}, {"no-check-prefix", 0, &config_checkprefix, 0}, @@ -265,7 +271,7 @@ int main(int argc, char *argv[]) }; /* Short options */ - while ((c = getopt_long(argc, argv, "hi:g:x:j:y:o:", + while ((c = getopt_long(argc, argv, "hi:g:x:j:y:o:b:", longopts, NULL)) != -1) { @@ -302,6 +308,15 @@ int main(int argc, char *argv[]) n_iter = atoi(optarg); break; + case 'b' : + beam = get_beam_parameters(optarg); + if ( beam == NULL ) { + ERROR("Failed to load beam parameters" + " from '%s'\n", optarg); + return 1; + } + break; + case 0 : break; @@ -348,6 +363,11 @@ int main(int argc, char *argv[]) } free(geomfile); + if ( beam == NULL ) { + ERROR("You must provide a beam parameters file.\n"); + return 1; + } + /* Prepare for iteration */ i_full = new_list_intensity(); obs = new_items(); @@ -388,13 +408,14 @@ int main(int argc, char *argv[]) snprintf(fnamereal, 1023, "%s%s", prefix, filename); images[i].filename = fnamereal; - images[i].div = 0.5e-3; - images[i].bw = 0.001; + images[i].div = beam->divergence; + images[i].bw = beam->bandwidth; 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; + images[i].beam = beam; /* Muppet proofing */ images[i].data = NULL; |