From 169e7c5677ffc9c296c0a7eeddb0b77e024a4a55 Mon Sep 17 00:00:00 2001 From: Thomas White Date: Mon, 18 Feb 2019 16:04:28 +0100 Subject: Be clear about whether functions need rational or integer transformations --- libcrystfel/src/cell-utils.c | 8 ++-- libcrystfel/src/cell.c | 88 ++++++++++++++++++++++++++++++++------------ libcrystfel/src/cell.h | 8 +++- libcrystfel/src/xgandalf.c | 2 +- src/cell_tool.c | 4 +- src/whirligig.c | 2 +- tests/transformation_check.c | 6 +-- 7 files changed, 81 insertions(+), 37 deletions(-) diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c index 507d4cc2..da1e6fd5 100644 --- a/libcrystfel/src/cell-utils.c +++ b/libcrystfel/src/cell-utils.c @@ -500,7 +500,7 @@ UnitCell *uncenter_cell(UnitCell *in, IntegerMatrix **t) tt = centering_transformation(in, &new_centering, &new_latt, &new_ua); if ( tt == NULL ) return NULL; - out = cell_transform_inverse(in, tt); + out = cell_transform_intmat_inverse(in, tt); if ( out == NULL ) return NULL; cell_set_lattice_type(out, new_latt); @@ -776,7 +776,7 @@ UnitCell *match_cell(UnitCell *cell_in, UnitCell *template_in, int verbose, /* Reverse the de-centering transformation */ if ( new_cell != NULL ) { - new_cell_trans = cell_transform(new_cell, centering); + new_cell_trans = cell_transform_intmat(new_cell, centering); cell_free(new_cell); cell_set_lattice_type(new_cell_trans, cell_get_lattice_type(template_in)); @@ -905,7 +905,7 @@ UnitCell *match_cell_ab(UnitCell *cell_in, UnitCell *template_in) new_cell = cell_new_from_direct_axes(real_a, real_b, real_c); /* Reverse the de-centering transformation */ - new_cell_trans = cell_transform_inverse(new_cell, to_given_cell); + new_cell_trans = cell_transform_intmat_inverse(new_cell, to_given_cell); cell_free(new_cell); cell_set_lattice_type(new_cell, cell_get_lattice_type(template_in)); cell_set_centering(new_cell, cell_get_centering(template_in)); @@ -1770,7 +1770,7 @@ int compare_reindexed_cell_parameters_and_orientation(UnitCell *a, UnitCell *b, if ( intmat_det(m) != +1 ) continue; - nc = cell_transform(b, m); + nc = cell_transform_intmat(b, m); if ( compare_cell_parameters_and_orientation(a, nc, ltl, atl) ) { if ( pmb != NULL ) *pmb = m; diff --git a/libcrystfel/src/cell.c b/libcrystfel/src/cell.c index 242094f6..5eb13349 100644 --- a/libcrystfel/src/cell.c +++ b/libcrystfel/src/cell.c @@ -688,16 +688,16 @@ UnitCell *cell_transform_gsl_reciprocal(UnitCell *in, gsl_matrix *m) /** - * cell_transform: + * cell_transform_rational: * @cell: A %UnitCell. - * @t: An %IntegerMatrix. + * @t: A %RationalMatrix. * * Applies @t to @cell. * * Returns: Transformed copy of @cell. * */ -UnitCell *cell_transform(UnitCell *cell, IntegerMatrix *m) +UnitCell *cell_transform_rational(UnitCell *cell, RationalMatrix *m) { UnitCell *out; gsl_matrix *tm; @@ -709,15 +709,15 @@ UnitCell *cell_transform(UnitCell *cell, IntegerMatrix *m) return NULL; } - gsl_matrix_set(tm, 0, 0, intmat_get(m, 0, 0)); - gsl_matrix_set(tm, 0, 1, intmat_get(m, 0, 1)); - gsl_matrix_set(tm, 0, 2, intmat_get(m, 0, 2)); - gsl_matrix_set(tm, 1, 0, intmat_get(m, 1, 0)); - gsl_matrix_set(tm, 1, 1, intmat_get(m, 1, 1)); - gsl_matrix_set(tm, 1, 2, intmat_get(m, 1, 2)); - gsl_matrix_set(tm, 2, 0, intmat_get(m, 2, 0)); - gsl_matrix_set(tm, 2, 1, intmat_get(m, 2, 1)); - gsl_matrix_set(tm, 2, 2, intmat_get(m, 2, 2)); + gsl_matrix_set(tm, 0, 0, rtnl_as_double(rtnl_mtx_get(m, 0, 0))); + gsl_matrix_set(tm, 0, 1, rtnl_as_double(rtnl_mtx_get(m, 0, 1))); + gsl_matrix_set(tm, 0, 2, rtnl_as_double(rtnl_mtx_get(m, 0, 2))); + gsl_matrix_set(tm, 1, 0, rtnl_as_double(rtnl_mtx_get(m, 1, 0))); + gsl_matrix_set(tm, 1, 1, rtnl_as_double(rtnl_mtx_get(m, 1, 1))); + gsl_matrix_set(tm, 1, 2, rtnl_as_double(rtnl_mtx_get(m, 1, 2))); + gsl_matrix_set(tm, 2, 0, rtnl_as_double(rtnl_mtx_get(m, 2, 0))); + gsl_matrix_set(tm, 2, 1, rtnl_as_double(rtnl_mtx_get(m, 2, 1))); + gsl_matrix_set(tm, 2, 2, rtnl_as_double(rtnl_mtx_get(m, 2, 2))); out = cell_transform_gsl_direct(cell, tm); @@ -730,16 +730,36 @@ UnitCell *cell_transform(UnitCell *cell, IntegerMatrix *m) /** - * cell_transform_inverse: + * cell_transform_intmat: * @cell: A %UnitCell. - * @m: An %IntegerMatrix + * @t: An %IntegerMatrix. + * + * Applies @t to @cell. + * + * Returns: Transformed copy of @cell. + * + */ +UnitCell *cell_transform_intmat(UnitCell *cell, IntegerMatrix *m) +{ + UnitCell *ans; + RationalMatrix *mtx = rtnl_mtx_from_intmat(m); + ans = cell_transform_rational(cell, mtx); + rtnl_mtx_free(mtx); + return ans; +} + + +/** + * cell_transform_rational_inverse: + * @cell: A %UnitCell. + * @m: A %RationalMatrix * * Applies the inverse of @m to @cell. * * Returns: Transformed copy of @cell. * */ -UnitCell *cell_transform_inverse(UnitCell *cell, IntegerMatrix *m) +UnitCell *cell_transform_rational_inverse(UnitCell *cell, RationalMatrix *m) { UnitCell *out; gsl_matrix *tm; @@ -754,15 +774,15 @@ UnitCell *cell_transform_inverse(UnitCell *cell, IntegerMatrix *m) return NULL; } - gsl_matrix_set(tm, 0, 0, intmat_get(m, 0, 0)); - gsl_matrix_set(tm, 0, 1, intmat_get(m, 0, 1)); - gsl_matrix_set(tm, 0, 2, intmat_get(m, 0, 2)); - gsl_matrix_set(tm, 1, 0, intmat_get(m, 1, 0)); - gsl_matrix_set(tm, 1, 1, intmat_get(m, 1, 1)); - gsl_matrix_set(tm, 1, 2, intmat_get(m, 1, 2)); - gsl_matrix_set(tm, 2, 0, intmat_get(m, 2, 0)); - gsl_matrix_set(tm, 2, 1, intmat_get(m, 2, 1)); - gsl_matrix_set(tm, 2, 2, intmat_get(m, 2, 2)); + gsl_matrix_set(tm, 0, 0, rtnl_as_double(rtnl_mtx_get(m, 0, 0))); + gsl_matrix_set(tm, 0, 1, rtnl_as_double(rtnl_mtx_get(m, 0, 1))); + gsl_matrix_set(tm, 0, 2, rtnl_as_double(rtnl_mtx_get(m, 0, 2))); + gsl_matrix_set(tm, 1, 0, rtnl_as_double(rtnl_mtx_get(m, 1, 0))); + gsl_matrix_set(tm, 1, 1, rtnl_as_double(rtnl_mtx_get(m, 1, 1))); + gsl_matrix_set(tm, 1, 2, rtnl_as_double(rtnl_mtx_get(m, 1, 2))); + gsl_matrix_set(tm, 2, 0, rtnl_as_double(rtnl_mtx_get(m, 2, 0))); + gsl_matrix_set(tm, 2, 1, rtnl_as_double(rtnl_mtx_get(m, 2, 1))); + gsl_matrix_set(tm, 2, 2, rtnl_as_double(rtnl_mtx_get(m, 2, 2))); perm = gsl_permutation_alloc(3); if ( perm == NULL ) { @@ -796,3 +816,23 @@ UnitCell *cell_transform_inverse(UnitCell *cell, IntegerMatrix *m) return out; } + + +/** + * cell_transform_intmat_inverse: + * @cell: A %UnitCell. + * @m: An %IntegerMatrix + * + * Applies the inverse of @m to @cell. + * + * Returns: Transformed copy of @cell. + * + */ +UnitCell *cell_transform_intmat_inverse(UnitCell *cell, IntegerMatrix *m) +{ + UnitCell *ans; + RationalMatrix *mtx = rtnl_mtx_from_intmat(m); + ans = cell_transform_rational_inverse(cell, mtx); + rtnl_mtx_free(mtx); + return ans; +} diff --git a/libcrystfel/src/cell.h b/libcrystfel/src/cell.h index e6ccbdf7..8e6809b5 100644 --- a/libcrystfel/src/cell.h +++ b/libcrystfel/src/cell.h @@ -150,8 +150,12 @@ extern const char *cell_rep(UnitCell *cell); extern UnitCell *cell_transform_gsl_direct(UnitCell *in, gsl_matrix *m); extern UnitCell *cell_transform_gsl_reciprocal(UnitCell *in, gsl_matrix *m); -extern UnitCell *cell_transform(UnitCell *cell, IntegerMatrix *m); -extern UnitCell *cell_transform_inverse(UnitCell *cell, IntegerMatrix *m); + +extern UnitCell *cell_transform_rational(UnitCell *cell, RationalMatrix *m); +extern UnitCell *cell_transform_rational_inverse(UnitCell *cell, RationalMatrix *m); + +extern UnitCell *cell_transform_intmat(UnitCell *cell, IntegerMatrix *m); +extern UnitCell *cell_transform_intmat_inverse(UnitCell *cell, IntegerMatrix *m); #ifdef __cplusplus } diff --git a/libcrystfel/src/xgandalf.c b/libcrystfel/src/xgandalf.c index e6618560..b122cb4c 100644 --- a/libcrystfel/src/xgandalf.c +++ b/libcrystfel/src/xgandalf.c @@ -114,7 +114,7 @@ int run_xgandalf(struct image *image, void *ipriv) if(xgandalf_private_data->cellTemplate != NULL){ restoreCell(uc, &xgandalf_private_data->latticeReductionTransform); - UnitCell *new_cell_trans = cell_transform(uc, xgandalf_private_data->centeringTransformation); + UnitCell *new_cell_trans = cell_transform_intmat(uc, xgandalf_private_data->centeringTransformation); cell_free(uc); uc = new_cell_trans; diff --git a/src/cell_tool.c b/src/cell_tool.c index 12a00192..20e156a1 100644 --- a/src/cell_tool.c +++ b/src/cell_tool.c @@ -119,7 +119,7 @@ static int comparecells(UnitCell *cell, const char *comparecell, if ( intmat_det(m) < 1 ) continue; - nc = cell_transform(cell, m); + nc = cell_transform_intmat(cell, m); if ( compare_cell_parameters(cell2, nc, ltl, atl) ) { STATUS("-----------------------------------------------" @@ -294,7 +294,7 @@ static int find_ambi(UnitCell *cell, SymOpList *sym, double ltl, double atl) if ( intmat_det(m) != +1 ) continue; - nc = cell_transform(cell, m); + nc = cell_transform_intmat(cell, m); if ( compare_cell_parameters(cell, nc, ltl, atl) ) { if ( !intmat_is_identity(m) ) add_symop(ops, m); diff --git a/src/whirligig.c b/src/whirligig.c index 7d50141a..209e7fbe 100644 --- a/src/whirligig.c +++ b/src/whirligig.c @@ -369,7 +369,7 @@ static int try_join(struct window *win, int sn) /* Get the appropriately transformed cell from the last crystal in this * series */ cr = win->img[sp].crystals[win->ser[sn][sp]]; - ref = cell_transform(crystal_get_cell(cr), win->mat[sn][sp]); + ref = cell_transform_intmat(crystal_get_cell(cr), win->mat[sn][sp]); for ( j=0; jimg[win->join_ptr].n_crystals; j++ ) { Crystal *cr2; diff --git a/tests/transformation_check.c b/tests/transformation_check.c index 2eb01aa6..c5dad16c 100644 --- a/tests/transformation_check.c +++ b/tests/transformation_check.c @@ -176,11 +176,11 @@ static int check_transformation(UnitCell *cell, IntegerMatrix *tfn, STATUS("-----------------------\n"); if ( ct == NULL ) { - cnew = cell_transform(cell, tfn); + cnew = cell_transform_intmat(cell, tfn); } else { cnew = ct; } - cback = cell_transform_inverse(cnew, tfn); + cback = cell_transform_intmat_inverse(cnew, tfn); STATUS("----> Before transformation:\n"); cell_print(cell); @@ -257,7 +257,7 @@ static int check_identity(UnitCell *cell, IntegerMatrix *tfn) STATUS("-----------------------\n"); - cnew = cell_transform(cell, tfn); + cnew = cell_transform_intmat(cell, tfn); STATUS("----> Before identity transformation:\n"); cell_print(cell); -- cgit v1.2.3