1 """@namespace IMP.EMageFit.imp_general.representation
2 Utility functions to handle representation.
12 log = logging.getLogger(
"representation")
17 Functions to deal with the representation of assemblies and managing
26 Builds the assembly setting the chains in the PDB file as components
29 hchains = IMP.atom.get_by_type(temp, IMP.atom.CHAIN_TYPE)
31 log.debug(
"Creating assembly from pdb %s,names: %s. Chains %s",
35 for i, h
in enumerate(hchains):
42 """ Read all the PDBs given in the list of names fn_pdbs and adds the
43 hierarchies to the model
46 for i, fn_pdb
in enumerate(fn_pdbs):
51 assembly.add_child(prot)
56 """ Read a PDB molecule, add atoms, and set a name
59 log.debug(
"reading component %s from %s", name, fn_pdb)
61 log.debug(
"reading component from %s", fn_pdb)
66 hierarchy.set_name(name)
72 """ set the children of a molecule type hierarchy as rigid bodies
73 In this case, all the children are the components of the complex.
74 I use the function create_rigid_body(), that creates a lot of
77 I have changed the function and now build the rigid body directly from
78 the leaves of each of the components. With this I guarantee that the
79 number of rigid members is going to be the same if the components have
80 the same number of atoms.
83 raise TypeError(
"create_rigid_bodies(): The argument is not a valid "
87 for c
in molecule.get_children():
98 """ Rename all the chains of an assembly so there are no conflicts with
99 the ids. The names are added sequentially.
102 raise TypeError(
"The argument is not a valid hierarchy")
105 letters = string.ascii_uppercase
106 n_chains = len(all_chains_as_hierarchies)
107 if(len(letters) < n_chains):
108 raise ValueError(
"There are more chains than letter ids")
109 ids = letters[0:n_chains]
110 for h, c_id
in zip(all_chains_as_hierarchies, ids):
113 chain.set_name(
"chain %s" % c_id)
117 """ Gets a hierarchy containing a molecule of DNA and simplifies it,
118 generating a coarse representation of spheres. The function returns
119 a hierarchy with the spheres.
120 n_res - Number of residues to use per sphere.
123 raise TypeError(
"create_simplified_dna: the hierarchy provided is not a "
126 model = dna_hierarchy.get_model()
131 residues = IMP.atom.get_by_type(dna_hierarchy, IMP.atom.RESIDUE_TYPE)
134 for i
in range(0, l, n_res):
136 equivalent_mass = 0.0
137 residues_numbers = []
138 for r
in residues[i: i + n_res]:
140 residues_numbers.append(rr.get_index())
143 for a
in rr.get_children()]
144 xyzrs += residue_xyzrs
146 equivalent_mass += get_residue_mass(r)
151 xyzr.set_radius(s.get_radius())
152 xyzr.set_coordinates(s.get_center())
154 fragment.set_residue_indexes(residues_numbers)
156 simplified_h.add_child(fragment)
157 simplified_h.set_name(
"DNA")
162 def get_residue_mass(residue):
164 raise TypeError(
"The argument is not a residue")
169 mass += ms.get_residue_mass()
174 """ Simplifies an assembly, by creating a hierarchy with one ball per
175 n_res residues. Each of the chains in the new hierarchy are added to
176 the rigid bodies for each of the components.
177 There must be correspondence between the children of the assembly
178 (components) and the rigid bodies. I check for the ids.
181 raise TypeError(
"The argument is not a valid hierarchy")
184 model = assembly.get_model()
185 n_children = molecule.get_number_of_children()
190 for i
in range(n_children):
191 component = molecule.get_child(i)
192 name = component.get_name()
193 rb = components_rbs[i]
195 raise ValueError(
"Rigid body and component do not match")
197 hchains = IMP.atom.get_by_type(component, IMP.atom.CHAIN_TYPE)
215 chain_rb.set_name(
"sub_rb" + name)
216 rb.add_member(chain_rb)
220 coarse_component_h.set_name(name)
221 simplified_hierarchy.add_child(coarse_component_h)
222 return simplified_hierarchy
226 """ Select a component of the assembly using the name """
227 for c
in assembly.get_children():
228 if (c.get_name() == name):
231 "The requested component %s is not in the assembly" %
236 """ Select a rigid body from the rigid_bodies using the name """
237 for rb
in rigid_bodies:
238 if (rb.get_name() == name):
240 raise ValueError(
"This rigid body is not in the set: %s" % name)
244 """ Name to use for the rigid body of a hierarch"""
249 """ Build the rigid body for all the particles in the selection S """
250 ps = S.get_selected_particles()
257 def get_selection_as_hierarchy(model, S):
260 for p
in S.get_selected_particles():
267 """ Gets a selection of particles and decorates them as Atoms.
268 Then all of them are put into a big residue. I have this to use
269 with the multifit.create_coarse_molecule_from_molecule() function
273 for p
in S.get_selected_particles():
279 """ The function returns the particles (fragments) in the coarse hierarchy
280 that were created by summarizing the residues_numbers.
282 Coarse hierarchy - Hierarchy formed by a bunch of
283 fragments. Each fragment must have the residue numbers that it contains
284 residue_numbers - list with the number of the residues that need to
286 The function returns the set of particles that are IMP.atom.Fragments
289 fragments = IMP.atom.get_by_type(coarse_h, IMP.atom.FRAGMENT_TYPE)
292 residues_in_f = ff.get_residue_indexes()
293 for number
in residues_in_f:
294 if number
in residues_numbers:
295 particles.append(ff.get_particle())
302 Rotates the reference frame of a rigid body around the centroid
304 c = rb.get_coordinates()
305 ref = rb.get_reference_frame()
306 R = ref.get_transformation_to().get_rotation()
310 rb.set_reference_frame(ref)
315 Applies a transformation around the centroid of a rigid body.
316 First does the rotation around the centroid and
317 then applies the transformation.
318 @param rb A IMP.core.RigidBody object
319 @param T a IMP.algebra.Transformation3D object
322 rb.set_coordinates(rb.get_coordinates() + T.get_translation())
327 Get the particle for a residue in a hierarchy.
328 @param h The hierarchy
329 @param chain_id If chain_id == False, just search for the residue
330 @param res Number of residue in the chain
337 return s.get_selected_particles()[0]
342 Get the coordinates of a residue (the coordinates of the first particle)
344 @param chain_id See help for get_residue_particle()
345 @param res See help for get_residue_particle()
352 hierarchy2, chain_id2, residue2):
354 Distance between two residues. See the help for get_residue_particle()
369 """ Gets all the chains in a set of hierarchies
370 @param hierarchies A set of IMP.atom.Hierarchy objects
373 for h
in hierarchies:
374 chains_in_h = IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE)
375 for ch
in chains_in_h:
380 def set_reference_frames(rigid_bodies, reference_frames):
381 for ref, rb
in zip(reference_frames, rigid_bodies):
382 rb.set_reference_frame(ref)
387 Returns the atoms in the backbone of the nucleic acid contained
389 backbone 'minimal' returns the atoms: ["P", "O5'", "C5'", "C4'", "C3'", "O3'"]
390 backbone 'trace' returns the atoms C4'
394 if backbone ==
'minimal':
395 backbone_atoms = [
"P",
"O5'",
"C5'",
"C4'",
"C3'",
"O3'"]
396 elif backbone ==
'trace':
397 backbone_atoms = [
"C4'"]
399 raise ValueError(
"Wrong value for the type of backbone")
401 h_chains = IMP.atom.get_by_type(hierarchy, IMP.atom.CHAIN_TYPE)
403 if len(h_chains) > 1:
404 raise ValueError(
"The hierarchy mas more than one chain")
405 h_residues = IMP.atom.get_by_type(hierarchy, IMP.atom.RESIDUE_TYPE)
406 for hr
in h_residues:
408 if not (res.get_is_dna()
or res.get_is_rna()):
409 raise ValueError(
"Residue is not part of a nucleic acid")
410 h_atoms = IMP.atom.get_by_type(hr, IMP.atom.ATOM_TYPE)
417 def get_calphas(chain_hierarchy):
418 h_residues = IMP.atom.get_by_type(chain_hierarchy, IMP.atom.RESIDUE_TYPE)
426 Get the backbone atoms for a hierarchy. It can be a protein or a
429 h_residues = IMP.atom.get_by_type(hierarchy, IMP.atom.RESIDUE_TYPE)
430 if len(h_residues) == 0:
431 raise ValueError(
"No residues!")
434 if res.get_is_dna()
or res.get_is_rna():
437 atoms = get_calphas(hierarchy)
443 Gets all the members of a set of rigid bodies, removing the subrigid
444 bodies. Returns all the plain atom or bead members
448 for rb
in rigid_bodies:
449 members += get_simple_members(rb)
453 def get_simple_members(rb):
455 members = [m
for m
in rb.get_members()
Select non water and non hydrogen atoms.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def get_selection_as_atom_hierarchy
Gets a selection of particles and decorates them as Atoms.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
A decorator to associate a particle with a part of a protein/DNA/RNA.
def get_nucleic_acid_backbone
Returns the atoms in the backbone of the nucleic acid contained in the hierarchy. ...
static Fragment setup_particle(Model *m, ParticleIndex pi)
static XYZR setup_particle(Model *m, ParticleIndex pi)
def create_simplified_assembly
Simplifies an assembly, by creating a hierarchy with one ball per n_res residues. ...
def get_residue_particle
Get the particle for a residue in a hierarchy.
def create_rigid_bodies
set the children of a molecule type hierarchy as rigid bodies In this case, all the children are the ...
void add_radii(Hierarchy d, const ForceFieldParameters *ffp=get_all_atom_CHARMM_parameters(), FloatKey radius_key=FloatKey("radius"))
def get_selection_rigid_body
Build the rigid body for all the particles in the selection S.
def get_component
Select a component of the assembly using the name.
def get_backbone
Get the backbone atoms for a hierarchy.
def apply_rotation_around_centroid
Rotates the reference frame of a rigid body around the centroid.
static Residue setup_particle(Model *m, ParticleIndex pi, ResidueType t, int index, int insertion_code)
Atom get_atom(Residue rd, AtomType at)
Return a particle atom from the residue.
def get_all_members
Gets all the members of a set of rigid bodies, removing the subrigid bodies.
def get_all_chains
Gets all the chains in a set of hierarchies.
void read_pdb(TextInput input, int model, Hierarchy h)
def get_rb_name
Name to use for the rigid body of a hierarch.
static Molecule setup_particle(Model *m, ParticleIndex pi)
def apply_transformation_around_centroid
Applies a transformation around the centroid of a rigid body.
def create_assembly_from_pdb
Builds the assembly setting the chains in the PDB file as components.
static Hierarchy setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor children=ParticleIndexesAdaptor())
Create a Hierarchy of level t by adding the needed attributes.
The standard decorator for manipulating molecular structures.
A decorator for a particle representing an atom.
static Mass setup_particle(Model *m, ParticleIndex pi, Float mass)
static Hierarchy setup_particle(Model *m, ParticleIndex pi, DecoratorTraits tr=get_default_decorator_traits())
def create_assembly
Read all the PDBs given in the list of names fn_pdbs and adds the hierarchies to the model...
def get_coarse_selection
The function returns the particles (fragments) in the coarse hierarchy that were created by summarizi...
A decorator for a particle with x,y,z coordinates.
def get_residues_distance
Distance between two residues.
Transformation3D compose(const Transformation3D &a, const Transformation3D &b)
Compose two transformations.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def rename_chains
Rename all the chains of an assembly so there are no conflicts with the ids.
algebra::Sphere3D get_enclosing_sphere(const XYZs &v)
Get a sphere enclosing the set of XYZRs.
A decorator for a residue.
Basic functionality that is expected to be used by a wide variety of IMP users.
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def get_rigid_body
Select a rigid body from the rigid_bodies using the name.
Hierarchy create_simplified_along_backbone(Chain input, const IntRanges &residue_segments, bool keep_detailed=false)
IMP::core::RigidBody create_rigid_body(Hierarchy h)
def get_residue_coordinates
Get the coordinates of a residue (the coordinates of the first particle)
Class to handle individual particles of a Model object.
double get_distance(const VectorD< D > &v1, const VectorD< D > &v2)
Compute the distance between two vectors.
Store info for a chain of a protein.
A decorator for a rigid body.
def create_simplified_dna
Gets a hierarchy containing a molecule of DNA and simplifies it, generating a coarse representation o...
Functionality for loading, creating, manipulating and scoring atomic structures.
def read_component
Read a PDB molecule, add atoms, and set a name.
static Chain setup_particle(Model *m, ParticleIndex pi, std::string id)
Hierarchies get_leaves(const Selection &h)
A decorator for a molecule.
Select hierarchy particles identified by the biological name.
static RigidBody setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor ps)
A decorator for a particle with x,y,z coordinates and a radius.