aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel/src
diff options
context:
space:
mode:
authorThomas White <taw@physics.org>2012-09-27 12:28:57 +0200
committerThomas White <taw@physics.org>2012-10-02 15:02:12 +0200
commit487f375b726fd555ff0c7eab40499201fc2e1d7d (patch)
tree28721dbbada5e27e770a60fbae2cb45d0ec1f8f1 /libcrystfel/src
parent7b31c8f41732b6ca5402afd8fe6cb4d6e185bc48 (diff)
Finishing tweaks for uncentering stuff
Today, I have mostly been having my life made difficult by the PDB's invention of "H centering".
Diffstat (limited to 'libcrystfel/src')
-rw-r--r--libcrystfel/src/cell-utils.c128
-rw-r--r--libcrystfel/src/cell-utils.h5
-rw-r--r--libcrystfel/src/cell.c11
3 files changed, 109 insertions, 35 deletions
diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c
index 2f966466..f854efe0 100644
--- a/libcrystfel/src/cell-utils.c
+++ b/libcrystfel/src/cell-utils.c
@@ -171,21 +171,19 @@ void cell_print(UnitCell *cell)
double a, b, c, alpha, beta, gamma;
double ax, ay, az, bx, by, bz, cx, cy, cz;
LatticeType lt;
+ char cen;
lt = cell_get_lattice_type(cell);
+ cen = cell_get_centering(cell);
- STATUS("%s %c", str_lattice(lt), cell_get_centering(cell));
+ STATUS("%s %c", str_lattice(lt), cen);
- switch ( lt )
+ if ( (lt==L_MONOCLINIC) || (lt==L_TETRAGONAL) || ( lt==L_HEXAGONAL)
+ || ( (lt==L_ORTHORHOMBIC) && (cen=='A') )
+ || ( (lt==L_ORTHORHOMBIC) && (cen=='B') )
+ || ( (lt==L_ORTHORHOMBIC) && (cen=='C') ) )
{
- case L_MONOCLINIC :
- case L_TETRAGONAL :
- case L_HEXAGONAL :
STATUS(", unique axis %c", cell_get_unique_axis(cell));
- break;
-
- default :
- break;
}
if ( right_handed(cell) ) {
@@ -263,6 +261,11 @@ int bravais_lattice(UnitCell *cell)
return 0;
case 'H' :
+ /* "Hexagonal H" is not a Bravais lattice, but rather something
+ * invented by the PDB to make life difficult for programmers.
+ * Accepting it as Bravais seems to be the least painful way to
+ * handle it correctly. Yuk. */
+ if ( ua != 'c' ) return 0;
if ( lattice == L_HEXAGONAL ) return 1;
return 0;
@@ -294,24 +297,28 @@ static UnitCellTransformation *uncentering_transformation(UnitCell *in,
t = tfn_identity();
if ( t == NULL ) return NULL;
- if ( (ua == 'a') || (cen == 'A') ) {
+ if ( ua == 'a' ) {
tfn_combine(t, tfn_vector(0,0,1),
tfn_vector(0,1,0),
tfn_vector(-1,0,0));
- if ( cen == 'A' ) cen = 'C';
}
- if ( (ua == 'b') || (cen == 'B') ) {
+ if ( ua == 'b' ) {
tfn_combine(t, tfn_vector(1,0,0),
tfn_vector(0,0,1),
tfn_vector(0,-1,0));
- if ( cen == 'B' ) cen = 'C';
}
switch ( cen ) {
case 'P' :
+ *new_latt = lt;
+ *new_centering = 'P';
+ break;
+
case 'R' :
+ *new_latt = L_RHOMBOHEDRAL;
+ *new_centering = 'R';
break;
case 'I' :
@@ -335,7 +342,6 @@ static UnitCellTransformation *uncentering_transformation(UnitCell *in,
if ( lt == L_CUBIC ) {
*new_latt = L_RHOMBOHEDRAL;
*new_centering = 'R';
-
} else {
assert(lt == L_ORTHORHOMBIC);
*new_latt = L_TRICLINIC;
@@ -343,6 +349,8 @@ static UnitCellTransformation *uncentering_transformation(UnitCell *in,
}
break;
+ case 'A' :
+ case 'B' :
case 'C' :
tfn_combine(t, tfn_vector(H,H,0),
tfn_vector(-H,H,0),
@@ -352,6 +360,7 @@ static UnitCellTransformation *uncentering_transformation(UnitCell *in,
break;
case 'H' :
+ /* Obverse setting */
tfn_combine(t, tfn_vector(TT,OT,OT),
tfn_vector(-OT,OT,OT),
tfn_vector(-OT,-TT,OT));
@@ -366,18 +375,20 @@ static UnitCellTransformation *uncentering_transformation(UnitCell *in,
}
- if ( ua == 'a' ) {
- tfn_combine(t, tfn_vector(0,0,-1),
- tfn_vector(0,1,0),
- tfn_vector(1,0,0));
- if ( cen == 'C' ) cen = 'A';
- }
+ /* Reverse the axis permutation, but only if this was not an H->R
+ * transformation */
+ if ( !((cen=='H') && (*new_latt == L_RHOMBOHEDRAL)) ) {
+ if ( ua == 'a' ) {
+ tfn_combine(t, tfn_vector(0,0,-1),
+ tfn_vector(0,1,0),
+ tfn_vector(1,0,0));
+ }
- if ( ua == 'b' ) {
- tfn_combine(t, tfn_vector(1,0,0),
- tfn_vector(0,0,-1),
- tfn_vector(0,1,0));
- if ( cen == 'C' ) cen = 'B';
+ if ( ua == 'b' ) {
+ tfn_combine(t, tfn_vector(1,0,0),
+ tfn_vector(0,0,-1),
+ tfn_vector(0,1,0));
+ }
}
return t;
@@ -406,6 +417,7 @@ UnitCell *uncenter_cell(UnitCell *in, UnitCellTransformation **t)
if ( !bravais_lattice(in) ) {
ERROR("Cannot uncenter: not a Bravais lattice.\n");
+ cell_print(in);
return NULL;
}
@@ -481,11 +493,11 @@ UnitCell *match_cell(UnitCell *cell_in, UnitCell *template_in, int verbose,
float angtol = deg2rad(tols[3]);
UnitCell *cell;
UnitCell *template;
- UnitCellTransformation *to_given_cell;
+ UnitCellTransformation *uncentering;
UnitCell *new_cell_trans;
/* "Un-center" the template unit cell to make the comparison easier */
- template = uncenter_cell(template_in, &to_given_cell);
+ template = uncenter_cell(template_in, &uncentering);
/* The candidate cell is also uncentered, because it might be centered
* if it came from (e.g.) MOSFLM */
@@ -671,7 +683,7 @@ UnitCell *match_cell(UnitCell *cell_in, UnitCell *template_in, int verbose,
free(cand[2]);
/* Reverse the de-centering transformation */
- new_cell_trans = cell_transform_inverse(new_cell, to_given_cell);
+ new_cell_trans = cell_transform_inverse(new_cell, uncentering);
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));
@@ -1115,10 +1127,13 @@ int cell_is_sensible(UnitCell *cell)
* lattice is a conventional Bravais lattice.
* Warnings are printied if any of the checks are failed.
*
+ * Returns: true if cell is invalid.
+ *
*/
-void validate_cell(UnitCell *cell)
+int validate_cell(UnitCell *cell)
{
int err = 0;
+ char cen, ua;
if ( !cell_is_sensible(cell) ) {
ERROR("Warning: Unit cell parameters are not sensible.\n");
@@ -1136,5 +1151,58 @@ void validate_cell(UnitCell *cell)
err = 1;
}
- if ( err ) cell_print(cell);
+ cen = cell_get_centering(cell);
+ ua = cell_get_unique_axis(cell);
+ if ( (cen == 'A') && (ua != 'a') ) {
+ ERROR("Warning: centering doesn't match unique axis.\n");
+ err = 1;
+ }
+ if ( (cen == 'B') && (ua != 'b') ) {
+ ERROR("Warning: centering doesn't match unique axis.\n");
+ err = 1;
+ }
+ if ( (cen == 'C') && (ua != 'c') ) {
+ ERROR("Warning: centering doesn't match unique axis.\n");
+ err = 1;
+ }
+
+ return err;
+}
+
+
+/**
+ * forbidden_reflection:
+ * @cell: A %UnitCell
+ * @h: h index to check
+ * @k: k index to check
+ * @l: l index to check
+ *
+ * Returns: true if this reflection is forbidden.
+ *
+ */
+int forbidden_reflection(UnitCell *cell,
+ signed int h, signed int k, signed int l)
+{
+ char cen;
+
+ cen = cell_get_centering(cell);
+
+ /* Reflection conditions here must match the transformation matrices
+ * in uncentering_transformation(). tests/centering_check verifies
+ * this (amongst other things). */
+
+ if ( cen == 'P' ) return 0;
+ if ( cen == 'R' ) return 0;
+
+ if ( cen == 'A' ) return (k+l) % 2;
+ if ( cen == 'B' ) return (h+l) % 2;
+ if ( cen == 'C' ) return (h+k) % 2;
+
+ if ( cen == 'I' ) return (h+k+l) % 2;
+ if ( cen == 'F' ) return ((h+k) % 2) || ((h+l) % 2) || ((k+l) % 2);
+
+ /* Obverse setting */
+ if ( cen == 'H' ) return (-h+k+l) % 3;
+
+ return 0;
}
diff --git a/libcrystfel/src/cell-utils.h b/libcrystfel/src/cell-utils.h
index 8acf2f85..f92ab22d 100644
--- a/libcrystfel/src/cell-utils.h
+++ b/libcrystfel/src/cell-utils.h
@@ -53,7 +53,7 @@ extern UnitCell *load_cell_from_pdb(const char *filename);
extern int cell_is_sensible(UnitCell *cell);
-extern void validate_cell(UnitCell *cell);
+extern int validate_cell(UnitCell *cell);
extern UnitCell *uncenter_cell(UnitCell *in, UnitCellTransformation **tr);
@@ -63,4 +63,7 @@ extern int right_handed(UnitCell *cell);
extern const char *str_lattice(LatticeType l);
+extern int forbidden_reflection(UnitCell *cell,
+ signed int h, signed int k, signed int l);
+
#endif /* CELL_UTILS_H */
diff --git a/libcrystfel/src/cell.c b/libcrystfel/src/cell.c
index 7a62f51c..e7c9dace 100644
--- a/libcrystfel/src/cell.c
+++ b/libcrystfel/src/cell.c
@@ -260,6 +260,9 @@ UnitCell *cell_new_from_cell(UnitCell *orig)
cell_get_cartesian(orig, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz);
cell_set_cartesian(new, ax, ay, az, bx, by, bz, cx, cy, cz);
cell_set_pointgroup(new, orig->pointgroup);
+ cell_set_lattice_type(new, orig->lattice_type);
+ cell_set_centering(new, orig->centering);
+ cell_set_unique_axis(new, orig->unique_axis);
return new;
}
@@ -691,12 +694,12 @@ UnitCell *cell_transform(UnitCell *cell, UnitCellTransformation *t)
if ( t == NULL ) return NULL;
- out = cell_new();
+ out = cell_new_from_cell(cell);
if ( out == NULL ) return NULL;
- cell_get_cartesian(cell, &ax, &ay, &az,
- &bx, &by, &bz,
- &cx, &cy, &cz);
+ cell_get_cartesian(out, &ax, &ay, &az,
+ &bx, &by, &bz,
+ &cx, &cy, &cz);
m = gsl_matrix_alloc(3,3);
a = gsl_matrix_calloc(3,3);