aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel/src/datatemplate.c
diff options
context:
space:
mode:
Diffstat (limited to 'libcrystfel/src/datatemplate.c')
-rw-r--r--libcrystfel/src/datatemplate.c362
1 files changed, 362 insertions, 0 deletions
diff --git a/libcrystfel/src/datatemplate.c b/libcrystfel/src/datatemplate.c
index 1ffa2777..bbcf1aa5 100644
--- a/libcrystfel/src/datatemplate.c
+++ b/libcrystfel/src/datatemplate.c
@@ -36,6 +36,7 @@
#include "utils.h"
#include "datatemplate.h"
+#include "image.h"
#include "datatemplate_priv.h"
@@ -1734,3 +1735,364 @@ double data_template_get_wavelength_if_possible(const DataTemplate *dt)
return NAN;
}
}
+
+
+static int separate_value_and_units(const char *from,
+ char **pvalue,
+ char **punits)
+{
+ char *sp;
+ char *fromcpy;
+ char *unitscpy;
+
+ if ( from == NULL ) return 1;
+
+ fromcpy = strdup(from);
+ if ( fromcpy == NULL ) return 1;
+
+ sp = strchr(fromcpy, ' ');
+ if ( sp == NULL ) {
+ unitscpy = NULL;
+ } else {
+ unitscpy = strdup(sp+1);
+ sp[0] = '\0';
+ }
+
+ *pvalue = fromcpy;
+ *punits = unitscpy;
+ return 0;
+}
+
+
+/* default_scale is a value to be used if both of the following
+ * conditions are met:
+ *
+ * 1. The value is a reference to image headers/metadata,
+ * rather than a literal number.
+ * 2. No units are specified in the number.
+ *
+ * This is totally horrible. Sorry. Blame history.
+ */
+static int im_get_length(struct image *image, const char *from,
+ double default_scale, double *pval)
+{
+ char *value_str;
+ char *units;
+
+ if ( from == NULL ) return 1;
+
+ if ( separate_value_and_units(from, &value_str, &units) ) return 1;
+
+ if ( units == NULL ) {
+
+ /* No units given */
+
+ if ( convert_float(value_str, pval) == 0 ) {
+
+ /* Literal value with no units */
+ free(value_str);
+ return 0;
+
+ } else {
+
+ int r;
+ r = image_read_header_float(image, value_str, pval);
+ free(value_str);
+
+ if ( r == 0 ) {
+ /* Value read from headers with no units */
+ *pval *= default_scale;
+ return 0;
+ } else {
+ /* Failed to read value from headers */
+ return 1;
+ }
+ }
+
+ } else {
+
+ /* Units are specified */
+
+ double scale;
+
+ if ( strcmp(units, "mm") == 0 ) {
+ scale = 1e-3;
+ } else if ( strcmp(units, "m") == 0 ) {
+ scale = 1.0;
+ } else {
+ ERROR("Invalid length unit '%s'\n", units);
+ free(value_str);
+ free(units);
+ return 1;
+ }
+
+ if ( convert_float(value_str, pval) == 0 ) {
+
+ /* Literal value, units specified */
+ free(value_str);
+ free(units);
+ *pval *= scale;
+ return 0;
+
+ } else {
+
+ int r;
+ r = image_read_header_float(image, value_str, pval);
+ free(value_str);
+
+ if ( r == 0 ) {
+ /* Value read from headers, units specified */
+ *pval *= scale;
+ return 0;
+ } else {
+ /* Failed to read value from headers */
+ return 1;
+ }
+ }
+ }
+}
+
+
+static int safe_strcmp(const char *a, const char *b)
+{
+ if ( (a==NULL) && (b==NULL) ) return 0;
+ if ( (a!=NULL) && (b!=NULL) ) return strcmp(a, b);
+ return 1;
+}
+
+
+static int all_panels_reference_same_clen(const DataTemplate *dtempl)
+{
+ int i;
+ char *first_val = NULL;
+ char *first_units = NULL;
+ int fail = 0;
+
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+ struct panel_template *p = &dtempl->panels[i];
+ char *val;
+ char *units;
+ if ( separate_value_and_units(p->cnz_from, &val, &units) ) {
+ /* Parse error */
+ return 0;
+ }
+ if ( i == 0 ) {
+ first_val = val;
+ first_units = units;
+ } else {
+ if ( safe_strcmp(val, first_val) != 0 ) fail = 1;
+ if ( safe_strcmp(units, first_units) != 0 ) fail = 1;
+ free(val);
+ free(units);
+ }
+ }
+
+ free(first_val);
+ free(first_units);
+ return fail;
+}
+
+
+static int all_coffsets_small(const DataTemplate *dtempl)
+{
+ int i;
+
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+ struct panel_template *p = &dtempl->panels[i];
+ if ( p->cnz_offset > 10.0*p->pixel_pitch ) return 0;
+ }
+
+ return 1;
+}
+
+
+static int all_panels_same_clen(const DataTemplate *dtempl)
+{
+ int i;
+ double *zvals;
+ double total = 0.0;
+ double mean;
+
+ zvals = malloc(sizeof(double)*dtempl->n_panels);
+ if ( zvals == NULL ) return 0;
+
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+ struct panel_template *p = &dtempl->panels[i];
+ if ( im_get_length(NULL, p->cnz_from, 1e-3, &zvals[i]) ) {
+ /* Can't get length because it used a header reference */
+ free(zvals);
+ return 0;
+ }
+ total += zvals[i];
+ }
+
+ mean = total/dtempl->n_panels;
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+ struct panel_template *p = &dtempl->panels[i];
+ if ( fabs(zvals[i] - mean) > 10.0*p->pixel_pitch ) return 0;
+ }
+
+ free(zvals);
+
+ return 1;
+}
+
+
+static int all_panels_perpendicular_to_beam(const DataTemplate *dtempl)
+{
+ int i;
+
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+ double z_diff;
+ struct panel_template *p = &dtempl->panels[i];
+ z_diff = p->fsz*PANEL_WIDTH(p) + p->ssz*PANEL_HEIGHT(p);
+ if ( z_diff > 10.0*p->pixel_pitch ) return 0;
+ }
+ return 1;
+}
+
+
+static int detector_flat(const DataTemplate *dtempl)
+{
+ return all_panels_perpendicular_to_beam(dtempl)
+ && ( (all_panels_reference_same_clen(dtempl) && all_coffsets_small(dtempl))
+ || all_panels_same_clen(dtempl) );
+}
+
+
+struct detgeom *create_detgeom(struct image *image,
+ const DataTemplate *dtempl,
+ int two_d_only)
+{
+ struct detgeom *detgeom;
+ int i;
+
+ if ( dtempl == NULL ) {
+ ERROR("NULL data template!\n");
+ return NULL;
+ }
+
+ detgeom = malloc(sizeof(struct detgeom));
+ if ( detgeom == NULL ) return NULL;
+
+ detgeom->panels = malloc(dtempl->n_panels*sizeof(struct detgeom_panel));
+ if ( detgeom->panels == NULL ) {
+ free(detgeom);
+ return NULL;
+ }
+
+ detgeom->n_panels = dtempl->n_panels;
+
+ if ( two_d_only ) {
+ if ( !detector_flat(dtempl) ) return NULL;
+ if ( dtempl->shift_x_from != NULL ) return NULL;
+ if ( dtempl->shift_y_from != NULL ) return NULL;
+ }
+
+ for ( i=0; i<dtempl->n_panels; i++ ) {
+
+ struct detgeom_panel *p = &detgeom->panels[i];
+ struct panel_template *tmpl = &dtempl->panels[i];
+ double shift_x, shift_y;
+
+ p->name = safe_strdup(tmpl->name);
+
+ p->pixel_pitch = tmpl->pixel_pitch;
+
+ /* NB cnx,cny are in pixels, cnz is in m */
+ p->cnx = tmpl->cnx;
+ p->cny = tmpl->cny;
+
+ if ( im_get_length(image, tmpl->cnz_from, 1e-3, &p->cnz) )
+ {
+ if ( two_d_only ) {
+ p->cnz = NAN;
+ } else {
+ ERROR("Failed to read length from '%s'\n", tmpl->cnz_from);
+ return NULL;
+ }
+ }
+
+ /* Apply offset (in m) and then convert cnz from
+ * m to pixels */
+ p->cnz += tmpl->cnz_offset;
+ p->cnz /= p->pixel_pitch;
+
+ /* Apply overall shift (already in m) */
+ if ( dtempl->shift_x_from != NULL ) {
+ if ( im_get_length(image, dtempl->shift_x_from, 1.0, &shift_x) ) {
+ ERROR("Failed to read length from '%s'\n",
+ dtempl->shift_x_from);
+ return NULL;
+ }
+ if ( im_get_length(image, dtempl->shift_y_from, 1.0, &shift_y) ) {
+ ERROR("Failed to read length from '%s'\n",
+ dtempl->shift_y_from);
+ return NULL;
+ }
+ } else {
+ shift_x = 0.0;
+ shift_y = 0.0;
+ }
+
+ if ( !isnan(shift_x) ) {
+ p->cnx += shift_x / p->pixel_pitch;
+ }
+ if ( !isnan(shift_y) ) {
+ p->cny += shift_y / p->pixel_pitch;
+ }
+
+ p->max_adu = tmpl->max_adu;
+
+ switch ( tmpl->adu_scale_unit ) {
+
+ case ADU_PER_PHOTON:
+ p->adu_per_photon = tmpl->adu_scale;
+ break;
+
+ case ADU_PER_EV:
+ p->adu_per_photon = tmpl->adu_scale
+ * ph_lambda_to_eV(image->lambda);
+ break;
+
+ default:
+ p->adu_per_photon = 1.0;
+ ERROR("Invalid ADU/ph scale unit (%i)\n",
+ tmpl->adu_scale_unit);
+ break;
+
+ }
+
+ p->w = tmpl->orig_max_fs - tmpl->orig_min_fs + 1;
+ p->h = tmpl->orig_max_ss - tmpl->orig_min_ss + 1;
+
+ p->fsx = tmpl->fsx;
+ p->fsy = tmpl->fsy;
+ p->fsz = tmpl->fsz;
+ p->ssx = tmpl->ssx;
+ p->ssy = tmpl->ssy;
+ p->ssz = tmpl->ssz;
+
+ }
+
+ return detgeom;
+}
+
+
+/**
+ * Create a detgeom structure from the DataTemplate, if possible, and ignoring
+ * 3D information.
+ *
+ * This procedure will create a detgeom structure provided that the detector
+ * is close to lying in a single flat plane perpendicular to the beam
+ * direction. If certain things (e.g. panel z-positions) refer to headers,
+ * it might not be possible to determine that the detector is really flat
+ * until an image is loaded. Therefore you must gracefully handle a NULL
+ * return value from this routine.
+ *
+ * \returns the detgeom structure, or NULL if impossible.
+ */
+struct detgeom *data_template_get_2d_detgeom_if_possible(const DataTemplate *dt)
+{
+ return create_detgeom(NULL, dt, 1);
+}
generated by cgit v1.2.3 (git 2.25.1) at 2024-06-28 00:38:21 +0000