aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2015-12-05 13:35:56 -0800
committerThomas White <taw@physics.org>2015-12-05 13:47:36 -0800
commitf094e70fd60899579a8452c1c08121f25e0be885 (patch)
tree379f0a5ec91694c2a4df26a76e3dbe600bd6fdf3 /libcrystfel
parent3fe4c04aacd623f76dcb35114a02b2b001900d4e (diff)
Add option for per-pixel saturation values
Diffstat (limited to 'libcrystfel')
-rw-r--r--libcrystfel/src/detector.c7
-rw-r--r--libcrystfel/src/detector.h2
-rw-r--r--libcrystfel/src/hdf5-file.c116
-rw-r--r--libcrystfel/src/image.h1
-rw-r--r--libcrystfel/src/integration.c7
5 files changed, 127 insertions, 6 deletions
diff --git a/libcrystfel/src/detector.c b/libcrystfel/src/detector.c
index cb0bff95..7689f5f8 100644
--- a/libcrystfel/src/detector.c
+++ b/libcrystfel/src/detector.c
@@ -956,6 +956,11 @@ static int parse_field_for_panel(struct panel *panel, const char *key,
} else if ( strcmp(key, "mask_file") == 0 ) {
panel->mask_file = strdup(val);
+ } else if ( strcmp(key, "saturation_map") == 0 ) {
+ panel->satmap = strdup(val);
+ } else if ( strcmp(key, "saturation_map_file") == 0 ) {
+ panel->satmap_file = strdup(val);
+
} else if ( strcmp(key, "coffset") == 0) {
panel->coffset = atof(val);
} else if ( strcmp(key, "res") == 0 ) {
@@ -1263,6 +1268,8 @@ struct detector *get_detector_geometry(const char *filename,
det->defaults.max_adu = +INFINITY;
det->defaults.mask = NULL;
det->defaults.mask_file = NULL;
+ det->defaults.satmap = NULL;
+ det->defaults.satmap_file = NULL;
det->defaults.data = NULL;
det->defaults.dim_structure = NULL;
strncpy(det->defaults.name, "", 1023);
diff --git a/libcrystfel/src/detector.h b/libcrystfel/src/detector.h
index 5e544b75..540db4b6 100644
--- a/libcrystfel/src/detector.h
+++ b/libcrystfel/src/detector.h
@@ -97,6 +97,8 @@ struct panel
char *clen_from;
char *mask;
char *mask_file;
+ char *satmap;
+ char *satmap_file;
double res; /* Resolution in pixels per metre */
char badrow; /* 'x' or 'y' */
int no_index; /* Don't index peaks in this panel if non-zero */
diff --git a/libcrystfel/src/hdf5-file.c b/libcrystfel/src/hdf5-file.c
index b159327e..09a34c31 100644
--- a/libcrystfel/src/hdf5-file.c
+++ b/libcrystfel/src/hdf5-file.c
@@ -1126,13 +1126,16 @@ static void debodge_saturation(struct hdfile *f, struct image *image)
static int unpack_panels(struct image *image, struct detector *det,
- float *data, uint16_t *flags)
+ float *data, uint16_t *flags, float *sat)
{
int pi;
image->dp = malloc(det->n_panels * sizeof(float *));
image->bad = malloc(det->n_panels * sizeof(int *));
- if ( (image->dp == NULL) || (image->bad == NULL) ) {
+ image->sat = malloc(det->n_panels * sizeof(float *));
+ if ( (image->dp == NULL) || (image->bad == NULL)
+ || (image->sat == NULL) )
+ {
ERROR("Failed to allocate panels.\n");
return 1;
}
@@ -1145,7 +1148,10 @@ static int unpack_panels(struct image *image, struct detector *det,
p = &det->panels[pi];
image->dp[pi] = malloc(p->w*p->h*sizeof(float));
image->bad[pi] = calloc(p->w*p->h, sizeof(int));
- if ( (image->dp[pi] == NULL) || (image->bad[pi] == NULL) ) {
+ image->sat[pi] = malloc(p->w*p->h*sizeof(float));
+ if ( (image->dp[pi] == NULL) || (image->bad[pi] == NULL)
+ || (image->sat[pi] == NULL) )
+ {
ERROR("Failed to allocate panel\n");
return 1;
}
@@ -1163,6 +1169,10 @@ static int unpack_panels(struct image *image, struct detector *det,
image->dp[pi][fs+p->w*ss] = data[idx];
+ if ( sat != NULL ) {
+ image->sat[pi][fs+p->w*ss] = sat[idx];
+ }
+
if ( p->no_index ) bad = 1;
if ( in_bad_region(det, cfs, css) ) {
@@ -1452,7 +1462,7 @@ int hdf5_read(struct hdfile *f, struct image *image, const char *element,
}
image->det = simple_geometry(image);
- unpack_panels(image, image->det, buf, NULL);
+ unpack_panels(image, image->det, buf, NULL, NULL);
if ( satcorr ) debodge_saturation(f, image);
if ( image->beam != NULL ) {
@@ -1551,15 +1561,95 @@ err:
}
+static int load_satmap(struct hdfile *f, struct event *ev, char *satmap,
+ const char *satmap_file,
+ const char *pname, struct image *image,
+ size_t p_w, size_t sum_p_h, float *smap,
+ hsize_t *f_offset, hsize_t *f_count,
+ hsize_t *m_offset, hsize_t *m_count)
+{
+ hid_t satmap_dataspace, satmap_dh;
+ int exists;
+ int check, r;
+ hid_t memspace;
+ hsize_t dimsm[2];
+ hid_t fh;
+
+ if ( satmap_file != NULL ) {
+ fh = H5Fopen(satmap_file, H5F_ACC_RDONLY, H5P_DEFAULT);
+ if ( fh < 0 ) {
+ ERROR("Couldn't open satmap file '%s'\n", satmap_file);
+ return 1;
+ }
+ } else {
+ fh = f->fh;
+ }
+
+ if ( ev != NULL ) {
+ satmap = retrieve_full_path(ev, satmap);
+ }
+
+ exists = check_path_existence(fh, satmap);
+ if ( !exists ) {
+ ERROR("Cannot find satmap for panel %s\n", pname);
+ goto err;
+ }
+
+ satmap_dh = H5Dopen2(fh, satmap, H5P_DEFAULT);
+ if ( satmap_dh <= 0 ) {
+ ERROR("Couldn't open satmap for panel %s\n", pname);
+ goto err;
+ }
+
+ satmap_dataspace = H5Dget_space(satmap_dh);
+ check = H5Sselect_hyperslab(satmap_dataspace, H5S_SELECT_SET,
+ f_offset, NULL, f_count, NULL);
+ if ( check < 0 ) {
+ ERROR("Error selecting satmap dataspace for panel %s\n", pname);
+ goto err;
+ }
+
+ dimsm[0] = sum_p_h;
+ dimsm[1] = p_w;
+ memspace = H5Screate_simple(2, dimsm, NULL);
+ check = H5Sselect_hyperslab(memspace, H5S_SELECT_SET,
+ m_offset, NULL, m_count, NULL);
+ if ( check < 0 ) {
+ ERROR("Error selecting memory dataspace for panel %s\n", pname);
+ goto err;
+ }
+
+ r = H5Dread(satmap_dh, H5T_NATIVE_FLOAT, memspace,
+ satmap_dataspace, H5P_DEFAULT, smap);
+ if ( r < 0 ) {
+ ERROR("Couldn't read satmap for panel %s\n", pname);
+ goto err;
+ }
+
+ H5Sclose(satmap_dataspace);
+ H5Dclose(satmap_dh);
+ if ( ev != NULL ) free(satmap);
+
+ return 0;
+
+err:
+ if ( satmap_file != NULL ) H5Fclose(fh);
+ if ( ev != NULL ) free(satmap);
+ return 1;
+}
+
+
int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
int satcorr)
{
herr_t r;
float *buf;
uint16_t *flags;
+ float *smap;
int sum_p_h;
int p_w;
int pi;
+ int i;
if ( image->det == NULL ) {
ERROR("Geometry not available\n");
@@ -1594,6 +1684,12 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
return 1;
}
+ smap = calloc(p_w*sum_p_h,sizeof(float));
+ if ( smap == NULL ) {
+ ERROR("Failed to allocate memory for satmap\n");
+ return 1;
+ }
+ for ( i=0; i<p_w*sum_p_h; i++ ) smap[i] = INFINITY;
for ( pi=0; pi<image->det->n_panels; pi++ ) {
@@ -1742,6 +1838,15 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
}
}
+ if ( p->satmap != NULL ) {
+ if ( load_satmap(f, ev, p->satmap, p->satmap_file,
+ p->name, image, p_w, sum_p_h, smap,
+ f_offset, f_count, m_offset, m_count) )
+ {
+ ERROR("Error loading saturation map!\n");
+ }
+ }
+
free(f_offset);
free(f_count);
@@ -1749,7 +1854,7 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
fill_in_values(image->det, f, ev);
- unpack_panels(image, image->det, buf, flags);
+ unpack_panels(image, image->det, buf, flags, smap);
if ( satcorr ) debodge_saturation(f, image);
if ( image->beam != NULL ) {
@@ -1768,6 +1873,7 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
free(buf);
free(flags);
+ free(smap);
return 0;
}
diff --git a/libcrystfel/src/image.h b/libcrystfel/src/image.h
index 52614d06..37904308 100644
--- a/libcrystfel/src/image.h
+++ b/libcrystfel/src/image.h
@@ -161,6 +161,7 @@ struct image {
float **dp; /* Data in panel */
int **bad; /* Bad pixels by panel */
+ float **sat; /* Per-pixel saturation values */
Crystal **crystals;
int n_crystals;
diff --git a/libcrystfel/src/integration.c b/libcrystfel/src/integration.c
index 730a14ae..b3bee120 100644
--- a/libcrystfel/src/integration.c
+++ b/libcrystfel/src/integration.c
@@ -804,6 +804,7 @@ static int check_box(struct intcontext *ic, struct peak_box *bx, int *sat)
double hd, kd, ld;
signed int h, k, l;
struct rvec dv;
+ float lsat;
fs = bx->cfs + p;
ss = bx->css + q;
@@ -844,9 +845,13 @@ static int check_box(struct intcontext *ic, struct peak_box *bx, int *sat)
}
}
+ /* Per-pixel saturation value */
+ lsat = ic->image->sat[bx->pn][fs + bx->p->w*ss];
if ( (bx->bm[p+ic->w*q] != BM_IG)
&& (bx->bm[p+ic->w*q] != BM_BH)
- && (boxi(ic, bx, p, q) > bx->p->max_adu) ) {
+ && ((boxi(ic, bx, p, q) > bx->p->max_adu)
+ || (boxi(ic, bx, p, q) > lsat)) )
+ {
if ( sat != NULL ) *sat = 1;
}