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 "
126 model = dna_hierarchy.get_model()
131 residues = IMP.atom.get_by_type(dna_hierarchy, IMP.atom.RESIDUE_TYPE)
133 for i
in range(0, lenr, n_res):
135 equivalent_mass = 0.0
136 residues_numbers = []
137 for r
in residues[i: i + n_res]:
139 residues_numbers.append(rr.get_index())
142 for a
in rr.get_children()]
143 xyzrs += residue_xyzrs
145 equivalent_mass += get_residue_mass(r)
150 xyzr.set_radius(s.get_radius())
151 xyzr.set_coordinates(s.get_center())
153 fragment.set_residue_indexes(residues_numbers)
155 simplified_h.add_child(fragment)
156 simplified_h.set_name(
"DNA")
161 def get_residue_mass(residue):
163 raise TypeError(
"The argument is not a residue")
168 mass += ms.get_residue_mass()
173 """ Simplifies an assembly, by creating a hierarchy with one ball per
174 n_res residues. Each of the chains in the new hierarchy are added to
175 the rigid bodies for each of the components.
176 There must be correspondence between the children of the assembly
177 (components) and the rigid bodies. I check for the ids.
180 raise TypeError(
"The argument is not a valid hierarchy")
183 model = assembly.get_model()
184 n_children = molecule.get_number_of_children()
189 for i
in range(n_children):
190 component = molecule.get_child(i)
191 name = component.get_name()
192 rb = components_rbs[i]
194 raise ValueError(
"Rigid body and component do not match")
196 hchains = IMP.atom.get_by_type(component, IMP.atom.CHAIN_TYPE)
212 chain_rb.set_name(
"sub_rb" + name)
213 rb.add_member(chain_rb)
217 coarse_component_h.set_name(name)
218 simplified_hierarchy.add_child(coarse_component_h)
219 return simplified_hierarchy
223 """ Select a component of the assembly using the name """
224 for c
in assembly.get_children():
225 if (c.get_name() == name):
228 "The requested component %s is not in the assembly" %
233 """ Select a rigid body from the rigid_bodies using the name """
234 for rb
in rigid_bodies:
235 if (rb.get_name() == name):
237 raise ValueError(
"This rigid body is not in the set: %s" % name)
241 """ Name to use for the rigid body of a hierarch"""
246 """ Build the rigid body for all the particles in the selection S """
247 ps = S.get_selected_particles()
254 def get_selection_as_hierarchy(model, S):
257 for p
in S.get_selected_particles():
264 """ Gets a selection of particles and decorates them as Atoms.
265 Then all of them are put into a big residue. I have this to use
266 with the multifit.create_coarse_molecule_from_molecule() function
270 for p
in S.get_selected_particles():
276 """ The function returns the particles (fragments) in the coarse hierarchy
277 that were created by summarizing the residues_numbers.
279 Coarse hierarchy - Hierarchy formed by a bunch of
280 fragments. Each fragment must have the residue numbers that it contains
281 residue_numbers - list with the number of the residues that need to
283 The function returns the set of particles that are IMP.atom.Fragments
286 fragments = IMP.atom.get_by_type(coarse_h, IMP.atom.FRAGMENT_TYPE)
289 residues_in_f = ff.get_residue_indexes()
290 for number
in residues_in_f:
291 if number
in residues_numbers:
292 particles.append(ff.get_particle())
299 Rotates the reference frame of a rigid body around the centroid
301 c = rb.get_coordinates()
302 ref = rb.get_reference_frame()
303 R = ref.get_transformation_to().get_rotation()
307 rb.set_reference_frame(ref)
312 Applies a transformation around the centroid of a rigid body.
313 First does the rotation around the centroid and
314 then applies the transformation.
315 @param rb A IMP.core.RigidBody object
316 @param T a IMP.algebra.Transformation3D object
319 rb.set_coordinates(rb.get_coordinates() + T.get_translation())
324 Get the particle for a residue in a hierarchy.
325 @param h The hierarchy
326 @param chain_id If chain_id == False, just search for the residue
327 @param res Number of residue in the chain
334 return s.get_selected_particles()[0]
339 Get the coordinates of a residue (the coordinates of the first
342 @param chain_id See help for get_residue_particle()
343 @param res See help for get_residue_particle()
350 hierarchy2, chain_id2, residue2):
352 Distance between two residues. See the help for get_residue_particle()
367 """ Gets all the chains in a set of hierarchies
368 @param hierarchies A set of IMP.atom.Hierarchy objects
371 for h
in hierarchies:
372 chains_in_h = IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE)
373 for ch
in chains_in_h:
378 def set_reference_frames(rigid_bodies, reference_frames):
379 for ref, rb
in zip(reference_frames, rigid_bodies):
380 rb.set_reference_frame(ref)
385 Returns the atoms in the backbone of the nucleic acid contained
387 backbone 'minimal' returns the atoms:
388 ["P", "O5'", "C5'", "C4'", "C3'", "O3'"]
389 backbone 'trace' returns the atoms C4'
393 if backbone ==
'minimal':
394 backbone_atoms = [
"P",
"O5'",
"C5'",
"C4'",
"C3'",
"O3'"]
395 elif backbone ==
'trace':
396 backbone_atoms = [
"C4'"]
398 raise ValueError(
"Wrong value for the type of backbone")
400 h_chains = IMP.atom.get_by_type(hierarchy, IMP.atom.CHAIN_TYPE)
402 if len(h_chains) > 1:
403 raise ValueError(
"The hierarchy mas more than one chain")
404 h_residues = IMP.atom.get_by_type(hierarchy, IMP.atom.RESIDUE_TYPE)
405 for hr
in h_residues:
407 if not (res.get_is_dna()
or res.get_is_rna()):
408 raise ValueError(
"Residue is not part of a nucleic acid")
409 h_atoms = IMP.atom.get_by_type(hr, IMP.atom.ATOM_TYPE)
416 def get_calphas(chain_hierarchy):
417 h_residues = IMP.atom.get_by_type(chain_hierarchy, IMP.atom.RESIDUE_TYPE)
425 Get the backbone atoms for a hierarchy. It can be a protein or a
428 h_residues = IMP.atom.get_by_type(hierarchy, IMP.atom.RESIDUE_TYPE)
429 if len(h_residues) == 0:
430 raise ValueError(
"No residues!")
433 if res.get_is_dna()
or res.get_is_rna():
436 atoms = get_calphas(hierarchy)
442 Gets all the members of a set of rigid bodies, removing the subrigid
443 bodies. Returns all the plain atom or bead members
447 for rb
in rigid_bodies:
448 members += get_simple_members(rb)
452 def get_simple_members(rb):
454 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"))
Add vdW radius from given force field.
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.