Commit 51cbfe09 authored by Thomas White's avatar Thomas White
Browse files

facetron: Placeholders for HRS scaling

parent fc4b1ab0
......@@ -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
......@@ -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@
......
......@@ -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);
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment