diff options
| -rw-r--r-- | libcrystfel/src/cell-utils.c | 114 | ||||
| -rw-r--r-- | libcrystfel/src/cell-utils.h | 16 | ||||
| -rw-r--r-- | libcrystfel/src/index.c | 6 | ||||
| -rw-r--r-- | src/cell_tool.c | 28 | ||||
| -rw-r--r-- | src/whirligig.c | 13 |
5 files changed, 116 insertions, 61 deletions
diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c index 7b1984bb..3d63249d 100644 --- a/libcrystfel/src/cell-utils.c +++ b/libcrystfel/src/cell-utils.c @@ -1694,6 +1694,48 @@ double cell_get_volume(UnitCell *cell) } +/** + * compare_cell_parameters: + * @a: A %UnitCell + * @b: Another %UnitCell + * @ltl: Maximum allowable fractional difference in axis lengths + * @atl: Maximum allowable difference in reciprocal angles (in radians) + * + * Compare the two unit cells. If the real space parameters match to within + * fractional difference @ltl, and the inter-axial angles match within @atl, + * and the centering matches, this function returns 1. Otherwise 0. + * + * This function considers the cell parameters and centering, but ignores the + * orientation of the cell. If you want to compare the orientation as well, + * use compare_cell_parameters_and_orientation() instead. + * + * Returns: non-zero if the cells match. + * + */ +int compare_cell_parameters(UnitCell *cell1, UnitCell *cell2, + float ltl, float atl) +{ + double a1, b1, c1, al1, be1, ga1; + double a2, b2, c2, al2, be2, ga2; + + /* Centering must match: we don't arbitrarte primitive vs centered, + * different cell choices etc */ + if ( cell_get_centering(cell1) != cell_get_centering(cell2) ) return 0; + + cell_get_parameters(cell1, &a1, &b1, &c1, &al1, &be1, &ga1); + cell_get_parameters(cell2, &a2, &b2, &c2, &al2, &be2, &ga2); + + if ( !within_tolerance(a1, a2, ltl) ) return 0; + if ( !within_tolerance(b1, b2, ltl) ) return 0; + if ( !within_tolerance(c1, c2, ltl) ) return 0; + if ( fabs(al1-al2) > atl ) return 0; + if ( fabs(be1-be2) > atl ) return 0; + if ( fabs(ga1-ga2) > atl ) return 0; + + return 1; +} + + static double moduli_check(double ax, double ay, double az, double bx, double by, double bz) { @@ -1703,28 +1745,42 @@ static double moduli_check(double ax, double ay, double az, } -static int cells_are_similar(UnitCell *cell1, UnitCell *cell2, - const double ltl, const double atl) +/** + * compare_cell_parameters_and_orientation: + * @a: A %UnitCell + * @b: Another %UnitCell + * @ltl: Maximum allowable fractional difference in reciprocal axis lengths + * @atl: Maximum allowable difference in reciprocal angles (in radians) + * + * Compare the two unit cells. If the axes match in length (to within + * fractional difference @ltl) and the axes are aligned to within @atl radians, + * this function returns non-zero. + * + * This function compares the orientation of the cell as well as the parameters. + * If you just want to see if the parameters are the same, use + * compare_cell_parameters() instead. + * + * The cells @a and @b must have the same centering. Otherwise, this function + * always returns zero. + * + * Returns: non-zero if the cells match. + * + */ +int compare_cell_parameters_and_orientation(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); + if ( cell_get_centering(cell1) != cell_get_centering(cell2) ) return 0; + cell_get_cartesian(cell1, &asx1, &asy1, &asz1, + &bsx1, &bsy1, &bsz1, + &csx1, &csy1, &csz1); - cell_free(pcell1); - cell_free(pcell2); + cell_get_cartesian(cell2, &asx2, &asy2, &asz2, + &bsx2, &bsy2, &bsz2, + &csx2, &csy2, &csz2); if ( angle_between(asx1, asy1, asz1, asx2, asy2, asz2) > atl ) return 0; if ( angle_between(bsx1, bsy1, bsz1, bsx2, bsy2, bsz2) > atl ) return 0; @@ -1738,28 +1794,38 @@ static int cells_are_similar(UnitCell *cell1, UnitCell *cell2, } - /** - * compare_cells: + * compare_reindexed_cell_parameters_and_orientation: * @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. + * Compare the two unit cells. If, using any permutation of the axes, the + * axes can be made to match in length (to within fractional difference @ltl) + * and the axes aligned to within @atl radians, this function returns non-zero + * and stores the transformation to map @b onto @a. + * + * Note that the orientations of the cells must match, not just the parameters. + * The comparison is done after reindexing using + * compare_cell_parameters_and_orientation(). + * + * The cells @a and @b must have the same centering. Otherwise, this function + * always returns zero. * * Returns: non-zero if the cells match. * */ -int compare_cells(UnitCell *a, UnitCell *b, double ltl, double atl, - IntegerMatrix **pmb) +int compare_reindexed_cell_parameters_and_orientation(UnitCell *a, UnitCell *b, + double ltl, double atl, + IntegerMatrix **pmb) { IntegerMatrix *m; int i[9]; + if ( cell_get_centering(a) != cell_get_centering(b) ) return 0; + m = intmat_new(3, 3); for ( i[0]=-1; i[0]<=+1; i[0]++ ) { @@ -1786,7 +1852,7 @@ int compare_cells(UnitCell *a, UnitCell *b, double ltl, double atl, tfn = tfn_from_intmat(m); nc = cell_transform(b, tfn); - if ( cells_are_similar(a, nc, ltl, atl) ) { + if ( compare_cell_parameters_and_orientation(a, nc, ltl, atl) ) { if ( pmb != NULL ) *pmb = m; tfn_free(tfn); cell_free(nc); diff --git a/libcrystfel/src/cell-utils.h b/libcrystfel/src/cell-utils.h index d0af717a..cf736e4a 100644 --- a/libcrystfel/src/cell-utils.h +++ b/libcrystfel/src/cell-utils.h @@ -80,8 +80,20 @@ 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); +extern int compare_cell_parameters(UnitCell *cell1, UnitCell *cell2, + float ltl, float atl); + + +extern int compare_cell_parameters_and_orientation(UnitCell *cell1, + UnitCell *cell2, + const double ltl, + const double atl); + +extern int compare_reindexed_cell_parameters_and_orientation(UnitCell *a, + UnitCell *b, + double ltl, + double atl, + IntegerMatrix **pmb); #ifdef __cplusplus } diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index f72334c4..ad9c0de3 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -698,9 +698,9 @@ static int try_indexer(struct image *image, IndexingMethod indm, /* Don't do similarity check against bad crystals */ if ( crystal_get_user_flag(that_cr) ) continue; - if ( compare_cells(crystal_get_cell(cr), - crystal_get_cell(that_cr), - 0.1, deg2rad(0.5), NULL) ) + if ( compare_reindexed_cell_parameters_and_orientation(crystal_get_cell(cr), + crystal_get_cell(that_cr), + 0.1, deg2rad(0.5), NULL) ) { crystal_set_user_flag(cr, 1); } diff --git a/src/cell_tool.c b/src/cell_tool.c index 1a2f33fb..ee82f332 100644 --- a/src/cell_tool.c +++ b/src/cell_tool.c @@ -68,30 +68,6 @@ static void show_help(const char *s) } -static int cells_the_same(UnitCell *cell1, UnitCell *cell2, float ltl, float atl) -{ - double a1, b1, c1, al1, be1, ga1; - double a2, b2, c2, al2, be2, ga2; - UnitCellTransformation *tfn1, *tfn2; - - /* Centering must match: we don't arbitrarte primitive vs centered, - * different cell choices etc */ - if ( cell_get_centering(cell1) != cell_get_centering(cell2) ) return 0; - - cell_get_parameters(cell1, &a1, &b1, &c1, &al1, &be1, &ga1); - cell_get_parameters(cell2, &a2, &b2, &c2, &al2, &be2, &ga2); - - if ( !within_tolerance(a1, a2, ltl) ) return 0; - if ( !within_tolerance(b1, b2, ltl) ) return 0; - if ( !within_tolerance(c1, c2, ltl) ) return 0; - if ( fabs(al1-al2) > atl ) return 0; - if ( fabs(be1-be2) > atl ) return 0; - if ( fabs(ga1-ga2) > atl ) return 0; - - return 1; -} - - static int comparecells(UnitCell *cell, const char *comparecell, double ltl, double atl) { @@ -146,7 +122,7 @@ static int comparecells(UnitCell *cell, const char *comparecell, tfn = tfn_from_intmat(m); nc = cell_transform(cell, tfn); - if ( cells_the_same(cell2, nc, ltl, atl) ) { + if ( compare_cell_parameters(cell2, nc, ltl, atl) ) { STATUS("-----------------------------------------------" "-------------------------------------------\n"); cell_print(nc); @@ -327,7 +303,7 @@ static int find_ambi(UnitCell *cell, SymOpList *sym, double ltl, double atl) tfn = tfn_from_intmat(m); nc = cell_transform(cell, tfn); - if ( cells_the_same(cell, nc, ltl, atl) ) { + if ( compare_cell_parameters(cell, nc, ltl, atl) ) { if ( !intmat_is_identity(m) ) add_symop(ops, m); STATUS("-----------------------------------------------" "-------------------------------------------\n"); diff --git a/src/whirligig.c b/src/whirligig.c index f5a435c2..54ade40a 100644 --- a/src/whirligig.c +++ b/src/whirligig.c @@ -309,9 +309,9 @@ static IntegerMatrix *try_all(struct window *win, int n1, int n2, for ( i=0; i<i1->n_crystals; i++ ) { for ( j=0; j<i2->n_crystals; j++ ) { - if ( compare_cells(crystal_get_cell(i1->crystals[i]), - crystal_get_cell(i2->crystals[j]), - 0.1, deg2rad(5.0), &m) ) + if ( compare_reindexed_cell_parameters_and_orientation(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) ) @@ -377,9 +377,10 @@ static int try_join(struct window *win, int sn) for ( j=0; j<win->img[win->join_ptr].n_crystals; j++ ) { Crystal *cr2; cr2 = win->img[win->join_ptr].crystals[j]; - if ( compare_cells(ref, crystal_get_cell(cr2), - 0.1, deg2rad(5.0), - &win->mat[sn][win->join_ptr]) ) { + if ( compare_reindexed_cell_parameters_and_orientation(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; |