diff options
| -rw-r--r-- | .gitignore | 1 | ||||
| -rw-r--r-- | Makefile.am | 10 | ||||
| -rw-r--r-- | Makefile.in | 46 | ||||
| -rw-r--r-- | src/partial_sim.c | 263 |
4 files changed, 308 insertions, 12 deletions
@@ -19,6 +19,7 @@ src/partialator src/cubeit src/check_hkl src/sum_stack +src/partial_sim src/.dirstamp *~ doc/reference/* diff --git a/Makefile.am b/Makefile.am index ba11eeed..f42e68b1 100644 --- a/Makefile.am +++ b/Makefile.am @@ -5,7 +5,7 @@ ACLOCAL_AMFLAGS = -I m4 bin_PROGRAMS = src/pattern_sim src/process_hkl src/get_hkl src/indexamajig \ src/compare_hkl src/powder_plot src/render_hkl \ src/calibrate_detector src/partialator \ - src/check_hkl src/sum_stack + src/check_hkl src/sum_stack src/partial_sim noinst_PROGRAMS = tests/list_check tests/integration_check @@ -30,6 +30,12 @@ AM_CFLAGS = -Wall AM_CPPFLAGS = -DDATADIR=\""$(datadir)"\" -I$(top_builddir)/lib -I$(top_srcdir)/lib LDADD = $(top_builddir)/lib/libgnu.a @IGNORE_UNUSED_LIBRARIES_CFLAGS@ +src_partial_sim_SOURCES = src/partial_sim.c src/cell.c src/detector.c \ + src/beam-parameters.c src/thread-pool.c src/utils.c \ + src/reflist-utils.c src/reflist.c src/symmetry.c \ + src/hdf5-file.c src/image.c src/geometry.c \ + src/peaks.c src/stream.c + src_pattern_sim_SOURCES = src/pattern_sim.c src/diffraction.c src/utils.c \ src/image.c src/cell.c src/hdf5-file.c \ src/detector.c src/peaks.c \ @@ -102,7 +108,7 @@ src_calibrate_detector_SOURCES = src/calibrate_detector.c src/utils.c \ src/hdf5-file.c src/image.c src/thread-pool.c \ src/detector.c src/stream.c src/cell.c \ src/reflist-utils.c src/reflist.c src/symmetry.c \ - src/peaks.c + src/peaks.c src_partialator_SOURCES = src/partialator.c src/cell.c src/hdf5-file.c \ src/utils.c src/detector.c src/peaks.c src/image.c \ diff --git a/Makefile.in b/Makefile.in index 1d3482bd..296e92e0 100644 --- a/Makefile.in +++ b/Makefile.in @@ -40,7 +40,8 @@ bin_PROGRAMS = src/pattern_sim$(EXEEXT) src/process_hkl$(EXEEXT) \ src/compare_hkl$(EXEEXT) src/powder_plot$(EXEEXT) \ src/render_hkl$(EXEEXT) src/calibrate_detector$(EXEEXT) \ src/partialator$(EXEEXT) src/check_hkl$(EXEEXT) \ - src/sum_stack$(EXEEXT) $(am__EXEEXT_1) $(am__EXEEXT_2) + src/sum_stack$(EXEEXT) src/partial_sim$(EXEEXT) \ + $(am__EXEEXT_1) $(am__EXEEXT_2) noinst_PROGRAMS = tests/list_check$(EXEEXT) \ tests/integration_check$(EXEEXT) $(am__EXEEXT_3) TESTS = tests/list_check$(EXEEXT) tests/first_merge_check \ @@ -166,6 +167,17 @@ am_src_indexamajig_OBJECTS = src/indexamajig.$(OBJEXT) \ src_indexamajig_OBJECTS = $(am_src_indexamajig_OBJECTS) src_indexamajig_LDADD = $(LDADD) src_indexamajig_DEPENDENCIES = $(top_builddir)/lib/libgnu.a +am_src_partial_sim_OBJECTS = src/partial_sim.$(OBJEXT) \ + src/cell.$(OBJEXT) src/detector.$(OBJEXT) \ + src/beam-parameters.$(OBJEXT) src/thread-pool.$(OBJEXT) \ + src/utils.$(OBJEXT) src/reflist-utils.$(OBJEXT) \ + src/reflist.$(OBJEXT) src/symmetry.$(OBJEXT) \ + src/hdf5-file.$(OBJEXT) src/image.$(OBJEXT) \ + src/geometry.$(OBJEXT) src/peaks.$(OBJEXT) \ + src/stream.$(OBJEXT) +src_partial_sim_OBJECTS = $(am_src_partial_sim_OBJECTS) +src_partial_sim_LDADD = $(LDADD) +src_partial_sim_DEPENDENCIES = $(top_builddir)/lib/libgnu.a am_src_partialator_OBJECTS = src/partialator.$(OBJEXT) \ src/cell.$(OBJEXT) src/hdf5-file.$(OBJEXT) src/utils.$(OBJEXT) \ src/detector.$(OBJEXT) src/peaks.$(OBJEXT) src/image.$(OBJEXT) \ @@ -283,19 +295,20 @@ am__v_GEN_0 = @echo " GEN " $@; SOURCES = $(src_calibrate_detector_SOURCES) $(src_check_hkl_SOURCES) \ $(src_compare_hkl_SOURCES) $(src_cubeit_SOURCES) \ $(src_get_hkl_SOURCES) $(src_hdfsee_SOURCES) \ - $(src_indexamajig_SOURCES) $(src_partialator_SOURCES) \ - $(src_pattern_sim_SOURCES) $(src_powder_plot_SOURCES) \ - $(src_process_hkl_SOURCES) $(src_render_hkl_SOURCES) \ - $(src_sum_stack_SOURCES) $(tests_gpu_sim_check_SOURCES) \ + $(src_indexamajig_SOURCES) $(src_partial_sim_SOURCES) \ + $(src_partialator_SOURCES) $(src_pattern_sim_SOURCES) \ + $(src_powder_plot_SOURCES) $(src_process_hkl_SOURCES) \ + $(src_render_hkl_SOURCES) $(src_sum_stack_SOURCES) \ + $(tests_gpu_sim_check_SOURCES) \ $(tests_integration_check_SOURCES) $(tests_list_check_SOURCES) DIST_SOURCES = $(src_calibrate_detector_SOURCES) \ $(src_check_hkl_SOURCES) $(src_compare_hkl_SOURCES) \ $(am__src_cubeit_SOURCES_DIST) $(src_get_hkl_SOURCES) \ $(am__src_hdfsee_SOURCES_DIST) \ - $(am__src_indexamajig_SOURCES_DIST) $(src_partialator_SOURCES) \ - $(am__src_pattern_sim_SOURCES_DIST) $(src_powder_plot_SOURCES) \ - $(src_process_hkl_SOURCES) $(src_render_hkl_SOURCES) \ - $(src_sum_stack_SOURCES) \ + $(am__src_indexamajig_SOURCES_DIST) $(src_partial_sim_SOURCES) \ + $(src_partialator_SOURCES) $(am__src_pattern_sim_SOURCES_DIST) \ + $(src_powder_plot_SOURCES) $(src_process_hkl_SOURCES) \ + $(src_render_hkl_SOURCES) $(src_sum_stack_SOURCES) \ $(am__tests_gpu_sim_check_SOURCES_DIST) \ $(tests_integration_check_SOURCES) $(tests_list_check_SOURCES) RECURSIVE_TARGETS = all-recursive check-recursive dvi-recursive \ @@ -608,6 +621,12 @@ ACLOCAL_AMFLAGS = -I m4 AM_CFLAGS = -Wall AM_CPPFLAGS = -DDATADIR=\""$(datadir)"\" -I$(top_builddir)/lib -I$(top_srcdir)/lib LDADD = $(top_builddir)/lib/libgnu.a @IGNORE_UNUSED_LIBRARIES_CFLAGS@ +src_partial_sim_SOURCES = src/partial_sim.c src/cell.c src/detector.c \ + src/beam-parameters.c src/thread-pool.c src/utils.c \ + src/reflist-utils.c src/reflist.c src/symmetry.c \ + src/hdf5-file.c src/image.c src/geometry.c \ + src/peaks.c src/stream.c + src_pattern_sim_SOURCES = src/pattern_sim.c src/diffraction.c \ src/utils.c src/image.c src/cell.c src/hdf5-file.c \ src/detector.c src/peaks.c src/reflist-utils.c \ @@ -664,7 +683,7 @@ src_calibrate_detector_SOURCES = src/calibrate_detector.c src/utils.c \ src/hdf5-file.c src/image.c src/thread-pool.c \ src/detector.c src/stream.c src/cell.c \ src/reflist-utils.c src/reflist.c src/symmetry.c \ - src/peaks.c + src/peaks.c src_partialator_SOURCES = src/partialator.c src/cell.c src/hdf5-file.c \ src/utils.c src/detector.c src/peaks.c src/image.c \ @@ -889,6 +908,11 @@ src/cl-utils.$(OBJEXT): src/$(am__dirstamp) \ src/indexamajig$(EXEEXT): $(src_indexamajig_OBJECTS) $(src_indexamajig_DEPENDENCIES) src/$(am__dirstamp) @rm -f src/indexamajig$(EXEEXT) $(AM_V_CCLD)$(LINK) $(src_indexamajig_OBJECTS) $(src_indexamajig_LDADD) $(LIBS) +src/partial_sim.$(OBJEXT): src/$(am__dirstamp) \ + src/$(DEPDIR)/$(am__dirstamp) +src/partial_sim$(EXEEXT): $(src_partial_sim_OBJECTS) $(src_partial_sim_DEPENDENCIES) src/$(am__dirstamp) + @rm -f src/partial_sim$(EXEEXT) + $(AM_V_CCLD)$(LINK) $(src_partial_sim_OBJECTS) $(src_partial_sim_LDADD) $(LIBS) src/partialator.$(OBJEXT): src/$(am__dirstamp) \ src/$(DEPDIR)/$(am__dirstamp) src/post-refinement.$(OBJEXT): src/$(am__dirstamp) \ @@ -971,6 +995,7 @@ mostlyclean-compile: -rm -f src/index.$(OBJEXT) -rm -f src/indexamajig.$(OBJEXT) -rm -f src/mosflm.$(OBJEXT) + -rm -f src/partial_sim.$(OBJEXT) -rm -f src/partialator.$(OBJEXT) -rm -f src/pattern_sim.$(OBJEXT) -rm -f src/peaks.$(OBJEXT) @@ -1017,6 +1042,7 @@ distclean-compile: @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/index.Po@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/indexamajig.Po@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/mosflm.Po@am__quote@ +@AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/partial_sim.Po@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/partialator.Po@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/pattern_sim.Po@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@src/$(DEPDIR)/peaks.Po@am__quote@ diff --git a/src/partial_sim.c b/src/partial_sim.c new file mode 100644 index 00000000..a023737f --- /dev/null +++ b/src/partial_sim.c @@ -0,0 +1,263 @@ +/* + * partial_sim.c + * + * Generate partials for testing scaling + * + * (c) 2006-2011 Thomas White <taw@physics.org> + * + * Part of CrystFEL - crystallography with a FEL + * + */ + + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <stdarg.h> +#include <stdlib.h> +#include <stdio.h> +#include <string.h> +#include <unistd.h> +#include <getopt.h> +#include <assert.h> + +#include "utils.h" +#include "reflist-utils.h" +#include "symmetry.h" +#include "beam-parameters.h" +#include "detector.h" +#include "geometry.h" +#include "stream.h" + + +/* For each reflection in "partial", fill in what the intensity would be + * according to "full" */ +static void calculate_partials(RefList *partial, RefList *full, const char *sym) +{ + Reflection *refl; + RefListIterator *iter; + + for ( refl = first_refl(partial, &iter); + refl != NULL; + refl = next_refl(refl, iter) ) { + + signed int h, k, l; + Reflection *rfull; + double p; + double Ip; + + get_indices(refl, &h, &k, &l); + get_asymm(h, k, l, &h, &k, &l, sym); + p = get_partiality(refl); + + rfull = find_refl(full, h, k, l); + if ( rfull == NULL ) { + set_redundancy(refl, 0); + } else { + Ip = p * get_intensity(rfull); + set_int(refl, Ip); + } + + } +} + + +static void show_help(const char *s) +{ + printf("Syntax: %s [options]\n\n", s); + printf( +"Generate a stream containing partials from a reflection list.\n" +"\n" +" -h, --help Display this help message.\n" +"\n" +"You need to provide the following basic options:\n" +" -i, --input=<file> Read reflections from <file>.\n" +" -o, --output=<file> Write partials in stream format to <file>.\n" +" -g. --geometry=<file> Get detector geometry from file.\n" +" -b, --beam=<file> Get beam parameters from file\n" +" -p, --pdb=<file> PDB file from which to get the unit cell.\n" +"\n" +" -y, --symmetry=<sym> Symmetry of the input reflection list.\n" +); +} + + +int main(int argc, char *argv[]) +{ + int c; + char *input_file = NULL; + char *output_file = NULL; + char *beamfile = NULL; + char *geomfile = NULL; + char *cellfile = NULL; + struct detector *det = NULL; + struct beam_params *beam = NULL; + RefList *full; + char *sym = NULL; + UnitCell *cell = NULL; + UnitCell *rcell; + struct quaternion orientation; + struct image image; + FILE *ofh; + + /* Long options */ + const struct option longopts[] = { + {"help", 0, NULL, 'h'}, + {"output", 1, NULL, 'o'}, + {"input", 1, NULL, 'i'}, + {"beam", 1, NULL, 'b'}, + {"pdb", 1, NULL, 'p'}, + {"geometry", 1, NULL, 'g'}, + {"symmetry", 1, NULL, 'y'}, + {0, 0, NULL, 0} + }; + + /* Short options */ + while ((c = getopt_long(argc, argv, "hi:o:b:p:g:y:", + longopts, NULL)) != -1) { + + switch (c) { + case 'h' : + show_help(argv[0]); + return 0; + + case 'o' : + output_file = strdup(optarg); + break; + + case 'i' : + input_file = strdup(optarg); + break; + + case 'b' : + beamfile = strdup(optarg); + break; + + case 'p' : + cellfile = strdup(optarg); + break; + + case 'g' : + geomfile = strdup(optarg); + break; + + case 'y' : + sym = strdup(optarg); + break; + + case 0 : + break; + + default : + return 1; + } + + } + + /* Load beam */ + if ( beamfile == NULL ) { + ERROR("You need to provide a beam parameters file.\n"); + return 1; + } + beam = get_beam_parameters(beamfile); + if ( beam == NULL ) { + ERROR("Failed to load beam parameters from '%s'\n", beamfile); + return 1; + } + free(beamfile); + + /* Load cell */ + if ( cellfile == NULL ) { + ERROR("You need to give a PDB file with the unit cell.\n"); + return 1; + } + cell = load_cell_from_pdb(cellfile); + if ( cell == NULL ) { + ERROR("Failed to get cell from '%s'\n", cellfile); + return 1; + } + free(cellfile); + + /* Load geometry */ + if ( geomfile == NULL ) { + ERROR("You need to give a geometry file.\n"); + return 1; + } + det = get_detector_geometry(geomfile); + if ( cell == NULL ) { + ERROR("Failed to read geometry from '%s'\n", geomfile); + return 1; + } + free(geomfile); + + /* Load (full) reflections */ + if ( input_file == NULL ) { + ERROR("You must provide a reflection list.\n"); + return 1; + } + full = read_reflections(input_file); + if ( full == NULL ) { + ERROR("Failed to read reflections from '%s'\n", input_file); + return 1; + } + free(input_file); + if ( check_list_symmetry(full, sym) ) { + ERROR("The input reflection list does not appear to" + " have symmetry %s\n", sym); + return 1; + } + + if ( output_file == NULL ) { + ERROR("You must pgive a filename for the output.\n"); + return 1; + } + ofh = fopen(output_file, "w"); + if ( ofh == NULL ) { + ERROR("Couldn't open output file '%s'\n", output_file); + return 1; + } + free(output_file); + write_stream_header(ofh, argc, argv); + + /* Set up a random orientation */ + orientation = random_quaternion(); + image.indexed_cell = cell_rotate(cell, orientation); + + image.det = det; + image.width = det->max_fs; + image.height = det->max_ss; + + image.lambda = ph_en_to_lambda(eV_to_J(beam->photon_energy)); + image.div = beam->divergence; + image.bw = beam->bandwidth; + image.profile_radius = 0.005e9; + image.i0_available = 0; + image.filename = "(simulated 1)"; + image.reflections = find_intersections(&image, image.indexed_cell, 0); + calculate_partials(image.reflections, full, sym); + write_chunk(ofh, &image, STREAM_INTEGRATED); + reflist_free(image.reflections); + + /* Rotate the cell by a tiny amount */ + rcell = rotate_cell(image.indexed_cell, + deg2rad(0.2), deg2rad(0.0), deg2rad(0.0)); + cell_free(image.indexed_cell); + image.indexed_cell = rcell; + image.filename = "(simulated 2)"; + + /* Write another chunk */ + image.reflections = find_intersections(&image, image.indexed_cell, 0); + calculate_partials(image.reflections, full, sym); + write_chunk(ofh, &image, STREAM_INTEGRATED); + reflist_free(image.reflections); + + fclose(ofh); + cell_free(cell); + free_detector_geometry(det); + free(beam); + free(sym); + reflist_free(full); + + return 0; +} |