POLARRFN (CCP4: Supported Program)
NAME
polarrfn
- fast rotation function that works in polar angles, updated to work with up to the 100th order of the spherical harmonics
SYNOPSIS
polarrfn hklin
foo_in_Fobs.mtz
hklin2
foo_in_Fmodel.mtz
mapin
bar.map
mapout
foo.map
PLOT
fred.plt
[Keyworded input]
DESCRIPTION
This is a fast rotation function that works in polar angles, written by
W. Kabsch. The program produces sections of constant rotation angle kappa for
different axis directions defined by omega (angle from pole) and phi (angle
around equator). The direction cosines of the axis is given by:
( sin omega cos phi )
( sin omega sin phi )
( cos omega )
These sections may be plotted as stereograms (PLOT option).
For example, the kappa = 180 section can be used to identify 2-fold axes in
a self-rotation function.
The program searches for peaks in the rotation function (FIND
option), and these are listed in the log output in both polar and Eulerian angles.
For each peak, peaks related by crystallographic symmetry are also listed.
If the peak, or one of its symmetry relations, corresponds to the identity
operation (i.e. zero rotation), then this is indicated by "Origin peak" in the
output. For a self-rotation function, peaks corresponding to genuine non-crystallographic
symmetry will appear after the origin peaks.
Options are provided in the program for writing the calculated rotation function
to a map file, that can be read by the program on another occasion for
plotting, peak searching etc (options MAP, READ).
The results of the first stage
of the calculation (ALMN) may also be saved for subsequent rotation function
summations on different grids (SAVE, SUM). A list of peaks may be read in,
and all symmetry related peaks generated (option PEAK).
Compared to the Crowther program (ALMN), this program is a little
slower but a self rotation function in polar angles is easier to understand.
This program can also be used to calculate a small part of the rotation
function on a fine grid (down to 1 degree steps). However, the crystal symmetry
is expressed less fully in polar angle space than in Eulerian angle space, so
for cross-rotation functions the Eulerian program is normally more appropriate.
This version of POLARRFN was modified for CCP4 5.0 to use up to the 100th
order of spherical harmonics (previously 30). This allows a larger radius
of integration and/or a higher resolution cutoff to be used.
KEYWORDED INPUT
Only the first 4 letters of each keyword are necessary. All input is free
format. Most data cards are optional. The various data control lines are
identified by keywords, those available being:
CROSS, CRYSTAL,
END, FIND,
JOIN, LABIN,
LIMITS, MAP,
MAXORD,
NOPRINT, PEAK,
PLOT, READ,
RESOLUTION, SAVE,
SELF, SUM, TITLE
The following keywords are compulsory:
SELF/CROSS, CRYSTAL, LABIN.
TITLE <title>
The rest of the line is taken as a title.
SELF <arad> <eps>
Set flag for self rotation (cf CROSS)
- <arad>
-
integration radius in Patterson space (default 20Å)
- <eps>
-
Radius (Angstrom) of sphere within which the
Patterson is removed. Default Eps set = 0.7*<resmax>.
CROSS <arad> <eps>
Set flag for cross rotation (cf SELF)
- <arad>
-
integration radius in Patterson space (default 20Å)
- <eps>
-
Radius (Angstrom) of sphere within which the
Patterson is removed. Default Eps set = 0.7*<resmax>.
RESOLUTION <resmin> <resmax>
Read resolution limits in Å, in either order. If only one treat
as high resolution limit. Note that because of the limit of 100 in the order of
the spherical harmonics in this version, the high resolution limit cannot be numerically
less than the integration radius (arad) / 17.4: if it is, the program resets <resmax>
to <arad> / 17.4 .
CRYSTAL
FILE <nfile> ORTH <ncode> BFAC <badd> FLIMITS <fmin> <fmax>
[CELL SYMMETRY]
This card is COMPULSORY for each hklin set (one for self-rotation,
two for cross-rotation). The first subkeyword FILE (or NUMBER) indicates
whether this CRYSTAL card refers to crystal 1 or to crystal 2.
- FILE
-
followed by crystal number <nfile>. This sets the crystal number, 1 or
2 (default 1). The syntax is any one of FILE 1 or FILE 2.
- ORTH
-
followed by orthogonalisation code <ncode> (default 1). If the PEAK
option is used ORTH should be given. Ncode defaults to 1 if
unspecified.
ncode orthogonalisation code for this crystal
= 1, orthogonal x y z along a,c*xa,c* (Brookhaven, default)
= 2 b,a*xb,a*
= 3 c,b*xc,b*
= 4 a+b,c*x(a+b),c*
= 5 a*,cxa*,c (Rollett)
- BFAC
-
followed by temperature factor <badd> (default = 0). This can be used
to sharpen the data (e.g. = -20). Fused = Finput * exp(-badd*ssq/lsq).
A better alternative is to use normalised amplitudes (from the program ECALC),
in which case the BFAC option must NOT be used.
- FLIMITS
-
followed by <fmin> and <fmax> (defaults to no cutoffs).
These tests are done after application of the temperature
factor, but before squaring to intensities.
- SYMMETRY
-
followed by symmetry operation or crystal space group name
or number. Default: take from MTZ file.
- CELL
-
followed by cell dimensions. Default: take from MTZ file.
LABIN <data assignments>
LABIN FILE <nfile> F=<label> SIGF=<label2>
File column data assignments for F, SIGF . If you use the output from ECALC
(recommended), set F=E and SIGF=SIGE .
FILE sets the crystal number, 1 or 2 (default 1).
The syntax is any one of FILE 1 or FILE 2.
LIMITS <iphi1> <iphi2> <nphi> <iomega1> <iomega2> <nomega> <ikappa1> <ikappa2> <nkappa>
i.e. limits and steps on phi, omega, kappa
iphi1 iomega1 ikappa1 start points in degrees
iphi2 iomega2 ikappa2 stop points in degrees
nphi nomega nkappa intervals in degrees
The default start & stop points are 0 & 180 which covers the basic asymmetric unit for the
self-rotation function; the default step size is set numerically equal to the high resolution
cutoff of the data, or that specified by the RESOLUTION keyword,
whichever is numerically the greater. It is recommended to omit specification of the angle
limits, at least for self-rotation functions.
[NO]PRINT <lprint>
Switch on or off printing of map as a figure-field in the logfile, and
set print threshold LPRINT. The rotation function is scaled to a
maximum of 100. Values less than 0 are printed as ' -' and values
between 0 and LPRINT are printed as ' .'. NOPRINT (= PRINT 0)
suppresses printing of the figure-field. The default is to print the
map (PRINT +1), except in the case of the READ option, for which the
default is NOPRINT.
MAP
Write rotation function to a map file. This may be read back later
for plotting, peak searching, etc. using the READ option.
MAXORD <maxord>
This is primarily for test purposes, and should not normally be set. If you set MAXORD 30,
you should get the same result as the original 30-order version, all else being equal. This
will have the effect of placing a numerically greater limit on the high resolution cutoff. Note
however that the old version does not set the default angle step equal to the high resolution
cutoff of the data, and it does not take the latter into account when determining the
maximum order to be used. What this means in practice is that to get the same result you
should set the resolution cutoffs and angular steps explicitly, as the default values used by
the two versions of the program are likely to be different.
PLOT <cstart> <cint> [ NOTITLE ]
Plot rotation function as a stereographic projection of each kappa
section. These sections have omega = 0 or 180 at the centre,
omega = 90 around the edge, and phi as marked around the periphery.
- <cstart>
-
contour start level (default = 10)
- <cint>
-
contour interval (default = 10)
It is recommended to set the interval equal to the RMS deviation of the map, and the start
level to twice this (this means you have to run the program twice, the first time to get
the value of the RMS deviation printed near the end of the log file and then again to get
the correctly contoured plot).
If the keyword NOTITLE is present, the maps are plotted with no
labels or titles: this may be useful for publication purposes.
The plots are not normally useful if a small part of the function is
being calculated.
JOIN <kappa> <omega1> <phi1> <omega2> <phi2> [DASHED [<repeat> <dash>]]
Draw arc connecting two peaks at (omega1,phi1) and (omega2,phi2) on
section kappa.
If the keyword DASHED is present, the arc will be dashed with repeat
and dash length as specified (in multiples of the radius of the stereogram).
Default values for <repeat> and <dash> are 0.05, 0.025.
Up to 30 JOIN commands may be given.
Warning: this option is primarily intended for self-rotation functions,
and probably will not work sensibly if the kappa sections are split
into two sections or half sections.
FIND <peak> <maxpeak> [RMS] [OUTPUT <filename>]
Read peak threshold and maximum number of peaks.
Peak searching may be suppressed by PEAK 0 0.
Note that although the program attempts to interpolate the positions
of peaks, this may not be very accurate in regions where the peaks
are very spread out (kappa near 0) or on the edges of the map.
- <peak>
-
threshold for peaks (default = 10).
- <maxpeak>
-
maximum number of peaks to find (default = 10)
Up to <maxpeak> peaks above <peak> will be found,
and all symmetry related peaks generated
If the keyword RMS is present, then the peak threshold is
<peak> * RMS density, otherwise <peak> is the absolute threshold
in the scaled map.
If the keyword OUTPUT is present, the unique peaks will be
output to a file whose name may also be given (default filename
or logical name PEAKS). The keyword RMS must precede OUTPUT if
both are present.
READ <kappa1> <kappa2>
Instead of calculating a rotation function, read a previously calculated
one (written out using the MAP flag). Sections <kappa1> to <kappa2> will
be processed (if no values are given, the whole map will be used),
printed (if PRINT is on, NOPRINT is default), plotted (if a PLOT card
is present) and searched for peaks (see FIND).
SAVE [ <fname> ]
save output from first stage of rotation (calculation of Almn
coefficients for each crystal) in files <fname>1.dat and <fname>2.dat.
<fname> defaults to 'COEFFS'. These files can be used to save time on a
later run of the program (see SUM card) to calculate a rotation
function on a different grid (e.g. on a fine grid around a peak).
SUM [ <fname> ]
read Almn coefficients from files generated in a previous run of the
program with the SAVE option. The files are <fname>1.dat and <fname>2.dat
<fname> defaults to 'COEFFS'
PEAK [ MATRICES ]
read a list of peak positions (omega, phi, kappa in degrees) from
subsequent lines, terminated by blank line, 'END' or end of file.
The program will not calculate a rotation function, but will generate
all symmetry related peaks from the peaks given, and print out the
corresponding matrices. This option requires cards SELF or CROSS,
CRYSTAL or CELL, and SYMMETRY (ie the number of crystals, their
symmetry and orthogonalisation codes must be defined). If the keyword
MATRICES is present, the rotation matrices corresponding to each
peak will be printed.
END
Terminates input. If present must be the last card.
A NOTE ON SYMMETRY AND ASYMMETRIC UNITS
The polar angles omega, phi, kappa have an intrinsic symmetry operation
180-omega, 180+phi, -kappa by definition. This gives an asymmetric unit in the
absence of other symmetry of omega 0 to 180, phi 0 to 360, kappa 0 to 180.
Additional symmetry relations between omega and phi for a given kappa are
created in self rotation functions, and by parallel dyad axes in the two
crystals (ie any dyad axis in a self rotation function). These are as follows
-
Self-rotation function, and also kappa=180 in all cases:
kappa = -kappa, hence 180-omega, 180+phi
-
dyad axes in both crystals parallel to orthogonal x axes (after any
permutation produced by the orthogonalisation code NCODE):
180 - omega, -phi
-
dyad axes parallel to orthogonal y axes:
180 - omega, 180 - phi
-
dyad axes parallel to orthogonal z axes:
omega, 180 + phi
If more than one of these symmetries is present, additional symmetry is
generate by their combination. The asymmetric units required are given in the
following table:
Symmetry Asymmetric unit
-------- omega phi kappa
no symmetry 180 360 180
1. 180 - omega, 180 + phi 90 360 180
or 180 180 180
2. 180 - omega, -phi 90 360 180
or 180 180 180
3. 180 - omega, 180 - phi 90 360 180
or 180 -90) 180
to +90)
4. omega, 180 + phi 180 180 180
In the program, the PLOT option will always try to plot a complete stereogram
of all kappa sections. In the absence of symmetry, each kappa section is split
into two hemispheres, omega 0 to 90 and omega 90 to 180. For symmetries 1, 2
and 3, only the hemisphere omega 0 to 90 is plotted. For symmetry 4, the two
quarter spheres phi 0 to 180 for omega 0 to 90 and omega 90 to 180 are plotted
together in the same circle. The PLOT option is not normally useful if only a
small part of the asymmetric unit is calculated.
INPUT AND OUTPUT FILES
- HKLIN
-
Input amplitudes for crystal 1, or for sole crystal in the case of
self-rotation function.
- HKLIN2
-
Input amplitudes for crystal 2, for cross-rotation function.
- MAPIN
-
Input map for READ option.
- MAPOUT
-
Output map containing rotation function, for MAP option.
- PLOT
-
Plotter output containing stereographic projection of each kappa
section, for PLOT option.
- COEFFS1
-
Almn coefficients for crystal 1. Default created as scratch file.
- COEFFS2
-
Almn coefficients for crystal 2. Default created as scratch file.
EXAMPLES
Self-rotation
polarrfn hklin tpfktrunc.mtz
mapout tself.map
plot self.plt
<< eof
TITLE T-PFK self-rotation R=29A, 5 - 7 A
SELF 29
RESOLUTION 7 5
CRYSTAL FILE 1 ORTH 1 BFAC -20 SYMMETRY P21
LABIN FILE 1 F=FP SIGF=SIGFP
MAP
PLOT 10 10
FIND 2 20 RMS OUTPUT peaks.dat
eof
#!/bin/csh -f
#
# Self-rotation function
#
set resolution = 15 3
set radius = 20
set u1af = {$scr0}/u1afobs_fc117
# Calculate whole map, & search for peaks
polarrfn hklin {$u1af} mapout u1aself << eof-1
title U1A Self-Rotation function, resolution ${resolution} radius ${radius}
self ${radius}
resolution ${resolution}
crystal file 1 bfac -20
labin file 1 F=FP SIGF=SIGFP
map
find 5 100
eof-1
# Now plot just the kappa = 180 deg section
plot:
polarrfn hklin {$u1af} mapin u1aself plot u1aself << eof-2
title U1A Self-Rotation function, resolution ${resolution} radius ${radius}
self ${radius}
resolution ${resolution}
crystal file 1 bfac -20
labin file 1 F=FP SIGF=SIGFP
read 180 180
plot 10 10
eof-2
Cross rotation
polarrfn hklin tpfktrunc
hklin2 uniqsfpfk.mtz
mapout tcross.map
<< eof
TITLE T-PFK cross-rotation R=29A, 5 - 7 A
CROSS 29
RESOLUTION 7 5
CRYSTAL FILE 1 ORTH 1 BFAC -20 SYMMETRY 4
LABIN FILE 1 F=FP SIGF=SIGFP
CRYSTAL FILE 2 ORTH 1 BFAC -20 SYMMETRY 1
LABIN FILE 2 F=FC SIGF=FC
LIMITS 0 180 5 0 360 5 0 180 5
MAP
FIND 10 20
eof
AUTHORS
Author: W. Kabsch
Contacts: Phil Evans (pre@mrc-lmb.cam.ac.uk)
Ian Tickle (i.tickle@astex-technology.com)