Density Map Instance Informations

class EMMap.Map(fullMap, origin, apix, filename, header=[])[source]

A class representing all information from a density map file. NOTE: Currently it can only read the CCP4/MRC format.

box_size()[source]
Return:
size of the map array, in ZYX format.
centre()[source]

Centre of the Map Instance

Return:
Vector instance of the centre of the map in Angstroms.
change_origin(x_origin, y_origin, z_origin)[source]

Change the origin of the map to a new origin.

Arguments:

x_origin, y_origin, z_origin
new co-ordinates of origin.
Return:
new Map instance
copy()[source]
Return:
copy of the Map.
fourier_transform()[source]

Apply Fourier transform on the density map.

Return:
new Map instance
getMap()[source]
Return:
array containing the map density data.
get_com()[source]

Retrieve the centre of mass of the map.

Return:
Vector instance
get_line_between_points(point1, point2)[source]

Return an array of float values representing a line of density values between two points on the map.

Arguments:

point1, point2
Vector instances of the end points co-ordinates of the line.
Return:
array of floating values
get_normal_vector(x_pos, y_pos, z_pos)[source]

Calculate the normal vector at the point specified. Point calculated using 3SOM algorithm used by Ceulemans H. & Russell R.B. (2004).

Arguments:
x_pos, y_pos, z_pos
pixel in map on which to calculate normal vector.
Returns:
Normal vector at the point specified
get_point_map(min_thr, percentage=0.2)[source]

Calculates the amount of point to use for the NV and CD score.

Arguments:

min_thr
run get_primary_boundary on target map.
percentage
percentage of the protein volume.

Return:

number points
get_pos(minDens, maxDens)[source]

Identify a set of voxels in the map whose density values fall between the specified minimum and maximum values.

Arguments:
minDens
minimum density value to include in array.
maxDens
maximum density value to include in array.
Return:
An array of 3-tuples (indices of the voxels in x,y,z format)
get_primary_boundary(molWeight, low, high, vol_factor=1.21)[source]

Calculates primary boundary density value. Volume of pixels with greater density than output is equivalent to volume given by molecular weight of protein. Uses recursive algorithm. NOTE: when used to calculated the NV score this value should be calculated for the experimental map.

Arguments:

molWeight
molecular weight of protein; use get_prot_mass_from_atoms() if your structure contains HETATOMS else use get_prot_mass_from_res().
low, high
minimum and maximum values between which the boundary will be taken. Initial values should be given by minimum and maximum density values in map.
vol_factor
in cubic Angstroms per Dalton. This is the approximate value for globular proteins used in Chimera (Petterson et al, 2004) from Harpaz 1994. Other recommended volume factor are 1.5 (1.1-1.9) cubic Angstroms per Dalton in EMAN Volume/mass conversions assume a density of 1.35 g/ml (0.81 Da/A3) (~1.23A3/Da)

Return:

primary boundary density value (float)
get_second_boundary(primary_boundary, noOfPoints, low, high, err_percent=1)[source]

Calculate the second bound density value. For a given boundary value, it calculates the second bound density value such that a specified number of points whose density values fall between the defined boundaries Uses recursive algorithm.

Arguments:
primary_boundary
primary threshold, normally given by get_primary_boundary method based on protein molecular weight.
noOfPoints
Number of points to use in the normal vector score - try first with 10% (maybe 5%better) of the number of points in the map ( round((self.map_size())*0.1)
low, high
minimum and maximum values between which the threshold will be taken. low should be equal to the value returned by the get_primary_boundary() method and high is the maximum density values in map.
err_percent
default value of 1. Allowed to find a secondary boundary that includes a 1% error.
Return:
secondary boundary density value
get_significant_points()[source]

Retrieve all points with a density greater than one standard deviation above the mean.

Return:
An array of 4-tuple (indices of the voxels in x,y,z format and density value)
get_vectors()[source]

Retrieve all non-zero density points in the form of Vector instances.

Return:
An array of 4-tuple (indices of the voxels in x,y,z format and density value)
header_to_binary()[source]

Returns binary version of map header data. For use in writing out density maps in MRC file format.

laplace_filtered()[source]

Apply Laplacian filter on density maps

Return:
new Map instance
makeKDTree(minDens, maxDens)[source]

Returns k-dimensional tree of points in the map with values between those chosen for quick nearest-neighbor lookup.

Arguments:
minDens
minimum density value to include in k-dimensional tree.
maxDens
maximum density value to include in k-dimensional tree.
Return:
index into a set of k-dimensional points.
make_bin_map(cutoff)[source]

Return a new Map instance that has been binarised. All voxel with densities above and below the specified cutoff value are assigned a value of 1 and 0 respectively.

Arguments:

cutoff
cutoff density value
Return:
new binarised Map instance
map_rotate_by_axis_angle(x, y, z, angle, CoM, rad=False)[source]

Return a new Map instance rotated around its centre.

Arguments:
angle
angle (in radians if rad == True, else in degrees) to rotate map.
x,y,z
axis to rotate about, ie. x,y,z = 0,0,1 rotates the Map round the xy-plane.
CoM
centre of mass around which map will be rotated.
Return:
Rotated new Map instance
map_size()[source]
Return:
size of the array fullMap.
matrix_transform(mat, x=0, y=0, z=0)[source]

Apply affine transform to the map.

Arguments:

mat
affine 3x3 transformation matrix
shape
new box dimensions
x, y, z
translation in angstroms.
Return:
new Map instance
max()[source]
Return:
maximum density value of the map.
mean()[source]
Return:
mean density value of map.
median()[source]
Return:
median density value of map.
min()[source]
Return:
minimum density value of the map.
normalise()[source]

Return a new Map instance with normalised density values.

Return:
new Map instance
origin_change_maps(MapRef)[source]

Return a new Map instance with origin changed accordingly to Reference Map

Arguments:
MapRef
Reference Map
Return:
new Map instance
pad_map(nx, ny, nz)[source]

Pad a map (in place) with specified increments along each dimension. Arguments:

nx,ny,nz
Number of slices to pad in each dimension.
Return:
new Map instance
pixel_centre()[source]
Return:
Vector instance of the centre of the map in pixels.
represent_normal_vectors(min_threshold, max_threshold)[source]

Create a Structure instance representing normal vectors of density points specified.

Arguments:
min_threshold, max_threshold
minimum/maximum values to include in normal vector representation.
Return:
Structure Instance
resample_by_apix(new_apix)[source]

Resample the map based on new_apix sampling.

Arguments:
new_apix
Angstrom per pixel sampling
Return:
new Map instance
resample_by_box_size(new_box_size)[source]

Resample the map based on new box size.

Arguments
new_box_size
An array containing box dimension in ZYX format
Return:
new Map instance
resize_map(new_size, centre=False)[source]

Resize Map instance.

Arguments:

new_size
3-tuple (x,y,z) giving the box size.
centre
default False
Return:
new Map instance
rotate_by_axis_angle(x, y, z, angle, CoM, rad=False)[source]

Rotate the map around its centre given an axis and angle.

Arguments:

angle
angle (in radians if rad == True, else in degrees) to rotate map.
x,y,z
axis to rotate about, ie. x,y,z = 0,0,1 rotates the Map round the xy-plane.
CoM
centre of mass around which map will be rotated.
Return:
new Map instance
rotate_by_euler(x, y, z, CoM, rad=False)[source]

Rotated map around pivot given by CoM using Euler angles.

Arguments:

x,y,z
Euler angles (in radians if rad == True, else in degrees) used to rotate map.
CoM
centre of mass around which map will be rotated.
x, y, z
translation in angstroms.
Return:
new Map instance
rotate_by_matrix(mat, CoM, rad=False)[source]

Rotated the map around pivot given by CoM using a rotation matrix

Arguments:
mat
3x3 matrix used to rotate map (in radians if rad == True, else in degrees).
CoM
rotation pivot, usually the centre of mass around which map will be rotated.
Return:
new Map instance
scale_map(scaling)[source]

Scaling Map by scaling factor

Return:
new Map instance
shift_origin(x_shift, y_shift, z_shift)[source]

Shift the Map origin.

Arguments:

x_origin, y_origin, z_origin
new co-ordinates of origin.
Return:
new Map instance
std()[source]
Return:
standard deviation of density values in map.
threshold_map(minDens, maxDens)[source]

Create a Map instance containing only density values between the specified min and max values.

Arguments:

minDens
minimum density threshold
maxDens
maximum density threshold
Return:
new Map instance
translate(x, y, z)[source]

Translate the map by changing origin

Arguments:
x,y,z
translation in angstroms
Return:
new Map instance
update_header()[source]

Update self.header to values currently relevant.

vectorise_point(x, y, z)[source]

Return a tuple of the Angstrom co-ordinates and density value of a particular density point in map. Transform the voxel specified by its indices (x,y,z) into a Vector object. The vector defines the position of the voxel with respect to the origin of the map. The magnitude of the vector is in Angstrom units.

Arguments:
x, y, z
co-ordinates of the density point to be vectorised.
Return:
Vector instance
write_to_MRC_file(mrcfilename)[source]

Write out a MRC file

Arguments:
mrcfilename
name of the output mrc file
x_origin()[source]
Return:
x-coordinate of the origin.
x_size()[source]
Return:
x size of the map array in x direction.
y_origin()[source]
Return:
y-coordinate of the origin.
y_size()[source]
Return:
y size of the map array in y direction.
z_origin()[source]
Return:
z-coordinate of the origin.
z_size()[source]
Return:
z size of the map array in z direction.