aboutsummaryrefslogtreecommitdiff
path: root/src/indexamajig.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/indexamajig.c')
-rw-r--r--src/indexamajig.c224
1 files changed, 133 insertions, 91 deletions
diff --git a/src/indexamajig.c b/src/indexamajig.c
index a1997b9e..ce87fde6 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -3,16 +3,17 @@
*
* Index patterns, output hkl+intensity etc.
*
- * Copyright © 2012-2015 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2017 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
* Copyright © 2012 Richard Kirian
* Copyright © 2012 Lorenzo Galli
*
* Authors:
- * 2010-2015 Thomas White <taw@physics.org>
+ * 2010-2017 Thomas White <taw@physics.org>
* 2011 Richard Kirian
* 2012 Lorenzo Galli
* 2012 Chunhong Yoon
+ * 2017 Valerio Mariani <valerio.mariani@desy.de>
*
* This file is part of CrystFEL.
*
@@ -88,6 +89,7 @@ static void show_help(const char *s)
" --peaks=<method> Use 'method' for finding peaks. Choose from:\n"
" zaef : Use Zaefferer (2000) gradient detection.\n"
" This is the default method.\n"
+" peakfinder8: Use Peakfinder8 algorithm.\n"
" hdf5 : Get from a table in HDF5 file.\n"
" cxi : Get from CXI format HDF5 file.\n"
" --hdf5-peaks=<p> Find peaks table in HDF5 file here.\n"
@@ -95,32 +97,52 @@ static void show_help(const char *s)
" --integration=<meth> Perform final pattern integration using <meth>.\n"
"\n\n"
"For more control over the process, you might need:\n\n"
-" --tolerance=<tol> Set the tolerances for cell comparison.\n"
-" Default: 5,5,5,1.5.\n"
-" --filter-noise Apply an aggressive noise filter which sets all\n"
-" pixels in each 3x3 region to zero if any of them\n"
-" have negative values. Intensity measurement will\n"
-" be performed on the image as it was before this.\n"
-" --median-filter=<n> Apply a median filter to the image data. Intensity\n"
-" measurement will be performed on the image as it\n"
-" was before this. The side length of the median\n"
-" filter box will be 2<n>+1. Default: 0 (no filter).\n"
-" --no-sat-corr Don't correct values of saturated peaks using a\n"
-" table included in the HDF5 file.\n"
-" --threshold=<n> Only accept peaks above <n> ADU. Default: 800.\n"
-" --min-gradient=<n> Minimum squared gradient for Zaefferer peak search.\n"
-" Default: 100,000.\n"
-" --min-snr=<n> Minimum signal-to-noise ratio for peaks.\n"
-" Default: 5.\n"
-" --check-hdf5-snr Check SNR for peaks from --peaks=hdf5.\n"
-" --peak-radius=<r> Integration radii for peak search.\n"
-" --int-radius=<r> Set the integration radii. Default: 4,5,7.\n"
-" --push-res=<n> Integrate higher than apparent resolution cutoff.\n"
-" --highres=<n> Absolute resolution cutoff in Angstroms.\n"
-" --fix-profile-radius Fix the reciprocal space profile radius for spot\n"
-" prediction (default: automatically determine.\n"
-" --fix-bandwidth Set the bandwidth for spot prediction.\n"
-" --fix-divergence Set the divergence (full angle) for spot prediction.\n"
+" --tolerance=<tol> Set the tolerances for cell comparison.\n"
+" Default: 5,5,5,1.5.\n"
+" --filter-noise Apply an aggressive noise filter which sets all\n"
+" pixels in each 3x3 region to zero if any of them\n"
+" have negative values. Intensity measurement will\n"
+" be performed on the image as it was before this.\n"
+" --median-filter=<n> Apply a median filter to the image data. Intensity\n"
+" measurement will be performed on the image as it\n"
+" was before this. The side length of the median\n"
+" filter box will be 2<n>+1.\n"
+" Default: 0 (no filter).\n"
+" --no-sat-corr Don't correct values of saturated peaks using a\n"
+" table included in the HDF5 file.\n"
+" --threshold=<n> Only accept peaks above <n> ADU in both the\n"
+" Zaefferer and Peakfinder8 algorithms.\n"
+" Default: 800.\n"
+" --min-gradient=<n> Minimum squared gradient for Zaefferer peak\n"
+" search. Default: 100,000.\n"
+" --min-snr=<n> Minimum signal-to-noise ratio for peaks, with both\n"
+" Zaefferer and Peakfinder8 algorithms.\n"
+" Default: 5.\n"
+" --min-pix-count=<n> Only accept peaks if they include more than <n>\n"
+" pixels, in the Peakfinder8 algorithm.\n"
+" Default: 2.\n"
+" --max-pix-count=<n> Only accept peaks if they include less than <n>\n"
+" pixels, in the Peakfinder8 algorithm.\n"
+" Default: 200.\n"
+" --min-peaks=<n> Minimum number of peaks for indexing.\n"
+" --local-bg-radius=<n> Radius (in pixels) to use for the estimation of\n"
+" local background in the Peakfinder8 algorithm.\n"
+" Default: 3.\n"
+" --min-res=<n> Only accept peaks if they lay at more than <n>\n"
+" pixels from the center of the detector, in the\n"
+" peakfinder8 algorithm. Default: 0.\n"
+" --max-res=<n> Only accept peaks if they lay at less than <n>\n"
+" pixels from the center of the detector, in the\n"
+" peakfinder8 algorithm. Default: 1200.\n"
+" --check-hdf5-snr Check SNR for peaks from --peaks=hdf5.\n"
+" --peak-radius=<r> Integration radii for peak search.\n"
+" --int-radius=<r> Set the integration radii. Default: 4,5,7.\n"
+" --push-res=<n> Integrate higher than apparent resolution cutoff.\n"
+" --highres=<n> Absolute resolution cutoff in Angstroms.\n"
+" --fix-profile-radius Fix the reciprocal space profile radius for spot\n"
+" prediction (default: automatically determine.\n"
+" --fix-bandwidth Set the bandwidth for spot prediction.\n"
+" --fix-divergence Set the divergence (full angle) for spot prediction.\n"
"\n"
"\nFor time-resolved stuff, you might want to use:\n\n"
" --copy-hdf5-field <f> Copy the value of field <f> into the stream. You\n"
@@ -131,16 +153,19 @@ static void show_help(const char *s)
" --temp-dir=<path> Put the temporary folder under <path>.\n"
"\n"
"\nOptions you probably won't need:\n\n"
-" --no-check-prefix Don't attempt to correct the --prefix.\n"
-" --no-use-saturated During the initial peak search, reject\n"
-" peaks which contain pixels above max_adu.\n"
-" --no-revalidate Don't re-integrate and check HDF5 peaks for\n"
-" validity.\n"
-" --no-peaks-in-stream Do not record peak search results in the stream.\n"
-" --no-refls-in-stream Do not record integrated reflections in the stream.\n"
-" --int-diag=<cond> Show debugging information about reflections.\n"
-" --no-refine Skip the prediction refinement step.\n"
-" --profile Show timing data for performance monitoring.\n"
+" --no-check-prefix Don't attempt to correct the --prefix.\n"
+" --no-use-saturated During the initial peak search, reject\n"
+" peaks which contain pixels above max_adu.\n"
+" --no-revalidate Don't re-integrate and check HDF5 peaks for\n"
+" validity.\n"
+" --no-peaks-in-stream Do not record peak search results in the stream.\n"
+" --no-refls-in-stream Do not record integrated reflections in the stream.\n"
+" --no-non-hits-in-stream Do not include non-hit frames in the stream.\n"
+" (see --min-peaks)\n"
+" --int-diag=<cond> Show debugging information about reflections.\n"
+" --no-refine Skip the prediction refinement step.\n"
+" --profile Show timing data for performance monitoring.\n"
+" --no-half-pixel-shift Don't offset the HDF5 peak locations by 0.5 px.\n"
"\nLow-level options for the felix indexer:\n\n"
" --felix-options Change the default arguments passed to the indexer.\n"
" Given as a list of comma separated list of \n"
@@ -150,9 +175,9 @@ static void show_help(const char *s)
}
-static void add_geom_beam_stuff_to_copy_hdf5(struct copy_hdf5_field *copyme,
- struct detector *det,
- struct beam_params *beam)
+static void add_geom_beam_stuff_to_field_list(struct imagefile_field_list *copyme,
+ struct detector *det,
+ struct beam_params *beam)
{
int i;
@@ -161,12 +186,12 @@ static void add_geom_beam_stuff_to_copy_hdf5(struct copy_hdf5_field *copyme,
struct panel *p = &det->panels[i];
if ( p->clen_from != NULL ) {
- add_copy_hdf5_field(copyme, p->clen_from);
+ add_imagefile_field(copyme, p->clen_from);
}
}
if ( beam->photon_energy_from != NULL ) {
- add_copy_hdf5_field(copyme, beam->photon_energy_from);
+ add_imagefile_field(copyme, beam->photon_energy_from);
}
}
@@ -181,8 +206,6 @@ int main(int argc, char *argv[])
int config_checkprefix = 1;
int config_basename = 0;
int integrate_saturated = 0;
- IndexingMethod *indm;
- IndexingPrivate **ipriv;
char *indm_str = NULL;
char *cellfile = NULL;
char *prefix = NULL;
@@ -214,11 +237,17 @@ int main(int argc, char *argv[])
iargs.threshold = 800.0;
iargs.min_gradient = 100000.0;
iargs.min_snr = 5.0;
+ iargs.min_pix_count = 2;
+ iargs.max_pix_count = 200;
+ iargs.min_res = 0;
+ iargs.max_res = 1200;
+ iargs.local_bg_radius = 3;
iargs.check_hdf5_snr = 0;
iargs.det = NULL;
iargs.peaks = PEAK_ZAEF;
iargs.beam = &beam;
iargs.hdf5_peak_path = NULL;
+ iargs.half_pixel_shift = 1;
iargs.copyme = NULL;
iargs.pk_inn = -1.0;
iargs.pk_mid = -1.0;
@@ -230,13 +259,14 @@ int main(int argc, char *argv[])
iargs.no_revalidate = 0;
iargs.stream_peaks = 1;
iargs.stream_refls = 1;
+ iargs.stream_nonhits = 1;
iargs.int_diag = INTDIAG_NONE;
- iargs.copyme = new_copy_hdf5_field_list();
+ iargs.copyme = new_imagefile_field_list();
+ iargs.min_peaks = 0;
if ( iargs.copyme == NULL ) {
ERROR("Couldn't allocate HDF5 field list.\n");
return 1;
}
- iargs.indm = NULL; /* No default */
iargs.ipriv = NULL; /* No default */
iargs.int_meth = integration_method("rings-nocen-nosat-nograd", NULL);
iargs.push_res = 0.0;
@@ -268,12 +298,14 @@ int main(int argc, char *argv[])
{"basename", 0, &config_basename, 1},
{"no-peaks-in-stream", 0, &iargs.stream_peaks, 0},
{"no-refls-in-stream", 0, &iargs.stream_refls, 0},
+ {"no-non-hits-in-stream", 0, &iargs.stream_nonhits, 0},
{"integrate-saturated",0, &integrate_saturated, 1},
{"no-use-saturated", 0, &iargs.use_saturated, 0},
{"no-revalidate", 0, &iargs.no_revalidate, 1},
{"check-hdf5-snr", 0, &iargs.check_hdf5_snr, 1},
{"no-refine", 0, &no_refine, 1},
{"profile", 0, &iargs.profile, 1},
+ {"no-half-pixel-shift",0, &iargs.half_pixel_shift, 0},
/* Long-only options which don't actually do anything */
{"no-sat-corr", 0, &iargs.satcorr, 0},
@@ -306,6 +338,12 @@ int main(int argc, char *argv[])
{"fix-bandwidth", 1, NULL, 23},
{"fix-divergence", 1, NULL, 24},
{"felix-options", 1, NULL, 25},
+ {"min-pix-count", 1, NULL, 26},
+ {"max-pix-count", 1, NULL, 27},
+ {"local-bg-radius", 1, NULL, 28},
+ {"min-res", 1, NULL, 29},
+ {"max-res", 1, NULL, 30},
+ {"min-peaks", 1, NULL, 31},
{0, 0, NULL, 0}
};
@@ -400,7 +438,7 @@ int main(int argc, char *argv[])
break;
case 10 :
- add_copy_hdf5_field(iargs.copyme, optarg);
+ add_imagefile_field(iargs.copyme, optarg);
break;
case 11 :
@@ -489,6 +527,30 @@ int main(int argc, char *argv[])
}
break;
+ case 26:
+ iargs.min_pix_count = atoi(optarg);
+ break;
+
+ case 27:
+ iargs.max_pix_count = atoi(optarg);
+ break;
+
+ case 28:
+ iargs.local_bg_radius = atoi(optarg);
+ break;
+
+ case 29:
+ iargs.min_res = atoi(optarg);
+ break;
+
+ case 30:
+ iargs.max_res = atoi(optarg);
+ break;
+
+ case 31:
+ iargs.min_peaks = atoi(optarg);
+ break;
+
case 0 :
break;
@@ -541,6 +603,8 @@ int main(int argc, char *argv[])
}
if ( strcmp(speaks, "zaef") == 0 ) {
iargs.peaks = PEAK_ZAEF;
+ } else if ( strcmp(speaks, "peakfinder8") == 0 ) {
+ iargs.peaks = PEAK_PEAKFINDER8;
} else if ( strcmp(speaks, "hdf5") == 0 ) {
iargs.peaks = PEAK_HDF5;
} else if ( strcmp(speaks, "cxi") == 0 ) {
@@ -574,7 +638,7 @@ int main(int argc, char *argv[])
geom_filename);
return 1;
}
- add_geom_beam_stuff_to_copy_hdf5(iargs.copyme, iargs.det, iargs.beam);
+ add_geom_beam_stuff_to_field_list(iargs.copyme, iargs.det, iargs.beam);
/* If no peak path from geometry file, use these (but see later) */
if ( iargs.hdf5_peak_path == NULL ) {
@@ -591,35 +655,6 @@ int main(int argc, char *argv[])
iargs.hdf5_peak_path = command_line_peak_path;
}
- /* Parse indexing methods */
- if ( indm_str == NULL ) {
-
- STATUS("You didn't specify an indexing method, so I won't try "
- " to index anything.\n"
- "If that isn't what you wanted, re-run with"
- " --indexing=<methods>.\n");
- indm = NULL;
-
- } else {
-
- int i = 0;
-
- indm = build_indexer_list(indm_str);
- if ( indm == NULL ) {
- ERROR("Invalid indexer list '%s'\n", indm_str);
- return 1;
- }
- free(indm_str);
-
- /* If --no-refine, unset the refinement flag on all methods */
- if ( no_refine ) {
- while ( indm[i] != INDEXING_NONE ) {
- indm[i] &= ~INDEXING_REFINE;
- i++;
- }
- }
- }
-
/* Parse integration method */
if ( int_str != NULL ) {
@@ -755,27 +790,34 @@ int main(int argc, char *argv[])
}
free(outfile);
- /* Prepare the indexer */
- if ( indm != NULL ) {
- ipriv = prepare_indexing(indm, iargs.cell, iargs.det,
- iargs.tols, iargs.felix_options);
- if ( ipriv == NULL ) {
- ERROR("Failed to prepare indexing.\n");
+ /* Prepare the indexing system */
+ if ( indm_str == NULL ) {
+
+ STATUS("You didn't specify an indexing method, so I won't try "
+ " to index anything.\n"
+ "If that isn't what you wanted, re-run with"
+ " --indexing=<methods>.\n");
+ iargs.ipriv = NULL;
+
+ } else {
+
+ iargs.ipriv = setup_indexing(indm_str, iargs.cell, iargs.det,
+ iargs.tols, no_refine,
+ iargs.felix_options);
+ if ( iargs.ipriv == NULL ) {
+ ERROR("Failed to set up indexing system\n");
return 1;
}
- } else {
- ipriv = NULL;
+ free(indm_str);
+
}
gsl_set_error_handler_off();
- iargs.indm = indm;
- iargs.ipriv = ipriv;
-
create_sandbox(&iargs, n_proc, prefix, config_basename, fh,
st, tempdir);
- free_copy_hdf5_field_list(iargs.copyme);
+ free_imagefile_field_list(iargs.copyme);
cell_free(iargs.cell);
free(iargs.beam->photon_energy_from);
free(prefix);
@@ -783,7 +825,7 @@ int main(int argc, char *argv[])
free_detector_geometry(iargs.det);
free(iargs.hdf5_peak_path);
close_stream(st);
- cleanup_indexing(indm, ipriv);
+ cleanup_indexing(iargs.ipriv);
return 0;
}