Loading...
Searching...
No Matches
ShapeInput.h
Go to the documentation of this file.
//
// Copyright (C) 2026 David Cosgrove and other RDKit contributors
//
// @@ All Rights Reserved @@
// This file is part of the RDKit.
// The contents are covered by the terms of the BSD license
// which is included in the file license.txt, found at the root
// of the RDKit source tree.
//
// Original author: David Cosgrove (CozChemIx Limited)
//
#ifndef RDKIT_SHAPEINPUT_GUARD
#define RDKIT_SHAPEINPUT_GUARD
#include <array>
#include <vector>
#include <RDGeneral/export.h>
#include <Geometry/Transform3D.h>
#include <RDGeneral/BoostStartInclude.h>
#include <boost/dynamic_bitset.hpp>
#ifdef RDK_USE_BOOST_SERIALIZATION
#include <boost/archive/text_oarchive.hpp>
#include <boost/archive/text_iarchive.hpp>
#include <boost/serialization/vector.hpp>
#include <boost/serialization/array.hpp>
#include <boost/serialization/unique_ptr.hpp>
#endif
#include <RDGeneral/BoostEndInclude.h>
#include <GraphMol/GaussianShape/ShapeOverlayOptions.h>
// The code below was provided by Claude (Sonnet 4.6).
// If first tried to get me to use boost/serialization/dynamic_bitset.hpp
// and then admitted that it had made that up.
namespace boost {
namespace serialization {
template <class Archive, typename Block, typename Allocator>
void serialize(Archive &ar, boost::dynamic_bitset<Block, Allocator> &bs,
const unsigned int /*version*/) {
size_t num_bits = bs.size();
ar & num_bits;
std::vector<Block> blocks;
if (Archive::is_saving::value) {
to_block_range(bs, std::back_inserter(blocks));
}
ar & blocks;
if (Archive::is_loading::value) {
bs.resize(num_bits);
from_block_range(blocks.begin(), blocks.end(), bs);
bs.resize(num_bits); // trim any excess bits
}
}
} // namespace serialization
} // namespace boost
namespace RDKit {
class ROMol;
class RWMol;
namespace GaussianShape {
// From Grant et al.
constexpr double P = 2.7;
constexpr double KAPPA = 2.41798793102;
using CustomFeatures =
std::vector<std::tuple<unsigned int, RDGeom::Point3D, double>>;
struct ShapeInputOptions {
ShapeInputOptions() = default;
ShapeInputOptions(const ShapeInputOptions &) = default;
ShapeInputOptions(ShapeInputOptions &&) = default;
ShapeInputOptions &operator=(const ShapeInputOptions &) = default;
ShapeInputOptions &operator=(ShapeInputOptions &&) = default;
~ShapeInputOptions() = default;
// By default, it will create features using the RDKit pharmacophore
// definitions.
bool useColors{
true}; //! Whether to build the color features. By default, it will
//! create features using the RDKit pharmacophore definitions.
CustomFeatures customFeatures; //! Custom color features used verbatim. A
//! vector of tuples of integer type, Point3D
//! coords, double radius.
std::vector<unsigned int>
atomSubset; //! If not empty, use just these atoms in the molecule to
//! form the ShapeInput object.
std::vector<std::pair<unsigned int, double>>
atomRadii; //! Use these non-standard radii for these atoms. The int is
//! for the atom index in the molecule, not the atomic number.
//! Not all atoms need be specified, just some radii can be
//! over-ridden, with the rest left as standard.
bool allCarbonRadii{
true}; //! Whether to use carbon radii for all atoms (which is quicker
//! but less accurate) or vdw radii appropriate for the elements.
};
// Data for shape alignment code
class RDKIT_GAUSSIANSHAPE_EXPORT ShapeInput {
public:
//! Create the ShapeInput object.
//! @param mol: The molecule of interest
//! @param confId: The conformer to use
//! @param opts: Options for setting up the shape
ShapeInput(const ROMol &mol, int confId = -1,
const ShapeInputOptions &opts = ShapeInputOptions(),
const ShapeOverlayOptions &overlayOpts = ShapeOverlayOptions());
ShapeInput(const std::string &str) {
#ifndef RDK_USE_BOOST_SERIALIZATION
PRECONDITION(0, "Boost SERIALIZATION is not enabled")
#else
std::stringstream ss(str);
boost::archive::text_iarchive ia(ss);
ia &*this;
#endif
}
ShapeInput(const ShapeInput &other);
ShapeInput(ShapeInput &&other) = default;
ShapeInput &operator=(const ShapeInput &other);
ShapeInput &operator=(ShapeInput &&other) = default;
virtual ~ShapeInput() = default;
std::string toString() const {
#ifndef RDK_USE_BOOST_SERIALIZATION
PRECONDITION(0, "Boost SERIALIZATION is not enabled")
#else
std::stringstream ss;
boost::archive::text_oarchive oa(ss);
oa &*this;
return ss.str();
#endif
}
// Note that the coords returned is a vector size 4*getNumAtoms()
// with the 4th value per atom being the alpha paramter.
const std::vector<double> &getCoords() const { return d_coords; }
//! Fetch the coordinates of the atoms and optionally features.
std::vector<RDGeom::Point3D> getAtomPoints(bool includeColors = false) const;
bool getNormalized() const { return d_normalized; }
const std::vector<int> &getTypes() const { return d_types; }
unsigned int getNumAtoms() const { return d_numAtoms; }
unsigned int getNumFeatures() const { return d_numFeats; }
double getShapeVolume() const { return d_selfOverlapVol; }
double getColorVolume() const { return d_selfOverlapColor; }
const boost::dynamic_bitset<> *getCarbonRadii() const {
return d_carbonRadii.get();
}
// These functions use cached values if available.
const std::array<double, 9> &calcCanonicalRotation();
const std::array<double, 3> &calcCanonicalTranslation();
const std::array<double, 3> &calcEigenValues();
const std::array<size_t, 6> &calcExtremes();
// Return the principal moments of inertia, if Eigen3 is available, and the
// eigenvalues of the canonical transformation if not.
std::array<double, 3> calcMomentsOfInertia(bool includeColors = false) const;
// Align the principal axes to the cartesian axes and centre on the origin.
// Doesn't require that the shape was created from a molecule. Creates
// the necessary transformation if not already done.
void normalizeCoords();
void transformCoords(RDGeom::Transform3D &xform);
// Mock a molecule up from the shape for visual inspection and sometimes
// calculation of the normalization matrices. No bonds.
// Atoms are C, features are N.
virtual std::unique_ptr<RWMol> shapeToMol(bool includeColors = true) const;
#ifdef RDK_USE_BOOST_SERIALIZATION
template <class Archive>
void serialize(Archive &ar, const unsigned int) {
ar & d_coords;
ar & d_types;
ar & d_numAtoms;
ar & d_numFeats;
ar & d_selfOverlapVol;
ar & d_selfOverlapColor;
ar & d_extremePoints;
ar & d_carbonRadii;
ar & d_normalized;
ar & d_normalizationOK;
ar & d_canonRot;
ar & d_canonTrans;
ar & d_eigenValues;
}
#endif
private:
void extractAtoms(const ROMol &mol, int confId,
const ShapeInputOptions &opts);
// Extract the features for the color scores, using RDKit pphore features
// for now. Other options to be added later.
void extractFeatures(const ROMol &mol, int confId,
const ShapeInputOptions &shapeOpts);
// Calculate the rotation and translation that will align the principal axes
// to the cartesian axes and centre on the origin.
void calcNormalization();
void calculateExtremes();
std::vector<double> d_coords; // The coordinates and alpha parameter for the
// atoms and features, packed as 4 floats per
// item - x, y, z and alpha. alpha is KAPPA / (r * r) where r is the radius
// of the atom. This is not used if using all_atoms_carbon mode.
std::vector<int> d_types; // The feature types. The size is the same
// as the number of coordinates, padded with 0
// for the atoms.
unsigned int d_numAtoms; // The number of atoms
unsigned int d_numFeats; // The number of features
double d_selfOverlapVol{0.0}; // Shape volume
double d_selfOverlapColor{0.0}; // Color volume
// These are the points at the extremes of the x, y and z axes.
// they are min_x, min_y, min_z and max_x, max_y, max_z.
std::array<size_t, 6> d_extremePoints;
std::unique_ptr<boost::dynamic_bitset<>>
d_carbonRadii; // Flags those atoms with a carbon radius, for faster
// calculation later.
// This is the rotation and translation to align the principal axes of the
// shape with cartesian axes. If d_normalized is true, it has been applied
// to the coordinates.
bool d_normalized{false};
// If the shape is moved, the normalization matrices are no longer valid.
// This flags that so it is re-computed as required.
bool d_normalizationOK{false};
std::array<double, 9> d_canonRot;
std::array<double, 3> d_canonTrans;
// The sorted eigenvalues of the principal axes.
std::array<double, 3> d_eigenValues;
};
// Calculate the mean position of the given atoms.
RDKIT_GAUSSIANSHAPE_EXPORT RDGeom::Point3D computeFeaturePos(
const ROMol &mol, int confId, const std::vector<unsigned int> &ats);
RDKIT_GAUSSIANSHAPE_EXPORT RDGeom::Transform3D quatTransToTransform(
const double *quat, const double *trans);
// Apply the transformation to the coordinates assumed to be in
// ShapeInput.d_coords form.
RDKIT_GAUSSIANSHAPE_EXPORT void applyTransformToShape(
std::vector<double> &shape, RDGeom::Transform3D &xform);
RDKIT_GAUSSIANSHAPE_EXPORT void applyTransformToShape(
const double *inShape, double *outShape, size_t numPoints,
RDGeom::Transform3D &xform);
RDKIT_GAUSSIANSHAPE_EXPORT void translateShape(
std::vector<double> &shape, const RDGeom::Point3D &translation);
RDKIT_GAUSSIANSHAPE_EXPORT void translateShape(
const double *inShape, double *outShape, size_t numPoints,
const RDGeom::Point3D &translation);
} // namespace GaussianShape
} // namespace RDKit
#endif // RDKIT_SHAPEINPUT_GUARD
BoostEndInclude.h
BoostStartInclude.h
PRECONDITION
#define PRECONDITION(expr, mess)
Definition Invariant.h:108
ShapeOverlayOptions.h
Transform3D.h
RDGeom::Point3D
Definition point.h:53
RDGeom::Transform3D
Definition Transform3D.h:23
RDKit::GaussianShape::ShapeInput::ShapeInput
ShapeInput(const std::string &str)
Definition ShapeInput.h:117
RDKit::GaussianShape::ShapeInput::ShapeInput
ShapeInput(const ROMol &mol, int confId=-1, const ShapeInputOptions &opts=ShapeInputOptions(), const ShapeOverlayOptions &overlayOpts=ShapeOverlayOptions())
RDKit::GaussianShape::ShapeInput::calcExtremes
const std::array< size_t, 6 > & calcExtremes()
RDKit::GaussianShape::ShapeInput::calcEigenValues
const std::array< double, 3 > & calcEigenValues()
RDKit::GaussianShape::ShapeInput::getAtomPoints
std::vector< RDGeom::Point3D > getAtomPoints(bool includeColors=false) const
Fetch the coordinates of the atoms and optionally features.
RDKit::GaussianShape::ShapeInput::operator=
ShapeInput & operator=(const ShapeInput &other)
RDKit::GaussianShape::ShapeInput::getNumAtoms
unsigned int getNumAtoms() const
Definition ShapeInput.h:150
RDKit::GaussianShape::ShapeInput::ShapeInput
ShapeInput(const ShapeInput &other)
RDKit::GaussianShape::ShapeInput::normalizeCoords
void normalizeCoords()
RDKit::GaussianShape::ShapeInput::toString
std::string toString() const
Definition ShapeInput.h:132
RDKit::GaussianShape::ShapeInput::getColorVolume
double getColorVolume() const
Definition ShapeInput.h:153
RDKit::GaussianShape::ShapeInput::operator=
ShapeInput & operator=(ShapeInput &&other)=default
RDKit::GaussianShape::ShapeInput::ShapeInput
ShapeInput(ShapeInput &&other)=default
RDKit::GaussianShape::ShapeInput::getTypes
const std::vector< int > & getTypes() const
Definition ShapeInput.h:149
RDKit::GaussianShape::ShapeInput::~ShapeInput
virtual ~ShapeInput()=default
RDKit::GaussianShape::ShapeInput::getNumFeatures
unsigned int getNumFeatures() const
Definition ShapeInput.h:151
RDKit::GaussianShape::ShapeInput::calcCanonicalTranslation
const std::array< double, 3 > & calcCanonicalTranslation()
RDKit::GaussianShape::ShapeInput::getCoords
const std::vector< double > & getCoords() const
Definition ShapeInput.h:145
RDKit::GaussianShape::ShapeInput::transformCoords
void transformCoords(RDGeom::Transform3D &xform)
RDKit::GaussianShape::ShapeInput::calcCanonicalRotation
const std::array< double, 9 > & calcCanonicalRotation()
RDKit::GaussianShape::ShapeInput::getShapeVolume
double getShapeVolume() const
Definition ShapeInput.h:152
RDKit::GaussianShape::ShapeInput::shapeToMol
virtual std::unique_ptr< RWMol > shapeToMol(bool includeColors=true) const
RDKit::GaussianShape::ShapeInput::getCarbonRadii
const boost::dynamic_bitset * getCarbonRadii() const
Definition ShapeInput.h:154
RDKit::GaussianShape::ShapeInput::calcMomentsOfInertia
std::array< double, 3 > calcMomentsOfInertia(bool includeColors=false) const
RDKit::GaussianShape::ShapeInput::getNormalized
bool getNormalized() const
Definition ShapeInput.h:148
RDKit::ROMol
Definition ROMol.h:285
RDKit::RWMol
RWMol is a molecule class that is intended to be edited.
Definition RWMol.h:32
export.h
RDKIT_GAUSSIANSHAPE_EXPORT
#define RDKIT_GAUSSIANSHAPE_EXPORT
Definition export.h:233
RDKit::GaussianShape
Definition GaussianShape.h:28
RDKit::GaussianShape::P
constexpr double P
Definition ShapeInput.h:71
RDKit::GaussianShape::CustomFeatures
std::vector< std::tuple< unsigned int, RDGeom::Point3D, double > > CustomFeatures
Definition ShapeInput.h:73
RDKit::GaussianShape::translateShape
RDKIT_GAUSSIANSHAPE_EXPORT void translateShape(std::vector< double > &shape, const RDGeom::Point3D &translation)
RDKit::GaussianShape::quatTransToTransform
RDKIT_GAUSSIANSHAPE_EXPORT RDGeom::Transform3D quatTransToTransform(const double *quat, const double *trans)
RDKit::GaussianShape::KAPPA
constexpr double KAPPA
Definition ShapeInput.h:72
RDKit::GaussianShape::applyTransformToShape
RDKIT_GAUSSIANSHAPE_EXPORT void applyTransformToShape(std::vector< double > &shape, RDGeom::Transform3D &xform)
RDKit::GaussianShape::computeFeaturePos
RDKIT_GAUSSIANSHAPE_EXPORT RDGeom::Point3D computeFeaturePos(const ROMol &mol, int confId, const std::vector< unsigned int > &ats)
RDKit
Std stuff.
Definition Abbreviations.h:19
boost::serialization
Definition ShapeInput.h:39
boost::serialization::serialize
void serialize(Archive &ar, boost::dynamic_bitset< Block, Allocator > &bs, const unsigned int)
Definition ShapeInput.h:42
boost
Definition ShapeInput.h:38
RDKit::GaussianShape::ShapeInputOptions
Definition ShapeInput.h:76
RDKit::GaussianShape::ShapeInputOptions::operator=
ShapeInputOptions & operator=(const ShapeInputOptions &)=default
RDKit::GaussianShape::ShapeInputOptions::useColors
bool useColors
Definition ShapeInput.h:87
RDKit::GaussianShape::ShapeInputOptions::customFeatures
CustomFeatures customFeatures
Definition ShapeInput.h:91
RDKit::GaussianShape::ShapeInputOptions::atomSubset
std::vector< unsigned int > atomSubset
Definition ShapeInput.h:95
RDKit::GaussianShape::ShapeInputOptions::atomRadii
std::vector< std::pair< unsigned int, double > > atomRadii
Definition ShapeInput.h:98
RDKit::GaussianShape::ShapeInputOptions::~ShapeInputOptions
~ShapeInputOptions()=default
RDKit::GaussianShape::ShapeInputOptions::allCarbonRadii
bool allCarbonRadii
Definition ShapeInput.h:102
RDKit::GaussianShape::ShapeInputOptions::ShapeInputOptions
ShapeInputOptions(ShapeInputOptions &&)=default
RDKit::GaussianShape::ShapeInputOptions::operator=
ShapeInputOptions & operator=(ShapeInputOptions &&)=default
RDKit::GaussianShape::ShapeInputOptions::ShapeInputOptions
ShapeInputOptions()=default
RDKit::GaussianShape::ShapeInputOptions::ShapeInputOptions
ShapeInputOptions(const ShapeInputOptions &)=default
RDKit::GaussianShape::ShapeOverlayOptions
Definition ShapeOverlayOptions.h:50