International
Tables for
Crystallography
Volume B
Reciprocal space
Edited by U. Shmueli

International Tables for Crystallography (2010). Vol. B, ch. 3.3, pp. 418-434   | 1 | 2 |

## Section 3.3.1. Graphics

R. Diamonda

### 3.3.1. Graphics

| top | pdf |

| top | pdf |

#### 3.3.1.1.1. Cartesian and crystallographic coordinates

| top | pdf |

It is usual, for purposes of molecular modelling and of computer graphics, to adopt a Cartesian coordinate system using mutually perpendicular axes in a right-handed system using the ångström unit or the nanometre as the unit of distance along such axes, and largely to ignore the existence of crystallographic coordinates expressed as fractions of unit-cell edges. Transformations between the two are thus associated, usually, with the input and output stages of any software concerned with modelling and graphics, and it will be assumed after this section that all coordinates are Cartesian using the chosen unit of distance as the unit of coordinates. For a discussion of coordinate transformations and rotations without making this assumption see Chapter 1.1 in which formulations using co- and contravariant forms are presented.

The relationship between these systems may be written in which X and x are position vectors in direct space, written as column vectors, with x expressed in crystallographic fractional coordinates (dimensionless) and X in Cartesian coordinates (dimension of length).

There are two forms of M in common use. The first of these sets the first component of X parallel to and the third parallel to c and is in which

ϕ is equal to the volume of the unit cell divided by abc, and is unchanged by cyclic permutation of α, β and γ and of , and . The Cartesian and crystallographic axes have the same chirality if the positive square root is taken.

The second form sets the first component of X parallel to a and the third component of X parallel to and is

A third form, suitable only for rhombohedral cells, is in which which preserves the equivalence of axes. Here the chiralities of the Cartesian and crystallographic axes are the same if p is chosen positive, and different otherwise, and the two sets of axes coincide in projection along the triad if q is chosen positive and are π out of phase otherwise.

#### 3.3.1.1.2. Homogeneous coordinates

| top | pdf |

Homogeneous coordinates have found wide application in computer graphics. For some equipment their use is essential, and they are of value analytically even if the available hardware does not require their use.

Homogeneous coordinates employ four quantities, X, Y, Z and W, to define the position of a point, rather than three. The fourth coordinate has a scaling function so that it is the quantity X/W (as delivered to the display hardware) which controls the left–right positioning of the point within the picture. A point with is in the picture, normally, and those with are outside it, but see Section 3.3.1.3.5.

There are many reasons why homogeneous coordinates may be adopted, among them the following:

 (i) X, Y, Z and W may be held as integers, thus enabling fast arithmetic whilst offering much of the flexibility of floating-point working. A single W value may be common to a whole array of X, Y, Z values. (ii) Perspective transformations can be implemented without the need for any division. Only high-speed matrix multiplication using integer arithmetic is necessary, provided only that the drawing hardware can provide displacements proportional to the ratio of two signals, X and W or Y and W. Rotation, translation, scaling and the application of perspective are all effected by operations of the same form, namely multiplication of a four-vector by a 4 × 4 matrix. The hardware may thus be kept relatively simple since only one type of operation needs to be provided for. (iii) Since kX, kY, kZ, kW represents the same point as X, Y, Z, W, the hardware may be arranged to maximize resolution without risk of integer overflow.

For analytical purposes it is convenient to regard homogeneous transformations in terms of partitioned matrices where M is a 3 × 3 matrix, V and X are three-element column vectors, U is a three-element row vector and N and W are scalars.

Matrices and vectors which are equivalent under the considerations of (iii) above will be related by the sign ≃ in what follows.

Hardware systems which use true floating-point representations have less need of homogeneous coordinates and for these N and W may normally be set to unity.

#### 3.3.1.1.3. Notation

| top | pdf |

In this chapter the conventions of matrix algebra will be adhered to except where it is convenient to show operations on elements of vectors, matrices and tensors, where a subscript notation will be used with a modified summation convention in which summation is over lower-case subscripts only. Thus the equation is to be read For any I, xI is AIjXj summed over j'.

Subscripts using the letter i or later in the alphabet will relate to the usual three dimensions and imply a three-term summation. Subscripts a to h are not necessarily so limited, and, in particular, the subscript a is used to imply summation over atoms of which there may be an arbitrary number.

We shall use the superscript T to denote a transpose, and also use the Kronecker delta, , which is 1 if I = J and zero otherwise, and the tensor which is 1 if I, J and K are a cyclic permutation of 1, 2, 3, −1 if an anticyclic permutation, and zero otherwise.

A useful identity is then Single modulus signs surrounding the symbol for a square matrix denote its determinant, and around a vector denote its length.

The symbol ≃ is defined in the previous section.

#### 3.3.1.1.4. Standards

| top | pdf |

The sections of this chapter concerned with graphics are primarily concerned with the mathematical aspects of graphics programming as they confront the applications programmer. The implementations outlined in Section 3.3.3 have all, so far as the author is aware, been developed ab initio by their inventors to deal with these aspects using their own and unrelated techniques and protocols. It is clear, however, that standards are now emerging, and it is to be hoped that future developments in applications software will handle the graphics aspects through one or other of these standards.

First among these standards is the Graphical Kernel System, GKS, defined in American National Standards Institute, American National Standard for Information Processing Systems – Computer Graphics – Graphical Kernel System (GKS) Functional Description (1985) and described and illustrated by Hopgood et al. (1986) and Enderle et al. (1984). GKS became a full International Standards Organization (ISO) standard in 1985, and its purpose is to standardize the interface between application software and the graphics system, thus enhancing portability of software. Specifications for Fortran, Pascal and Ada formulations are at an advanced stage of development. Its value to crystallographers is limited by the fact that it is only two-dimensional. A three-dimensional extension known as GKS-3D, defined in International Standards Organisation, International Standard Information Processing Systems – Computer Graphics – Graphical Kernel System for Three Dimensions (GKS-3D), Functional Description (1988) became an ISO standard in 1988. Perhaps of greatest interest to crystallographers, however, is the Programmers' Hierarchical Interactive Graphics System (PHIGS) (Brown, 1985; Abi-Ezzi & Bunshaft, 1986) since this allows hierarchical segmentation of picture content to exist in both the applications software and the graphics device in a related manner, which GKS does not. Some graphics devices now available support this type of working and its exploitation indicates the choice of PHIGS. Furthermore, Fortran implementations of GKS and GKS-3D require points to be stored in arrays dimensioned as X(N), Y(N), Z(N) which may be equivalenced (in the Fortran sense) to XYZ(N, 3) but not to XYZ(3, N), which may not be convenient. PHIGS also became an International Standard in 1988: American National Standards Institute, American National Standard for Information Processing Systems – Computer Graphics – Programmer's Hierarchical Graphics System (PHIGS) Functional Description, Archive File Format, Clear-Text Encoding of Archive File (1988). PHIGS has also been extended to support the capability of raster-graphics machines to represent reflections, shadows, see-through effects etc. in a version known as PHIGS+ (van Dam, 1988).

Increasingly, manufacturers of graphics equipment are orienting their products towards one or other of these standards. While these standards are not the subject of this chapter it is recommended that they be studied before investing in equipment.

In addition to these standards, related standards are evolving under the auspices of the ISO for defining images in a file-storage, or metafile, form, and for the interface between the device-independent and device-dependent parts of a graphics package. Arnold & Bono (1988) describe the ANSI and ISO Computer Graphics Metafile standard which provides for the definition of (two-dimensional) images. The definition of three-dimensional scenes requires the use of (PHIGS) archive files.

#### 3.3.1.2. Orthogonal (or rotation) matrices

| top | pdf |

It is a basic requirement for any graphics or molecular-modelling system to be able to control and manipulate the orientation of the structures involved and this is achieved using orthogonal matrices which are the subject of these sections.

#### 3.3.1.2.1. General form

| top | pdf |

If a vector v is expressed in terms of its components resolved onto an axial set of vectors X, Y, Z which are of unit length and mutually perpendicular and right handed in the sense that , and if these components are , and if a second set of axes X′, Y′, Z′ is similarly established, with the same origin and chirality, and if v has components on these axes then in which is the cosine of the angle between the ith primed axis and the jth unprimed axis. Evidently the elements comprise a matrix R, such that any row represents one of the primed axial vectors, such as X′, expressed as components on the unprimed axes, and each column represents one of the unprimed axial vectors expressed as components on the primed axes. It follows that since elements of the product are scalar products among perpendicular unit vectors.

A real matrix whose transpose equals its inverse is said to be orthogonal.

Since X, Y and Z can simultaneously be superimposed on X′, Y′ and Z′ without deformation or change of scale the relationship is one of rotation, and orthogonal matrices are often referred to as rotation matrices. The operation of replacing the vector v by Rv corresponds to rotating the axes from the unprimed to the primed set with v itself unchanged. Equally, the same operation corresponds to retaining fixed axes and rotating the vector in the opposite sense. The second interpretation is the one more frequently helpful since conceptually it corresponds more closely to rotational operations on objects, and it is primarily in this sense that the following is written.

If three vectors u, v and w form the edges of a parallelepiped, then its volume V is and if these vectors are transformed by the matrix R as above, then the transformed volume V′ is But the determinant of R is given by so that and the determinant of R must therefore be +1 for a transformation which is a pure rotation. Nevertheless orthogonal matrices with determinant −1 exist though these do not describe a pure rotation. They may always be described as the product of a pure rotation and inversion through the origin and are referred to here as improper rotations. In what follows all references to orthogonal matrices refer to those with positive determinant only, unless stated otherwise.

An important general form of an orthogonal matrix in three dimensions was derived as equation (1.1.4.32 ) and is or in which l, m and n are the direction cosines of the axis of rotation (which are the same when referred to either set of axes under either interpretation) and θ is the angle of rotation. In this form, and with R operating on column vectors on the right, the sign of θ is such that, when viewed along the rotation axis from the origin towards the point lmn, the object is rotated clockwise for positive θ with a fixed right-handed axial system. If, under the same viewing conditions, the axes are to be rotated clockwise through θ with the object fixed then the components of vectors in the object, on the new axes, are given by R with the same lmn and with θ negated. This is the transpose of R, and if R is constructed from a product, as below, then each factor matrix in the product must be transposed and their order reversed to achieve this. Note that if, for a given rotation, the viewing direction from the origin is reversed, l, m, n and θ are all reversed and the matrix is unchanged.

Any rotation about a reference axis such that two of the direction cosines are zero is termed a primitive rotation, and it is frequently a requirement to generate or to interpret a general rotation as a product of primitive rotations.

A second important general form is based on Eulerian angles and is the product of three such primitives. It is which is commonly employed in four-circle diffractometers for which , and . In terms of the fixed-axes–moving-object conceptualization this corresponds to a rotation about Z followed by about Y followed by about Z. In the familiar diffractometer example, when the ϕ and ω axes are both vertical and equivalent. If ϕ is altered first, then the χ axis is still in the direction of a fixed Y axis, but if ω is altered first it is not. Since all angles are to be rotations about fixed axes to describe a rotating object it follows that it is ϕ rather than ω which corresponds to . In general, when rotating parts are mounted on rotating parts the rotation closest to the moved object must be applied first, forming the right-most factor in any multiple transformation, with the rotation closest to the fixed part as the left-most factor, assuming data supplied as column vectors on the right.

Given an orthogonal matrix, in either numerical or analytical form, it may be required to discover θ and the axis of rotation, or to factorize it as a product of primitives. From the first form we see that the vector consisting of the antisymmetric part of R, has elements times the direction cosines l, m, n, which establishes the direction immediately, and normalization using determines . Furthermore, the trace is so that the quadrant of θ is also fixed. This method fails, however, if the matrix is symmetrical, which occurs if . In this case only the direction of the axis is required, which is given by for nonzero elements, or etc., with the signs chosen to satisfy etc.

The Eulerian form may be factorized by noting that . There is then freedom to choose the sign of , but the choice then fixes the quadrants of and through the elements in the last row and column, and the primitives may then be constructed. These expressions for and fail if , in which case the rotation reduces to a primitive rotation about Z with angle , or .

Eulerian angles are unlikely to be the best choice of primitive angles unless they are directly related to the parameters of a system, as with the diffractometer. It is often more important that the changes to primitive angles should be quasi-linearly related to θ for any small rotations, which is not the case with Eulerian angles when the required rotation axis is close to the X axis. In such a case linearized techniques for solving for the primitive angles will fail. Furthermore, if the required rotation is about Z only is determinate.

Quasi-linear relationships between θ and the primitive rotations arise if the primitives are one each about X, Y and Z. Any order of the three factors may be chosen, but the choice must then be adhered to since these factors do not commute. For sufficiently small rotations the primitive rotations are then , and , whilst for larger θ linearized iterative techniques for finding the primitive rotations are likely to be convergent and well conditioned.

The three-dimensional space of the angles , and in either case is nonlinearly related to θ. In the Eulerian case the worst nonlinearities occur at the origin of ϕ-space. Equally severe nonlinearities occur in the quasi-linear case also but are 90° away from the origin and less likely to be troublesome.

Neither of the foregoing general forms of orthogonal matrix has ideally convenient properties. The first is inconvenient because it uses four nonequivalent variables l, m, n and θ, with a linking equation involving l, m and n, so that they cannot be treated as independent variables for analytical purposes. The second form (the product of primitives) is not ideal because the three angles, though independent, are not equivalent, the non­equivalence arising from the noncommutation of the primitive factors. In the remainder of this section we give two further forms of orthogonal matrix which each use three variables which are independent and strictly equivalent, and a third form using four whose squares sum to unity.

The first of these is based on the diagonal and uses the three independent variables p, q, r, from which we construct the auxiliary variables then is orthogonal with positive determinant for any of the sixteen sign combinations. The signs of P, Q, R and S are, respectively, the signs of the direction cosines of the rotation axis and of . Using also , which may be deemed positive without loss of generality,

Although p, q and r are independent, the point [pqr] is bound, by the requirement that P, Q, R and S be real, to lie within a tetrahedron whose vertices are the points [111], , and , corresponding to the identity and to 180° rotations about each of the axes. The facts that the identity occurs at a vertex of the feasible region and that , rather than , is linear on p, q and r in this vicinity make this form suitable only for substantial rotations.

The second form consists in defining a rotation vector r with components u, v, w such that , , with and . Then the matrix is orthogonal and the variables u, v, w are independent, equivalent and unbounded, and, unlike the previous form, small rotations are quasi-linear on these variables. As examples, r = [100] gives 90° about X, r = [111] gives 120° about [111].

R then transforms a vector d according to

Multiplying two such matrices together allows us to establish the manner in which the rotation vectors r1 and r2 combine. for a rotation r1 followed by r2, so that rotations expressed in terms of rotation angles and axes may be compounded into a single such rotation without the need to form and decompose a product matrix.

Note that if r1 and r2 are parallel this reduces to the formula for the tangent of the sum of two angles, and that if the combined rotation is always 180°. Note, too, that reversing the order of application of the rotations reverses only the vector product.

If three rotations r1, r2 and r3 are applied successively, r1 first, then their combined rotation is

Note the irregular pattern of signs in the numerator.

Similar ideas, using a vector of magnitude , are developed in Aharonov et al. (1977).

The third form of orthogonal matrix uses four variables, λ, μ, ν and σ, which comprise a four-dimensional vector , such that , , with and . In terms of these variables Two further matrices S and T may be defined (Diamond, 1988), which are themselves orthogonal (though S has determinant −1) and which have the property that so that, for example, if homogeneous coordinates are being employed (Section 3.3.1.1.2) is a rotation of (x, y, z, w) through the angle θ about the axis (l, m, n). With suitably pipelined hardware this forms an efficient means of applying rotations since the overhead' of establishing S is so trivial.

T has the property that the rotation vector arising from a concatenation of n rotations is in which is the vector (0, 0, 0, 1) which defines a null rotation. This equation may be used as a basis for factorizing a given rotation into a concatenation of rotations about designated axes (Diamond, 1990a).

Finally, an exact rotation of the vector d may be obtained without using matrices at all by writing in which and is the initial position which is to be rotated. Here is a vector with direction cosines l, m and n, and magnitude equal to the required rotation angle in radians (Diamond, 1966). This method is particularly efficient when or when the number of vectors to be transformed is small since the overhead of establishing R is eliminated and the process is simple to program. It is the three-dimensional analogue of the power series for sin θ and cos θ and has the same convergence properties.

#### 3.3.1.2.2. Measurement of rotations and strains from coordinates

| top | pdf |

Given the coordinates of a molecular fragment it is often a requirement to relate the fragment to its image in some standard orientation by a transformation which may be required to be a pure rotation, or may be required to be a combination of rotation and strain. Of the methods reviewed in this section all except (iv) are concerned with pure rotation, ignoring any strain that may be present, and give the best rigid-body superposition. In all these methods, unless inhomogeneous strain is being considered, the best possible superposition is obtained if the centroids of the two images are first brought into coincidence by translation and treated as the origin.

Methods (i) to (v) seek transformations which perform the superposition and impose on these, in various ways, the requirements of orthogonality for the rotational part. All these methods therefore need some defence against indeterminacy that arises in the general transformation if one or both of the fragments is planar, and, if improper rotations are to be excluded, need a defence against these also if the fragment and its image are of opposite chirality. Methods (vi) and (vii) pay no attention to the general transformation (which defines the general superposition) and work with variables which are intrinsically rotational in character, and always produce an orthogonal transformation with positive determinant, with no degeneracy arising from planar fragments which need no special attention. Even collinear atoms cause no problem, the superposition being performed correctly but with an arbitrary rotation about the length of the line being present in the result. These methods are therefore to be preferred over the earlier ones unless the purpose of the operation is to detect differences of chirality, although this, too, can be detected with a simple test.

In this review we adopt the same notation for all the methods which, unavoidably, means that symbols are used in ways which differ from the original publications. We use the symbol x for the vector set which is to be rotated and X for the vector set whose orientation is not to be altered, and write the residuals as and, by choice of origin, for weights W. The quadratic residual to be minimized is and we define the matrix and use l for the direction cosines of the rotation axis.

• (i) McLachlan's first method (McLachlan, 1972, 1982) is iterative and conceptually the simplest. It sets in which A and R are both orthogonal with R being a current estimate of D and A being an adjustment which, at the beginning of each cycle, has a zero angle associated with it. One iterative cycle estimates a nontrivial A, after which the product AR replaces R. and therefore For this to vanish for all possible rotation axes l the vector must vanish, i.e. at the end of the iteration R must be such that the matrix is symmetrical. The vector g represents the couple exerted on the rotating body by forces acting at the atoms. Choosing gives the greatest and vanishes when in which N is constructed from the current R matrix. A is then constructed from l and this θ and AR replaces R. The process is iterative because a couple about some new axis can appear when rotation about g eliminates the couple about g.

Note that for each rotation axis l there are two values of θ, differing by π, which reduce to zero, corresponding to maximum and minimum values of E. The minimum is that which makes positive. Adding π to θ alters R and N and negates this quantity.

Note, too, that the process is essentially characterized as that which makes the product RM symmetrical with R orthogonal. We return to this point in (iii).

• (ii) Kabsch's method (Kabsch, 1976, 1978) minimizes E with respect to the nine elements of D, subject to the six constraints by using an auxiliary function in which L is symmetric containing six Lagrange multipliers. The Lagrangian functionthen has minima with respect to the elements of D at locations which are dependent, inter alia, on the elements of L. By suitably choosing L a minimum of G may be brought into coincidence with the constrained minimum of E. A minimum of G occurs where and the matrix is positive definite, block diagonal, and has which is symmetrical. Thus L must be chosen so as to make the symmetric matrix such that with D orthogonal, or with R replacing D since we are now confined to the orthogonal case, and N is symmetric and positive definite.

• (iii) Comparison of the Kabsch and McLachlan methods. Using the initials of these authors as subscripts, we have seen that the Kabsch solution involves solving for an orthogonal matrix given that is symmetrical and positive definite. Thus and

By comparison, the McLachlan treatment leads to an orthogonal R matrix satisfying in which is also symmetric and positive definite, which similarly leads to

These seemingly different expressions for and are, in fact, equal, as the following shows therefore Multiplying on both sides by gives and since both N matrices are positive definite and conversely therefore

• (iv) Diamond's first method. This method (Diamond, 1976a) differs from the previous ones in that the transformation D is allowed to be a general transformation which is then factorized into the product of an orthogonal matrix R and a symmetrical matrix T. The transformation of x to fit X is thus interpreted as the combination of homogeneous strain and pure rotation in which x is subjected to strain and the result is rotated. Furthermore, the solution for D is (in the notation of Kabsch), so that which may be compared with the results of the previous paragraph.

Although this R matrix by itself (i.e. applied without T) does not produce the best rotational superposition (i.e. smallest E), it is the one which exactly superposes the only three vectors in x whose mutual dispositions are conserved, on their equivalents in X, so that the rotation so found is arguably the best defined one.

Alternatives based on , , are all easily developed, and these ideas are developed by Diamond (1976a) to include nonhomogeneous strains also.

• (v) McLachlan's second method. This method (McLachlan, 1979) is based on the properties of the matrix and is immune to singularity of M. If p and q are three-dimensional vectors such that is an eigenvector of this matrix then

If q is negated the second equality is maintained provided λ is also negated. Therefore an orthogonal matrix (consisting of partitions) exists for which in which is diagonal and contains non-negative eigenvalues. The reverse transformation shows that and multiplying the eigenvectors together gives Therefore but is orthogonal and is symmetrical, therefore [by paragraphs (i) and (iii) above] is the required rotation. Similarly, forming corresponds to the Kabsch formulation [paragraphs (ii) and (iii)] since is symmetrical and the same rotation, , appears.

Note that the determinant of the orthogonal matrix so found is twice the product of the determinants of H and of K, and since the positive eigenvalues are collected into it follows that the sign of the determinant of M is the same as the sign of the determinant of the resulting orthogonal matrix. If this is negative it means that the best superposition is obtained if one vector set is inverted and that x and X are of opposite chirality.

Expanding the expression for E, the weighted sum of squares of errors, for an orthogonal transformation shows that this is least when the trace of the product RM is greatest. In this treatment Hence, if the eigenvalues in and − are arranged in decreasing order of modulus, and if the determinant of M is negative, then exchanging the third and sixth columns of produces a product with positive determinant (i.e. a proper rotation) at minimum cost in residual. Similarly, if M is singular and one or more eigenvalues in vanishes it is necessary only to complete an orthonormal set of eigenvectors such that the determinants of H and K have the same sign.

• (vi) MacKay's method. MacKay (1984) was the first to consider the rotational superposition problem in terms of the vector r of Section 3.3.1.2.1. Using quaternion algebra he showed that if a vector x is rotated to then where and the direction of r is the axis of rotation, as may also be shown from elementary considerations. MacKay then solves this for the vector r by least squares given the vector pairs X and x. The individual errors are and Setting gives which reduces to in which Thus a direct solution for r is obtained, the elements of which are u, v and w, and may be used to construct the orthogonal matrix as in Section 3.3.1.2.1. may be obtained directly from .

If the requisite rotation is 180°, is singular and cannot be inverted. In this case any row or column of the adjoint of is a vector in the direction of the axis. Normalizing this vector to unity, giving l, gives the requisite orthogonal matrix as

Note that MacKay's residual E is quadratic in r. E therefore has one minimum and no maximum, and the minimum is reached on the first cycle of least squares. By contrast, the objective function E that is minimized in methods (i), (ii), (v) and (vii) has one minimum, one maximum and two saddle points in the space of the vector r, as shown in (vii).

It may be shown (Diamond, 1989) that if MacKay's solution vector r is denoted by and that given by the other methods [except (iv)] by then in which A and B are real symmetric, positive semi-definite. A is positive definite unless all the individual vector sums are parallel, as can happen when the best rotation is 180°. Thus the MacKay method only gives the same result as the other methods if:

 (a) the initial orientation is optimal, for then , or (b) perfect fitting is possible, for then , or (c) all the residual vectors (after fitting by ) are parallel to , for then B is singular such that . In general, . may be found by iterating , but x must be replaced by Rx on each iteration.

• (vii) Diamond's second method. This is closely related to MacKay's method, but uses a four-dimensional vector with components λ, μ, ν and σ in which λ, μ and ν are the direction cosines of the rotation axis multiplied by and σ is . In terms of such a vector Diamond (1988) showed that in which E is the weighted sum of squares of coordinate differences, as before, is its value before any rotation is applied and P is the matrix The rotation matrix R corresponding to the vector is then the last of the forms for R given in Section 3.3.1.2.1.

The minimum E is therefore minus twice the largest eigenvalue of P since , and a stationary value of E occurs when is any of the four eigenvectors of P. E thus has a maximum, a minimum and two saddle points, in general, and its value may be determined before any coordinates are transformed. Diamond also showed that the orientations giving these stationary values are related by the operations of 222 symmetry. Equivalent results have also been obtained by Kearsley (1989).

As an alternative to solving a eigenproblem, Diamond also showed that the vector r, as in MacKay's solution, may be obtained by iterating which has the property that if X and x are exactly superposable then is the exact solution and no iteration is necessary. If X and x are similar but not exactly superposable then a small number of iterations may be required to reach a stable r vector, though the matrix is not required. As in MacKay's solution, is singular at the end of the iteration if the required rotation is 180°, but the MacKay and Diamond methods both have the advantage that improper rotations are never generated by these means, and methods based on P and rather than Q and r are trouble-free for 180° rotations. The iterative loop in this method does not require Rx to be redetermined on each cycle.

Finally, it may be shown that if are the eigenvalues of P arranged in descending order and is negative, then a closer superposition may be obtained by reversing the chirality of one of the vector sets, and the R matrix constructed from optimally superimposes Rx onto − X, the enantiomer of X (Diamond, 1990b).

These ideas lend themselves to the problem of the multiple simultaneous superposition of many structures, which is a nonlinear problem with the possibility of multiple solutions [Diamond (1992); see also Kearsley (1990) and Shapiro et al. (1992)], and they can provide for cluster analysis based on structural similarity to find subsets of similar structures within an ensemble of structures, such as may arise from NMR calculations (Diamond, 1995).

#### 3.3.1.2.3. Orthogonalization of impure rotations

| top | pdf |

There are several ways of deriving a strictly orthogonal matrix from a given approximately orthogonal matrix, among them the following.

 (i) The Gram–Schmidt process. This is probably the simplest and the easiest to compute. If the given matrix consists of three column vectors v1, v2 and v3 (later referred to as primers) which are to be replaced by three column vectors u1, u2 and u3 then the process is As successive vectors are established, each vector v has subtracted from it its components in the directions of established vectors, and the remainder is normalized. The method will fail at the normalization step if the vectors v are not linearly independent. Otherwise, the process may be extended to any number of dimensions. The weakness of the method is that, though u1 differs from v1 only in scale, uN may differ grossly from vN as the various columns are not treated equivalently. (ii) A preferable method which treats all vectors equivalently is to iteratively replace the matrix M by . Defining the residual matrix E as then on each iteration E is replaced by and convergence necessarily ensues. (iii) A third method resolves M into its symmetric and antisymmetric parts and constructs an orthogonal matrix for which only S is altered. A determines l, m, n and θ as shown in Section 3.3.1.2.1, and from these a new S may be constructed. (iv) A fourth method is to treat the general matrix M as a combination of pure strain and pure rotation. Setting with R orthogonal and T symmetrical gives The rotation so found is the one which exactly superposes those three mutually perpendicular directions which remain mutually perpendicular under the transformation M. is then the strain tensor of an unrotated body. Writing , , may also be useful, in which is the strain tensor of a rotated body. See also Section 3.3.1.2.2 (iv).

#### 3.3.1.2.4. Eigenvalues and eigenvectors of orthogonal matrices

| top | pdf |

If R is the orthogonal matrix given in Section 3.3.1.2.1 in terms of the direction cosines l, m and n of the axis of rotation, then it is clear that (l, m, n) is an eigenvector of R with eigenvalue unity because

Consideration of the determinant shows that the sum of the three eigenvalues is and that their product is unity. Hence the three eigenvalues are 1, and . Since R is real, its product with any real vector is also real, yet its product with an eigenvector must, in general, be complex. Thus the eigenvectors must themselves be complex.

The remaining two eigenvectors u may be found using the results of Section 3.3.1.2.1 (q.v.) according to which is solved by any vector of the form for any real vector v, where l is the normalized axis vector, , , . Eigenvectors for the two eigenvalues may have unrelated v vectors though the sign choices are coupled. If the vector v is rotated about l through an angle ϕ the corresponding vector u is multiplied by and remains an eigenvector. Using superscript signs to denote the sign of θ in the eigenvalue with which each vector is associated, the matrix has the properties that and which places restrictions on v if this is to be the identity. Note that the 23 element vanishes even in the absence of any relationship between and .

A convenient form for U, symmetrical in the elements of l, is obtained by setting and is in which the normalizing denominator is given by

#### 3.3.1.3. Projection transformations and spaces

| top | pdf |

In the following section we address the question of the relationship between the coordinates of a molecular model and the corresponding coordinates on the screen of the graphics device. A good introduction to this topic is given by Newman & Sproull (1973), and Foley et al. (1990) give a comprehensive account of the field, including recent developments, especially those arising from the development of raster-graphics technologies.

#### 3.3.1.3.1. Definitions

| top | pdf |

Typically, the coordinates, X, of points in an object to be drawn are held in homogeneous Cartesian form as described in Section 3.3.1.1.2. Such coordinates are said to be in data space or world coordinates and this coordinate system is generally a constant aspect of the problem.

In order to view these data in convenient ways such coordinates may be subjected to a viewing transformation T, affecting orientation, scale etc., the resulting coordinates TX being then in display space. Here, and throughout what follows, we treat position vectors as columns with transformation matrices as factors on the left, though some writers do the reverse.

In general, only some portion of display space which lies inside a certain frustum of a pyramid is required to fall within the picture. The pyramid may be thought of as having the observer's eye at its vertex, with a rectangular base corresponding to the picture area. This volume is called a window. A transformation, U, which takes display-space coordinates as input and generates vectors (X, Y, Z, W) for which and for points on the left, right, top and bottom boundaries of the window and for which takes particular values on the front and back planes of the window, is said to be a windowing transformation. In machines for which controls intensity depth cueing, the range of corresponding to the window is likely to be 0 to 1 rather than −1 to 1. Coordinates obtained by multiplying display-space coordinates by U are termed picture-space coordinates. Mathematically, U is a matrix like any other, but functionally it is special. Declaring a transformation to be a windowing transformation implies that only resulting points having and positive are to be plotted. Machines with clipping hardware to truncate lines which run out of the picture perform clipping on the output from the windowing transformation.

Finally, the picture has to be drawn in some rectangular portion of the screen which is allocated for the purpose. Such an area is termed a viewport and is defined in terms of screen coordinates which are defined absolutely for the hardware in question as ±n for full-screen deflection, where n is declared by the manufacturer. Screen coordinates are obtained from picture coordinates with a viewport transformation, V.1

To summarize, screen coordinates, S, are given by

#### 3.3.1.3.2. Translation

| top | pdf |

The transformation evidently corresponds to the addition of the vector to the components of X or of to the components of . (I is the identity.) Displacements may thus be affected by expressing the required displacement vector in homogeneous coordinates with any suitable choice of N (commonly, ), with V scaled to correspond to this choice, and loading the transformation matrix as indicated above.

#### 3.3.1.3.3. Rotation

| top | pdf |

Rotation about the origin is achieved by in which R is an orthogonal matrix. R necessarily has elements not exceeding one in modulus. For machines using integer arithmetic, therefore, N would be chosen large enough (usually half the largest possible integer) for the product NR to be well represented in the available word length. Characteristically, N affects resolution but not scale.

#### 3.3.1.3.4. Scale

| top | pdf |

The transformation scales the vector (X, W) by the factor S. For integer working and , N should be set to the largest representable integer. For the product SN should be the largest representable integer, N being reduced accordingly.

#### 3.3.1.3.5. Windowing and perspective

| top | pdf |

It is necessary at this point to relate the discussion to the axial system inherent in the graphics device employed. One common system adopts X horizontal and to the right when viewing the screen, Y vertically upwards in the plane of the screen, and Z normal to X and Y with +Z into the screen. This is, unfortunately, a left-handed system in that is negative. Since it is usual in crystallographic work to use right-handed axial systems it is necessary to incorporate a matrix of the form either as the left-most factor in the matrix T or as the right-most factor in the windowing transformation U (see Section 3.3.1.3.1). The latter choice is to be preferred and is adopted here. The former choice leads to complications if transformations in display space will be required. Display-space coordinates are necessarily referred to this axial system.

Let L, R, T, B, N and F be the left, right, top, bottom, near and far boundaries of the windowed volume , S be the Z coordinate of the screen, and C, D and E be the coordinates of the observer's eye position, all ten of these parameters being referred to the origin of display space as origin, which may be anywhere in relation to the hardware. L, R, T and B are to be evaluated in the screen plane. All ten parameters may be referred to their own fourth coordinate, V, meaning that the point (X, Y, Z, W) in display space will be on the left boundary of the picture if when . V may be freely chosen so that all eleven quantities and all elements of U suit the word length of the machine. These relationships are illustrated in Fig. 3.3.1.1.

 Figure 3.3.1.1 | top | pdf |The relationship between display-space coordinates (X, Y, Z, W) and picture-space coordinates (x, y, z, w) derived from them by the window transformation, U. (a) Display space (in X, Z projection) showing a square object P, Q, R, S for display viewed from the position (C, D, E, V). The bold trapezium is the window (volume) and the bold line is the viewport portion of the screen. The points P, Q, R and S must be plotted at p, q, r and s to give the correct impression of the object. (b) Picture space (in x, z projection). The window is mapped to a rectangle and all sight lines are parallel to the z axis, but the object P, Q, R, S is no longer square. The distribution of p, q, r and s is identical in the two cases. Note that z/w values are not linear on Z/W, and that the origin of picture space arises at the midpoint of the near clipping plane, regardless of the location of the origin of display space. The figure is accurately to scale for coincident viewport positions. The words Left clipping plane', if part of the scene in display space, would currently be obscured, but would come into view if the eye moved to the right, increasing C, as the left clipping plane would pivot about the point L/V in the screen plane.

Since is a display-space coordinate on the same scale as the window parameters. This must be plotted on the screen at an X coordinate (on the scale of the window parameters) which is the weighted mean of and C, the weights being and , respectively. This provides perspective because the weighted mean is at the point where the straight line from to the eye intersects the screen. This then has to be mapped into the L-to-R interval, so that picture-space coordinates are given by which provides for and to be unity on the picture boundaries, which is usually a requirement of the clipping hardware, and for , zero being for the near-plane boundary. Even though z/w is not linear on Z/W, straight lines and planes in display space transform to straight lines and planes in picture space, the nonlinearity affecting only distances. Thus vector-drawing machines are not disadvantaged by the introduction of perspective.

Note that the dimensionality of X/W must equal that of S/V and that this may be regarded as length or as a pure number, but that in either case x/w is dimensionless, consistent with the stipulation that the picture boundaries be defined by the pure number ±1.

The above matrix is U and is suited to left-handed hardware systems. Note that only the last column of U (the translational part) is sensitive to the location of the origin of display space and that if the eye is on the normal to the picture centre then , and simplifications result. If C, D and E can be continuously monitored then dynamic parallax as well as perspective may be obtained (Diamond et al., 1982).

If data space is referred to right-handed axes, the viewing transformation T involves only proper rotations and the hardware uses a left-handed axial system then elements in the third column of U should be negated, as explained in the opening paragraph.

To provide for orthographic projection, multiply every element of U by and then let , choosing some positive K to suit the word length of the machine [see Section 3.3.1.1.2 (iii)]. The result is which is the orthographic window.

It may be convenient in some applications to separate the functions of windowing and the application of perspective, and to write where U and U′ are as above and P is a perspective transformation given by which involves F and N but not R, L, T or B. In this form the action of P may be thought of as compressing distant parts of display space prior to an orthographic projection by U′ into picture space.

Other factorizations of U are possible, for example with which renders P′ independent of all six boundary planes, but U′′ is no longer independent of E. It is not possible to factorize U so that the left factor is a function only of the boundary planes and the right factor a function only of eye and screen positions.

Note that as , , P and .

#### 3.3.1.3.6. Stereoviews

| top | pdf |

Assuming that left- and right-eye views are to be presented through the same viewport (next section) or that their viewports are to be superimposed by an external optical system, e.g. Ortony mirrors, then stereopairs are obtained by using appropriate eye coordinates in the U matrix of the previous section. However, U may be factorized according to in which U′′′ is the matrix U obtained by setting to correspond to the point midway between the viewer's eyes and in which (c, 0, 0, V) is the position of the right eye relative to the mean eye position, and the left-eye view is obtained by negating c.

Stereo is often approximated by introducing a rotation about the Y axis of to the views or to one of them. The first corresponds to with . The main difference is in the resulting Z value, which only affects depth cueing and z clipping. The X translation which arises if is also suppressed, but this is not likely to be noticeable. σ is often treated as a constant, such as .

The distinction in principle between the true S and the rotational approximation is that with the true S the eye moves relative to the screen and the displayed object, whereas with the approximation the eye and the screen are moved relative to the displayed object, in going from one view to the other.

Strobing of left and right images may conveniently be accomplished with an electro-optic liquid-crystal shutter as described by Harris et al. (1985). The shutter is switched by the display itself, thus solving the synchronization problem in a manner free of inertia.

A further discussion of stereopairs may be found in Johnson (1970) and in Thomas (1993), the second of which generalizes the treatment to allow for the possible presence of an optical system.

#### 3.3.1.3.7. Viewports

| top | pdf |

The window transformation of the previous two sections has been constructed to yield picture coordinates (X, Y, Z, W) (formerly called x, y, z, w) such that a point having or is on the boundary of the picture, and the clipping hardware operates on this basis. However, the edges of the picture need not be at the edges of the screen and a viewport transformation, V, is therefore needed to position the picture in the requisite part of the screen. where r, l, t and b are now the right, left, top and bottom boundaries of the picture area, or viewport, expressed in screen coordinates, and n is the full-screen deflection value. Thus a point with in picture space plots on the screen with an X coordinate which is a fraction of full-screen deflection to the right. is unchanged by V and is used only to control intensity in a technique known as depth cueing.

It is necessary, of course, to arrange for the aspect ratio of the viewport, , to equal that of the window otherwise distortions are introduced.

#### 3.3.1.3.8. Compound transformations

| top | pdf |

In this section we consider the viewing transformation T of Section 3.3.1.3.1 and its construction in terms of translation, rotation and scaling, Sections 3.3.1.3.2–4. We use T′ to denote a new transformation in terms of the prevailing transformation T.

We note first that any matrix of the form with U a scalar, may be factorized according to and also that multiplying by an isotropic scaling matrix, a rotation, or a translation, either on the left or on the right, yields a product matrix of the same form, and its inverse is also of this form, i.e. any combination of these three operations in any order may be reduced by the above factorization to a rotation about the original origin, a translation (which defines a new origin) and an expansion or contraction about the new origin, applied in that order.

If is a rotation matrix as in Section 3.3.1.3.3, its application produces a rotation about an axis through the origin defined only in the space in which it is applied. For example, if rotates the image about the z axis of data space, whatever the prevailing viewing transformation, T.

Forming rotates the image about the z axis of display space, i.e. the normal to the tube face under the usual conventions, whatever the prevailing T. Furthermore, if this rotation is to appear to be about some chosen position in the picture, e.g. the centre, then the window transformation U, Section 3.3.1.3.5, must place the origin of display space there by setting , in the notation of that section.

If a rotation is to be about a point then or according to whether R and V are both defined in display space or both in data space. If the rotation is defined in display space and the position of the centre of rotation is defined in data space, then the first form of T′ must be used, in which V is first computed from for a rotation centre at in data space.

For continuous rotations defined in display space it is usually worthwhile to bring the centre of rotation to the origin of display space and leave it there, i.e. to omit the left-most factor in the first expression for T′. Incremental rotations can then be made by further rotational factors on the left without further attention to V. When continuous rotations are implemented by repeated multiplication of T by a rotation matrix, say thirty times a second for a minute or so, the orthogonality of the top-left partition of T may become degraded by accumulation of round-off error and this should be corrected occasionally by one of the methods of Section 3.3.1.2.3.

It is sometimes a requirement, depending on hardware capabilities, to effect a transformation in display space when access to data space is all that is readily available. In such a case where is the required alteration to the prevailing viewing transformation T and is the data-space equivalent,

An important special case is when is to effect a rotation about the origin of display space without change of scale, so that , for then

If r is the required axis of rotation of in display space then the axis of rotation of in data space is since . This gives a particularly simple result if is to be a primitive rotation for then s is the relevant row of R, and can be constructed directly from this and the required angle of rotation.

#### 3.3.1.3.9. Inverse transformations

| top | pdf |

It is frequently a requirement to be able to identify a feature or position in data space from its position on the screen. Facilities for identifying an existing feature on the screen are in many instances provided by the manufacturer as a hit' function which correlates the position indicated on the screen by the user (with a tablet or light pen) with the action of drawing and flags the corresponding item in the drawing internally as having been hit. In other instances it may be necessary to be able to indicate a position in data space independently of any drawn feature and this may be done by setting two or more nonparallel sight lines through the displayed volume and finding their best point of intersection in data space.

In Section 3.3.1.3.1 the relationship between data-space coordinates and screen-space coordinates was given as hence data-space coordinates are given by A line of sight through the displayed volume passing through the point on the screen is the line joining the two position vectors in screen-space coordinates, as in Section 3.3.1.3.7, from which the corresponding two points in data space may be obtained using and in the notation of Section 3.3.1.3.5, and was given in Section 3.3.1.3.8. If orthographic projection is being used then simplifies to Each of these inverse matrices may be suitably scaled to suit the word length of the machine [Section 3.3.1.1.2 (iii)].

Having determined the end points of one sight line in data space the viewing transformation T may then be changed and the required position marked again through the screen in the new orientation. Each such operation generates a pair of points in data space, expressed in homogeneous form, with a variety of values for the fourth coordinate. Each such point must then be converted to three dimensions in the form (X/W, Y/W, Z/W), and for each sight line any (three-dimensional) point pA on the line and the direction qA of the line are established. For each sight line a rank 2 projector matrix MA of order 3 is formed as and the best point of intersection of the sight lines is given by to which three-vector a fourth coordinate of unity may be applied.

#### 3.3.1.3.10. The three-axis joystick

| top | pdf |

The three-axis joystick is a device which depends on compound transformations for its exploitation. As it is usually mounted it consists of a vertical shaft, mounted at its lower end, which can rotate about its own length (the Y axis of display space, Section 3.3.1.3.1), its angular setting, ϕ, being measured by a shaft encoder in its mounting. At the top of this shaft is a knee-joint coupling to a second shaft. The first angle ϕ is set to zero when the second shaft is in some selected direction, e.g. normal to the screen and towards the viewer, and goes positive if the second shaft is moved clockwise when seen from above. The knee joint itself contains a shaft encoder, providing an angle, ψ, which may be set to zero when the second shaft is horizontal and goes positive when its free end is raised. A knob on the tip of the second shaft can then rotate about an axis along the second shaft, driving a third shaft encoder providing an angle θ. The device may then be used to produce a rotation of the object on the screen about an axis parallel to the second shaft through an angle given by the knob. The necessary transformation is then which is in which cos and sin are abbreviated to c and s, which is the standard form with , , .

#### 3.3.1.3.11. Other useful rotations

| top | pdf |

If rotations in display space are to be controlled by trackerball or tablet then there are two measures available, an x and a y, which can define an axis of rotation in the plane of the screen and an angle θ. If x and y are suitably scaled coordinates of a pen on a tablet then the rotation with is about an axis in the xy plane (i.e. the screen face) normal to and with . Applied repetitively this gives a quadratic velocity characteristic. Similarly, if an atom at in display space is to be brought onto the z axis by a rotation with its axis in the xy plane the necessary matrix, in homogeneous form, is with .

#### 3.3.1.3.12. Symmetry

| top | pdf |

In Section 3.3.1.1.1 it was pointed out that it is usual to express coordinates for graphical purposes in Cartesian coordinates in ångström units or nanometres. Symmetry, however, is best expressed in crystallographic fractional coordinates. If a molecule, with Cartesian coordinates, is being displayed, and a symmetry-related neighbour is also to be displayed, then the data-space coordinates must be multiplied by where are the data-space coordinates of the crystallographic origin, M and are as in Section 3.3.1.1.1 and is a crystallographic symmetry operator in homogeneous coordinates, expressed relative to the same crystallographic origin.

For example, in P21 with the origin on the screw dyad along b, and

comprises a proper or improper rotational partition, S, in the upper-left 3 × 3 in the sense that is orthogonal, and with the associated fractional lattice translation in the last column, with the last row always consisting of three zeros and 1 at the 4, 4 position. See IT A (2005, Chapters 5.2 and 8.1 ) for a fuller discussion of symmetry using augmented (i.e. 4 × 4) matrices.

#### 3.3.1.4. Modelling transformations

| top | pdf |

The two sections under this heading are concerned only with the graphical aspects of conformational changes. Determination of such changes is considered under Section 3.3.2.2.

#### 3.3.1.4.1. Rotation about a bond

| top | pdf |

It is a common requirement in molecular modelling to be able to rotate part of a molecule relative to the remainder about a bond between two atoms.

If four atoms are bonded 1–2–3–4 then the dihedral angle in the bond 2–3 is zero if the four atoms are cis planar, and a rotation in the 2–3 bond is, by convention (IUPAC–IUB Commission on Biochemical Nomenclature, 1970), positive if, when looking along the 2–3 bond, the far end rotates clockwise relative to the near end. This is valid for either viewing direction. This sign convention, when applied to the R matrix of Section 3.3.1.2.1, leads to the following statement.

If one of the two atoms is selected as the near atom and the direction cosines are those of the vector from the near atom to the far atom, and if the matrix is to rotate material attached to the far atom (with the reference axes fixed), then a positive rotation in the foregoing sense is generated by a positive θ.

Rotation about a bond normally involves compounding R with translations in the manner of Section 3.3.1.3.8.

#### 3.3.1.4.2. Stacked transformations

| top | pdf |

A flexible molecule may require to be drawn in any of a number of conformations which are related to one another by, for example, rotations about single bonds, changes of bond angles or changes of bond lengths, all of which changes may be brought about by the application of suitable homogeneous transformations during the drawing of the molecule (Section 3.3.1.3.8). With suitable organization, this may be done without necessarily altering the coordinates of the atoms in the coordinate list, only the transformations being manipulated during drawing.

The use of transformations in the manner shown below is straightforward for simply connected structures or structures containing only rigid rings. Flexible rings may be similarly handled provided that the matrices employed are consistent with the consequential constraints as described in Section 3.3.2.2.1, though this requirement may make real-time folding of flexible rings difficult.

Any simply connected structure may be organized as a tree with a node at each branch point and with an arbitrary number of sites of conformational change between one node and the next. We shall call such sites and their associated matrices conformons'. The technique then depends on the stacking technique in which matrices are stored and later recovered in the reverse order of their storage.

One begins at some reference point deemed to be fixed in data space and at this point one stacks the prevailing viewing transformation. From this reference point one advances through the molecule along the structural tree and as each conformon is encountered its matrix is calculated. The product of the prevailing matrix with the conformon matrix is formed and stacked, and this product becomes the prevailing matrix. This product is constructed with the conformon matrix as a factor on the right, i.e. in data space as defined in Section 3.3.1.3.1, and is calculated using the coordinates of the molecule in their unmod­ified form, i.e. before any shape changes are brought about.

This progression leads eventually to an extremity of the tree. At this point drawing is commenced using the prevailing matrix and working backwards towards the fixed root, unstacking (or popping') a matrix as each conformon is passed until a node is reached, which, in general, will occur only part way back to the root. On reaching such a node drawing is suspended and one advances along the newly found branch as before, stacking matrices, until another extremity is reached when drawing towards the root is resumed. This alternation of stacking matrices while moving away from the root and drawing and unstacking matrices while moving towards the root is continued until the whole tree is traversed.

This process is illustrated schematically in Fig. 3.3.1.2 for a simple tree with one node, numbered 1, and three conformons at a, b and c. One enters the tree with a current viewing transformation T and progresses upwards from the fixed lower extremity. When the conformon at a is encountered, T is stacked and the product is formed. Continuing up the tree, at node 1 either branch may be chosen; we choose the left and, on reaching b, is stacked and is formed. On reaching the tip drawing down to b is done with this transformation, is then unstacked and drawing continues with this matrix until node 1 is reached. The other branch is then followed to c whereupon is again stacked and the product is formed. From the tip down as far as c is drawn with this matrix, whereupon is unstacked and drawing continues down to a, where T is unstacked before drawing the section nearest the root.

 Figure 3.3.1.2 | top | pdf |Schematic representation of a simple branched-chain molecule with a stationary root and two extremities. The positions marked a, b and c are the loci of possible conformational change, here called conformons, and there is a single, numbered branch point.

With this organization the matrices associated with b and c are unaffected by changes in the conformation at a, notwithstanding the fact that changes at a alter the direction of the axis of rotation at b or c.

Two other approaches are also possible. One of these is to start at the tip of the left branch, replace the coordinates of atoms between b and the tip by , and later replace all coordinates between the tip and a by , with a similar treatment for the other branch. The advantage of this is that no storage is required for stacked matrices, but the disadvantage is that atoms near the tips of the tree have to be reprocessed for every conformon. It also modifies the stored coordinates, which may or may not be desirable.

The second alternative is to draw upwards from the root using T until a is reached, then using until b is reached, then using to the tip, but in this formulation must be based on the geometry that exists at b after the transformation has been applied to this region of the molecule, i.e. is characteristic of the final conformation rather than the initial one.

| top | pdf |

#### 3.3.1.5.1. Types of hardware

| top | pdf |

There are two main types of graphical hardware in use for interactive work, in addition to plotters used for batch work. These main types are raster and vector. In raster equipment the screen is scanned as in television, with a grid of points, called pixels, addressed sequentially as the scan proceeds. Associated with each pixel is a word of memory, usually containing something in the range of 1 to 24 bits per pixel, which controls the colour and intensity to be displayed. Many computer terminals have one bit per pixel (said to be single-plane' systems) and these are essentially monochrome and have no grey scale. Four-plane systems are cheap and popular and commonly provide 4-bit by 4-bit look-up tables between the pixel memory and the monitor with one such table for each of the colours red, green and blue. If these tables are each loaded identically then 16 levels of monochrome grey scale are available, but if they are loaded differently 16 different colours are available simultaneously chosen from a total of 4096 possibilities. Four-plane systems are adequate for many applications where colour is used for coding, but are inadequate if colour is intended also to provide realism, where brilliance and saturation must be varied as well as hue. For these applications eight-plane systems are commonly used which permit 256 colours chosen from 16 million using three look-up tables, though the limitations of these can also be felt and full colour is only regarded as being available in 24-plane systems.

Raster-graphics devices are ideal for drawing objects represented by opaque surfaces which can be endowed with realistic reflecting properties (Max, 1984) and they have been successfully used to give effects of transparency. They are also capable of representing shadows, though these are generally difficult to calculate (see Section 3.3.1.5.5). Many devices of this type provide vectorization, area fill and anti-aliasing. Vectorization provides automatically for the loading of relevant pixels on a straight line between specified points. Area fill automatically fills any irregular pre-defined polygon on the screen with a uniform colour with the user specifying only the colour and one point within the polygon. Anti-aliasing is the term used for a technique which softens the staircase effect that may be seen on a line which runs at a small angle to a vertical or horizontal row of pixels.

The main drawback with this type of equipment is that it is slow compared to vector machines. Only relatively simple objects can be displayed with smooth rotation in real time as transformed coordinates have to be converted to pixel addresses and the previous frame needs to be deleted with each new frame unless it is known that each new frame will specify every pixel. However, the technology is advancing rapidly and these restrictions are already disappearing.

Vector machines, on the other hand, are specialized to drawing straight lines between specified points by driving the electron beam along such lines. No time is wasted on blank areas of the screen. Dots may be drawn with arbitrary coordinates, in any order, but areas, if they are to be filled, must be done with a ruling technique which is very seldom done. Images produced by vector machines are naturally transparent in that foreground does not obscure background, which makes them ideal for seeing into representations of molecular structure.

#### 3.3.1.5.2. Optimization of line drawings

| top | pdf |

A line drawing consisting of n line segments may be specified by anything from to 2n position vectors depending on whether the lines are end-to-end connected or independent. Appreciable gains in both processing time and storage requirements may be made in complicated drawings by arranging for line segments to be end-to-end connected as far as possible, and an algorithm for doing this is outlined below. For further details see Diamond (1984a).

Supposing that a list of nodal points (atoms if a covalent skeleton is being drawn) exists within a computer with each node appearing only once and that the line segments to be drawn between them are already determined, then at each node there are, generally, both forward and backward connections, forward connections being those to nodes further down the list. A quantity D is calculated at each node which is the number of forward connections minus the number of backward connections. At the commencement of drawing, the first connected node in the list must have a positive D, the last must have a negative D, the sum of all D values must be zero and the sum of the positive ones is the number of strokes required to draw the drawing, a stroke' being a sequence of end-to-end connected line segments drawn without interruption. The total number of position vectors required to specify the drawing is then the number of nodes plus the number of strokes plus the number of rings minus one.

Drawing should then be done by scanning the list of nodes from the top looking for a positive D (usually found at the first node), commencing a stroke at this node and decrementing its D value by 1. This stroke is continued from node to node using the specified connections until a negative D is encountered, at which point the stroke is terminated and the D value at the terminating node is incremented by 1. This is done even though this terminating node may also possess some forward connections, as the total number of strokes required is not minimized by keeping a stroke going as far as possible, but by terminating a stroke as soon as it reaches a node at which some stroke is bound to terminate.

The next stroke is initiated by resuming the scan for positive D values at the point in the node list where the previous stroke began. If this scan encounters a zero D value at a node which has not hitherto been drawn to, or drawn from, then the node concerned is isolated and not connected to any other, and such nodes may require to be drawn with some special symbol. The expression already given for the number of vectors required is valid in the presence of isolated nodes if drawing an isolated node is allowed one position vector, this vector not being counted as a stroke.

The number of strokes generated by this algorithm is sensitive to the order in which the nodes are listed, but if this resembles a natural order then the number of strokes generated is usually close to the minimum, which is half the number of nodes having an odd number of connections. For example, the letter has six nodes, four of which have an odd number of connections, so it may be drawn with two strokes.

#### 3.3.1.5.3. Representation of surfaces by lines

| top | pdf |

The commonest means of representing surfaces, especially contour surfaces, is to consider evenly spaced serial sections and to perform two-dimensional contouring on each section. Repeating this on serial sections in two other orientations then provides a good representation of the surface in three dimensions when all such contours are displayed. The density is normally cited on a grid with submultiples of a, b and c as grid vectors, inverse linear interpolation being used between adjacent grid points to locate points on the contour. For vector-graphics applications it is expedient to connect such points with straight lines; some equipment may be capable of connecting them with splines though this is burdensome or impossible if real-time rotation of the scene is required. Precalculation of splines stored as short vectors is always possible if the proliferation of vectors is acceptable. For efficient drawing it is necessary for the line segments of a contour to be end-to-end connected, which means that it is necessary to contour by following contours wherever they go and not by scanning the grid. Algorithms which function in this way have been given by Heap & Pink (1969) and Diamond (1982a). Contouring by grid scanning followed by line connection by the methods of the previous section would be possible but less efficient. Further contouring methods are described by Sutcliffe (1980) and Cockrell (1983).

For raster-graphics devices there is little disadvantage in using curved contours though many raster devices now have vectorizing hardware for loading a line of pixels given only the end points. For these devices well shaped contours may be computed readily, using only linear arithmetic and a grid-scanning approach (Gossling, 1967). Others have colour-coded each pixel according to the density, which provides a contoured visual impression without performing contouring (Hubbard, 1983).

#### 3.3.1.5.4. Representation of surfaces by dots

| top | pdf |

Connolly (Langridge et al., 1981; Connolly, 1983a,b) represents surfaces by placing dots on the surface with an approximately uniform superficial density. Connolly's algorithm was developed to display solvent-accessible surfaces of macromolecules and provides for curved concave portions where surface atoms meet. Pearl & Honegger (1983) have developed a similar algorithm, based on a grid, which generates only convex portions which meet in cusps, but is faster to compute than the Connolly surface. Bash et al. (1983) have produced a van der Waals surface algorithm fast enough to permit real-time changes to the structure without tearing the surface.

It has become customary to use a dot representation to display computed surfaces, such as the surface at a van der Waals radius from atomic centres, and to use lines to represent experimentally determined surfaces, especially density contours.

#### 3.3.1.5.5. Representation of surfaces by shading

| top | pdf |

Many techniques have been developed, mainly for raster-graphics devices, for representing molecular surfaces and these have been very well reviewed by Max (1984).

The simplest technique in this class consists in representing each atom by a uniform disc, or high polygon, which can be colour-coded and area-filled by the firmware of the device. If such atoms are sorted on their z coordinate and drawn in order, furthest ones first, so that nearer ones partly or completely overwrite the further ones then the result is a simple representation of the molecule as seen from the front. This technique is fast and has its uses when a rapid schematic is all that is required. In one sense it is wasteful to process distant atoms when they are going to be overwritten by foreground atoms, but front-to-back processing requires the boundaries of visible parts of partially obscured atoms near the front to be determined before they can be painted or, alternatively, every pixel must be tested before loading to see if it is already loaded. Not only does this approach give a uniform rendering over the whole area of one atom, it also gives a boundary between overlapping atoms with almost equal z values which completes the circle of the nearer atom, though it should be an arc of an ellipse when the atoms are drawn with radii exceeding their covalent radii.

Greater realism is achieved by establishing a z buffer, which is an additional area of memory with one word per pixel, in which is stored the z value of the currently loaded feature in each pixel. Treatments which take account of the sphericity are then possible and correct arcs of intersection for interpenetrating spheres and more complicated entities arise naturally through loading a colour value into a pixel only if the z coordinate is less than that of the currently loaded value. This z buffer and the associated x, y coordinates should be in picture space or screen space rather than display space since only after the application of perspective can points with the same x/w and y/w coordinates obscure one another.

It is usual in such systems to vary the intensity of colour within one atom by darkening it towards the circumference on the basis of the z coordinate. Some systems augment this impression of sphericity by highlighting. The simplest form of highlighting is an extension of the uniform disc treatment in which additional, brighter discs, possibly off centre, are associated with each atom. More general highlighting (Phong, 1975) is computed from four unit vectors, these being the normal to the surface, the direction to a light source, the direction to the viewer and the normalized vector sum of these last two. Intensity levels may then be set as the sum of three terms: a constant, a term proportional to the scalar product of the first two vectors (if positive) and a term proportional to a high power of the scalar product of the first and last vectors; the higher the power the glossier the surface appears to be. This final term normally adds a white term, rather than the surface colour, supposing the light source to be white.

Shadows may also be rendered to give even greater realism. In addition to the z buffer and (x, y) frame buffer a second z buffer for z′ values associated with x′ and y′ is also required. These coordinates are then related by , , z′ = z. The second buffer is a ray buffer since xy′ are the coordinates with which an illuminating ray passing through (xyz) passes through the z = 0 plane, and z′, stored at x′, y′, records the depth at which this ray encounters material. Thus any two pixels and are on the same illuminating ray if their x′ and y′ values are equal and the one with smaller z′ shadows the other. Processing a pixel at therefore involves first determining its visibility on the basis of the z buffer, as before, then, whether or not it is visible, setting and considering the value of z′ currently stored at xy′, which we call .

If then is in light and must be loaded accordingly. From we find the previously processed pixel which is now in shade and which was in light when originally processed, so that the colour value stored at needs to be altered unless the pixel at is now with , in which case the pixel which has now become shadowed by has, in the meantime, been obscured by which is not shadowed by and no change is therefore needed. In either event then replaces .

If then , if visible, is in shade and must be coloured accordingly, and in this case is not superseded.

This shadowing scheme corresponds to illumination by a light source at infinity in picture space or, equivalently, with a z coordinate equal to that of the eye in display space. For its implementation x, y and z may be in any convenient coordinate system, e.g. pixel addresses, but if x and y are expressed with the range −1 to 1 and z with the range 0 to 1 corresponding to the window then they may be identified as the quantities x/w, y/w and z/w of picture space (Section 3.3.1.3.1).

If, in the notation of Section 3.3.1.3.5, the light source is placed at (P, Q, E, V) in display space and a ray leaves it in the direction (p, q, r, V) thenwhich varies only with beam direction, and similarly for y′ and β.

#### 3.3.1.5.6. Advanced hidden-line and hidden-surface algorithms

| top | pdf |

Hidden surfaces may be handled quite generally with the z-buffer technique described in the previous section but this technique becomes very inefficient with very complicated scenes. Faster techniques have been developed to handle computations in real time (e.g. 25 frames s−1) on raster machines when both the viewpoint and parts of the environment are moving and substantial complexity is required. These techniques generally represent surfaces by a number of points in the surface, connected by lines to form panels. Many algorithms require these panels to be planar and some require them to be triangular. Of those that permit polygonal panels, most require the polygons to be convex with no re-entrant angles. Yet others are limited to cases where the objects themselves are convex. Some can handle interpenetrating surfaces, others exclude these. Some make enormous gains in efficiency if the objects in the scene are separable by the insertion of planes between them and degrade to lower efficiency if required, for example, to draw a chain. Some are especially suited to vector machines and others to raster machines, the latter capitalizing on the finite resolution of such systems. In all of these the basic entities for consideration are entire panels or edges, and in some cases vertices, point-by-point treatment of the entire surface being avoided until after all decisions are made concerning what is or is not visible.

All of these algorithms strive to derive economies from the notion of coherence'. The fact that, in a cine context, one frame is likely to be similar to the previous frame is referred to as frame coherence'. In raster scans line coherence also exists, and other kinds of coherence can also be identified. The presence of any form of coherence may enable the computation to be concerned primarily with changes in the situation, rather than with the totality of the situation so that, for example, computation is required where one edge crosses in front of another, but only trivial actions are involved so long as scan lines encounter the projections of edges in the same order.

The choice of technique from among many possibilities may even depend on the viewpoint if the scene has a statistical anisotropy. For example, the depiction of a city seen from a viewpoint near ground level involves many hidden surfaces. Distant buildings may be hidden many times over. The same scene depicted from an aerial viewpoint shows many more surfaces and fewer overlaps. This difference may swing the balance of advantage between an algorithm which sorts first on z or one which leaves that till last.

These advanced techniques have, so far, found little application in crystallography, but this may change. Ten such techniques are critically reviewed and compared by Sutherland et al. (1974), and three of these are described in detail by Newman & Sproull (1973).

### References

American National Standards Institute, American National Standard for Information Processing Systems – Computer Graphics – Graphical Kernel System (GKS) Functional Description (1985). ISO 7942, ISO Central Secretariat, Geneva, Switzerland.
American National Standards Institute, American National Standard for Information Processing Systems – Computer Graphics – Programmer's Hierarchical Graphics System (PHIGS) Functional Description, Archive File Format, Clear-Text Encoding of Archive File (1988). ANSI X3.144–1988. ANSI, New York, USA.
International Standards Organisation, International Standard Information Processing Systems – Computer Graphics – Graphical Kernel System for Three Dimensions (GKS-3D), Functional Description (1988). ISO Document No. 8805:1988(E). American National Standards Institute, New York, USA.
International Tables for Crystallography (2005). Vol. A. Space-Group Symmetry, edited by Th. Hahn. Heidelberg: Springer.
IUPAC–IUB Commission on Biochemical Nomenclature (1970). Abbreviations and symbols for the description of the conformation of polypeptide chains. J. Biol. Chem. 245, 6489–6497.
Abi-Ezzi, S. S. & Bunshaft, A. J. (1986). An implementer's view of PHIGS. IEEE Comput. Graphics Appl. Vol. 6, Part 2.
Aharonov, Y., Farach, H. A. & Poole, C. P. (1977). Non-linear vector product to describe rotations. Am. J. Phys. 45, 451–454.
Arnold, D. B. & Bono, P. R. (1988). CGM and CGI: Metafile and Interface Standards for Computer Graphics. Berlin: Springer-Verlag.
Bash, P. A., Pattabiraman, N., Huang, C., Ferrin, T. E. & Langridge, R. (1983). Van der Waals surfaces in molecular modelling: implementation with real-time computer graphics. Science, 222, 1325–1327.
Brown, M. D. (1985). Understanding PHIGS. Template, Megatek Corp., San Diego, California, USA.
Cockrell, P. R. (1983). A new general purpose method for large volume production of contour charts. Comput. Graphics Forum, 2, 35–47.
Connolly, M. L. (1983a). Solvent-accessible surfaces of proteins and nucleic acids. Science, 221, 709–713.
Connolly, M. L. (1983b). Analytical molecular surface calculation. J. Appl. Cryst. 16, 548–558.
Dam, A. van (1988). PHIGS+ functional description, revision 3.0. Comput. Graphics, 22, 125–218.
Diamond, R. (1966). A mathematical model-building procedure for proteins. Acta Cryst. 21, 253–266.
Diamond, R. (1976a). On the comparison of conformations using linear and quadratic transformations. Acta Cryst. A32, 1–10.
Diamond, R. (1982a). Two contouring algorithms. In Computational Crystallography, edited by D. Sayre, pp. 266–272. Oxford University Press.
Diamond, R. (1984a). Applications of computer graphics in molecular biology. Comput. Graphics Forum, 3, 3–11.
Diamond, R. (1988). A note on the rotational superposition problem. Acta Cryst. A44, 211–216.
Diamond, R. (1989). A comparison of three recently published methods for superimposing vector sets by pure rotation. Acta Cryst. A45, 657.
Diamond, R. (1990a). On the factorisation of rotations with special reference to diffractometry. Proc. R. Soc. London Ser. A, 428, 451–472.
Diamond, R. (1990b). Chirality in rotational superposition. Acta Cryst. A46, 423.
Diamond, R. (1992). On the multiple simultaneous superposition of molecular structures by rigid body transformations. Protein Sci. 1, 1279–1287.
Diamond, R. (1995). Coordinate based cluster analysis. Acta Cryst. D51, 127–135.
Diamond, R., Wynn, A., Thomsen, K. & Turner, J. (1982). Three-dimensional perception for one-eyed guys, or, the use of dynamic parallax. In Computational Crystallography, edited by D. Sayre, pp. 286–293. Oxford University Press.
Enderle, G., Kansy, K. & Pfaff, G. (1984). Computer Graphics Programming, GKS – the Graphics Standard. Berlin: Springer-Verlag.
Foley, J. D., van Dam, A., Feiner, S. K. & Hughes, J. F. (1990). Computer Graphics Principles and Practice, 2nd ed. New York: Addison Wesley.
Gossling, T. H. (1967). Two methods of presentation of electron-density maps using a small-store computer. Acta Cryst. 22, 465–468.
Harris, M. R., Geddes, A. J. & North, A. C. T. (1985). A liquid crystal stereo-viewer for molecular graphics. J. Mol. Graphics, 3, 121–122.
Heap, B. R. & Pink, M. G. (1969). Three contouring algorithms, DNAM Rep. 81. National Physical Laboratory, Teddington, England.
Hopgood, F. R. A., Duce, D. A., Gallop, J. R. & Sutcliffe, D. C. (1986). Introduction to the Graphical Kernel System, 2nd ed. London: Academic Press.
Hubbard, R. E. (1983). Colour molecular graphics on a microcomputer. J. Mol. Graphics, 1, 13–16, C3–C4.
Johnson, C. K. (1970). Drawing crystal structures by computer. In Crystallographic Computing, edited by F. R. Ahmed, pp. 227–230. Copenhagen: Munksgaard.
Kabsch, W. (1976). A solution for the best rotation to relate two sets of vectors. Acta Cryst. A32, 922–923.
Kabsch, W. (1978). A discussion of the solution for the best rotation to relate two sets of vectors. Acta Cryst. A34, 827–828.
Kearsley, S. K. (1989). On the orthogonal transformation used for structural comparisons. Acta Cryst. A45, 208–210.
Kearsley, S. K. (1990). An algorithm for the simultaneous superposition of a structural series. J. Comput. Chem. 11, 1187–1192.
Langridge, R., Ferrin, T. E., Kuntz, I. D. & Connolly, M. L. (1981). Real-time color graphics in studies of molecular interactions. Science, 211, 661–666.
Mackay, A. L. (1984). Quaternion transformation of molecular orientation. Acta Cryst. A40, 165–166.
Max, N. L. (1984). Computer representation of molecular surfaces. J. Mol. Graphics, 2, 8–13, C2–C4.
McLachlan, A. D. (1972). A mathematical procedure for superimposing atomic coordinates of proteins. Acta Cryst. A28, 656–657.
McLachlan, A. D. (1979). Gene duplications in the structural evolution of chymotrypsin. Appendix: Least squares fitting of two structures. J. Mol. Biol. 128, 49–79.
McLachlan, A. D. (1982). Rapid comparison of protein structures. Acta Cryst. A38, 871–873.
Newman, W. M. & Sproull, R. F. (1973). Principles of Inter-active Computer Graphics. New York: McGraw-Hill.
Pearl, L. H. & Honegger, A. (1983). Generation of molecular surfaces for graphic display. J. Mol. Graphics, 1, 9–12, C2.
Phong, B. T. (1975). Illumination for computer generated images. Commun. ACM, 18, 311–317.
Shapiro, A., Botha, J. D., Pastore, A. & Lesk, A. M. (1992). A method for multiple superposition of structures. Acta Cryst. A48, 11–14.
Sutcliffe, D. C. (1980). Contouring over rectangular and skewed rectangular grids – an introduction. In Mathematical Methods in Computer Graphics and Design, edited by K. W. Brodie, pp. 39–62. London: Academic Press.
Sutherland, I. E., Sproull, R. F. & Schumacker, R. A. (1974). A characterization of ten hidden surface algorithms. Comput. Surv. 6, 1–55.
Thomas, D. J. (1993). Toward more reliable printed stereo. J. Mol. Graphics, 11, 15–22.