IMP logo
IMP Reference Guide  2.13.0
The Integrative Modeling Platform
charmm_topology.h
Go to the documentation of this file.
1 /**
2  * \file IMP/atom/charmm_topology.h
3  * \brief Classes for handling CHARMM-style topology.
4  *
5  * Copyright 2007-2020 IMP Inventors. All rights reserved.
6  *
7  */
8 
9 #ifndef IMPATOM_CHARMM_TOPOLOGY_H
10 #define IMPATOM_CHARMM_TOPOLOGY_H
11 
12 #include "IMP/Object.h"
13 #include "Hierarchy.h"
14 #include "Atom.h"
15 #include <IMP/atom/atom_config.h>
16 
17 #include <string>
18 #include <vector>
19 
20 IMPATOM_BEGIN_NAMESPACE
21 
22 //! A single atom in a CHARMM topology.
23 /** Each atom has a name (unique to the residue), a CHARMM type (used to
24  look up parameters such as radii and bond lengths in the parameter file)
25  and an electrostatic charge.
26  \see CHARMMAtom
27  */
29  std::string name_;
30  std::string charmm_type_;
31  double charge_;
32 
33  public:
34  CHARMMAtomTopology(std::string name) : name_(name) {}
35 
36  CHARMMAtomTopology(std::string name, const CHARMMAtomTopology &other)
37  : name_(name), charmm_type_(other.charmm_type_), charge_(other.charge_) {}
38 
39  std::string get_name() const { return name_; }
40  std::string get_charmm_type() const { return charmm_type_; }
41  double get_charge() const { return charge_; }
42  void set_charmm_type(std::string charmm_type) { charmm_type_ = charmm_type; }
43  void set_charge(double charge) { charge_ = charge; }
45  out << "name: " << name_ << "; CHARMM type: " << charmm_type_
46  << "; charge: " << charge_;
47  });
48 };
49 
51 
53 
54 //! The end of a bond, angle, dihedral, improper, or internal coordinate.
55 /** An endpoint is an atom in a residue; bonds have two endpoints,
56  angles three, and dihedrals, impropers and internal coordinates, four.
57 
58  If residue_ is not nullptr, the endpoint is the named atom in the
59  pointed-to residue. (This is used for bonds that can span multiple
60  residues, perhaps in different chains, and is typically only created
61  by applying a two-residue patch such as DISU or LINK).
62 
63  If residue_ is nullptr, the endpoint is the named atom in the "current"
64  residue, unless the atom name has a +, -, 1: or 2: prefix.
65  Atoms prefixed with + are found in the next residue in the chain.
66  Atoms prefixed with - are found in the previous residue in the chain.
67  1: or 2: prefixes are used by two-residue patches to refer to atoms in
68  each residue patched by the patch.
69  */
70 class IMPATOMEXPORT CHARMMBondEndpoint {
71  std::string atom_name_;
72  Pointer<Object> residue_;
73 
74  public:
75  CHARMMBondEndpoint(std::string atom_name,
76  CHARMMResidueTopology *residue = nullptr);
77 
78  std::string get_atom_name() const { return atom_name_; }
79 
80  //! Map the endpoint to an Atom particle.
81  Atom get_atom(
82  const CHARMMResidueTopology *current_residue,
83  const CHARMMResidueTopology *previous_residue,
84  const CHARMMResidueTopology *next_residue,
85  const std::map<const CHARMMResidueTopology *, Hierarchy> &resmap) const;
86 
87  IMP_SHOWABLE_INLINE(CHARMMBondEndpoint, { out << atom_name_; });
88 };
89 
91 
92 //! A connection (bond, angle, dihedral) between some number of endpoints.
93 template <unsigned int D>
95  protected:
96  Vector<CHARMMBondEndpoint> endpoints_;
97 
98  public:
99  CHARMMConnection(const IMP::Strings &atoms) {
100  IMP_INTERNAL_CHECK(atoms.size() == D, "wrong number of bond endpoints");
101  for (Strings::const_iterator it = atoms.begin(); it != atoms.end(); ++it) {
102  endpoints_.push_back(CHARMMBondEndpoint(*it));
103  }
104  }
105 
106 #ifndef SWIG
108  : endpoints_(endpoints) {
109  IMP_INTERNAL_CHECK(endpoints.size() == D, "wrong number of bond endpoints");
110  }
111 #endif
112 
113  const IMP::atom::CHARMMBondEndpoint &get_endpoint(unsigned int i) const {
114  return endpoints_[i];
115  }
116 
117  //! Return true if the bond contains the named atom.
118  bool get_contains_atom(std::string name) const {
120  endpoints_.begin();
121  it != endpoints_.end(); ++it) {
122  if (it->get_atom_name() == name) {
123  return true;
124  }
125  }
126  return false;
127  }
128 
129 #ifndef SWIG
130  //! Map the bond to a list of Atom particles.
132  const CHARMMResidueTopology *current_residue,
133  const CHARMMResidueTopology *previous_residue,
134  const CHARMMResidueTopology *next_residue,
135  const std::map<const CHARMMResidueTopology *, Hierarchy> &resmap) const {
136  Atoms as;
138  endpoints_.begin();
139  it != endpoints_.end(); ++it) {
140  Atom a =
141  it->get_atom(current_residue, previous_residue, next_residue, resmap);
142  if (a) {
143  as.push_back(a);
144  } else {
145  return Atoms();
146  }
147  }
148  return as;
149  }
150 #endif
153  endpoints_.begin();
154  it != endpoints_.end(); ++it) {
155  if (it != endpoints_.begin()) {
156  out << "-";
157  }
158  out << it->get_atom_name();
159  }
160  });
161 };
162 
163 typedef CHARMMConnection<2> CHARMMBond;
164 typedef CHARMMConnection<3> CHARMMAngle;
165 typedef CHARMMConnection<4> CHARMMDihedral;
166 
167 IMP_VALUES(CHARMMBond, CHARMMBonds);
168 
169 IMP_VALUES(CHARMMAngle, CHARMMAngles);
170 
171 IMP_VALUES(CHARMMDihedral, CHARMMDihedrals);
172 
173 #ifdef SWIG
174 // Ugly, but SWIG needs a template instantiation before it is used as
175 // a base class
176 %template(CHARMMDihedral) CHARMMConnection<4>;
177 #endif
178 
179 //! A geometric relationship between four atoms.
180 /** The atoms (denoted i,j,k,l here) are uniquely positioned in 3D space
181  relative to each other by means of two distances, two angles, and
182  a dihedral.
183 
184  A regular internal coordinate stores the distances between ij and
185  kl respectively, and the angles between ijk and jkl.
186 
187  An improper internal coordinate stores the distances between ik and
188  kl respectively, and the angles between ikj and jkl.
189 
190  In both cases the dihedral is the angle between the two planes formed
191  by ijk and jkl.
192  */
194  float first_distance_, second_distance_, first_angle_, second_angle_,
195  dihedral_;
196  bool improper_;
197 
198  public:
199  CHARMMInternalCoordinate(const IMP::Strings &atoms, float first_distance,
200  float first_angle, float dihedral,
201  float second_angle, float second_distance,
202  bool improper)
203  : CHARMMConnection<4>(atoms),
204  first_distance_(first_distance),
205  second_distance_(second_distance),
206  first_angle_(first_angle),
207  second_angle_(second_angle),
208  dihedral_(dihedral),
209  improper_(improper) {}
211  float first_distance, float first_angle,
212  float dihedral, float second_angle,
213  float second_distance, bool improper)
214  : CHARMMConnection<4>(endpoints),
215  first_distance_(first_distance),
216  second_distance_(second_distance),
217  first_angle_(first_angle),
218  second_angle_(second_angle),
219  dihedral_(dihedral),
220  improper_(improper) {}
221 
222  float get_first_distance() const { return first_distance_; }
223  float get_second_distance() const { return second_distance_; }
224  float get_first_angle() const { return first_angle_; }
225  float get_second_angle() const { return second_angle_; }
226  float get_dihedral() const { return dihedral_; }
227  bool get_improper() const { return improper_; }
228 
231  out << "; distances: " << first_distance_ << ", " << second_distance_
232  << "; angles: " << first_angle_ << ", " << second_angle_
233  << "; dihedral: " << dihedral_;
234  if (improper_) {
235  out << "; improper";
236  }
237  });
238 };
240 
241 //! Base class for all CHARMM residue-based topology
242 class IMPATOMEXPORT CHARMMResidueTopologyBase : public IMP::Object {
243  std::string type_;
244 
245  protected:
247  CHARMMBonds bonds_;
248  CHARMMAngles angles_;
249  CHARMMDihedrals dihedrals_;
250  CHARMMDihedrals impropers_;
251  CHARMMInternalCoordinates internal_coordinates_;
252 
253  CHARMMResidueTopologyBase(std::string type)
254  : Object("CharmmResidueTopoBase%1%"), type_(type) {
255  set_name(std::string("CHARMM residue ") + type);
256  }
257 
258  public:
259  std::string get_type() const { return type_; }
260 
261  unsigned int get_number_of_atoms() const { return atoms_.size(); }
262  const CHARMMAtomTopology &get_atom(unsigned int i) const { return atoms_[i]; }
263 
264  void add_atom(const CHARMMAtomTopology &atom);
265 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
266  CHARMMAtomTopology &get_atom(std::string name);
267 #endif
268  const CHARMMAtomTopology &get_atom(AtomType type) const {
269  return get_atom(type.get_string());
270  }
271 
272  const CHARMMAtomTopology &get_atom(std::string name) const;
273 
274  unsigned int get_number_of_bonds() const { return bonds_.size(); }
275  void add_bond(const CHARMMBond &bond) { bonds_.push_back(bond); }
276 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
277  // cannot return non const ref to swig
278  CHARMMBond &get_bond(unsigned int index) { return bonds_[index]; }
279 #endif
280 
281  unsigned int get_number_of_angles() const { return angles_.size(); }
282  void add_angle(const CHARMMAngle &bond) { angles_.push_back(bond); }
283 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
284  // cannot return non const ref to swig
285  CHARMMAngle &get_angle(unsigned int index) { return angles_[index]; }
286 #endif
287  unsigned int get_number_of_dihedrals() const { return dihedrals_.size(); }
288  void add_dihedral(const CHARMMDihedral &bond) { dihedrals_.push_back(bond); }
289 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
290  // cannot return non const ref to swig
291  CHARMMDihedral &get_dihedral(unsigned int index) { return dihedrals_[index]; }
292 #endif
293  unsigned int get_number_of_impropers() const { return impropers_.size(); }
294  void add_improper(const CHARMMDihedral &bond) { impropers_.push_back(bond); }
295 #if !defined(SWIG) && !defined(IMP_DOXYGEN)
296  CHARMMDihedral &get_improper(unsigned int index) { return impropers_[index]; }
297 #endif
298 
299  unsigned int get_number_of_internal_coordinates() const {
300  return internal_coordinates_.size();
301  }
302  void add_internal_coordinate(const CHARMMInternalCoordinate &ic) {
303  internal_coordinates_.push_back(ic);
304  }
305  const CHARMMInternalCoordinate &get_internal_coordinate(unsigned int index)
306  const {
307  return internal_coordinates_[index];
308  }
309 
310  const CHARMMBond &get_bond(unsigned int index) const { return bonds_[index]; }
311  const CHARMMAngle &get_angle(unsigned int index) const {
312  return angles_[index];
313  }
314  const CHARMMDihedral &get_dihedral(unsigned int index) const {
315  return dihedrals_[index];
316  }
317  const CHARMMDihedral &get_improper(unsigned int index) const {
318  return impropers_[index];
319  }
320 
322 };
323 
325 
326 //! The ideal topology of a single residue
327 /** These residue topologies can be constructed manually (by adding
328  CHARMMAtomTopology objects with add_atom()). However, they are more
329  commonly populated automatically by reading a CHARMM topology file
330  (CHARMMParameters::get_residue_topology()).
331  */
332 class IMPATOMEXPORT CHARMMIdealResidueTopology
333  : public CHARMMResidueTopologyBase {
334  std::string default_first_patch_, default_last_patch_;
335 
336  public:
337  CHARMMIdealResidueTopology(std::string type)
338  : CHARMMResidueTopologyBase(type) {}
340  : CHARMMResidueTopologyBase(type.get_string()) {}
341 
342  //! Remove the named atom
343  /** Any bonds/angles that involve this atom are also removed.
344  */
345  void remove_atom(std::string name);
346 
347  void set_default_first_patch(std::string patch) {
348  default_first_patch_ = patch;
349  }
350  void set_default_last_patch(std::string patch) {
351  default_last_patch_ = patch;
352  }
353  std::string get_default_first_patch() const { return default_first_patch_; }
354  std::string get_default_last_patch() const { return default_last_patch_; }
356 };
357 
359 
361 
362 //! A CHARMM patch residue.
363 /** Patch residues are similar to regular residues, except that they are
364  used to modify an existing residue. Any atoms they contain replace or
365  add to those in the residue; they can also remove atoms.
366  Atom names are as for regular residues, except for patches that affect
367  two residues (e.g. DISU, LINK) in which each atom contains a 1 or 2
368  prefix (plus an optional colon) to identify the residue which will be
369  patched; for example, both 1:CA and 1CA identify the CA atom in the
370  first residue.
371  */
372 class IMPATOMEXPORT CHARMMPatch : public CHARMMResidueTopologyBase {
373  Vector<std::string> deleted_atoms_;
374 
375  public:
376  //! Construct a new, empty patch residue.
377  /** To get an existing patch, use CHARMMParameters::get_patch() instead.
378  */
379  CHARMMPatch(std::string type) : CHARMMResidueTopologyBase(type) {
380  set_name(std::string("CHARMM patching residue ") + type);
381  }
382 
383  void add_removed_atom(std::string name) { deleted_atoms_.push_back(name); }
384 
385  unsigned int get_number_of_removed_atoms() const {
386  return deleted_atoms_.size();
387  }
388 
389  std::string get_removed_atom(unsigned int i) const {
390  return deleted_atoms_[i];
391  }
392 
393  //! Apply the patch to the residue, modifying its topology accordingly.
394  /** \note Most CHARMM patches are designed to be applied in isolation;
395  it is usually an error to try to apply two different patches
396  to the same residue. Thus, by default \imp prohibits this.
397  To allow an already-patched residue to be re-patched, first
398  call CHARMMResidueTopology::set_patched(false).
399  */
400  void apply(CHARMMResidueTopology *res) const;
401 
402  //! Apply the patch to the given pair of residues.
403  /** This can only be used for special two-residue patches, such as
404  DISU or LINK. In a two-residue patch, each atom has a 1: or 2: prefix
405  to identify the residue it refers to.
406 
407  \throws ValueException if the patch is not a two-residue patch.
408  */
409  void apply(CHARMMResidueTopology *res1, CHARMMResidueTopology *res2) const;
411 };
412 
414 
415 //! The topology of a single residue in a model.
416 /** Each CHARMMResidueTopology object can represent an 'unknown' residue
417  (containing no atoms or bonds) or is a copy of an existing
418  CHARMMIdealResidueTopology. Optionally, patches can be applied to
419  residues (see CHARMMPatch::apply()) to add N- or C-termini, disulfide
420  bridges, sidechain modifications, etc.
421  */
422 class IMPATOMEXPORT CHARMMResidueTopology : public CHARMMIdealResidueTopology {
423  bool patched_;
424 
425  public:
426  //! Create an empty topology, containing no atoms or bonds.
428  : CHARMMIdealResidueTopology(type), patched_(false) {}
429 
430  //! Construct residue topology as a copy of an existing topology.
432  : CHARMMIdealResidueTopology(ideal->get_type()), patched_(false) {
433  set_default_first_patch(ideal->get_default_first_patch());
434  set_default_last_patch(ideal->get_default_last_patch());
435 
436  // Add atoms from existing topology
437  for (unsigned int i = 0; i < ideal->get_number_of_atoms(); ++i) {
438  add_atom(ideal->get_atom(i));
439  }
440  // Add angles/bonds/dihedrals/impropers
441  for (unsigned int i = 0; i < ideal->get_number_of_bonds(); ++i) {
442  add_bond(ideal->get_bond(i));
443  }
444  for (unsigned int i = 0; i < ideal->get_number_of_angles(); ++i) {
445  add_angle(ideal->get_angle(i));
446  }
447  for (unsigned int i = 0; i < ideal->get_number_of_dihedrals(); ++i) {
448  add_dihedral(ideal->get_dihedral(i));
449  }
450  for (unsigned int i = 0; i < ideal->get_number_of_impropers(); ++i) {
451  add_improper(ideal->get_improper(i));
452  }
453  for (unsigned int i = 0; i < ideal->get_number_of_internal_coordinates();
454  ++i) {
455  add_internal_coordinate(ideal->get_internal_coordinate(i));
456  }
457  }
458 
459  bool get_patched() const { return patched_; }
460  void set_patched(bool patched) { patched_ = patched; }
461 
463 };
464 
466 
467 IMPATOM_END_NAMESPACE
468 
469 #endif /* IMPATOM_CHARMM_TOPOLOGY_H */
CHARMMResidueTopology(CHARMMIdealResidueTopology *ideal)
Construct residue topology as a copy of an existing topology.
The ideal topology of a single residue.
#define IMP_SHOWABLE_INLINE(Name, how_to_show)
Declare the methods needed by an object that can be printed.
A connection (bond, angle, dihedral) between some number of endpoints.
A single atom in a CHARMM topology.
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Definition: object_macros.h:25
double get_dihedral(XYZ a, XYZ b, XYZ c, XYZ d)
Compute the dihedral angle (in radians) between the four particles.
Simple atom decorator.
A geometric relationship between four atoms.
The type of an atom.
Atom get_atom(Residue rd, AtomType at)
Return a particle atom from the residue.
A more IMP-like version of the std::vector.
Definition: Vector.h:39
The end of a bond, angle, dihedral, improper, or internal coordinate.
double get_angle(const Line3D &a, const Line3D &b)
Get angle in radians between two lines around their closest points.
The topology of a single residue in a model.
A smart pointer to a reference counted object.
Definition: Pointer.h:87
#define IMP_INTERNAL_CHECK(expr, message)
An assertion to check for internal errors in IMP. An IMP::ErrorException will be thrown.
Definition: check_macros.h:139
Decorator for helping deal with a hierarchy of molecules.
#define IMP_VALUES(Name, PluralName)
Define the type for storing sets of values.
Definition: value_macros.h:23
Common base class for heavy weight IMP objects.
Definition: Object.h:106
Base class for all CHARMM residue-based topology.
A decorator for a particle representing an atom.
Definition: atom/Atom.h:234
The type for a residue.
CHARMMPatch(std::string type)
Construct a new, empty patch residue.
bool get_contains_atom(std::string name) const
Return true if the bond contains the named atom.
A CHARMM patch residue.
std::ostream & show(Hierarchy h, std::ostream &out=std::cout)
Print the hierarchy using a given decorator to display each node.
#define IMP_OBJECTS(Name, PluralName)
Define the types for storing lists of object pointers.
Definition: object_macros.h:44
A shared base class to help in debugging and things.
Bond get_bond(Bonded a, Bonded b)
Get the bond between two particles.
CHARMMResidueTopology(ResidueType type)
Create an empty topology, containing no atoms or bonds.
Atoms get_atoms(const CHARMMResidueTopology *current_residue, const CHARMMResidueTopology *previous_residue, const CHARMMResidueTopology *next_residue, const std::map< const CHARMMResidueTopology *, Hierarchy > &resmap) const
Map the bond to a list of Atom particles.