From 141fb2cd3af6ea14098cc65ba3c64f93f72c9cf3 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 6 Mar 2017 17:14:54 +0100 Subject: Add compare_cells() (and use it in whirligig) --- libcrystfel/src/cell-utils.c | 121 ++++++++++++++++++++++++++++++++++++++++++- libcrystfel/src/cell-utils.h | 9 ++-- src/whirligig.c | 111 +++------------------------------------ 3 files changed, 131 insertions(+), 110 deletions(-) diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c index 6610abc3..41436bb4 100644 --- a/libcrystfel/src/cell-utils.c +++ b/libcrystfel/src/cell-utils.c @@ -3,12 +3,12 @@ * * Unit Cell utility functions * - * Copyright © 2012-2015 Deutsches Elektronen-Synchrotron DESY, + * Copyright © 2012-2017 Deutsches Elektronen-Synchrotron DESY, * a research centre of the Helmholtz Association. * Copyright © 2012 Lorenzo Galli * * Authors: - * 2009-2012,2014-2015 Thomas White + * 2009-2012,2014-2017 Thomas White * 2012 Lorenzo Galli * * This file is part of CrystFEL. @@ -1686,3 +1686,120 @@ double cell_get_volume(UnitCell *cell) return 1/rec_volume; } + + +static double moduli_check(double ax, double ay, double az, + double bx, double by, double bz) +{ + double ma = modulus(ax, ay, az); + double mb = modulus(bx, by, bz); + return fabs(ma-mb)/ma; +} + + +static int cells_are_similar(UnitCell *cell1, UnitCell *cell2, + const double ltl, const double atl) +{ + double asx1, asy1, asz1, bsx1, bsy1, bsz1, csx1, csy1, csz1; + double asx2, asy2, asz2, bsx2, bsy2, bsz2, csx2, csy2, csz2; + UnitCell *pcell1, *pcell2; + + /* Compare primitive cells, not centered */ + pcell1 = uncenter_cell(cell1, NULL); + pcell2 = uncenter_cell(cell2, NULL); + + cell_get_reciprocal(pcell1, &asx1, &asy1, &asz1, + &bsx1, &bsy1, &bsz1, + &csx1, &csy1, &csz1); + + cell_get_reciprocal(pcell2, &asx2, &asy2, &asz2, + &bsx2, &bsy2, &bsz2, + &csx2, &csy2, &csz2); + + + cell_free(pcell1); + cell_free(pcell2); + + if ( angle_between(asx1, asy1, asz1, asx2, asy2, asz2) > atl ) return 0; + if ( angle_between(bsx1, bsy1, bsz1, bsx2, bsy2, bsz2) > atl ) return 0; + if ( angle_between(csx1, csy1, csz1, csx2, csy2, csz2) > atl ) return 0; + + if ( moduli_check(asx1, asy1, asz1, asx2, asy2, asz2) > ltl ) return 0; + if ( moduli_check(bsx1, bsy1, bsz1, bsx2, bsy2, bsz2) > ltl ) return 0; + if ( moduli_check(csx1, csy1, csz1, csx2, csy2, csz2) > ltl ) return 0; + + return 1; +} + + + +/** + * compare_cells: + * @a: A %UnitCell + * @b: Another %UnitCell + * @ltl: Maximum allowable fractional difference in reciprocal axis lengths + * @atl: Maximum allowable difference in reciprocal angles (in radians) + * @pmb: Place to store pointer to matrix + * + * Compare the two units cells. If they agree to within @ltl and @atl, using + * any change of axes, returns non-zero and stores the transformation to map @b + * onto @a. + * + * Returns: non-zero if the cells match. + * + */ +int compare_cells(UnitCell *a, UnitCell *b, double ltl, double atl, + IntegerMatrix **pmb) +{ + IntegerMatrix *m; + int i[9]; + + m = intmat_new(3, 3); + + for ( i[0]=-1; i[0]<=+1; i[0]++ ) { + for ( i[1]=-1; i[1]<=+1; i[1]++ ) { + for ( i[2]=-1; i[2]<=+1; i[2]++ ) { + for ( i[3]=-1; i[3]<=+1; i[3]++ ) { + for ( i[4]=-1; i[4]<=+1; i[4]++ ) { + for ( i[5]=-1; i[5]<=+1; i[5]++ ) { + for ( i[6]=-1; i[6]<=+1; i[6]++ ) { + for ( i[7]=-1; i[7]<=+1; i[7]++ ) { + for ( i[8]=-1; i[8]<=+1; i[8]++ ) { + + UnitCellTransformation *tfn; + UnitCell *nc; + int j, k; + int l = 0; + + for ( j=0; j<3; j++ ) + for ( k=0; k<3; k++ ) + intmat_set(m, j, k, i[l++]); + + if ( intmat_det(m) != +1 ) continue; + + tfn = tfn_from_intmat(m); + nc = cell_transform(b, tfn); + + if ( cells_are_similar(a, nc, ltl, atl) ) { + if ( pmb != NULL ) *pmb = m; + tfn_free(tfn); + cell_free(nc); + return 1; + } + + tfn_free(tfn); + cell_free(nc); + + } + } + } + } + } + } + } + } + } + + intmat_free(m); + return 0; +} diff --git a/libcrystfel/src/cell-utils.h b/libcrystfel/src/cell-utils.h index efb4b25b..5e2b2825 100644 --- a/libcrystfel/src/cell-utils.h +++ b/libcrystfel/src/cell-utils.h @@ -3,13 +3,13 @@ * * Unit Cell utility functions * - * Copyright © 2012-2014 Deutsches Elektronen-Synchrotron DESY, + * Copyright © 2012-2017 Deutsches Elektronen-Synchrotron DESY, * a research centre of the Helmholtz Association. * Copyright © 2012 Lorenzo Galli * * Authors: - * 2009-2013,2014 Thomas White - * 2012 Lorenzo Galli + * 2009-2013,2014,2017 Thomas White + * 2012 Lorenzo Galli * * This file is part of CrystFEL. * @@ -80,6 +80,9 @@ extern int forbidden_reflection(UnitCell *cell, extern double cell_get_volume(UnitCell *cell); +extern int compare_cells(UnitCell *a, UnitCell *b, double ltl, double atl, + IntegerMatrix **pmb); + #ifdef __cplusplus } #endif diff --git a/src/whirligig.c b/src/whirligig.c index c5e8084b..e8829880 100644 --- a/src/whirligig.c +++ b/src/whirligig.c @@ -279,107 +279,6 @@ static void find_and_process_series(struct window *win, int is_last_frame, } -static double moduli_check(double ax, double ay, double az, - double bx, double by, double bz) -{ - double ma = modulus(ax, ay, az); - double mb = modulus(bx, by, bz); - return fabs(ma-mb)/ma; -} - - -static int cells_are_similar(UnitCell *cell1, UnitCell *cell2) -{ - double asx1, asy1, asz1, bsx1, bsy1, bsz1, csx1, csy1, csz1; - double asx2, asy2, asz2, bsx2, bsy2, bsz2, csx2, csy2, csz2; - UnitCell *pcell1, *pcell2; - const double atl = deg2rad(5.0); - const double ltl = 0.1; - - /* Compare primitive cells, not centered */ - pcell1 = uncenter_cell(cell1, NULL); - pcell2 = uncenter_cell(cell2, NULL); - - cell_get_reciprocal(pcell1, &asx1, &asy1, &asz1, - &bsx1, &bsy1, &bsz1, - &csx1, &csy1, &csz1); - - cell_get_reciprocal(pcell2, &asx2, &asy2, &asz2, - &bsx2, &bsy2, &bsz2, - &csx2, &csy2, &csz2); - - - cell_free(pcell1); - cell_free(pcell2); - - if ( angle_between(asx1, asy1, asz1, asx2, asy2, asz2) > atl ) return 0; - if ( angle_between(bsx1, bsy1, bsz1, bsx2, bsy2, bsz2) > atl ) return 0; - if ( angle_between(csx1, csy1, csz1, csx2, csy2, csz2) > atl ) return 0; - - if ( moduli_check(asx1, asy1, asz1, asx2, asy2, asz2) > ltl ) return 0; - if ( moduli_check(bsx1, bsy1, bsz1, bsx2, bsy2, bsz2) > ltl ) return 0; - if ( moduli_check(csx1, csy1, csz1, csx2, csy2, csz2) > ltl ) return 0; - - return 1; -} - - -static int gatinator(UnitCell *a, UnitCell *b, IntegerMatrix **pmb) -{ - IntegerMatrix *m; - int i[9]; - - m = intmat_new(3, 3); - - for ( i[0]=-1; i[0]<=+1; i[0]++ ) { - for ( i[1]=-1; i[1]<=+1; i[1]++ ) { - for ( i[2]=-1; i[2]<=+1; i[2]++ ) { - for ( i[3]=-1; i[3]<=+1; i[3]++ ) { - for ( i[4]=-1; i[4]<=+1; i[4]++ ) { - for ( i[5]=-1; i[5]<=+1; i[5]++ ) { - for ( i[6]=-1; i[6]<=+1; i[6]++ ) { - for ( i[7]=-1; i[7]<=+1; i[7]++ ) { - for ( i[8]=-1; i[8]<=+1; i[8]++ ) { - - UnitCellTransformation *tfn; - UnitCell *nc; - int j, k; - int l = 0; - - for ( j=0; j<3; j++ ) - for ( k=0; k<3; k++ ) - intmat_set(m, j, k, i[l++]); - - if ( intmat_det(m) != +1 ) continue; - - tfn = tfn_from_intmat(m); - nc = cell_transform(b, tfn); - - if ( cells_are_similar(a, nc) ) { - *pmb = m; - tfn_free(tfn); - cell_free(nc); - return 1; - } - - tfn_free(tfn); - cell_free(nc); - - } - } - } - } - } - } - } - } - } - - intmat_free(m); - return 0; -} - - static int crystal_used(struct window *win, int pos, int cn) { int i; @@ -411,8 +310,9 @@ static IntegerMatrix *try_all(struct window *win, int n1, int n2, for ( i=0; in_crystals; i++ ) { for ( j=0; jn_crystals; j++ ) { - if ( gatinator(crystal_get_cell(i1->crystals[i]), - crystal_get_cell(i2->crystals[j]), &m) ) + if ( compare_cells(crystal_get_cell(i1->crystals[i]), + crystal_get_cell(i2->crystals[j]), + 0.1, deg2rad(5.0), &m) ) { if ( !crystal_used(win, n1, i) && !crystal_used(win, n2, j) ) @@ -478,8 +378,9 @@ static int try_join(struct window *win, int sn) for ( j=0; jimg[win->join_ptr].n_crystals; j++ ) { Crystal *cr2; cr2 = win->img[win->join_ptr].crystals[j]; - if ( gatinator(ref, crystal_get_cell(cr2), - &win->mat[sn][win->join_ptr]) ) { + if ( compare_cells(ref, crystal_get_cell(cr2), + 0.1, deg2rad(5.0), + &win->mat[sn][win->join_ptr]) ) { win->ser[sn][win->join_ptr] = j; cell_free(ref); return 1; -- cgit v1.2.3 From 2ea5c535ac15aa0e429eb985a82f326cd2bce7c2 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 6 Mar 2017 17:15:49 +0100 Subject: remove_flagged_crystals(): Return the number of crystals removed --- libcrystfel/src/image.c | 5 ++++- libcrystfel/src/image.h | 2 +- 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/libcrystfel/src/image.c b/libcrystfel/src/image.c index aad5017c..8de16bd5 100644 --- a/libcrystfel/src/image.c +++ b/libcrystfel/src/image.c @@ -288,9 +288,10 @@ void image_add_crystal(struct image *image, Crystal *cryst) } -void remove_flagged_crystals(struct image *image) +int remove_flagged_crystals(struct image *image) { int i; + int n_bad = 0; for ( i=0; in_crystals; i++ ) { if ( crystal_get_user_flag(image->crystals[i]) ) { @@ -302,9 +303,11 @@ void remove_flagged_crystals(struct image *image) image->crystals[j] = image->crystals[j+1]; } image->n_crystals--; + n_bad++; } } + return n_bad; } diff --git a/libcrystfel/src/image.h b/libcrystfel/src/image.h index d7884e4c..c900bd29 100644 --- a/libcrystfel/src/image.h +++ b/libcrystfel/src/image.h @@ -257,7 +257,7 @@ extern struct imagefeature *image_get_feature(ImageFeatureList *flist, int idx); extern ImageFeatureList *sort_peaks(ImageFeatureList *flist); extern void image_add_crystal(struct image *image, Crystal *cryst); -extern void remove_flagged_crystals(struct image *image); +extern int remove_flagged_crystals(struct image *image); extern void free_all_crystals(struct image *image); /* Image files */ -- cgit v1.2.3 From 71659b553db6b4edda325bbfd2b39a59fd0fff0f Mon Sep 17 00:00:00 2001 From: Thomas White Date: Tue, 7 Mar 2017 13:03:06 +0100 Subject: Use compare_cells() for multi-lattice indexing --- libcrystfel/src/index.c | 42 ++++++++++++++++++++++++++++++++---------- 1 file changed, 32 insertions(+), 10 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 644ea31f..36f99995 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -354,9 +354,15 @@ static int try_indexer(struct image *image, IndexingMethod indm, } + /* For all the crystals found this time ... */ for ( i=0; icrystals[image->n_crystals-i-1]; + int j; + int this_crystal = image->n_crystals - i - 1; + + /* ... starting at the end of the (complete) list ... */ + Crystal *cr = image->crystals[this_crystal]; + crystal_set_image(cr, image); crystal_set_user_flag(cr, 0); crystal_set_profile_radius(cr, 0.02e9); @@ -385,7 +391,6 @@ static int try_indexer(struct image *image, IndexingMethod indm, if ( out == NULL ) { crystal_set_user_flag(cr, 1); - n_bad++; } cell_free(out); @@ -394,13 +399,32 @@ static int try_indexer(struct image *image, IndexingMethod indm, } - } + /* Don't do similarity check if this crystal is bad */ + if ( crystal_get_user_flag(cr) ) continue; + + /* Check if cell is too similar to existing ones */ + for ( j=0; jcrystals[j]; + + /* Don't do similarity check against bad crystals */ + if ( crystal_get_user_flag(that_cr) ) continue; - remove_flagged_crystals(image); + if ( compare_cells(crystal_get_cell(cr), + crystal_get_cell(that_cr), + 0.1, deg2rad(5.0), NULL) ) + { + crystal_set_user_flag(cr, 1); + ERROR("Excluding duplicate cell\n"); + } + } - if ( n_bad == r ) return 0; + } + + n_bad = remove_flagged_crystals(image); + assert(r >= n_bad); - return r; + return r - n_bad; } @@ -511,12 +535,9 @@ static int finished_retry(IndexingMethod indm, int r, struct image *image) * after deleting the explained peaks */ if ( indm & INDEXING_MULTI ) { - /* Remove "used" spots and try for - * another lattice */ + /* Remove "used" spots and try for another lattice */ Crystal *cr; cr = image->crystals[image->n_crystals-1]; - STATUS("WARNING: Multi-lattice indexing does not work" - " well in this version.\n"); return delete_explained_peaks(image, cr); } else { return 1; @@ -525,6 +546,7 @@ static int finished_retry(IndexingMethod indm, int r, struct image *image) } } + void index_pattern(struct image *image, IndexingPrivate *ipriv) { index_pattern_2(image, ipriv, NULL); -- cgit v1.2.3 From 3f1c6a1fbc241721aafeb9e1e6c0cbdb1eb810b5 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Thu, 7 Sep 2017 16:49:42 +0200 Subject: Fix uninitialised peak counts (for info in stream only) --- libcrystfel/src/hdf5-file.c | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/libcrystfel/src/hdf5-file.c b/libcrystfel/src/hdf5-file.c index 13339506..6fc10357 100644 --- a/libcrystfel/src/hdf5-file.c +++ b/libcrystfel/src/hdf5-file.c @@ -409,6 +409,8 @@ int get_peaks_cxi_2(struct image *image, struct hdfile *f, const char *p, } image->features = image_feature_list_new(); + image->num_peaks = 0; + image->num_saturated_peaks = 0; for ( pk=0; pkorig_min_ss; image_add_feature(image->features, fs, ss, p, image, val, NULL); + image->num_peaks++; } @@ -554,6 +557,8 @@ int get_peaks_2(struct image *image, struct hdfile *f, const char *p, } image->features = image_feature_list_new(); + image->num_peaks = 0; + image->num_saturated_peaks = 0; for ( i=0; ifeatures, fs, ss, p, image, val, NULL); + image->num_peaks++; } -- cgit v1.2.3 From 60c48f3876b0d3c30bf729ce691dbd3f56665c4c Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 8 Sep 2017 11:41:51 +0200 Subject: Simplify indexing options mosflm-cell-nolatt-multi-refine-bad-retry.... it's all getting a bit too much. This commit removes all the options from the indexing methods except for the prior information flags (cell/nocell and latt/nolatt). All the other options are now given at the top level with new indexamajig options, e.g. --no-multi and --no-retry. They are applied to all the methods, with no option to (say) use retry for one method and not for another. This also makes the "peak alignment test" default to OFF, and multi-lattice indexing ON. I think that the prediction refinement does a good enough job of checking the validity of indexing solutions, allowing us to turn off the alignment check which filters out most multi-lattice solutions. --- libcrystfel/src/asdf.c | 42 +--------- libcrystfel/src/dirax.c | 31 +------- libcrystfel/src/felix.c | 16 +--- libcrystfel/src/index.c | 198 +++++++++++++++++++--------------------------- libcrystfel/src/index.h | 77 +++++++----------- libcrystfel/src/mosflm.c | 63 +-------------- libcrystfel/src/taketwo.c | 19 +---- libcrystfel/src/xds.c | 52 +++--------- src/indexamajig.c | 38 ++++++++- 9 files changed, 168 insertions(+), 368 deletions(-) diff --git a/libcrystfel/src/asdf.c b/libcrystfel/src/asdf.c index 0c43fe7a..bbfd1a56 100644 --- a/libcrystfel/src/asdf.c +++ b/libcrystfel/src/asdf.c @@ -272,39 +272,14 @@ static int calc_reciprocal(gsl_vector **direct, gsl_vector **reciprocal) static int check_cell(struct asdf_private *dp, struct image *image, UnitCell *cell) { - UnitCell *out; Crystal *cr; - if ( dp->indm & INDEXING_CHECK_CELL_COMBINATIONS ) { - - out = match_cell(cell, dp->template, 0, dp->ltl, 1); - if ( out == NULL ) return 0; - - } else if ( dp->indm & INDEXING_CHECK_CELL_AXES ) { - - out = match_cell(cell, dp->template, 0, dp->ltl, 0); - if ( out == NULL ) return 0; - - } else { - out = cell_new_from_cell(cell); - } - cr = crystal_new(); if ( cr == NULL ) { ERROR("Failed to allocate crystal.\n"); return 0; } - - crystal_set_cell(cr, out); - - if ( dp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, &cr, 1) ) { - crystal_free(cr); /* Frees the cell as well */ - cell_free(out); - return 0; - } - } - + crystal_set_cell(cr, cell); image_add_crystal(image, cr); return 1; @@ -1206,22 +1181,9 @@ void *asdf_prepare(IndexingMethod *indm, UnitCell *cell, struct detector *det, float *ltl) { struct asdf_private *dp; - int need_cell = 0; - - if ( *indm & INDEXING_CHECK_CELL_COMBINATIONS ) need_cell = 1; - if ( *indm & INDEXING_CHECK_CELL_AXES ) need_cell = 1; - - if ( need_cell && !cell_has_parameters(cell) ) { - ERROR("Altering your asdf flags because cell parameters were" - " not provided.\n"); - *indm &= ~INDEXING_CHECK_CELL_COMBINATIONS; - *indm &= ~INDEXING_CHECK_CELL_AXES; - } /* Flags that asdf knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_CELL_COMBINATIONS - | INDEXING_CHECK_CELL_AXES | INDEXING_CHECK_PEAKS - | INDEXING_CONTROL_FLAGS; + *indm &= INDEXING_METHOD_MASK; dp = malloc(sizeof(struct asdf_private)); if ( dp == NULL ) return NULL; diff --git a/libcrystfel/src/dirax.c b/libcrystfel/src/dirax.c index b263eea0..29b735ad 100644 --- a/libcrystfel/src/dirax.c +++ b/libcrystfel/src/dirax.c @@ -108,39 +108,14 @@ struct dirax_data { static int check_cell(struct dirax_private *dp, struct image *image, UnitCell *cell) { - UnitCell *out; Crystal *cr; - if ( dp->indm & INDEXING_CHECK_CELL_COMBINATIONS ) { - - out = match_cell(cell, dp->template, 0, dp->ltl, 1); - if ( out == NULL ) return 0; - - } else if ( dp->indm & INDEXING_CHECK_CELL_AXES ) { - - out = match_cell(cell, dp->template, 0, dp->ltl, 0); - if ( out == NULL ) return 0; - - } else { - out = cell_new_from_cell(cell); - } - cr = crystal_new(); if ( cr == NULL ) { ERROR("Failed to allocate crystal.\n"); return 0; } - - crystal_set_cell(cr, out); - - if ( dp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, &cr, 1) ) { - crystal_free(cr); /* Frees the cell as well */ - cell_free(out); - return 0; - } - } - + crystal_set_cell(cr, cell); image_add_crystal(image, cr); return 1; @@ -655,9 +630,7 @@ void *dirax_prepare(IndexingMethod *indm, UnitCell *cell, } /* Flags that DirAx knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_CELL_COMBINATIONS - | INDEXING_CHECK_CELL_AXES | INDEXING_CHECK_PEAKS - | INDEXING_CONTROL_FLAGS; + *indm &= INDEXING_METHOD_MASK; dp = malloc(sizeof(struct dirax_private)); if ( dp == NULL ) return NULL; diff --git a/libcrystfel/src/felix.c b/libcrystfel/src/felix.c index 6daf4408..4c12778b 100644 --- a/libcrystfel/src/felix.c +++ b/libcrystfel/src/felix.c @@ -222,16 +222,7 @@ static int read_felix(struct felix_private *gp, struct image *image, fclose(fh); - if ( gp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, image->crystals, - image->n_crystals) ) - { - free_all_crystals(image); - return 0; - } - } - - return n_crystals; + return n_crystals; } @@ -672,9 +663,8 @@ void *felix_prepare(IndexingMethod *indm, UnitCell *cell, if ( gp == NULL ) return NULL; /* Flags that Felix knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_PEAKS - | INDEXING_USE_LATTICE_TYPE | INDEXING_USE_CELL_PARAMETERS - | INDEXING_CONTROL_FLAGS; + *indm &= INDEXING_METHOD_MASK + | INDEXING_USE_LATTICE_TYPE | INDEXING_USE_CELL_PARAMETERS; gp->cell = cell; gp->indm = *indm; diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 36f99995..e7ae7a37 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -61,6 +61,7 @@ struct _indexingprivate { + IndexingFlags flags; UnitCell *target_cell; float tolerance[4]; @@ -72,6 +73,30 @@ struct _indexingprivate }; +static const char *onoff(int a) +{ + if ( a ) return "on"; + return "off"; +} + + +static void show_indexing_flags(IndexingFlags flags) +{ + STATUS(" Check unit cell (combinations): %s\n", + onoff(flags & INDEXING_CHECK_CELL_COMBINATIONS)); + STATUS(" Check unit cell (axis permutations only): %s\n", + onoff(flags & INDEXING_CHECK_CELL_AXES)); + STATUS(" Check peak alignment: %s\n", + onoff(flags & INDEXING_CHECK_PEAKS)); + STATUS(" Refine indexing solutions: %s\n", + onoff(flags & INDEXING_REFINE)); + STATUS(" Multi-lattice indexing (\"delete and retry\"): %s\n", + onoff(flags & INDEXING_MULTI)); + STATUS(" Retry indexing: %s\n", + onoff(flags & INDEXING_RETRY)); +} + + static int debug_index(struct image *image) { Crystal *cr = crystal_new(); @@ -159,7 +184,7 @@ static void *prepare_method(IndexingMethod *m, UnitCell *cell, IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, struct detector *det, float *ltl, - int no_refine, const char *options, + IndexingFlags flags, const char *options, struct taketwo_options *ttopts) { int i, n; @@ -178,10 +203,30 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, for ( i=0; imethods = methods; ipriv->n_methods = n; + ipriv->flags = flags; if ( cell != NULL ) { ipriv->target_cell = cell_new_from_cell(cell); @@ -222,6 +268,8 @@ IndexingPrivate *setup_indexing(const char *method_list, UnitCell *cell, ipriv->ttopts = ttopts; + show_indexing_flags(flags); + return ipriv; } @@ -369,34 +417,43 @@ static int try_indexer(struct image *image, IndexingMethod indm, crystal_set_mosaicity(cr, 0.0); /* Prediction refinement if requested */ - if ( indm & INDEXING_REFINE ) { - - UnitCell *out; - + if ( ipriv->flags & INDEXING_REFINE ) + { if ( refine_prediction(image, cr) ) { crystal_set_user_flag(cr, 1); n_bad++; continue; } + } - if ( (indm & INDEXING_CHECK_CELL_COMBINATIONS) - || (indm & INDEXING_CHECK_CELL_AXES) ) - { - - /* Check that the cell parameters are still - * within the tolerance */ - out = match_cell(crystal_get_cell(cr), - ipriv->target_cell, 0, - ipriv->tolerance, 0); + /* Unit cell check if requested */ + assert( !((ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS) + && (ipriv->flags & INDEXING_CHECK_CELL_AXES)) ); + if ( (ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS) + || (ipriv->flags & INDEXING_CHECK_CELL_AXES) ) + { + UnitCell *out; + int reduce; - if ( out == NULL ) { - crystal_set_user_flag(cr, 1); - } + if ( ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS ) + { + reduce = 1; + } else { + reduce = 0; + } - cell_free(out); + out = match_cell(crystal_get_cell(cr), + ipriv->target_cell, 0, + ipriv->tolerance, reduce); + if ( out == NULL ) { + crystal_set_user_flag(cr, 1); + n_bad++; + continue; } + cell_free(crystal_get_cell(cr)); + crystal_set_cell(cr, out); } /* Don't do similarity check if this crystal is bad */ @@ -514,14 +571,15 @@ static int delete_explained_peaks(struct image *image, Crystal *cr) * * Returns false for "try again", true for "no, stop now" */ -static int finished_retry(IndexingMethod indm, int r, struct image *image) +static int finished_retry(IndexingMethod indm, IndexingFlags flags, + int r, struct image *image) { if ( r == 0 ) { /* Indexing failed on the previous attempt. Maybe try again * after poking the peak list a bit */ - if ( indm & INDEXING_RETRY ) { + if ( flags & INDEXING_RETRY ) { /* Retry with fewer peaks */ return delete_weakest_peaks(image->features); } else { @@ -534,7 +592,7 @@ static int finished_retry(IndexingMethod indm, int r, struct image *image) /* Indexing succeeded on previous attempt. Maybe try again * after deleting the explained peaks */ - if ( indm & INDEXING_MULTI ) { + if ( flags & INDEXING_MULTI ) { /* Remove "used" spots and try for another lattice */ Crystal *cr; cr = image->crystals[image->n_crystals-1]; @@ -581,7 +639,8 @@ void index_pattern_2(struct image *image, IndexingPrivate *ipriv, int *ping) ipriv, ipriv->engine_private[n]); success += r; ntry++; - done = finished_retry(ipriv->methods[n], r, image); + done = finished_retry(ipriv->methods[n], ipriv->flags, + r, image); if ( ntry > 5 ) done = 1; if ( ping != NULL ) (*ping)++; @@ -605,39 +664,6 @@ void index_pattern_2(struct image *image, IndexingPrivate *ipriv, int *ping) } -/* Set the indexer flags for "raw mode" ("--cell-reduction=none") */ -static IndexingMethod set_raw(IndexingMethod a) -{ - /* Disable all unit cell checks */ - a &= ~(INDEXING_CHECK_CELL_COMBINATIONS | INDEXING_CHECK_CELL_AXES); - return a; -} - - -/* Set the indexer flags for "bad mode" ("--insane) */ -static IndexingMethod set_bad(IndexingMethod a) -{ - /* Disable the peak check */ - return a & ~INDEXING_CHECK_PEAKS; -} - - -/* Set the indexer flags for "axes mode" ("--cell-reduction=compare") */ -static IndexingMethod set_axes(IndexingMethod a) -{ - return (a & ~INDEXING_CHECK_CELL_COMBINATIONS) - | INDEXING_CHECK_CELL_AXES; -} - - -/* Set the indexer flags for "combination mode" ("--cell-reduction=reduce") */ -static IndexingMethod set_comb(IndexingMethod a) -{ - return (a & ~INDEXING_CHECK_CELL_AXES) - | INDEXING_CHECK_CELL_COMBINATIONS; -} - - /* Set the indexer flags for "use no lattice type information" */ static IndexingMethod set_nolattice(IndexingMethod a) { @@ -725,18 +751,6 @@ char *indexer_str(IndexingMethod indm) if ( (indm & INDEXING_METHOD_MASK) == INDEXING_SIMULATION ) return str; - if ( indm & INDEXING_CHECK_CELL_COMBINATIONS ) { - strcat(str, "-comb"); - } else if ( indm & INDEXING_CHECK_CELL_AXES ) { - strcat(str, "-axes"); - } else { - strcat(str, "-raw"); - } - - if ( !(indm & INDEXING_CHECK_PEAKS) ) { - strcat(str, "-bad"); - } - if ( indm & INDEXING_USE_LATTICE_TYPE ) { strcat(str, "-latt"); } else { @@ -749,24 +763,6 @@ char *indexer_str(IndexingMethod indm) strcat(str, "-nocell"); } - if ( indm & INDEXING_RETRY ) { - strcat(str, "-retry"); - } else { - strcat(str, "-noretry"); - } - - if ( indm & INDEXING_MULTI ) { - strcat(str, "-multi"); - } else { - strcat(str, "-nomulti"); - } - - if ( indm & INDEXING_REFINE ) { - strcat(str, "-refine"); - } else { - strcat(str, "-norefine"); - } - return str; } @@ -835,18 +831,6 @@ IndexingMethod get_indm_from_string(const char *str) method = INDEXING_DEBUG; return method; - } else if ( strcmp(bits[i], "raw") == 0) { - method = set_raw(method); - - } else if ( strcmp(bits[i], "bad") == 0) { - method = set_bad(method); - - } else if ( strcmp(bits[i], "comb") == 0) { - method = set_comb(method); /* Default */ - - } else if ( strcmp(bits[i], "axes") == 0) { - method = set_axes(method); - } else if ( strcmp(bits[i], "latt") == 0) { method = set_lattice(method); @@ -859,24 +843,6 @@ IndexingMethod get_indm_from_string(const char *str) } else if ( strcmp(bits[i], "nocell") == 0) { method = set_nocellparams(method); - } else if ( strcmp(bits[i], "retry") == 0) { - method |= INDEXING_RETRY; - - } else if ( strcmp(bits[i], "noretry") == 0) { - method &= ~INDEXING_RETRY; - - } else if ( strcmp(bits[i], "multi") == 0) { - method |= INDEXING_MULTI; - - } else if ( strcmp(bits[i], "nomulti") == 0) { - method &= ~INDEXING_MULTI; - - } else if ( strcmp(bits[i], "refine") == 0) { - method |= INDEXING_REFINE; - - } else if ( strcmp(bits[i], "norefine") == 0) { - method &= ~INDEXING_REFINE; - } else { ERROR("Bad list of indexing methods: '%s'\n", str); return INDEXING_ERROR; diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h index 8b3aadeb..2aad377b 100644 --- a/libcrystfel/src/index.h +++ b/libcrystfel/src/index.h @@ -39,36 +39,22 @@ #endif -#define INDEXING_DEFAULTS_DIRAX (INDEXING_DIRAX | INDEXING_CHECK_PEAKS \ - | INDEXING_CHECK_CELL_COMBINATIONS \ - | INDEXING_RETRY | INDEXING_REFINE) - -#define INDEXING_DEFAULTS_ASDF (INDEXING_ASDF | INDEXING_CHECK_PEAKS \ - | INDEXING_CHECK_CELL_COMBINATIONS \ - | INDEXING_RETRY | INDEXING_REFINE) - -#define INDEXING_DEFAULTS_MOSFLM (INDEXING_MOSFLM | INDEXING_CHECK_PEAKS \ - | INDEXING_CHECK_CELL_COMBINATIONS \ - | INDEXING_USE_LATTICE_TYPE \ - | INDEXING_USE_CELL_PARAMETERS \ - | INDEXING_RETRY | INDEXING_REFINE) - -#define INDEXING_DEFAULTS_FELIX (INDEXING_FELIX \ - | INDEXING_USE_LATTICE_TYPE \ - | INDEXING_USE_CELL_PARAMETERS \ - | INDEXING_RETRY | INDEXING_REFINE) - -#define INDEXING_DEFAULTS_TAKETWO (INDEXING_TAKETWO | INDEXING_CHECK_PEAKS \ - | INDEXING_USE_CELL_PARAMETERS \ - | INDEXING_USE_LATTICE_TYPE \ - | INDEXING_RETRY | INDEXING_REFINE) - -/* Axis check is needed for XDS, because it likes to permute the axes */ -#define INDEXING_DEFAULTS_XDS (INDEXING_XDS | INDEXING_USE_LATTICE_TYPE \ - | INDEXING_USE_CELL_PARAMETERS \ - | INDEXING_CHECK_CELL_AXES \ - | INDEXING_CHECK_PEAKS \ - | INDEXING_RETRY | INDEXING_REFINE) +#define INDEXING_DEFAULTS_DIRAX (INDEXING_DIRAX) + +#define INDEXING_DEFAULTS_ASDF (INDEXING_ASDF) + +#define INDEXING_DEFAULTS_MOSFLM (INDEXING_MOSFLM | INDEXING_USE_LATTICE_TYPE \ + | INDEXING_USE_CELL_PARAMETERS) + +#define INDEXING_DEFAULTS_FELIX (INDEXING_FELIX | INDEXING_USE_LATTICE_TYPE \ + | INDEXING_USE_CELL_PARAMETERS) + +#define INDEXING_DEFAULTS_TAKETWO (INDEXING_TAKETWO \ + | INDEXING_USE_CELL_PARAMETERS \ + | INDEXING_USE_LATTICE_TYPE) + +#define INDEXING_DEFAULTS_XDS (INDEXING_XDS | INDEXING_USE_LATTICE_TYPE \ + | INDEXING_USE_CELL_PARAMETERS) /** * IndexingMethod: @@ -82,21 +68,10 @@ * @INDEXING_ASDF: Use in-built "asdf" indexer * @INDEXING_TAKETWO: Use in-built "taketwo" indexer * @INDEXING_ERROR: Special value for unrecognised indexing engine name - * @INDEXING_CHECK_CELL_COMBINATIONS: Check linear combinations of unit cell - * axes for agreement with given cell. - * @INDEXING_CHECK_CELL_AXES: Check unit cell axes for agreement with given - * cell, and permute them if necessary. - * @INDEXING_CHECK_PEAKS: Check that the peaks can be explained by the indexing - * result. * @INDEXING_USE_LATTICE_TYPE: Use lattice type and centering information to * guide the indexing process. * @INDEXING_USE_CELL_PARAMETERS: Use the unit cell parameters to guide the * indexing process. - * @INDEXING_RETRY: If the indexer doesn't succeed, delete the weakest peaks - * and try again. - * @INDEXING_MULTI: If the indexer succeeds, delete the peaks explained by the - * lattice and try again in the hope of finding another crystal. - * @INDEXING_REFINE: Perform "prediction refinement" after indexing. * * An enumeration of all the available indexing methods. The dummy value * @INDEXING_SIMULATION is used by partial_sim to indicate that no indexing was @@ -120,14 +95,8 @@ typedef enum { /* Bits at the top of the IndexingMethod are flags which modify the * behaviour of the indexer. */ - INDEXING_CHECK_CELL_COMBINATIONS = 256, - INDEXING_CHECK_CELL_AXES = 512, - INDEXING_CHECK_PEAKS = 1024, INDEXING_USE_LATTICE_TYPE = 2048, INDEXING_USE_CELL_PARAMETERS = 4096, - INDEXING_RETRY = 8192, - INDEXING_MULTI = 16384, - INDEXING_REFINE = 32768, } IndexingMethod; @@ -135,8 +104,16 @@ typedef enum { * core of the indexing method */ #define INDEXING_METHOD_MASK (0xff) -/* Indexing flags which the indexing method does not need to know about */ -#define INDEXING_CONTROL_FLAGS (INDEXING_RETRY | INDEXING_MULTI | INDEXING_REFINE) +typedef enum { + + INDEXING_RETRY = 1, + INDEXING_MULTI = 2, + INDEXING_REFINE = 4, + INDEXING_CHECK_CELL_COMBINATIONS = 8, + INDEXING_CHECK_CELL_AXES = 16, + INDEXING_CHECK_PEAKS = 32, + +} IndexingFlags; #ifdef __cplusplus extern "C" { @@ -162,7 +139,7 @@ extern IndexingMethod get_indm_from_string(const char *method); extern IndexingPrivate *setup_indexing(const char *methods, UnitCell *cell, struct detector *det, float *ltl, - int no_refine, const char *options, + IndexingFlags flags, const char *options, struct taketwo_options *ttopts); extern void index_pattern(struct image *image, IndexingPrivate *ipriv); diff --git a/libcrystfel/src/mosflm.c b/libcrystfel/src/mosflm.c index edf8524e..110243cb 100644 --- a/libcrystfel/src/mosflm.c +++ b/libcrystfel/src/mosflm.c @@ -128,7 +128,6 @@ struct mosflm_data { static int check_cell(struct mosflm_private *mp, struct image *image, UnitCell *cell) { - UnitCell *out; Crystal *cr; /* If we sent lattice information, make sure that we got back what we @@ -163,35 +162,13 @@ static int check_cell(struct mosflm_private *mp, struct image *image, } - if ( mp->indm & INDEXING_CHECK_CELL_COMBINATIONS ) { - - out = match_cell(cell, mp->template, 0, mp->ltl, 1); - if ( out == NULL ) return 0; - - } else if ( mp->indm & INDEXING_CHECK_CELL_AXES ) { - - out = match_cell(cell, mp->template, 0, mp->ltl, 0); - if ( out == NULL ) return 0; - - } else { - out = cell_new_from_cell(cell); - } - cr = crystal_new(); if ( cr == NULL ) { ERROR("Failed to allocate crystal.\n"); return 0; } - crystal_set_cell(cr, out); - - if ( mp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, &cr, 1) ) { - cell_free(out); - crystal_free(cr); - return 0; - } - } + crystal_set_cell(cr, cell); image_add_crystal(image, cr); @@ -340,8 +317,6 @@ static int read_newmat(struct mosflm_data *mosflm, const char *filename, mosflm->done = 1; } - cell_free(cell); - return 0; } @@ -850,14 +825,10 @@ void *mosflm_prepare(IndexingMethod *indm, UnitCell *cell, int need_cell = 0; /* Check if cell parameters are needed/provided */ - if ( *indm & INDEXING_CHECK_CELL_COMBINATIONS ) need_cell = 1; - if ( *indm & INDEXING_CHECK_CELL_AXES ) need_cell = 1; if ( *indm & INDEXING_USE_CELL_PARAMETERS ) need_cell = 1; if ( need_cell && !cell_has_parameters(cell) ) { ERROR("Altering your MOSFLM flags because cell parameters were" " not provided.\n"); - *indm &= ~INDEXING_CHECK_CELL_COMBINATIONS; - *indm &= ~INDEXING_CHECK_CELL_AXES; *indm &= ~INDEXING_USE_CELL_PARAMETERS; } @@ -869,36 +840,8 @@ void *mosflm_prepare(IndexingMethod *indm, UnitCell *cell, } /* Flags that MOSFLM knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_CELL_COMBINATIONS - | INDEXING_CHECK_CELL_AXES | INDEXING_CHECK_PEAKS - | INDEXING_USE_LATTICE_TYPE| INDEXING_USE_CELL_PARAMETERS - | INDEXING_CONTROL_FLAGS; - - if ( (*indm & INDEXING_USE_LATTICE_TYPE) - && !((*indm & INDEXING_CHECK_CELL_COMBINATIONS) - || (*indm & INDEXING_CHECK_CELL_AXES)) - && (cell_has_parameters(cell) - || (cell_get_unique_axis(cell) == 'a') - || (cell_get_unique_axis(cell) == 'b') - || (cell_get_unique_axis(cell) == 'c')) ) - { - ERROR("WARNING: The unit cell from %s might have had " - "its axes permuted from the unit cell you gave.\n" - "If this is a problem, consider using " - "mosflm-axes-latt or mosflm-comb-latt instead of " - "mosflm-raw-latt.\n", indexer_str(*indm)); - /* It'll be OK if INDEXING_CHECK_CELL_PARAMETERS and - * MOSFLM 7.2.0 or later, but that's getting a bit too - * complicated to explain to the user. */ - } - - if ( (*indm & INDEXING_USE_CELL_PARAMETERS) - && !((*indm & INDEXING_CHECK_CELL_COMBINATIONS) - || (*indm & INDEXING_CHECK_CELL_AXES)) ) - { - ERROR("WARNING: Using 'mosflm-raw' but providing cell " - "parameters as prior information to mosflm.\n"); - } + *indm &= INDEXING_METHOD_MASK + | INDEXING_USE_LATTICE_TYPE | INDEXING_USE_CELL_PARAMETERS; mp = malloc(sizeof(struct mosflm_private)); if ( mp == NULL ) return NULL; diff --git a/libcrystfel/src/taketwo.c b/libcrystfel/src/taketwo.c index 87da44a1..f7a2fda9 100644 --- a/libcrystfel/src/taketwo.c +++ b/libcrystfel/src/taketwo.c @@ -1460,7 +1460,7 @@ static void cleanup_taketwo_cell(struct TakeTwoCell *ttCell) * successful. **/ static UnitCell *run_taketwo(UnitCell *cell, const struct taketwo_options *opts, - struct rvec *rlps, int rlp_count) + struct rvec *rlps, int rlp_count) { int cell_vec_count = 0; int asym_vec_count = 0; @@ -1581,18 +1581,6 @@ int taketwo_index(struct image *image, const struct taketwo_options *opts, crystal_set_cell(cr, cell); - if ( tp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, &cr, 1) ) { - cell_free(cell); - crystal_free(cr); - // STATUS("Rubbish!!\n"); - - return 0; - } else { - // STATUS("That's good!\n"); - } - } - image_add_crystal(image, cr); return 1; @@ -1605,9 +1593,8 @@ void *taketwo_prepare(IndexingMethod *indm, UnitCell *cell, struct taketwo_private *tp; /* Flags that TakeTwo knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_PEAKS - | INDEXING_USE_LATTICE_TYPE | INDEXING_USE_CELL_PARAMETERS - | INDEXING_CONTROL_FLAGS; + *indm &= INDEXING_METHOD_MASK | INDEXING_USE_LATTICE_TYPE + | INDEXING_USE_CELL_PARAMETERS; if ( !( (*indm & INDEXING_USE_LATTICE_TYPE) && (*indm & INDEXING_USE_CELL_PARAMETERS)) ) diff --git a/libcrystfel/src/xds.c b/libcrystfel/src/xds.c index 093b9160..39578c60 100644 --- a/libcrystfel/src/xds.c +++ b/libcrystfel/src/xds.c @@ -187,39 +187,14 @@ static int xds_readable(struct image *image, struct xds_data *xds) static int check_cell(struct xds_private *xp, struct image *image, UnitCell *cell) { - UnitCell *out; Crystal *cr; - if ( xp->indm & INDEXING_CHECK_CELL_COMBINATIONS ) { - - out = match_cell(cell, xp->cell, 0, xp->ltl, 1); - if ( out == NULL ) return 0; - - } else if ( xp->indm & INDEXING_CHECK_CELL_AXES ) { - - out = match_cell(cell, xp->cell, 0, xp->ltl, 0); - if ( out == NULL ) return 0; - - } else { - out = cell_new_from_cell(cell); - } - cr = crystal_new(); if ( cr == NULL ) { ERROR("Failed to allocate crystal.\n"); return 0; } - - crystal_set_cell(cr, out); - - if ( xp->indm & INDEXING_CHECK_PEAKS ) { - if ( !peak_sanity_check(image, &cr, 1) ) { - cell_free(out); - crystal_free(cr); - return 0; - } - } - + crystal_set_cell(cr, cell); image_add_crystal(image, cr); return 1; @@ -626,14 +601,10 @@ void *xds_prepare(IndexingMethod *indm, UnitCell *cell, int need_cell = 0; /* Check if cell parameters are needed/provided */ - if ( *indm & INDEXING_CHECK_CELL_COMBINATIONS ) need_cell = 1; - if ( *indm & INDEXING_CHECK_CELL_AXES ) need_cell = 1; if ( *indm & INDEXING_USE_CELL_PARAMETERS ) need_cell = 1; if ( need_cell && !cell_has_parameters(cell) ) { ERROR("Altering your XDS flags because cell parameters were" " not provided.\n"); - *indm &= ~INDEXING_CHECK_CELL_COMBINATIONS; - *indm &= ~INDEXING_CHECK_CELL_AXES; *indm &= ~INDEXING_USE_CELL_PARAMETERS; } @@ -652,36 +623,35 @@ void *xds_prepare(IndexingMethod *indm, UnitCell *cell, } if ( (*indm & INDEXING_USE_LATTICE_TYPE) - && !(*indm & INDEXING_USE_CELL_PARAMETERS) ) { + && !(*indm & INDEXING_USE_CELL_PARAMETERS) ) + { ERROR("Invalid XDS options (-latt-nocell): " "try xds-nolatt-nocell.\n"); return NULL; } if ( (*indm & INDEXING_USE_CELL_PARAMETERS) - && !(*indm & INDEXING_USE_LATTICE_TYPE) ) { + && !(*indm & INDEXING_USE_LATTICE_TYPE) ) + { ERROR("Invalid XDS options (-cell-nolatt): " "try xds-nolatt-nocell.\n"); return NULL; } if ( ((*indm & INDEXING_USE_CELL_PARAMETERS) - || (*indm & INDEXING_USE_LATTICE_TYPE)) - && !(*indm & INDEXING_CHECK_CELL_AXES) - && !(*indm & INDEXING_CHECK_CELL_COMBINATIONS) ) { - ERROR("The cell from xds-raw-cell or xds-raw-latt may have had" + || (*indm & INDEXING_USE_LATTICE_TYPE))) + { + ERROR("The cell from xds-cell or xds-latt may have had" " its axes permuted from the cell you provided. If this" - " is a problem, consider using xds-axes-cell.\n"); + " is a problem, consider using xds-cell.\n"); } xp = calloc(1, sizeof(*xp)); if ( xp == NULL ) return NULL; /* Flags that XDS knows about */ - *indm &= INDEXING_METHOD_MASK | INDEXING_CHECK_CELL_COMBINATIONS - | INDEXING_CHECK_CELL_AXES | INDEXING_USE_LATTICE_TYPE - | INDEXING_CHECK_PEAKS | INDEXING_USE_CELL_PARAMETERS - | INDEXING_CONTROL_FLAGS; + *indm &= INDEXING_METHOD_MASK | INDEXING_USE_LATTICE_TYPE + | INDEXING_USE_CELL_PARAMETERS; xp->ltl = ltl; xp->cell = cell; diff --git a/src/indexamajig.c b/src/indexamajig.c index a79234d3..44b91338 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -227,8 +227,13 @@ int main(int argc, char *argv[]) struct beam_params beam; int have_push_res = 0; int len; - int no_refine = 0; char *command_line_peak_path = NULL; + int if_refine = 1; + int if_comb = 1; + int if_axes = 1; + int if_peaks = 0; + int if_multi = 1; + int if_retry = 1; /* Defaults */ iargs.cell = NULL; @@ -312,9 +317,14 @@ int main(int argc, char *argv[]) {"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}, + {"no-refine", 0, &if_refine, 0}, + {"no-cell-combinations",0,&if_comb, 0}, + {"no-check-cell", 0, &if_axes, 0}, + {"check-peaks", 0, &if_peaks, 1}, + {"no-retry", 0, &if_retry, 0}, + {"no-multi", 0, &if_multi, 0}, /* Long-only options which don't actually do anything */ {"no-sat-corr", 0, &iargs.satcorr, 0}, @@ -856,8 +866,30 @@ int main(int argc, char *argv[]) } else { + IndexingFlags flags = 0; + + if ( if_axes ) { + flags |= INDEXING_CHECK_CELL_AXES; + } + if ( if_comb ) { + flags |= INDEXING_CHECK_CELL_COMBINATIONS; + flags &= ~INDEXING_CHECK_CELL_AXES; /* Not needed */ + } + if ( if_refine ) { + flags |= INDEXING_REFINE; + } + if ( if_peaks ) { + flags |= INDEXING_CHECK_PEAKS; + } + if ( if_multi ) { + flags |= INDEXING_MULTI; + } + if ( if_retry ) { + flags |= INDEXING_RETRY; + } + iargs.ipriv = setup_indexing(indm_str, iargs.cell, iargs.det, - iargs.tols, no_refine, + iargs.tols, flags, iargs.felix_options, &iargs.taketwo_opts); if ( iargs.ipriv == NULL ) { -- cgit v1.2.3 From 9e1d57633c8e788acd9a966857f112386aad52f1 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 11 Sep 2017 13:32:19 +0200 Subject: Factorise cell parameter check --- libcrystfel/src/index.c | 75 +++++++++++++++++++++++++++++++------------------ 1 file changed, 48 insertions(+), 27 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index e7ae7a37..d8810716 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -356,6 +356,36 @@ void map_all_peaks(struct image *image) } +static int check_cell(IndexingFlags flags, Crystal *cr, UnitCell *target, + float *tolerance) +{ + if ( (flags & INDEXING_CHECK_CELL_COMBINATIONS) + || (flags & INDEXING_CHECK_CELL_AXES) ) + { + UnitCell *out; + int reduce; + + if ( flags & INDEXING_CHECK_CELL_COMBINATIONS ) + { + reduce = 1; + } else { + reduce = 0; + } + + out = match_cell(crystal_get_cell(cr), + target, 0, tolerance, reduce); + + if ( out == NULL ) { + return 1; + } + + cell_free(crystal_get_cell(cr)); + crystal_set_cell(cr, out); + } + return 0; +} + + /* Return non-zero for "success" */ static int try_indexer(struct image *image, IndexingMethod indm, IndexingPrivate *ipriv, void *mpriv) @@ -416,6 +446,18 @@ static int try_indexer(struct image *image, IndexingMethod indm, crystal_set_profile_radius(cr, 0.02e9); crystal_set_mosaicity(cr, 0.0); + assert( !((ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS) + && (ipriv->flags & INDEXING_CHECK_CELL_AXES)) ); + + /* Pre-refinement unit cell check if requested */ + if ( check_cell(ipriv->flags, cr, ipriv->target_cell, + ipriv->tolerance) ) + { + crystal_set_user_flag(cr, 1); + n_bad++; + continue; + } + /* Prediction refinement if requested */ if ( ipriv->flags & INDEXING_REFINE ) { @@ -426,34 +468,13 @@ static int try_indexer(struct image *image, IndexingMethod indm, } } - /* Unit cell check if requested */ - assert( !((ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS) - && (ipriv->flags & INDEXING_CHECK_CELL_AXES)) ); - if ( (ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS) - || (ipriv->flags & INDEXING_CHECK_CELL_AXES) ) + /* After refinement unit cell check if requested */ + if ( check_cell(ipriv->flags, cr, ipriv->target_cell, + ipriv->tolerance) ) { - UnitCell *out; - int reduce; - - if ( ipriv->flags & INDEXING_CHECK_CELL_COMBINATIONS ) - { - reduce = 1; - } else { - reduce = 0; - } - - out = match_cell(crystal_get_cell(cr), - ipriv->target_cell, 0, - ipriv->tolerance, reduce); - - if ( out == NULL ) { - crystal_set_user_flag(cr, 1); - n_bad++; - continue; - } - - cell_free(crystal_get_cell(cr)); - crystal_set_cell(cr, out); + crystal_set_user_flag(cr, 1); + n_bad++; + continue; } /* Don't do similarity check if this crystal is bad */ -- cgit v1.2.3 From 318a10ca2aded24327456d0fb302be61f870ec02 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Wed, 13 Sep 2017 15:52:16 +0200 Subject: Record indices in write_pairs() (debugging function) --- libcrystfel/src/predict-refine.c | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c index 61ce3629..015fd94d 100644 --- a/libcrystfel/src/predict-refine.c +++ b/libcrystfel/src/predict-refine.c @@ -138,10 +138,12 @@ static void UNUSED write_pairs(const char *filename, struct reflpeak *rps, double write_fs, write_ss; double fs, ss; struct panel *p; + signed int h, k, l; fs = rps[i].peak->fs; ss = rps[i].peak->ss; p = rps[i].panel; + get_indices(rps[i].refl, &h, &k, &l); write_fs = fs + p->orig_min_fs; write_ss = ss + p->orig_min_ss; @@ -152,6 +154,9 @@ static void UNUSED write_pairs(const char *filename, struct reflpeak *rps, x_dev(&rps[i], det), y_dev(&rps[i], det)); + //fprintf(fh, "%4i %4i %4i 0.0 - 0.0 1 %7.2f %7.2f %s\n", + // h, k, l, write_fs, write_ss, p->name); + } fclose(fh); -- cgit v1.2.3 From 0101b82caa3517be533a31f5b6bf5636c25a94da Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 15 Sep 2017 14:28:10 +0200 Subject: Better way of displaying cell check option --- libcrystfel/src/index.c | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index d8810716..258c4bc1 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -82,10 +82,18 @@ static const char *onoff(int a) static void show_indexing_flags(IndexingFlags flags) { - STATUS(" Check unit cell (combinations): %s\n", - onoff(flags & INDEXING_CHECK_CELL_COMBINATIONS)); - STATUS(" Check unit cell (axis permutations only): %s\n", - onoff(flags & INDEXING_CHECK_CELL_AXES)); + char check[64]; + + assert( !((flags & INDEXING_CHECK_CELL_COMBINATIONS) + && (flags & INDEXING_CHECK_CELL_AXES)) ); + strcpy(check, onoff(flags & (INDEXING_CHECK_CELL_COMBINATIONS | INDEXING_CHECK_CELL_AXES))); + if ( flags & INDEXING_CHECK_CELL_AXES ) { + strcat(check, " (axis permutations only)"); + } + if ( flags & INDEXING_CHECK_CELL_COMBINATIONS ) { + strcat(check, " (axis combinations)"); + } + STATUS(" Check unit cell parameters: %s\n", check); STATUS(" Check peak alignment: %s\n", onoff(flags & INDEXING_CHECK_PEAKS)); STATUS(" Refine indexing solutions: %s\n", -- cgit v1.2.3 From f8ec4432d198e00998576108199bf0ab0e4b08d4 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 15 Sep 2017 14:28:49 +0200 Subject: Remove debugging statement --- libcrystfel/src/index.c | 1 - 1 file changed, 1 deletion(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 258c4bc1..94a66566 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -501,7 +501,6 @@ static int try_indexer(struct image *image, IndexingMethod indm, 0.1, deg2rad(5.0), NULL) ) { crystal_set_user_flag(cr, 1); - ERROR("Excluding duplicate cell\n"); } } -- cgit v1.2.3 From 228f4da5d30418f6176557d07d0b063f5ee6ac00 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 15 Sep 2017 14:29:02 +0200 Subject: Add restraints to prediction refinement --- libcrystfel/src/predict-refine.c | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c index 015fd94d..e702cf33 100644 --- a/libcrystfel/src/predict-refine.c +++ b/libcrystfel/src/predict-refine.c @@ -488,6 +488,18 @@ static int iterate(struct reflpeak *rps, int n, UnitCell *cell, } + int k; + for ( k=0; k Date: Fri, 15 Sep 2017 15:42:30 +0200 Subject: Perform peak alignment check if requested --- libcrystfel/src/index.c | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 94a66566..8c1f52e5 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -485,6 +485,16 @@ static int try_indexer(struct image *image, IndexingMethod indm, continue; } + /* Peak alignment check if requested */ + if ( ipriv->flags & INDEXING_CHECK_PEAKS ) + { + if ( !peak_sanity_check(image, &cr, 1) ) { + crystal_set_user_flag(cr, 1); + n_bad++; + continue; + } + } + /* Don't do similarity check if this crystal is bad */ if ( crystal_get_user_flag(cr) ) continue; -- cgit v1.2.3