From 65deebdcc51030002d1f210cc121a06255055e67 Mon Sep 17 00:00:00 2001 From: Yaroslav Gevorkov Date: Fri, 26 Apr 2019 16:23:38 +0200 Subject: Add pinkIndexer interface --- doc/man/indexamajig.1 | 37 +++++ libcrystfel/CMakeLists.txt | 3 + libcrystfel/config.h.cmake.in | 2 + libcrystfel/src/index.c | 31 +++- libcrystfel/src/index.h | 6 +- libcrystfel/src/pinkindexer.c | 361 ++++++++++++++++++++++++++++++++++++++++++ libcrystfel/src/pinkindexer.h | 46 ++++++ src/indexamajig.c | 128 +++++++++++++++ src/process_image.h | 2 + 9 files changed, 613 insertions(+), 3 deletions(-) create mode 100644 libcrystfel/src/pinkindexer.c create mode 100644 libcrystfel/src/pinkindexer.h diff --git a/doc/man/indexamajig.1 b/doc/man/indexamajig.1 index 5b51566f..3d96cf0b 100644 --- a/doc/man/indexamajig.1 +++ b/doc/man/indexamajig.1 @@ -90,6 +90,11 @@ Use the TakeTwo algorithm. See Ginn et al., Acta Cryst. (2016). D72, 956-965. .PD Invoke XGANDALF - eXtended GrAdieNt Descent Algorithm for Lattice Finding. Xgandalf must be installed in order to be able to use it. +.IP \fBpinkIndexer\fR +.PD +Invoke pinkIndexer. pinkIndexer must be installed in order to be able to use it. The geometry file should contain the value photon_energy_bandwidth, which sets the bandwidth as (delta energy)/(mean energy). + + .PP You can add one or more of the following to the above indexing methods, to control what information should be provided to them. Note that indexamajig performs a series of checks on the indexing results, including checking that the result is consistent with the target unit cell parameters. To get completely "raw" indexing, you need to disable these checks (see below) \fBand\fR not provide prior information. @@ -437,6 +442,38 @@ These set low-level parameters for the XGANDALF indexing algorithm. .IP \fB--xgandalf-fast-execution\fR Shortcut to set --xgandalf-sampling-pitch=2 --xgandalf-grad-desc-iterations=3 +.PD 0 +.IP \fB--pinkIndexer-considered-peaks-count=\fIn\fR +.IP \fB--pinkIndexer-angle-resolution=\fIn\fR +.IP \fB--pinkIndexer-refinement-type=\fIn\fR +.IP \fB--pinkIndexer-tolerance=\fIn\fR +.IP \fB--pinkIndexer-reflection-radius=\fIn\fR +.IP \fB--pinkIndexer-max-resolution-for-indexing=\fIn\fR +.IP \fB--pinkIndexer-multi=\fIn\fR +.IP \fB--pinkIndexer-thread-count=\fIn\fR +.IP \fB--pinkIndexer-no-check-indexed=\fIn\fR + +.PD +These set low-level parameters for the XGANDALF indexing algorithm. +.IP +\fB--pinkIndexer-considered-peaks-count\fR selects how many peaks are considered for indexing. [0-4] (veryFew to manyMany). Default is 4 (manyMany). +.IP +\fB--pinkIndexer-angle-resolution\fR selects how dense the orientation angles of the sample lattice are sampled. [0-4] (extremelyLoose to extremelyDense). Default is 2 (normal). +.IP +\fB--pinkIndexer-refinement-type\fR selects the refinement type. 0 = none, 1 = fixedLatticeParameters, 2 = variableLatticeParameters, 3 = firstFixedThenVariableLatticeParameters, 4 = firstFixedThenVariableLatticeParametersMultiSeed, 5 = firstFixedThenVariableLatticeParametersCenterAdjustmentMultiSeed. +.IP +\fB--pinkIndexer-tolerance\fR selects the tolerance of the pinkIndexer (relative tolerance of the lattice vectors). Default is 0.06. For bad geometrys or cell parameters use a high tolerance. For a well known geometry and cell use a small tolerance. Only important for refinement and indexed/not indexed identificaton. Too small tolerance will lead to refining to only a fraction of the peaks and possibly discarding of correctly indexed images. Too high tolerance will lead to bad fitting in presence of multiples or noise and can mark wrongly-indexed patterns as indexed. +.IP +\fB--pinkIndexer-reflection-radius\fR sets radius of the reflections in reciprocal space in 1/A. Default is 2%% of a* (which works quiet well for X-rays). Should be chosen much bigger for electrons (~0.002). +.IP +\fB--pinkIndexer-max-resolution-for-indexing\fR sets the maximum resolition in 1/A used for indexing. Peaks at high resolution don't add much information, but they add a lot of computation time. Default is infinity. Does not influence the refinement. +.IP +\fB--pinkIndexer-multi\fR Use pinkIndexers own multi indexing. Should be combined with the --no-multi flag. +.IP +\fB--pinkIndexer-thread-count\fR sets the thread count for internal parallelization. Default is 1. Very useful for small datasets (e.g. for screening). Internal parallelization does not significantly increase the amount of RAM needed, whereas CrystFELs parallelization does. For HPCs typically a mixture of both parallelizations leads to best results. +.IP +\fB--pinkIndexer-no-check-indexed\fR Leave the check whether a pattern is indexed completely to CrystFEL. Useful for monochromatic (since CrystFELs prediction model is smarter than the one of pinkIndexer) or in combnation with --no-check-peaks for geometry optimization. This flag is meant to eventually disappear, when the full pink pipeline is implemented. + .SH INTEGRATION OPTIONS .PD 0 .IP \fB--integration=\fR\fImethod\fR diff --git a/libcrystfel/CMakeLists.txt b/libcrystfel/CMakeLists.txt index 67a20ee5..fdba725e 100644 --- a/libcrystfel/CMakeLists.txt +++ b/libcrystfel/CMakeLists.txt @@ -14,6 +14,7 @@ pkg_search_module(FFTW fftw3) set(HAVE_CURSES ${CURSES_FOUND}) set(HAVE_FFTW ${FFTW_FOUND}) set(HAVE_XGANDALF ${XGANDALF_FOUND}) +set(HAVE_PINKINDEXER ${PINKINDEXER_FOUND}) set(HAVE_FDIP ${FDIP_FOUND}) # Find out where forkpty() is declared @@ -67,6 +68,7 @@ set(LIBCRYSTFEL_SOURCES src/peakfinder8.c src/taketwo.c src/xgandalf.c + src/pinkindexer.c src/rational.c src/spectrum.c ${BISON_symopp_OUTPUTS} @@ -108,6 +110,7 @@ set(LIBCRYSTFEL_HEADERS src/peakfinder8.h src/taketwo.h src/xgandalf.h + src/pinkindexer.h src/rational.h src/spectrum.h ) diff --git a/libcrystfel/config.h.cmake.in b/libcrystfel/config.h.cmake.in index ea26c1a2..b7ad4c7a 100644 --- a/libcrystfel/config.h.cmake.in +++ b/libcrystfel/config.h.cmake.in @@ -4,6 +4,8 @@ #cmakedefine HAVE_CPU_AFFINITY #cmakedefine HAVE_FFTW #cmakedefine HAVE_XGANDALF +#cmakedefine HAVE_NBP +#cmakedefine HAVE_PINKINDEXER #cmakedefine HAVE_FDIP #cmakedefine HAVE_CURSES diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 4c82e868..cc004d9a 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -58,6 +58,7 @@ #include "predict-refine.h" #include "taketwo.h" #include "xgandalf.h" +#include "pinkindexer.h" /** \file index.h */ @@ -71,6 +72,7 @@ struct _indexingprivate struct taketwo_options *ttopts; struct xgandalf_options *xgandalf_opts; + struct pinkIndexer_options *pinkIndexer_opts; int n_methods; IndexingMethod *methods; @@ -195,6 +197,10 @@ static char *base_indexer_str(IndexingMethod indm) strcpy(str, "xgandalf"); break; + case INDEXING_PINKINDEXER: + strcpy(str, "pinkIndexer"); + break; + case INDEXING_SIMULATION : strcpy(str, "simulation"); break; @@ -233,6 +239,7 @@ static char *friendly_indexer_name(IndexingMethod m) static void *prepare_method(IndexingMethod *m, UnitCell *cell, struct xgandalf_options *xgandalf_opts, + struct pinkIndexer_options* pinkIndexer_opts, struct felix_options *felix_opts) { char *str; @@ -277,6 +284,10 @@ static void *prepare_method(IndexingMethod *m, UnitCell *cell, priv = xgandalf_prepare(m, cell, xgandalf_opts); break; + case INDEXING_PINKINDEXER : + priv = pinkIndexer_prepare(m, cell, pinkIndexer_opts); + break; + default : ERROR("Don't know how to prepare indexing method %i\n", *m); break; @@ -309,6 +320,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, IndexingFlags flags, struct taketwo_options *ttopts, struct xgandalf_options *xgandalf_opts, + struct pinkIndexer_options *pinkIndexer_opts, struct felix_options *felix_opts) { int i, n; @@ -403,7 +415,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, int j; ipriv->engine_private[i] = prepare_method(&methods[i], cell, - xgandalf_opts, + xgandalf_opts, pinkIndexer_opts, felix_opts); if ( ipriv->engine_private[i] == NULL ) return NULL; @@ -430,6 +442,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, ipriv->ttopts = ttopts; ipriv->xgandalf_opts = xgandalf_opts; + ipriv->pinkIndexer_opts = pinkIndexer_opts; STATUS("List of indexing methods:\n"); for ( i=0; iengine_private[n]); break; + case INDEXING_PINKINDEXER : + pinkIndexer_cleanup(ipriv->engine_private[n]); + break; + default : ERROR("Don't know how to clean up indexing method %i\n", ipriv->methods[n]); @@ -601,6 +618,10 @@ static int try_indexer(struct image *image, IndexingMethod indm, r = taketwo_index(image, ipriv->ttopts, mpriv); break; + case INDEXING_PINKINDEXER : + r = run_pinkIndexer(image, mpriv); + break; + case INDEXING_XGANDALF : set_last_task(last_task, "indexing:xgandalf"); r = run_xgandalf(image, mpriv); @@ -999,6 +1020,11 @@ IndexingMethod get_indm_from_string_2(const char *str, int *err) method = INDEXING_DEFAULTS_XGANDALF; have_method = 1; + } else if ( strcmp(bits[i], "pinkIndexer") == 0) { + if ( have_method ) return warn_method(str); + method = INDEXING_DEFAULTS_PINKINDEXER; + have_method = 1; + } else if ( strcmp(bits[i], "none") == 0) { if ( have_method ) return warn_method(str); method = INDEXING_NONE; @@ -1097,9 +1123,10 @@ char *detect_indexing_methods(UnitCell *cell) do_probe(asdf_probe, cell, methods); do_probe(xds_probe, cell, methods); do_probe(xgandalf_probe, cell, methods); - /* Don't automatically use TakeTwo or Felix (yet) */ + /* Don't automatically use TakeTwo, Felix or PinkIndexer (yet) */ //do_probe(taketwo_probe, cell, methods); //do_probe(felix_probe, cell, methods); + //do_probe(pinkIndexer_probe, cell, methods); if ( strlen(methods) == 0 ) { free(methods); diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h index ec7ed61b..ab38016f 100644 --- a/libcrystfel/src/index.h +++ b/libcrystfel/src/index.h @@ -63,6 +63,8 @@ #define INDEXING_DEFAULTS_XGANDALF (INDEXING_XGANDALF | INDEXING_USE_CELL_PARAMETERS) +#define INDEXING_DEFAULTS_PINKINDEXER (INDEXING_PINKINDEXER | INDEXING_USE_CELL_PARAMETERS) + /** * An enumeration of all the available indexing methods. **/ @@ -79,6 +81,7 @@ typedef enum { INDEXING_ASDF = 8, /**< Use built-in ASDF algorithm */ INDEXING_TAKETWO = 9, /**< Use built-in TakeTwo algorithm */ INDEXING_XGANDALF = 10, /**< Use XGANDALF (via optional library) */ + INDEXING_PINKINDEXER = 11,/**< Use PinkIndexer (via optional library) */ INDEXING_ERROR = 255, /**< Special value for unrecognised indexing * engine */ @@ -146,14 +149,15 @@ extern IndexingMethod get_indm_from_string_2(const char *method, int *err); #include "image.h" #include "taketwo.h" #include "xgandalf.h" +#include "pinkindexer.h" #include "felix.h" - extern IndexingPrivate *setup_indexing(const char *methods, UnitCell *cell, struct detector *det, float *ltl, IndexingFlags flags, struct taketwo_options *ttopts, struct xgandalf_options *xgandalf_opts, + struct pinkIndexer_options *pinkIndexer_opts, struct felix_options *felix_opts); extern char *detect_indexing_methods(UnitCell *cell); diff --git a/libcrystfel/src/pinkindexer.c b/libcrystfel/src/pinkindexer.c new file mode 100644 index 00000000..2002533f --- /dev/null +++ b/libcrystfel/src/pinkindexer.c @@ -0,0 +1,361 @@ +/* + * pinkIndexer.c + * + * Created on: Nov 27, 2017 + * Author: gevorkov + */ + +#include "pinkindexer.h" + +#ifdef HAVE_PINKINDEXER + +#include + +#include "utils.h" +#include "cell-utils.h" +#include "peaks.h" + +#include "pinkIndexer/adaptions/crystfel/Lattice.h" +#include "pinkIndexer/adaptions/crystfel/ExperimentSettings.h" +#include "pinkIndexer/adaptions/crystfel/PinkIndexer.h" + +#define MAX_MULTI_LATTICE_COUNT 8 + +struct pinkIndexer_private_data { + PinkIndexer *pinkIndexer; + reciprocalPeaks_1_per_A_t reciprocalPeaks_1_per_A; + float *intensities; + + IndexingMethod indm; + UnitCell *cellTemplate; + int threadCount; + int multi; + int min_peaks; + + int no_check_indexed; + + IntegerMatrix *centeringTransformation; + LatticeTransform_t latticeReductionTransform; +}; + +//static void reduceCell(UnitCell* cell, LatticeTransform_t* appliedReductionTransform); +//static void restoreCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform); +static void reduceReciprocalCell(UnitCell* cell, LatticeTransform_t* appliedReductionTransform); +static void restoreReciprocalCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform); +static void makeRightHanded(UnitCell* cell); +static void update_detector(struct detector *det, double xoffs, double yoffs); + +int run_pinkIndexer(struct image *image, void *ipriv) +{ + struct pinkIndexer_private_data* pinkIndexer_private_data = (struct pinkIndexer_private_data*) ipriv; + reciprocalPeaks_1_per_A_t* reciprocalPeaks_1_per_A = &(pinkIndexer_private_data->reciprocalPeaks_1_per_A); + float *intensities = pinkIndexer_private_data->intensities; + + int peakCountMax = image_feature_count(image->features); + if (peakCountMax < 5) { + int goodLatticesCount = 0; + return goodLatticesCount; + } + reciprocalPeaks_1_per_A->peakCount = 0; + for (int i = 0; i < peakCountMax && i < MAX_PEAK_COUNT_FOR_INDEXER; i++) { + struct imagefeature *f; + f = image_get_feature(image->features, i); + if (f == NULL) { + continue; + } + + reciprocalPeaks_1_per_A->coordinates_x[reciprocalPeaks_1_per_A->peakCount] = f->rz * 1e-10; + reciprocalPeaks_1_per_A->coordinates_y[reciprocalPeaks_1_per_A->peakCount] = f->rx * 1e-10; + reciprocalPeaks_1_per_A->coordinates_z[reciprocalPeaks_1_per_A->peakCount] = f->ry * 1e-10; + intensities[reciprocalPeaks_1_per_A->peakCount] = (float) (f->intensity); + reciprocalPeaks_1_per_A->peakCount++; + } + int indexed = 0; + Lattice_t indexedLattice[MAX_MULTI_LATTICE_COUNT]; + float center_shift[MAX_MULTI_LATTICE_COUNT][2]; + + float maxRefinementDisbalance = 0.4; + + do { + int peakCount = reciprocalPeaks_1_per_A->peakCount; + printf("\ntotal peaks count %i\n", peakCount); + int matchedPeaksCount = PinkIndexer_indexPattern(pinkIndexer_private_data->pinkIndexer, + &(indexedLattice[indexed]), center_shift[indexed], reciprocalPeaks_1_per_A, intensities, + maxRefinementDisbalance, + pinkIndexer_private_data->threadCount); + + printf("center shift = %f %f", center_shift[indexed][0], center_shift[indexed][1]); + printf("\nmatchedPeaksCount %i from %i\n", matchedPeaksCount, peakCount); + + if ((matchedPeaksCount >= 25 && matchedPeaksCount >= peakCount * 0.30) + || matchedPeaksCount >= peakCount * 0.4 + || matchedPeaksCount >= 70 + || pinkIndexer_private_data->no_check_indexed == 1) + { + UnitCell *uc; + uc = cell_new(); + + Lattice_t *l = &(indexedLattice[indexed]); + + cell_set_reciprocal(uc, l->ay * 1e10, l->az * 1e10, l->ax * 1e10, + l->by * 1e10, l->bz * 1e10, l->bx * 1e10, + l->cy * 1e10, l->cz * 1e10, l->cx * 1e10); + printf("before restoration\n"); + cell_print(uc); + + restoreReciprocalCell(uc, &pinkIndexer_private_data->latticeReductionTransform); + + UnitCell *new_cell_trans = cell_transform_intmat(uc, pinkIndexer_private_data->centeringTransformation); + cell_free(uc); + uc = new_cell_trans; + printf("after restoration\n"); + cell_print(uc); + + cell_set_lattice_type(new_cell_trans, cell_get_lattice_type(pinkIndexer_private_data->cellTemplate)); + cell_set_centering(new_cell_trans, cell_get_centering(pinkIndexer_private_data->cellTemplate)); + cell_set_unique_axis(new_cell_trans, cell_get_unique_axis(pinkIndexer_private_data->cellTemplate)); + + if (validate_cell(uc)) { + printf("pinkIndexer: problem with returned cell!\n"); + } + + Crystal * cr = crystal_new(); + if (cr == NULL) { + ERROR("Failed to allocate crystal.\n"); + return 0; + } + crystal_set_cell(cr, uc); + crystal_set_det_shift(cr, center_shift[indexed][0], center_shift[indexed][1]); + update_detector(image->det, center_shift[indexed][0], center_shift[indexed][1]); + image_add_crystal(image, cr); + indexed++; + + printf("crystal %i in image added\n", indexed); + } + else { + break; + } + } while (pinkIndexer_private_data->multi + && indexed <= MAX_MULTI_LATTICE_COUNT + && reciprocalPeaks_1_per_A->peakCount >= pinkIndexer_private_data->min_peaks); + + printf("\n\nfound %i crystals in image\n", indexed); + return indexed; +} + +void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, + struct pinkIndexer_options *pinkIndexer_opts) +{ + printf("preparing pink indexer\n"); + + if (pinkIndexer_opts->beamEnergy == 0.0 || pinkIndexer_opts->detectorDistance <= 0) { + ERROR("ERROR!!!!!! photon_energy and photon_energy_bandwidth must be defined as constants " + "in the geometry file for the pinkIndexer!!!!!!"); + } + if (pinkIndexer_opts->beamBandwidth == 0.0) { + printf("using default bandwidth of 0.01 for pinkIndexer!"); + pinkIndexer_opts->beamBandwidth = 0.01; + } + if (pinkIndexer_opts->detectorDistance == 0.0 && pinkIndexer_opts->refinement_type == + REFINEMENT_TYPE_firstFixedThenVariableLatticeParametersCenterAdjustmentMultiSeed) { + ERROR("Using center refinement makes it necessary to have the detector distance fixed in the geometry file!"); + } + + struct pinkIndexer_private_data* pinkIndexer_private_data = malloc(sizeof(struct pinkIndexer_private_data)); + allocReciprocalPeaks(&(pinkIndexer_private_data->reciprocalPeaks_1_per_A)); + pinkIndexer_private_data->intensities = malloc(MAX_PEAK_COUNT_FOR_INDEXER * sizeof(float)); + pinkIndexer_private_data->indm = *indm; + pinkIndexer_private_data->cellTemplate = cell; + pinkIndexer_private_data->threadCount = pinkIndexer_opts->thread_count; + pinkIndexer_private_data->multi = pinkIndexer_opts->multi; + pinkIndexer_private_data->min_peaks = pinkIndexer_opts->min_peaks; + pinkIndexer_private_data->no_check_indexed = pinkIndexer_opts->no_check_indexed; + + UnitCell* primitiveCell = uncenter_cell(cell, &pinkIndexer_private_data->centeringTransformation, NULL); + + //reduceCell(primitiveCell, &pinkIndexer_private_data->latticeReductionTransform); + reduceReciprocalCell(primitiveCell, &pinkIndexer_private_data->latticeReductionTransform); + + printf("reduced cell:\n"); + cell_print(primitiveCell); + + double asx, asy, asz, bsx, bsy, bsz, csx, csy, csz; + int ret = cell_get_reciprocal(primitiveCell, &asx, &asy, &asz, &bsx, &bsy, &bsz, &csx, &csy, &csz); + if (ret != 0) { + ERROR("cell_get_reciprocal did not finish properly!"); + } + + Lattice_t lattice = { .ax = asz * 1e-10, .ay = asx * 1e-10, .az = asy * 1e-10, + .bx = bsz * 1e-10, .by = bsx * 1e-10, .bz = bsy * 1e-10, + .cx = csz * 1e-10, .cy = csx * 1e-10, .cz = csy * 1e-10 }; + + float detectorDistance_m = pinkIndexer_opts->detectorDistance; + float beamEenergy_eV = pinkIndexer_opts->beamEnergy; + float nonMonochromaticity = pinkIndexer_opts->beamBandwidth; + float reflectionRadius_1_per_A; + if (pinkIndexer_opts->reflectionRadius < 0) { + reflectionRadius_1_per_A = 0.02 + * sqrt(lattice.ax * lattice.ax + lattice.ay * lattice.ay + lattice.az * lattice.az); + } + else { + reflectionRadius_1_per_A = pinkIndexer_opts->reflectionRadius; + } + + float divergenceAngle_deg = 0.01; //fake + + float tolerance = pinkIndexer_opts->tolerance; + Lattice_t sampleReciprocalLattice_1_per_A = lattice; + float detectorRadius_m = 0.03; //fake, only for prediction + ExperimentSettings* experimentSettings = ExperimentSettings_new(beamEenergy_eV, detectorDistance_m, + detectorRadius_m, divergenceAngle_deg, nonMonochromaticity, sampleReciprocalLattice_1_per_A, tolerance, + reflectionRadius_1_per_A); + + consideredPeaksCount_t consideredPeaksCount = pinkIndexer_opts->considered_peaks_count; + angleResolution_t angleResolution = pinkIndexer_opts->angle_resolution; + refinementType_t refinementType = pinkIndexer_opts->refinement_type; + float maxResolutionForIndexing_1_per_A = pinkIndexer_opts->maxResolutionForIndexing_1_per_A; + pinkIndexer_private_data->pinkIndexer = PinkIndexer_new(experimentSettings, consideredPeaksCount, angleResolution, + refinementType, + maxResolutionForIndexing_1_per_A); + + ExperimentSettings_delete(experimentSettings); + cell_free(primitiveCell); + + /* Flags that pinkIndexer knows about */ + *indm &= INDEXING_METHOD_MASK + | INDEXING_USE_CELL_PARAMETERS; + + return pinkIndexer_private_data; +} + +//static void reduceCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform) +//{ +// double ax, ay, az, bx, by, bz, cx, cy, cz; +// cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); +// +// Lattice_t l = { ax, ay, az, bx, by, bz, cx, cy, cz }; +// +// reduceLattice(&l, appliedReductionTransform); +// +// cell_set_cartesian(cell, l.ax, l.ay, l.az, +// l.bx, l.by, l.bz, +// l.cx, l.cy, l.cz); +// +// makeRightHanded(cell); +//} +// +//static void restoreCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform) +//{ +// +// double ax, ay, az, bx, by, bz, cx, cy, cz; +// cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); +// +// Lattice_t l = { ax, ay, az, bx, by, bz, cx, cy, cz }; +// +// restoreLattice(&l, appliedReductionTransform); +// +// cell_set_cartesian(cell, l.ax, l.ay, l.az, +// l.bx, l.by, l.bz, +// l.cx, l.cy, l.cz); +// +// makeRightHanded(cell); +//} + +static void reduceReciprocalCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform) +{ + double ax, ay, az, bx, by, bz, cx, cy, cz; + cell_get_reciprocal(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); + + Lattice_t l = { ax, ay, az, bx, by, bz, cx, cy, cz }; + + reduceLattice(&l, appliedReductionTransform); + + cell_set_reciprocal(cell, l.ax, l.ay, l.az, + l.bx, l.by, l.bz, + l.cx, l.cy, l.cz); + + makeRightHanded(cell); +} + +static void restoreReciprocalCell(UnitCell *cell, LatticeTransform_t* appliedReductionTransform) +{ + + double ax, ay, az, bx, by, bz, cx, cy, cz; + cell_get_reciprocal(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); + + Lattice_t l = { ax, ay, az, bx, by, bz, cx, cy, cz }; + + restoreLattice(&l, appliedReductionTransform); + + cell_set_reciprocal(cell, l.ax, l.ay, l.az, + l.bx, l.by, l.bz, + l.cx, l.cy, l.cz); + + makeRightHanded(cell); +} + +static void makeRightHanded(UnitCell *cell) +{ + double ax, ay, az, bx, by, bz, cx, cy, cz; + cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); + + if (!right_handed(cell)) { + cell_set_cartesian(cell, -ax, -ay, -az, -bx, -by, -bz, -cx, -cy, -cz); + } +} + +//hack for electron crystallography while crystal_set_det_shift is not working approprietly +static void update_detector(struct detector *det, double xoffs, double yoffs) +{ + int i; + + for (i = 0; i < det->n_panels; i++) { + struct panel *p = &det->panels[i]; + p->cnx += xoffs * p->res; + p->cny += yoffs * p->res; + } +} + +void pinkIndexer_cleanup(void *pp) +{ + printf("cleaning up pink indexer\n"); + + struct pinkIndexer_private_data* pinkIndexer_private_data = (struct pinkIndexer_private_data*) pp; + + freeReciprocalPeaks(pinkIndexer_private_data->reciprocalPeaks_1_per_A); + free(pinkIndexer_private_data->intensities); + intmat_free(pinkIndexer_private_data->centeringTransformation); + PinkIndexer_delete(pinkIndexer_private_data->pinkIndexer); +} + +const char *pinkIndexer_probe(UnitCell *cell) +{ + return "pinkIndexer"; +} + +#else /* HAVE_PINKINDEXER */ + +int run_pinkIndexer(struct image *image, void *ipriv) +{ + ERROR("This copy of CrystFEL was compiled without PINKINDEXER support.\n"); + return 0; +} + +extern void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, + struct pinkIndexer_options *pinkIndexer_opts) +{ + ERROR("This copy of CrystFEL was compiled without PINKINDEXER support.\n"); + ERROR("To use PINKINDEXER indexing, recompile with PINKINDEXER.\n"); + return NULL; +} + +void pinkIndexer_cleanup(void *pp) +{ +} + +const char *pinkIndexer_probe(UnitCell *cell) +{ + return NULL; +} + +#endif /* HAVE_PINKINDEXER */ diff --git a/libcrystfel/src/pinkindexer.h b/libcrystfel/src/pinkindexer.h new file mode 100644 index 00000000..cd2aaba2 --- /dev/null +++ b/libcrystfel/src/pinkindexer.h @@ -0,0 +1,46 @@ +/* + * pinkIndexer.h + * + * Created on: Nov 27, 2017 + * Author: gevorkov + */ + +#ifndef LIBCRYSTFEL_SRC_PINKINDEXER_H_ +#define LIBCRYSTFEL_SRC_PINKINDEXER_H_ + +#ifdef HAVE_CONFIG_H +#include +#endif + +struct pinkIndexer_options { + unsigned int considered_peaks_count; + unsigned int angle_resolution; + unsigned int refinement_type; + float maxResolutionForIndexing_1_per_A; + float tolerance; + int multi; + int thread_count; + int min_peaks; + + int no_check_indexed; + + float beamEnergy; //in eV + float beamBandwidth; //(delta lambda)/lambda + float detectorDistance; //in m + + float reflectionRadius; //in 1/A +}; + +#include +#include "index.h" + +extern int run_pinkIndexer(struct image *image, void *ipriv); + +extern void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, + struct pinkIndexer_options *pinkIndexer_opts); + +extern void pinkIndexer_cleanup(void *pp); + +extern const char *pinkIndexer_probe(UnitCell *cell); + +#endif /* LIBCRYSTFEL_SRC_PINKINDEXER_H_ */ diff --git a/src/indexamajig.c b/src/indexamajig.c index bd159a9a..5770051a 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -204,6 +204,32 @@ static void show_help(const char *s) " All peaks are used for refinement.\n" " Default: 250\n" "\n" +" --pinkIndexer-considered-peaks-count Considered peaks count selector \n" +" [0-4] (veryFew to manyMany)\n" +" Default is 4 (manyMany)\n" +" --pinkIndexer-angle-resolution Angle resolution selector \n" +" [0-4] (extremelyLoose to extremelyDense)\n" +" Default is 2 (normal)\n" +" --pinkIndexer-refinement-type Refinement type \n" +" Default is 1\n" +" 0 = none\n" +" 1 = fixedLatticeParameters\n" +" 2 = variableLatticeParameters\n" +" 3 = firstFixedThenVariableLatticeParameters\n" +" 4 = firstFixedThenVariableLatticeParametersMultiSeed\n" +" 5 = firstFixedThenVariableLatticeParametersCenterAdjustmentMultiSeed\n" +" --pinkIndexer-tolerance Relative tolerance of the lattice vectors.\n" +" Default is 0.06\n" +" --pinkIndexer-reflection-radius radius of the reflections in reciprocal space in 1/A.\n" +" Default is 2%% of a*.\n" +" Should be chosen ~0.002 for electrons.\n" +" --pinkIndexer-max-resolution-for-indexing Measured in 1/A\n" +" --pinkIndexer-multi Use pinkIndexers own multi indexing.\n" +" --pinkIndexer-thread-count thread count for internal parallelization \n" +" Default is 1\n" +" --pinkIndexer-no-check-indexed Leave the check whether a pattern is indexed completely to CrystFEL\n" +" useful for monochromatic and tweaking." +"\n" "\nIntegration options:\n\n" " --integration= Integration method (rings,prof2d)-(cen,nocen)\n" " Default: rings-nocen\n" @@ -253,6 +279,7 @@ static void add_geom_beam_stuff_to_field_list(struct imagefile_field_list *copym int main(int argc, char *argv[]) { int c; + unsigned int tmp_enum; char *filename = NULL; char *outfile = NULL; FILE *fh; @@ -360,6 +387,16 @@ int main(int argc, char *argv[]) iargs.xgandalf_opts.minLatticeVectorLength_A = 30; iargs.xgandalf_opts.maxLatticeVectorLength_A = 250; iargs.xgandalf_opts.maxPeaksForIndexing = 250; + iargs.pinkIndexer_opts.considered_peaks_count = 4; + iargs.pinkIndexer_opts.angle_resolution = 2; + iargs.pinkIndexer_opts.refinement_type = 1; + iargs.pinkIndexer_opts.tolerance = 0.06; + iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A = FLT_MAX; + iargs.pinkIndexer_opts.thread_count = 1; + iargs.pinkIndexer_opts.multi = 0; + iargs.pinkIndexer_opts.no_check_indexed = 0; + iargs.pinkIndexer_opts.min_peaks = 2; + iargs.pinkIndexer_opts.reflectionRadius = -1; iargs.felix_opts.ttmin = -1.0; iargs.felix_opts.ttmax = -1.0; iargs.felix_opts.min_visits = 0; @@ -491,6 +528,21 @@ int main(int argc, char *argv[]) {"min-sq-gradient", 1, NULL, 359}, /* compat */ {"xgandalf-fast-execution", 0, NULL, 360}, {"xgandalf-max-peaks", 1, NULL, 361}, + {"pinkIndexer-considered-peaks-count", 1, NULL, 362}, + {"pinkIndexer-cpc", 1, NULL, 362}, + {"pinkIndexer-angle-resolution", 1, NULL, 363}, + {"pinkIndexer-ar", 1, NULL, 363}, + {"pinkIndexer-refinement-type", 1, NULL, 364}, + {"pinkIndexer-rt", 1, NULL, 364}, + {"pinkIndexer-thread-count", 1, NULL, 365}, + {"pinkIndexer-tc", 1, NULL, 365}, + {"pinkIndexer-max-resolution-for-indexing", 1, NULL, 366}, + {"pinkIndexer-mrfi", 1, NULL, 366}, + {"pinkIndexer-tolerance", 1, NULL, 367}, + {"pinkIndexer-tol", 1, NULL, 367}, + {"pinkIndexer-multi", 0, NULL, 368}, + {"pinkIndexer-no-check-indexed", 0, NULL, 369}, + {"pinkIndexer-reflection-radius", 1, NULL, 370}, {0, 0, NULL, 0} }; @@ -689,6 +741,7 @@ int main(int argc, char *argv[]) case 331: iargs.min_peaks = atoi(optarg); + iargs.pinkIndexer_opts.min_peaks = iargs.min_peaks; break; case 332: @@ -922,6 +975,74 @@ int main(int argc, char *argv[]) } break; + case 362: + if (sscanf(optarg, "%u", &tmp_enum) != 1) + { + ERROR("Invalid value for " + "--pinkIndexer-considered-peaks-count\n"); + return 1; + } + iargs.pinkIndexer_opts.considered_peaks_count = tmp_enum; + break; + + case 363: + if (sscanf(optarg, "%u", &tmp_enum) != 1) + { + ERROR("Invalid value for --pinkIndexer-angle-resolution \n"); + return 1; + } + iargs.pinkIndexer_opts.angle_resolution = tmp_enum; + break; + + case 364: + if (sscanf(optarg, "%u", &tmp_enum) != 1) + { + ERROR("Invalid value for --pinkIndexer-refinement-type \n"); + return 1; + } + iargs.pinkIndexer_opts.refinement_type = tmp_enum; + break; + + case 365: + if (sscanf(optarg, "%d", &iargs.pinkIndexer_opts.thread_count) != 1) + { + ERROR("Invalid value for --pinkIndexer-thread-count \n"); + return 1; + } + break; + + case 366: + if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A) != 1) + { + ERROR("Invalid value for --pinkIndexer-max-resolution-for-indexing \n"); + return 1; + } + break; + + case 367: + if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.tolerance) != 1) + { + ERROR("Invalid value for --pinkIndexer-tolerance \n"); + return 1; + } + break; + + case 368: + iargs.pinkIndexer_opts.multi = 1; + break; + + case 369: + iargs.pinkIndexer_opts.no_check_indexed = 1; + break; + + case 370: + if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.reflectionRadius) != 1) + { + ERROR("Invalid value for --pinkIndexer-reflection-radius \n"); + return 1; + } + break; + case 0 : break; @@ -1021,6 +1142,12 @@ int main(int argc, char *argv[]) return 1; } add_geom_beam_stuff_to_field_list(iargs.copyme, iargs.det, iargs.beam); + iargs.pinkIndexer_opts.beamEnergy = iargs.beam->photon_energy; + iargs.pinkIndexer_opts.beamBandwidth = iargs.beam->photon_energy_bandwidth; + iargs.pinkIndexer_opts.detectorDistance = iargs.det->panels[0].clen; + if(iargs.det->panels[0].clen_from != NULL){ + iargs.pinkIndexer_opts.detectorDistance = 0; + } /* If no peak path from geometry file, use these (but see later) */ if ( iargs.hdf5_peak_path == NULL ) { @@ -1263,6 +1390,7 @@ int main(int argc, char *argv[]) iargs.tols, flags, &iargs.taketwo_opts, &iargs.xgandalf_opts, + &iargs.pinkIndexer_opts, &iargs.felix_opts); if ( iargs.ipriv == NULL ) { ERROR("Failed to set up indexing system\n"); diff --git a/src/process_image.h b/src/process_image.h index 9a58a64f..ced3bcf2 100644 --- a/src/process_image.h +++ b/src/process_image.h @@ -46,6 +46,7 @@ struct index_args; #include "time-accounts.h" #include "taketwo.h" #include "xgandalf.h" +#include "../libcrystfel/src/pinkindexer.h" #include "felix.h" @@ -117,6 +118,7 @@ struct index_args int profile; /* Whether or not to do wall clock profiling */ struct taketwo_options taketwo_opts; struct xgandalf_options xgandalf_opts; + struct pinkIndexer_options pinkIndexer_opts; struct felix_options felix_opts; Spectrum *spectrum; signed int wait_for_file; /* -1 means wait forever */ -- cgit v1.2.3 From 4e1d9e63c06169baf5c1966aeddc958591890d9a Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 16:27:50 +0200 Subject: Use +INFINITY instead of FLT_MAX --- src/indexamajig.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indexamajig.c b/src/indexamajig.c index 5770051a..a67c1a5d 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -391,7 +391,7 @@ int main(int argc, char *argv[]) iargs.pinkIndexer_opts.angle_resolution = 2; iargs.pinkIndexer_opts.refinement_type = 1; iargs.pinkIndexer_opts.tolerance = 0.06; - iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A = FLT_MAX; + iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A = +INFINITY; iargs.pinkIndexer_opts.thread_count = 1; iargs.pinkIndexer_opts.multi = 0; iargs.pinkIndexer_opts.no_check_indexed = 0; -- cgit v1.2.3 From 70ce07879e5fb6b72681f86fcf7545f7f4a95f98 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 16:28:00 +0200 Subject: Use hard-coded default bandwidth for pinkIndexer (Temporarily) --- src/indexamajig.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indexamajig.c b/src/indexamajig.c index a67c1a5d..cded5f34 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -1143,7 +1143,7 @@ int main(int argc, char *argv[]) } add_geom_beam_stuff_to_field_list(iargs.copyme, iargs.det, iargs.beam); iargs.pinkIndexer_opts.beamEnergy = iargs.beam->photon_energy; - iargs.pinkIndexer_opts.beamBandwidth = iargs.beam->photon_energy_bandwidth; + iargs.pinkIndexer_opts.beamBandwidth = 0.01; iargs.pinkIndexer_opts.detectorDistance = iargs.det->panels[0].clen; if(iargs.det->panels[0].clen_from != NULL){ iargs.pinkIndexer_opts.detectorDistance = 0; -- cgit v1.2.3 From 91a16342b42900292cf828f22516a932a2c598eb Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 16:33:05 +0200 Subject: pinkIndexer: Formatting, and get rid of tmp_enum --- src/indexamajig.c | 51 ++++++++++++++++++++++++++++++--------------------- 1 file changed, 30 insertions(+), 21 deletions(-) diff --git a/src/indexamajig.c b/src/indexamajig.c index cded5f34..da5c3f7e 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -279,7 +279,6 @@ static void add_geom_beam_stuff_to_field_list(struct imagefile_field_list *copym int main(int argc, char *argv[]) { int c; - unsigned int tmp_enum; char *filename = NULL; char *outfile = NULL; FILE *fh; @@ -940,7 +939,7 @@ int main(int argc, char *argv[]) &iargs.xgandalf_opts.maxLatticeVectorLength_A) != 1) { ERROR("Invalid value for " - "--xgandalf-max-lattice-vector-length\n"); + "--xgandalf-max-lattice-vector-length\n"); return 1; } break; @@ -976,69 +975,79 @@ int main(int argc, char *argv[]) break; case 362: - if (sscanf(optarg, "%u", &tmp_enum) != 1) + if (sscanf(optarg, "%u", + &iargs.pinkIndexer_opts.considered_peaks_count) != 1) { ERROR("Invalid value for " - "--pinkIndexer-considered-peaks-count\n"); + "--pinkIndexer-considered-peaks-count\n"); return 1; } - iargs.pinkIndexer_opts.considered_peaks_count = tmp_enum; break; case 363: - if (sscanf(optarg, "%u", &tmp_enum) != 1) + if (sscanf(optarg, "%u", + &iargs.pinkIndexer_opts.angle_resolution) != 1) { - ERROR("Invalid value for --pinkIndexer-angle-resolution \n"); + ERROR("Invalid value for " + "--pinkIndexer-angle-resolution\n"); return 1; } - iargs.pinkIndexer_opts.angle_resolution = tmp_enum; break; case 364: - if (sscanf(optarg, "%u", &tmp_enum) != 1) + if (sscanf(optarg, "%u", + &iargs.pinkIndexer_opts.refinement_type) != 1) { - ERROR("Invalid value for --pinkIndexer-refinement-type \n"); + ERROR("Invalid value for " + "--pinkIndexer-refinement-type\n"); return 1; } - iargs.pinkIndexer_opts.refinement_type = tmp_enum; break; case 365: - if (sscanf(optarg, "%d", &iargs.pinkIndexer_opts.thread_count) != 1) + if (sscanf(optarg, "%d", + &iargs.pinkIndexer_opts.thread_count) != 1) { - ERROR("Invalid value for --pinkIndexer-thread-count \n"); + ERROR("Invalid value for " + "--pinkIndexer-thread-count\n"); return 1; } break; case 366: - if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A) != 1) + if (sscanf(optarg, "%f", + &iargs.pinkIndexer_opts.maxResolutionForIndexing_1_per_A) != 1) { - ERROR("Invalid value for --pinkIndexer-max-resolution-for-indexing \n"); + ERROR("Invalid value for " + "--pinkIndexer-max-resolution-for-indexing\n"); return 1; } break; case 367: - if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.tolerance) != 1) + if (sscanf(optarg, "%f", + &iargs.pinkIndexer_opts.tolerance) != 1) { - ERROR("Invalid value for --pinkIndexer-tolerance \n"); + ERROR("Invalid value for " + "--pinkIndexer-tolerance\n"); return 1; } break; case 368: - iargs.pinkIndexer_opts.multi = 1; + iargs.pinkIndexer_opts.multi = 1; break; case 369: - iargs.pinkIndexer_opts.no_check_indexed = 1; + iargs.pinkIndexer_opts.no_check_indexed = 1; break; case 370: - if (sscanf(optarg, "%f", &iargs.pinkIndexer_opts.reflectionRadius) != 1) + if (sscanf(optarg, "%f", + &iargs.pinkIndexer_opts.reflectionRadius) != 1) { - ERROR("Invalid value for --pinkIndexer-reflection-radius \n"); + ERROR("Invalid value for " + "--pinkIndexer-reflection-radius\n"); return 1; } break; -- cgit v1.2.3 From 8ccf4f93017ce6637ded0763693ddfe39daa568f Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 17:20:47 +0200 Subject: Whitespace --- libcrystfel/src/index.c | 4 ++-- src/process_image.h | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index cc004d9a..89c0667a 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -504,8 +504,8 @@ void cleanup_indexing(IndexingPrivate *ipriv) break; case INDEXING_PINKINDEXER : - pinkIndexer_cleanup(ipriv->engine_private[n]); - break; + pinkIndexer_cleanup(ipriv->engine_private[n]); + break; default : ERROR("Don't know how to clean up indexing method %i\n", diff --git a/src/process_image.h b/src/process_image.h index ced3bcf2..0b266de8 100644 --- a/src/process_image.h +++ b/src/process_image.h @@ -118,7 +118,7 @@ struct index_args int profile; /* Whether or not to do wall clock profiling */ struct taketwo_options taketwo_opts; struct xgandalf_options xgandalf_opts; - struct pinkIndexer_options pinkIndexer_opts; + struct pinkIndexer_options pinkIndexer_opts; struct felix_options felix_opts; Spectrum *spectrum; signed int wait_for_file; /* -1 means wait forever */ -- cgit v1.2.3 From fa5d0f9561642e7dd5edc94297b313cf6e55b35b Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 17:20:56 +0200 Subject: Add GPL boilerplate --- libcrystfel/src/pinkindexer.c | 27 ++++++++++++++++++++++++--- libcrystfel/src/pinkindexer.h | 27 ++++++++++++++++++++++++--- 2 files changed, 48 insertions(+), 6 deletions(-) diff --git a/libcrystfel/src/pinkindexer.c b/libcrystfel/src/pinkindexer.c index 2002533f..3653a355 100644 --- a/libcrystfel/src/pinkindexer.c +++ b/libcrystfel/src/pinkindexer.c @@ -1,8 +1,29 @@ /* - * pinkIndexer.c + * pinkindexer.c + * + * Interface to PinkIndexer + * + * Copyright © 2017-2019 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2017-2019 Yaroslav Gevorkov + * + * 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 . * - * Created on: Nov 27, 2017 - * Author: gevorkov */ #include "pinkindexer.h" diff --git a/libcrystfel/src/pinkindexer.h b/libcrystfel/src/pinkindexer.h index cd2aaba2..5be2c1e7 100644 --- a/libcrystfel/src/pinkindexer.h +++ b/libcrystfel/src/pinkindexer.h @@ -1,8 +1,29 @@ /* - * pinkIndexer.h + * pinkindexer.h + * + * Interface to PinkIndexer + * + * Copyright © 2017-2019 Deutsches Elektronen-Synchrotron DESY, + * a research centre of the Helmholtz Association. + * + * Authors: + * 2017-2019 Yaroslav Gevorkov + * + * 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 . * - * Created on: Nov 27, 2017 - * Author: gevorkov */ #ifndef LIBCRYSTFEL_SRC_PINKINDEXER_H_ -- cgit v1.2.3 From 7f5a04465249744d4ba1d0a9f50c350ec7aaaab2 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 17:21:05 +0200 Subject: Format help items for pinkIndexer --- src/indexamajig.c | 50 ++++++++++++++++++++++++-------------------------- 1 file changed, 24 insertions(+), 26 deletions(-) diff --git a/src/indexamajig.c b/src/indexamajig.c index da5c3f7e..9aef3f2e 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -204,32 +204,30 @@ static void show_help(const char *s) " All peaks are used for refinement.\n" " Default: 250\n" "\n" -" --pinkIndexer-considered-peaks-count Considered peaks count selector \n" -" [0-4] (veryFew to manyMany)\n" -" Default is 4 (manyMany)\n" -" --pinkIndexer-angle-resolution Angle resolution selector \n" -" [0-4] (extremelyLoose to extremelyDense)\n" -" Default is 2 (normal)\n" -" --pinkIndexer-refinement-type Refinement type \n" -" Default is 1\n" -" 0 = none\n" -" 1 = fixedLatticeParameters\n" -" 2 = variableLatticeParameters\n" -" 3 = firstFixedThenVariableLatticeParameters\n" -" 4 = firstFixedThenVariableLatticeParametersMultiSeed\n" -" 5 = firstFixedThenVariableLatticeParametersCenterAdjustmentMultiSeed\n" -" --pinkIndexer-tolerance Relative tolerance of the lattice vectors.\n" -" Default is 0.06\n" -" --pinkIndexer-reflection-radius radius of the reflections in reciprocal space in 1/A.\n" -" Default is 2%% of a*.\n" -" Should be chosen ~0.002 for electrons.\n" -" --pinkIndexer-max-resolution-for-indexing Measured in 1/A\n" -" --pinkIndexer-multi Use pinkIndexers own multi indexing.\n" -" --pinkIndexer-thread-count thread count for internal parallelization \n" -" Default is 1\n" -" --pinkIndexer-no-check-indexed Leave the check whether a pattern is indexed completely to CrystFEL\n" -" useful for monochromatic and tweaking." -"\n" +" --pinkIndexer-considered-peaks-count\n" +" Considered peaks count, 0 (fewest) to 4 (most)\n" +" Default: 4\n" +" --pinkIndexer-angle-resolution\n" +" Angle resolution, 0 (loosest) to 4 (most dense)\n" +" Default: 2\n" +" --pinkIndexer-refinement-type\n" +" Refinement type, 0 (none) to 5 (most accurate)\n" +" Default: 1\n" +" --pinkIndexer-tolerance\n" +" Relative tolerance of the lattice vectors.\n" +" Default 0.06\n" +" --pinkIndexer-reflection-radius\n" +" Radius of the reflections in reciprocal space.\n" +" Specified in 1/A. Default is 2%% of a*.\n" +" --pinkIndexer-max-resolution-for-indexing\n" +" Measured in 1/A\n" +" --pinkIndexer-multi Use pinkIndexers own multi indexing.\n" +" --pinkIndexer-thread-count\n" +" Thread count for internal parallelization \n" +" Default: 1\n" +" --pinkIndexer-no-check-indexed\n" +" Disable internal check for correct indexing\n" +" solutions\n" "\nIntegration options:\n\n" " --integration= Integration method (rings,prof2d)-(cen,nocen)\n" " Default: rings-nocen\n" -- cgit v1.2.3 From 3a12faaa0b22762c218eb1a71b3144fea8c56d64 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 26 Apr 2019 17:21:18 +0200 Subject: Fix header include path --- src/process_image.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/process_image.h b/src/process_image.h index 0b266de8..5691a0e4 100644 --- a/src/process_image.h +++ b/src/process_image.h @@ -46,7 +46,7 @@ struct index_args; #include "time-accounts.h" #include "taketwo.h" #include "xgandalf.h" -#include "../libcrystfel/src/pinkindexer.h" +#include "pinkindexer.h" #include "felix.h" -- cgit v1.2.3 From 9805a1aa7475de72bd9da251949490345a93987c Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 29 Apr 2019 13:09:04 +0200 Subject: Set task name for PinkIndexer --- libcrystfel/src/index.c | 1 + 1 file changed, 1 insertion(+) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 89c0667a..2653ec74 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -619,6 +619,7 @@ static int try_indexer(struct image *image, IndexingMethod indm, break; case INDEXING_PINKINDEXER : + set_last_task(last_task, "indexing:pinkindexer"); r = run_pinkIndexer(image, mpriv); break; -- cgit v1.2.3 From a700f8e170a00a58d1d78351bc4297ac41f01aae Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 29 Apr 2019 13:09:14 +0200 Subject: Accept --indexing=pinkindexer (lower case) --- libcrystfel/src/index.c | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 2653ec74..80a773ab 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -1021,7 +1021,9 @@ IndexingMethod get_indm_from_string_2(const char *str, int *err) method = INDEXING_DEFAULTS_XGANDALF; have_method = 1; - } else if ( strcmp(bits[i], "pinkIndexer") == 0) { + } else if ( (strcmp(bits[i], "pinkIndexer") == 0) + || (strcmp(bits[i], "pinkindexer") == 0) ) + { if ( have_method ) return warn_method(str); method = INDEXING_DEFAULTS_PINKINDEXER; have_method = 1; -- cgit v1.2.3 From 4620e00fd1fb75392fe516b3358fa2eaec3e189c Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 29 Apr 2019 13:09:35 +0200 Subject: pinkIndexer: Remove debugging messages --- libcrystfel/src/pinkindexer.c | 30 +++++++----------------------- 1 file changed, 7 insertions(+), 23 deletions(-) diff --git a/libcrystfel/src/pinkindexer.c b/libcrystfel/src/pinkindexer.c index 3653a355..8195002b 100644 --- a/libcrystfel/src/pinkindexer.c +++ b/libcrystfel/src/pinkindexer.c @@ -99,15 +99,11 @@ int run_pinkIndexer(struct image *image, void *ipriv) do { int peakCount = reciprocalPeaks_1_per_A->peakCount; - printf("\ntotal peaks count %i\n", peakCount); int matchedPeaksCount = PinkIndexer_indexPattern(pinkIndexer_private_data->pinkIndexer, &(indexedLattice[indexed]), center_shift[indexed], reciprocalPeaks_1_per_A, intensities, maxRefinementDisbalance, pinkIndexer_private_data->threadCount); - printf("center shift = %f %f", center_shift[indexed][0], center_shift[indexed][1]); - printf("\nmatchedPeaksCount %i from %i\n", matchedPeaksCount, peakCount); - if ((matchedPeaksCount >= 25 && matchedPeaksCount >= peakCount * 0.30) || matchedPeaksCount >= peakCount * 0.4 || matchedPeaksCount >= 70 @@ -121,23 +117,19 @@ int run_pinkIndexer(struct image *image, void *ipriv) cell_set_reciprocal(uc, l->ay * 1e10, l->az * 1e10, l->ax * 1e10, l->by * 1e10, l->bz * 1e10, l->bx * 1e10, l->cy * 1e10, l->cz * 1e10, l->cx * 1e10); - printf("before restoration\n"); - cell_print(uc); restoreReciprocalCell(uc, &pinkIndexer_private_data->latticeReductionTransform); UnitCell *new_cell_trans = cell_transform_intmat(uc, pinkIndexer_private_data->centeringTransformation); cell_free(uc); uc = new_cell_trans; - printf("after restoration\n"); - cell_print(uc); cell_set_lattice_type(new_cell_trans, cell_get_lattice_type(pinkIndexer_private_data->cellTemplate)); cell_set_centering(new_cell_trans, cell_get_centering(pinkIndexer_private_data->cellTemplate)); cell_set_unique_axis(new_cell_trans, cell_get_unique_axis(pinkIndexer_private_data->cellTemplate)); if (validate_cell(uc)) { - printf("pinkIndexer: problem with returned cell!\n"); + ERROR("pinkIndexer: problem with returned cell!\n"); } Crystal * cr = crystal_new(); @@ -151,30 +143,27 @@ int run_pinkIndexer(struct image *image, void *ipriv) image_add_crystal(image, cr); indexed++; - printf("crystal %i in image added\n", indexed); - } - else { + } else { break; } } while (pinkIndexer_private_data->multi && indexed <= MAX_MULTI_LATTICE_COUNT && reciprocalPeaks_1_per_A->peakCount >= pinkIndexer_private_data->min_peaks); - printf("\n\nfound %i crystals in image\n", indexed); return indexed; } void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, struct pinkIndexer_options *pinkIndexer_opts) { - printf("preparing pink indexer\n"); - if (pinkIndexer_opts->beamEnergy == 0.0 || pinkIndexer_opts->detectorDistance <= 0) { - ERROR("ERROR!!!!!! photon_energy and photon_energy_bandwidth must be defined as constants " - "in the geometry file for the pinkIndexer!!!!!!"); + ERROR("For pinkIndexer, the photon_energy and " + "photon_energy_bandwidth must be defined as constants in " + "the geometry file\n"); + return NULL; } if (pinkIndexer_opts->beamBandwidth == 0.0) { - printf("using default bandwidth of 0.01 for pinkIndexer!"); + STATUS("Using default bandwidth of 0.01 for pinkIndexer\n"); pinkIndexer_opts->beamBandwidth = 0.01; } if (pinkIndexer_opts->detectorDistance == 0.0 && pinkIndexer_opts->refinement_type == @@ -197,9 +186,6 @@ void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, //reduceCell(primitiveCell, &pinkIndexer_private_data->latticeReductionTransform); reduceReciprocalCell(primitiveCell, &pinkIndexer_private_data->latticeReductionTransform); - printf("reduced cell:\n"); - cell_print(primitiveCell); - double asx, asy, asz, bsx, bsy, bsz, csx, csy, csz; int ret = cell_get_reciprocal(primitiveCell, &asx, &asy, &asz, &bsx, &bsy, &bsz, &csx, &csy, &csz); if (ret != 0) { @@ -339,8 +325,6 @@ static void update_detector(struct detector *det, double xoffs, double yoffs) void pinkIndexer_cleanup(void *pp) { - printf("cleaning up pink indexer\n"); - struct pinkIndexer_private_data* pinkIndexer_private_data = (struct pinkIndexer_private_data*) pp; freeReciprocalPeaks(pinkIndexer_private_data->reciprocalPeaks_1_per_A); -- cgit v1.2.3 From 10f9432595f5425ff2f8605ea181f7c643ab81ab Mon Sep 17 00:00:00 2001 From: Yaroslav Gevorkov Date: Mon, 29 Apr 2019 21:38:47 +0200 Subject: bug fix: asking for detector distance != 0 instad of beamBandwidth != 0 --- libcrystfel/src/pinkindexer.c | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/libcrystfel/src/pinkindexer.c b/libcrystfel/src/pinkindexer.c index 8195002b..fe3b1625 100644 --- a/libcrystfel/src/pinkindexer.c +++ b/libcrystfel/src/pinkindexer.c @@ -156,10 +156,9 @@ int run_pinkIndexer(struct image *image, void *ipriv) void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, struct pinkIndexer_options *pinkIndexer_opts) { - if (pinkIndexer_opts->beamEnergy == 0.0 || pinkIndexer_opts->detectorDistance <= 0) { - ERROR("For pinkIndexer, the photon_energy and " - "photon_energy_bandwidth must be defined as constants in " - "the geometry file\n"); + if (pinkIndexer_opts->beamEnergy == 0.0) { + ERROR("For pinkIndexer, the photon_energy must be defined as a " + "constant in the geometry file\n"); return NULL; } if (pinkIndexer_opts->beamBandwidth == 0.0) { @@ -171,6 +170,10 @@ void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, ERROR("Using center refinement makes it necessary to have the detector distance fixed in the geometry file!"); } + if(pinkIndexer_opts->detectorDistance <= 0.0){ + pinkIndexer_opts->detectorDistance = 0.25; //fake value + } + struct pinkIndexer_private_data* pinkIndexer_private_data = malloc(sizeof(struct pinkIndexer_private_data)); allocReciprocalPeaks(&(pinkIndexer_private_data->reciprocalPeaks_1_per_A)); pinkIndexer_private_data->intensities = malloc(MAX_PEAK_COUNT_FOR_INDEXER * sizeof(float)); -- cgit v1.2.3 From 117779c5fd9321b4a262f144798f081321a65a38 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Tue, 30 Apr 2019 14:43:06 +0200 Subject: Avoid duplication of beam and detector parameters for pinkIndexer --- libcrystfel/src/detector.c | 14 ++++++++++++++ libcrystfel/src/image.h | 2 ++ libcrystfel/src/index.c | 12 ++++++++---- libcrystfel/src/index.h | 3 ++- libcrystfel/src/pinkindexer.c | 30 +++++++++++++++--------------- libcrystfel/src/pinkindexer.h | 11 +++-------- src/indexamajig.c | 10 +++------- 7 files changed, 47 insertions(+), 35 deletions(-) diff --git a/libcrystfel/src/detector.c b/libcrystfel/src/detector.c index 201b9041..aab72341 100644 --- a/libcrystfel/src/detector.c +++ b/libcrystfel/src/detector.c @@ -1070,6 +1070,19 @@ static void parse_toplevel(struct detector *det, struct beam_params *beam, } } + } else if ( strcmp(key, "photon_energy_bandwidth") == 0 ) { + if ( beam != NULL ) { + double v; + char *end; + v = strtod(val, &end); + if ( (val[0] != '\0') && (end[0] == '\0') ) { + beam->bandwidth = v; + } else { + ERROR("Invalid value for " + "photon_energy_bandwidth\n"); + } + } + } else if ( strcmp(key, "photon_energy_scale") == 0 ) { if ( beam != NULL ) { beam->photon_energy_scale = atof(val); @@ -1202,6 +1215,7 @@ struct detector *get_detector_geometry_from_string(const char *string, beam->photon_energy = 0.0; beam->photon_energy_from = NULL; beam->photon_energy_scale = 1.0; + beam->bandwidth = 0.0; } det->n_panels = 0; diff --git a/libcrystfel/src/image.h b/libcrystfel/src/image.h index 706b06b1..04fb18a0 100644 --- a/libcrystfel/src/image.h +++ b/libcrystfel/src/image.h @@ -110,6 +110,8 @@ struct beam_params char *photon_energy_from; /**< HDF5 dataset name */ double photon_energy_scale; /**< Scale factor for photon energy, if it * comes from an image header */ + double bandwidth; /**< FWHM bandwidth as a fraction of + * wavelength */ }; diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 80a773ab..d48a6634 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -238,6 +238,7 @@ static char *friendly_indexer_name(IndexingMethod m) static void *prepare_method(IndexingMethod *m, UnitCell *cell, + struct detector *det, struct beam_params *beam, struct xgandalf_options *xgandalf_opts, struct pinkIndexer_options* pinkIndexer_opts, struct felix_options *felix_opts) @@ -285,7 +286,8 @@ static void *prepare_method(IndexingMethod *m, UnitCell *cell, break; case INDEXING_PINKINDEXER : - priv = pinkIndexer_prepare(m, cell, pinkIndexer_opts); + priv = pinkIndexer_prepare(m, cell, pinkIndexer_opts, + det, beam); break; default : @@ -316,8 +318,8 @@ static void *prepare_method(IndexingMethod *m, UnitCell *cell, IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, - struct detector *det, float *tols, - IndexingFlags flags, + struct detector *det, struct beam_params *beam, + float *tols, IndexingFlags flags, struct taketwo_options *ttopts, struct xgandalf_options *xgandalf_opts, struct pinkIndexer_options *pinkIndexer_opts, @@ -415,7 +417,9 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, int j; ipriv->engine_private[i] = prepare_method(&methods[i], cell, - xgandalf_opts, pinkIndexer_opts, + det, beam, + xgandalf_opts, + pinkIndexer_opts, felix_opts); if ( ipriv->engine_private[i] == NULL ) return NULL; diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h index ab38016f..883df4a3 100644 --- a/libcrystfel/src/index.h +++ b/libcrystfel/src/index.h @@ -153,7 +153,8 @@ extern IndexingMethod get_indm_from_string_2(const char *method, int *err); #include "felix.h" extern IndexingPrivate *setup_indexing(const char *methods, UnitCell *cell, - struct detector *det, float *ltl, + struct detector *det, + struct beam_params *beam, float *ltl, IndexingFlags flags, struct taketwo_options *ttopts, struct xgandalf_options *xgandalf_opts, diff --git a/libcrystfel/src/pinkindexer.c b/libcrystfel/src/pinkindexer.c index fe3b1625..eb8ec05e 100644 --- a/libcrystfel/src/pinkindexer.c +++ b/libcrystfel/src/pinkindexer.c @@ -154,26 +154,19 @@ int run_pinkIndexer(struct image *image, void *ipriv) } void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, - struct pinkIndexer_options *pinkIndexer_opts) + struct pinkIndexer_options *pinkIndexer_opts, + struct detector *det, struct beam_params *beam) { - if (pinkIndexer_opts->beamEnergy == 0.0) { + if ( beam->photon_energy_from != NULL ) { ERROR("For pinkIndexer, the photon_energy must be defined as a " "constant in the geometry file\n"); return NULL; } - if (pinkIndexer_opts->beamBandwidth == 0.0) { - STATUS("Using default bandwidth of 0.01 for pinkIndexer\n"); - pinkIndexer_opts->beamBandwidth = 0.01; - } - if (pinkIndexer_opts->detectorDistance == 0.0 && pinkIndexer_opts->refinement_type == + if ( (det->panels[0].clen_from != NULL) && pinkIndexer_opts->refinement_type == REFINEMENT_TYPE_firstFixedThenVariableLatticeParametersCenterAdjustmentMultiSeed) { ERROR("Using center refinement makes it necessary to have the detector distance fixed in the geometry file!"); } - if(pinkIndexer_opts->detectorDistance <= 0.0){ - pinkIndexer_opts->detectorDistance = 0.25; //fake value - } - struct pinkIndexer_private_data* pinkIndexer_private_data = malloc(sizeof(struct pinkIndexer_private_data)); allocReciprocalPeaks(&(pinkIndexer_private_data->reciprocalPeaks_1_per_A)); pinkIndexer_private_data->intensities = malloc(MAX_PEAK_COUNT_FOR_INDEXER * sizeof(float)); @@ -199,9 +192,15 @@ void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, .bx = bsz * 1e-10, .by = bsx * 1e-10, .bz = bsy * 1e-10, .cx = csz * 1e-10, .cy = csx * 1e-10, .cz = csy * 1e-10 }; - float detectorDistance_m = pinkIndexer_opts->detectorDistance; - float beamEenergy_eV = pinkIndexer_opts->beamEnergy; - float nonMonochromaticity = pinkIndexer_opts->beamBandwidth; + float detectorDistance_m; + if ( det->panels[0].clen_from != NULL ) { + detectorDistance_m = 0.25; /* fake value */ + } else { + detectorDistance_m = det->panels[0].clen + det->panels[0].coffset; + } + + float beamEenergy_eV = beam->photon_energy; + float nonMonochromaticity = beam->bandwidth; float reflectionRadius_1_per_A; if (pinkIndexer_opts->reflectionRadius < 0) { reflectionRadius_1_per_A = 0.02 @@ -350,7 +349,8 @@ int run_pinkIndexer(struct image *image, void *ipriv) } extern void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, - struct pinkIndexer_options *pinkIndexer_opts) + struct pinkIndexer_options *pinkIndexer_opts, + struct detector *det, struct beam_params *beam) { ERROR("This copy of CrystFEL was compiled without PINKINDEXER support.\n"); ERROR("To use PINKINDEXER indexing, recompile with PINKINDEXER.\n"); diff --git a/libcrystfel/src/pinkindexer.h b/libcrystfel/src/pinkindexer.h index 5be2c1e7..1fcb3711 100644 --- a/libcrystfel/src/pinkindexer.h +++ b/libcrystfel/src/pinkindexer.h @@ -42,14 +42,8 @@ struct pinkIndexer_options { int multi; int thread_count; int min_peaks; - int no_check_indexed; - - float beamEnergy; //in eV - float beamBandwidth; //(delta lambda)/lambda - float detectorDistance; //in m - - float reflectionRadius; //in 1/A + float reflectionRadius; /* In m^-1 */ }; #include @@ -58,7 +52,8 @@ struct pinkIndexer_options { extern int run_pinkIndexer(struct image *image, void *ipriv); extern void *pinkIndexer_prepare(IndexingMethod *indm, UnitCell *cell, - struct pinkIndexer_options *pinkIndexer_opts); + struct pinkIndexer_options *pinkIndexer_opts, + struct detector *det, struct beam_params *beam); extern void pinkIndexer_cleanup(void *pp); diff --git a/src/indexamajig.c b/src/indexamajig.c index 9aef3f2e..9fe14c65 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -1048,6 +1048,8 @@ int main(int argc, char *argv[]) "--pinkIndexer-reflection-radius\n"); return 1; } + /* A^-1 to m^-1 */ + iargs.pinkIndexer_opts.reflectionRadius /= 1e10; break; case 0 : @@ -1149,12 +1151,6 @@ int main(int argc, char *argv[]) return 1; } add_geom_beam_stuff_to_field_list(iargs.copyme, iargs.det, iargs.beam); - iargs.pinkIndexer_opts.beamEnergy = iargs.beam->photon_energy; - iargs.pinkIndexer_opts.beamBandwidth = 0.01; - iargs.pinkIndexer_opts.detectorDistance = iargs.det->panels[0].clen; - if(iargs.det->panels[0].clen_from != NULL){ - iargs.pinkIndexer_opts.detectorDistance = 0; - } /* If no peak path from geometry file, use these (but see later) */ if ( iargs.hdf5_peak_path == NULL ) { @@ -1394,7 +1390,7 @@ int main(int argc, char *argv[]) } iargs.ipriv = setup_indexing(indm_str, iargs.cell, iargs.det, - iargs.tols, flags, + iargs.beam, iargs.tols, flags, &iargs.taketwo_opts, &iargs.xgandalf_opts, &iargs.pinkIndexer_opts, -- cgit v1.2.3 From a249da0dae1ab6fa5d708f3ba4fa412c7c715ce1 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 2 Sep 2019 14:56:29 +0200 Subject: Remove message directing to old option --- libcrystfel/src/index.c | 1 - 1 file changed, 1 deletion(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index d48a6634..cc9f0248 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -347,7 +347,6 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, ERROR("The indexing method should contain only the method itself and "); ERROR("prior information modifiers ('cell' or 'latt')\n"); ERROR("To disable prediction refinement ('norefine'), use --no-refine.\n"); - ERROR("To check cell axes only ('axes'), use --no-cell-combinations.\n"); ERROR("To disable all unit cell checks ('raw'), use --no-check-cell.\n"); ERROR("To disable peak alignment check ('bad'), use --no-check-peaks.\n"); ERROR("To disable indexing retry ('noretry'), use --no-retry.\n"); -- cgit v1.2.3 From 6e713e3a78167d421acd7b511d07675ca0032d5a Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 2 Sep 2019 14:56:53 +0200 Subject: Add indexing_methods() to examine the actual list of methods --- libcrystfel/src/index.c | 6 ++++++ libcrystfel/src/index.h | 2 ++ 2 files changed, 8 insertions(+) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index cc9f0248..05d669a1 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -461,6 +461,12 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, } +const IndexingMethod *indexing_methods(IndexingPrivate *p, int *n) +{ + return p->methods; +} + + void cleanup_indexing(IndexingPrivate *ipriv) { int n; diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h index 883df4a3..4fc95b1d 100644 --- a/libcrystfel/src/index.h +++ b/libcrystfel/src/index.h @@ -161,6 +161,8 @@ extern IndexingPrivate *setup_indexing(const char *methods, UnitCell *cell, struct pinkIndexer_options *pinkIndexer_opts, struct felix_options *felix_opts); +extern const IndexingMethod *indexing_methods(IndexingPrivate *p, int *n); + extern char *detect_indexing_methods(UnitCell *cell); extern void index_pattern(struct image *image, IndexingPrivate *ipriv); -- cgit v1.2.3 From b0ceb28b75b7b0642cc40fa5dc4a48dcf8ec92ec Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 2 Sep 2019 14:57:16 +0200 Subject: Increase indexing timeout when using pinkIndexer --- libcrystfel/src/index.c | 1 + src/im-sandbox.c | 20 ++++++++++++++------ src/im-sandbox.h | 2 +- src/indexamajig.c | 14 +++++++++++++- 4 files changed, 29 insertions(+), 8 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 05d669a1..5a1673c3 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -463,6 +463,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, const IndexingMethod *indexing_methods(IndexingPrivate *p, int *n) { + *n = p->n_methods; return p->methods; } diff --git a/src/im-sandbox.c b/src/im-sandbox.c index dacc2dbb..d9348b37 100644 --- a/src/im-sandbox.c +++ b/src/im-sandbox.c @@ -73,6 +73,12 @@ struct sandbox int n_processed_last_stats; time_t t_last_stats; + /* Processing timeout in seconds. After this long without responding + * to a ping, the worker will be killed. After 3 times this long + * working on one image, even with ping responses, a warning will be + * shown to the user. */ + int timeout; + struct index_args *iargs; /* Worker processes */ @@ -157,17 +163,18 @@ static void check_hung_workers(struct sandbox *sb) stamp_response(sb, i); } - if ( tnow - sb->last_response[i] > 240 ) { - STATUS("Worker %i did not respond for 4 minutes - " - "sending it SIGKILL.\n", i); + if ( tnow - sb->last_response[i] > sb->timeout ) { + STATUS("Worker %i did not respond for %i seconds - " + "sending it SIGKILL.\n", i, sb->timeout); kill(sb->pids[i], SIGKILL); stamp_response(sb, i); } - if ( tnow - sb->shared->time_last_start[i] > 600 ) { + if ( tnow - sb->shared->time_last_start[i] > sb->timeout*3 ) { if ( !sb->shared->warned_long_running[i] ) { STATUS("Worker %i has been working on one " - "frame for more than 10 minutes.\n", i); + "frame for more than %i seconds (just " + "for info).\n", i, sb->timeout); STATUS("Event ID is: %s\n", sb->shared->last_ev[i]); STATUS("Task ID is: %s\n", @@ -1044,7 +1051,7 @@ char *create_tempdir(const char *temp_location) int create_sandbox(struct index_args *iargs, int n_proc, char *prefix, int config_basename, FILE *fh, Stream *stream, const char *tmpdir, int serial_start, - const char *zmq_address) + const char *zmq_address, int timeout) { int i; struct sandbox *sb; @@ -1073,6 +1080,7 @@ int create_sandbox(struct index_args *iargs, int n_proc, char *prefix, sb->iargs = iargs; sb->serial = serial_start; sb->tmpdir = tmpdir; + sb->timeout = timeout; if ( zmq_address != NULL ) { sb->zmq = 1; sb->zmq_address = zmq_address; diff --git a/src/im-sandbox.h b/src/im-sandbox.h index 9da11526..b1f5502a 100644 --- a/src/im-sandbox.h +++ b/src/im-sandbox.h @@ -84,6 +84,6 @@ extern void set_last_task(char *lt, const char *task); extern int create_sandbox(struct index_args *iargs, int n_proc, char *prefix, int config_basename, FILE *fh, Stream *stream, const char *tempdir, int serial_start, - const char *zmq_address); + const char *zmq_address, int timeout); #endif /* IM_SANDBOX_H */ diff --git a/src/indexamajig.c b/src/indexamajig.c index 9fe14c65..13a44471 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -313,6 +313,7 @@ int main(int argc, char *argv[]) char *spectrum_fn = NULL; int zmq = 0; char *zmq_address = NULL; + int timeout = 240; /* Defaults */ iargs.cell = NULL; @@ -1363,6 +1364,8 @@ int main(int argc, char *argv[]) } else { + int i, n; + const IndexingMethod *methods; IndexingFlags flags = 0; if ( iargs.cell != NULL ) { @@ -1400,6 +1403,15 @@ int main(int argc, char *argv[]) return 1; } + methods = indexing_methods(iargs.ipriv, &n); + for ( i=0; i Date: Tue, 3 Sep 2019 17:18:24 +0200 Subject: Force --no-check-peaks and --no-refine if pinkIndexer is used, and add a warning --- libcrystfel/src/index.c | 23 +++++++++++++++++++++++ 1 file changed, 23 insertions(+) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 5a1673c3..a3479a6c 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -423,6 +423,29 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, if ( ipriv->engine_private[i] == NULL ) return NULL; + if ( methods[i] & INDEXING_PINKINDEXER ) { + if ( n > 1 ) { + ERROR("WARNING: Using PinkIndexer at the same " + "time as other indexers is not " + "recommended.\n"); + } + + if ( flags & INDEXING_CHECK_PEAKS ) { + ERROR("WARNING: Setting --no-check-peaks " + "because PinkIndexer is in use.\n"); + } + flags |= INDEXING_CHECK_PEAKS; + flags ^= INDEXING_CHECK_PEAKS; + + if ( flags & INDEXING_REFINE ) { + ERROR("WARNING: Setting --no-refine because " + "PinkIndexer is in use.\n"); + } + flags |= INDEXING_REFINE; + flags ^= INDEXING_REFINE; + } + + for ( j=0; j Date: Tue, 10 Sep 2019 15:57:42 +0200 Subject: indexamajig: Get rid of --fix-bandwidth / use bandwidth from geometry file --- doc/man/indexamajig.1 | 3 +-- src/indexamajig.c | 12 ++++-------- src/process_image.c | 13 +++++-------- src/process_image.h | 1 - 4 files changed, 10 insertions(+), 19 deletions(-) diff --git a/doc/man/indexamajig.1 b/doc/man/indexamajig.1 index 3d96cf0b..9dfdd70d 100644 --- a/doc/man/indexamajig.1 +++ b/doc/man/indexamajig.1 @@ -482,10 +482,9 @@ Integrate the reflections using \fImethod\fR. See the section titled \fBPEAK IN .PD 0 .IP \fB--fix-profile-radius=\fIn\fR -.IP \fB--fix-bandwidth=\fIn\fR .IP \fB--fix-divergence=\fIn\fR .PD -Fix the beam and crystal paramters to the given values. The profile radius is given in m^-1, the bandwidth as a decimal fraction and the divergence in radians (full angle). The default is to set the divergence to zero, the bandwidth to a very small value, and then to automatically determine the profile radius. +Fix the beam and crystal paramters to the given values. The profile radius is given in m^-1 and the divergence in radians (full angle). The default is to set the divergence to zero, and then to automatically determine the profile radius. .IP You do not have to use all three of these options together. For example, if the automatic profile radius determination is not working well for your data set, you could fix that alone and continue using the default values for the other parameters (which might be automatically determined in future versions of CrystFEL, but are not currently). diff --git a/src/indexamajig.c b/src/indexamajig.c index 13a44471..0077c7fc 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -232,8 +232,7 @@ static void show_help(const char *s) " --integration= Integration method (rings,prof2d)-(cen,nocen)\n" " Default: rings-nocen\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" +" prediction (default: automatically determine)\n" " --fix-divergence Set the divergence (full angle) for spot prediction\n" " --int-radius= Set the integration radii. Default: 4,5,7.\n" " --int-diag= Show debugging information about reflections\n" @@ -370,7 +369,6 @@ int main(int argc, char *argv[]) iargs.push_res = 0.0; iargs.highres = +INFINITY; iargs.fix_profile_r = -1.0; - iargs.fix_bandwidth = -1.0; iargs.fix_divergence = -1.0; iargs.profile = 0; iargs.no_image_data = 0; @@ -699,11 +697,9 @@ int main(int argc, char *argv[]) break; case 323 : - if ( sscanf(optarg, "%f", &iargs.fix_bandwidth) != 1 ) { - ERROR("Invalid value for --fix-bandwidth\n"); - return 1; - } - break; + ERROR("Option --fix-bandwidth is no longer used.\n"); + ERROR("Set the bandwidth in the geometry file instead.\n"); + return 1; case 324 : if ( sscanf(optarg, "%f", &iargs.fix_divergence) != 1 ) { diff --git a/src/process_image.c b/src/process_image.c index 654ffd46..f28952bc 100644 --- a/src/process_image.c +++ b/src/process_image.c @@ -226,6 +226,8 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs, &imfile) ) return; } + image.bw = iargs->beam->bandwidth; + /* Take snapshot of image before applying horrible noise filters */ time_accounts_set(taccs, TACC_FILTER); set_last_task(last_task, "image filter"); @@ -369,17 +371,12 @@ void process_image(const struct index_args *iargs, struct pattern_args *pargs, } else { image.div = 0.0; } - if ( iargs->fix_bandwidth >= 0.0 ) { - image.bw = iargs->fix_bandwidth; - } else { - image.bw = 0.00000001; - } - - /* Set beam spectrum for pink beam data */ + /* Set beam spectrum */ if ( iargs->spectrum != NULL ) { image.spectrum = iargs->spectrum; } else { - image.spectrum = spectrum_generate_gaussian(image.lambda, image.bw); + image.spectrum = spectrum_generate_gaussian(image.lambda, + image.bw); } if ( image_feature_count(image.features) < iargs->min_peaks ) { diff --git a/src/process_image.h b/src/process_image.h index 5691a0e4..2e03900e 100644 --- a/src/process_image.h +++ b/src/process_image.h @@ -112,7 +112,6 @@ struct index_args float push_res; float highres; float fix_profile_r; - float fix_bandwidth; float fix_divergence; int overpredict; int profile; /* Whether or not to do wall clock profiling */ -- cgit v1.2.3 From 8d8ce4a3c7b340c1ce7d538c29fd08c59190d288 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 11 Sep 2019 11:12:48 +0200 Subject: Fix logic --- libcrystfel/src/index.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index a3479a6c..440a5b97 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -423,7 +423,7 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, if ( ipriv->engine_private[i] == NULL ) return NULL; - if ( methods[i] & INDEXING_PINKINDEXER ) { + if ( (methods[i] & INDEXING_METHOD_MASK) == INDEXING_PINKINDEXER ) { if ( n > 1 ) { ERROR("WARNING: Using PinkIndexer at the same " "time as other indexers is not " -- cgit v1.2.3 From 8695e135823b6b6708713761dd41a0871551522c Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 11 Sep 2019 11:28:21 +0200 Subject: Restore default bandwidth (if not specified in geometry file) --- libcrystfel/src/detector.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/libcrystfel/src/detector.c b/libcrystfel/src/detector.c index aab72341..4d008e65 100644 --- a/libcrystfel/src/detector.c +++ b/libcrystfel/src/detector.c @@ -1215,7 +1215,7 @@ struct detector *get_detector_geometry_from_string(const char *string, beam->photon_energy = 0.0; beam->photon_energy_from = NULL; beam->photon_energy_scale = 1.0; - beam->bandwidth = 0.0; + beam->bandwidth = 0.00000001; } det->n_panels = 0; -- cgit v1.2.3 From aad71902c94b3bc21ab34b0e2c92972400921574 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 11 Sep 2019 17:05:33 +0200 Subject: Fix sign of zl comparison --- libcrystfel/src/geometry.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/libcrystfel/src/geometry.c b/libcrystfel/src/geometry.c index 5b2107d1..fc42a8de 100644 --- a/libcrystfel/src/geometry.c +++ b/libcrystfel/src/geometry.c @@ -332,7 +332,7 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst, /* Closest point on Ewald sphere. * Project zl to 0, bit of a hack... */ - const double zlp0 = zl>0?zl:0; + const double zlp0 = zl<0?zl:0; exerr2 = (x-xl)*(x-xl) + (y-yl)*(y-yl) + (z-zl)*(z-zl); /* Weighted average between projected lattice point -- cgit v1.2.3