aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel/src
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2019-08-21 12:23:39 +0200
committerThomas White <taw@physics.org>2019-08-22 17:03:28 +0200
commit0b1db945a12c0e68491412baf322b761511eb016 (patch)
treec725147e0e1a68814c2eb24754847c2181402064 /libcrystfel/src
parent5ab96056a6d1972a248abb154ea5df0ffda44d06 (diff)
Tidy up comparison function definitions
Especially, remove the last ltl/atl tolerance values.
Diffstat (limited to 'libcrystfel/src')
-rw-r--r--libcrystfel/src/cell-utils.c120
-rw-r--r--libcrystfel/src/cell-utils.h20
-rw-r--r--libcrystfel/src/index.c9
3 files changed, 87 insertions, 62 deletions
diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c
index 097ee463..961a6554 100644
--- a/libcrystfel/src/cell-utils.c
+++ b/libcrystfel/src/cell-utils.c
@@ -1671,7 +1671,7 @@ double cell_get_volume(UnitCell *cell)
/**
* \param cell: A UnitCell
* \param reference: Another UnitCell
- * \param tolerance: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
+ * \param tols: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
*
* Compare the two unit cells. If the real space parameters match to within
* the specified tolerances, and the centering matches, this function returns 1.
@@ -1684,7 +1684,8 @@ double cell_get_volume(UnitCell *cell)
* \returns non-zero if the cells match.
*
*/
-int compare_cell_parameters(UnitCell *cell, UnitCell *reference, double *tolerance)
+int compare_cell_parameters(UnitCell *cell, UnitCell *reference,
+ const double *tols)
{
double a1, b1, c1, al1, be1, ga1;
double a2, b2, c2, al2, be2, ga2;
@@ -1696,12 +1697,12 @@ int compare_cell_parameters(UnitCell *cell, UnitCell *reference, double *toleran
cell_get_parameters(cell, &a1, &b1, &c1, &al1, &be1, &ga1);
cell_get_parameters(reference, &a2, &b2, &c2, &al2, &be2, &ga2);
- if ( !within_tolerance(a1, a2, tolerance[0]*100.0) ) return 0;
- if ( !within_tolerance(b1, b2, tolerance[1]*100.0) ) return 0;
- if ( !within_tolerance(c1, c2, tolerance[2]*100.0) ) return 0;
- if ( fabs(al1-al2) > tolerance[3] ) return 0;
- if ( fabs(be1-be2) > tolerance[4] ) return 0;
- if ( fabs(ga1-ga2) > tolerance[5] ) return 0;
+ if ( !within_tolerance(a1, a2, tols[0]*100.0) ) return 0;
+ if ( !within_tolerance(b1, b2, tols[1]*100.0) ) return 0;
+ if ( !within_tolerance(c1, c2, tols[2]*100.0) ) return 0;
+ if ( fabs(al1-al2) > tols[3] ) return 0;
+ if ( fabs(be1-be2) > tols[4] ) return 0;
+ if ( fabs(ga1-ga2) > tols[5] ) return 0;
return 1;
}
@@ -1719,17 +1720,22 @@ static double moduli_check(double ax, double ay, double az,
/**
* \param cell: A UnitCell
* \param reference: Another UnitCell
- * \param ltl: Maximum allowable fractional difference in axis lengths
- * \param atl: Maximum allowable difference in angles (in radians)
+ * \param tols: Pointer to six tolerance values (see below)
*
* Compare the two unit cells. If the axes match in length (to within
- * fractional difference \p ltl) and the axes are aligned to within \p atl radians,
- * this function returns non-zero.
+ * the specified tolerances), 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 comparison is done by checking that the lengths of the unit cell axes are
+ * the same between the two cells, and that the axes have similar directions in
+ * 3D space. The first three tolerance values are the maximum allowable fractional
+ * differences between the a,b,c axis lengths (respectively) of the two cells.
+ * The last three tolerance values are the maximum allowable angles, in radians,
+ * between the directions of the a,b,c axes of the two cells.
+ *
* \p cell and \p reference must have the same centering. Otherwise, this
* function always returns zero.
*
@@ -1737,7 +1743,7 @@ static double moduli_check(double ax, double ay, double az,
*
*/
int compare_cell_parameters_and_orientation(UnitCell *cell, UnitCell *reference,
- const double ltl, const double atl)
+ const double *tols)
{
double ax1, ay1, az1, bx1, by1, bz1, cx1, cy1, cz1;
double ax2, ay2, az2, bx2, by2, bz2, cx2, cy2, cz2;
@@ -1752,13 +1758,13 @@ int compare_cell_parameters_and_orientation(UnitCell *cell, UnitCell *reference,
&bx2, &by2, &bz2,
&cx2, &cy2, &cz2);
- if ( angle_between(ax1, ay1, az1, ax2, ay2, az2) > atl ) return 0;
- if ( angle_between(bx1, by1, bz1, bx2, by2, bz2) > atl ) return 0;
- if ( angle_between(cx1, cy1, cz1, cx2, cy2, cz2) > atl ) return 0;
+ if ( angle_between(ax1, ay1, az1, ax2, ay2, az2) > tols[3] ) return 0;
+ if ( angle_between(bx1, by1, bz1, bx2, by2, bz2) > tols[4] ) return 0;
+ if ( angle_between(cx1, cy1, cz1, cx2, cy2, cz2) > tols[5] ) return 0;
- if ( moduli_check(ax1, ay1, az1, ax2, ay2, az2) > ltl ) return 0;
- if ( moduli_check(bx1, by1, bz1, bx2, by2, bz2) > ltl ) return 0;
- if ( moduli_check(cx1, cy1, cz1, cx2, cy2, cz2) > ltl ) return 0;
+ if ( moduli_check(ax1, ay1, az1, ax2, ay2, az2) > tols[0] ) return 0;
+ if ( moduli_check(bx1, by1, bz1, bx2, by2, bz2) > tols[1] ) return 0;
+ if ( moduli_check(cx1, cy1, cz1, cx2, cy2, cz2) > tols[2] ) return 0;
return 1;
}
@@ -1767,18 +1773,23 @@ int compare_cell_parameters_and_orientation(UnitCell *cell, UnitCell *reference,
/**
* \param cell: A UnitCell
* \param reference: Another UnitCell
- * \param ltl: Maximum allowable fractional difference in reciprocal axis lengths
- * \param atl: Maximum allowable difference in reciprocal angles (in radians)
+ * \param tols: Pointer to six tolerance values (see below)
* \param pmb: Place to store pointer to matrix
*
* 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 \p ltl)
- * and the axes aligned to within \p atl radians, this function returns non-zero
- * and stores the transformation which must be applied to \p cell at \p pmb.
+ * axes can be made to match in length and the axes aligned in space, this
+ * function returns non-zero and stores the transformation which must be applied
+ * to \p cell at \p pmb.
+ *
+ * A "permutation" means a transformation represented by a matrix with all
+ * elements equal to +1, 0 or -1, having determinant +1 or -1. That means that
+ * this function will find the relationship between a left-handed and a right-
+ * handed basis.
*
* 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().
+ * compare_cell_parameters_and_orientation(). See that function for more
+ * details.
*
* \p cell and \p reference must have the same centering. Otherwise, this
* function always returns zero.
@@ -1788,7 +1799,7 @@ int compare_cell_parameters_and_orientation(UnitCell *cell, UnitCell *reference,
*/
int compare_permuted_cell_parameters_and_orientation(UnitCell *cell,
UnitCell *reference,
- double ltl, double atl,
+ const double *tols,
IntegerMatrix **pmb)
{
IntegerMatrix *m;
@@ -1823,7 +1834,7 @@ int compare_permuted_cell_parameters_and_orientation(UnitCell *cell,
nc = cell_transform_intmat(cell, m);
if ( compare_cell_parameters_and_orientation(nc, reference,
- ltl, atl) )
+ tols) )
{
if ( pmb != NULL ) *pmb = m;
cell_free(nc);
@@ -1950,13 +1961,13 @@ static double g6_distance(UnitCell *cell1, UnitCell *cell2)
/**
* \param cell_in: A UnitCell
* \param reference_in: Another UnitCell
- * \param tolerance: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
+ * \param tols: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
* \param csl: Non-zero to look for coincidence site lattice relationships
* \param pmb: Place to store pointer to matrix
*
* Compare the \p cell_in with \p reference_in. If \p cell is a derivative lattice
- * of \p reference, within fractional axis length differences \p tolerance[0..2]
- * and absolute angle difference \p atl[3..5] (in radians), this function returns
+ * of \p reference, within fractional axis length differences \p tols[0..2]
+ * and absolute angle difference \p tols[3..5] (in radians), this function returns
* non-zero and stores the transformation which needs to be applied to
* \p cell_in at \p pmb.
*
@@ -1975,12 +1986,16 @@ static double g6_distance(UnitCell *cell1, UnitCell *cell2)
* meaning that the lattice points of the transformed version of \p cell_in
* might not coincide with lattice points of \p reference_in.
*
+ * This function is used by CrystFEL's cell_tool program to find non-obvious
+ * relationships between crystal lattices. For most routine comparisons, this
+ * function is probably not the one you need!
+ *
* \returns non-zero if the cells match, zero for no match or error.
*
*/
-int compare_reindexed_cell_parameters(UnitCell *cell_in, UnitCell *reference_in,
- double *tolerance, int csl,
- RationalMatrix **pmb)
+int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in,
+ const double *tols, int csl,
+ RationalMatrix **pmb)
{
UnitCell *cell;
UnitCell *reference;
@@ -2016,9 +2031,9 @@ int compare_reindexed_cell_parameters(UnitCell *cell_in, UnitCell *reference_in,
&cv[0], &cv[1], &cv[2]);
/* Find vectors in 'cell' with lengths close to a, b and c */
- cand_a = find_candidates(a, av, bv, cv, tolerance[0], csl, &ncand_a);
- cand_b = find_candidates(b, av, bv, cv, tolerance[1], csl, &ncand_b);
- cand_c = find_candidates(c, av, bv, cv, tolerance[2], csl, &ncand_c);
+ cand_a = find_candidates(a, av, bv, cv, tols[0], csl, &ncand_a);
+ cand_b = find_candidates(b, av, bv, cv, tols[1], csl, &ncand_b);
+ cand_c = find_candidates(c, av, bv, cv, tols[2], csl, &ncand_c);
if ( (ncand_a==0) || (ncand_b==0) || (ncand_c==0) ) {
*pmb = NULL;
@@ -2055,7 +2070,7 @@ int compare_reindexed_cell_parameters(UnitCell *cell_in, UnitCell *reference_in,
test = cell_transform_rational(cell, M);
cell_get_parameters(test, &at, &bt, &ct, &alt, &bet, &gat);
cell_free(test);
- if ( fabs(gat - ga) > tolerance[5] ) continue;
+ if ( fabs(gat - ga) > tols[5] ) continue;
/* Gamma OK, now look for place for c axis */
for ( ic=0; ic<ncand_c; ic++ ) {
@@ -2081,11 +2096,11 @@ int compare_reindexed_cell_parameters(UnitCell *cell_in, UnitCell *reference_in,
cell_free(test);
continue;
}
- if ( fabs(alt - al) > tolerance[3] ) {
+ if ( fabs(alt - al) > tols[3] ) {
cell_free(test);
continue;
}
- if ( fabs(bet - be) > tolerance[4] ) {
+ if ( fabs(bet - be) > tols[4] ) {
cell_free(test);
continue;
}
@@ -2434,25 +2449,30 @@ IntegerMatrix *reduce_g6(struct g6 g, double epsrel)
/**
* \param cell_in: A UnitCell
* \param reference_in: Another UnitCell
- * \param tolerance: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
- * \param pmb: Place to store pointer to matrix
+ * \param tols: Pointer to tolerances for a,b,c (fractional), al,be,ga (radians)
+ * \param pmb: Place to store pointer to matrix, or NULL if not needed
*
* Compare the \p cell_in with \p reference_in. If they represent the same
- * lattice, this function returns non-zero and stores the transformation which
- * needs to be applied to \p cell_in at \p pmb.
+ * lattice, this function returns a copy of \p cell_in transformed to look
+ * similar to \p reference_in. Otherwise, it returns NULL.
*
- * Subject to the tolerances, this function will find the transformation which
- * gives the best match to the reference cell, using the Euclidian norm in
- * G6 [see e.g. Andrews and Bernstein, Acta Cryst. A44 (1988) p1009].
+ * If \pmb is non-NULL, the transformation which needs to be applied to
+ * \p cell_in will be stored there.
*
* Only the cell parameters will be compared. The relative orientations are
- * irrelevant.
+ * irrelevant. The tolerances will be applied to the transformed copy of
+ * \p cell_in.
*
- * \returns non-zero if the cells match, zero for no match or error.
+ * This is the right function to use for deciding if an indexing solution
+ * matches a reference cell or not.
+ *
+ * \returns A newly allocated UnitCell, or NULL.
*
*/
-int compare_lattices(UnitCell *cell_in, UnitCell *reference_in,
- double *tolerance, RationalMatrix **pmb)
+UnitCell *compare_reindexed_cell_parameters(UnitCell *cell_in,
+ UnitCell *reference_in,
+ const double *tols,
+ RationalMatrix **pmb)
{
UnitCell *cell;
UnitCell *reference;
diff --git a/libcrystfel/src/cell-utils.h b/libcrystfel/src/cell-utils.h
index e46d860b..bb1e5cd0 100644
--- a/libcrystfel/src/cell-utils.h
+++ b/libcrystfel/src/cell-utils.h
@@ -88,25 +88,25 @@ extern int forbidden_reflection(UnitCell *cell,
extern double cell_get_volume(UnitCell *cell);
extern int compare_cell_parameters(UnitCell *cell, UnitCell *reference,
- double *tolerance);
+ const double *tols);
extern int compare_cell_parameters_and_orientation(UnitCell *cell,
UnitCell *reference,
- const double ltl,
- const double atl);
+ const double *tols);
extern int compare_permuted_cell_parameters_and_orientation(UnitCell *cell,
UnitCell *reference,
- double ltl,
- double atl,
+ const double *tols,
IntegerMatrix **pmb);
-extern int compare_reindexed_cell_parameters(UnitCell *cell, UnitCell *reference,
- double *tolerance, int csl,
- RationalMatrix **pmb);
+extern int compare_derivative_cell_parameters(UnitCell *cell, UnitCell *reference,
+ const double *tols, int csl,
+ RationalMatrix **pmb);
-extern int compare_lattices(UnitCell *cell_in, UnitCell *reference_in,
- double *tolerance, RationalMatrix **pmb);
+extern UnitCell *compare_reindexed_cell_parameters(UnitCell *cell_in,
+ UnitCell *reference_in,
+ const double *tols,
+ RationalMatrix **pmb);
#ifdef __cplusplus
}
diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c
index aa0cfb6a..45b24afb 100644
--- a/libcrystfel/src/index.c
+++ b/libcrystfel/src/index.c
@@ -553,7 +553,8 @@ static int check_cell(IndexingFlags flags, Crystal *cr, UnitCell *target,
if ( ! ((flags & INDEXING_CHECK_CELL_COMBINATIONS)
|| (flags & INDEXING_CHECK_CELL_AXES)) ) return 0;
- if ( compare_lattices(crystal_get_cell(cr), target, tolerance, &rm) )
+ if ( compare_reindexed_cell_parameters(crystal_get_cell(cr), target,
+ tolerance, &rm) )
{
out = cell_transform_rational(crystal_get_cell(cr), rm);
cell_free(crystal_get_cell(cr));
@@ -691,13 +692,17 @@ static int try_indexer(struct image *image, IndexingMethod indm,
for ( j=0; j<this_crystal; j++ ) {
Crystal *that_cr = image->crystals[j];
+ const double tols[] = {0.1, 0.1, 0.1,
+ deg2rad(5.0),
+ deg2rad(5.0),
+ deg2rad(5.0)};
/* Don't do similarity check against bad crystals */
if ( crystal_get_user_flag(that_cr) ) continue;
if ( compare_permuted_cell_parameters_and_orientation(crystal_get_cell(cr),
crystal_get_cell(that_cr),
- 0.1, deg2rad(0.5), NULL) )
+ tols, NULL) )
{
crystal_set_user_flag(cr, 1);
}