aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2011-04-13 17:31:13 +0200
committerThomas White <taw@physics.org>2012-02-22 15:27:24 +0100
commitde215b2a46b3452e8d540a536ecaf53dd839e4b0 (patch)
tree1e0264cf4982a9bb9bfb34b7922c772b024177a1
parent36b436334514124e444fc5ccb42e1dc1d1bb25bc (diff)
Add partial_sim, for generating test streams by geometrical methods
-rw-r--r--.gitignore1
-rw-r--r--Makefile.am10
-rw-r--r--Makefile.in46
-rw-r--r--src/partial_sim.c263
4 files changed, 308 insertions, 12 deletions
diff --git a/.gitignore b/.gitignore
index 5a1f78d5..906b0c75 100644
--- a/.gitignore
+++ b/.gitignore
@@ -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;
+}