1 """@namespace IMP.sampcon.precision_rmsd
2 Calculation of precision and RMSD."""
4 from __future__
import print_function
9 import pyRMSD.RMSDCalculator
13 def parse_custom_ranges(ranges_file):
16 with open(ranges_file)
as fh:
19 return d[
'density_custom_ranges']
22 def get_particles_from_superposed(
23 cluster_conform_i, cluster_conform_0, align, ps, trans,
25 def _to_vector3ds(numpy_array):
31 calculator_name =
"QCP_SERIAL_CALCULATOR"
33 calculator_name =
"NOSUP_SERIAL_CALCULATOR"
35 conforms = numpy.array([cluster_conform_0, cluster_conform_i])
37 if symm_groups
is None:
38 calculator = pyRMSD.RMSDCalculator.RMSDCalculator(
42 if calculator_name ==
'NOSUP_SERIAL_CALCULATOR':
44 s1 = parse_symm_groups_for_pyrmsd(symm_groups)
45 calculator = pyRMSD.RMSDCalculator.RMSDCalculator(
47 fittingCoordsets=conforms,
48 calculationCoordsets=conforms,
49 calcSymmetryGroups=s1,
52 calculator = pyRMSD.RMSDCalculator.RMSDCalculator(
54 fittingCoordsets=conforms,
55 calcSymmetryGroups=[],
56 fitSymmetryGroups=symm_groups)
58 check1 = (calculator_name ==
'NOSUP_SERIAL_CALCULATOR')
59 check2 = symm_groups
is not None
63 rmsd, superposed_fit, _calc_fit = calculator.pairwise(
64 0, 1, get_superposed_coordinates=
True)
66 rmsd, superposed_fit = calculator.pairwise(
67 0, 1, get_superposed_coordinates=
True)
74 _to_vector3ds(superposed_fit[0]), _to_vector3ds(cluster_conform_0))
76 for particle_index
in range(len(superposed_fit[1])):
81 return rmsd, ps, trans
85 """Compute mean density maps from structures.
86 Keeps a dictionary of density maps,
87 keys are in the custom ranges. When you call add_subunits_density, it adds
88 particle coordinates to the existing density maps.
91 def __init__(self, custom_ranges=None, resolution=20.0, voxel=5.0,
94 @param list of particles decorated with mass, radius, and XYZ
95 @param resolution The MRC resolution of the output map
97 @param voxel The voxel size for the output map (lower is slower)
100 self.MRCresolution = resolution
102 self.count_models = 0.0
104 self.bead_names = bead_names
105 self.custom_ranges = custom_ranges
109 self.particle_indices_in_custom_ranges = {}
111 for density_name
in self.custom_ranges:
112 self.particle_indices_in_custom_ranges[density_name] = []
115 for index, beadname
in enumerate(self.bead_names):
116 for density_name
in self.custom_ranges:
118 for domain
in self.custom_ranges[density_name]:
119 if self._is_contained(beadname, domain):
120 self.particle_indices_in_custom_ranges[
121 density_name].append(index)
124 def normalize_density(self):
127 def _create_density_from_particles(self, ps, name,
128 kernel_type=
'GAUSSIAN'):
129 '''Internal function for adding to densities.
130 pass XYZR particles with mass and create a density from them.
131 kernel type options are GAUSSIAN, BINARIZED_SPHERE, and SPHERE.'''
134 dmap.set_was_used(
True)
136 if name
not in self.densities:
137 self.densities[name] = dmap
143 dmap3.set_was_used(
True)
145 dmap3.add(self.densities[name])
146 self.densities[name] = dmap3
148 def _is_contained(self, bead_name, domain):
149 """ domain can be the name of a single protein or a tuple
150 (start_residue,end_residue,protein_name)
151 bead is a string of type moleculeName_startResidue_endResidue
154 (bead_protein, bead_res_start,
155 bead_res_end, bead_copy) = bead_name.rsplit(
"_", 3)
158 if isinstance(domain, tuple):
159 domain_protein = domain[2]
161 domain_protein = domain
163 if "." in domain_protein:
164 spl = domain_protein.split(
".")
165 domain_protein = spl[0]
166 domain_copy = int(spl[1])
168 domain_copy = bead_copy = -1
170 if bead_protein != domain_protein
or int(bead_copy) != domain_copy:
174 if isinstance(domain, tuple):
175 bead_residues = set(range(int(bead_res_start),
176 int(bead_res_end)+1))
177 domain_residues = set(range(int(domain[0]),
179 return not domain_residues.isdisjoint(bead_residues)
184 """Add a frame to the densities.
185 @param ps List of particles decorated with XYZR and Mass.
187 self.count_models += 1.0
189 particles_custom_ranges = {}
190 for density_name
in self.custom_ranges:
191 particles_custom_ranges[density_name] = []
194 for density_name
in self.custom_ranges:
196 in self.particle_indices_in_custom_ranges[density_name]:
197 particles_custom_ranges[density_name].append(
201 for density_name
in self.custom_ranges:
202 self._create_density_from_particles(
203 particles_custom_ranges[density_name], density_name)
205 def get_density_keys(self):
206 return list(self.densities.keys())
209 """Get the current density for some component name"""
210 if name
not in self.densities:
213 return self.densities[name]
215 def write_mrc(self, path=".", file_prefix=""):
217 for density_name
in self.densities:
218 mrc = os.path.join(path, file_prefix +
"_" + density_name +
".mrc")
219 self.densities[density_name].
multiply(1. / self.count_models)
220 densmap[density_name] = mrc
222 self.densities[density_name], mrc,
def get_density
Get the current density for some component name.
Compute mean density maps from structures.
Class for sampling a density map from particles.
DensityMap * multiply(const DensityMap *m1, const DensityMap *m2)
Return a density map for which voxel i contains the result of m1[i]*m2[i].
DensityMap * create_density_map(const IMP::algebra::GridD< 3, S, V, E > &arg)
Create a density map from an arbitrary IMP::algebra::GridD.
def add_subunits_density
Add a frame to the densities.
Basic utilities for handling cryo-electron microscopy 3D density maps.
A decorator for a particle with x,y,z coordinates.
algebra::BoundingBoxD< 3 > get_bounding_box(const DensityMap *m)
Utilities to help with RMSD calculation.
Transformation3D get_transformation_aligning_first_to_second(Vector3Ds a, Vector3Ds b)