IMP logo
IMP Reference Guide  2.5.0
The Integrative Modeling Platform
EMageFit/imp_general/representation.py
1 """@namespace IMP.EMageFit.imp_general.representation
2  Utility functions to handle representation.
3 """
4 
5 import IMP
6 import IMP.atom as atom
7 import IMP.core as core
8 import IMP.algebra as alg
9 import IMP.display as display
10 import IMP.container as container
11 import string
12 from math import pi
13 import logging
14 
15 log = logging.getLogger("representation")
16 
17 #
18 """
19 
20  Functions to deal with the representation of assemblies and managing
21  rigid bodies.
22 
23 """
24 #
25 
26 
27 def create_assembly_from_pdb(model, fn_pdb, names=False):
28  """
29  Builds the assembly setting the chains in the PDB file as components
30  """
31  temp = read_component(model, fn_pdb)
32  hchains = atom.get_by_type(temp, atom.CHAIN_TYPE)
33  ids = [atom.Chain(h).get_id() for h in hchains]
34  log.debug("Creating assembly from pdb %s,names: %s. Chains %s",
35  fn_pdb, names, ids)
36  atom.add_radii(temp)
37  if(names):
38  for i, h in enumerate(hchains):
39  h.set_name(names[i])
40  assembly = IMP.atom.Molecule.setup_particle(temp)
41  return assembly
42 
43 
44 def create_assembly(model, fn_pdbs, names=False):
45  """ Read all the PDBs given in the list of names fn_pdbs and adds the
46  hierarchies to the model
47  """
49  for i, fn_pdb in enumerate(fn_pdbs):
50  if(names):
51  prot = read_component(model, fn_pdb, names[i])
52  else:
53  prot = read_component(model, fn_pdb)
54  assembly.add_child(prot)
55  return assembly
56 
57 
58 def read_component(model, fn_pdb, name=False):
59  """ Read a PDB molecule, add atoms, and set a name
60  """
61  if name:
62  log.debug("reading component %s from %s", name, fn_pdb)
63  else:
64  log.debug("reading component from %s", fn_pdb)
65 
66  hierarchy = atom.read_pdb(fn_pdb, model,
68  if name:
69  hierarchy.set_name(name)
70  atom.add_radii(hierarchy)
71  return hierarchy
72 
73 
74 def create_rigid_bodies(assembly):
75  """ set the children of a molecule type hierarchy as rigid bodies
76  In this case, all the children are the components of the complex.
77  I use the function create_rigid_body(), that creates a lot of
78  sub-rigid bodies.
79 
80  I have changed the function and now build the rigid body directly from
81  the leaves of each of the components. With this I guarantee that the
82  number of rigid members is going to be the same if the components have
83  the same number of atoms.
84  """
85  molecule = assembly.get_as_molecule()
86  if(not molecule.get_is_valid(True)):
87  raise TypeError("create_rigid_bodies(): The argument is not a valid "
88  "hierarchy")
89  rbs = []
90  for c in molecule.get_children():
91  p = IMP.Particle(c.get_model())
92  core.RigidBody.setup_particle(p, atom.get_leaves(c))
93  rb = core.RigidBody(p)
94 # rb = atom.create_rigid_body(c)
95  rb.set_name(get_rb_name(c.get_name()))
96  rbs.append(rb)
97  return rbs
98 
99 
100 def rename_chains(assembly):
101  """ Rename all the chains of an assembly so there are no conflicts with
102  the ids. The names are added sequentially.
103  """
104  m = assembly.get_as_molecule()
105  if(not m.get_is_valid(True)):
106  raise TypeError("The argument is not a valid hierarchy")
107  all_chains_as_hierarchies = get_all_chains(m.get_children())
108  letters = string.ascii_uppercase
109  n_chains = len(all_chains_as_hierarchies)
110  if(len(letters) < n_chains):
111  raise ValueError("There are more chains than letter ids")
112  ids = letters[0:n_chains]
113  for h, c_id in zip(all_chains_as_hierarchies, ids):
114  chain = atom.Chain(h.get_particle())
115  chain.set_id(c_id)
116  chain.set_name("chain %s" % c_id)
117 
118 
119 def create_simplified_dna(dna_hierarchy, n_res):
120  """ Gets a hierarchy containing a molecule of DNA and simplifies it,
121  generating a coarse representation of spheres. The function returns
122  a hierarchy with the spheres.
123  n_res - Number of residues to use per sphere.
124  """
125  chain = dna_hierarchy.get_as_chain()
126  if(not chain.get_is_valid(True)):
127  raise TypeError("create_simplified_dna: the hierarchy provided is not a "
128  "chain.")
129 
130  model = dna_hierarchy.get_model()
131  ph = IMP.Particle(model)
132  simplified_h = atom.Hierarchy.setup_particle(ph)
133  atom.Chain.setup_particle(ph, "0")
134 
135  residues = atom.get_by_type(dna_hierarchy, atom.RESIDUE_TYPE)
136  l = len(residues)
137  # print "the DNA has ",l,"residues"
138  for i in range(0, l, n_res):
139  xyzrs = []
140  equivalent_mass = 0.0
141  residues_numbers = []
142  for r in residues[i: i + n_res]:
143  rr = atom.Residue(r)
144  residues_numbers.append(rr.get_index())
145  # print "residue",rr.get_name(),rr.get_index()
146  residue_xyzrs = [core.XYZ(a.get_particle())
147  for a in rr.get_children()]
148  xyzrs += residue_xyzrs
149 # print "residue",r,"mass",get_residue_mass(r)
150  equivalent_mass += get_residue_mass(r)
151 
152  s = core.get_enclosing_sphere(xyzrs)
153  p = IMP.Particle(model)
154  xyzr = core.XYZR.setup_particle(p)
155  xyzr.set_radius(s.get_radius())
156  xyzr.set_coordinates(s.get_center())
157  fragment = atom.Fragment.setup_particle(p)
158  fragment.set_residue_indexes(residues_numbers)
159  atom.Mass.setup_particle(p, equivalent_mass)
160  simplified_h.add_child(fragment)
161  simplified_h.set_name("DNA")
162 # print "simplified_h is valid:",simplified_h.get_is_valid(True)
163  return simplified_h
164 
165 
166 def get_residue_mass(residue):
167  r = residue.get_as_residue()
168  if(r.get_is_valid(True) == False):
169  raise TypeError("The argument is not a residue")
170  mass = 0.0
171  for l in atom.get_leaves(r):
172  ms = atom.Mass(l)
173  mass += ms.get_residue_mass()
174  return mass
175 
176 
177 def create_simplified_assembly(assembly, components_rbs, n_res):
178  """ Simplifies an assembly, by creating a hierarchy with one ball per
179  n_res residues. Each of the chains in the new hierarchy are added to
180  the rigid bodies for each of the components.
181  There must be correspondence between the children of the assembly
182  (components) and the rigid bodies. I check for the ids.
183  """
184  molecule = assembly.get_as_molecule()
185  if(not molecule.get_is_valid(True)):
186  raise TypeError("The argument is not a valid hierarchy")
187 
188  model = assembly.get_model()
189  n_children = molecule.get_number_of_children()
190 
191  sh = IMP.Particle(model)
192  simplified_hierarchy = atom.Molecule.setup_particle(sh)
193 
194  for i in range(n_children): # for all members of the assembly
195  component = molecule.get_child(i)
196  name = component.get_name()
197  rb = components_rbs[i]
198  if(rb.get_name() != get_rb_name(name)):
199  raise ValueError("Rigid body and component do not match")
200 
201  hchains = atom.get_by_type(component, atom.CHAIN_TYPE)
202  ch = IMP.Particle(model)
203  coarse_component_h = atom.Molecule.setup_particle(ch)
204  # simplify all the chains in the member
205  for h in hchains:
206  chain = atom.Chain(h.get_particle())
207  coarse_h = None
208  if(name == "DNA"):
209  # print "simplifying DNA"
210  coarse_h_particle = create_simplified_dna(h, n_res)
211  coarse_h = atom.Hierarchy(coarse_h_particle)
212  else:
213  coarse_h = atom.create_simplified_along_backbone(chain, n_res)
214 
215  # does not work for DNA
216  chain_rb = atom.create_rigid_body(coarse_h)
217  # chain_rb = atom.setup_as_rigid_body(coarse_h) # works with DNA
218  chain_rb.set_name("sub_rb" + name)
219  rb.add_member(chain_rb)
220 
221  # are required to have excluded volume
222  coarse_component_h.add_child(atom.Chain(coarse_h))
223  coarse_component_h.set_name(name)
224  simplified_hierarchy.add_child(coarse_component_h)
225  return simplified_hierarchy
226 
227 
228 def get_component(assembly, name):
229  """ Select a component of the assembly using the name """
230  for c in assembly.get_children():
231  if (c.get_name() == name):
232  return c
233  raise ValueError(
234  "The requested component %s is not in the assembly" %
235  name)
236 
237 
238 def get_rigid_body(rigid_bodies, name):
239  """ Select a rigid body from the rigid_bodies using the name """
240  for rb in rigid_bodies:
241  if (rb.get_name() == name):
242  return rb
243  raise ValueError("This rigid body is not in the set: %s" % name)
244 
245 
246 def get_rb_name(name):
247  """ Name to use for the rigid body of a hierarch"""
248  return "rb_" + name
249 
250 
252  """ Build the rigid body for all the particles in the selection S """
253  ps = S.get_selected_particles()
254  xyzrs = [core.XYZR(p) for p in ps]
255  p_rbS = IMP.Particle(model)
256  rbS = core.RigidBody.setup_particle(p_rbS, xyzrs)
257  return rbS
258 
259 
260 def get_selection_as_hierarchy(model, S):
261  ph = IMP.Particle(model)
262  h = core.Hierarchy.setup_particle(ph)
263  for p in S.get_selected_particles():
264  x = core.Hierarchy.setup_particle(p)
265  h.add_child(x)
266  return h
267 
268 
270  """ Gets a selection of particles and decorates them as Atoms.
271  Then all of them are put into a big residue. I have this to use
272  with the multifit.create_coarse_molecule_from_molecule() function
273  """
274  ph = IMP.Particle(model)
275  h = atom.Residue.setup_particle(ph)
276  for p in S.get_selected_particles():
277  h.add_child(atom.Atom(p))
278  return h
279 
280 
281 def get_coarse_selection(coarse_h, residues_numbers):
282  """ The function returns the particles (fragments) in the coarse hierarchy
283  that were created by summarizing the residues_numbers.
284 
285  Coarse hierarchy - Hierarchy formed by a bunch of
286  fragments. Each fragment must have the residue numbers that it contains
287  residue_numbers - list with the number of the residues that need to
288  be recovered.
289  The function returns the set of particles that are atom.Fragments
290  """
291  particles = []
292  fragments = atom.get_by_type(coarse_h, atom.FRAGMENT_TYPE)
293  for f in fragments:
294  ff = atom.Fragment(f)
295  residues_in_f = ff.get_residue_indexes()
296  for number in residues_in_f:
297  if number in residues_numbers:
298  particles.append(ff.get_particle())
299  break
300  return particles
301 
302 
304  """
305  Rotates the reference frame of a rigid body around the centroid
306  """
307  c = rb.get_coordinates()
308  ref = rb.get_reference_frame()
309  R = ref.get_transformation_to().get_rotation()
310  R2 = alg.compose(rot, R)
311  T = alg.Transformation3D(R2, c)
312  ref = alg.ReferenceFrame3D(T)
313  rb.set_reference_frame(ref)
314 
315 
317  """
318  Aplies a transformation around the centroid of a rigid body.
319  First does the rotation around the centroid and
320  then applies the transformation.
321  @param rb A IMP.core.RigidBody object
322  @param T a IMP.algebra.Transformation3D object
323  """
324  apply_rotation_around_centroid(rb, T.get_rotation())
325  rb.set_coordinates(rb.get_coordinates() + T.get_translation())
326 
327 
328 def get_residue_particle(h, chain_id=False, res=1):
329  """
330  Get the particle for a residue in a hierarchy.
331  @param h The hierarchy
332  @param chain_id If chain_id == False, just search for the residue
333  @param res Number of residue in the chain
334  """
335 # log.debug("get_residue_particle: chain_id %s, res %s",chain_id, res)
336  if(chain_id):
337  s = IMP.atom.Selection(h, chain=chain_id, residue_index=res)
338  else:
339  s = IMP.atom.Selection(h, residue_index=res)
340  return s.get_selected_particles()[0]
341 
342 
343 def get_residue_coordinates(h, chain_id=False, res=1):
344  """
345  Get the coordinates of a residue (the coordinates of the first particle)
346  @param h Hierarchy
347  @param chain_id See help for get_residue_particle()
348  @param res See help for get_residue_particle()
349  """
350  p = get_residue_particle(h, chain_id, res)
351  return core.XYZ(p).get_coordinates()
352 
353 
354 def get_residues_distance(hierarchy1, chain_id1, residue1,
355  hierarchy2, chain_id2, residue2):
356  """
357  Distance between two residues. See the help for get_residue_particle()
358  @param hierarchy1
359  @param chain_id1
360  @param residue1
361  @param hierarchy2
362  @param chain_id2
363  @param residue2
364  """
365  coords1 = get_residue_coordinates(hierarchy1, chain_id1, residue1)
366  coords2 = get_residue_coordinates(hierarchy2, chain_id2, residue2)
367  d = alg.get_distance(coords1, coords2)
368  return d
369 
370 
371 def get_all_chains(hierarchies):
372  """ Gets all the chains in a set of hierarchies
373  @param hierarchies A set of IMP.atom.Hierarchy objects
374  """
375  chains = []
376  for h in hierarchies:
377  chains_in_h = atom.get_by_type(h, atom.CHAIN_TYPE)
378  for ch in chains_in_h:
379  chains.append(ch.get_as_chain())
380  return chains
381 
382 
383 def set_reference_frames(rigid_bodies, reference_frames):
384  for ref, rb in zip(reference_frames, rigid_bodies):
385  rb.set_reference_frame(ref)
386 
387 
388 def get_nucleic_acid_backbone(hierarchy, backbone='minimal'):
389  """
390  Returns the atoms in the backbone of the nucleic acid contained
391  in the hierarchy.
392  backbone 'minimal' returns the atoms: ["P", "O5'", "C5'", "C4'", "C3'", "O3'"]
393  backbone 'trace' returns the atoms C4'
394  """
395 # log.debug("get_nucleic_acid_backbone")
396  backbone_atoms = []
397  if backbone == 'minimal':
398  backbone_atoms = ["P", "O5'", "C5'", "C4'", "C3'", "O3'"]
399  elif backbone == 'trace':
400  backbone_atoms = ["C4'"]
401  else:
402  raise ValueError("Wrong value for the type of backbone")
403  backbone_atom_types = [atom.AtomType(t) for t in backbone_atoms]
404  h_chains = atom.get_by_type(hierarchy, atom.CHAIN_TYPE)
405  backbone = []
406  if len(h_chains) > 1:
407  raise ValueError("The hierarchy mas more than one chain")
408  h_residues = atom.get_by_type(hierarchy, atom.RESIDUE_TYPE)
409  for hr in h_residues:
410  res = atom.Residue(hr)
411  if not (res.get_is_dna() or res.get_is_rna()):
412  raise ValueError("Residue is not part of a nucleic acid")
413  h_atoms = atom.get_by_type(hr, atom.ATOM_TYPE)
414  for at in h_atoms:
415  if atom.Atom(at).get_atom_type() in backbone_atom_types:
416  backbone.append(at)
417  return backbone
418 
419 
420 def get_calphas(chain_hierarchy):
421  h_residues = atom.get_by_type(chain_hierarchy, atom.RESIDUE_TYPE)
422  cas = [atom.get_atom(atom.Residue(r), atom.AtomType("CA"))
423  for r in h_residues]
424  return cas
425 
426 
427 def get_backbone(hierarchy):
428  """
429  Get the backbone atoms for a hierarchy. It can be a protein or a
430  nucleic acid
431  """
432  h_residues = atom.get_by_type(hierarchy, atom.RESIDUE_TYPE)
433  if len(h_residues) == 0:
434  raise ValueError("No residues!")
435  atoms = []
436  res = atom.Residue(h_residues[0])
437  if res.get_is_dna() or res.get_is_rna():
438  atoms = get_nucleic_acid_backbone(hierarchy, 'trace')
439  else:
440  atoms = get_calphas(hierarchy)
441  return atoms
442 
443 
444 def get_all_members(rigid_bodies):
445  """
446  Gets all the members of a set of rigid bodies, removing the subrigid
447  bodies. Returns all the plain atom or bead members
448  """
449 
450  members = []
451  for rb in rigid_bodies:
452  members += get_simple_members(rb)
453  return members
454 
455 
456 def get_simple_members(rb):
457  # Add members if they are not sub-rigid bodies
458  members = [m for m in rb.get_members()
459  if not core.RigidBody.get_is_setup(m.get_particle())]
460  return members
Select non water and non hydrogen atoms.
Definition: pdb.h:197
Add mass to a particle.
Definition: Mass.h:23
def get_selection_as_atom_hierarchy
Gets a selection of particles and decorates them as Atoms.
A decorator to associate a particle with a part of a protein/DNA/RNA.
Definition: Fragment.h:20
def get_nucleic_acid_backbone
Returns the atoms in the backbone of the nucleic acid contained in the hierarchy. ...
Various classes to hold sets of particles.
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 ...
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.
The type of an atom.
def apply_rotation_around_centroid
Rotates the reference frame of a rigid body around the centroid.
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.
def get_rb_name
Name to use for the rigid body of a hierarch.
static Molecule setup_particle(Model *m, ParticleIndex pi)
Definition: Molecule.h:37
def apply_transformation_around_centroid
Aplies 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.
The standard decorator for manipulating molecular structures.
A decorator for a particle representing an atom.
Definition: atom/Atom.h:234
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.
Definition: XYZ.h:30
def rename_chains
Rename all the chains of an assembly so there are no conflicts with the ids.
A decorator for a residue.
Definition: Residue.h:134
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...
def get_rigid_body
Select a rigid body from the rigid_bodies using the name.
def get_residue_coordinates
Get the coordinates of a residue (the coordinates of the first particle)
Class to handle individual model particles.
Definition: Particle.h:37
Store info for a chain of a protein.
Definition: Chain.h:21
A decorator for a rigid body.
Definition: rigid_bodies.h:75
def create_simplified_dna
Gets a hierarchy containing a molecule of DNA and simplifies it, generating a coarse representation o...
Output IMP model data in various file formats.
Functionality for loading, creating, manipulating and scoring atomic structures.
def read_component
Read a PDB molecule, add atoms, and set a name.
Select hierarchy particles identified by the biological name.
Definition: Selection.h:65
A decorator for a particle with x,y,z coordinates and a radius.
Definition: XYZR.h:27