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