Mark Dow



P22 virion structure

Test renderings of the P22 virion EMD-1220 data set

EMD-1220 rendering notes

"For a physicist, it is indeed a great joy to learn how we can use beautiful mathematics to understand the real world."  In Stephen Weinberg's "Without God" (9/25/08, The New York Review of Books)

    The EMD_1220 volume data set was reconstructed from cryo-electron microscopy images. It represents, roughly speaking, a three dimensional estimate of electron densities of an average P22 virion particle.

    These notes outline post-processing that I performed in order to do renderings of the structure and sub-structures represented in the EMD_1220 volume.

Data source

Processing notes

Convert from CCP4 map (.map/.xml) to Space volume (.vol) file format
Allign the virion's raw data longitudinal axis with the volume's z-axis accurately
Coregister over-sampled original volume data to axially symmetric target
Rotational symmetry decomposition
Spectral density segmentation

Programs and code
Space Software

Other structural questions

The homogeneous (coordinate) transformation matrix for 3D points

2008-09-08 I just took a look at the EMD_1220 volume data set and it is very nice. My goal is to fully segment and color render, with cutaways, this volume data.

2008-09-18  The post-processing has morphed into a somewhat separate project: make an average over axial n-fold rotations, where n is point dependent and point-wise derived. I considered this approach to "cleaning" the data (averaging out noise that doesn't have a rotational symmetry) while rendering EMD-1222 [To Do: link]. I've tested the idea (see rotational averages) and now can write the MATLAB code that implements it with precision.

2009-10-08  These notes are reconstructed, fixed up and reformatted from a mangled version that was dumped from the previous site, Currently I'm reinstalling FSL so I can make more progress on this.

Data source

2008-08-03  mcman12 wrote a comment on a short rendered animation loop of mine posted on YouTube, Virus structure, which was rendered from the EMD-1222 data set [To Do link].
"for a higher resolution (better) reconstruction check out EMD-1220"

    This EMD-1220 data set (a cryo-EM assymetric reconstruction, 3-D grayscale) is available (accessed 2008-09-08) at, in the European Bioinformatics Institute's Protein Data Bank Europe (PDBe). Also see the EM Data Bank, which provides three-dimensional density maps for non-viruses.

From the data set web page description:

    The structure of an infectious P22 virion shows the signal for headful DNA packaging

    Gabriel C. Lander, Liang Tang, Sherwood R. Casjens, Eddie B. Gilcrease, Peter Prevelige, Anton Poliakov, Clinton S. Potter, Bridget Carragher, and John E. Johnson:
    The structure of an infectious p22 virion shows the signal for headful DNA packaging J Mol Biol (2006) 312, pp. 1791-1795 [PubMed entry 16709746]


    singleParticle, 17 Å (resolution determined by FSC at 0.5 cut-off


Processing notes

Convert from CCP4 map (.map/.xml) to Space volume (.vol) file format

    MATLAB program convert_map_volume.m does this (see attached). First read the .xml file (as text) to find the volume dimensions (the is 256 x 256 x 256). A hardcoded black and white point must be chosen for each data set conversion.  A histogram of final values is displayed that can help with this.

    I'm not sure I captured the full dynamic range, particularly at the lowest values. Add gamma conversion control? Or redo with 32-bit NIFTI format?

MATLAB command:
>> volData = convert_map_volume( '', 'emd_1220_8-bit', [ 256 256 256 ] );

    Rewrote convert_map_volume.m for NIFTI 32-bit output.

MATLAB command:
>> volData = convert_map_volume( '', 'emd_1220', [ 256 256 256 ] );  [To Do: link]

Allign the virion's raw data longitudinal axis with the volume's z-axis accurately

    The approach I chose was to find a best fit longitudinal axis took advantage of the high radial 5-fold symmetry of the capsid. The capsid was first roughly alligned with the volume's z-axis (within about an angle of 1 px./80 px., just an eyeball resampling at this orientation in Space), and roughly segmented (a binary mask that included almost all the capsid, but not the interior or tail structures). The masked volume (just the capsid) was then averaged with its four roughly equivalent 2*pi/5 rotations. This average volume is axially symmetric, and served as a target for coregistration of the raw data.

Rough (low resolution) segmentation of capsule head and contents (DNA)

    Cubic resampled (x2) and linearly resampled in an approximately axial orientation. Binary mask for substructures (capsid, inner head, outer head, DNA shells 1, 2, 3, and 4) using 2 and 3-D fill on the full binay image. The capsid has the primary orientation information, at least for any 5-fold symmetry structure. All using Space Software.

5-fold rotational averaging

    After masking the high-value 8-bit volume with the rough capsid binary segmentation (using Space), FSL was used to make 5 copies, at n*2*pi/5 rotations.
MATLAB calculation for xt and yt coordinate transformation matrix:

>> xt = xc + ( sqrt( xc*xc + yc*yc ) )*sin( (2*pi)/5 - atan( xc/yc ) )
>> yt = yc - ( sqrt( xc*xc + yc*yc ) )*cos( (2*pi)/5 - atan( xc/yc ) )   ,
where xc = 165 and yc = 163 are the center coordinates with respect to the corner origin.

The full transformation matrix, in homogenous coordinates:

0.309017  -0.951057   0    269.0934
0.951057   0.309017   0    -44.2941
0          0          1      0
0          0          0      1

FLIRT | Utils | Apply FLIRT transforms (ApplyXFM)
flirt -applyxfm -init -out -ref

FSL merge of the five capsid rotations to 4-D volume and collapse across the "time" (n=5) dimension:

fslmerge -t capsid_avg.nii.gz emd_1220_rot_0 emd_1220_rot_72 emd_1220_rot_144 emd_1220_rot_216 emd_1220_rot_288

fslmaths capsid_avg -Tmean capsid_mean

Note that this took a long time, apparently due to memory limitations. I think this time could be reduced using fslmaths -add, but would require a conversion to a floating point format, and adding each one sequentially. Adding in the 8-bit format caused value wrap-around (8-bit overflow).

I think that adding sequentially will be required for the the alligned image rotation averages, to maintain the bit depth without memory limitations.

                    [To Do: Screenshot of capsid and capsid_mean section.]

An assymetric capsid slice (left, high contrast image of an angled slice, resampled using linear interpolation) compared with the coregisitration target (right), a 5-fold mean image, on orthogonal slice planes:
Capsid coresitration target comparison
Comparison of the nearly alligned roughly segmented capsid (left) and the 5-fold rotational average (right). At top is an axial slice through the vertex at top of capsid (magnified x2 relative to bottom images), where the mis-allignment is most pronounced.

Coregister over-sampled original volume data to axially symmetric target

    Raw data was cropped (to limit memory footprint) and over-sampled (tri-cubic) by a factor of x2, in Space Software.

    The raw 256x256x256x4 byte = 64 MB volume is too large to resample x23. So need to crop (in Space after Volume | Volume Info | Center, exact), both this raw data and the capsid weighting volume, emd_1220_8-bit_high-125-230_x2_capsid_weight.nii.gz.

    The raw data crop at native resolution:

x: [ 41 : 212 ]  (even width s.t. the resampled dimensions will be odd)
y: [ 35 : 218 ] 
z: [ 39 : 256 ]

    A binary coregistration weighting volume was created (emd_1220_8-bit_high-125-230_x2_capsid_weight_cropped.nii) by segmenting the capsid at the original orientation, in a procedure identical to that above (rough segmentation of capsule head).

    This was using a high values (125-230) 8-bit volume, but at the original orientation. This tilted orientation made it more difficult to segment the tail/capsid boundary, but it doesn't need to be perfect -- even a majority portion of the capsid would coregister well.

    Weighting origin tweaked to exactly overlap raw data:

Capsid binary wieghting on EMD-1220 raw data at x2 resampling
Note the value wrap flaw on the high luminance "needle" center at bottom-left. This might be due to a flaw in the Space tri-cubic resampling algorithm, or in the MATLAB conversion code. [ 2008-09-15 I can minimize the effect of this, but it seems to be a problem with Space's storage or interpretation of NIFTI range information.]

FSL coregistration:

Target: capsid_mean.nii
Object: emd_1220_cropped_x2.nii
Weighting of object: emd_1220_8-bit_high-125-230_x2_capsid_weight_cropped.nii
Coregistered output: emd_1220_cropped_x2_alligned.nii
Normalized correlation method. Straight correlation method doesn't work, because the raw values (32-bit) are a very different range than the mean capsid (8-bit).
Tri-linear resampling.

flirt -in /emd_1220_cropped_x2.nii.gz -ref capsid_mean.nii.gz -out /emd_1220_cropped_x2_alligned.nii.gz -omat /emd_1220_cropped_x2_alligned.mat -bins 256 -cost normcorr -searchrx 25 35 -searchry -5 5 -searchrz -5 5 -dof 9  -inweight /emd_1220_8-bit_high-125-230_x2_capsid_weight_cropped.nii.gz -interp trilinear

    The final alligned volume, emd_1220_cropped_x2_alligned.nii, is 32-bit, 333 x 339 x 435 vx. (after post-cropping), origin at (164, 169, 217), 188 MB:

EMD-1220 orthogonal slices of alligned volume
EMD-1220, orthogonal slices of alligned volume, high contrast.
Averaged 5-fold rotation of capsid/DNA and 6-fold rotation of tail structures.

Rotational averages

    Write a specialized MATLAB program for radial symmetries analysis and n-fold rotational averaging. Started rotation_average_volume.m.

    Here's a proof of concept image using a test version of rotation_average_image.m. The concentric DNA rings near the axial capsid center are nearly circular. By rotational averaging the central radial regularities are highlighted:

DNA nested ring rings compared with 2-fold symmetry average
The original image (left) and an average with its pi/2 radian rotation about the center (right). The ring definition near the center is enhanced. Note that the capsid (pink at left) does not have a 2-fold symmetry on this plane, so its rotation has only a small overlap (pink at right).

DNA rings averaged over 5- and 90-fold symmetries
DNA rings averaged over 5-fold (left) and 90-fold (right) rotations. The 5-fold rotation  matches the capsid symmetry. Nine DNA rings are clear, but the central features are likely artifacts (noise, or different rotational symmetry).

    Tentative procedure to find center accurately and average over 2*pi/5 and 2*pi/6 rotations. I've eperimented with this procedure on 2-D slices, but only minor modification need to be made to do 3-D serial slices:

    Did this process on emd_1220_cropped_x2_alligned.nii. Note that it is not needed for the rotational averaging algorithms (below). But for finding the center in single slices, it worked well to remove the contribution from noise, only taking into account the standard deviation of histogram values that were not on the peripheral noise.

>> rotation_average_image( 'test_remove_background-zeroback.png', [169 169], 5, true, 'test_' );

    To Do: Ignore zero values (background), and require a minimum number of low deviation points. Ignore background point in net deviation calculation. Set map points to zero if not enough forground points on rotation.

    To Do: Accept an n-fold rotation map and only calculate deviation of rotation points for each corresponding n-fold rotation.

>>rotation_nfold_image( 'test_remove_background-zeroback.png', [169.2383  170.0625], [5 6 10 12 15 18 24 25 30], true, '' );

    To Do: Jitter the center point for each point rotation. Save the maximal point for the n-fold rotation that was most stable with jitter. Also save the jitter maps.

    Experimented with this a lot. This point by point method (evaluating best n-fold fit at each point) is not the way to go -- to noisy and not easy to integrate neighborhood and global n-fold information. See more simple approach below.

    Repeat center finding, using the new n-fold rotation map. Different centers for different n?

    Rough 2 part segmentation of alligned volume by symmetry (5-fold for capsid and DNA, 6 or 12-fold for tail). Include some overlap at intersections. Find best fit n-fold rotation and center for all of one segmented slice, checking for other possible symmetries (circular symmetry, nearby odd symmetries). Repeat on all slices. Trim and blend two (or more) symmetry averaged volumes.

Averaged rotation test, 5-fold and 6-fold composite
5-fold rotational averaging of capsid (white on right) and 6-fold averaging of tail (pink on right). Some (unknown) parts of the 12-fold symmetric "bumps" at the tail periphery are artifacts of the rotational averaging. Ameliorate by excluding oulliers on annular rings before averaging.

    Used rotational center finding and averaging on several planes of the 16-bit alligned data with exterior noise zeroed (see Exterior mask above). It is beautiful, both 5 and 6-fold components on the critical boundary planes. The calculated centers were dead on for all capsid only planes, and only a possible +-.25 pixel shift of tail only planes, with little or no axis skewing.
    Currently used this program, with only one plane of interest and one symmetry hardcoded, to center find and rotationally average:

>> rotation_average_volume( 'emd_1220_cropped_x2_alligned_background_masked.nii', [170 165] );

[To Do: Example average images, derived center estimates.]

Rotational symmetry decomposition


    Testing radial_components_image.m

5 and 6-fold decomposition of EMD-1220 capsid tail boundary
5-fold symmetric capsid structures (blue) overlaid with 12-fold injectosome structure. This axial plane is near where these two symmetries come in contact, or overlap.


>> fPwr = radial_components_volume( 'emd_1220_cropped_x2_alligned.nii', [170  165] );
--> calls:     fpwower = radial_components_image( imPlane, xyCenter, false, num2str( iP ) );

% fPwr = radial_components_volume( strInputFilename, xyCenter )
% ------------------------------------------------------------------------------
%   Fourier decomposition of annuli about the center of a selected serial
%   grayscale z-planes of a NIFTI volume.

% Hardcoded information:

bShow  = false;   % [true, false] If true, displays each raw slice.
bClear = false;   % [true, false] Unload volume from memory while analyzing each plane.

iPlanes = 1:435    % [ n : m ] z-coordinates of planes to be analyzed.
rRes = 2;


% fPwr = radial_components_image( strInputFilename, xyCenter, bShow, strOutputFilePrefix )
% ------------------------------------------------------------------------------
%   Fourier decomposition of annuli about the center of a grayscale image.
%   This version is specialized for analysis of the 5 and 6-fold symmetry
%   of the P22 virion (particularly the EMD-1220 data set), but it could be
%   modified for other annular symmetries.

% Hardcoded information:

bWrite = true;
bDoReconstruct = true;
rRes = 2;
% To Do: Calculate rmax from the given center and image size, use this
% parameter as a second hard maximum.
rmax = 163;

% Annular sampling rate; the number of sample points on each annulus.
% Lower integer factors have some computational speed advantage.

% This is longer than needed, and takes up considerable memory.
np = 2*2*3*3*5*7*11     % = 13860, excludes 8 as a factor


%     fPwr           - total Fourier power of annuli at each radius.
%                       1  - 1 cycle/annulus
%                       2  - 2*n, excluding mod 5 and mod 6
%                       3  - 3*n, excluding mod 5 and mod 6
%                       4  - 4*n, excluding mod 5 and mod 6
%                       5  - 5*n cycles/annulus, excluding 30, 60, ...
%                       6  - 6*n cycles/annulus, excluding 30, 60, ...
%                       7  - 7*n, excluding mod 5 and mod 6
%                       8  - all nonzero and not in 5 or 6-fold
%                       7  - [not used, empty]
%                       10 - zeroth (DC, or constant) component
%                       11 - 11*n, excluding mod 5 and mod 6
%                       12 - 12
%                       13 - 13*n, excluding mod 5 and mod 6
%                       14 - 2
%                       15 - 3
%                       16 - 4
%                       17 - 6
%                       18 - 8
%                       19 - 30*n
%                       20 - 60*n

>> serial_mat_to_nii( 'test_', 'Test_out.nii' )

% serial_mat_to_nii( strInputFilestem, strOuputFilename )
% ------------------------------------------------------------------------------
%   Read a series of 2-D matrices from .mat format into a 3-D volume and
%   save as a NIFTI volume.
%   Currently it is used for assembling the planes output by
%   radial_components_image.m.
%   This version is specific to a set of .mat files that
%   contained an image in a structure field .im56fold.

% Hardcoded information:

bShow   = false;   % [true, false] If true, displays each raw slice.
iPlanes = 1:435  % [ n : m ] z-coordinates of planes to be analyzed.


Before and after rotational symmetry filtering
Before (left) and after (right) rotational symmetry filtering. High contrast to emphasize DNA packing structures. The 5 and 6-fold filtering parameters are just a rough guess here. There are some flaws, for example the blobs at far right are artifacts because I neglected to include overtones of 5-fold symmetries that are also overtones of 6-fold symetries, even though there are negligable 6-fold components at that location.

    Spectral powers above a threshold, for several sets of frequencies using  spectral_segmentation_EMD_1220.m:

MATLAB commands:

>> load  emd-1220_RPwr
>> spectral_segmentation_EMD_1220( 'Cross_section_half.png', fPwr )

            where emd-1220_RPwr.mat was generated by radial_components_volume.m (see above), containing the fPwr 3-D array.

% Hardcoded information:

bShow   = false;   % [true, false] If true, displays each raw slice.
bWrite  = true;    % [true, false] If true, write segmentation image.
strOutFilestem = 'test_seg_fundamental_6'

iChannels = [ 21:33 ]
%   iChannels = [ 32 ]
%ndev = 5.0;       % Segmentation threshold, n standard deviations from median cut.
ndev = 6.0;        % Segmentation threshold, n standard deviations from median cut.

Annular spectral segmentation of P22 volume
Reperesentative cross section (left) compared with spectral segmentation (right).

Red: f = 5 or 10 cycles/annulus
Green and cyan: f = 6 cycles/annulus
Blue: f = 12 cycles/annulus
Magenta: f = 5 and 12 cycles/annulus 
Yellow: f = 6 and 12 cycles/annulus 
White: Other frequencies, primarily 7 and 13 cycles/annulus, and overlap of 5, 6 and 12  cycles/annulus.
The white in the bottom right corner is an artifact due to constant (zero) filling where there was no image data.

Spectral density segmentation

    For second itteration of spectral density segmentation:
    Done and tested. Need to add logic for overlapping regions -- what ratio of overtones to use.
    Made imPwrSigMap. Find clusters with at least one value == 1 and remove all other values (set to zero).
    It is apparent that a large fraction of the non-5-6-12-fold power spikes (1-7-13-fold, larger white blobs in the spectral segmentation map above) are due to mis-estimation of the center(s) of rotation. Remember that the original (single) estimate was based on the capsid segmentation, all 5-fold.
    Re-estimate the centers for the 5-fold, 6-12-fold and central (about circularly symmetric) regions independently. Use hill climbing and spectral decomposition. For which center is there the least residual power?

    Currently testing  rotation_center_volume/image.m for center estimation by n-fold spectral power peak , about the longitudinal axis.

rotation_center_image.m calls radial_components_image.m

with the hardcoded parameter:
nCriterion = 2

To Do

    Tested a 90 degree (volume corner) wedge, it works well. In Space Software, simpy create an all white volume and adjust its x/y origin (Volume|Volume Info) until its corner is near the origin of the virion volume object. Apply Overlays|Mask with this overlay and Overlays|Merge All.

Programs and code

    Space Software

    My free Windows GUI for display, navigation, rendering, editing, and measurement of 3-D data sets.

    Space was used for some of the processing (resampling, value scaling, etc), most of the image slices, and all of the rendering.

    FMRIB Software Library (FSL)

    Free linux based analysis tools intended for for FMRI, MRI and DTI brain imaging data, but the coregistration and volume math tools are suitable for use with other volume data. GUI and command line tools. 

    Used for coregistration, application of some rotation transforms, and volume averaging.

    To avoid writing code that allows finding the axis of rotation accurately (see above, iterated least squares self-coregistration fit about capsule's 5-fold symmetric axis, and head's six-fold rotational axis), I decided to try a well tested coregistration program. The FLIRT tool in the Linux FSL toolbox does volume coregistration with a variety of options, and also has functions that can resample, transform and average volumes.

    As I'm using a Windows operating system, I needed to install VMware Player to run the FSL virtual machine. This was a huge hassle, taking a lot of time to get everything working about right.  [To Do: link to FSL instructions, outline installation.]

    Space Software still has some compatibility issues (and a bug or two) with respect to the NIFTI volume format. In particular, the first value of the dim[8] header item gets written as 1 (short integer), while it should be 3 (the number of dimensions of the data set. The work-around is to manually set this value in the Space Volume | Volume Information | View/Edit Header dialog box. Also, the header length appears to by off by 4 bytes, but this didn't crash FSLView (v3.0) or FLIRT

    4x4 transforms can by applied to volumes using the FLIRT utility (flirt -applyxfm, -init and -out). This transform program can be accessed through the FLIRT GUI using the Util | Apply FLIRT Transform  (ApplyXFM) button (lower right).

    MATLAB code

    MATLAB is an expensive high-level language and interactive environment that enables fast coding with matrix based data. I wrote these MATLAB programs for several specific purposes.


    Reads 32-bit floating point volume data from a CCP4 format file, and writes the 3-D array into an 8-bit unsigned Space Volume Type 1 or a 32-bit single file format NIFTI(.nii) file. Both formats are readable by Space Software, and FSL and other programs read the NIFTI format.

    [To Do: Add argument to convert_map_volume.m for desired output volume format. Document and include the NIFTI MATLAB tools I used: make_nii.m and save_nii.m

"Part of this file is copied and modified under GNU license from MRI_TOOLBOX developed by CNSP in Flinders University, Australia NIFTI data format can be found on:
                - Jimmy Shen ( "

    Rewrote convert_map_volume.m to include NIFTI 32-bit output option.

[To Do: Document and attach this new code.]


    Average n-fold rotations of an image about a near-central point. The location of the center is refined using a hill climbing algorithm that maximizes the standard deviation of values in the average image.

[To Do: Example usage.]




    See MATLAB code file, remove_background_peak_image.m.

    [ vPeak fFWHM_range imZeroBackground ] = remove_background_peak_image( strInputFilename, fFWHM, bShow )

    USAGE:  [ vPeak fFWHM_range imZeroBackground ] = remove_background_peak_image( 'test_remove_background.png', 1.5, true );

        For an image with a relatively uniform noisy background (resulting in a luminance histogram peak), find the peak luminance histogram hump and the range of values above  a given fraction of the full-width at half-maximum FWHM. For example, with fFWHM = 1.5 the width of the range of values that is set to zero is 1.5 times the derived FWHM.
        For an image with a relatively uniform noisy background (composing a luminance histogram peak), find the peak luminance histogram hump and the range of values above  a given fraction of the full-width at half-maximum FWHM. For example, with fFWHM = 1.5 the width of the range of values that is set to zero is 1.5 times the derived FWHM.

    I'm not going to use this for this project, although this is why wrote it. I think it will be easier to mask the backgound on the whole volume using the Space Software 3-D Fill tool. This avoids dealing with interior pixels/voxels that this program removes from the interior.
Example usage:

    This example uses a hardcoded parameter used to remove the uniform image area (at bottom of original image below) before finding the histogram peak:

     valueIgnore = 125; % This values will be set to zero before processing and in output.

MATLAB command and command line output:
>> remove_background_peak_image( 'test_remove_background.png', 1.5, true );

    szInput =

       339   333

    Peak value:        128.5
    Low range value:   118.3749
    High range value:  137.8805

    Original and final images, and intermediate (image and histogram) outputs:

original (input) image

 original (input) image

 final (output) image

 final (output) image


histogram of original after ignored value removed
histogram of original after ignored value removed

 histogram of final image
histogram of final image

 image after ignore value set to zero
image after ignore value set to zero



Other structural questions

    Central and below axial center rings don't fit the rest of the DNA spool pattern. What are they and how do they just float, apparently disconnected?

    Some of the peripheral "stacked rings" of DNA are very robustly detected, indicative of "spooling". But the non-uniform breaks and bridges between the rings also appear to be real. How could this be true for a rotationally averaged image? There doesn't seem to much rhyme or reason - not correlated across concentric rings or capsid structures.

    Central axial structures appear to be real but not "lead DNA". There's a hyper-density plug, and a hypo-density cap.

    What do small Fourier components that are not multiples of 5 or 6 mean? Could there be admixtures of, say, 11-fold tail structures, and would these show up as 11-fold components? Would their phase be coherent if there were an admixture? Can non-harmonic distortions be distinguished from real n-fold regularities?

    The 12-fold portion of the tail is fit into a 5-fold vertex of the capsid, like a square peg in a round hole. The Fourier decomposition shows the boundary has a narrow region where both annular frequencies are significant. If the 5-fold and 6-fold symmetries are independently reconstructed, are the overlapping features due to distortion of one or both structures, or are the "electron densities" (or whatever is measured by EM) of both structures are interpenetrating? Can these two mechanisms be distinguished in some ideal case?

    I got distracted fom this project in making an animation that I naimed "Mixed symmetries", because I combined some 3-D periodic, 3-D perceptual rotation, and 2-D periodic symmetries. After I posted it I googled "mixed symmetries", yielding some physics energy level papers. But the top hit on the "images" search with that phrase was a rendering of a Caudovirus, Epsilon 15, a virus that infects the bacterium Salmonella based on Cryo-EM reconstruction! There's a good 2006 paper [To Do: Refind this link in "Virus_structure_and_mechanics_information"] with very nice illustrations that explores just this question. Quite a coincidence as I had not looked into this close relative of P22, and it has obvious distortions of molecular configuration at the 6/5-fold boundary. I though that there would be a long list of topics that used the phrase, but the top two happen to be obscure side topics I've worked on.

    Is the axis of the capsid alligned with that of the tail?
    The allignment was based on the capsid symmetry. I expect the capsid center derived from each axial slice will be close to parallel with the volume's z-axis. But the tail structures might be skewed due to the 5 to 6/12-fold interlocking.

    Is there a distortion of the icosohedral capsid tail vertex due to the tail?
    After the full n-fold segmentation is done, a difference with one (or an average of several) vertex will demonstrate any distortion.

    Rotional symmetry program, for EMD-1220:

        I've looked closely at how to symmetrically average planes normal to the axis, particularly near the boundary of 5 and 6-fold symmetry. This tail-capsid boundary has components of both symmetries, and the image can be partitioned with a simple linear model.
        The tail-capsid interface is particularly interesting with respect to molecular machines. [2008-10-15 1950's molecular machine book by ???von Hippel] The dsDNA flows through the axis of the tail in both directions, while being packaged and while being inserted in a new host cell. The mechanism of both processes is not known. It is likely that the central tail structure is the basis of a molecular motor. It has been suggested that its incommensurate 5 and 6-fold boundary at the capsid vertex might form a kind of ratchet.

        "Viral DNA packaging motors inject viral genomic DNA into capsids as part of their replication cycle, packing it very tightly. [5]"

The homogeneous (coordinate) transformation matrix for 3D points

    FSL, and many other 3D programs use a 4x4 spatial transfomation matrix of the following form.  This 4x4 matrix operator performs the rotation given by
$ R(\alpha,\beta,\gamma)$
followed by a translation given by
$ x_t,y_t,z_t$

The result is
homogenous 3D transform matrix
The order of operations is critical. The matrix T represents the following sequence of transformations:

1. Roll by gamma
2. Pitch by beta
3. Yaw by alpha
4. Translate by $ x_t,y_t,z_t$


 Public Domain DedicationThis work is dedicated to the Public Domain.
There are no restrictions on use of the images or code on this page. Claiming to be the originator of the material, explicitly or implicitly, is bad
karma. A link (if appropriate), a note to dow[at], and credit are appreciated but not required.

Comments are welcome (dow[at]