RDKit
Open-source cheminformatics and machine learning.
MolTransforms.h
Go to the documentation of this file.
1 //
2 // Copyright (C) 2003-2006 Rational Discovery LLC
3 //
4 // @@ All Rights Reserved @@
5 // This file is part of the RDKit.
6 // The contents are covered by the terms of the BSD license
7 // which is included in the file license.txt, found at the root
8 // of the RDKit source tree.
9 //
10 #include <RDGeneral/export.h>
11 #ifndef _RD_MOLTRANSFORMS_H_
12 #define _RD_MOLTRANSFORMS_H_
13 
14 #include <Geometry/point.h>
15 #include <Numerics/SymmMatrix.h>
16 
17 #ifdef RDK_HAS_EIGEN3
18 #include <Eigen/Dense>
19 #endif
20 
21 namespace RDKit {
22 class ROMol;
23 class Atom;
24 class Conformer;
25 } // namespace RDKit
26 
27 namespace RDGeom {
28 class Transform3D;
29 }
30 
31 namespace MolTransforms {
33  RDGeom::Transform3D &tform);
35  RDGeom::Transform3D &tform);
36 
37 //! Compute the centroid of a conformer
38 /*!
39  This is simple the average of the heavy atom locations in the conformer,
40  not attention is paid to hydrogens or the differences in atomic radii
41 
42  \param conf Conformer of interest
43  \param ignoreHs If true, ignore hydrogen atoms
44 */
46  const RDKit::Conformer &conf, bool ignoreHs = true);
47 
48 #ifdef RDK_HAS_EIGEN3
49 //! Compute principal axes and moments of inertia for a conformer
50 /*!
51  These values are calculated from the inertia tensor:
52  Iij = - sum_{s=1..N}(w_s * r_{si} * r_{sj}) i != j
53  Iii = sum_{s=1..N} sum_{j!=i} (w_s * r_{sj} * r_{sj})
54  where the coordinates are relative to the center of mass.
55 
56 
57  \param conf Conformer of interest
58  \param axes used to return the principal axes
59  \param moments used to return the principal moments
60  \param ignoreHs If true, ignore hydrogen atoms
61  \param force If true, the calculation will be carried out even if a
62  cached value is present
63  \param weights If present used to weight the atomic coordinates
64 
65  \returns whether or not the calculation was successful
66 */
67 RDKIT_MOLTRANSFORMS_EXPORT bool computePrincipalAxesAndMoments(
68  const RDKit::Conformer &conf, Eigen::Matrix3d &axes,
69  Eigen::Vector3d &moments, bool ignoreHs = false, bool force = false,
70  const std::vector<double> *weights = NULL);
71 //! Compute principal axes and moments from the gyration matrix of a conformer
72 /*!
73 
74  These values are calculated from the gyration matrix/tensor:
75  Iij = sum_{s=1..N}(w_s * r_{si} * r_{sj}) i != j
76  Iii = sum_{s=1..N} sum_{t!=s}(w_s * r_{si} * r_{ti})
77  where the coordinates are relative to the center of mass.
78 
79  \param conf Conformer of interest
80  \param axes used to return the principal axes
81  \param moments used to return the principal moments
82  \param ignoreHs If true, ignore hydrogen atoms
83  \param force If true, the calculation will be carried out even if a
84  cached value is present
85  \param weights If present used to weight the atomic coordinates
86 
87  \returns whether or not the calculation was successful
88 */
90 computePrincipalAxesAndMomentsFromGyrationMatrix(
91  const RDKit::Conformer &conf, Eigen::Matrix3d &axes,
92  Eigen::Vector3d &moments, bool ignoreHs = false, bool force = false,
93  const std::vector<double> *weights = NULL);
94 #endif
95 
96 //! Compute the transformation require to orient the conformation
97 //! along the principal axes about the center; i.e. center is made to coincide
98 // with the
99 //! origin, the largest princiapl axis with the x-axis, the next largest with
100 // the y-axis
101 //! and the smallest with the z-axis
102 /*!
103  If center is not specified the centroid of the conformer will be used
104  \param conf Conformer of interest
105  \param center Center to be used for canonicalization, defaults to
106  the centroid of the
107  conformation
108  \param normalizeCovar Normalize the covariance matrix with the number of
109  atoms
110  \param ignoreHs Optionally ignore hydrogens
111 */
113  const RDKit::Conformer &conf, const RDGeom::Point3D *center = 0,
114  bool normalizeCovar = false, bool ignoreHs = true);
115 
116 //! Transform the conformation using the specified transformation
118  RDKit::Conformer &conf, const RDGeom::Transform3D &trans);
119 
120 //! Canonicalize the orientation of a conformer so that its principal axes
121 //! around the specified center point coincide with the x, y, z axes
122 /*!
123  \param conf The conformer of interest
124  \param center Optional center point about which the principal axes are
125  computed
126  if not specified the centroid of the conformer will be
127  used
128  \param normalizeCovar Optionally normalize the covariance matrix by the number
129  of atoms
130  \param ignoreHs If true, ignore hydrogen atoms
131 
132 */
134  RDKit::Conformer &conf, const RDGeom::Point3D *center = 0,
135  bool normalizeCovar = false, bool ignoreHs = true);
136 
137 //! Canonicalize all the conformations in a molecule
138 /*!
139  \param mol the molecule of interest
140  \param normalizeCovar Optionally normalize the covariance matrix by the number
141  of atoms
142  \param ignoreHs If true, ignore hydrogens
143 */
145  bool normalizeCovar = false,
146  bool ignoreHs = true);
147 
148 //! Get the bond length between the specified atoms i, j
150  unsigned int iAtomId,
151  unsigned int jAtomId);
152 
153 //! Set the bond length between the specified atoms i, j
154 //! (all atoms bonded to atom j are moved)
156  unsigned int iAtomId,
157  unsigned int jAtomId,
158  double value);
159 
160 //! Get the angle in radians among the specified atoms i, j, k
162  unsigned int iAtomId,
163  unsigned int jAtomId,
164  unsigned int kAtomId);
165 
166 //! Get the angle in degrees among the specified atoms i, j, k
167 inline double getAngleDeg(const RDKit::Conformer &conf, unsigned int iAtomId,
168  unsigned int jAtomId, unsigned int kAtomId) {
169  return (180. / M_PI * getAngleRad(conf, iAtomId, jAtomId, kAtomId));
170 }
171 
172 //! Set the angle in radians among the specified atoms i, j, k
173 //! (all atoms bonded to atom k are moved)
175  unsigned int iAtomId,
176  unsigned int jAtomId,
177  unsigned int kAtomId, double value);
178 
179 //! Set the angle in degrees among the specified atoms i, j, k
180 //! (all atoms bonded to atom k are moved)
181 inline void setAngleDeg(RDKit::Conformer &conf, unsigned int iAtomId,
182  unsigned int jAtomId, unsigned int kAtomId,
183  double value) {
184  setAngleRad(conf, iAtomId, jAtomId, kAtomId, value / 180. * M_PI);
185 }
186 
187 //! Get the dihedral angle in radians among the specified atoms i, j, k, l
189  unsigned int iAtomId,
190  unsigned int jAtomId,
191  unsigned int kAtomId,
192  unsigned int lAtomId);
193 
194 //! Get the dihedral angle in degrees among the specified atoms i, j, k, l
195 inline double getDihedralDeg(const RDKit::Conformer &conf, unsigned int iAtomId,
196  unsigned int jAtomId, unsigned int kAtomId,
197  unsigned int lAtomId) {
198  return (180. / M_PI *
199  getDihedralRad(conf, iAtomId, jAtomId, kAtomId, lAtomId));
200 }
201 
202 //! Set the dihedral angle in radians among the specified atoms i, j, k, l
203 //! (all atoms bonded to atom l are moved)
205  RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId,
206  unsigned int kAtomId, unsigned int lAtomId, double value);
207 
208 //! Set the dihedral angle in degrees among the specified atoms i, j, k, l
209 //! (all atoms bonded to atom l are moved)
210 inline void setDihedralDeg(RDKit::Conformer &conf, unsigned int iAtomId,
211  unsigned int jAtomId, unsigned int kAtomId,
212  unsigned int lAtomId, double value) {
213  setDihedralRad(conf, iAtomId, jAtomId, kAtomId, lAtomId, value / 180. * M_PI);
214 }
215 } // namespace MolTransforms
216 #endif
point.h
MolTransforms::getAngleDeg
double getAngleDeg(const RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId)
Get the angle in degrees among the specified atoms i, j, k.
Definition: MolTransforms.h:167
RDGeom
Definition: TorsionAngleM6.h:20
MolTransforms::setAngleRad
RDKIT_MOLTRANSFORMS_EXPORT void setAngleRad(RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, double value)
RDGeom::Point3D
Definition: point.h:46
MolTransforms::transformMolsAtoms
RDKIT_MOLTRANSFORMS_EXPORT void transformMolsAtoms(RDKit::ROMol *mol, RDGeom::Transform3D &tform)
MolTransforms::getAngleRad
RDKIT_MOLTRANSFORMS_EXPORT double getAngleRad(const RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId)
Get the angle in radians among the specified atoms i, j, k.
MolTransforms
Definition: MolTransforms.h:31
RDKit::Atom
The class for representing atoms.
Definition: Atom.h:69
MolTransforms::computeCentroid
RDKIT_MOLTRANSFORMS_EXPORT RDGeom::Point3D computeCentroid(const RDKit::Conformer &conf, bool ignoreHs=true)
Compute the centroid of a conformer.
M_PI
#define M_PI
Definition: MMFF/Params.h:27
RDKit::ROMol
Definition: ROMol.h:171
MolTransforms::setAngleDeg
void setAngleDeg(RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, double value)
Definition: MolTransforms.h:181
MolTransforms::setBondLength
RDKIT_MOLTRANSFORMS_EXPORT void setBondLength(RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, double value)
RDKit::Conformer
The class for representing 2D or 3D conformation of a molecule.
Definition: Conformer.h:43
MolTransforms::getDihedralDeg
double getDihedralDeg(const RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, unsigned int lAtomId)
Get the dihedral angle in degrees among the specified atoms i, j, k, l.
Definition: MolTransforms.h:195
MolTransforms::transformConformer
RDKIT_MOLTRANSFORMS_EXPORT void transformConformer(RDKit::Conformer &conf, const RDGeom::Transform3D &trans)
Transform the conformation using the specified transformation.
MolTransforms::canonicalizeConformer
RDKIT_MOLTRANSFORMS_EXPORT void canonicalizeConformer(RDKit::Conformer &conf, const RDGeom::Point3D *center=0, bool normalizeCovar=false, bool ignoreHs=true)
RDKIT_MOLTRANSFORMS_EXPORT
#define RDKIT_MOLTRANSFORMS_EXPORT
Definition: export.h:450
MolTransforms::getDihedralRad
RDKIT_MOLTRANSFORMS_EXPORT double getDihedralRad(const RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, unsigned int lAtomId)
Get the dihedral angle in radians among the specified atoms i, j, k, l.
RDKit
Std stuff.
Definition: Atom.h:30
MolTransforms::setDihedralRad
RDKIT_MOLTRANSFORMS_EXPORT void setDihedralRad(RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, unsigned int lAtomId, double value)
RDGeom::Transform3D
Definition: Transform3D.h:22
MolTransforms::getBondLength
RDKIT_MOLTRANSFORMS_EXPORT double getBondLength(const RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId)
Get the bond length between the specified atoms i, j.
MolTransforms::computeCanonicalTransform
RDKIT_MOLTRANSFORMS_EXPORT RDGeom::Transform3D * computeCanonicalTransform(const RDKit::Conformer &conf, const RDGeom::Point3D *center=0, bool normalizeCovar=false, bool ignoreHs=true)
origin, the largest princiapl axis with the x-axis, the next largest with
MolTransforms::setDihedralDeg
void setDihedralDeg(RDKit::Conformer &conf, unsigned int iAtomId, unsigned int jAtomId, unsigned int kAtomId, unsigned int lAtomId, double value)
Definition: MolTransforms.h:210
MolTransforms::canonicalizeMol
RDKIT_MOLTRANSFORMS_EXPORT void canonicalizeMol(RDKit::ROMol &mol, bool normalizeCovar=false, bool ignoreHs=true)
Canonicalize all the conformations in a molecule.
MolTransforms::transformAtom
RDKIT_MOLTRANSFORMS_EXPORT void transformAtom(RDKit::Atom *atom, RDGeom::Transform3D &tform)
SymmMatrix.h
export.h