diff options
| -rw-r--r-- | Makefile.am | 8 | ||||
| -rw-r--r-- | Makefile.in | 24 | ||||
| -rw-r--r-- | README | 7 | ||||
| -rw-r--r-- | src/hdf5-file.c | 18 | ||||
| -rw-r--r-- | src/image.h | 4 | ||||
| -rw-r--r-- | src/indexamajig.c | 37 | ||||
| -rw-r--r-- | src/pattern_sim.c | 4 | ||||
| -rw-r--r-- | src/reflist.c | 1 | ||||
| -rw-r--r-- | src/stream.c | 152 | ||||
| -rw-r--r-- | src/stream.h | 14 |
10 files changed, 197 insertions, 72 deletions
diff --git a/Makefile.am b/Makefile.am index 5a309f76..ce2d06d7 100644 --- a/Makefile.am +++ b/Makefile.am @@ -45,7 +45,7 @@ src_process_hkl_SOURCES = src/process_hkl.c src/sfac.c src/statistics.c \ src/cell.c src/utils.c src/reflections.c \ src/symmetry.c src/stream.c src/beam-parameters.c \ src/thread-pool.c src/image.c src/detector.c \ - src/hdf5-file.c + src/hdf5-file.c src/reflist.c src_indexamajig_SOURCES = src/indexamajig.c src/hdf5-file.c src/utils.c \ src/cell.c src/image.c src/peaks.c src/index.c \ @@ -104,11 +104,13 @@ src_partialator_SOURCES = src/partialator.c src/cell.c src/hdf5-file.c \ if BUILD_CUBEIT src_cubeit_SOURCES = src/cubeit.c src/cell.c src/hdf5-file.c src/utils.c \ src/detector.c src/render.c src/filters.c src/image.c \ - src/symmetry.c src/stream.c src/thread-pool.c + src/symmetry.c src/stream.c src/thread-pool.c src/reflist.c endif src_estimate_background_SOURCES = src/estimate_background.c src/stream.c \ - src/utils.c src/cell.c src/thread-pool.c + src/utils.c src/cell.c src/thread-pool.c \ + src/image.c src/detector.c src/hdf5-file.c \ + src/reflist.c tests_list_check_SOURCES = tests/list_check.c src/reflist.c src/thread-pool.c \ src/utils.c diff --git a/Makefile.in b/Makefile.in index b76f0796..b65b60a3 100644 --- a/Makefile.in +++ b/Makefile.in @@ -106,21 +106,25 @@ src_compare_hkl_LDADD = $(LDADD) src_compare_hkl_DEPENDENCIES = $(top_builddir)/lib/libgnu.a am__src_cubeit_SOURCES_DIST = src/cubeit.c src/cell.c src/hdf5-file.c \ src/utils.c src/detector.c src/render.c src/filters.c \ - src/image.c src/symmetry.c src/stream.c src/thread-pool.c + src/image.c src/symmetry.c src/stream.c src/thread-pool.c \ + src/reflist.c @BUILD_CUBEIT_TRUE@am_src_cubeit_OBJECTS = src/cubeit.$(OBJEXT) \ @BUILD_CUBEIT_TRUE@ src/cell.$(OBJEXT) src/hdf5-file.$(OBJEXT) \ @BUILD_CUBEIT_TRUE@ src/utils.$(OBJEXT) src/detector.$(OBJEXT) \ @BUILD_CUBEIT_TRUE@ src/render.$(OBJEXT) src/filters.$(OBJEXT) \ @BUILD_CUBEIT_TRUE@ src/image.$(OBJEXT) src/symmetry.$(OBJEXT) \ @BUILD_CUBEIT_TRUE@ src/stream.$(OBJEXT) \ -@BUILD_CUBEIT_TRUE@ src/thread-pool.$(OBJEXT) +@BUILD_CUBEIT_TRUE@ src/thread-pool.$(OBJEXT) \ +@BUILD_CUBEIT_TRUE@ src/reflist.$(OBJEXT) src_cubeit_OBJECTS = $(am_src_cubeit_OBJECTS) src_cubeit_LDADD = $(LDADD) src_cubeit_DEPENDENCIES = $(top_builddir)/lib/libgnu.a am_src_estimate_background_OBJECTS = \ src/estimate_background.$(OBJEXT) src/stream.$(OBJEXT) \ src/utils.$(OBJEXT) src/cell.$(OBJEXT) \ - src/thread-pool.$(OBJEXT) + src/thread-pool.$(OBJEXT) src/image.$(OBJEXT) \ + src/detector.$(OBJEXT) src/hdf5-file.$(OBJEXT) \ + src/reflist.$(OBJEXT) src_estimate_background_OBJECTS = \ $(am_src_estimate_background_OBJECTS) src_estimate_background_LDADD = $(LDADD) @@ -205,7 +209,7 @@ am_src_process_hkl_OBJECTS = src/process_hkl.$(OBJEXT) \ src/symmetry.$(OBJEXT) src/stream.$(OBJEXT) \ src/beam-parameters.$(OBJEXT) src/thread-pool.$(OBJEXT) \ src/image.$(OBJEXT) src/detector.$(OBJEXT) \ - src/hdf5-file.$(OBJEXT) + src/hdf5-file.$(OBJEXT) src/reflist.$(OBJEXT) src_process_hkl_OBJECTS = $(am_src_process_hkl_OBJECTS) src_process_hkl_LDADD = $(LDADD) src_process_hkl_DEPENDENCIES = $(top_builddir)/lib/libgnu.a @@ -611,7 +615,7 @@ src_process_hkl_SOURCES = src/process_hkl.c src/sfac.c src/statistics.c \ src/cell.c src/utils.c src/reflections.c \ src/symmetry.c src/stream.c src/beam-parameters.c \ src/thread-pool.c src/image.c src/detector.c \ - src/hdf5-file.c + src/hdf5-file.c src/reflist.c src_indexamajig_SOURCES = src/indexamajig.c src/hdf5-file.c \ src/utils.c src/cell.c src/image.c src/peaks.c src/index.c \ @@ -661,10 +665,12 @@ src_partialator_SOURCES = src/partialator.c src/cell.c src/hdf5-file.c \ @BUILD_CUBEIT_TRUE@src_cubeit_SOURCES = src/cubeit.c src/cell.c src/hdf5-file.c src/utils.c \ @BUILD_CUBEIT_TRUE@ src/detector.c src/render.c src/filters.c src/image.c \ -@BUILD_CUBEIT_TRUE@ src/symmetry.c src/stream.c src/thread-pool.c +@BUILD_CUBEIT_TRUE@ src/symmetry.c src/stream.c src/thread-pool.c src/reflist.c src_estimate_background_SOURCES = src/estimate_background.c src/stream.c \ - src/utils.c src/cell.c src/thread-pool.c + src/utils.c src/cell.c src/thread-pool.c \ + src/image.c src/detector.c src/hdf5-file.c \ + src/reflist.c tests_list_check_SOURCES = tests/list_check.c src/reflist.c src/thread-pool.c \ src/utils.c @@ -843,6 +849,8 @@ src/filters.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) src/stream.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) +src/reflist.$(OBJEXT): src/$(am__dirstamp) \ + src/$(DEPDIR)/$(am__dirstamp) src/cubeit$(EXEEXT): $(src_cubeit_OBJECTS) $(src_cubeit_DEPENDENCIES) src/$(am__dirstamp) @rm -f src/cubeit$(EXEEXT) $(AM_V_CCLD)$(LINK) $(src_cubeit_OBJECTS) $(src_cubeit_LDADD) $(LIBS) @@ -876,8 +884,6 @@ src/mosflm.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) src/geometry.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) -src/reflist.$(OBJEXT): src/$(am__dirstamp) \ - src/$(DEPDIR)/$(am__dirstamp) src/diffraction-gpu.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) src/cl-utils.$(OBJEXT): src/$(am__dirstamp) \ @@ -61,12 +61,7 @@ In addition, there is also: - estimate_background, for calculating signal to noise ratios from the indexed peaks. - - calibrate_detector, which splits - - - reintegrate, which is like "indexamajig" but without the indexing - step, instead getting the orientation matrix from the - output of a previous run of either "indexamajig" or - "reintegrate". + - calibrate_detector, which does nothing useful at the moment. Included at no extra cost are: diff --git a/src/hdf5-file.c b/src/hdf5-file.c index 20141522..eef85591 100644 --- a/src/hdf5-file.c +++ b/src/hdf5-file.c @@ -292,21 +292,21 @@ static double get_wavelength(struct hdfile *f) } -static double get_f0(struct hdfile *f) +static double get_i0(struct hdfile *f) { herr_t r; hid_t dh; - double f0; + double i0; dh = H5Dopen2(f->fh, "/LCLS/f_11_ENRC", H5P_DEFAULT); if ( dh < 0 ) return -1.0; r = H5Dread(dh, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, - H5P_DEFAULT, &f0); + H5P_DEFAULT, &i0); H5Dclose(dh); if ( r < 0 ) return -1.0; - return f0; + return i0; } @@ -434,13 +434,13 @@ int hdf5_read(struct hdfile *f, struct image *image, int satcorr) /* Read wavelength from file */ image->lambda = get_wavelength(f); - image->f0 = get_f0(f); - if ( image->f0 < 0.0 ) { + image->i0 = get_i0(f); + if ( image->i0 < 0.0 ) { ERROR("Couldn't read incident intensity - using 1.0.\n"); - image->f0 = 1.0; - image->f0_available = 0; + image->i0 = 1.0; + image->i0_available = 0; } else { - image->f0_available = 1; + image->i0_available = 1; } if ( satcorr ) debodge_saturation(f, image); diff --git a/src/image.h b/src/image.h index 69c36de8..692ef7f8 100644 --- a/src/image.h +++ b/src/image.h @@ -76,8 +76,8 @@ struct image { double lambda; /* Wavelength in m */ double div; /* Divergence in radians */ double bw; /* Bandwidth as a fraction */ - double f0; /* Incident intensity */ - int f0_available; /* 0 if f0 wasn't available + double i0; /* Incident intensity */ + int i0_available; /* 0 if f0 wasn't available * from the input. */ double osf; /* Overall scaling factor */ double profile_radius; /* Radius of reflection */ diff --git a/src/indexamajig.c b/src/indexamajig.c index 72491450..320439b3 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -52,7 +52,7 @@ struct static_index_args int config_cmfilter; int config_noisefilter; int config_verbose; - int stream_flags; + int stream_flags; /* What goes into the output? */ int config_polar; int config_satcorr; int config_closer; @@ -62,7 +62,7 @@ struct static_index_args struct detector *det; IndexingMethod *indm; IndexingPrivate **ipriv; - int peaks; + int peaks; /* Peak detection method */ int cellr; struct beam_params *beam; const char *element; @@ -132,15 +132,20 @@ static void show_help(const char *s) " in the HDF5 file.\n" "\n\n" "You can control what information is included in the output stream using\n" -"' --record=<flags>'. Possible flags are:\n" -"pixels Include a list of sums of pixel values within the\n" -" integration domain, correcting for individual pixel\n" -" solid angles.\n" -"integrated Include a list of reflection intensities, produced by\n" -" integrating around predicted peak locations.\n" -" The flags 'pixels' and 'integrated' are mutually exclusive.\n" -"peaks Include peak locations and intensities from the peak search.\n" -"peaksifindexed Include peak locations only if the pattern could be indexed.\n" +"' --record=<flags>'. Possible flags are:\n\n" +" pixels Include a list of sums of pixel values within the\n" +" integration domain, correcting for individual pixel\n" +" solid angles.\n" +"\n" +" integrated Include a list of reflection intensities, produced by\n" +" integrating around predicted peak locations.\n" +"\n" +" peaks Include peak locations and intensities from the peak search.\n" +"\n" +" peaksifindexed As 'peaks', but only if the pattern could be indexed.\n\n" +"\n" +"The default is '--record=integrated'. The flags 'pixels' and 'integrated'\n" +"are mutually exclusive, as are the flags 'peaks' and 'peaksifindexed'.\n" "\n\n" "For more control over the process, you might need:\n\n" " --cell-reduction=<m> Use <m> as the cell reduction method. Choose from:\n" @@ -162,7 +167,7 @@ static void show_help(const char *s) " --min-gradient=<n> Minimum gradient for Zaefferer peak search.\n" " Default: 100,000.\n" " -e, --image=<element> Use this image from the HDF5 file.\n" -" Example: /data/data0." +" Example: /data/data0.\n" " Default: The first one found.\n" "\n" "\nOptions for greater performance or verbosity:\n\n" @@ -431,6 +436,7 @@ int main(int argc, char *argv[]) struct beam_params *beam = NULL; char *element = NULL; double nominal_photon_energy; + int stream_flags = STREAM_INTEGRATED; /* Long options */ const struct option longopts[] = { @@ -460,6 +466,7 @@ int main(int argc, char *argv[]) {"insane", 0, &config_insane, 1}, {"image", 1, NULL, 'e'}, {"basename", 0, &config_basename, 1}, + {"record", 1, NULL, 5}, {0, 0, NULL, 0} }; @@ -529,6 +536,11 @@ int main(int argc, char *argv[]) element = strdup(optarg); break; + case 5 : + stream_flags = parse_stream_flags(optarg); + if ( stream_flags < 0 ) return 1; + break; + case 0 : break; @@ -673,6 +685,7 @@ int main(int argc, char *argv[]) free(pdb); /* Start by writing the entire command line to stdout */ + fprintf(ofh, "CrystFEL stream format 2.0\n"); fprintf(ofh, "Command line:"); for ( i=0; i<argc; i++ ) { fprintf(ofh, " %s", argv[i]); diff --git a/src/pattern_sim.c b/src/pattern_sim.c index 948458bc..c0318176 100644 --- a/src/pattern_sim.c +++ b/src/pattern_sim.c @@ -472,8 +472,8 @@ int main(int argc, char *argv[]) image.filename = NULL; image.features = NULL; image.flags = NULL; - image.f0 = 1.0; - image.f0_available = 1; + image.i0 = 1.0; + image.i0_available = 1; powder = calloc(image.width*image.height, sizeof(*powder)); diff --git a/src/reflist.c b/src/reflist.c index e7d072c0..3a53330c 100644 --- a/src/reflist.c +++ b/src/reflist.c @@ -44,6 +44,7 @@ struct _refldata { /* Intensity */ double intensity; + double esd_i; }; diff --git a/src/stream.c b/src/stream.c index fbe0abad..4e34a48c 100644 --- a/src/stream.c +++ b/src/stream.c @@ -21,6 +21,63 @@ #include "cell.h" #include "utils.h" #include "image.h" +#include "stream.h" + + +static void exclusive(const char *a, const char *b) +{ + ERROR("The stream options '%s' and '%s' are mutually exclusive.\n", + a, b); +} + + +int parse_stream_flags(const char *a) +{ + int n, i; + int ret = STREAM_NONE; + char **flags; + + n = assplode(a, ",", &flags, ASSPLODE_NONE); + + for ( i=0; i<n; i++ ) { + + if ( strcmp(flags[i], "pixels") == 0) { + if ( ret & STREAM_INTEGRATED ) { + exclusive("pixels", "integrated"); + return -1; + } + ret |= STREAM_PIXELS; + } else if ( strcmp(flags[i], "integrated") == 0) { + if ( ret & STREAM_PIXELS ) { + exclusive("pixels", "integrated"); + return -1; + } + ret |= STREAM_INTEGRATED; + } else if ( strcmp(flags[i], "peaks") == 0) { + if ( ret & STREAM_PEAKS_IF_INDEXED ) { + exclusive("peaks", "peaksifindexed"); + return -1; + } + ret |= STREAM_PEAKS; + } else if ( strcmp(flags[i], "peaksifindexed") == 0) { + if ( ret & STREAM_PEAKS ) { + exclusive("peaks", "peaksifindexed"); + return -1; + } + ret |= STREAM_PEAKS_IF_INDEXED; + } else { + ERROR("Unrecognised stream flag '%s'\n", flags[i]); + return 0; + } + + free(flags[i]); + + } + free(flags); + + return ret; + +} int count_patterns(FILE *fh) @@ -76,6 +133,31 @@ static UnitCell *read_orientation_matrix(FILE *fh) static void write_reflections(struct image *image, FILE *ofh) { + Reflection *refl; + RefListIterator *iter; + + fprintf(ofh, "Reflections measured after indexing\n"); + /* FIXME: Unify this with write_reflections() over in reflections.c */ + fprintf(ofh, " h k l I phase sigma(I) " + " 1/d(nm^-1) counts\n"); + + for ( refl = first_refl(image->reflections, &iter); + refl != NULL; + refl = next_refl(refl, iter) ) { + + signed int h, k, l; + double intensity, esd_i, s; + + get_indices(refl, &h, &k, &l); + intensity = get_intensity(refl); + esd_i = 0.0; /* FIXME! */ + s = 0.0; /* FIXME! */ + + /* h, k, l, I, sigma(I), s */ + fprintf(ofh, "%3i %3i %3i %10.2f %s %10.2f %10.2f %7i\n", + h, k, l, intensity, " -", esd_i, s/1.0e9, 1); + + } } @@ -83,7 +165,7 @@ static void write_peaks(struct image *image, FILE *ofh) { int i; - fprintf(ofh, "Peaks from peak search in %s\n", image->filename); + fprintf(ofh, "Peaks from peak search\n"); fprintf(ofh, " x/px y/px (1/d)/nm^-1 Intensity\n"); for ( i=0; i<image_feature_count(image->features); i++ ) { @@ -102,13 +184,10 @@ static void write_peaks(struct image *image, FILE *ofh) f->x, f->y, q/1.0e9, f->intensity); } - - fprintf(ofh, "\n"); } - -void write_chunk(FILE *ofh, struct image *image, int flags) +void write_chunk(FILE *ofh, struct image *i, int f) { double asx, asy, asz; double bsx, bsy, bsz; @@ -117,36 +196,51 @@ void write_chunk(FILE *ofh, struct image *image, int flags) fprintf(ofh, "----- Begin chunk -----\n"); - fprintf(ofh, "Image filename: %s\n", image->filename); - - cell_get_parameters(image->indexed_cell, &a, &b, &c, &al, &be, &ga); - fprintf(ofh, "Cell parameters %7.5f %7.5f %7.5f nm," - " %7.5f %7.5f %7.5f deg\n", - a*1.0e9, b*1.0e9, c*1.0e9, - rad2deg(al), rad2deg(be), rad2deg(ga)); - - cell_get_reciprocal(image->indexed_cell, &asx, &asy, &asz, - &bsx, &bsy, &bsz, - &csx, &csy, &csz); - fprintf(ofh, "astar = %+9.7f %+9.7f %+9.7f nm^-1\n", - asx/1e9, asy/1e9, asz/1e9); - fprintf(ofh, "bstar = %+9.7f %+9.7f %+9.7f nm^-1\n", - bsx/1e9, bsy/1e9, bsz/1e9); - fprintf(ofh, "cstar = %+9.7f %+9.7f %+9.7f nm^-1\n", - csx/1e9, csy/1e9, csz/1e9); - - if ( image->f0_available ) { - fprintf(ofh, "I0 = %7.5f (arbitrary gas detector units)\n", - image->f0); + fprintf(ofh, "Image filename: %s\n", i->filename); + + if ( i->indexed_cell != NULL ) { + cell_get_parameters(i->indexed_cell, &a, &b, &c, + &al, &be, &ga); + fprintf(ofh, "Cell parameters %7.5f %7.5f %7.5f nm," + " %7.5f %7.5f %7.5f deg\n", + a*1.0e9, b*1.0e9, c*1.0e9, + rad2deg(al), rad2deg(be), rad2deg(ga)); + + cell_get_reciprocal(i->indexed_cell, &asx, &asy, &asz, + &bsx, &bsy, &bsz, + &csx, &csy, &csz); + fprintf(ofh, "astar = %+9.7f %+9.7f %+9.7f nm^-1\n", + asx/1e9, asy/1e9, asz/1e9); + fprintf(ofh, "bstar = %+9.7f %+9.7f %+9.7f nm^-1\n", + bsx/1e9, bsy/1e9, bsz/1e9); + fprintf(ofh, "cstar = %+9.7f %+9.7f %+9.7f nm^-1\n", + csx/1e9, csy/1e9, csz/1e9); + + } else { + + fprintf(ofh, "No unit cell from indexing.\n"); + + } + + if ( i->i0_available ) { + fprintf(ofh, "I0 = %7.5f (arbitrary units)\n", i->i0); } else { fprintf(ofh, "I0 = invalid\n"); } fprintf(ofh, "photon_energy_eV = %f\n", - J_to_eV(ph_lambda_to_en(image->lambda))); + J_to_eV(ph_lambda_to_en(i->lambda))); + + if ( (f & STREAM_PEAKS) + || ((f & STREAM_PEAKS_IF_INDEXED) && (i->indexed_cell != NULL)) ) { + fprintf(ofh, "\n"); + write_peaks(i, ofh); + } - write_peaks(image, ofh); - write_reflections(image, ofh); + if ( (f & STREAM_PIXELS) || (f & STREAM_INTEGRATED) ) { + fprintf(ofh, "\n"); + write_reflections(i, ofh); + } fprintf(ofh, "----- End chunk -----\n\n"); } diff --git a/src/stream.h b/src/stream.h index b5f05b19..71ab4b79 100644 --- a/src/stream.h +++ b/src/stream.h @@ -19,9 +19,23 @@ struct image; +/* Possible options dictating what goes into the output stream */ +enum +{ + STREAM_NONE = 0, + STREAM_INTEGRATED = 1<<0, + STREAM_PIXELS = 1<<1, + STREAM_PEAKS = 1<<2, + STREAM_PEAKS_IF_INDEXED = 1<<3, +}; + + extern int count_patterns(FILE *fh); + extern int find_chunk(FILE *fh, UnitCell **cell, char **filename, double *ev); + extern void write_chunk(FILE *ofh, struct image *image, int flags); +extern int parse_stream_flags(const char *a); #endif /* STREAM_H */ |