Merge branch 'tom/cbf'

8 files changed, 776 insertions, 64 deletions

diff --git a/libcrystfel/src/detector.c b/libcrystfel/src/detector.c
index e160d1e1..e3b351ff 100644
--- a/libcrystfel/src/detector.c
+++ b/libcrystfel/src/detector.c
@@ -528,37 +528,15 @@ int panel_number(struct detector *det, struct panel *p)
 }
 
 
-void fill_in_values(struct detector *det, struct hdfile *f, struct event* ev)
+void adjust_centering_for_rail(struct panel *p)
 {
-	int i;
-
-	for ( i=0; i<det->n_panels; i++ ) {
-
-		double offs;
-		struct panel *p = &det->panels[i];
+	double offs;
 
-		if ( p->clen_from != NULL ) {
-
-			double val;
-			int r;
-
-			r = hdfile_get_value(f, p->clen_from, ev, &val,
-			                     H5T_NATIVE_DOUBLE);
-			if ( r ) {
-				ERROR("Failed to read '%s'\n", p->clen_from);
-			} else {
-				p->clen = val * 1.0e-3;
-			}
-
-		}
-
-                /* Offset in +z direction from calibrated clen to actual */
-		offs = p->clen - p->clen_for_centering;
-		p->cnx += p->rail_x * offs;
-		p->cny += p->rail_y * offs;
-		p->clen = p->clen_for_centering + p->coffset + p->rail_z * offs;
-
-	}
+	/* Offset in +z direction from calibrated clen to actual */
+	offs = p->clen - p->clen_for_centering;
+	p->cnx += p->rail_x * offs;
+	p->cny += p->rail_y * offs;
+	p->clen = p->clen_for_centering + p->coffset + p->rail_z * offs;
 }
 
 
@@ -1076,12 +1054,15 @@ static void parse_toplevel(struct detector *det, struct beam_params *beam,
 
 	} else if ( strcmp(key, "photon_energy") == 0 ) {
 		if ( beam != NULL ) {
-			if ( strncmp(val, "/", 1) == 0 ) {
+			double v;
+			char *end;
+			v = strtod(val, &end);
+			if ( (val[0] != '\0') && (end[0] == '\0') ) {
+				beam->photon_energy = v;
+				beam->photon_energy_from = NULL;
+			} else {
 				beam->photon_energy = 0.0;
 				beam->photon_energy_from = strdup(val);
-			} else {
-				beam->photon_energy = atof(val);
-				beam->photon_energy_from = NULL;
 			}
 		}
 
diff --git a/libcrystfel/src/detector.h b/libcrystfel/src/detector.h
index 04e3c1ec..e9dd1154 100644
--- a/libcrystfel/src/detector.h
+++ b/libcrystfel/src/detector.h
@@ -224,9 +224,8 @@ extern void get_pixel_extents(struct detector *det,
                               double *min_x, double *min_y,
                               double *max_x, double *max_y);
 
-extern void fill_in_values(struct detector *det, struct hdfile *f,
-                           struct event* ev);
 extern void fill_in_adu(struct image *image);
+extern void adjust_centering_for_rail(struct panel *p);
 
 extern int panel_is_in_rigid_group(const struct rigid_group *rg,
                                    struct panel *p);
diff --git a/libcrystfel/src/hdf5-file.c b/libcrystfel/src/hdf5-file.c
index c36a58b3..dc86d641 100644
--- a/libcrystfel/src/hdf5-file.c
+++ b/libcrystfel/src/hdf5-file.c
@@ -134,8 +134,7 @@ struct hdfile *hdfile_open(const char *filename)
 }
 
 
-int hdfile_set_image(struct hdfile *f, const char *path,
-                     struct panel *p)
+int hdfile_set_image(struct hdfile *f, const char *path)
 {
 	f->dh = H5Dopen2(f->fh, path, H5P_DEFAULT);
 	if ( f->dh < 0 ) {
@@ -1391,8 +1390,10 @@ int hdfile_get_value(struct hdfile *f, const char *name,
 }
 
 
-void fill_in_beam_parameters(struct beam_params *beam, struct hdfile *f,
-                             struct event *ev, struct image *image)
+static void hdfile_fill_in_beam_parameters(struct beam_params *beam,
+                                           struct hdfile *f,
+                                           struct event *ev,
+                                           struct image *image)
 {
 	double eV;
 
@@ -1405,8 +1406,8 @@ void fill_in_beam_parameters(struct beam_params *beam, struct hdfile *f,
 
 		int r;
 
-		r = hdfile_get_value(f, beam->photon_energy_from, ev, &eV,
-		                     H5T_NATIVE_DOUBLE);
+		r = hdfile_get_value(f, beam->photon_energy_from,
+		                     ev, &eV, H5T_NATIVE_DOUBLE);
 		if ( r ) {
 			ERROR("Failed to read '%s'\n",
 			      beam->photon_energy_from);
@@ -1418,6 +1419,36 @@ void fill_in_beam_parameters(struct beam_params *beam, struct hdfile *f,
 }
 
 
+static void hdfile_fill_in_clen(struct detector *det, struct hdfile *f,
+                                struct event *ev)
+{
+	int i;
+
+	for ( i=0; i<det->n_panels; i++ ) {
+
+		struct panel *p = &det->panels[i];
+
+		if ( p->clen_from != NULL ) {
+
+			double val;
+			int r;
+
+			r = hdfile_get_value(f, p->clen_from, ev, &val,
+			                     H5T_NATIVE_DOUBLE);
+			if ( r ) {
+				ERROR("Failed to read '%s'\n", p->clen_from);
+			} else {
+				p->clen = val * 1.0e-3;
+			}
+
+		}
+
+		adjust_centering_for_rail(p);
+
+	}
+}
+
+
 int hdf5_read(struct hdfile *f, struct image *image, const char *element,
               int satcorr)
 {
@@ -1433,7 +1464,7 @@ int hdf5_read(struct hdfile *f, struct image *image, const char *element,
 	if ( element == NULL ) {
 		fail = hdfile_set_first_image(f, "/");
 	} else {
-		fail = hdfile_set_image(f, element, NULL);
+		fail = hdfile_set_image(f, element);
 	}
 
 	if ( fail ) {
@@ -1482,7 +1513,7 @@ int hdf5_read(struct hdfile *f, struct image *image, const char *element,
 
 	if ( image->beam != NULL ) {
 
-		fill_in_beam_parameters(image->beam, f, NULL, image);
+		hdfile_fill_in_beam_parameters(image->beam, f, NULL, image);
 
 		if ( image->lambda > 1000 ) {
 			/* Error message covers a silly value in the beam file
@@ -1747,7 +1778,7 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
 				return 1;
 			}
 
-			fail = hdfile_set_image(f, panel_full_path, p);
+			fail = hdfile_set_image(f, panel_full_path);
 
 			free(panel_full_path);
 
@@ -1769,7 +1800,7 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
 					free(image->sat);
 					return 1;
 				}
-				fail = hdfile_set_image(f, p->data, p);
+				fail = hdfile_set_image(f, p->data);
 
 			}
 
@@ -1886,13 +1917,13 @@ int hdf5_read2(struct hdfile *f, struct image *image, struct event *ev,
 
 	H5Dclose(f->dh);
 	f->data_open = 0;
-	fill_in_values(image->det, f, ev);
+	hdfile_fill_in_clen(image->det, f, ev);
 
 	if ( satcorr ) debodge_saturation(f, image);
 
 	if ( image->beam != NULL ) {
 
-		fill_in_beam_parameters(image->beam, f, ev, image);
+		hdfile_fill_in_beam_parameters(image->beam, f, ev, image);
 
 		if ( (image->lambda > 1.0) || (image->lambda < 1e-20) ) {
 
@@ -2299,7 +2330,7 @@ int hdfile_set_first_image(struct hdfile *f, const char *group)
 	for ( i=0; i<n; i++ ) {
 
 		if ( is_image[i] ) {
-			hdfile_set_image(f, names[i], NULL);
+			hdfile_set_image(f, names[i]);
 			for ( j=0; j<n; j++ ) free(names[j]);
 			free(is_image);
 			free(is_group);
diff --git a/libcrystfel/src/hdf5-file.h b/libcrystfel/src/hdf5-file.h
index f22d4d54..ab53dd2e 100644
--- a/libcrystfel/src/hdf5-file.h
+++ b/libcrystfel/src/hdf5-file.h
@@ -67,8 +67,8 @@ extern int hdf5_read2(struct hdfile *f, struct image *image,
 extern int check_path_existence(hid_t fh, const char *path);
 
 extern struct hdfile *hdfile_open(const char *filename);
-int hdfile_set_image(struct hdfile *f, const char *path,
-                     struct panel *p);
+int hdfile_set_image(struct hdfile *f, const char *path);
+
 extern int16_t *hdfile_get_image_binned(struct hdfile *hdfile,
                                          int binning, int16_t *maxp);
 extern char **hdfile_read_group(struct hdfile *f, int *n, const char *parent,
diff --git a/libcrystfel/src/image.c b/libcrystfel/src/image.c
index 09f4958d..4e1a8b86 100644
--- a/libcrystfel/src/image.c
+++ b/libcrystfel/src/image.c
@@ -3,12 +3,12 @@
  *
  * Handle images and image features
  *
- * Copyright © 2012-2016 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2017 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
  *   2014      Kenneth Beyerlein <kenneth.beyerlein@desy.de>
- *   2011-2016 Thomas White <taw@physics.org>
+ *   2011-2017 Thomas White <taw@physics.org>
  *
  * This file is part of CrystFEL.
  *
@@ -31,9 +31,14 @@
 #include <assert.h>
 #include <math.h>
 #include <stdio.h>
+#include <hdf5.h>
+#include <cbflib/cbf.h>
 
 #include "image.h"
 #include "utils.h"
+#include "events.h"
+#include "hdf5-file.h"
+#include "detector.h"
 
 /**
  * SECTION:image
@@ -50,6 +55,14 @@
  */
 
 
+struct imagefile
+{
+	enum imagefile_type type;
+	char *filename;
+	struct hdfile *hdfile;
+};
+
+
 struct _imagefeaturelist
 {
 	struct imagefeature	*features;
@@ -299,3 +312,654 @@ void free_all_crystals(struct image *image)
 	free(image->crystals);
 	image->n_crystals = 0;
 }
+
+
+/**************************** Image field lists *******************************/
+
+struct imagefile_field_list
+{
+	char **fields;
+	int n_fields;
+	int max_fields;
+};
+
+
+struct imagefile_field_list *new_imagefile_field_list()
+{
+	struct imagefile_field_list *n;
+
+	n = calloc(1, sizeof(struct imagefile_field_list));
+	if ( n == NULL ) return NULL;
+
+	n->max_fields = 32;
+	n->fields = malloc(n->max_fields*sizeof(char *));
+	if ( n->fields == NULL ) {
+		free(n);
+		return NULL;
+	}
+
+	return n;
+}
+
+
+void free_imagefile_field_list(struct imagefile_field_list *n)
+{
+	int i;
+	for ( i=0; i<n->n_fields; i++ ) {
+		free(n->fields[i]);
+	}
+	free(n->fields);
+	free(n);
+}
+
+
+void add_imagefile_field(struct imagefile_field_list *copyme, const char *name)
+{
+	int i;
+
+	/* Already on the list?   Don't re-add if so. */
+	for ( i=0; i<copyme->n_fields; i++ ) {
+		if ( strcmp(copyme->fields[i], name) == 0 ) return;
+	}
+
+	/* Need more space? */
+	if ( copyme->n_fields == copyme->max_fields ) {
+
+		char **nfields;
+		int nmax = copyme->max_fields + 32;
+
+		nfields = realloc(copyme->fields, nmax*sizeof(char *));
+		if ( nfields == NULL ) {
+			ERROR("Failed to allocate space for new HDF5 field.\n");
+			return;
+		}
+
+		copyme->max_fields = nmax;
+		copyme->fields = nfields;
+
+	}
+
+	copyme->fields[copyme->n_fields] = strdup(name);
+	if ( copyme->fields[copyme->n_fields] == NULL ) {
+		ERROR("Failed to add field for copying '%s'\n", name);
+		return;
+	}
+
+	copyme->n_fields++;
+}
+
+
+/******************************* CBF files ************************************/
+
+
+static char *cbf_strerr(int e)
+{
+	char *err;
+
+	err = malloc(1024);
+	if ( err == NULL ) return NULL;
+
+	err[0] = '\0';
+
+	/* NB Sum of lengths of all strings must be less than 1024 */
+	if ( e & CBF_FORMAT ) strcat(err, "Invalid format");
+	if ( e & CBF_ALLOC ) strcat(err, "Memory allocation failed");
+	if ( e & CBF_ARGUMENT ) strcat(err, "Invalid argument");
+	if ( e & CBF_ASCII ) strcat(err, "Value is ASCII");
+	if ( e & CBF_BINARY ) strcat(err, "Value is binary");
+	if ( e & CBF_BITCOUNT ) strcat(err, "Wrong number of bits");
+	if ( e & CBF_ENDOFDATA ) strcat(err, "End of data");
+	if ( e & CBF_FILECLOSE ) strcat(err, "File close error");
+	if ( e & CBF_FILEOPEN ) strcat(err, "File open error");
+	if ( e & CBF_FILEREAD ) strcat(err, "File read error");
+	if ( e & CBF_FILETELL ) strcat(err, "File tell error");
+	if ( e & CBF_FILEWRITE ) strcat(err, "File write error");
+	if ( e & CBF_IDENTICAL ) strcat(err, "Name already exists");
+	if ( e & CBF_NOTFOUND ) strcat(err, "Not found");
+	if ( e & CBF_OVERFLOW ) strcat(err, "Overflow");
+	if ( e & CBF_UNDEFINED ) strcat(err, "Number undefined");
+	if ( e & CBF_NOTIMPLEMENTED ) strcat(err, "Not implemented");
+
+	return err;
+}
+
+
+static int unpack_panels(struct image *image, signed int *data, int data_width,
+                         int data_height)
+{
+	int pi;
+
+	/* FIXME: Load these masks from an HDF5 file, if filenames are
+	 * given in the geometry file */
+	uint16_t *flags = NULL;
+	float *sat = NULL;
+
+	image->dp = malloc(image->det->n_panels * sizeof(float *));
+	image->bad = malloc(image->det->n_panels * sizeof(int *));
+	image->sat = malloc(image->det->n_panels * sizeof(float *));
+	if ( (image->dp == NULL) || (image->bad == NULL)
+	  || (image->sat == NULL) )
+	{
+		ERROR("Failed to allocate panels.\n");
+		return 1;
+	}
+
+	for ( pi=0; pi<image->det->n_panels; pi++ ) {
+
+		struct panel *p;
+		int fs, ss;
+
+		p = &image->det->panels[pi];
+		image->dp[pi] = malloc(p->w*p->h*sizeof(float));
+		image->bad[pi] = calloc(p->w*p->h, sizeof(int));
+		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;
+		}
+
+		if ( p->mask != NULL ) {
+			ERROR("WARNING: Bad pixel masks do not currently work "
+			      "with CBF files\n");
+			ERROR(" (bad pixel regions specified in the geometry "
+			      "file will be used, however)\n");
+		}
+
+		if ( p->satmap != NULL ) {
+			ERROR("WARNING: Saturation maps do not currently work "
+			      "with CBF files\n");
+		}
+
+		if ( (p->orig_min_fs + p->w > data_width)
+		  || (p->orig_min_ss + p->h > data_height) )
+		{
+			ERROR("Panel %s is outside range of data in CBF file\n",
+			      p->name);
+			return 1;
+		}
+
+		for ( ss=0; ss<p->h; ss++ ) {
+		for ( fs=0; fs<p->w; fs++ ) {
+
+			int idx;
+			int cfs, css;
+			int bad = 0;
+
+			cfs = fs+p->orig_min_fs;
+			css = ss+p->orig_min_ss;
+			idx = cfs + css*data_width;
+
+			image->dp[pi][fs+p->w*ss] = data[idx];
+
+			if ( sat != NULL ) {
+				image->sat[pi][fs+p->w*ss] = sat[idx];
+			} else {
+				image->sat[pi][fs+p->w*ss] = INFINITY;
+			}
+
+			if ( p->no_index ) bad = 1;
+
+			if ( in_bad_region(image->det, p, cfs, css) ) {
+				bad = 1;
+			}
+
+			if ( flags != NULL ) {
+
+				int f;
+
+				f = flags[idx];
+
+				/* Bad if it's missing any of the "good" bits */
+				if ( (f & image->det->mask_good)
+				       != image->det->mask_good ) bad = 1;
+
+				/* Bad if it has any of the "bad" bits. */
+				if ( f & image->det->mask_bad ) bad = 1;
+
+			}
+			image->bad[pi][fs+p->w*ss] = bad;
+
+		}
+		}
+
+	}
+
+	return 0;
+}
+
+
+static void cbf_fill_in_beam_parameters(struct beam_params *beam,
+                                        struct imagefile *f,
+                                        struct image *image)
+{
+	double eV;
+
+	if ( beam->photon_energy_from == NULL ) {
+
+		/* Explicit value given */
+		eV = beam->photon_energy;
+
+	} else {
+
+		ERROR("Can't get photon energy from CBF yet.\n");
+		eV = 0.0;
+
+	}
+
+	image->lambda = ph_en_to_lambda(eV_to_J(eV))*beam->photon_energy_scale;
+}
+
+
+static void cbf_fill_in_clen(struct detector *det, struct imagefile *f)
+{
+	int i;
+
+	for ( i=0; i<det->n_panels; i++ ) {
+
+		struct panel *p = &det->panels[i];
+
+		if ( p->clen_from != NULL ) {
+
+			ERROR("Can't get clen from CBF yet.\n");
+
+		}
+
+		adjust_centering_for_rail(p);
+
+	}
+}
+
+
+static int read_cbf(struct imagefile *f, struct image *image)
+{
+	cbf_handle cbfh;
+	FILE *fh;
+	int r;
+	unsigned int compression;
+	int binary_id, minelement, maxelement, elsigned, elunsigned;
+	size_t elsize, elements, elread, dimfast, dimmid, dimslow, padding;
+	const char *byteorder;
+	signed int *data;
+
+	if ( image->det == NULL ) {
+		ERROR("read_cbf() needs a geometry\n");
+		return 1;
+	}
+
+	if ( cbf_make_handle(&cbfh) ) {
+		ERROR("Failed to allocate CBF handle\n");
+		return 1;
+	}
+
+	fh = fopen(f->filename, "rb");
+	if ( fh == NULL ) {
+		ERROR("Failed to open '%s'\n", f->filename);
+		return 1;
+	}
+	/* CBFlib calls fclose(fh) when it's ready */
+
+	if ( cbf_read_widefile(cbfh, fh, 0) ) {
+		ERROR("Failed to read CBF file '%s'\n", f->filename);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	/* Select row 0 in data column inside array_data */
+	cbf_find_category(cbfh, "array_data");
+	cbf_find_column(cbfh, "data");
+	cbf_select_row(cbfh, 0);
+
+	/* Get parameters for array read */
+	r = cbf_get_integerarrayparameters_wdims(cbfh, &compression, &binary_id,
+	                                         &elsize, &elsigned, &elunsigned,
+	                                         &elements,
+	                                         &minelement, &maxelement,
+	                                         &byteorder,
+	                                         &dimfast, &dimmid, &dimslow,
+	                                         &padding);
+	if ( r ) {
+		char *err = cbf_strerr(r);
+		ERROR("Failed to read CBF array parameters: %s\n", err);
+		free(err);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( dimslow != 0 ) {
+		ERROR("CBF data array is 3D - don't know what to do with it\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( dimfast*dimmid*elsize > 10e9 ) {
+		ERROR("CBF data is far too big (%i x %i x %i bytes).\n",
+		      (int)dimfast, (int)dimmid, (int)elsize);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( elsize != 4 ) {
+		STATUS("Don't know what to do with element size %i\n",
+		       (int)elsize);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( strcmp(byteorder, "little_endian") != 0 ) {
+		STATUS("Don't know what to do with non-little-endian datan\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	data = malloc(elsize*dimfast*dimmid);
+	if ( data == NULL ) {
+		ERROR("Failed to allocate memory for CBF data\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	r = cbf_get_integerarray(cbfh, &binary_id, data, elsize, 1,
+	                         elsize*dimfast*dimmid, &elread);
+	if ( r ) {
+		char *err = cbf_strerr(r);
+		ERROR("Failed to read CBF array: %s\n", err);
+		free(err);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	unpack_panels(image, data, dimfast, dimmid);
+	free(data);
+
+	if ( image->beam != NULL ) {
+		cbf_fill_in_beam_parameters(image->beam, f, image);
+		if ( image->lambda > 1000 ) {
+			ERROR("WARNING: Missing or nonsensical wavelength "
+			      "(%e m) for %s.\n",
+			      image->lambda, image->filename);
+		}
+	}
+	cbf_fill_in_clen(image->det, f);
+	fill_in_adu(image);
+
+	cbf_free_handle(cbfh);
+	return 0;
+}
+
+
+static float *convert_float(signed int *data, int w, int h)
+{
+	float *df;
+	long int i;
+
+	df = malloc(sizeof(float)*w*h);
+	if ( df == NULL ) return NULL;
+
+	for ( i=0; i<w*h; i++ ) {
+		df[i] = data[i];
+	}
+
+	return df;
+}
+
+
+static int read_cbf_simple(struct imagefile *f, struct image *image)
+{
+	cbf_handle cbfh;
+	FILE *fh;
+	int r;
+	unsigned int compression;
+	int binary_id, minelement, maxelement, elsigned, elunsigned;
+	size_t elsize, elements, elread, dimfast, dimmid, dimslow, padding;
+	const char *byteorder;
+	signed int *data;
+
+	if ( cbf_make_handle(&cbfh) ) {
+		ERROR("Failed to allocate CBF handle\n");
+		return 1;
+	}
+
+	fh = fopen(f->filename, "rb");
+	if ( fh == NULL ) {
+		ERROR("Failed to open '%s'\n", f->filename);
+		return 1;
+	}
+	/* CBFlib calls fclose(fh) when it's ready */
+
+	if ( cbf_read_widefile(cbfh, fh, 0) ) {
+		ERROR("Failed to read CBF file '%s'\n", f->filename);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	/* Select row 0 in data column inside array_data */
+	cbf_find_category(cbfh, "array_data");
+	cbf_find_column(cbfh, "data");
+	cbf_select_row(cbfh, 0);
+
+	/* Get parameters for array read */
+	r = cbf_get_integerarrayparameters_wdims(cbfh, &compression, &binary_id,
+	                                         &elsize, &elsigned, &elunsigned,
+	                                         &elements,
+	                                         &minelement, &maxelement,
+	                                         &byteorder,
+	                                         &dimfast, &dimmid, &dimslow,
+	                                         &padding);
+	if ( r ) {
+		char *err = cbf_strerr(r);
+		ERROR("Failed to read CBF array parameters: %s\n", err);
+		free(err);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( dimslow != 0 ) {
+		ERROR("CBF data array is 3D - don't know what to do with it\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( dimfast*dimmid*elsize > 10e9 ) {
+		ERROR("CBF data is far too big (%i x %i x %i bytes).\n",
+		      (int)dimfast, (int)dimmid, (int)elsize);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( elsize != 4 ) {
+		STATUS("Don't know what to do with element size %i\n",
+		       (int)elsize);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	if ( strcmp(byteorder, "little_endian") != 0 ) {
+		STATUS("Don't know what to do with non-little-endian datan\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	data = malloc(elsize*dimfast*dimmid);
+	if ( data == NULL ) {
+		ERROR("Failed to allocate memory for CBF data\n");
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	r = cbf_get_integerarray(cbfh, &binary_id, data, elsize, 1,
+	                         elsize*dimfast*dimmid, &elread);
+	if ( r ) {
+		char *err = cbf_strerr(r);
+		ERROR("Failed to read CBF array: %s\n", err);
+		free(err);
+		cbf_free_handle(cbfh);
+		return 1;
+	}
+
+	image->det = simple_geometry(image, dimfast, dimmid);
+	image->dp = malloc(sizeof(float *));
+	if ( image->dp == NULL ) {
+		ERROR("Failed to allocate dp array\n");
+		return 1;
+	}
+	image->dp[0] = convert_float(data, dimfast, dimmid);
+	if ( image->dp[0] == NULL ) {
+		ERROR("Failed to allocate dp array\n");
+		return 1;
+	}
+
+	if ( image->beam != NULL ) {
+		cbf_fill_in_beam_parameters(image->beam, f, image);
+		if ( image->lambda > 1000 ) {
+			ERROR("WARNING: Missing or nonsensical wavelength "
+			      "(%e m) for %s.\n",
+			      image->lambda, image->filename);
+		}
+	}
+	cbf_fill_in_clen(image->det, f);
+	fill_in_adu(image);
+
+	cbf_free_handle(cbfh);
+	return 0;
+}
+
+
+/****************************** Image files ***********************************/
+
+
+static signed int is_cbf_file(const char *filename)
+{
+	FILE *fh;
+	char line[1024];
+
+	fh = fopen(filename, "r");
+	if ( fh == NULL ) return -1;
+
+	if ( fgets(line, 1024, fh) == NULL ) return -1;
+	fclose(fh);
+
+	if ( strstr(line, "CBF") == NULL ) {
+		return 0;
+	}
+
+	return 1;
+}
+
+
+struct imagefile *imagefile_open(const char *filename)
+{
+	struct imagefile *f;
+
+	f = malloc(sizeof(struct imagefile));
+	if ( f == NULL ) return NULL;
+
+	if ( H5Fis_hdf5(filename) > 0 ) {
+
+		/* This is an HDF5, pass through to HDF5 layer */
+		f->type = IMAGEFILE_HDF5;
+		f->hdfile = hdfile_open(filename);
+
+		if ( f->hdfile == NULL ) {
+			free(f);
+			return NULL;
+		}
+
+	} else if ( is_cbf_file(filename) > 0 ) {
+
+		f->type = IMAGEFILE_CBF;
+
+	}
+
+	f->filename = strdup(filename);
+	return f;
+}
+
+
+int imagefile_read(struct imagefile *f, struct image *image,
+                   struct event *event)
+{
+	if ( f->type == IMAGEFILE_HDF5 ) {
+		return hdf5_read2(f->hdfile, image, event, 0);
+	} else if ( f->type == IMAGEFILE_CBF ) {
+		return read_cbf(f, image);
+	} else {
+		ERROR("Unknown file type %i\n", f->type);
+		return 1;
+	}
+}
+
+
+/* Read a simple file, no multi-event, no prior geometry etc, and
+ * generate a geometry for it */
+int imagefile_read_simple(struct imagefile *f, struct image *image)
+{
+	if ( f->type == IMAGEFILE_HDF5 ) {
+		return hdf5_read(f->hdfile, image, NULL, 0);
+	} else {
+		return read_cbf_simple(f, image);
+	}
+}
+
+
+enum imagefile_type imagefile_get_type(struct imagefile *f)
+{
+	assert(f != NULL);
+	return f->type;
+}
+
+
+struct hdfile *imagefile_get_hdfile(struct imagefile *f)
+{
+	if ( f->type != IMAGEFILE_HDF5 ) {
+		ERROR("Not an HDF5 file!\n");
+		return NULL;
+	}
+
+	return f->hdfile;
+}
+
+
+void imagefile_copy_fields(struct imagefile *f,
+                           const struct imagefile_field_list *copyme,
+                           FILE *fh, struct event *ev)
+{
+	int i;
+
+	if ( copyme == NULL ) return;
+
+	for ( i=0; i<copyme->n_fields; i++ ) {
+
+		char *val;
+		char *field;
+
+		field = copyme->fields[i];
+
+		if ( f->type == IMAGEFILE_HDF5 ) {
+			val = hdfile_get_string_value(f->hdfile, field, ev);
+			if ( field[0] == '/' ) {
+				fprintf(fh, "hdf5%s = %s\n", field, val);
+			} else {
+				fprintf(fh, "hdf5/%s = %s\n", field, val);
+			}
+			free(val);
+
+		} else {
+			STATUS("Mock CBF variable\n");
+			fprintf(fh, "cbf/%s = %s\n", field, "(FIXME)");
+		}
+
+
+	}
+}
+
+
+void imagefile_close(struct imagefile *f)
+{
+	if ( f->type == IMAGEFILE_HDF5 ) {
+		hdfile_close(f->hdfile);
+	}
+	free(f->filename);
+	free(f);
+}
diff --git a/libcrystfel/src/image.h b/libcrystfel/src/image.h
index 9fd9b495..ead5cd4d 100644
--- a/libcrystfel/src/image.h
+++ b/libcrystfel/src/image.h
@@ -3,11 +3,11 @@
  *
  * Handle images and image features
  *
- * Copyright © 2012-2016 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2017 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
- *   2009-2016 Thomas White <taw@physics.org>
+ *   2009-2017 Thomas White <taw@physics.org>
  *   2014      Valerio Mariani
  *
  *
@@ -44,6 +44,8 @@ struct detector;
 struct imagefeature;
 struct sample;
 struct image;
+struct imagefile;
+struct imagefile_field_list;
 
 #include "utils.h"
 #include "cell.h"
@@ -51,6 +53,7 @@ struct image;
 #include "reflist.h"
 #include "crystal.h"
 #include "index.h"
+#include "events.h"
 
 /**
  * SpectrumType:
@@ -87,6 +90,15 @@ struct imagefeature {
 	const char                      *name;
 };
 
+
+/* An enum representing the image file formats we can handle */
+enum imagefile_type
+{
+	IMAGEFILE_HDF5,
+	IMAGEFILE_CBF
+};
+
+
 /* An opaque type representing a list of image features */
 typedef struct _imagefeaturelist ImageFeatureList;
 
@@ -120,7 +132,7 @@ struct beam_params
  *    struct detector         *det;
  *    struct beam_params      *beam;
  *    char                    *filename;
- *    const struct copy_hdf5_field *copyme;
+ *    const struct imagefile_field_list *copyme;
  *
  *    int                     id;
  *
@@ -152,8 +164,8 @@ struct beam_params
  * returned by the low-level indexing system. <structfield>n_crystals</structfield>
  * is the number of crystals which were found in the image.
  *
- * <structfield>copyme</structfield> represents a list of HDF5 fields to copy
- * to the output stream.
+ * <structfield>copyme</structfield> represents a list of fields in the image
+ * file (e.g. HDF5 fields or CBF headers) to copy to the output stream.
  **/
 struct image;
 
@@ -171,7 +183,7 @@ struct image {
 	struct beam_params      *beam;  /* The nominal beam parameters */
 	char                    *filename;
 	struct event            *event;
-	const struct copy_hdf5_field *copyme;
+	const struct imagefile_field_list *copyme;
 	struct stuff_from_stream *stuff_from_stream;
 
 	double                  avg_clen;  /* Average camera length extracted
@@ -233,6 +245,25 @@ extern void image_add_crystal(struct image *image, Crystal *cryst);
 extern void remove_flagged_crystals(struct image *image);
 extern void free_all_crystals(struct image *image);
 
+/* Image files */
+extern struct imagefile *imagefile_open(const char *filename);
+extern int imagefile_read(struct imagefile *f, struct image *image,
+                          struct event *event);
+extern int imagefile_read_simple(struct imagefile *f, struct image *image);
+extern struct hdfile *imagefile_get_hdfile(struct imagefile *f);
+extern enum imagefile_type imagefile_get_type(struct imagefile *f);
+extern void imagefile_copy_fields(struct imagefile *f,
+                                  const struct imagefile_field_list *copyme,
+                                  FILE *fh, struct event *ev);
+extern void imagefile_close(struct imagefile *f);
+
+/* Field lists */
+extern struct imagefile_field_list *new_imagefile_field_list(void);
+extern void free_imagefile_field_list(struct imagefile_field_list *f);
+
+extern void add_imagefile_field(struct imagefile_field_list *copyme,
+                                const char *name);
+
 #ifdef __cplusplus
 }
 #endif
diff --git a/libcrystfel/src/stream.c b/libcrystfel/src/stream.c
index 17da74b2..fb4b0c70 100644
--- a/libcrystfel/src/stream.c
+++ b/libcrystfel/src/stream.c
@@ -806,8 +806,8 @@ static int write_crystal(Stream *st, Crystal *cr, int include_reflections)
 }
 
 
-int write_chunk(Stream *st, struct image *i, struct hdfile *hdfile,
-                int include_peaks, int include_reflections, struct event* ev)
+int write_chunk(Stream *st, struct image *i, struct imagefile *imfile,
+                int include_peaks, int include_reflections, struct event *ev)
 {
 	int j;
 	char *indexer;
@@ -832,7 +832,7 @@ int write_chunk(Stream *st, struct image *i, struct hdfile *hdfile,
 	fprintf(st->fh, "beam_divergence = %.2e rad\n", i->div);
 	fprintf(st->fh, "beam_bandwidth = %.2e (fraction)\n", i->bw);
 
-	copy_hdf5_fields(hdfile, i->copyme, st->fh, ev);
+	imagefile_copy_fields(imfile, i->copyme, st->fh, ev);
 
 	if ( i->det != NULL ) {
 
diff --git a/libcrystfel/src/stream.h b/libcrystfel/src/stream.h
index ff8628c0..764e3e36 100644
--- a/libcrystfel/src/stream.h
+++ b/libcrystfel/src/stream.h
@@ -3,11 +3,11 @@
  *
  * Stream tools
  *
- * Copyright © 2013-2014 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2013-2017 Deutsches Elektronen-Synchrotron DESY,
  *                       a research centre of the Helmholtz Association.
  *
  * Authors:
- *   2010-2014 Thomas White <taw@physics.org>
+ *   2010-2017 Thomas White <taw@physics.org>
  *   2014      Valerio Mariani
  *   2011      Andrew Aquila
  *
@@ -39,6 +39,7 @@
 struct image;
 struct hdfile;
 struct event;
+struct imagefile;
 #include "cell.h"
 
 #define GEOM_START_MARKER "----- Begin geometry file -----"
@@ -106,10 +107,15 @@ extern int read_chunk(Stream *st, struct image *image);
 extern int read_chunk_2(Stream *st, struct image *image,
                            StreamReadFlags srf);
 
-extern int write_chunk(Stream *st, struct image *image, struct hdfile *hdfile,
+extern int write_chunk(Stream *st, struct image *image, struct imagefile *imfile,
                        int include_peaks, int include_reflections,
                        struct event *ev);
 
+extern int write_chunk_2(Stream *st, struct image *image,
+                         struct imagefile *imfile,
+                         int include_peaks, int include_reflections,
+                         struct event *ev);
+
 extern void write_command(Stream *st, int argc, char *argv[]);
 extern void write_geometry_file(Stream *st, const char *geom_filename);
 extern int rewind_stream(Stream *st);