CHAPTER 12

Unit Cell and Coordinate Systems

 

Although we mentioned slightly about the volume of a unit cell and polar orientation, we have not seriously involved crystallographic topics yet.  With the tools introduced in the previous chapters, this chapter provides an ¡®orientation¡¯ for the rest of the crystallographic concepts.  Several coordinate systems are used in crystallography, three of which are directly related to the concept of unit cell.  Fractional coordinate system is often used for symmetry operation (See Chapter 13).  Cartesian or orthogonal coordinate system services as an interface for crystallography to communicate with other scientific fields.  The coordinate system in the reciprocal space is heavily used without a name.  Miller indices are the coordinates in the discrete reciprocal space.  See Chapter 13 for details.

 

12.1 ccl::cell Class

 

The class ccl::cell implements the concept of unit cell.  It can be either a direct or a reciprocal cell.  This class itself contains the information about the type, and provides functions for the inter-conversion between direct and reciprocal cells.  This class also provides matrices for conversions between various coordinate systems.

 

12.1.1 Constructors, assignment and parameter access

 

cell(const double a     = 0,

     const double b     = 0,

     const double c     = 0,

     const double alpha = 0,

     const double beta  = 0,

     const double gamma = 0, const bool d = 1);

cell(const mac::vector3D<double> l,

     const mac::vector3D<double> a, const bool d = 1);

The general constructors take either two mac::vector3D<double>s and a Boolean, or six doubles and a Boolean as the cell lengths a, b, c and angles a, b, g.  The Boolean indicates whether the cell is direct.  The default constructor constructs a cell with 0-lengths and 0-angles, which will cause assertion error if further calculation is carried out.  This will prevent unintended result.

 

cell(const cell& c);

A copy constructor is available.

 

const cell& operator=(const cell& c);

The assignment operator = is defined.

 

void setA(const double a);

void setB(const double b);

void setC(const double c);

void setAlpha(const double a);

void setBeta(const double b);

void setGamma(const double g);

void setLengths(const mac::vector3D<double> l);

void setAngles(const  mac::vector3D<double> a);

void setCell(const double a,

             const double b,

             const double c,

             const double alpha,

             const double beta,

             const double gamma);

void setCell(const mac::vector3D<double> l,

             const mac::vector3D<double> a);

Each individual cell parameter, the cell lengths, angles, or all six parameters can be set by access functions.

 

double getA()     const;

double getB()     const;

double getC()     const;

double getAlpha() const;

double getBeta()  const;

double getGamma() const;

mac::vector3D<double> getLengths() const;

mac::vector3D<double> getAngles()  const;

Return individual cell parameter, the cell lengths, and angles.

 

isDirect()

isReciprocal()

Return Boolean that tests the cell type.

 

12.1.2 Cell volume and comparison

 

static double volume(const cell&);

double        volume() const;

The cell volume is calculated by a static and a member function volume().

 

bool operator==(const cell&) const;

bool operator!=(const cell&) const;

bool operator&=(const cell&) const;

The operators ==, &= and != test absolutely equal, approximately equal, and not equal between two cells, respectively.

 

bool operator< (const cell&) const;

bool operator<=(const cell&) const;

bool operator> (const cell&) const;

bool operator>=(const cell&) const;

The comparison operator < is defined by comparing the cell volumes.  If two cell volumes are equal, use the comparison of cell lengths.  If the cell lengths are all equal, use the comparison of cell angles.  The other comparison operators >, <=, and >= are defined accordingly.

 

12.1.3 Conversion between direct and reciprocal cells

 

Every cell has its reciprocal cell.  The reciprocal cell of a direct cell is, of course, called reciprocal cell.  The reciprocal cell of a reciprocal cell becomes the direct cell again.

 

static cell reciprocal(const cell&,

                       mac::matrix3D<double>* = NULL);

static cell reciprocalScalar(const cell&);

Two static functions reciprocal() and reciprocalScalar() are defined to convert a cell to its reciprocal cell.  The former is implemented in a vector form.  The latter is a scalar form and obsolete, but it serves as an example to show a scalar implementation.  The vector form reciprocal() may also modify a mac::matrix3D<double>, if a valid pointer points to it when the vector form is called.  The matrix will store three vectors by row, which can be accessed by using the operator [] of mac::matrix3D<double>.  These vectors are a*, b*, and c* in the orthogonal coordinate system of the cell.  See below for details of coordinate systems.

 

cell reciprocal(mac::matrix3D<double>* = NULL);

The member function reciprocal() converts the cell itself to its reciprocal cell, and optionally modifies a matrix if a valid pointer is given.

 

cell operator*() const;

cell operator~() const;

Two unary operators * and ~ have the exactly same definition that returns the reciprocal cell, and leaves the original cell unchanged.

 

12.1.4 Coordinate systems and orthogonalization

 

The three cell edges a, b, c naturally form a coordinate system, often called crystallographic coordinate system.  This system could be orthogonal for some crystal systems, but it is not in general.  Either the absolute length in Å or the three cell lengths a, b, c can be the units.  If the cell lengths are used as units for their axes, this crystallographic coordinate system is a fractional coordinate system.  Symmetry operation from one atomic position to its symmetry-related positions is very convenient in the fractional system.

 

However, arbitrary rotation is difficult using fractional coordinate system.  An orthogonal system is therefore required.  It is entirely a convention in the field to choose one unique orthogonal coordinate system respect to the crystallographic coordinate systems.  The unique orthogonal system also serves as an interface when structural results are communicated to other scientific fields.  We therefore use the definition of orthogonal system from the Protein Data Bank (PDB; http://www.rcsb.org/pdb).  The PDB convention is to select a-axis as the x-axis, a ´ b as the z-axis, and of course z ´ x as y-axis.  Considering the relationship between direct and reciprocal cells, this is also to say that c* is selected as z-axis, b* ¡Á c* as x-axis, and z ¡Á x as y-axis.  Sometimes other convention may also be used.  For example, Rossmann convention selects a as x-axis, c ¡Á a as y-axis, and x ¡Á y as z-axis, that is equivalent to b* as y-axis, b* ¡Á c* as x-axis, and x ¡Á y as z-axis.

 

static mac::matrix3D<double>   orthogonal(

   const cell&, const bool = 1);

static mac::matrix3D<double> deorthogonal(

   const cell&, const bool = 1);

Two static functions orthogonal() and deorthogonal() take a cell and an optional Boolean, and return the orthogonalization and deorthogonalization matrices, respectively.  The Boolean specifies whether PDB convention is used.  The default is certainly the PDB convention.  If the Boolean is unset, Rossmann convention will be used.  The orthogonalization matrix pre-multiplies a vector in crystallographic coordinate system, and results in a vector in the orthogonal system.  The deorthogonalization matrix is the inverse of orthogonalization matrix.  The required argument of these functions is ccl::cell, which can be either a direct or reciprocal cell.  If it is a direct cell, the returned matrices can ONLY pre-multiply a vector in the direct coordinate system.  It is an error to pre-multiply vectors in the reciprocal system.  USE CAUTION.

 

mac::matrix3D<double>   orthogonal(

   const bool = 1) const;

mac::matrix3D<double> deorthogonal(

   const bool = 1) const;

Two member functions calculate the orthogonalization and deorthogonalization matrices, respectively.  If the optional Boolean is set, PDB convention will used; otherwise, use Rossmann convention.

 

mac::matrix3D<double> operator+() const;

mac::matrix3D<double> operator-() const;

The unary operators + and ¨C also return these matrices, but under the PDB convention only.

 

12.1.5 Transformation between direct and reciprocal coordinate systems

 

The previous section is about the transformation between a crystallographic coordinate system, direct or reciprocal, and an orthogonal system.  This section describes the transformation between two crystallographic coordinate systems.  Transformation of one location in a system to the other does not make much sense, since these two systems are reciprocal to each other and do not share the same unit.  However, the relative orientation between these two systems are fixed.  Transformation of one orientation in one system to the other is useful in the process of indexing diffraction patterns, for example.

 

static mac::matrix3D<double> directToReciprocal(cell);

static mac::matrix3D<double> reciprocalToDirect(cell);

mac::matrix3D<double> directToReciprocal() const;

mac::matrix3D<double> reciprocalToDirect() const;

Static and member functions directToReciprocal() and reciprocalToDirect() are available.  They return matrices that can be pre-multiplied to vectors in direct and reciprocal systems.  The static functions require an argument of direct or reciprocal cell.  The type of the cell does not affect the result.

 

12.1.6 Output

 

The operator << is overloaded to print the cell parameters to an ostream.

 

The following program tests the class ccl::cell.

 

#include <assert.h>

#include <CCL/cell.h>

 

main()

{

   ccl::cell c(40, 30, 80, 90, 90, 120); // Construct a direct cell

   ccl::cell d = c; // Assignment

   assert (d == c); // Absolutely equal

   assert (d &= c); // Approximately equal

   c.setA(50);      // Set a

   assert (d != c); // Not equal

 

   assert (c.isDirect());     // c is a direct cell.

   cout << c << endl << endl;

   c.reciprocal();

   assert (c.isReciprocal()); // c becomes a reciprocal cell now.

   cout << c << endl << endl;

 

   assert (ccl::cell::reciprocal(d) &= ccl::cell::reciprocalScalar(d));

      // The vector and scalar versions.  The scalar version is obsolete.

   assert (*d == ~d); // The unary operators * and ~ are the same.

 

   c = ~d; // c is the reciprocal cell of d.

   mac::vector3D<double> v(1, 1);

 

   // Orthogonalization

   assert (+d * v != +c * v);

      // Orthogonalization matrix +d can only pre-multiply a vector in direct

      // space.  Orthogonalization matrix +c pre-multiplies a numerically same

      // vector v, but now v is considered in reciprocal space, since c is the

      // reciprocal cell.

 

   assert (d.reciprocalToDirect() * (d.directToReciprocal() * v) &= v);

   assert (c.reciprocalToDirect() * (d.directToReciprocal() * v) &= v);

   assert (d.reciprocalToDirect() * (c.directToReciprocal() * v) &= v);

   assert (c.reciprocalToDirect() * (c.directToReciprocal() * v) &= v);

   assert (d.reciprocalToDirect() * (c.directToReciprocal() * v) &= v);

   assert (c.reciprocalToDirect() * (d.directToReciprocal() * v) &= v);

      // Transform vector v from direct space to reciprocal and back to

      // direct.  It does not make a difference whether a direct or reciprocal

      // cell is used.

 

   return 0;

}

 

Listing 12.1.0.0.1 Testing ccl:cell.

 

12.2 The Shadow Class cpl.ccl.symm.cell

 

The shadow class cpl.ccl.symm.cell has almost all functions of the original class ccl::cell.

 

12.2.1 Constructors and parameter access

 

cell(a = 0, b = 0, c = 0,

     alpha = 0, beta = 0, gamma = 0, direct = 1)

Constructs a cell object.  direct indicates whether the cell is direct cell.  If no arguments are given, the default constructor returns a cell with 0-lengths and 0-angles.  Further calculation based on the default cell may cause assertion error.

 

cell(len, ang, direct = 1)

Constructs a cell object with two cpl.mac.vm.vector3DDBLs being cell lengths and angles.  If the Boolean direct is set, returns a direct cell; otherwise a reciprocal cell.

 

If cl is a cell object, the following functions are defined:

 

cell(cl)

Constructs a new cell object identical to cl.

 

cl.setA(a)

cl.setB(b)

cl.setC(c)

cl.setAlpha(alpha)

cl.setBeta(beta)

cl.setGamma(gamma)

cl.setLengths(len)

cl.setAngles(ang)

cl.setCell(a, b, c, alpha, beta, gamma)

cl.setCell(len, ang)

Set or change the cell parameters of cl.  len and ang are cpl.mac.vm.vector3DDBL, the others are all floats.

 

cl.getA()

cl.getB()

cl.getC()

cl.getAlpha()

cl.getBeta()

cl.getGamma()

cl.getLengths()

cl.getAngles()

Return the cell parameters of cl.  getLengths and getAngles return cpl.mac.vm.vector3DDBL, the others return floats.

 

cl.isDirect()

cl.isReciprocal()

Return Booleans that test whether cl is direct or reciprocal cell.

 

12.2.2 Cell volume and comparison

 

cell.volume(cl)

cl.volume()

The static and member functions compute the cell volume.

 

cl == c

cl.ap_eq(c)

cl != c

Return true if cl is absolutely, approximately, or not equal to another cell c, respectively.

 

cl < c

Returns true if the volume of cl is less than that of c.  If these cells are equal in volume, returns true if the vector of cl cell lengths is less than that of c.  If the cell length vectors are also equal, returns true if the vector of cl cell angles is less than that of c.  Returns false otherwise.

 

cl <= c

Returns true if cl < c or cl == c is true.

 

cl > c

Returns true if cl <= c is false.

 

cl >= c

Returns true if cl > c or cl == c is true.

 

12.2.3 Conversion between direct and reciprocal cells

 

cell_reciprocalScalar(cl)

cell.reciprocalScalar(cl)

Computes the reciprocal cell of cl by a scalar implementation.  This static function is obsolete.

 

cell_recip(cl, m = Null)

cell.recip(cl, m = Null)

Computes the reciprocal cell of cl by a vector implementation.  If m is a valid cpl.mac.vm.matrix3DDBL, it is modified to store vectors a*, b*, c* by rows in the orthogonal coordinate system of the cell cl.

 

cl.reciprocal(m = Null)

Converts cl to its reciprocal cell and returns the reciprocal cell.  See above for m.

 

~cl

Returns the reciprocal cell of cl.  cl remains unchanged.

 

12.2.4 Orthogonalization

 

cell.orthogonal(cl, pdb = 1)

cell.deorthogonal(cl, pdb = 1)

Return the orthogonalization and deorthogonalization matrices, respectively.  cl can be either a direct or a reciprocal cell.  If cl is a direct cell, the returned matrices can be used to pre-multiply vectors in the direct space only.  It will result in unintended solution if vectors in reciprocal space are pre-multiplied by these matrices.  Use a reciprocal cell cl to generate matrices for orthogonalization or deorthogonalization of vectors in reciprocal space.  pdb indicates whether the PDB convention is used.  If it is unset, Rossmann convention will be used.

 

cl.orthogonal(pdb = 1)

cl.deorthogonal(pdb = 1)

Return the orthogonalization and deorthogonalization matrices, respectively.  See above.

 

+cl

-cl

Return the orthogonalization and deorthogonalization matrices, respectively.  Always use the PDB convention.  See above.

 

12.2.5 Transformation between direct and reciprocal coordinate systems

 

cell_directToReciprocal(cl)

cell_reciprocalToDirect(cl)

Return matrices that can transform between direct and reciprocal coordinate systems.  cl can be either direct or reciprocal cell, and does not affect the results.

 

cl.directToReciprocal()

cl.reciprocalToDirect()

Same as above.

 

12.2.6 Output

 

cl.__str__()

The special function __str__ is added, so that print cl in Python will print out the cell parameters.

 

This example tests the shadow class cpl.ccl.symm.cell.  This program shall print a direct cell and its reciprocal cell, as well as the vectors a*, b*, c* in the PDB orthogonal coordinate system.

 

#!/usr/bin/python

 

import cpl

 

C = cpl.ccl.symm.cell

 

c = C(30, 50, 40, 90, 110, 90)  # Constructs a direct cell

d = cpl.ccl.symm.cell(c)        # Copy constructor

assert d == c                   # Absolutely equal

assert d.ap_eq(c)               # Approximately equal

c.setA(60)                      # Set a

assert d != c                   # Not equal

assert c.isDirect()             # c is a direct cell.

print c                         # Output

c.reciprocal()                  # Reciprocalize

print c

assert c.isReciprocal()         # c becomes a reciprocal cell now.

d = ~c # The unary operator ~.  Reciprocal of reciprocal.

m = cpl.mac.vm.matrix3DDBL()

assert cpl.ccl.symm.cell_recip(d, m).ap_eq(\

       cpl.ccl.symm.cell_reciprocalScalar(d))

assert cpl.ccl.symm.cell.recip(d, m).ap_eq(\

       cpl.ccl.symm.cell.reciprocalScalar(d))

   # The vector and scalar versions.  The scalar version is obsolete.

print 'a* =', m(1)

print 'b* =', m(2)

print 'c* =', m(3)

 

# Orthogonalization

v = cpl.mac.vm.vector3DDBL(1, 1)

assert +d * v != +c * v

   # Orthogonalization matrix +d can only pre-multiply a vector in direct

   # space.  Orthogonalization matrix +c pre-multiplies a numerically same

   # vector v, but now v is considered in reciprocal space, since c is the

   # reciprocal cell.

assert d.orthogonal() * v != c.orthogonal() * v

assert d.deorthogonal() * v != c.deorthogonal() * v

assert cpl.ccl.symm.cell.orthogonal(d) * v !=\

       cpl.ccl.symm.cell.orthogonal(c) * v

assert cpl.ccl.symm.cell.deorthogonal(d) * v !=\

       cpl.ccl.symm.cell.deorthogonal(c) * v

 

assert v.ap_eq(d.reciprocalToDirect() * (d.directToReciprocal() * v))

assert v.ap_eq(c.reciprocalToDirect() * (d.directToReciprocal() * v))

assert v.ap_eq(d.reciprocalToDirect() * (c.directToReciprocal() * v))

assert v.ap_eq(c.reciprocalToDirect() * (c.directToReciprocal() * v))

assert v.ap_eq(d.reciprocalToDirect() * (c.directToReciprocal() * v))

assert v.ap_eq(c.reciprocalToDirect() * (d.directToReciprocal() * v))

   # Transform vector v from direct space to reciprocal and back to

   # direct.  It does not make a difference whether a direct or reciprocal

   # cell is used.

 

Listing 12.2.0.0.1 Testing cpl.ccl.symm.cell.

 

12.3 Coordinates

 

A location in 3-dimensional space is described by three coordinates.  The numerical values of the coordinates depend on coordinate system and unit, so that one location may have different sets of coordinates.  Conversion between these sets is therefore required.  CCL implements three types of the most frequently-used coordinates and the inter-conversion between them: grid coordinates, fractional coordinates, and Cartesian coordinates.

 

12.3.1 Grid coordinates

 

If a unit cell length is divided into a number of grids, any location in the unit cell can be represented by the closest grid point.  The grid numbers along three cell edges are the grid coordinates.  Electron density map usually has the density values only on the grids, since it is often calculated by discrete Fourier transform.  Grid coordinates can only be integer.  A negative grid can be translated one or several unit cells to be converted into a non-negative grid.

 

12.3.1.1 ccl::gridCoordinates

The class ccl::gridCoordinates is derived from mac::vector3D<short>.

 

Constructors and destructor

gridCoordinates(const int x = 0,

                const int y = 0,

                const int z = 0,

                const unsigned short u = 1,

                const unsigned short v = 1,

                const unsigned short w = 1,

                const bool unity = 1);

The general constructor takes six integers.  The first three are the grid coordinates, the last three are the total grid numbers that all three cell edges are divided into.  The total grid numbers are also the cell dimensions in grid instead of in Angstrom.  By default, the constructor returns the origin.  However, the cell dimensions are set to 1, which is usually not intended.  It is not a good idea to use the default constructor to generate the origin.  If the Boolean unity is set, this constructor returns grid coordinates in the first unit cell.

 

gridCoordinates(

   const mac::vector3D<int>&,

   const mac::vector3D<unsigned short>&,

   const bool unity = 1);

Obviously, the first mac::vector3D<int> is the coordinates, and the mac::vector3D<unsigned short> is the cell dimensions.

 

gridCoordinates(

   const gridCoordinates&,

   const bool = 0);

gridCoordinates(

   const fractionalCoordinates&,

   const mac::vector3D<unsigned short>&,

   const bool = 1);

gridCoordinates(

   const CartesianCoordinates&,

   const mac::vector3D<unsigned short>&,

   const cell&,

   const bool = 1);

These three are the copy constructors.  The first is the ordinary one that copies one gridCoordinates to another.  The other two copy fractional coordinates and Cartesian coordinates to grid coordinates, respectively.  Both of these require the total grid dimensions, and the last copy constructor requires the unit cell.  See 12.3.2 and 12.3.3 for more detail.

 

All constructors including the copy constructor may accept an optional Boolean unity.  All defaults of this argument are true, except for the copy constructor.  By default, all constructors return gridCoordinates in the first unit cell, while the copy constructor does not translate the location of the coordinates unless explicitly indicated.

 

const gridCoordinates& operator=(

   const gridCoordinates&);

Assignment operator = is defined.

 

operator mac::vector3D<int>() const;

Converts to mac::vector3D<int>.

 

~gridCoordinates();

There is nothing to do in destructor.

 

Equivalent positions

gridCoordinates equivalent(const spaceGroup&,

                           const int i) const;

Returns the i-th equivalent position, where i is a 1-offset subscript.  The function requires the space group information.  i must be from 1 to the number of symmetry operators of the space group.  See the next chapter for space group.

 

vector<gridCoordinates> equivalent(

   const spaceGroup&,

   int* multiplicity = NULL) const;

Returns all equivalent positions in a vector<gridCoordinates>.  If the last pointer is valid, the multiplicity of this grid coordinates will also be returned.  This function obviously requires the symmetry information of the space group.  See the next chapter for space group.

 

vector<gridCoordinates> equivalentAU(

   const spaceGroup&,

   int* multiplicity = NULL);

This function has the exact same interface as the previous one, except this is a non-const version.  This function returns the same as the previous one, but after return, the original grid coordinates will be replaced by the equivalent position in the default asymmetric unit.  The default asymmetric unit is defined by the volume in which a gridCoordinates is less than any of its equivalent ones.  See 3.1.1.19 for vector comparison.

 

int multiplicity(const spaceGroup&) const;

Returns the multiplicity of this grid coordinates, given a space group.

 

Cell dimensions and grid distance

mac::vector3D<unsigned short>

   getDimensions() const;

Returns the cell dimensions in grid.

 

double gridDistance(

   const gridCoordinates& gc = gridCoordinates())

   const;

Returns grid distance from this gridCoordinates to the given one gc.  If no one is given, returns the grid distance to the origin.  Please note that grid distance is not the normal sense distance.  Its unit is grid.  Greater value may not mean longer distance.  Avoid using this function if possible or use caution.  This function may be removed or put into protected members in future release.

 

Output

friend ostream &operator<<(

   ostream&, const gridCoordinates&);

Operator << is defined to return an ostream.

 

#include <assert.h>

#include <ccl/coordinates.h>

 

main()

{

   typedef ccl::gridCoordinates       GC;

   typedef ccl::fractionalCoordinates FC;

   typedef ccl::CartesianCoordinates  CC;

   typedef mac::vector3D<unsigned short> V;

 

   GC gc; // Construct a default grid coordinates

   assert (!gc);

   assert (gc.getDimensions() == V(1, 1, 1));

   gc = GC(9, -28, 20, 30, 40, 50);

   assert (gc.gridDistance() == 25);

      // Grid distance is not a distance in the normal sense.

 

   const ccl::spaceGroup c2("c2");

   int m;

   const vector<GC> eq = gc.equivalent(c2, &m);

   assert (GC(V(9, 12, 20), V(30, 40, 50)) == gc); // Didn't change gc.

   assert (m == 1);                                // Multiplicity is 1.

   for (int i = 0; i < eq.size(); i++) cout << eq[i] << endl;

 

   gc.equivalentAU(c2, &m);

   assert (GC(V(9, 12, 20), V(30, 40, 50)) != gc); // gc is changed in AU.

 

   GC(0, 12, 25, 30, 40, 50).equivalent(c2, &m);

   assert (m == 2);                                // Multiplicity is 2.

   GC(15, 12, 25, 30, 40, 50).equivalent(c2, &m);

   assert (m == 2);                                // Multiplicity is 2.

   GC(0, 0, 0, 30, 40, 50).equivalent(c2, &m);

   assert (m == 2);                                // Multiplicity is 2.

 

   const ccl::cell cl(62, 80, 95, 90, 100, 90);

   FC fc = FC(0.5, 0.3, 0.8);

   assert (GC(fc, V(30, 40, 50)) == GC(15, 12, 40, 30, 40, 50)); // FC -> GC

   assert (FC(GC(15, 12, 40, 30, 40, 50)) &= fc);                // GC -> FC

   assert (GC(CC(gc, cl), V(30, 40, 50), cl) == gc);             // GC <-> CC

 

   return 0;

}

 

Listing 12.3.1.1.1 Testing ccl::gridCoordinates.

 

12.3.1.2 cpl.ccl.symm.gridCoordinates

gc, fc, and cc are grid, fractional, and Cartesian coordinates in this section, respectively.

 

Constructors

gridCoordinates(x = 0, y = 0, z = 0,

                u = 1, v = 1, w = 1, unity = 1)

Constructs a gridCoordinates from three integers x, y, z, and three unsigned shorts u, v, w as the total grid numbers in three cell dimensions.  There is no numerical type unsigned short in Python, so that u, v, w appear to be integer in Python, but they will be unsigned short when mapped into C++.  For example, if u is mistakenly assigned ¨C30, it will be 65506 in C++.  If u is assigned 70000, it will be 4464.  No error checking is provided in the current release.  The default constructor sets all grid dimensions to 1, which does not represent the origin properly.  Do not use the default as the origin.

 

gridCoordinates(gc, unity = 0)

Copies gc to another gridCoordinates.

 

gridCoordinates(u, v, unity = 1)

gridCoordinates(fc, v, unity = 1)

gridCoordinates(cc, v, cl, unity = 1)

Convert vector3DDBL u, fc and cc to gridCoordinates, respectively.  v is a cpl.mac.vm.vector3DSUI that contains the grid dimensions.  cl is a cell.  If unity is true, the returned gridCoordinates is kept in the first unit cell.

 

gc.vector3DINT()

Converts gc to cpl.mac.vm.vector3DINT.

 

Equivalent positions

gc.equivalent(sg, i)

Given a space group sg, returns the i-th equivalent position, where i is a 1-offset subscript between 1 and the number of symmetry operators of the space group.

 

gc.euivalent(sg)

Given a space group sg, returns a tuple of all equivalent positions in a unit cell.

 

gc.equivalentAU(sg)

Given a space group sg, replaces gc by its equivalent position in the default asymmetric unit.  See 12.3.1.1 for definition of default asymmetric unit.  This function also returns all equivalent positions in a unit cell.

 

gc.multiplicity(sg)

Returns the multiplicity of the gridCoordinates gc, if the symmetry is given by space group sg.

 

Cell dimensions and grid distance

gc.getDimensions()

Returns the total grid numbers in three cell dimensions as a cpl.mac.vm.vector3DSUI.

 

gc.gridDistance(p = gridCoordinates())

Returns the grid distance from gc to the given position p.  If nothing is given, returns the grid distance to the origin.  Please note that grid distance is not the distance in normal sense.  Avoid using this function if possible.

 

Output

gc.__str__()

The special function __str__ is defined to print gridCoordinates.

 

#!/usr/bin/python

 

import cpl

 

GC = cpl.ccl.symm.gridCoordinates

FC = cpl.ccl.symm.fractionalCoordinates

CC = cpl.ccl.symm.CartesianCoordinates

V  = cpl.mac.vm.vector3DSUI

c2 = cpl.ccl.symm.spaceGroup("c2")

cl = cpl.ccl.symm.cell(62, 80, 95, 90, 100, 90)

 

gc = GC() # Construct a default grid coordinates

assert gc.is_zero()

assert gc.getDimensions() == V(1, 1, 1)

gc = GC(9, -28, 20, 30, 40, 50)

assert gc.gridDistance() == 25

   # Grid distance is not a distance in the normal sense.

for i in range(c2.operatorSize()): print gc.equivalent(c2, i + 1)

assert GC(9, 12, 20, 30, 40, 50) == gc # Didn't change gc.

gc.equivalentAU(c2)

assert GC(9, 12, 20, 30, 40, 50) != gc # gc is changed in AU.

assert gc.multiplicity(c2) == 1                         # Multiplicity is 1.

assert GC(0, 12, 25, 30, 40, 50).multiplicity(c2)  == 2 # Multiplicity is 2.

assert GC(15, 12, 25, 30, 40, 50).multiplicity(c2) == 2 # Multiplicity is 2.

assert GC(0, 0, 0, 30, 40, 50).multiplicity(c2)    == 2 # Multiplicity is 2.

 

fc = FC(0.5, 0.3, 0.8)

assert GC(fc, V(30, 40, 50)) == GC(15, 12, 40, 30, 40, 50) # FC -> GC

assert FC(GC(15, 12, 40, 30, 40, 50)).ap_eq(fc)            # GC -> FC

assert GC(CC(gc, cl), V(30, 40, 50), cl) == gc             # GC <-> CC

 

Listing 12.3.1.2.1 Testing cpl.ccl.symm.gridCoordinates.

 

12.3.1.3 Container of grid coordinates

std::vector<ccl::gridCoordinates> is instantiated, can be used as cpl.ccl.symm.vectorGDC.

 

12.3.2 Fractional coordinates

 

Grid coordinates, by definition, are only approximate coordinates.  Fractional coordinates could represent any location in a unit cell precisely.  A fractional coordinate is a floating point value from 0 to 1 with no unit.  All cell lengths are 1.

 

12.3.2.1 ccl::fractionalCoordinates

The class ccl::fractionalCoordinates is derived from mac::vector3D<double>.

 

Constructors and destructor

fractionalCoordinates(const double x = 0,

                      const double y = 0,

                      const double z = 0,

                      const int unity = 1);

Constructs a fractionalCoordinates from x, y, z.  If unity is 1, all coordinates will be translated to within 0 to 1.  If unity is -1, they will be translated to (-0.5, 0.5].  The default is the origin.

 

fractionalCoordinates(

   const mac::vector3D<double>& v,

   const int unity = 1);

fractionalCoordinates can also be constructed from a mac::vector3D<double>.

 

fractionalCoordinates(

   const fractionalCoordinates&,

   const int unity = 0);

fractionalCoordinates(

   const gridCoordinates&,

   const int unity = 1);

fractionalCoordinates(

   const CartesianCoordinates&,

   const cell&,

   const int unity = 1);

The three copy constructors copy fractionalCoordinates, gridCoordinates, and CartesianCoordinates to fractionalCoordinates, respectively.  The last one requires a cell.

 

All constructors above may accept an optional int unity.  They all have 1 as default, except the copy constructor from fractionalCoordinates to fractionalCoordinates.

 

operator mac::vector3D<float>() const;

The conversion operator converts to mac::vector3D<float> when needed.

 

~fractionalCoordinates();

There is nothing to do in destructor.

 

Equivalent positions

fractionalCoordinates equivalent(

   const spaceGroup&, const int i) const;

Returns the i-th equivalent position, given a space group, where i is a 1-offset subscript.  i must be from 1 to the number of symmetry operators of the space group.

 

vector<fractionalCoordinates> equivalent(

   const spaceGroup&) const;

vector<fractionalCoordinates> equivalentAU(

   const spaceGroup&);

Return the equivalent positions in a unit cell as a vector.  These functions obviously need to know the symmetry information about the space group.  See the next chapter for space group.  The first function equivalent() is a const form that does not change the fractionalCoordinates itself, but the second equivalentAU() is a non-const form.  The fractionalCoordinates will be replaced by its equivalent position in the default asymmetric unit.  The default asymmetric unit is defined by the volume in which a fractionalCoordinates is less than any of its equivalent ones.  This definition may not be exactly same as the one for gridCoordinates.  See 3.1.1.19 for vector comparison.

 

Fractional distance

double fractionalDistance(

   const fractionalCoordinates& =

         fractionalCoordinates()) const;

Returns the fractional distance to the given fractionalCoordinates.  If no one is given, the fractional distance to the origin is returned.  Please note fractional distance is not the distance in normal sense.  Greater value may not reflect greater distance.  Avoid this function if possible.  This function may be protected in future releases.

 

Output

friend ostream &operator<<(ostream&,

   const fractionalCoordinates&);

Operator << is defined to return an ostream.

 

#include <assert.h>

#include <CCL/coordinates.h>

 

main()

{

   typedef ccl::gridCoordinates       GC;

   typedef ccl::fractionalCoordinates FC;

   typedef ccl::CartesianCoordinates  CC;

   typedef mac::vector3D<unsigned short> V;

 

   const ccl::spaceGroup c2("c2");

   const ccl::cell cl(62, 80, 95, 90, 100, 90);

 

   FC fc;        // Construct a default fractional coordinates

   assert (!fc); // The default is the origin.

   fc = FC(0.5, 0.3, 0.8);

   assert (fc &= FC(0.5, -0.7, 1.8));

   assert (fc != FC(0.5, -0.7, 1.8, 0));

   assert (fabs(FC(0.09, 0.12, 0.20).fractionalDistance() - 0.25)

           < FLT_EPSILON);

 

   assert (GC(fc, V(30, 40, 50)) == GC(15, 12, 40, 30, 40, 50)); // FC -> GC

   assert (FC(GC(15, 12, 40, 30, 40, 50)) &= fc);                // GC -> FC

   assert (FC(CC(fc, cl), cl) &= fc);                            // FC <-> CC

 

   FC cp = fc;

   const vector<FC> ef = fc.equivalent(c2);

   assert (fc == cp); // fc is not changed.

   for (int i = 0; i < ef.size(); i++) cout << ef[i] << endl;

   fc.equivalentAU(c2);

   assert (fc != cp); // fc is changed into AU.

 

   return 0;

}

 

Listing 12.3.2.1.1 Testing ccl::fractionalCoordinates.

 

12.3.2.2 cpl.ccl.symm.fractionalCoordinates

Constructors

fractionalCoordinates(x = 0, y = 0, z = 0,

                      unity = 1)

Constructs a fractionalCoordinates from three doubles.  If unity is 1, all coordinates will be translated to between 0 and 1.  If it is , all coordinates will be translated to (-0.5, 0.5].  The default is the origin.

 

fractionalCoordinates(fc, unity = 0)

Copies fc to another fractionalCoordinates.

 

fractionalCoordinates(v, unity = 1)

Converts a vector3DDBL to fractionalCoordinates.

 

fractionalCoordinates(gc, unity = 1)

Converts a gridCoordinates gc to a fractionalCoordinates.

 

fractionalCoordinates(cc, cl, unity = 1)

Given a cell cl, converts a CartesianCoordinates cc to a fractionalCoordinates.

 

If unity is 1, all coordinates will translated to between 0 and 1.  If unity is -1, they will be translated to .  If it is 0, no translation will be done.  This applies to all constructors above.  Notice all defaults are 1, except the copy constructor.

 

Equivalent positions

fc.equivalent(sg, i)

Given a space group sg, returns the i-th equivalent position, where i is a 1-offset subscript.

 

fc.euivalent(sg)

Given a space group sg, returns a tuple of all equivalent positions in a unit cell.

 

fc.equivalentAU(sg)

Given a space group sg, replaces fc by its equivalent position in the default asymmetric unit.  See 12.3.2.1 for definition of default asymmetric unit.  This function also returns all equivalent positions in a unit cell.

 

Fractional distance

fc.fractionalDistance(

   p = fractionalCoordinates())

Returns the fractional distance from fc to the given p.  If nothing is given, returns the fractional distance to the origin.  Fractional distance is not the distance in normal sense.  Avoid this function if possible.

 

Output

fc.__str__()

The special function __str__ is defined to print fractionalCoordinates.

 

#!/usr/bin/python

 

import cpl

 

GC = cpl.ccl.symm.gridCoordinates

FC = cpl.ccl.symm.fractionalCoordinates

CC = cpl.ccl.symm.CartesianCoordinates

V  = cpl.mac.vm.vector3DSUI

c2 = cpl.ccl.symm.spaceGroup("c2")

cl = cpl.ccl.symm.cell(62, 80, 95, 90, 100, 90)

 

fc = FC()           # Construct a default fractional coordinates

assert fc.is_zero() # The default is the origin.

fc = FC(0.5, 0.3, 0.8);

assert fc.ap_eq(FC(0.5, -0.7, 1.8))

assert fc != FC(0.5, -0.7, 1.8, 0)

assert abs(FC(0.09, 0.12, 0.20).fractionalDistance() - 0.25) <\

       cpl.mac.utility.TINY

 

assert GC(fc, V(30, 40, 50)) == GC(15, 12, 40, 30, 40, 50) # FC -> GC

assert FC(GC(15, 12, 40, 30, 40, 50)).ap_eq(fc)            # GC -> FC

assert FC(CC(fc, cl), cl).ap_eq(fc)                        # FC <-> CC

 

cp = FC(fc)

for i in range(c2.operatorSize()): print fc.equivalent(c2, i + 1)

assert fc == cp # fc is not changed.

fc.equivalentAU(c2);

assert fc != cp # fc is changed into AU.

 

Listing 12.3.2.2.1 Testing cpl.ccl.symm.fractionalCoordinates.

 

12.3.2.3 Container of fractional coordinates

std::vector<ccl::fractionalCoordinates> is instantiated, can be used as cpl.ccl.symm.vectorFRC.

 

12.3.3 Cartesian coordinates

 

Both grid and fractional coordinates are in crystallographic coordinate system, an oblique system.  Rotation in oblique system is difficult.  In addition, oblique coordinate system is not common in other scientific fields.  The coordinates in an orthogonal system is Cartesian coordinates.  Unlike the other two types of coordinates, Calculation of equivalent positions from Cartesian coordinates in an orthogonal system is not as convenient.  In the current release of CCL, no function is implemented to return equivalent positions.  These may be added in the future.  Nevertheless, one may first convert a CartesianCoordinates to a fractionalCoordinates, calculate its equivalent positions, and then convert them back.

 

12.3.3.1 ccl::CartesianCoordinates

Constructors and destructor

CartesianCoordinates(const double x = 0,

                     const double y = 0,

                     const double z = 0);

CartesianCoordinates(

   const mac::vector3D<double>&);

Constructs a CartesianCoordinates from three doubles or a mac::vector3D<double>.  The default is the origin.  The second constructor services also as the copy constructor.

 

CartesianCoordinates(const gridCoordinates&,

                     const cell&);

CartesianCoordinates(

   const fractionalCoordinates&, const cell&);

The other two copy constructors copy gridCoordinates and fractionalCoordinates to CartesianCoordinates, respectively.  They require the cell information since orthogonalization is needed.

 

operator mac::vector3D<float>() const;

Converts to mac::vector3D<float>.

 

~CartesianCoordinates();

Nothing to do in destructor.

 

Distance

double distance(

   const CartesianCoordinates& =

         CartesianCoordinates()) const;

Returns the distance to the given CartesianCoordinates, or to the origin if nothing given.

 

Output

friend ostream &operator<<(ostream&,

   const CartesianCoordinates&);

Operator << is defined.

 

#include <assert.h>

#include <CCL/coordinates.h>

 

main()

{

   typedef ccl::gridCoordinates       GC;

   typedef ccl::fractionalCoordinates FC;

   typedef ccl::CartesianCoordinates  CC;

   typedef mac::vector3D<unsigned short> V;

 

   const ccl::spaceGroup c2("c2");

   const ccl::cell cl(62, 80, 95, 90, 100, 90);

 

   CC cc;        // Construct a default Cartesian coordinates

   assert (!cc); // The default is the origin.

   cc = CC(9, 12, 20);

   assert (fabs(cc.distance() - 25) < FLT_EPSILON); // Distance to origin

 

   // assert (CC(GC(cc, V(30, 40, 50), cl), cl) &= cc);

   // This is not guaranteed, since GC is only approximate.

   GC gc(cc, V(30, 40, 50), cl);      // CC -> GC

   cc = CC(gc, cl);                   // GC -> CC

   assert (CC(GC(cc, V(30, 40, 50), cl), cl) &= cc); // cc is on grid now.

   assert (CC(FC(cc, cl), cl) &= cc); // FC <-> CC

 

   return 0;

}

 

Listing 12.3.3.1.1 Testing ccl::CartesianCoordinates.

 

12.3.3.2 cpl.ccl.symm.CartesianCoordinates

Constructors

CartesianCoordinates(x = 0, y = 0, z = 0)

Constructs a CartesianCoordinates from three doubles.  The default is the origin.

 

CartesianCoordinates(cc)

Copies cc to another CartesianCoordinates, where cc could also be a vector3DDBL.

 

CartesianCoordinates(gc, cl)

CartesianCoordinates(fc, cl)

Given a cell cl, convert gridCoordinates and fractionalCoordinates to CartesianCoordinates.

 

Distance

cc.distance(p = CartesianCoordinates())

Returns the distance from cc to the given position p.  If nothing is given, returns the distance to the origin.

 

Output

cc.__str__()

The special function __str__ is defined to print CartesianCoordinates.

 

#!/usr/bin/python

 

import cpl

 

GC = cpl.ccl.symm.gridCoordinates

FC = cpl.ccl.symm.fractionalCoordinates

CC = cpl.ccl.symm.CartesianCoordinates

V  = cpl.mac.vm.vector3DSUI

c2 = cpl.ccl.symm.spaceGroup("c2")

cl = cpl.ccl.symm.cell(62, 80, 95, 90, 100, 90)

 

cc = CC()           # Construct a default Cartesian coordinates

assert cc.is_zero() # The default is the origin.

cc = CC(9, 12, 20)

assert abs(cc.distance() - 25) < cpl.mac.utility.TINY # Distance to origin

 

# assert CC(GC(cc, V(30, 40, 50), cl), cl).ap_eq(cc)

# This is not guaranteed, since GC is only approximate.

gc = GC(cc, V(30, 40, 50), cl)                     # CC -> GC

cc = CC(gc, cl)                                    # GC -> CC

assert CC(GC(cc, V(30, 40, 50), cl), cl).ap_eq(cc) # cc is on grid now.

assert CC(FC(cc, cl), cl).ap_eq(cc)                # FC <-> CC

 

Listing 12.3.3.2.1 Testing cpl.ccl.symm.CartesianCoordinates.

 

12.3.3.3 Container of Cartesian coordinates

std::vector<ccl::CartesianCoordinates> is instantiated, can be used as cpl.ccl.symm.vectorCTC.