IMP  2.0.0
The Integrative Modeling Platform
atom/Hierarchy.h
Go to the documentation of this file.
1 /**
2  * \file IMP/atom/Hierarchy.h
3  * \brief Decorator for helping deal with a hierarchy of molecules.
4  *
5  * Copyright 2007-2013 IMP Inventors. All rights reserved.
6  *
7  */
8 
9 #ifndef IMPATOM_HIERARCHY_H
10 #define IMPATOM_HIERARCHY_H
11 
12 #include <IMP/atom/atom_config.h>
13 #include <IMP/core/utility.h>
14 #include <IMP/core/Hierarchy.h>
15 #include "bond_decorators.h"
16 #include "atom_macros.h"
17 #include <IMP/core/XYZR.h>
18 #include <IMP/core/rigid_bodies.h>
19 
20 #include <IMP/Particle.h>
21 #include <IMP/Model.h>
22 
23 #include <vector>
24 #include <deque>
25 
26 
27 #define IMP_GET_AS_DECL(UCName, lcname, CAPSNAME) \
28  UCName get_as_##lcname() const;
29 
30 // figure out how to inline
31 #define IMP_GET_AS_DEF(UCName, lcname, CAPSNAME) \
32  UCName Hierarchy::get_as_##lcname() const { \
33  if (UCName::particle_is_instance(get_particle())) { \
34  return UCName(get_particle()); \
35  } else { \
36  return UCName(); \
37  } \
38  }
39 
40 // DOMAIN is defined to be 1 by a fedora math header
41 #define IMP_FOREACH_HIERARCHY_TYPE_LIST(macro) \
42  macro(Atom, atom, ATOM_TYPE), \
43  macro(Residue, residue, RESIDUE_TYPE), \
44  macro(Chain, chain, CHAIN_TYPE), \
45  macro(Molecule, molecule, MOLECULE_TYPE), \
46  macro(Domain, domain, DOMAIN_TYPE), \
47  macro(Fragment, fragment, FRAGMENT_TYPE), \
48  macro(core::XYZ, xyz, XYZ_TYPE), \
49  macro(core::XYZR, xyzr, XYZR_TYPE), \
50  macro(Mass, mass, MASS_TYPE)
51 
52 // DOMAIN is defined to be 1 by a fedora math header
53 #define IMP_FOREACH_HIERARCHY_TYPE_STATEMENTS(macro) \
54  macro(Atom, atom, ATOM_TYPE); \
55  macro(Residue, residue, RESIDUE_TYPE); \
56  macro(Chain, chain, CHAIN_TYPE); \
57  macro(Molecule, molecule, MOLECULE_TYPE); \
58  macro(Domain, domain, DOMAIN_TYPE); \
59  macro(Fragment, fragment, FRAGMENT_TYPE); \
60  macro(core::XYZ, xyz, XYZ_TYPE); \
61  macro(core::XYZR, xyzr, XYZR_TYPE); \
62  macro(Mass, mass, MASS_TYPE)
63 
64 // DOMAIN is defined to be 1 by a fedora math header
65 #define IMP_FOREACH_HIERARCHY_TYPE_FUNCTIONS(macro) \
66  macro(Atom, atom, ATOM_TYPE) \
67  macro(Residue, residue, RESIDUE_TYPE) \
68  macro(Chain, chain, CHAIN_TYPE) \
69  macro(Molecule, molecule, MOLECULE_TYPE) \
70  macro(Domain, domain, DOMAIN_TYPE) \
71  macro(Fragment, fragment, FRAGMENT_TYPE) \
72  macro(core::XYZ, xyz, XYZ_TYPE) \
73  macro(core::XYZR, xyzr, XYZR_TYPE) \
74  macro(Mass, mass, MASS_TYPE) \
75  IMP_REQUIRE_SEMICOLON_NAMESPACE
76 
77 
78 #define IMP_CAPS_NAME(UCName, lcname, CAPSNAME) \
79  CAPSNAME
80 
81 
82 IMPATOM_BEGIN_NAMESPACE
83 class Atom;
84 class Residue;
85 class Domain;
86 class Fragment;
87 class Chain;
88 class Molecule;
89 class Mass;
90 
91 IMP_DECORATORS_DECL(Hierarchy, Hierarchies);
92 
93 
94 //! The standard decorator for manipulating molecular structures.
95 /** \imp represents molecular structures using the Hierachy decorator.
96  Molecules and collections of molecules each are stored as a
97  hierarchy (or tree) where the resolution of the representation increases
98  as you move further from the root. That is, if a parent has
99  some particular property (eg, marks out a volume by having
100  a x,y,z coordinates and a radius), then the children should have
101  a higher resolution form of that information (eg, mark out a more
102  detailed excluded volume by defining a set of balls which having
103  approximately the same total volume).
104 
105  \section tree_basics Tree Basics
106  In a tree you have a set of nodes, represented by Hierarchy particles.
107  Each node can have a node can have at most one parent. The node with no
108  parent is known as the root of the tree.
109 
110  Here is a simple example with a protein with three residues. Two of the
111  residues have atoms, where as the third is coarse grained.
112  \dotgraph{\dot
113  digraph example {
114  node [shape=record\, fontname= Helvetica\, fontsize=10]
115  a [label="Protein A (the root)"];
116  b [label="Residue 0"\, URL="Residue"];
117  c [label="Residue 1"];
118  cp [label="Residue 2"];
119  d0 [label="CA"];
120  e0 [label="CA"];
121  d1 [label="C"];
122  e1 [label="C"];
123  d2 [label="N"];
124  e2 [label="N"];
125  a -> b [arrowhead="open"];
126  a -> c [arrowhead="open"]
127  a -> cp [arrowhead="open"];
128  b -> d0 [arrowhead="open"];
129  c -> e0 [arrowhead="open"];
130  b -> d1 [arrowhead="open"];
131  c -> e1 [arrowhead="open"];
132  b -> d2 [arrowhead="open"];
133  c -> e2 [arrowhead="open"];
134  }
135  \enddot
136  }
137 
138 
139  The nodes in the hierarchy can correspond to arbitrary bits of a
140  molecule and do not need to have any biological significant. For
141  example we could introduce a fragment containing residues 0 and 1:
142  \dotgraph{\dot
143  digraph example {
144  node [shape=record\, fontname= Helvetica\, fontsize=10]
145  a [label="Protein A (the root)"\, URL="\ref B"];
146  aa [label="Fragment 0"];
147  b [label="Residue 0"];
148  c [label="Residue 1"];
149  cp [label="Residue 2"];
150  d0 [label="CA"];
151  e0 [label="CA"];
152  d1 [label="C"];
153  e1 [label="C"];
154  d2 [label="N"];
155  e2 [label="N"];
156  a -> aa [arrowhead="open"];
157  aa -> b [arrowhead="open"];
158  aa -> c [arrowhead="open"]
159  a -> cp [arrowhead="open"];
160  b -> d0 [arrowhead="open"];
161  c -> e0 [arrowhead="open"];
162  b -> d1 [arrowhead="open"];
163  c -> e1 [arrowhead="open"];
164  b -> d2 [arrowhead="open"];
165  c -> e2 [arrowhead="open"];
166  }
167  \enddot}
168 
169 
170  A hierarchy can have any tree structure as long as:
171  - the type of the parent makes sense for the child: eg a Residue
172  cannot be the parent of a Chain.
173  - the leaves always have coordinates, radius and mass
174  - all particles in hierarchy are from the same model
175  - all Atoms has a Residue for as parent
176  - any Atom with a non-heterogen atom type is part of a protein,
177  DNA or RNA molecule.
178  - all Residue children of a particle appear in order based
179  on their index
180  - all Atom children in of a particle appear in order of their
181  AtomType
182 
183  The get_is_valid() method checks some of these properties. Any
184  method taking a hierarchy as an argument should do
185  \code
186  IMP_USAGE_CHECK(h.get_is_valid(), "Invalid hierarchy as input");
187  \endcode
188  to make sure the hierarchy makes sense.
189 
190  A number of decorator types are associated with the Hierarchy
191  to store the information associated with that node in the
192  hierarchy. Examples include Residue, Atom, XYZ, Chain, XYZR,
193  Mass, Domain, Molecule etc.
194  We provide a get_as_x() function for each such decorator which
195  returns either X() (a null type) if the node is not a particle
196  of type x, or an X decorator wrapping the current particle if
197  it is.
198 
199  \ingroup hierarchy
200  \ingroup decorators
201  \see Atom
202  \see Residue
203  \see Chain
204  \see Molecule
205  \see Domain
206  \see Fragment
207  \see Mass
208  */
209 class IMPATOMEXPORT Hierarchy: public core::Hierarchy
210 {
211  typedef core::Hierarchy H;
212 public:
213  IMP_NO_DOXYGEN(typedef boost::false_type DecoratorHasTraits);
214  explicit Hierarchy(Particle *p): H(p, get_traits()) {
215  }
216 
217  Hierarchy(Model *m, ParticleIndex pi): H(m, pi, get_traits()) {
218  }
219 
220  //! null constructor
222 
223  //! cast a particle which has the needed attributes
225  H::decorate_particle(p, get_traits());
226  return Hierarchy(p);
227  }
228 
229  //! The traits must match
230  explicit Hierarchy(IMP::core::Hierarchy h): H(h) {
232  || h.get_traits() == get_traits(),
233  "Cannot construct a IMP.atom.Hierarchy from a general "
234  " IMP.core.Hierarchy");
235  }
236 
238  const ParticlesTemp &children) {
239  H::setup_particle(p, get_traits());
240  Hierarchy ret(p);
241  for (unsigned int i=0; i< children.size(); ++i) {
242  if (!particle_is_instance(children[i])) {
243  setup_particle(children[i]);
244  }
245  ret.add_child(Hierarchy(children[i]));
246  }
247  return ret;
248  }
249 
250  /** Create a Hierarchy of level t by adding the needed
251  attributes. */
253  const ParticlesTemp &children
254  = ParticlesTemp()) {
255  Particle *p= m->get_particle(pi);
256  H::setup_particle(p, get_traits());
257  Hierarchy ret(p);
258  for (unsigned int i=0; i< children.size(); ++i) {
259  if (!particle_is_instance(children[i])) {
260  setup_particle(children[i]);
261  }
262  ret.add_child(Hierarchy(children[i]));
263  }
264  return ret;
265  }
266 
267  static Hierarchy setup_particle(Particle *p) {
268  return setup_particle(p->get_model(),
269  p->get_index());
270  }
271 
272  /** Check if the particle has the needed attributes for a
273  cast to succeed */
275  return H::particle_is_instance(p, get_traits());
276  }
277 
278  static bool particle_is_instance(Model *m, ParticleIndex p){
279  return H::particle_is_instance(m->get_particle(p), get_traits());
280  }
281 
282  //! Return true if the hierarchy is valid.
283  /** Print information about the hierarchy if print_info is
284  true and things are invalid.
285  \note Returning true only means that no problems were
286  found, it can't check everything.*/
287  bool get_is_valid(bool print_info) const;
288  //! Add a child and check that the types are appropriate
289  /** A child must have a type that is listed before the parent in the
290  Type enum list.
291  */
293  IMP_USAGE_CHECK(o != *this, "Can't add something as its own child");
294  H::add_child(o);
295  }
296 
297  /** Get the ith child based on the order they were added. */
298  Hierarchy get_child(unsigned int i) const {
299  H hd= H::get_child(i);
300  return Hierarchy(hd);
301  }
302  //! Return the children in the order they were added
303  Hierarchies get_children() const {
304  Hierarchies ret(get_number_of_children());
305  for (unsigned int i=0; i< get_number_of_children(); ++i) {
306  ret[i]= get_child(i);
307  }
308  return ret;
309  }
310 
311  //! Get the children in a container of your choosing, eg ParticlesTemp
312  template <class C>
313  C get_children() const {
314  C ret(get_number_of_children());
315  for (unsigned int i=0; i< get_number_of_children(); ++i) {
316  ret[i]= get_child(i);
317  }
318  return ret;
319  }
320 
321  /** Get the parent particle. */
323  H hd= H::get_parent();
324  if (hd == H()) {
325  return Hierarchy();
326  } else {
327  return Hierarchy(hd);
328  }
329  }
330 
331  /** \name Methods to get associated decorators
332 
333  We provide a number of helper methods to get associated
334  decorators for the current node in the hierarchy. As an
335  example, if the particle decorated by this decorator is
336  a Residue particle, then get_as_residue() return
337  Residue(get_particle()), if not it returns Residue().
338  @{
339  */
340  IMP_FOREACH_HIERARCHY_TYPE_FUNCTIONS(IMP_GET_AS_DECL);
341  /** @} */
342 
343  //! Get the molecular hierarchy HierararchyTraits.
344  static const IMP::core::HierarchyTraits& get_traits();
345 
346  // swig overwrites __repr__ if it is inherited
348 };
349 
350 IMP_DECORATORS_DEF(Hierarchy, Hierarchies);
351 
352 
353 #ifdef IMP_DOXYGEN
354 /** The different types which can be passed to get_by_type()
355  */
356 enum GetByType {ATOM_TYPE, RESIDUE_TYPE, CHAIN_TYPE, MOLECULE_TYPE.
357  DOMAIN_TYPE, FRAGMENT_TYPE,
358  XYZ_TYPE,XYZR_TYPE,MASS_TYPE};
359 #else
360 enum GetByType {
361  IMP_FOREACH_HIERARCHY_TYPE_LIST(IMP_CAPS_NAME)
362 };
363 #endif
364 
365 /**
366  Gather all the molecular particles of a certain level
367  in the molecular hierarchy.
368  \ingroup hierarchy
369  \relatesalso Hierarchy
370 */
371 IMPATOMEXPORT Hierarchies
372 get_by_type(Hierarchy mhd, GetByType t);
373 
374 
375 //! Get the residue with the specified index
376 /** Find the leaf containing the residue with the appropriate index.
377  This is the PDB index, not the offset in the chain (if they are different).
378 
379  The function returns a Hierarchy, rather than a Residue since the
380  residue may not be explicitly represented and may just be part of some
381  fragment.
382 
383  \throw ValueException if mhd's type is not one of CHAIN, PROTEIN, NUCLEOTIDE
384  \return Hierarchy() if that residue is not found.
385 
386  \ingroup hierarchy
387  \relatesalso Hierarchy
388  */
389 IMPATOMEXPORT Hierarchy
390 get_residue(Hierarchy mhd, unsigned int index);
391 
392 
393 //! Create a fragment containing the specified nodes
394 /** A particle representing the fragment is created and initialized.
395 
396  The Fragment is inserted as a child of the parent (and the particles are
397  removed). The particles become children of the fragment.
398 
399  \throw ValueException If all the particles do not have the same parent.
400  \ingroup hierarchy
401  \relatesalso Hierarchy
402  */
403 IMPATOMEXPORT Hierarchy
404 create_fragment(const Hierarchies &ps);
405 
406 //! Get the bonds internal to this tree
407 /** \relatesalso Hierarchy
408  \see Bond
409  \relatesalso Bond
410  */
411 IMPATOMEXPORT Bonds
412 get_internal_bonds(Hierarchy mhd);
413 
414 
415 //! Return the root of the hierarchy
416 /** \relatesalso Hierarchy */
418  while (h.get_parent()) {
419  h= h.get_parent();
420  }
421  return h;
422 }
423 
424 /** \relatesalso Hierarchy */
425 inline Hierarchies get_leaves(Hierarchy h) {
426  return Hierarchies(IMP::core::get_leaves(h));
427 }
428 
429 /** \relatesalso Hierarchy */
430 inline Hierarchies get_leaves(const Hierarchies& h) {
431  ParticlesTemp ret;
432  for (unsigned int i=0; i< h.size(); ++i) {
433  core::GenericHierarchies cur=IMP::core::get_leaves(h[i]);
434  ret.insert(ret.end(), cur.begin(), cur.end());
435  }
436  return get_as<Hierarchies>(ret);
437 }
438 
439 //! Print out a molecular hierarchy
440 /** \relatesalso Hierarchy
441  */
442 inline void show(Hierarchy h, std::ostream &out=std::cout) {
443  IMP::core::show<Hierarchy>(h, out);
444 }
445 
446 //! Rigidify a molecule or collection of molecules.
447 /** The rigid body created has all the leaves as members and a
448  member rigid body for each internal node in the tree. The
449  particle created to be the rigid body is returned.
450 
451  A name can be passed as it is not easy to automatically pick
452  a decent name.
453  \see create_aligned_rigid_body()
454  \relatesalso Hierarchy
455  \relatesalso IMP::core::RigidBody
456 */
457 IMPATOMEXPORT IMP::core::RigidBody create_rigid_body(const Hierarchies& h,
458  std::string name=std::string("created rigid body"));
459 
460 /** \see create_rigid_body(const Hierarchies&)
461  */
462 IMPATOMEXPORT IMP::core::RigidBody create_rigid_body(Hierarchy h);
463 
464 //! Rigidify a molecule or collection of molecules.
465 /** This method is identical to create_rigid_body() except that
466  the chosen reference frame is aligned with that of reference
467  (which must have exactly the same set of particles). This allows
468  one to make sure the rigid body is equivalent when you have several
469  copies of the same molecule.
470 
471  \relatesalso Hierarchy
472  \relatesalso IMP::core::RigidBody
473 */
474 IMPATOMEXPORT IMP::core::RigidBody create_compatible_rigid_body(Hierarchy h,
475  Hierarchy reference);
476 
477 
478 #ifndef IMP_DOXYGEN
479 IMPATOMEXPORT IMP::core::RigidBody setup_as_rigid_body(Hierarchy h);
480 #endif
481 
482 //! Return true if the piece of hierarchy should be classified as a heterogen
483 /** For the purposes of classification, a heterogen is anything that
484  - is a heterogen atom (one whose name starts with HET:)
485  - is or is part of a Residue that is not a normal protein, rna or
486  dna residue
487  - or is not part of a Chain
488  For the moment, this can only be called on residues or atoms.
489 */
490 IMPATOMEXPORT bool get_is_heterogen(Hierarchy h);
491 
492 //! Clone the Hierarchy
493 /** This method copies the Bond, Bonded, Atom,
494  Residue, and Domain data and the particle name to the
495  new copies in addition to the Hierarchy relationships.
496 
497  \relatesalso Hierarchy
498 */
499 IMPATOMEXPORT
500 Hierarchy create_clone(Hierarchy d);
501 
502 //! Clone the node in the Hierarchy
503 /** This method copies the Atom,
504  Residue, Chain and Domain data and the particle name.
505 
506  \relatesalso Hierarchy
507 */
508 IMPATOMEXPORT
509 Hierarchy create_clone_one(Hierarchy d);
510 
511 
512 //! Delete the Hierarchy
513 /** All bonds connecting to these atoms are destroyed as are
514  hierarchy links in the Hierarchy and the particles are
515  removed from the Model. If this particle has a parent, it is
516  removed from the parent.
517  \relatesalso Hierarchy
518 */
519 IMPATOMEXPORT
520 void destroy(Hierarchy d);
521 
522 
523 
524 //! Get a bounding box for the Hierarchy
525 /** This bounding box is that of the highest (in the CS sense of a tree
526  growing down from the root) cut
527  through the tree where each node in the cut has x,y,z, and r.
528  That is, if the root has x,y,z,r then it is the bounding box
529  of that sphere. If only the leaves have radii, it is the bounding
530  box of the leaves. If no such cut exists, the behavior is undefined.
531  \relatesalso Hierarchy
532  \relatesalso IMP::algebra::BoundingBoxD
533  */
534 IMPATOMEXPORT
535 algebra::BoundingBoxD<3> get_bounding_box(const Hierarchy &h);
536 
537 
538 /** See get_bounding_box() for more details.
539  \relatesalso Hierarchy
540  */
541 IMPATOMEXPORT
542 algebra::Sphere3D get_bounding_sphere(const Hierarchy &h);
543 
544 
545 IMPATOM_END_NAMESPACE
546 
547 
548 #endif /* IMPATOM_HIERARCHY_H */