aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2010-11-30 16:39:15 +0100
committerThomas White <taw@physics.org>2012-02-22 15:27:07 +0100
commit51cbfe0995cb1c3206841f8a33b75bd64447cb4e (patch)
tree870d247dd2d4e472d3bcbf035eb526f8a181d220
parentfc4b1ab0a554a21374b03b502c1c5a5958429be9 (diff)
facetron: Placeholders for HRS scaling
-rw-r--r--src/Makefile.am5
-rw-r--r--src/Makefile.in8
-rw-r--r--src/facetron.c141
3 files changed, 43 insertions, 111 deletions
diff --git a/src/Makefile.am b/src/Makefile.am
index 5599da98..6031e441 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -64,7 +64,7 @@ calibrate_detector_LDADD = @LIBS@
facetron_SOURCES = facetron.c cell.c hdf5-file.c utils.c detector.c peaks.c \
image.c geometry.c reflections.c stream.c thread-pool.c \
- beam-parameters.c symmetry.c post-refinement.c
+ beam-parameters.c symmetry.c post-refinement.c hrs-scaling.c
facetron_LDADD = @LIBS@
cubeit_SOURCES = cubeit.c cell.c hdf5-file.c utils.c detector.c render.c \
@@ -87,4 +87,5 @@ EXTRA_DIST = cell.h hdf5-file.h image.h utils.h diffraction.h detector.h \
render.h hdfsee.h dirax.h peaks.h index.h filters.h \
diffraction-gpu.h cl-utils.h symmetry.h \
povray.h index-priv.h geometry.h templates.h render_hkl.h \
- stream.h thread-pool.h beam-parameters.h post-refinement.h
+ stream.h thread-pool.h beam-parameters.h post-refinement.h \
+ hrs-scaling.h
diff --git a/src/Makefile.in b/src/Makefile.in
index fbdcd53c..a427a4ef 100644
--- a/src/Makefile.in
+++ b/src/Makefile.in
@@ -85,7 +85,7 @@ am_facetron_OBJECTS = facetron.$(OBJEXT) cell.$(OBJEXT) \
peaks.$(OBJEXT) image.$(OBJEXT) geometry.$(OBJEXT) \
reflections.$(OBJEXT) stream.$(OBJEXT) thread-pool.$(OBJEXT) \
beam-parameters.$(OBJEXT) symmetry.$(OBJEXT) \
- post-refinement.$(OBJEXT)
+ post-refinement.$(OBJEXT) hrs-scaling.$(OBJEXT)
facetron_OBJECTS = $(am_facetron_OBJECTS)
facetron_DEPENDENCIES =
am_get_hkl_OBJECTS = get_hkl.$(OBJEXT) sfac.$(OBJEXT) cell.$(OBJEXT) \
@@ -325,7 +325,7 @@ calibrate_detector_SOURCES = calibrate_detector.c utils.c hdf5-file.c image.c \
calibrate_detector_LDADD = @LIBS@
facetron_SOURCES = facetron.c cell.c hdf5-file.c utils.c detector.c peaks.c \
image.c geometry.c reflections.c stream.c thread-pool.c \
- beam-parameters.c symmetry.c post-refinement.c
+ beam-parameters.c symmetry.c post-refinement.c hrs-scaling.c
facetron_LDADD = @LIBS@
cubeit_SOURCES = cubeit.c cell.c hdf5-file.c utils.c detector.c render.c \
@@ -346,7 +346,8 @@ EXTRA_DIST = cell.h hdf5-file.h image.h utils.h diffraction.h detector.h \
render.h hdfsee.h dirax.h peaks.h index.h filters.h \
diffraction-gpu.h cl-utils.h symmetry.h \
povray.h index-priv.h geometry.h templates.h render_hkl.h \
- stream.h thread-pool.h beam-parameters.h post-refinement.h
+ stream.h thread-pool.h beam-parameters.h post-refinement.h \
+ hrs-scaling.h
all: all-am
@@ -487,6 +488,7 @@ distclean-compile:
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/get_hkl.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/hdf5-file.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/hdfsee.Po@am__quote@
+@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/hrs-scaling.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/image.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/index.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/indexamajig.Po@am__quote@
diff --git a/src/facetron.c b/src/facetron.c
index d14141c9..38d701e6 100644
--- a/src/facetron.c
+++ b/src/facetron.c
@@ -33,6 +33,7 @@
#include "thread-pool.h"
#include "beam-parameters.h"
#include "post-refinement.h"
+#include "hrs-scaling.h"
/* Maximum number of iterations of NLSq to do for each image per macrocycle. */
@@ -134,26 +135,39 @@ static void refine_image(int mytask, void *tasks)
}
-struct integrate_args
+static void refine_all(struct image *images, int n_total_patterns,
+ struct detector *det, const char *sym,
+ ReflItemList *obs, double *i_full, int nthreads,
+ FILE *graph, FILE *pgraph)
{
- const char *sym;
- ReflItemList *obs;
- double *i_full;
- unsigned int *cts;
- pthread_mutex_t *list_lock;
- struct image *image;
-};
+ struct refine_args *tasks;
+ int i;
+ tasks = malloc(n_total_patterns * sizeof(struct refine_args));
+ for ( i=0; i<n_total_patterns; i++ ) {
+
+ tasks[i].sym = sym;
+ tasks[i].obs = obs;
+ tasks[i].i_full = i_full;
+ tasks[i].image = &images[i];
+ tasks[i].graph = graph;
+ tasks[i].pgraph = pgraph;
+
+ }
-static void integrate_image(int mytask, void *tasks)
+ run_thread_range(n_total_patterns, nthreads, "Refining",
+ refine_image, tasks);
+
+ free(tasks);
+}
+
+
+static void integrate_image(struct image *image)
{
- struct integrate_args *all_args = tasks;
- struct integrate_args *pargs = &all_args[mytask];
struct cpeak *spots;
int j, n;
struct hdfile *hdfile;
- struct image *image = pargs->image;
- double nominal_photon_energy = pargs->image->beam->photon_energy;
+ double nominal_photon_energy = image->beam->photon_energy;
hdfile = hdfile_open(image->filename);
if ( hdfile == NULL ) {
@@ -165,7 +179,7 @@ static void integrate_image(int mytask, void *tasks)
return;
}
- if ( hdf5_read(hdfile, pargs->image, 0, nominal_photon_energy) ) {
+ if ( hdf5_read(hdfile, image, 0, nominal_photon_energy) ) {
ERROR("Couldn't read '%s'\n", image->filename);
hdfile_close(hdfile);
return;
@@ -178,10 +192,8 @@ static void integrate_image(int mytask, void *tasks)
for ( j=0; j<n; j++ ) {
signed int h, k, l;
- signed int ha, ka, la;
float i_partial;
float xc, yc;
- float i_full_est;
h = spots[j].h;
k = spots[j].k;
@@ -198,26 +210,13 @@ static void integrate_image(int mytask, void *tasks)
&xc, &yc, &i_partial, NULL, NULL, 1, 1) ) {
continue;
}
- i_partial *= image->osf;
-
- i_full_est = i_partial / spots[j].p;
-
- get_asymm(h, k, l, &ha, &ka, &la, pargs->sym);
-
- pthread_mutex_lock(pargs->list_lock);
- integrate_intensity(pargs->i_full, ha, ka, la, i_full_est);
- integrate_count(pargs->cts, ha, ka, la, 1);
- if ( !find_item(pargs->obs, ha, ka, la) ) {
- add_item(pargs->obs, ha, ka, la);
- }
- pthread_mutex_unlock(pargs->list_lock);
}
free(image->data);
if ( image->flags != NULL ) free(image->flags);
hdfile_close(hdfile);
- free(spots);
+ image->cpeaks = spots;
/* Muppet proofing */
image->data = NULL;
@@ -225,78 +224,6 @@ static void integrate_image(int mytask, void *tasks)
}
-static void refine_all(struct image *images, int n_total_patterns,
- struct detector *det, const char *sym,
- ReflItemList *obs, double *i_full, int nthreads,
- FILE *graph, FILE *pgraph)
-{
- struct refine_args *tasks;
- int i;
-
- tasks = malloc(n_total_patterns * sizeof(struct refine_args));
- for ( i=0; i<n_total_patterns; i++ ) {
-
- tasks[i].sym = sym;
- tasks[i].obs = obs;
- tasks[i].i_full = i_full;
- tasks[i].image = &images[i];
- tasks[i].graph = graph;
- tasks[i].pgraph = pgraph;
-
- }
-
- run_thread_range(n_total_patterns, nthreads, "Refining",
- refine_image, tasks);
-
- free(tasks);
-}
-
-
-static void estimate_full(struct image *images, int n_total_patterns,
- struct detector *det, const char *sym,
- ReflItemList *obs, double *i_full, unsigned int *cts,
- int nthreads)
-{
- int i;
- struct integrate_args *tasks;
- pthread_mutex_t list_lock = PTHREAD_MUTEX_INITIALIZER;
-
- clear_items(obs);
- memset(i_full, 0, LIST_SIZE*sizeof(double));
- memset(cts, 0, LIST_SIZE*sizeof(unsigned int));
-
- tasks = malloc(n_total_patterns * sizeof(struct integrate_args));
- for ( i=0; i<n_total_patterns; i++ ) {
-
- tasks[i].sym = sym;
- tasks[i].obs = obs;
- tasks[i].i_full = i_full;
- tasks[i].cts = cts;
- tasks[i].list_lock = &list_lock;
- tasks[i].image = &images[i];
-
- }
-
- run_thread_range(n_total_patterns, nthreads, "Integrating",
- integrate_image, tasks);
-
- free(tasks);
-
- /* Divide the totals to get the means */
- for ( i=0; i<num_items(obs); i++ ) {
-
- struct refl_item *it;
- double total;
-
- it = get_item(obs, i);
- total = lookup_intensity(i_full, it->h, 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);
-
- }
-}
-
-
int main(int argc, char *argv[])
{
int c;
@@ -485,6 +412,9 @@ int main(int argc, char *argv[])
images[i].data = NULL;
images[i].flags = NULL;
+ /* Get reflections from this image */
+ integrate_image(&images[i]);
+
free(filename);
progress_bar(i, n_total_patterns-1, "Loading pattern data");
@@ -496,8 +426,8 @@ int main(int argc, char *argv[])
cts = new_list_count();
/* Make initial estimates */
- estimate_full(images, n_total_patterns, det, sym, obs, i_full, cts,
- nthreads);
+ STATUS("Performing initial scaling.\n");
+ scale_intensities(images, n_total_patterns, sym);
/* Iterate */
for ( i=0; i<n_iter; i++ ) {
@@ -527,8 +457,7 @@ int main(int argc, char *argv[])
nthreads, fhg, fhp);
/* Re-estimate all the full intensities */
- estimate_full(images, n_total_patterns, det, sym, obs, i_full,
- cts, nthreads);
+ scale_intensities(images, n_total_patterns, sym);
fclose(fhg);
fclose(fhp);