diff options
author | Thomas White <taw@physics.org> | 2020-06-02 17:56:00 +0200 |
---|---|---|
committer | Thomas White <taw@physics.org> | 2020-07-29 18:53:33 +0200 |
commit | 6b9dc7ada924a69b24e0525dabf96e32fbc45ab9 (patch) | |
tree | 3cf05671a8d508065e3dd50891dc1216fa5fffc0 /libcrystfel | |
parent | df2d6ef008187db7a8df6e65e507de7a98d9b071 (diff) |
Move HDF5 and CBF stuff into separate files
Diffstat (limited to 'libcrystfel')
-rw-r--r-- | libcrystfel/CMakeLists.txt | 2 | ||||
-rw-r--r-- | libcrystfel/src/image-cbf.c | 636 | ||||
-rw-r--r-- | libcrystfel/src/image-cbf.h | 58 | ||||
-rw-r--r-- | libcrystfel/src/image-hdf5.c | 1223 | ||||
-rw-r--r-- | libcrystfel/src/image-hdf5.h | 67 | ||||
-rw-r--r-- | libcrystfel/src/image.c | 1888 |
6 files changed, 2043 insertions, 1831 deletions
diff --git a/libcrystfel/CMakeLists.txt b/libcrystfel/CMakeLists.txt index 24278009..32bf6fe3 100644 --- a/libcrystfel/CMakeLists.txt +++ b/libcrystfel/CMakeLists.txt @@ -65,6 +65,8 @@ set(LIBCRYSTFEL_SOURCES src/datatemplate.c src/colscale.c src/detgeom.c + src/image-cbf.c + src/image-hdf5.c ${BISON_symopp_OUTPUTS} ${FLEX_symopl_OUTPUTS} ) diff --git a/libcrystfel/src/image-cbf.c b/libcrystfel/src/image-cbf.c new file mode 100644 index 00000000..19c4d5c3 --- /dev/null +++ b/libcrystfel/src/image-cbf.c @@ -0,0 +1,636 @@ +/* + * image-cbf.c + * + * Image loading, CBF parts + * + * Copyright © 2012-2020 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2020 Thomas White <taw@physics.org> + * + * This file is part of CrystFEL. + * + * CrystFEL is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * CrystFEL is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with CrystFEL. If not, see <http://www.gnu.org/licenses/>. + * + */ + +#include <config.h> + +#include <stdlib.h> +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <hdf5.h> +#include <zlib.h> + +#include "image.h" +#include "utils.h" +#include "events.h" +#include "detgeom.h" + +#include "datatemplate.h" +#include "datatemplate_priv.h" + +static void add_out(float val, float *data_out, int nmemb_out, + int *outpos, int *nrej) +{ + if ( *outpos < nmemb_out ) { + data_out[(*outpos)++] = val; + } else { + (*nrej)++; + } +} + + +/* Reverses byte offset compression and converts to single precision float. + * Note that this compression scheme specifies the data format of the input + * data, therefore the X-Binary-Element-Type is completely ignored. */ +static void decode_cbf_byte_offset(float *data_out, int nmemb_out, + const int8_t *data_in, const size_t n) +{ + int inpos = 0; + int outpos = 0; + int nrej = 0; + float val = 0.0; + + while ( inpos < n ) { + + int64_t delta = data_in[inpos++]; + + if ( (delta >= -127) && (delta <= 127) ) { + val += delta; + add_out(val, data_out, nmemb_out, &outpos, &nrej); + continue; + } + + delta = *(int16_t *)(data_in+inpos); + inpos += 2; + + if ( (delta >= -32767) && (delta <= 32767) ) { + val += delta; + add_out(val, data_out, nmemb_out, &outpos, &nrej); + continue; + } + + delta = *(int32_t *)(data_in+inpos); + inpos += 4; + + if ( (delta >= -2147483647) && (delta <= 2147483647) ) { + val += delta; + add_out(val, data_out, nmemb_out, &outpos, &nrej); + continue; + } + + delta = *(int64_t *)(data_in+inpos); + inpos += 8; + val += delta; + add_out(val, data_out, nmemb_out, &outpos, &nrej); + + } + + if ( nrej > 0 ) { + STATUS("%i elements rejected\n", nrej); + } +} + + +static int binary_start(char *data) +{ + char *datac = data; + if ( (datac[0] == (char)0x0c) && (datac[1] == (char)0x1a) + && (datac[2] == (char)0x04) && (datac[3] == (char)0xd5) ) return 1; + return 0; +} + + +enum cbf_data_conversion +{ + CBF_NO_CONVERSION, + CBF_BYTE_OFFSET, + CBF_PACKED, + CBF_CANONICAL +}; + +enum cbf_data_type +{ + CBF_NO_TYPE, + CBF_ELEMENT_U8, + CBF_ELEMENT_S8, + CBF_ELEMENT_U16, + CBF_ELEMENT_S16, + CBF_ELEMENT_U32, + CBF_ELEMENT_S32, + CBF_ELEMENT_F32, + CBF_ELEMENT_F64, +}; + + +static enum cbf_data_type parse_element_type(const char *t) +{ + if ( strstr(t, "signed 8-bit integer") != NULL ) + { + return CBF_ELEMENT_S8; + } + + if ( strstr(t, "unsigned 8-bit integer") != NULL ) + { + return CBF_ELEMENT_U8; + } + + if ( strstr(t, "signed 16-bit integer") != NULL ) + { + return CBF_ELEMENT_S16; + } + + if ( strstr(t, "unsigned 16-bit integer") != NULL ) + { + return CBF_ELEMENT_U16; + } + + if ( strstr(t, "signed 32-bit integer") != NULL ) + { + return CBF_ELEMENT_S32; + } + + if ( strstr(t, "unsigned 32-bit integer") != NULL ) + { + return CBF_ELEMENT_U32; + } + + if ( strstr(t, "signed 32-bit real IEEE") != NULL ) + { + return CBF_ELEMENT_F32; + } + + if ( strstr(t, "signed 64-bit real IEEE") != NULL ) + { + return CBF_ELEMENT_F64; + } + + /* complex type is unsupported */ + + return CBF_NO_TYPE; +} + + +static size_t element_size(enum cbf_data_type t) +{ + switch ( t ) { + case CBF_ELEMENT_S8 : return 1; + case CBF_ELEMENT_U8 : return 1; + case CBF_ELEMENT_S16 : return 2; + case CBF_ELEMENT_U16 : return 2; + case CBF_ELEMENT_S32 : return 4; + case CBF_ELEMENT_U32 : return 4; + case CBF_ELEMENT_F32 : return 4; + case CBF_ELEMENT_F64 : return 8; + default : return 0; + } +} + + + +static int convert_type(float *data_out, long nmemb_exp, + enum cbf_data_type eltype, + void *data_in, size_t data_in_len) +{ + long int i; + long int o = 0; + size_t elsize = element_size(eltype); + + if ( elsize == 0 ) return 1; + + if ( nmemb_exp * elsize > data_in_len ) { + ERROR("Not enough CBF data for image size/type!\n"); + return 1; + } + + for ( i=0; i<nmemb_exp; i++ ) { + switch ( eltype ) { + + case CBF_ELEMENT_S8: + data_out[o++] = ((int8_t *)data_in)[i]; + break; + + case CBF_ELEMENT_U8: + data_out[o++] = ((uint8_t *)data_in)[i]; + break; + + case CBF_ELEMENT_S16: + data_out[o++] = ((int16_t *)data_in)[i]; + break; + + case CBF_ELEMENT_U16: + data_out[o++] = ((uint16_t *)data_in)[i]; + break; + + case CBF_ELEMENT_S32: + data_out[o++] = ((int32_t *)data_in)[i]; + break; + + case CBF_ELEMENT_U32: + data_out[o++] = ((uint32_t *)data_in)[i]; + break; + + case CBF_ELEMENT_F32: + data_out[o++] = ((float *)data_in)[i]; + break; + + case CBF_ELEMENT_F64: + data_out[o++] = ((double *)data_in)[i]; + break; + + case CBF_NO_TYPE: + break; + } + } + + return 0; +} + + +static float *read_cbf_data(const char *filename, int gz, int *w, int *h) +{ + FILE *fh; + void *buf = NULL; + char *rval; + size_t data_compressed_len = 0; + float *data_out = NULL; + enum cbf_data_conversion data_conversion = CBF_NO_CONVERSION; + enum cbf_data_type data_type = CBF_ELEMENT_U32; /* ITG (2006) 2.3.3.3 */ + int in_binary_section = 0; + + *w = 0; + *h = 0; + + if ( !gz ) { + + fh = fopen(filename, "rb"); + if ( fh == NULL ) { + ERROR("Failed to open '%s'\n", filename); + return NULL; + } + + } else { + + gzFile gzfh; + size_t len, len_read; + const size_t bufinc = 8*1024*1024; /* Allocate buffer in 8Mb chunks */ + size_t bufsz = bufinc; + + gzfh = gzopen(filename, "rb"); + if ( gzfh == NULL ) return NULL; + + /* Set larger buffer size for hopefully faster uncompression */ + gzbuffer(gzfh, 128*1024); + + buf = malloc(bufsz); + if ( buf == NULL ) return NULL; + + len = 0; + do { + + len_read = gzread(gzfh, buf+len, bufinc); + if ( len_read == -1 ) return NULL; + len += len_read; + + if ( len_read == bufinc ) { + bufsz += bufinc; + buf = realloc(buf, bufsz); + if ( buf == NULL ) return NULL; + } + + } while ( len_read == bufinc ); + + fh = fmemopen(buf, len, "rb"); + if ( fh == NULL ) return NULL; + + gzclose(gzfh); + + } + + /* This is really horrible, but there are at least three different types + * of header mingled together (CIF, MIME, DECTRIS), so a real parser + * would be very complicated and much more likely to have weird bugs. */ + do { + + char line[1024]; + long line_start; + + line_start = ftell(fh); + rval = fgets(line, 1023, fh); + if ( rval == NULL ) break; + chomp(line); + + if ( strcmp(line, "--CIF-BINARY-FORMAT-SECTION--") == 0 ) { + in_binary_section = 1; + } + + if ( strcmp(line, "--CIF-BINARY-FORMAT-SECTION----") == 0 ) { + in_binary_section = 0; + } + + if ( in_binary_section ) { + + if ( strncmp(line, "X-Binary-Size: ", 15) == 0 ) { + data_compressed_len = atoi(line+15); + } + + if ( strncmp(line, "X-Binary-Element-Byte-Order: ", 29) == 0 ) { + const char *elbo = line+29; + if ( strcmp(elbo, "LITTLE_ENDIAN") != 0 ) { + ERROR("Unsupported endianness: %s\n", elbo); + free(buf); + fclose(fh); + return NULL; + } + } + + /* Try to spot compression algorithm */ + if ( strstr(line, "conversions=\"x-CBF_BYTE_OFFSET\"") != NULL ) { + data_conversion = CBF_BYTE_OFFSET; + } else if ( strstr(line, "conversions=\"x-CBF_CANONICAL\"") != NULL ) { + data_conversion = CBF_CANONICAL; + } else if ( strstr(line, "conversions=\"x-CBF_PACKED\"") != NULL ) { + data_conversion = CBF_PACKED; + } else if ( strstr(line, "conversions=") != NULL ) { + ERROR("Unrecognised CBF content conversion: %s\n", line); + free(buf); + fclose(fh); + return NULL; + } + + /* Likewise, element type */ + if ( strncmp(line, "X-Binary-Element-Type: ", 23) == 0 ) + { + const char *eltype = (line+23); + data_type = parse_element_type(eltype); + if ( data_type == CBF_NO_TYPE ) { + ERROR("Unrecognised element type: %s\n", + eltype); + free(buf); + fclose(fh); + return NULL; + } + } + + if ( strncmp(line, "X-Binary-Size-Fastest-Dimension: ", 33) == 0 ) { + *w = atoi(line+33); + } + + if ( strncmp(line, "X-Binary-Size-Second-Dimension: ", 32) == 0 ) { + *h = atoi(line+32); + } + + } + + if ( in_binary_section && binary_start(line) ) { + + size_t len_read; + int nmemb_exp; + void *data_compressed; + int r = 0; + + if ( data_compressed_len == 0 ) { + ERROR("Found CBF data before X-Binary-Size!\n"); + free(buf); + fclose(fh); + return NULL; + } + + if ( (*w == 0) || (*h == 0) ) { + ERROR("Found CBF data before dimensions!\n"); + free(buf); + fclose(fh); + return NULL; + } + + if ( data_compressed_len > 100*1024*1024 ) { + ERROR("Stated CBF data size too big\n"); + free(buf); + fclose(fh); + return NULL; + } + + data_compressed = malloc(data_compressed_len); + if ( data_compressed == NULL ) { + ERROR("Failed to allocate memory for CBF data\n"); + free(buf); + fclose(fh); + return NULL; + } + + fseek(fh, line_start+4, SEEK_SET); + len_read = fread(data_compressed, 1, data_compressed_len, fh); + if ( len_read < data_compressed_len ) { + ERROR("Couldn't read entire CBF data\n"); + free(buf); + free(data_compressed); + fclose(fh); + return NULL; + } + + nmemb_exp = (*w) * (*h); + data_out = malloc(nmemb_exp*sizeof(float)); + if ( data_out == NULL ) { + ERROR("Failed to allocate memory for CBF data\n"); + free(buf); + free(data_compressed); + fclose(fh); + return NULL; + } + + switch ( data_conversion ) { + + case CBF_NO_CONVERSION: + r = convert_type(data_out, nmemb_exp, data_type, + data_compressed, + data_compressed_len); + break; + + case CBF_BYTE_OFFSET: + decode_cbf_byte_offset(data_out, nmemb_exp, + data_compressed, + data_compressed_len); + break; + + case CBF_PACKED: + case CBF_CANONICAL: + ERROR("Don't yet know how to decompress " + "CBF_PACKED or CBF_CANONICAL\n"); + free(buf); + free(data_compressed); + fclose(fh); + return NULL; + + } + + free(data_compressed); + + if ( r ) { + free(buf); + free(data_out); + fclose(fh); + return NULL; + } + + free(buf); + fclose(fh); + return data_out; + + } + + } while ( rval != NULL ); + + ERROR("Reached end of CBF file before finding data.\n"); + free(buf); /* might be NULL */ + return NULL; +} + + +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; +} + + +signed int is_cbfgz_file(const char *filename) +{ + gzFile gzfh; + char line[1024]; + + gzfh = gzopen(filename, "rb"); + if ( gzfh == NULL ) return -1; + if ( gzgets(gzfh, line, 1024) == NULL ) return -1; + gzclose(gzfh); + + if ( strstr(line, "CBF") == NULL ) { + return 0; + } + + return 1; +} + + +int image_cbf_read_mask(struct panel_template *p, + const char *filename, const char *event, + int gz, int *bad, int mask_good, int mask_bad) +{ + ERROR("Mask loading from CBF not yet supported\n"); + return 1; +} + + +static int unpack_panels(struct image *image, DataTemplate *dtempl, + float *data, int data_width, int data_height) +{ + int pi; + + image->dp = malloc(dtempl->n_panels * sizeof(float *)); + if ( image->dp == NULL ) { + ERROR("Failed to allocate panels.\n"); + return 1; + } + + for ( pi=0; pi<dtempl->n_panels; pi++ ) { + + struct panel_template *p; + int fs, ss; + int p_w, p_h; + + p = &dtempl->panels[pi]; + p_w = p->orig_max_fs - p->orig_min_fs + 1; + p_h = p->orig_max_ss - p->orig_min_ss + 1; + image->dp[pi] = malloc(p_w*p_h*sizeof(float)); + if ( image->dp[pi] == NULL ) { + ERROR("Failed to allocate panel\n"); + return 1; + } + + 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; + + 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]; + + } + } + + } + + return 0; +} + + +struct image *image_cbf_read(DataTemplate *dtempl, const char *filename, + const char *event, int gz) +{ + struct image *image; + float *data; + int w, h; + + if ( access(filename, R_OK) == -1 ) { + ERROR("File does not exist or cannot be read: %s\n", filename); + return NULL; + } + + image = image_new(); + if ( image == NULL ) { + ERROR("Couldn't allocate image structure.\n"); + return NULL; + } + + data = read_cbf_data(filename, gz, &w, &h); + if ( data == NULL ) { + ERROR("Failed to read CBF data\n"); + return NULL; + } + + unpack_panels(image, dtempl, data, w, h); + free(data); + + //cbf_fill_in_beam_parameters(image->beam, f, image); + //cbf_fill_in_clen(image->det, f); + //fill_in_adu(image); + + return image; +} diff --git a/libcrystfel/src/image-cbf.h b/libcrystfel/src/image-cbf.h new file mode 100644 index 00000000..400d5132 --- /dev/null +++ b/libcrystfel/src/image-cbf.h @@ -0,0 +1,58 @@ +/* + * image-cbf.h + * + * Image loading, CBF parts + * + * Copyright © 2012-2020 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2020 Thomas White <taw@physics.org> + * + * This file is part of CrystFEL. + * + * CrystFEL is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * CrystFEL is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with CrystFEL. If not, see <http://www.gnu.org/licenses/>. + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +/* NB This file is NOT part of the public API, and should NOT + * be installed, but rather stays in the libcrystfel source folder. */ + +#ifndef IMAGE_CBF_H +#define IMAGE_CBF_H + +#include "datatemplate_priv.h" + +extern signed int is_cbf_file(const char *filename); + +extern signed int is_cbfgz_file(const char *filename); + +extern int load_mask_cbf(struct panel_template *p, + const char *filename, const char *event, + int gz, int *bad, int mask_good, int mask_bad); + +extern struct image *image_cbf_read(DataTemplate *dtempl, + const char *filename, + const char *event, int gz); + +extern int image_cbf_read_mask(struct panel_template *p, + const char *filename, const char *event, + int gz, int *bad, + int mask_good, int mask_bad); + +#endif /* IMAGE_CBF_H */ diff --git a/libcrystfel/src/image-hdf5.c b/libcrystfel/src/image-hdf5.c new file mode 100644 index 00000000..6b6b70c5 --- /dev/null +++ b/libcrystfel/src/image-hdf5.c @@ -0,0 +1,1223 @@ +/* + * image-hdf5.c + * + * Image loading, HDF5 parts + * + * Copyright © 2012-2020 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2020 Thomas White <taw@physics.org> + * + * This file is part of CrystFEL. + * + * CrystFEL is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * CrystFEL is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with CrystFEL. If not, see <http://www.gnu.org/licenses/>. + * + */ + +#include <config.h> + +#include <stdlib.h> +#include <assert.h> +#include <math.h> +#include <stdio.h> +#include <hdf5.h> +#include <unistd.h> + +#include "image.h" +#include "utils.h" +#include "events.h" +#include "detgeom.h" + +#include "datatemplate.h" +#include "datatemplate_priv.h" + + +static int load_hdf5_hyperslab(struct panel_template *p, + const char *filename, + const char *event, + void **pdata, + hid_t el_type, size_t el_size, + int skip_placeholders_ok, + const char *path_spec) +{ + struct event *ev; + hid_t fh; + herr_t r; + hsize_t *f_offset, *f_count; + hid_t dh; + int hsi; + herr_t check; + hid_t dataspace, memspace; + hsize_t dims[2]; + char *panel_full_path; + void *data; + int ndims, dpos; + int skip_placeholders = 0; + + if ( access(filename, R_OK) == -1 ) { + ERROR("File does not exist or cannot be read: %s\n", + filename); + return 1; + } + + fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); + if ( fh < 0 ) { + ERROR("Couldn't open file: %s\n", filename); + return 1; + } + + ev = get_event_from_event_string(event); + if ( (ev == NULL) && (event != NULL) ) { + ERROR("Invalid event identifier '%s'\n", event); + H5Fclose(fh); + return 1; + } + + panel_full_path = retrieve_full_path(ev, path_spec); + + if ( !check_path_existence(fh, panel_full_path) ) { + ERROR("Cannot find data for panel %s (%s)\n", + p->name, panel_full_path); + free_event(ev); + H5Fclose(fh); + return 1; + } + + dh = H5Dopen2(fh, panel_full_path, H5P_DEFAULT); + if ( dh < 0 ) { + ERROR("Cannot open data for panel %s (%s)\n", + p->name, panel_full_path); + free(panel_full_path); + free_event(ev); + H5Fclose(fh); + return 1; + } + + free(panel_full_path); + + /* Set up dataspace for file + * (determine where to read the data from) */ + dataspace = H5Dget_space(dh); + ndims = H5Sget_simple_extent_ndims(dataspace); + if ( ndims < 0 ) { + ERROR("Failed to get number of dimensions for panel %s\n", + p->name); + free_event(ev); + H5Fclose(fh); + return 1; + } + + if ( ndims != p->dim_structure->num_dims ) { + /* Dimensionality doesn't match */ + int n_nonplaceholder = 0; + for ( hsi=0; hsi<p->dim_structure->num_dims; hsi++ ) { + if ( p->dim_structure->dims[hsi] != HYSL_PLACEHOLDER ) { + n_nonplaceholder++; + } + } + + /* If the dimensions match after excluding + * placeholders, it's OK - probably a static mask + * in a multi-event file. */ + if ( ndims == n_nonplaceholder ) { + skip_placeholders = 1; + } + } + + f_offset = malloc(ndims*sizeof(hsize_t)); + f_count = malloc(ndims*sizeof(hsize_t)); + if ( (f_offset == NULL) || (f_count == NULL ) ) { + ERROR("Failed to allocate offset or count.\n"); + free_event(ev); + H5Fclose(fh); + return 1; + } + + dpos = 0; + for ( hsi=0; hsi<p->dim_structure->num_dims; hsi++ ) { + + switch ( p->dim_structure->dims[hsi] ) { + + case HYSL_FS: + f_offset[dpos] = p->orig_min_fs; + f_count[dpos] = p->orig_max_fs - p->orig_min_fs+1; + dpos++; + break; + + case HYSL_SS: + f_offset[dpos] = p->orig_min_ss; + f_count[dpos] = p->orig_max_ss - p->orig_min_ss+1; + dpos++; + break; + + case HYSL_PLACEHOLDER: + if ( !skip_placeholders ) { + f_offset[dpos] = ev->dim_entries[0]; + f_count[dpos] = 1; + dpos++; + } + break; + + default: + f_offset[dpos] = p->dim_structure->dims[hsi]; + f_count[dpos] = 1; + dpos++; + break; + + } + + } + + check = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, + f_offset, NULL, f_count, NULL); + if ( check < 0 ) { + ERROR("Error selecting file dataspace for panel %s\n", + p->name); + free(f_offset); + free(f_count); + free_event(ev); + H5Fclose(fh); + return 1; + } + + dims[0] = p->orig_max_ss - p->orig_min_ss + 1; + dims[1] = p->orig_max_fs - p->orig_min_fs + 1; + memspace = H5Screate_simple(2, dims, NULL); + + data = malloc(dims[0]*dims[1]*el_size); + if ( data == NULL ) { + ERROR("Failed to allocate panel %s\n", p->name); + free(f_offset); + free(f_count); + free_event(ev); + H5Fclose(fh); + return 1; + } + + r = H5Dread(dh, el_type, memspace, dataspace, + H5P_DEFAULT, data); + if ( r < 0 ) { + ERROR("Couldn't read data for panel %s\n", + p->name); + free(f_offset); + free(f_count); + free(data); + free_event(ev); + H5Fclose(fh); + return 1; + } + + H5Dclose(dh); + H5Sclose(dataspace); + free(f_offset); + free(f_count); + free_event(ev); + H5Fclose(fh); + + *pdata = data; + return 0; +} + + +struct image *image_hdf5_read(DataTemplate *dtempl, + const char *filename, const char *event) +{ + struct image *image; + int i; + + image = image_new(); + if ( image == NULL ) { + ERROR("Couldn't allocate image structure.\n"); + return NULL; + } + + image->dp = malloc(dtempl->n_panels*sizeof(float *)); + if ( image->dp == NULL ) { + ERROR("Failed to allocate data array.\n"); + image_free(image); + return NULL; + } + + /* Set all pointers to NULL for easier clean-up */ + for ( i=0; i<dtempl->n_panels; i++ ) image->dp[i] = NULL; + + for ( i=0; i<dtempl->n_panels; i++ ) { + if ( load_hdf5_hyperslab(&dtempl->panels[i], filename, + event, (void *)&image->dp[i], + H5T_NATIVE_FLOAT, + sizeof(float), 0, + dtempl->panels[i].data) ) + { + ERROR("Failed to load panel data\n"); + image_free(image); + return NULL; + } + } + + image->filename = strdup(filename); + image->ev = safe_strdup(event); + + return image; +} + + +int image_hdf5_read_mask(struct panel_template *p, + const char *filename, const char *event, + int *bad, int mask_good, int mask_bad) +{ + int p_w, p_h; + int *mask = NULL; + long unsigned int j; + + p_w = p->orig_max_fs - p->orig_min_fs + 1; + p_h = p->orig_max_ss - p->orig_min_ss + 1; + + if ( load_hdf5_hyperslab(p, filename, event, + (void *)&mask, H5T_NATIVE_INT, + sizeof(int), 1, p->mask) ) + { + ERROR("Failed to load mask data\n"); + free(mask); + return 1; + } + + for ( j=0; j<p_w*p_h; j++ ) { + + /* Bad if it's missing any of the "good" bits */ + if ( (mask[j] & mask_good) != mask_good ) bad[j] = 1; + + /* Bad if it has any of the "bad" bits. */ + if ( mask[j] & mask_bad ) bad[j] = 1; + + } + + free(mask); + return 0; +} + + +double image_hdf5_get_value(const char *name, const char *filename, + const char *event) +{ + hid_t dh; + hid_t type; + hid_t class; + hid_t sh; + hid_t ms; + hsize_t *f_offset = NULL; + hsize_t *f_count = NULL; + hsize_t m_offset[1]; + hsize_t m_count[1]; + hsize_t msdims[1]; + hsize_t size[64]; + herr_t r; + herr_t check; + int check_pe; + int dim_flag; + int ndims; + int i; + char *subst_name = NULL; + struct event *ev; + hid_t fh; + double val; + + if ( access(filename, R_OK) == -1 ) { + ERROR("File does not exist or cannot be read: %s\n", filename); + return NAN; + } + + fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); + if ( fh < 0 ) { + ERROR("Couldn't open file: %s\n", filename); + return NAN; + } + ev = get_event_from_event_string(event); + if ( (ev == NULL) && (event != NULL) ) { + ERROR("Invalid event identifier '%s'\n", event); + H5Fclose(fh); + return NAN; + } + + subst_name = retrieve_full_path(ev, name); + + check_pe = check_path_existence(fh, subst_name); + if ( check_pe == 0 ) { + ERROR("No such event-based numeric field '%s'\n", subst_name); + return NAN; + } + + dh = H5Dopen2(fh, subst_name, H5P_DEFAULT); + type = H5Dget_type(dh); + class = H5Tget_class(type); + + if ( (class != H5T_FLOAT) && (class != H5T_INTEGER) ) { + ERROR("Not a floating point or integer value.\n"); + H5Tclose(type); + H5Dclose(dh); + return NAN; + } + + /* Get the dimensionality. We have to cope with scalars expressed as + * arrays with all dimensions 1, as well as zero-d arrays. */ + sh = H5Dget_space(dh); + ndims = H5Sget_simple_extent_ndims(sh); + if ( ndims > 64 ) { + ERROR("Too many dimensions for numeric value\n"); + H5Tclose(type); + H5Dclose(dh); + return NAN; + } + H5Sget_simple_extent_dims(sh, size, NULL); + + m_offset[0] = 0; + m_count[0] = 1; + msdims[0] = 1; + ms = H5Screate_simple(1,msdims,NULL); + + /* Check that the size in all dimensions is 1 + * or that one of the dimensions has the same + * size as the hyperplane events */ + + dim_flag = 0; + + for ( i=0; i<ndims; i++ ) { + if ( size[i] == 1 ) continue; + if ( ( i==0 ) && (ev != NULL) && (size[i] > ev->dim_entries[0]) ) { + dim_flag = 1; + } else { + H5Tclose(type); + H5Dclose(dh); + return NAN; + } + } + + if ( dim_flag == 0 ) { + + if ( H5Dread(dh, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, &val) < 0 ) + { + ERROR("Couldn't read value.\n"); + H5Tclose(type); + H5Dclose(dh); + return NAN; + } + + } else { + + f_offset = malloc(ndims*sizeof(hsize_t)); + f_count = malloc(ndims*sizeof(hsize_t)); + + for ( i=0; i<ndims; i++ ) { + + if ( i == 0 ) { + f_offset[i] = ev->dim_entries[0]; + f_count[i] = 1; + } else { + f_offset[i] = 0; + f_count[i] = 0; + } + + } + + check = H5Sselect_hyperslab(sh, H5S_SELECT_SET, + f_offset, NULL, f_count, NULL); + if ( check <0 ) { + ERROR("Error selecting dataspace for float value"); + free(f_offset); + free(f_count); + return NAN; + } + + ms = H5Screate_simple(1,msdims,NULL); + check = H5Sselect_hyperslab(ms, H5S_SELECT_SET, + m_offset, NULL, m_count, NULL); + if ( check < 0 ) { + ERROR("Error selecting memory dataspace for float value"); + free(f_offset); + free(f_count); + return NAN; + } + + r = H5Dread(dh, H5T_NATIVE_DOUBLE, ms, sh, H5P_DEFAULT, &val); + if ( r < 0 ) { + ERROR("Couldn't read value.\n"); + H5Tclose(type); + H5Dclose(dh); + return NAN; + } + + } + + free_event(ev); + free(subst_name); + H5Fclose(fh); + + return val; +} + + +static int read_peak_count(hid_t fh, char *path, int line, + int *num_peaks) +{ + + hid_t dh, sh, mh; + hsize_t size[1]; + hsize_t max_size[1]; + hsize_t offset[1], count[1]; + hsize_t m_offset[1], m_count[1], dimmh[1]; + int tw, r; + + dh = H5Dopen2(fh, path, H5P_DEFAULT); + if ( dh < 0 ) { + ERROR("Data block %s not found.\n", path); + return 1; + } + + sh = H5Dget_space(dh); + if ( sh < 0 ) { + H5Dclose(dh); + ERROR("Couldn't get dataspace for data.\n"); + return 1; + } + + if ( H5Sget_simple_extent_ndims(sh) != 1 ) { + ERROR("Data block %s has the wrong dimensionality (%i).\n", + path, H5Sget_simple_extent_ndims(sh)); + H5Sclose(sh); + H5Dclose(dh); + return 1; + } + + H5Sget_simple_extent_dims(sh, size, max_size); + + tw = size[0]; + + if ( line > tw-1 ) { + H5Sclose(sh); + H5Dclose(dh); + ERROR("Data block %s does not contain data for required event.\n", + path); + return 1; + } + + offset[0] = line; + count[0] = 1; + + r = H5Sselect_hyperslab(sh, H5S_SELECT_SET, + offset, NULL, count, NULL); + if ( r < 0 ) { + ERROR("Error selecting file dataspace " + "for data block %s\n", path); + H5Dclose(dh); + H5Sclose(sh); + return 1; + } + + m_offset[0] = 0; + m_count[0] = 1; + dimmh[0] = 1; + mh = H5Screate_simple(1, dimmh, NULL); + r = H5Sselect_hyperslab(mh, H5S_SELECT_SET, + m_offset, NULL, m_count, NULL); + if ( r < 0 ) { + ERROR("Error selecting memory dataspace " + "for data block %s\n", path); + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return 1; + } + + r = H5Dread(dh, H5T_NATIVE_INT, mh, + sh, H5P_DEFAULT, num_peaks); + if ( r < 0 ) { + ERROR("Couldn't read data for block %s, line %i\n", path, line); + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return 1; + } + + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return 0; +} + + +static float *read_hdf5_data(hid_t fh, char *path, int line) +{ + + hid_t dh, sh, mh; + hsize_t size[2]; + hsize_t max_size[2]; + hsize_t offset[2], count[2]; + hsize_t m_offset[2], m_count[2], dimmh[2]; + float *buf; + int tw, r; + + dh = H5Dopen2(fh, path, H5P_DEFAULT); + if ( dh < 0 ) { + ERROR("Data block (%s) not found.\n", path); + return NULL; + } + + sh = H5Dget_space(dh); + if ( sh < 0 ) { + H5Dclose(dh); + ERROR("Couldn't get dataspace for data.\n"); + return NULL; + } + + if ( H5Sget_simple_extent_ndims(sh) != 2 ) { + ERROR("Data block %s has the wrong dimensionality (%i).\n", + path, H5Sget_simple_extent_ndims(sh)); + H5Sclose(sh); + H5Dclose(dh); + return NULL; + } + + H5Sget_simple_extent_dims(sh, size, max_size); + + tw = size[0]; + if ( line> tw-1 ) { + H5Sclose(sh); + H5Dclose(dh); + ERROR("Data block %s does not contain data for required event.\n", + path); + return NULL; + } + + offset[0] = line; + offset[1] = 0; + count[0] = 1; + count[1] = size[1]; + + r = H5Sselect_hyperslab(sh, H5S_SELECT_SET, offset, NULL, count, NULL); + if ( r < 0 ) { + ERROR("Error selecting file dataspace " + "for data block %s\n", path); + H5Dclose(dh); + H5Sclose(sh); + return NULL; + } + + m_offset[0] = 0; + m_offset[1] = 0; + m_count[0] = 1; + m_count[1] = size[1]; + dimmh[0] = 1; + dimmh[1] = size[1]; + + mh = H5Screate_simple(2, dimmh, NULL); + r = H5Sselect_hyperslab(mh, H5S_SELECT_SET, + m_offset, NULL, m_count, NULL); + if ( r < 0 ) { + ERROR("Error selecting memory dataspace " + "for data block %s\n", path); + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return NULL; + } + + buf = malloc(size[1]*sizeof(float)); + if ( buf == NULL ) return NULL; + r = H5Dread(dh, H5T_NATIVE_FLOAT, mh, sh, H5P_DEFAULT, buf); + if ( r < 0 ) { + ERROR("Couldn't read data for block %s, line %i\n", path, line); + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return NULL; + } + + H5Dclose(dh); + H5Sclose(sh); + H5Sclose(mh); + return buf; +} + + +ImageFeatureList *image_hdf5_read_peaks_cxi(const DataTemplate *dtempl, + const char *filename, + const char *event, + int half_pixel_shift) +{ + ImageFeatureList *features; + hid_t fh; + char path_n[1024]; + char path_x[1024]; + char path_y[1024]; + char path_i[1024]; + int r; + int pk; + struct event *ev; + char *subst_name; + + int line = 0; + int num_peaks; + + float *buf_x; + float *buf_y; + float *buf_i; + + double peak_offset = half_pixel_shift ? 0.5 : 0.0; + + if ( access(filename, R_OK) == -1 ) { + ERROR("File does not exist or cannot be read: %s\n", + filename); + return NULL; + } + + fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); + if ( fh < 0 ) { + ERROR("Couldn't open file: %s\n", filename); + return NULL; + } + + ev = get_event_from_event_string(event); + if ( (ev == NULL) && (event != NULL) ) { + ERROR("Invalid event identifier '%s'\n", event); + H5Fclose(fh); + return NULL; + } + + if ( ev->dim_entries == NULL ) { + ERROR("CXI format peak list format selected," + "but file has no event structure"); + return NULL; + } + line = ev->dim_entries[0]; + + subst_name = retrieve_full_path(ev, dtempl->peak_list); + free_event(ev); + if ( subst_name == NULL ) { + ERROR("Invalid peak path %s\n", subst_name); + H5Fclose(fh); + return NULL; + } + + if ( check_path_existence(fh, subst_name) == 0 ) { + ERROR("Peak path not found: %s:%s", + filename, subst_name); + free(subst_name); + H5Fclose(fh); + return NULL; + } + + snprintf(path_n, 1024, "%s/nPeaks", subst_name); + snprintf(path_x, 1024, "%s/peakXPosRaw", subst_name); + snprintf(path_y, 1024, "%s/peakYPosRaw", subst_name); + snprintf(path_i, 1024, "%s/peakTotalIntensity", subst_name); + + r = read_peak_count(fh, path_n, line, &num_peaks); + if ( r != 0 ) return NULL; + + buf_x = read_hdf5_data(fh, path_x, line); + if ( r != 0 ) return NULL; + + buf_y = read_hdf5_data(fh, path_y, line); + if ( r != 0 ) return NULL; + + buf_i = read_hdf5_data(fh, path_i, line); + if ( r != 0 ) return NULL; + + features = image_feature_list_new(); + + for ( pk=0; pk<num_peaks; pk++ ) { + + float fs, ss, val; + int pn; + + fs = buf_x[pk] + peak_offset; + ss = buf_y[pk] + peak_offset; + val = buf_i[pk]; + + if ( data_template_file_to_panel_coords(dtempl, + &fs, &ss, + &pn) ) + { + ERROR("Failed to convert %i,%i to " + "panel-relative coordinates\n", fs, ss); + } else { + image_add_feature(features, fs, ss, pn, + NULL, val, NULL); + } + + } + + return features; +} + + +ImageFeatureList *image_hdf5_read_peaks_hdf5(const DataTemplate *dtempl, + const char *filename, + const char *event, + int half_pixel_shift) +{ + hid_t fh, dh, sh; + hsize_t size[2]; + hsize_t max_size[2]; + int i; + float *buf; + herr_t r; + int tw; + struct event *ev; + char *subst_name; + ImageFeatureList *features; + double peak_offset = half_pixel_shift ? 0.5 : 0.0; + + if ( access(filename, R_OK) == -1 ) { + ERROR("File does not exist or cannot be read: %s\n", + filename); + return NULL; + } + + fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); + if ( fh < 0 ) { + ERROR("Couldn't open file: %s\n", filename); + return NULL; + } + + ev = get_event_from_event_string(event); + if ( (ev == NULL) && (event != NULL) ) { + ERROR("Invalid event identifier '%s'\n", event); + H5Fclose(fh); + return NULL; + } + + subst_name = retrieve_full_path(ev, dtempl->peak_list); + free_event(ev); + if ( subst_name == NULL ) { + ERROR("Invalid peak path %s\n", subst_name); + free_event(ev); + H5Fclose(fh); + return NULL; + } + + if ( check_path_existence(fh, subst_name) == 0 ) { + ERROR("Peak path not found: %s:%s", + filename, subst_name); + free(subst_name); + H5Fclose(fh); + return NULL; + } + + dh = H5Dopen2(fh, subst_name, H5P_DEFAULT); + free(subst_name); + if ( dh < 0 ) { + ERROR("Peak list (%s) not found.\n", subst_name); + H5Fclose(fh); + return NULL; + } + + sh = H5Dget_space(dh); + if ( sh < 0 ) { + ERROR("Couldn't get dataspace for peak list.\n"); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + + if ( H5Sget_simple_extent_ndims(sh) != 2 ) { + ERROR("Peak list has the wrong dimensionality (%i).\n", + H5Sget_simple_extent_ndims(sh)); + H5Sclose(sh); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + + H5Sget_simple_extent_dims(sh, size, max_size); + H5Sclose(sh); + + tw = size[1]; + if ( (tw != 3) && (tw != 4) ) { + ERROR("Peak list has the wrong dimensions.\n"); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + + buf = malloc(sizeof(float)*size[0]*size[1]); + if ( buf == NULL ) { + ERROR("Couldn't reserve memory for the peak list.\n"); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + r = H5Dread(dh, H5T_NATIVE_FLOAT, H5S_ALL, H5S_ALL, + H5P_DEFAULT, buf); + if ( r < 0 ) { + ERROR("Couldn't read peak list.\n"); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + + features = image_feature_list_new(); + if ( features == NULL ) { + ERROR("Failed to allocate peak list\n"); + H5Dclose(dh); + H5Fclose(fh); + return NULL; + } + + for ( i=0; i<size[0]; i++ ) { + + float fs, ss, val; + int pn; + + fs = buf[tw*i+0] + peak_offset; + ss = buf[tw*i+1] + peak_offset; + val = buf[tw*i+2]; + + if ( data_template_file_to_panel_coords(dtempl, + &fs, &ss, + &pn) ) + { + ERROR("Failed to convert %i,%i to " + "panel-relative coordinates\n", fs, ss); + } else { + image_add_feature(features, fs, ss, pn, + NULL, val, NULL); + } + + } + + free(buf); + H5Dclose(dh); + H5Fclose(fh); + + return features; +} + + +struct parse_params { + hid_t fh; + int path_dim; + const char *path; + struct event *curr_event; + struct event_list *ev_list; + int top_level; +}; + + +static herr_t parse_file_event_structure(hid_t loc_id, char *name, + const H5L_info_t *info, + struct parse_params *pp) + +{ + char *substituted_path; + char *ph_loc; + char *truncated_path; + htri_t check; + herr_t herrt_iterate, herrt_info; + struct H5O_info_t object_info; + + if ( !pp->top_level ) { + + int fail_push; + + fail_push = push_path_entry_to_event(pp->curr_event, name); + if ( fail_push ) { + return -1; + } + + substituted_path = event_path_placeholder_subst(name, pp->path); + + } else { + substituted_path = strdup(pp->path); + } + + if ( pp->top_level == 1 ) { + pp->top_level = 0; + } + + truncated_path = strdup(substituted_path); + ph_loc = strstr(substituted_path,"%"); + if ( ph_loc != NULL) { + truncated_path[ph_loc-substituted_path] = '\0'; + } + + herrt_iterate = 0; + herrt_info = 0; + + check = check_path_existence(pp->fh, truncated_path); + if ( check == 0 ) { + pop_path_entry_from_event(pp->curr_event); + return 0; + } else { + + herrt_info = H5Oget_info_by_name(pp->fh, truncated_path, + &object_info, H5P_DEFAULT); + if ( herrt_info < 0 ) { + free(truncated_path); + free(substituted_path); + return -1; + } + + if ( pp->curr_event->path_length == pp->path_dim + && object_info.type == H5O_TYPE_DATASET ) + { + + int fail_append; + + fail_append = append_event_to_event_list(pp->ev_list, + pp->curr_event); + if ( fail_append ) { + free(truncated_path); + free(substituted_path); + return -1; + } + + pop_path_entry_from_event(pp->curr_event); + return 0; + + } else { + + pp->path = substituted_path; + + if ( object_info.type == H5O_TYPE_GROUP ) { + + herrt_iterate = H5Literate_by_name(pp->fh, + truncated_path, H5_INDEX_NAME, + H5_ITER_NATIVE, NULL, + (H5L_iterate_t)parse_file_event_structure, + (void *)pp, H5P_DEFAULT); + } + } + } + + pop_path_entry_from_event(pp->curr_event); + + free(truncated_path); + free(substituted_path); + + return herrt_iterate; +} + + +static int fill_paths(hid_t fh, const DataTemplate *dtempl, int pi, + struct event_list *master_el) +{ + struct parse_params pparams; + struct event *empty_event; + struct event_list *panel_ev_list; + int ei; + int check; + + empty_event = initialize_event(); + panel_ev_list = initialize_event_list(); + if ( (empty_event == NULL) || (panel_ev_list == NULL) ) + { + ERROR("Failed to allocate memory for event list.\n"); + return 1; + } + + pparams.path = dtempl->panels[pi].data; + pparams.fh = fh; + pparams.path_dim = dtempl->path_dim; + pparams.curr_event = empty_event; + pparams.top_level = 1; + pparams.ev_list = panel_ev_list; + + check = parse_file_event_structure(fh, NULL, NULL, &pparams); + if ( check < 0 ) { + free_event(empty_event); + free_event_list(panel_ev_list); + return 1; + } + + for ( ei=0; ei<panel_ev_list->num_events; ei++ ) { + + int fail_add; + + fail_add = add_non_existing_event_to_event_list(master_el, + panel_ev_list->events[ei]); + if ( fail_add ) { + free_event(empty_event); + free_event_list(panel_ev_list); + return 1; + } + + } + + free_event(empty_event); + free_event_list(panel_ev_list); + + return 0; +} + + +static int check_dims(hid_t fh, struct panel_template *p, + struct event *ev, + struct event_list *events, int *global_path_dim) +{ + char *full_panel_path; + hid_t dh; + hid_t sh; + int dims; + hsize_t *size; + hsize_t *max_size; + int hsdi; + int panel_path_dim = 0; + struct dim_structure *panel_dim_structure; + + /* Get the full path for this panel in this event */ + full_panel_path = retrieve_full_path(ev, p->data); + + dh = H5Dopen2(fh, full_panel_path, H5P_DEFAULT); + if ( dh < 0 ) { + ERROR("Error opening '%s'\n", full_panel_path); + ERROR("Failed to enumerate events. " + "Check your geometry file.\n"); + return 1; + } + + sh = H5Dget_space(dh); + dims = H5Sget_simple_extent_ndims(sh); + size = malloc(dims*sizeof(hsize_t)); + max_size = malloc(dims*sizeof(hsize_t)); + if ( (size==NULL) || (max_size==NULL) ) { + ERROR("Failed to allocate memory for dimensions\n"); + return 1; + } + + dims = H5Sget_simple_extent_dims(sh, size, max_size); + + panel_dim_structure = p->dim_structure; + for ( hsdi=0; hsdi<panel_dim_structure->num_dims; hsdi++ ) { + if ( panel_dim_structure->dims[hsdi] == HYSL_PLACEHOLDER ) { + panel_path_dim = size[hsdi]; + break; + } + } + + if ( *global_path_dim == -1 ) { + + *global_path_dim = panel_path_dim; + + } else if ( panel_path_dim != *global_path_dim ) { + + ERROR("All panels must have the same number of frames\n"); + ERROR("Panel %s has %i frames in one dimension, but the first " + "panel has %i.\n", + p->name, panel_path_dim, *global_path_dim); + free(size); + free(max_size); + return 1; + } + + H5Sclose(sh); + H5Dclose(dh); + + return 0; +} + + +struct event_list *image_hdf5_expand_frames(const DataTemplate *dtempl, + const char *filename) +{ + struct event_list *master_el; + hid_t fh; + + fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); + if ( fh < 0 ) { + ERROR("Couldn't open file '%s'\n", filename); + return NULL; + } + + master_el = initialize_event_list(); + if ( master_el == NULL ) { + ERROR("Failed to allocate event list.\n"); + H5Fclose(fh); + return NULL; + } + + /* First expand any placeholders in the HDF5 paths */ + if ( dtempl->path_dim != 0 ) { + int pi; + for ( pi=0; pi<dtempl->n_panels; pi++ ) { + if ( fill_paths(fh, dtempl, pi, master_el) ) { + ERROR("Failed to enumerate paths.\n"); + H5Fclose(fh); + return NULL; + } + } + } + + /* Now enumerate the placeholder dimensions */ + if ( dtempl->dim_dim > 0 ) { + + struct event_list *master_el_with_dims; + int evi; + + /* If there were no HDF5 path placeholders, add a dummy event */ + if ( master_el->num_events == 0 ) { + struct event *empty_ev; + empty_ev = initialize_event(); + append_event_to_event_list(master_el, empty_ev); + free(empty_ev); + } + + master_el_with_dims = initialize_event_list(); + + /* For each event so far, expand the dimensions */ + for ( evi=0; evi<master_el->num_events; evi++ ) { + + int pi; + int global_path_dim = -1; + int mlwd; + + /* Check the dimensionality of each panel */ + for ( pi=0; pi<dtempl->n_panels; pi++ ) { + if ( check_dims(fh, &dtempl->panels[pi], + master_el->events[evi], + master_el_with_dims, + &global_path_dim) ) + { + ERROR("Failed to enumerate dims.\n"); + H5Fclose(fh); + return NULL; + } + } + + /* Add this dimensionality to all events */ + for ( mlwd=0; mlwd<global_path_dim; mlwd++ ) { + + struct event *mlwd_ev; + + mlwd_ev = copy_event(master_el->events[evi]); + push_dim_entry_to_event(mlwd_ev, mlwd); + append_event_to_event_list(master_el_with_dims, + mlwd_ev); + free_event(mlwd_ev); + } + + } + + free_event_list(master_el); + H5Fclose(fh); + return master_el_with_dims; + + } else { + + H5Fclose(fh); + return master_el; + + } +} diff --git a/libcrystfel/src/image-hdf5.h b/libcrystfel/src/image-hdf5.h new file mode 100644 index 00000000..93f887e6 --- /dev/null +++ b/libcrystfel/src/image-hdf5.h @@ -0,0 +1,67 @@ +/* + * image-hdf5.h + * + * Image loading, HDF5 parts + * + * Copyright © 2012-2020 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2020 Thomas White <taw@physics.org> + * + * This file is part of CrystFEL. + * + * CrystFEL is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * CrystFEL is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with CrystFEL. If not, see <http://www.gnu.org/licenses/>. + * + */ + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +/* NB This file is NOT part of the public API, and should NOT + * be installed, but rather stays in the libcrystfel source folder. */ + +#ifndef IMAGE_HDF5_H +#define IMAGE_HDF5_H + +#include "datatemplate_priv.h" + +extern double image_hdf5_get_value(const char *from, + const char *filename, + const char *ev); + +extern struct image *image_hdf5_read(DataTemplate *dtempl, + const char *filename, + const char *event); + +extern int image_hdf5_read_mask(struct panel_template *p, + const char *filename, + const char *event, int *bad, + int mask_good, int mask_bad); + +extern ImageFeatureList *image_hdf5_read_peaks_cxi(const DataTemplate *dtempl, + const char *filename, + const char *event, + int half_pixel_shift); + +extern ImageFeatureList *image_hdf5_read_peaks_hdf5(const DataTemplate *dtempl, + const char *filename, + const char *event, + int half_pixel_shift); + +extern struct event_list *image_hdf5_expand_frames(const DataTemplate *dtempl, + const char *filename); + +#endif /* IMAGE_HDF5_H */ diff --git a/libcrystfel/src/image.c b/libcrystfel/src/image.c index 898f4863..46b110ac 100644 --- a/libcrystfel/src/image.c +++ b/libcrystfel/src/image.c @@ -33,15 +33,13 @@ #include <assert.h> #include <math.h> #include <stdio.h> -#include <hdf5.h> -#include <zlib.h> #include "image.h" #include "utils.h" #include "events.h" -#include "hdf5-file.h" -#include "detector.h" #include "detgeom.h" +#include "image-hdf5.h" +#include "image-cbf.h" #include "datatemplate.h" #include "datatemplate_priv.h" @@ -280,1068 +278,6 @@ void free_all_crystals(struct image *image) } -static void add_out(float val, float *data_out, int nmemb_out, - int *outpos, int *nrej) -{ - if ( *outpos < nmemb_out ) { - data_out[(*outpos)++] = val; - } else { - (*nrej)++; - } -} - - -/* Reverses byte offset compression and converts to single precision float. - * Note that this compression scheme specifies the data format of the input - * data, therefore the X-Binary-Element-Type is completely ignored. */ -static void decode_cbf_byte_offset(float *data_out, int nmemb_out, - const int8_t *data_in, const size_t n) -{ - int inpos = 0; - int outpos = 0; - int nrej = 0; - float val = 0.0; - - while ( inpos < n ) { - - int64_t delta = data_in[inpos++]; - - if ( (delta >= -127) && (delta <= 127) ) { - val += delta; - add_out(val, data_out, nmemb_out, &outpos, &nrej); - continue; - } - - delta = *(int16_t *)(data_in+inpos); - inpos += 2; - - if ( (delta >= -32767) && (delta <= 32767) ) { - val += delta; - add_out(val, data_out, nmemb_out, &outpos, &nrej); - continue; - } - - delta = *(int32_t *)(data_in+inpos); - inpos += 4; - - if ( (delta >= -2147483647) && (delta <= 2147483647) ) { - val += delta; - add_out(val, data_out, nmemb_out, &outpos, &nrej); - continue; - } - - delta = *(int64_t *)(data_in+inpos); - inpos += 8; - val += delta; - add_out(val, data_out, nmemb_out, &outpos, &nrej); - - } - - if ( nrej > 0 ) { - STATUS("%i elements rejected\n", nrej); - } -} - - -static int binary_start(char *data) -{ - char *datac = data; - if ( (datac[0] == (char)0x0c) && (datac[1] == (char)0x1a) - && (datac[2] == (char)0x04) && (datac[3] == (char)0xd5) ) return 1; - return 0; -} - - -enum cbf_data_conversion -{ - CBF_NO_CONVERSION, - CBF_BYTE_OFFSET, - CBF_PACKED, - CBF_CANONICAL -}; - -enum cbf_data_type -{ - CBF_NO_TYPE, - CBF_ELEMENT_U8, - CBF_ELEMENT_S8, - CBF_ELEMENT_U16, - CBF_ELEMENT_S16, - CBF_ELEMENT_U32, - CBF_ELEMENT_S32, - CBF_ELEMENT_F32, - CBF_ELEMENT_F64, -}; - - -static enum cbf_data_type parse_element_type(const char *t) -{ - if ( strstr(t, "signed 8-bit integer") != NULL ) - { - return CBF_ELEMENT_S8; - } - - if ( strstr(t, "unsigned 8-bit integer") != NULL ) - { - return CBF_ELEMENT_U8; - } - - if ( strstr(t, "signed 16-bit integer") != NULL ) - { - return CBF_ELEMENT_S16; - } - - if ( strstr(t, "unsigned 16-bit integer") != NULL ) - { - return CBF_ELEMENT_U16; - } - - if ( strstr(t, "signed 32-bit integer") != NULL ) - { - return CBF_ELEMENT_S32; - } - - if ( strstr(t, "unsigned 32-bit integer") != NULL ) - { - return CBF_ELEMENT_U32; - } - - if ( strstr(t, "signed 32-bit real IEEE") != NULL ) - { - return CBF_ELEMENT_F32; - } - - if ( strstr(t, "signed 64-bit real IEEE") != NULL ) - { - return CBF_ELEMENT_F64; - } - - /* complex type is unsupported */ - - return CBF_NO_TYPE; -} - - -static size_t element_size(enum cbf_data_type t) -{ - switch ( t ) { - case CBF_ELEMENT_S8 : return 1; - case CBF_ELEMENT_U8 : return 1; - case CBF_ELEMENT_S16 : return 2; - case CBF_ELEMENT_U16 : return 2; - case CBF_ELEMENT_S32 : return 4; - case CBF_ELEMENT_U32 : return 4; - case CBF_ELEMENT_F32 : return 4; - case CBF_ELEMENT_F64 : return 8; - default : return 0; - } -} - - - -static int convert_type(float *data_out, long nmemb_exp, - enum cbf_data_type eltype, - void *data_in, size_t data_in_len) -{ - long int i; - long int o = 0; - size_t elsize = element_size(eltype); - - if ( elsize == 0 ) return 1; - - if ( nmemb_exp * elsize > data_in_len ) { - ERROR("Not enough CBF data for image size/type!\n"); - return 1; - } - - for ( i=0; i<nmemb_exp; i++ ) { - switch ( eltype ) { - - case CBF_ELEMENT_S8: - data_out[o++] = ((int8_t *)data_in)[i]; - break; - - case CBF_ELEMENT_U8: - data_out[o++] = ((uint8_t *)data_in)[i]; - break; - - case CBF_ELEMENT_S16: - data_out[o++] = ((int16_t *)data_in)[i]; - break; - - case CBF_ELEMENT_U16: - data_out[o++] = ((uint16_t *)data_in)[i]; - break; - - case CBF_ELEMENT_S32: - data_out[o++] = ((int32_t *)data_in)[i]; - break; - - case CBF_ELEMENT_U32: - data_out[o++] = ((uint32_t *)data_in)[i]; - break; - - case CBF_ELEMENT_F32: - data_out[o++] = ((float *)data_in)[i]; - break; - - case CBF_ELEMENT_F64: - data_out[o++] = ((double *)data_in)[i]; - break; - - case CBF_NO_TYPE: - break; - } - } - - return 0; -} - - -static float *read_cbf_data(const char *filename, int gz, int *w, int *h) -{ - FILE *fh; - void *buf = NULL; - char *rval; - size_t data_compressed_len = 0; - float *data_out = NULL; - enum cbf_data_conversion data_conversion = CBF_NO_CONVERSION; - enum cbf_data_type data_type = CBF_ELEMENT_U32; /* ITG (2006) 2.3.3.3 */ - int in_binary_section = 0; - - *w = 0; - *h = 0; - - if ( !gz ) { - - fh = fopen(filename, "rb"); - if ( fh == NULL ) { - ERROR("Failed to open '%s'\n", filename); - return NULL; - } - - } else { - - gzFile gzfh; - size_t len, len_read; - const size_t bufinc = 8*1024*1024; /* Allocate buffer in 8Mb chunks */ - size_t bufsz = bufinc; - - gzfh = gzopen(filename, "rb"); - if ( gzfh == NULL ) return NULL; - - #ifdef HAVE_GZBUFFER - /* Set larger buffer size for hopefully faster uncompression */ - gzbuffer(gzfh, 128*1024); - #endif - - buf = malloc(bufsz); - if ( buf == NULL ) return NULL; - - len = 0; - do { - - len_read = gzread(gzfh, buf+len, bufinc); - if ( len_read == -1 ) return NULL; - len += len_read; - - if ( len_read == bufinc ) { - bufsz += bufinc; - buf = realloc(buf, bufsz); - if ( buf == NULL ) return NULL; - } - - } while ( len_read == bufinc ); - - fh = fmemopen(buf, len, "rb"); - if ( fh == NULL ) return NULL; - - gzclose(gzfh); - - } - - /* This is really horrible, but there are at least three different types - * of header mingled together (CIF, MIME, DECTRIS), so a real parser - * would be very complicated and much more likely to have weird bugs. */ - do { - - char line[1024]; - long line_start; - - line_start = ftell(fh); - rval = fgets(line, 1023, fh); - if ( rval == NULL ) break; - chomp(line); - - if ( strcmp(line, "--CIF-BINARY-FORMAT-SECTION--") == 0 ) { - in_binary_section = 1; - } - - if ( strcmp(line, "--CIF-BINARY-FORMAT-SECTION----") == 0 ) { - in_binary_section = 0; - } - - if ( in_binary_section ) { - - if ( strncmp(line, "X-Binary-Size: ", 15) == 0 ) { - data_compressed_len = atoi(line+15); - } - - if ( strncmp(line, "X-Binary-Element-Byte-Order: ", 29) == 0 ) { - const char *elbo = line+29; - if ( strcmp(elbo, "LITTLE_ENDIAN") != 0 ) { - ERROR("Unsupported endianness: %s\n", elbo); - free(buf); - fclose(fh); - return NULL; - } - } - - /* Try to spot compression algorithm */ - if ( strstr(line, "conversions=\"x-CBF_BYTE_OFFSET\"") != NULL ) { - data_conversion = CBF_BYTE_OFFSET; - } else if ( strstr(line, "conversions=\"x-CBF_CANONICAL\"") != NULL ) { - data_conversion = CBF_CANONICAL; - } else if ( strstr(line, "conversions=\"x-CBF_PACKED\"") != NULL ) { - data_conversion = CBF_PACKED; - } else if ( strstr(line, "conversions=") != NULL ) { - ERROR("Unrecognised CBF content conversion: %s\n", line); - free(buf); - fclose(fh); - return NULL; - } - - /* Likewise, element type */ - if ( strncmp(line, "X-Binary-Element-Type: ", 23) == 0 ) - { - const char *eltype = (line+23); - data_type = parse_element_type(eltype); - if ( data_type == CBF_NO_TYPE ) { - ERROR("Unrecognised element type: %s\n", - eltype); - free(buf); - fclose(fh); - return NULL; - } - } - - if ( strncmp(line, "X-Binary-Size-Fastest-Dimension: ", 33) == 0 ) { - *w = atoi(line+33); - } - - if ( strncmp(line, "X-Binary-Size-Second-Dimension: ", 32) == 0 ) { - *h = atoi(line+32); - } - - } - - if ( in_binary_section && binary_start(line) ) { - - size_t len_read; - int nmemb_exp; - void *data_compressed; - int r = 0; - - if ( data_compressed_len == 0 ) { - ERROR("Found CBF data before X-Binary-Size!\n"); - free(buf); - fclose(fh); - return NULL; - } - - if ( (*w == 0) || (*h == 0) ) { - ERROR("Found CBF data before dimensions!\n"); - free(buf); - fclose(fh); - return NULL; - } - - if ( data_compressed_len > 100*1024*1024 ) { - ERROR("Stated CBF data size too big\n"); - free(buf); - fclose(fh); - return NULL; - } - - data_compressed = malloc(data_compressed_len); - if ( data_compressed == NULL ) { - ERROR("Failed to allocate memory for CBF data\n"); - free(buf); - fclose(fh); - return NULL; - } - - fseek(fh, line_start+4, SEEK_SET); - len_read = fread(data_compressed, 1, data_compressed_len, fh); - if ( len_read < data_compressed_len ) { - ERROR("Couldn't read entire CBF data\n"); - free(buf); - free(data_compressed); - fclose(fh); - return NULL; - } - - nmemb_exp = (*w) * (*h); - data_out = malloc(nmemb_exp*sizeof(float)); - if ( data_out == NULL ) { - ERROR("Failed to allocate memory for CBF data\n"); - free(buf); - free(data_compressed); - fclose(fh); - return NULL; - } - - switch ( data_conversion ) { - - case CBF_NO_CONVERSION: - r = convert_type(data_out, nmemb_exp, data_type, - data_compressed, - data_compressed_len); - break; - - case CBF_BYTE_OFFSET: - decode_cbf_byte_offset(data_out, nmemb_exp, - data_compressed, - data_compressed_len); - break; - - case CBF_PACKED: - case CBF_CANONICAL: - ERROR("Don't yet know how to decompress " - "CBF_PACKED or CBF_CANONICAL\n"); - free(buf); - free(data_compressed); - fclose(fh); - return NULL; - - } - - free(data_compressed); - - if ( r ) { - free(buf); - free(data_out); - fclose(fh); - return NULL; - } - - free(buf); - fclose(fh); - return data_out; - - } - - } while ( rval != NULL ); - - ERROR("Reached end of CBF file before finding data.\n"); - free(buf); /* might be NULL */ - return NULL; -} - - -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; -} - - -signed int is_cbfgz_file(const char *filename) -{ - gzFile gzfh; - char line[1024]; - - gzfh = gzopen(filename, "rb"); - if ( gzfh == NULL ) return -1; - if ( gzgets(gzfh, line, 1024) == NULL ) return -1; - gzclose(gzfh); - - if ( strstr(line, "CBF") == NULL ) { - return 0; - } - - return 1; -} - - -/************************** New API (DataTemplate) ****************************/ - -struct image *image_new() -{ - struct image *image; - - image = malloc(sizeof(struct image)); - if ( image == NULL ) return NULL; - - image->dp = NULL; - image->bad = NULL; - image->sat = NULL; - image->hit = 0; - image->crystals = NULL; - image->n_crystals = 0; - image->indexed_by = INDEXING_NONE; - image->detgeom = NULL; - image->filename = NULL; - image->ev = NULL; - image->copied_headers = NULL; - image->avg_clen = -1.0; - image->id = 0; - image->serial = 0; - image->spectrum = NULL; - image->lambda = -1.0; - image->div = -1.0; - image->bw = -1.0; - image->peak_resolution = -1.0; - image->features = NULL; - - /* Deprecated stuff */ - image->beam = NULL; - image->det = NULL; - - return image; -} - - -static int unpack_panels_dtempl(struct image *image, - DataTemplate *dtempl, - float *data, - int data_width, int data_height) -{ - int pi; - - image->dp = malloc(dtempl->n_panels * sizeof(float *)); - if ( image->dp == NULL ) { - ERROR("Failed to allocate panels.\n"); - return 1; - } - - for ( pi=0; pi<dtempl->n_panels; pi++ ) { - - struct panel_template *p; - int fs, ss; - int p_w, p_h; - - p = &dtempl->panels[pi]; - p_w = p->orig_max_fs - p->orig_min_fs + 1; - p_h = p->orig_max_ss - p->orig_min_ss + 1; - image->dp[pi] = malloc(p_w*p_h*sizeof(float)); - if ( image->dp[pi] == NULL ) { - ERROR("Failed to allocate panel\n"); - return 1; - } - - 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; - - 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]; - - } - } - - } - - return 0; -} - - -static int load_hdf5_hyperslab(struct panel_template *p, - const char *filename, - const char *event, - void **pdata, - hid_t el_type, size_t el_size, - int skip_placeholders_ok, - const char *path_spec) -{ - struct event *ev; - hid_t fh; - herr_t r; - hsize_t *f_offset, *f_count; - hid_t dh; - int hsi; - herr_t check; - hid_t dataspace, memspace; - hsize_t dims[2]; - char *panel_full_path; - void *data; - int ndims, dpos; - int skip_placeholders = 0; - - if ( access(filename, R_OK) == -1 ) { - ERROR("File does not exist or cannot be read: %s\n", - filename); - return 1; - } - - fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); - if ( fh < 0 ) { - ERROR("Couldn't open file: %s\n", filename); - return 1; - } - - ev = get_event_from_event_string(event); - if ( (ev == NULL) && (event != NULL) ) { - ERROR("Invalid event identifier '%s'\n", event); - H5Fclose(fh); - return 1; - } - - panel_full_path = retrieve_full_path(ev, path_spec); - - if ( !check_path_existence(fh, panel_full_path) ) { - ERROR("Cannot find data for panel %s (%s)\n", - p->name, panel_full_path); - free_event(ev); - H5Fclose(fh); - return 1; - } - - dh = H5Dopen2(fh, panel_full_path, H5P_DEFAULT); - if ( dh < 0 ) { - ERROR("Cannot open data for panel %s (%s)\n", - p->name, panel_full_path); - free(panel_full_path); - free_event(ev); - H5Fclose(fh); - return 1; - } - - free(panel_full_path); - - /* Set up dataspace for file - * (determine where to read the data from) */ - dataspace = H5Dget_space(dh); - ndims = H5Sget_simple_extent_ndims(dataspace); - if ( ndims < 0 ) { - ERROR("Failed to get number of dimensions for panel %s\n", - p->name); - free_event(ev); - H5Fclose(fh); - return 1; - } - - if ( ndims != p->dim_structure->num_dims ) { - /* Dimensionality doesn't match */ - int n_nonplaceholder = 0; - for ( hsi=0; hsi<p->dim_structure->num_dims; hsi++ ) { - if ( p->dim_structure->dims[hsi] != HYSL_PLACEHOLDER ) { - n_nonplaceholder++; - } - } - - /* If the dimensions match after excluding - * placeholders, it's OK - probably a static mask - * in a multi-event file. */ - if ( ndims == n_nonplaceholder ) { - skip_placeholders = 1; - } - } - - f_offset = malloc(ndims*sizeof(hsize_t)); - f_count = malloc(ndims*sizeof(hsize_t)); - if ( (f_offset == NULL) || (f_count == NULL ) ) { - ERROR("Failed to allocate offset or count.\n"); - free_event(ev); - H5Fclose(fh); - return 1; - } - - dpos = 0; - for ( hsi=0; hsi<p->dim_structure->num_dims; hsi++ ) { - - switch ( p->dim_structure->dims[hsi] ) { - - case HYSL_FS: - f_offset[dpos] = p->orig_min_fs; - f_count[dpos] = p->orig_max_fs - p->orig_min_fs+1; - dpos++; - break; - - case HYSL_SS: - f_offset[dpos] = p->orig_min_ss; - f_count[dpos] = p->orig_max_ss - p->orig_min_ss+1; - dpos++; - break; - - case HYSL_PLACEHOLDER: - if ( !skip_placeholders ) { - f_offset[dpos] = ev->dim_entries[0]; - f_count[dpos] = 1; - dpos++; - } - break; - - default: - f_offset[dpos] = p->dim_structure->dims[hsi]; - f_count[dpos] = 1; - dpos++; - break; - - } - - } - - check = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, - f_offset, NULL, f_count, NULL); - if ( check < 0 ) { - ERROR("Error selecting file dataspace for panel %s\n", - p->name); - free(f_offset); - free(f_count); - free_event(ev); - H5Fclose(fh); - return 1; - } - - dims[0] = p->orig_max_ss - p->orig_min_ss + 1; - dims[1] = p->orig_max_fs - p->orig_min_fs + 1; - memspace = H5Screate_simple(2, dims, NULL); - - data = malloc(dims[0]*dims[1]*el_size); - if ( data == NULL ) { - ERROR("Failed to allocate panel %s\n", p->name); - free(f_offset); - free(f_count); - free_event(ev); - H5Fclose(fh); - return 1; - } - - r = H5Dread(dh, el_type, memspace, dataspace, - H5P_DEFAULT, data); - if ( r < 0 ) { - ERROR("Couldn't read data for panel %s\n", - p->name); - free(f_offset); - free(f_count); - free(data); - free_event(ev); - H5Fclose(fh); - return 1; - } - - H5Dclose(dh); - H5Sclose(dataspace); - free(f_offset); - free(f_count); - free_event(ev); - H5Fclose(fh); - - *pdata = data; - return 0; -} - - -static struct image *image_read_hdf5(DataTemplate *dtempl, - const char *filename, - const char *event) -{ - struct image *image; - int i; - - image = image_new(); - if ( image == NULL ) { - ERROR("Couldn't allocate image structure.\n"); - return NULL; - } - - image->dp = malloc(dtempl->n_panels*sizeof(float *)); - if ( image->dp == NULL ) { - ERROR("Failed to allocate data array.\n"); - image_free(image); - return NULL; - } - - /* Set all pointers to NULL for easier clean-up */ - for ( i=0; i<dtempl->n_panels; i++ ) image->dp[i] = NULL; - - for ( i=0; i<dtempl->n_panels; i++ ) { - if ( load_hdf5_hyperslab(&dtempl->panels[i], filename, - event, (void *)&image->dp[i], - H5T_NATIVE_FLOAT, - sizeof(float), 0, - dtempl->panels[i].data) ) - { - ERROR("Failed to load panel data\n"); - image_free(image); - return NULL; - } - } - - image->filename = strdup(filename); - image->ev = safe_strdup(event); - - return image; -} - - -static int load_mask_cbf(struct panel_template *p, - const char *filename, const char *event, - int gz, int *bad, int mask_good, int mask_bad) -{ - ERROR("Mask loading from CBF not yet supported\n"); - return 1; -} - - -/* Load bad pixels for this panel from given filename/event, and merge - * with (already allocated/initialised) mask "bad" */ -static int load_mask_hdf5(struct panel_template *p, - const char *filename, const char *event, - int *bad, int mask_good, int mask_bad) -{ - int p_w, p_h; - int *mask = NULL; - long unsigned int j; - - p_w = p->orig_max_fs - p->orig_min_fs + 1; - p_h = p->orig_max_ss - p->orig_min_ss + 1; - - if ( load_hdf5_hyperslab(p, filename, event, - (void *)&mask, H5T_NATIVE_INT, - sizeof(int), 1, p->mask) ) - { - ERROR("Failed to load mask data\n"); - free(mask); - return 1; - } - - for ( j=0; j<p_w*p_h; j++ ) { - - /* Bad if it's missing any of the "good" bits */ - if ( (mask[j] & mask_good) != mask_good ) bad[j] = 1; - - /* Bad if it has any of the "bad" bits. */ - if ( mask[j] & mask_bad ) bad[j] = 1; - - } - - free(mask); - return 0; -} - - -struct image *image_read_cbf(DataTemplate *dtempl, const char *filename, - const char *event, int gz) -{ - struct image *image; - float *data; - int w, h; - - if ( access(filename, R_OK) == -1 ) { - ERROR("File does not exist or cannot be read: %s\n", filename); - return NULL; - } - - image = image_new(); - if ( image == NULL ) { - ERROR("Couldn't allocate image structure.\n"); - return NULL; - } - - data = read_cbf_data(filename, gz, &w, &h); - if ( data == NULL ) { - ERROR("Failed to read CBF data\n"); - return NULL; - } - - unpack_panels_dtempl(image, dtempl, data, w, h); - free(data); - - //cbf_fill_in_beam_parameters(image->beam, f, image); - //cbf_fill_in_clen(image->det, f); - //fill_in_adu(image); - - return image; -} - - -static double get_value_hdf5(const char *name, const char *filename, - const char *event) -{ - hid_t dh; - hid_t type; - hid_t class; - hid_t sh; - hid_t ms; - hsize_t *f_offset = NULL; - hsize_t *f_count = NULL; - hsize_t m_offset[1]; - hsize_t m_count[1]; - hsize_t msdims[1]; - hsize_t size[64]; - herr_t r; - herr_t check; - int check_pe; - int dim_flag; - int ndims; - int i; - char *subst_name = NULL; - struct event *ev; - hid_t fh; - double val; - - if ( access(filename, R_OK) == -1 ) { - ERROR("File does not exist or cannot be read: %s\n", filename); - return NAN; - } - - fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); - if ( fh < 0 ) { - ERROR("Couldn't open file: %s\n", filename); - return NAN; - } - ev = get_event_from_event_string(event); - if ( (ev == NULL) && (event != NULL) ) { - ERROR("Invalid event identifier '%s'\n", event); - H5Fclose(fh); - return NAN; - } - - subst_name = retrieve_full_path(ev, name); - - check_pe = check_path_existence(fh, subst_name); - if ( check_pe == 0 ) { - ERROR("No such event-based numeric field '%s'\n", subst_name); - return NAN; - } - - dh = H5Dopen2(fh, subst_name, H5P_DEFAULT); - type = H5Dget_type(dh); - class = H5Tget_class(type); - - if ( (class != H5T_FLOAT) && (class != H5T_INTEGER) ) { - ERROR("Not a floating point or integer value.\n"); - H5Tclose(type); - H5Dclose(dh); - return NAN; - } - - /* Get the dimensionality. We have to cope with scalars expressed as - * arrays with all dimensions 1, as well as zero-d arrays. */ - sh = H5Dget_space(dh); - ndims = H5Sget_simple_extent_ndims(sh); - if ( ndims > 64 ) { - ERROR("Too many dimensions for numeric value\n"); - H5Tclose(type); - H5Dclose(dh); - return NAN; - } - H5Sget_simple_extent_dims(sh, size, NULL); - - m_offset[0] = 0; - m_count[0] = 1; - msdims[0] = 1; - ms = H5Screate_simple(1,msdims,NULL); - - /* Check that the size in all dimensions is 1 - * or that one of the dimensions has the same - * size as the hyperplane events */ - - dim_flag = 0; - - for ( i=0; i<ndims; i++ ) { - if ( size[i] == 1 ) continue; - if ( ( i==0 ) && (ev != NULL) && (size[i] > ev->dim_entries[0]) ) { - dim_flag = 1; - } else { - H5Tclose(type); - H5Dclose(dh); - return NAN; - } - } - - if ( dim_flag == 0 ) { - - if ( H5Dread(dh, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, - H5P_DEFAULT, &val) < 0 ) - { - ERROR("Couldn't read value.\n"); - H5Tclose(type); - H5Dclose(dh); - return NAN; - } - - } else { - - f_offset = malloc(ndims*sizeof(hsize_t)); - f_count = malloc(ndims*sizeof(hsize_t)); - - for ( i=0; i<ndims; i++ ) { - - if ( i == 0 ) { - f_offset[i] = ev->dim_entries[0]; - f_count[i] = 1; - } else { - f_offset[i] = 0; - f_count[i] = 0; - } - - } - - check = H5Sselect_hyperslab(sh, H5S_SELECT_SET, - f_offset, NULL, f_count, NULL); - if ( check <0 ) { - ERROR("Error selecting dataspace for float value"); - free(f_offset); - free(f_count); - return NAN; - } - - ms = H5Screate_simple(1,msdims,NULL); - check = H5Sselect_hyperslab(ms, H5S_SELECT_SET, - m_offset, NULL, m_count, NULL); - if ( check < 0 ) { - ERROR("Error selecting memory dataspace for float value"); - free(f_offset); - free(f_count); - return NAN; - } - - r = H5Dread(dh, H5T_NATIVE_DOUBLE, ms, sh, H5P_DEFAULT, &val); - if ( r < 0 ) { - ERROR("Couldn't read value.\n"); - H5Tclose(type); - H5Dclose(dh); - return NAN; - } - - } - - free_event(ev); - free(subst_name); - H5Fclose(fh); - - return val; -} - - static double get_value(struct image *image, const char *from) { double val; @@ -1353,7 +289,9 @@ static double get_value(struct image *image, const char *from) if ( *rval == '\0' ) return val; if ( H5Fis_hdf5(image->filename) > 0 ) { - return get_value_hdf5(from, image->filename, image->ev); + return image_hdf5_get_value(from, + image->filename, + image->ev); } else if ( is_cbf_file(image->filename) > 0 ) { /* FIXME: From headers */ @@ -1489,13 +427,13 @@ struct image *image_read(DataTemplate *dtempl, const char *filename, } if ( H5Fis_hdf5(filename) > 0 ) { - image = image_read_hdf5(dtempl, filename, event); + image = image_hdf5_read(dtempl, filename, event); } else if ( is_cbf_file(filename) > 0 ) { - image = image_read_cbf(dtempl, filename, event, 0); + image = image_cbf_read(dtempl, filename, event, 0); } else if ( is_cbfgz_file(filename) ) { - image = image_read_cbf(dtempl, filename, event, 1); + image = image_cbf_read(dtempl, filename, event, 1); } else { ERROR("Unrecognised file type: %s\n", filename); @@ -1557,22 +495,22 @@ struct image *image_read(DataTemplate *dtempl, const char *filename, mask_fn = p->mask_file; } if ( H5Fis_hdf5(mask_fn) > 0 ) { - load_mask_hdf5(p, mask_fn, event, - image->bad[i], - dtempl->mask_good, - dtempl->mask_bad); + image_hdf5_read_mask(p, mask_fn, event, + image->bad[i], + dtempl->mask_good, + dtempl->mask_bad); } else if ( is_cbf_file(filename) > 0 ) { - load_mask_cbf(p, mask_fn, event, - 0, image->bad[i], - dtempl->mask_good, - dtempl->mask_bad); + image_cbf_read_mask(p, mask_fn, event, + 0, image->bad[i], + dtempl->mask_good, + dtempl->mask_bad); } else if ( is_cbfgz_file(filename) ) { - load_mask_cbf(p, mask_fn, event, - 1, image->bad[i], - dtempl->mask_good, - dtempl->mask_bad); + image_cbf_read_mask(p, mask_fn, event, + 1, image->bad[i], + dtempl->mask_good, + dtempl->mask_bad); } else { ERROR("Unrecognised mask file type" @@ -1621,441 +559,39 @@ void image_free(struct image *image) } -static int read_peak_count(hid_t fh, char *path, int line, - int *num_peaks) -{ - - hid_t dh, sh, mh; - hsize_t size[1]; - hsize_t max_size[1]; - hsize_t offset[1], count[1]; - hsize_t m_offset[1], m_count[1], dimmh[1]; - int tw, r; - - dh = H5Dopen2(fh, path, H5P_DEFAULT); - if ( dh < 0 ) { - ERROR("Data block %s not found.\n", path); - return 1; - } - - sh = H5Dget_space(dh); - if ( sh < 0 ) { - H5Dclose(dh); - ERROR("Couldn't get dataspace for data.\n"); - return 1; - } - - if ( H5Sget_simple_extent_ndims(sh) != 1 ) { - ERROR("Data block %s has the wrong dimensionality (%i).\n", - path, H5Sget_simple_extent_ndims(sh)); - H5Sclose(sh); - H5Dclose(dh); - return 1; - } - - H5Sget_simple_extent_dims(sh, size, max_size); - - tw = size[0]; - - if ( line > tw-1 ) { - H5Sclose(sh); - H5Dclose(dh); - ERROR("Data block %s does not contain data for required event.\n", - path); - return 1; - } - - offset[0] = line; - count[0] = 1; - - r = H5Sselect_hyperslab(sh, H5S_SELECT_SET, - offset, NULL, count, NULL); - if ( r < 0 ) { - ERROR("Error selecting file dataspace " - "for data block %s\n", path); - H5Dclose(dh); - H5Sclose(sh); - return 1; - } - - m_offset[0] = 0; - m_count[0] = 1; - dimmh[0] = 1; - mh = H5Screate_simple(1, dimmh, NULL); - r = H5Sselect_hyperslab(mh, H5S_SELECT_SET, - m_offset, NULL, m_count, NULL); - if ( r < 0 ) { - ERROR("Error selecting memory dataspace " - "for data block %s\n", path); - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return 1; - } - - r = H5Dread(dh, H5T_NATIVE_INT, mh, - sh, H5P_DEFAULT, num_peaks); - if ( r < 0 ) { - ERROR("Couldn't read data for block %s, line %i\n", path, line); - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return 1; - } - - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return 0; -} - - -static float *read_hdf5_data(hid_t fh, char *path, int line) -{ - - hid_t dh, sh, mh; - hsize_t size[2]; - hsize_t max_size[2]; - hsize_t offset[2], count[2]; - hsize_t m_offset[2], m_count[2], dimmh[2]; - float *buf; - int tw, r; - - dh = H5Dopen2(fh, path, H5P_DEFAULT); - if ( dh < 0 ) { - ERROR("Data block (%s) not found.\n", path); - return NULL; - } - - sh = H5Dget_space(dh); - if ( sh < 0 ) { - H5Dclose(dh); - ERROR("Couldn't get dataspace for data.\n"); - return NULL; - } - - if ( H5Sget_simple_extent_ndims(sh) != 2 ) { - ERROR("Data block %s has the wrong dimensionality (%i).\n", - path, H5Sget_simple_extent_ndims(sh)); - H5Sclose(sh); - H5Dclose(dh); - return NULL; - } - - H5Sget_simple_extent_dims(sh, size, max_size); - - tw = size[0]; - if ( line> tw-1 ) { - H5Sclose(sh); - H5Dclose(dh); - ERROR("Data block %s does not contain data for required event.\n", - path); - return NULL; - } - - offset[0] = line; - offset[1] = 0; - count[0] = 1; - count[1] = size[1]; - - r = H5Sselect_hyperslab(sh, H5S_SELECT_SET, offset, NULL, count, NULL); - if ( r < 0 ) { - ERROR("Error selecting file dataspace " - "for data block %s\n", path); - H5Dclose(dh); - H5Sclose(sh); - return NULL; - } - - m_offset[0] = 0; - m_offset[1] = 0; - m_count[0] = 1; - m_count[1] = size[1]; - dimmh[0] = 1; - dimmh[1] = size[1]; - - mh = H5Screate_simple(2, dimmh, NULL); - r = H5Sselect_hyperslab(mh, H5S_SELECT_SET, - m_offset, NULL, m_count, NULL); - if ( r < 0 ) { - ERROR("Error selecting memory dataspace " - "for data block %s\n", path); - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return NULL; - } - - buf = malloc(size[1]*sizeof(float)); - if ( buf == NULL ) return NULL; - r = H5Dread(dh, H5T_NATIVE_FLOAT, mh, sh, H5P_DEFAULT, buf); - if ( r < 0 ) { - ERROR("Couldn't read data for block %s, line %i\n", path, line); - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return NULL; - } - - H5Dclose(dh); - H5Sclose(sh); - H5Sclose(mh); - return buf; -} - - -ImageFeatureList *get_peaks_cxi_dtempl(const DataTemplate *dtempl, - const char *filename, - const char *event, - int half_pixel_shift) -{ - ImageFeatureList *features; - hid_t fh; - char path_n[1024]; - char path_x[1024]; - char path_y[1024]; - char path_i[1024]; - int r; - int pk; - struct event *ev; - char *subst_name; - - int line = 0; - int num_peaks; - - float *buf_x; - float *buf_y; - float *buf_i; - - double peak_offset = half_pixel_shift ? 0.5 : 0.0; - - if ( access(filename, R_OK) == -1 ) { - ERROR("File does not exist or cannot be read: %s\n", - filename); - return NULL; - } - - fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); - if ( fh < 0 ) { - ERROR("Couldn't open file: %s\n", filename); - return NULL; - } - - ev = get_event_from_event_string(event); - if ( (ev == NULL) && (event != NULL) ) { - ERROR("Invalid event identifier '%s'\n", event); - H5Fclose(fh); - return NULL; - } - - if ( ev->dim_entries == NULL ) { - ERROR("CXI format peak list format selected," - "but file has no event structure"); - return NULL; - } - line = ev->dim_entries[0]; - - subst_name = retrieve_full_path(ev, dtempl->peak_list); - free_event(ev); - if ( subst_name == NULL ) { - ERROR("Invalid peak path %s\n", subst_name); - H5Fclose(fh); - return NULL; - } - - if ( check_path_existence(fh, subst_name) == 0 ) { - ERROR("Peak path not found: %s:%s", - filename, subst_name); - free(subst_name); - H5Fclose(fh); - return NULL; - } - - snprintf(path_n, 1024, "%s/nPeaks", subst_name); - snprintf(path_x, 1024, "%s/peakXPosRaw", subst_name); - snprintf(path_y, 1024, "%s/peakYPosRaw", subst_name); - snprintf(path_i, 1024, "%s/peakTotalIntensity", subst_name); - - r = read_peak_count(fh, path_n, line, &num_peaks); - if ( r != 0 ) return NULL; - - buf_x = read_hdf5_data(fh, path_x, line); - if ( r != 0 ) return NULL; - - buf_y = read_hdf5_data(fh, path_y, line); - if ( r != 0 ) return NULL; - - buf_i = read_hdf5_data(fh, path_i, line); - if ( r != 0 ) return NULL; - - features = image_feature_list_new(); - - for ( pk=0; pk<num_peaks; pk++ ) { - - float fs, ss, val; - int pn; - - fs = buf_x[pk] + peak_offset; - ss = buf_y[pk] + peak_offset; - val = buf_i[pk]; - - if ( data_template_file_to_panel_coords(dtempl, - &fs, &ss, - &pn) ) - { - ERROR("Failed to convert %i,%i to " - "panel-relative coordinates\n", fs, ss); - } else { - image_add_feature(features, fs, ss, pn, - NULL, val, NULL); - } - - } - - return features; -} - - -ImageFeatureList *get_peaks_hdf5_dtempl(const DataTemplate *dtempl, - const char *filename, - const char *event, - int half_pixel_shift) +struct image *image_new() { - hid_t fh, dh, sh; - hsize_t size[2]; - hsize_t max_size[2]; - int i; - float *buf; - herr_t r; - int tw; - struct event *ev; - char *subst_name; - ImageFeatureList *features; - double peak_offset = half_pixel_shift ? 0.5 : 0.0; - - if ( access(filename, R_OK) == -1 ) { - ERROR("File does not exist or cannot be read: %s\n", - filename); - return NULL; - } - - fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); - if ( fh < 0 ) { - ERROR("Couldn't open file: %s\n", filename); - return NULL; - } - - ev = get_event_from_event_string(event); - if ( (ev == NULL) && (event != NULL) ) { - ERROR("Invalid event identifier '%s'\n", event); - H5Fclose(fh); - return NULL; - } - - subst_name = retrieve_full_path(ev, dtempl->peak_list); - free_event(ev); - if ( subst_name == NULL ) { - ERROR("Invalid peak path %s\n", subst_name); - free_event(ev); - H5Fclose(fh); - return NULL; - } - - if ( check_path_existence(fh, subst_name) == 0 ) { - ERROR("Peak path not found: %s:%s", - filename, subst_name); - free(subst_name); - H5Fclose(fh); - return NULL; - } - - dh = H5Dopen2(fh, subst_name, H5P_DEFAULT); - free(subst_name); - if ( dh < 0 ) { - ERROR("Peak list (%s) not found.\n", subst_name); - H5Fclose(fh); - return NULL; - } - - sh = H5Dget_space(dh); - if ( sh < 0 ) { - ERROR("Couldn't get dataspace for peak list.\n"); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - - if ( H5Sget_simple_extent_ndims(sh) != 2 ) { - ERROR("Peak list has the wrong dimensionality (%i).\n", - H5Sget_simple_extent_ndims(sh)); - H5Sclose(sh); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - - H5Sget_simple_extent_dims(sh, size, max_size); - H5Sclose(sh); - - tw = size[1]; - if ( (tw != 3) && (tw != 4) ) { - ERROR("Peak list has the wrong dimensions.\n"); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - - buf = malloc(sizeof(float)*size[0]*size[1]); - if ( buf == NULL ) { - ERROR("Couldn't reserve memory for the peak list.\n"); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - r = H5Dread(dh, H5T_NATIVE_FLOAT, H5S_ALL, H5S_ALL, - H5P_DEFAULT, buf); - if ( r < 0 ) { - ERROR("Couldn't read peak list.\n"); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - - features = image_feature_list_new(); - if ( features == NULL ) { - ERROR("Failed to allocate peak list\n"); - H5Dclose(dh); - H5Fclose(fh); - return NULL; - } - - for ( i=0; i<size[0]; i++ ) { + struct image *image; - float fs, ss, val; - int pn; + image = malloc(sizeof(struct image)); + if ( image == NULL ) return NULL; - fs = buf[tw*i+0] + peak_offset; - ss = buf[tw*i+1] + peak_offset; - val = buf[tw*i+2]; + image->dp = NULL; + image->bad = NULL; + image->sat = NULL; + image->hit = 0; + image->crystals = NULL; + image->n_crystals = 0; + image->indexed_by = INDEXING_NONE; + image->detgeom = NULL; + image->filename = NULL; + image->ev = NULL; + image->copied_headers = NULL; + image->avg_clen = -1.0; + image->id = 0; + image->serial = 0; + image->spectrum = NULL; + image->lambda = -1.0; + image->div = -1.0; + image->bw = -1.0; + image->peak_resolution = -1.0; + image->features = NULL; - if ( data_template_file_to_panel_coords(dtempl, - &fs, &ss, - &pn) ) - { - ERROR("Failed to convert %i,%i to " - "panel-relative coordinates\n", fs, ss); - } else { - image_add_feature(features, fs, ss, pn, - NULL, val, NULL); - } - - } - - free(buf); - H5Dclose(dh); - H5Fclose(fh); + /* Deprecated stuff */ + image->beam = NULL; + image->det = NULL; - return features; + return image; } @@ -2069,16 +605,16 @@ ImageFeatureList *image_read_peaks(const DataTemplate *dtempl, const char *ext; ext = filename_extension(filename, NULL); if ( strcmp(ext, ".cxi") == 0 ) { - return get_peaks_cxi_dtempl(dtempl, - filename, - event, - half_pixel_shift); + return image_hdf5_read_peaks_cxi(dtempl, + filename, + event, + half_pixel_shift); } else { - return get_peaks_hdf5_dtempl(dtempl, - filename, - event, - half_pixel_shift); + return image_hdf5_read_peaks_hdf5(dtempl, + filename, + event, + half_pixel_shift); } } else { @@ -2088,318 +624,8 @@ ImageFeatureList *image_read_peaks(const DataTemplate *dtempl, } -struct parse_params { - hid_t fh; - int path_dim; - const char *path; - struct event *curr_event; - struct event_list *ev_list; - int top_level; -}; - - -static herr_t parse_file_event_structure(hid_t loc_id, char *name, - const H5L_info_t *info, - struct parse_params *pp) - -{ - char *substituted_path; - char *ph_loc; - char *truncated_path; - htri_t check; - herr_t herrt_iterate, herrt_info; - struct H5O_info_t object_info; - - if ( !pp->top_level ) { - - int fail_push; - - fail_push = push_path_entry_to_event(pp->curr_event, name); - if ( fail_push ) { - return -1; - } - - substituted_path = event_path_placeholder_subst(name, pp->path); - - } else { - substituted_path = strdup(pp->path); - } - - if ( pp->top_level == 1 ) { - pp->top_level = 0; - } - - truncated_path = strdup(substituted_path); - ph_loc = strstr(substituted_path,"%"); - if ( ph_loc != NULL) { - truncated_path[ph_loc-substituted_path] = '\0'; - } - - herrt_iterate = 0; - herrt_info = 0; - - check = check_path_existence(pp->fh, truncated_path); - if ( check == 0 ) { - pop_path_entry_from_event(pp->curr_event); - return 0; - } else { - - herrt_info = H5Oget_info_by_name(pp->fh, truncated_path, - &object_info, H5P_DEFAULT); - if ( herrt_info < 0 ) { - free(truncated_path); - free(substituted_path); - return -1; - } - - if ( pp->curr_event->path_length == pp->path_dim - && object_info.type == H5O_TYPE_DATASET ) - { - - int fail_append; - - fail_append = append_event_to_event_list(pp->ev_list, - pp->curr_event); - if ( fail_append ) { - free(truncated_path); - free(substituted_path); - return -1; - } - - pop_path_entry_from_event(pp->curr_event); - return 0; - - } else { - - pp->path = substituted_path; - - if ( object_info.type == H5O_TYPE_GROUP ) { - - herrt_iterate = H5Literate_by_name(pp->fh, - truncated_path, H5_INDEX_NAME, - H5_ITER_NATIVE, NULL, - (H5L_iterate_t)parse_file_event_structure, - (void *)pp, H5P_DEFAULT); - } - } - } - - pop_path_entry_from_event(pp->curr_event); - - free(truncated_path); - free(substituted_path); - - return herrt_iterate; -} - - -static int fill_paths(hid_t fh, const DataTemplate *dtempl, int pi, - struct event_list *master_el) -{ - struct parse_params pparams; - struct event *empty_event; - struct event_list *panel_ev_list; - int ei; - int check; - - empty_event = initialize_event(); - panel_ev_list = initialize_event_list(); - if ( (empty_event == NULL) || (panel_ev_list == NULL) ) - { - ERROR("Failed to allocate memory for event list.\n"); - return 1; - } - - pparams.path = dtempl->panels[pi].data; - pparams.fh = fh; - pparams.path_dim = dtempl->path_dim; - pparams.curr_event = empty_event; - pparams.top_level = 1; - pparams.ev_list = panel_ev_list; - - check = parse_file_event_structure(fh, NULL, NULL, &pparams); - if ( check < 0 ) { - free_event(empty_event); - free_event_list(panel_ev_list); - return 1; - } - - for ( ei=0; ei<panel_ev_list->num_events; ei++ ) { - - int fail_add; - - fail_add = add_non_existing_event_to_event_list(master_el, - panel_ev_list->events[ei]); - if ( fail_add ) { - free_event(empty_event); - free_event_list(panel_ev_list); - return 1; - } - - } - - free_event(empty_event); - free_event_list(panel_ev_list); - - return 0; -} - - -static int check_dims(hid_t fh, struct panel_template *p, - struct event *ev, - struct event_list *events, int *global_path_dim) -{ - char *full_panel_path; - hid_t dh; - hid_t sh; - int dims; - hsize_t *size; - hsize_t *max_size; - int hsdi; - int panel_path_dim = 0; - struct dim_structure *panel_dim_structure; - - /* Get the full path for this panel in this event */ - full_panel_path = retrieve_full_path(ev, p->data); - - dh = H5Dopen2(fh, full_panel_path, H5P_DEFAULT); - if ( dh < 0 ) { - ERROR("Error opening '%s'\n", full_panel_path); - ERROR("Failed to enumerate events. " - "Check your geometry file.\n"); - return 1; - } - - sh = H5Dget_space(dh); - dims = H5Sget_simple_extent_ndims(sh); - size = malloc(dims*sizeof(hsize_t)); - max_size = malloc(dims*sizeof(hsize_t)); - if ( (size==NULL) || (max_size==NULL) ) { - ERROR("Failed to allocate memory for dimensions\n"); - return 1; - } - - dims = H5Sget_simple_extent_dims(sh, size, max_size); - - panel_dim_structure = p->dim_structure; - for ( hsdi=0; hsdi<panel_dim_structure->num_dims; hsdi++ ) { - if ( panel_dim_structure->dims[hsdi] == HYSL_PLACEHOLDER ) { - panel_path_dim = size[hsdi]; - break; - } - } - - if ( *global_path_dim == -1 ) { - - *global_path_dim = panel_path_dim; - - } else if ( panel_path_dim != *global_path_dim ) { - - ERROR("All panels must have the same number of frames\n"); - ERROR("Panel %s has %i frames in one dimension, but the first " - "panel has %i.\n", - p->name, panel_path_dim, *global_path_dim); - free(size); - free(max_size); - return 1; - } - - H5Sclose(sh); - H5Dclose(dh); - - return 0; -} - - struct event_list *image_expand_frames(const DataTemplate *dtempl, const char *filename) { - struct event_list *master_el; - hid_t fh; - - fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT); - if ( fh < 0 ) { - ERROR("Couldn't open file '%s'\n", filename); - return NULL; - } - - master_el = initialize_event_list(); - if ( master_el == NULL ) { - ERROR("Failed to allocate event list.\n"); - H5Fclose(fh); - return NULL; - } - - /* First expand any placeholders in the HDF5 paths */ - if ( dtempl->path_dim != 0 ) { - int pi; - for ( pi=0; pi<dtempl->n_panels; pi++ ) { - if ( fill_paths(fh, dtempl, pi, master_el) ) { - ERROR("Failed to enumerate paths.\n"); - H5Fclose(fh); - return NULL; - } - } - } - - /* Now enumerate the placeholder dimensions */ - if ( dtempl->dim_dim > 0 ) { - - struct event_list *master_el_with_dims; - int evi; - - /* If there were no HDF5 path placeholders, add a dummy event */ - if ( master_el->num_events == 0 ) { - struct event *empty_ev; - empty_ev = initialize_event(); - append_event_to_event_list(master_el, empty_ev); - free(empty_ev); - } - - master_el_with_dims = initialize_event_list(); - - /* For each event so far, expand the dimensions */ - for ( evi=0; evi<master_el->num_events; evi++ ) { - - int pi; - int global_path_dim = -1; - int mlwd; - - /* Check the dimensionality of each panel */ - for ( pi=0; pi<dtempl->n_panels; pi++ ) { - if ( check_dims(fh, &dtempl->panels[pi], - master_el->events[evi], - master_el_with_dims, - &global_path_dim) ) - { - ERROR("Failed to enumerate dims.\n"); - H5Fclose(fh); - return NULL; - } - } - - /* Add this dimensionality to all events */ - for ( mlwd=0; mlwd<global_path_dim; mlwd++ ) { - - struct event *mlwd_ev; - - mlwd_ev = copy_event(master_el->events[evi]); - push_dim_entry_to_event(mlwd_ev, mlwd); - append_event_to_event_list(master_el_with_dims, - mlwd_ev); - free_event(mlwd_ev); - } - - } - - free_event_list(master_el); - H5Fclose(fh); - return master_el_with_dims; - - } else { - - H5Fclose(fh); - return master_el; - - } + return image_hdf5_expand_frames(dtempl, filename); } |