PyTom: Align projections and reconstruct a tomogram

Prior to a 3D reconstruction projections must be aligned to a common coordinate system because despite the best efforts to keep the feature of interest in the center during acquisition movement is inevitable. Moreover, the electron optical lense system causes image rotations and even subtle focus-dependent magnification changes, which need to be compensated for.

Currently, PyTom only supports alignment by fiducials (gold beads). The positions
of markers in the projection images need to be provided in an EM-file. This
file can be obtained by interactive marker localization in the `TOM/AV3`

or `EM`

packages. There is also a function in AV3 to convert IMOD
marker files to an EM file (av3_wimp2em.m).

The PyTom alignment function simply minimizes the residual error of the
coordinates projected from a 3D model and the oberserved coordinates in the
projections. In more detail, 3D coordinates of the marker model, the projection
shifts, and the image rotation is first computed using an analytical approach,
in which a constant rotation (=tilt axis) and magnification is assumed
(`TiltAlignment.computeCoarseAlignment`

). This
function determines a tilt axis (to x-axis) between 0 and 180 degrees; by definition
for a tilt axis x there is a second solution, x+180 deg, that would satisfy the
marker equations equally well. You need prior knowledge to decide, which tilt
axis is correct - assuming the wrong tilt axis will yield reconstructions
with the inverted handedness. In the alignment it can be specified if the solution
180<x<360 deg is supposed to be correct rather than the default solution
0<x<180 (handflip=True).

To determine the 3D coordinate system of the markers the coordinates of one
marker point need to be defined. Thus, a reference marker needs to be chosen
(`irefmark`

). To this marker reference coordinates (`r`

)
can be assigned or default values can be chosen. The default value is
`r=(x_iref-center,y_iref-center,0)`

. `x_iref`

and
`y_iref`

are the coordinates of the marker in a reference projection
`iref`

.

After approximate determination of the image rotations, translations, and 3D
coordinates of the markers an optimization algorithm
(`TiltAlignment.alignFromFiducials`

) determines the tilt-angle specific image
rotation, tilt-angle specific magnification, and again translations and 3D coordinates.
For the reference projection `iref`

the magnification is defined as 1.

The projections are aligned by inverse application of the determined translations, rotations, and magnifications to each projection. The transformation is carried out as a single interpolation in real space (default interpolation method: 3rd order spline).

If a tomogram were reconstructed from the unprocessed projections the low-frequency information would be artificially enhanced as illustrated in Fig.1. To obtain a truthful 3D image the information below the so-called Crowther frequency need to be attenuated. This is accomplished by weighting the projections in Fourier space proportional to its frequency perpendicular to the tilt axis. This can be approximated as a 'ramp' function (analytical weighting). A more accurate weighting scheme would be to explicitly compute the overlapping information in Fourier space (exact weighting, not yet implemented in PyTom). In many cases it makes sense to furthermore apply a low-pass filter to the data to eliminate noise, which is dominant in the higher frequencies.

**Fig. 1.**Sampling of Fourier space by projections in 2D. The
Fourier transform of each projection samples a slice of Fourier space
(black lines); its normal is determined by the tilt angle. The signal
‘leak’ of each sample point into the neighboring area is reciprocally
proportional to the diameter D of the object. For reconstruction, the
Fourier coefficients on the 3D grid (intersected lines) need to be
approximated from the sampling points of the projections by an
appropriate algorithm. In principle, the 3D signal can be
reconstructed without any gaps to a resolution k_C if the entire tilt
range of +/-90 deg was accessible. Restricting the tilting to +/-90 deg
gives rise to a ‘missing wedge’.

In PyTom reconstructions are performed using backprojection. The algorithm can use different interpolation kernels. By default we use a 3rd order spline.

`reconstructTomogram.py`

In our tutorial data-set you will find a directory `reconstructTomo`

where you can perform reconstruction of a full tomogram by weighted backprojection.
Projections will be aligned and weighted during this step so that the data you start with are unaligned and unweighted projections. Essentially raw projections only sorted according to the tilt angle.
Running this script will generate a full `tomogram.em`

file, either 2x binned (4x downscaled) (`512,512,128`

) or in original size (`2048,2048,512`

).

The `reconstructTomo.sh`

is essentially the only script you have to run to obtain the aligned, weighted projections and the 2x binned tomogram.

```
```reconstructTomogram.py --tiltSeriesPath ../projections/ \
--tiltSeriesPrefix tomo01_sorted \
--firstIndex 1 \
--lastIndex 41 \
--referenceIndex 21 \
--markerFile ../projections/markfile_temp.mark \
--referenceMarkerIndex 1 \
--projectionTargets ./alignedProjections/tomo01 \
--projectionBinning 4 \
--lowpassFilter 0.5 \
--tomogramFile tomogram.em \
--fileType em \
--tomogramSizeX 512 \
--tomogramSizeY 512 \
--tomogramSizeZ 128 \
--reconstructionCenterX 0 \
--reconstructionCenterY 0 \
--reconstructionCenterZ 0

`reconstructTomogram.py`

, simply run
```
reconstructTomogram.py --help
```

**tiltSeriesPath**: Relative or absolute path to tilt series`EM`

**tiltSeriesPrefix**: Prefix of projections before running index (projections in`EM`

or`MRC`

format)**firstIndex**: Index of first projection**lastIndex**: Index of last projection**referenceIndex**: Index of reference projection used for alignment**markerFile**: Name of markerfile (in`EM`

format, .e.g., generated using TOM)**referenceMarkerIndex**: Index of reference marker to set up coordinate system**handflip**: Is your tilt series outside of 0-180deg**projectionTargets**: Relative or absolute path to the aligned projections that will be generated + file prefix**fineAlignFile**: Relative or absolute path to the file with fineAlign parameters (type should be *.dat)**projectionBinning**: Binning of projections during read**lowpassFilter**: Lowpass filter after binning**tomogramFile**: Relative or absolute path to final tomogram (no tomogram written if not specified)**fileType**: File type (can be`EM`

or`MRC`

- no tomogram written if not specified)**tomogramSizeX**: Size of tomogram in x (no tomogram written if not specified)**tomogramSizeY**: Size of tomogram in y (no tomogram written if not specified)**tomogramSizeZ**: Size of tomogram in z (no tomogram written if not specified)**reconstructionCenterX**: Center where tomogram will be reconstructed (no tomogram written if not specified)**reconstructionCenterY**: Center where tomogram will be reconstructed (no tomogram written if not specified)**reconstructionCenterZ**: Center where tomogram will be reconstructed (no tomogram written if not specified)**weightingType**: type of weighting. -1: r-weighting (default), 0: no weighting.**verbose**: verbose mode

Now PyTom also includes the support of a new reconstruction algorithm called
Iterative Nonuniform fast Fourier transform based Reconstruction method
(*INFR*). For more details about this algorithm, please check:
Iterative reconstruction of cryo-electron tomograms using nonuniform fast
Fourier transforms, Y. Chen et al., JSB 2014.

The way to use it is simple:

`pytom PathToPytom/reconstruction/reconstruct_INFR.py`

**-d**: The directory name, in which all the aligned and unweighted projections are stored. Note that the projections should be EM format and the tilt angle information must be specified in the headers!**-i**: Optional, the number of iterations to run. 10 by default.**-o**: Output filename.

```
``` pytom PathToPytom/bin/reconstructTomogram.py --tiltSeriesPath MyTiltSeriesDir \
--tiltSeriesPrefix NameOfTiltSeriesFiles --firstIndex 1 --lastIndex MyLastIndex \
--referenceIndex MyReferenceProjectionIndex --markerFile MyMarkerFile \
--referenceMarkerIndex MyReferenceMarkerIndex --projectionTargets \
MyFileNameOfAlignedProjections --projectionBinning MyBinning --weightingType 0

Sometimes you might want to do simple things on your reconstructed tomogram like
reducing its size or low-pass filtering it without reconstructing it again.
Here is an example how to reduce the size in *z*:

```
$ipytom
```

from pytom_volume import read, subvolume

v = read('tomogram.em')

#this function will cut out a subvolume the size 512,512,256 out of the original volume starting at position 0,0,128

s = subvolume(v,0,0,128,512,512,256)

s.write('subVolume.em')

Make sure you store the position at which the subvolume was cut. This will coordinate is required for reconstructing subtomograms at full resolution (see Reconstruct Subtomograms Tutorial).

If you want to lowpass filter the tomogram / subvolume to a specified resolution you can do that interactively as well:

```
$ipytom
volume = read('subvolume.em')
```

from pytom.basic.filter import lowpassFilter

filteredVolume = lowpassFilter(volume,100,10)[0] #filter the volume to the 100th
pixel in Fourier space with 10 pixel smoothing

filteredVolume.write('filteredVolume.em')