diff --git a/Code/GraphMol/Descriptors/CMakeLists.txt b/Code/GraphMol/Descriptors/CMakeLists.txt index 7c737dc8c..24540e455 100644 --- a/Code/GraphMol/Descriptors/CMakeLists.txt +++ b/Code/GraphMol/Descriptors/CMakeLists.txt @@ -1,7 +1,7 @@ if(RDK_BUILD_DESCRIPTORS3D) - set(DESC3D_HDRS MolDescriptors3D.h EEM.h PBF.h PMI.h AUTOCORR3D.h RDF.h MORSE.h GETAWAY.h WHIM.h) - set(DESC3D_SOURCES EEM.cpp PBF.cpp PMI.cpp AUTOCORR3D.cpp RDF.cpp MORSE.cpp GETAWAY.cpp WHIM.cpp) + set(DESC3D_HDRS MolDescriptors3D.h EEM.h PBF.h PMI.h AUTOCORR3D.h RDF.h MORSE.h GETAWAY.h WHIM.h CoulombMat.h) + set(DESC3D_SOURCES EEM.cpp PBF.cpp PMI.cpp AUTOCORR3D.cpp RDF.cpp MORSE.cpp GETAWAY.cpp WHIM.cpp CoulombMat.cpp) endif(RDK_BUILD_DESCRIPTORS3D) @@ -62,6 +62,8 @@ LINK_LIBRARIES Descriptors FileParsers SmilesParse MolTransforms PartialCharges LINK_LIBRARIES Descriptors FileParsers SmilesParse MolTransforms PartialCharges GraphMol DataStructs EigenSolvers RDGeneral RDGeometryLib ${RDKit_THREAD_LIBS} ) rdkit_test(testAUTOCORR3D testAUTOCORR3D.cpp LINK_LIBRARIES Descriptors FileParsers SmilesParse MolTransforms PartialCharges GraphMol DataStructs EigenSolvers RDGeneral RDGeometryLib ${RDKit_THREAD_LIBS} ) + rdkit_test(testCoulombMat testCoulombMat.cpp +LINK_LIBRARIES Descriptors FileParsers SmilesParse MolTransforms PartialCharges GraphMol DataStructs EigenSolvers RDGeneral RDGeometryLib ForceField ForceFieldHelpers Optimizer DistGeomHelpers ${RDKit_THREAD_LIBS} ) endif(RDK_BUILD_DESCRIPTORS3D) diff --git a/Code/GraphMol/Descriptors/CoulombMat.cpp b/Code/GraphMol/Descriptors/CoulombMat.cpp new file mode 100644 index 000000000..8c9302f5c --- /dev/null +++ b/Code/GraphMol/Descriptors/CoulombMat.cpp @@ -0,0 +1,61 @@ +// +// Copyright (c) 2018, Guillaume GODIN +// All rights reserved. +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are +// met: +// +// * Redistributions of source code must retain the above copyright +// notice, this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above +// copyright notice, this list of conditions and the following +// disclaimer in the documentation and/or other materials provided +// with the distribution. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +// OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +// +// Modified by Greg Landrum, March 2020 +// + +#include +#include +#include "CoulombMat.h" + +namespace RDKit { +namespace Descriptors { + +void CoulombMat(const ROMol &mol, std::vector> &res, int confId) { + PRECONDITION(mol.getNumConformers() >= 1, "molecule has no conformers"); + + unsigned int numAtoms = mol.getNumAtoms(); + + const auto conf = mol.getConformer(confId); + + res.resize(numAtoms); + for(unsigned int i=0;igetAtomicNum(); + res[i][i] = 0.5 * pow(Zi,2.4); + + const auto Pi = conf.getAtomPos(i); + for(unsigned int j=0;jgetAtomicNum(); + res[i][j] = res[j][i] = Zi*Zj / (Pi-Pj).length(); + } + } +} +} // namespace Descriptors +} // namespace RDKit \ No newline at end of file diff --git a/Code/GraphMol/Descriptors/CoulombMat.h b/Code/GraphMol/Descriptors/CoulombMat.h new file mode 100644 index 000000000..d6e8feb68 --- /dev/null +++ b/Code/GraphMol/Descriptors/CoulombMat.h @@ -0,0 +1,42 @@ +// +// Copyright (c) 2018, Guillaume GODIN +// All rights reserved. +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are +// met: +// +// * Redistributions of source code must retain the above copyright +// notice, this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above +// copyright notice, this list of conditions and the following +// disclaimer in the documentation and/or other materials provided +// with the distribution. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +// OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + + +#ifndef CoulombMatRDKIT_H_MAY2018 +#define CoulombMatRDKIT_H_MAY2018 + +#ifdef RDK_BUILD_DESCRIPTORS3D +namespace RDKit { +class ROMol; +namespace Descriptors { +const std::string CoulombMatVersion = "1.0.0"; +RDKIT_DESCRIPTORS_EXPORT void CoulombMat(const ROMol &mol, std::vector> &res, + int confId=-1); +} +} +#endif +#endif diff --git a/Code/GraphMol/Descriptors/MolDescriptors3D.h b/Code/GraphMol/Descriptors/MolDescriptors3D.h index 4640768e6..b2673991c 100644 --- a/Code/GraphMol/Descriptors/MolDescriptors3D.h +++ b/Code/GraphMol/Descriptors/MolDescriptors3D.h @@ -12,6 +12,7 @@ #ifndef RD_MOLDESCRIPTORS3D_H_SEPT2016 #define RD_MOLDESCRIPTORS3D_H_SEPT2016 +#include #include #include #include diff --git a/Code/GraphMol/Descriptors/Wrap/rdMolDescriptors.cpp b/Code/GraphMol/Descriptors/Wrap/rdMolDescriptors.cpp index ce471fc2c..0173e4782 100644 --- a/Code/GraphMol/Descriptors/Wrap/rdMolDescriptors.cpp +++ b/Code/GraphMol/Descriptors/Wrap/rdMolDescriptors.cpp @@ -118,6 +118,16 @@ python::tuple calcCrippenDescriptors(const RDKit::ROMol &mol, #ifdef RDK_BUILD_DESCRIPTORS3D +python::tuple calcCoulombMat(const RDKit::ROMol &mol, int confId) { + std::vector> results; + RDKit::Descriptors::CoulombMat(mol, results, confId); + python::list result; + for (auto &res : results) { + result.append(res); + } + return python::tuple(result); +} + python::list calcEEMcharges(RDKit::ROMol &mol, int confId) { std::vector res; RDKit::Descriptors::EEM(mol, res, confId); @@ -332,7 +342,8 @@ double hkAlphaHelper(const RDKit::ROMol &mol, python::object atomContribs) { RDKit::SparseIntVect *MorganFingerprintHelper( const RDKit::ROMol &mol, int radius, int nBits, python::object invariants, python::object fromAtoms, bool useChirality, bool useBondTypes, - bool useFeatures, bool useCounts, python::object bitInfo, bool includeRedundantEnvironments) { + bool useFeatures, bool useCounts, python::object bitInfo, + bool includeRedundantEnvironments) { std::vector *invars = nullptr; if (invariants) { unsigned int nInvar = @@ -371,7 +382,8 @@ RDKit::SparseIntVect *MorganFingerprintHelper( if (nBits < 0) { res = RDKit::MorganFingerprints::getFingerprint( mol, static_cast(radius), invars, froms, useChirality, - useBondTypes, useCounts, false, bitInfoMap, includeRedundantEnvironments); + useBondTypes, useCounts, false, bitInfoMap, + includeRedundantEnvironments); } else { res = RDKit::MorganFingerprints::getHashedFingerprint( mol, static_cast(radius), @@ -402,22 +414,20 @@ RDKit::SparseIntVect *MorganFingerprintHelper( return res; } -#ifdef RDK_HAS_EIGEN3 -std::pair BCUT2D_list(const RDKit::ROMol &m, python::list atomprops) -{ +#ifdef RDK_HAS_EIGEN3 +std::pair BCUT2D_list(const RDKit::ROMol &m, + python::list atomprops) { std::vector dvec; - for (int i = 0; i < len(atomprops); ++i) - { + for (int i = 0; i < len(atomprops); ++i) { dvec.push_back(boost::python::extract(atomprops[i])); } return RDKit::Descriptors::BCUT2D(m, dvec); } -std::pair BCUT2D_tuple(const RDKit::ROMol &m, python::tuple atomprops) -{ +std::pair BCUT2D_tuple(const RDKit::ROMol &m, + python::tuple atomprops) { std::vector dvec; - for (int i = 0; i < len(atomprops); ++i) - { + for (int i = 0; i < len(atomprops); ++i) { dvec.push_back(boost::python::extract(atomprops[i])); } return RDKit::Descriptors::BCUT2D(m, dvec); @@ -426,43 +436,40 @@ std::pair BCUT2D_tuple(const RDKit::ROMol &m, python::tuple atomp // From boost::python examples // Converts a std::pair instance to a Python tuple. template -struct std_pair_to_tuple -{ - static PyObject* convert(std::pair const& p) - { +struct std_pair_to_tuple { + static PyObject *convert(std::pair const &p) { return boost::python::incref( - boost::python::make_tuple(p.first, p.second).ptr()); + boost::python::make_tuple(p.first, p.second).ptr()); } - static PyTypeObject const *get_pytype () {return &PyTuple_Type; } + static PyTypeObject const *get_pytype() { return &PyTuple_Type; } }; - + // Helper for convenience. template -struct std_pair_to_python_converter -{ - std_pair_to_python_converter() - { - boost::python::to_python_converter< - std::pair, - std_pair_to_tuple, - true //std_pair_to_tuple has get_pytype - >(); +struct std_pair_to_python_converter { + std_pair_to_python_converter() { + boost::python::to_python_converter, + std_pair_to_tuple, + true // std_pair_to_tuple has get_pytype + >(); } }; -#endif +#endif } // namespace RDKit::SparseIntVect *GetMorganFingerprint( const RDKit::ROMol &mol, int radius, python::object invariants, python::object fromAtoms, bool useChirality, bool useBondTypes, - bool useFeatures, bool useCounts, python::object bitInfo, bool includeRedundantEnvironments) { - return MorganFingerprintHelper(mol, radius, -1, invariants, fromAtoms, - useChirality, useBondTypes, useFeatures, - useCounts, bitInfo, includeRedundantEnvironments); + bool useFeatures, bool useCounts, python::object bitInfo, + bool includeRedundantEnvironments) { + return MorganFingerprintHelper( + mol, radius, -1, invariants, fromAtoms, useChirality, useBondTypes, + useFeatures, useCounts, bitInfo, includeRedundantEnvironments); } RDKit::SparseIntVect *GetHashedMorganFingerprint( const RDKit::ROMol &mol, int radius, int nBits, python::object invariants, python::object fromAtoms, bool useChirality, bool useBondTypes, - bool useFeatures, python::object bitInfo, bool includeRedundantEnvironments) { + bool useFeatures, python::object bitInfo, + bool includeRedundantEnvironments) { return MorganFingerprintHelper(mol, radius, nBits, invariants, fromAtoms, useChirality, useBondTypes, useFeatures, true, bitInfo, includeRedundantEnvironments); @@ -471,7 +478,8 @@ RDKit::SparseIntVect *GetHashedMorganFingerprint( ExplicitBitVect *GetMorganFingerprintBV( const RDKit::ROMol &mol, int radius, unsigned int nBits, python::object invariants, python::object fromAtoms, bool useChirality, - bool useBondTypes, bool useFeatures, python::object bitInfo, bool includeRedundantEnvironments) { + bool useBondTypes, bool useFeatures, python::object bitInfo, + bool includeRedundantEnvironments) { std::vector *invars = nullptr; if (invariants) { unsigned int nInvar = @@ -501,7 +509,8 @@ ExplicitBitVect *GetMorganFingerprintBV( ExplicitBitVect *res; res = RDKit::MorganFingerprints::getFingerprintAsBitVect( mol, static_cast(radius), nBits, invars, froms.get(), - useChirality, useBondTypes, false, bitInfoMap, includeRedundantEnvironments); + useChirality, useBondTypes, false, bitInfoMap, + includeRedundantEnvironments); if (bitInfoMap) { bitInfo.attr("clear")(); for (RDKit::MorganFingerprints::BitInfoMap::const_iterator iter = @@ -521,15 +530,13 @@ ExplicitBitVect *GetMorganFingerprintBV( return res; } -python::list GetAtomFeatures(const RDKit::ROMol &mol, - int atomid, - bool addchiral) { - +python::list GetAtomFeatures(const RDKit::ROMol &mol, int atomid, + bool addchiral) { std::vector res; - RDKit::Descriptors::AtomFeatVect(mol, res, atomid, addchiral ); + RDKit::Descriptors::AtomFeatVect(mol, res, atomid, addchiral); python::list pyres; - for( auto iv : res ) { - pyres.append(iv); + for (auto iv : res) { + pyres.append(iv); } return pyres; } @@ -1413,12 +1420,13 @@ BOOST_PYTHON_MODULE(rdMolDescriptors) { (python::arg("mol")), docString.c_str(), python::return_value_policy()); - python::scope().attr("_GetAtomFeatures_version") = RDKit::Descriptors::AtomFeatVersion; + python::scope().attr("_GetAtomFeatures_version") = + RDKit::Descriptors::AtomFeatVersion; docString = "Returns the Atom Features vector"; - python::def( - "GetAtomFeatures", GetAtomFeatures, - (python::arg("mol"), python::arg("atomid"), python::arg("addchiral") = false), - docString.c_str()); + python::def("GetAtomFeatures", GetAtomFeatures, + (python::arg("mol"), python::arg("atomid"), + python::arg("addchiral") = false), + docString.c_str()); python::scope().attr("_CalcNumSpiroAtoms_version") = RDKit::Descriptors::NumSpiroAtomsVersion; @@ -1548,6 +1556,13 @@ BOOST_PYTHON_MODULE(rdMolDescriptors) { python::return_value_policy()); #ifdef RDK_BUILD_DESCRIPTORS3D + python::scope().attr("_CalcCoulombMat_version") = + RDKit::Descriptors::CoulombMatVersion; + docString = "Returns severals Coulomb randomized matrices"; + python::def("CalcCoulombMat", calcCoulombMat, + (python::arg("mol"), python::arg("confId") = -1), + docString.c_str()); + python::scope().attr("_CalcEMMcharges_version") = RDKit::Descriptors::EEMVersion; docString = "Returns EEM atomic partial charges"; @@ -1688,31 +1703,34 @@ BOOST_PYTHON_MODULE(rdMolDescriptors) { docString.c_str()); #endif - + #ifdef RDK_HAS_EIGEN3 - python::scope().attr("_BCUT2D_version") = - RDKit::Descriptors::BCUT2DVersion; - std::vector (*BCUT)(const RDKit::ROMol&) = & RDKit::Descriptors::BCUT2D; - std::pair (*BCUT_atomprops)(const RDKit::ROMol&, const std::string&) = & RDKit::Descriptors::BCUT2D; - docString = \ - "Implements BCUT descriptors From J. Chem. Inf. Comput. Sci., Vol. 39, No. 1, 1999" \ - "Diagonal elements are (currently) atomic mass, gasteiger charge,"\ - "crippen logP and crippen MRReturns the 2D BCUT2D descriptors vector as described in\n"\ - "returns [mass eigen value high, mass eigen value low,\n"\ - " gasteiger charge eigenvalue high, gasteiger charge low,\n"\ - " crippen lowgp eigenvalue high, crippen lowgp low,\n"\ - " crippen mr eigenvalue high, crippen mr low]\n"\ - ""; - - python::def("BCUT2D", BCUT, - (python::arg("mol")), - docString.c_str()); + python::scope().attr("_BCUT2D_version") = RDKit::Descriptors::BCUT2DVersion; + std::vector (*BCUT)(const RDKit::ROMol &) = + &RDKit::Descriptors::BCUT2D; + std::pair (*BCUT_atomprops)( + const RDKit::ROMol &, const std::string &) = &RDKit::Descriptors::BCUT2D; + docString = + "Implements BCUT descriptors From J. Chem. Inf. Comput. Sci., Vol. 39, " + "No. 1, 1999" + "Diagonal elements are (currently) atomic mass, gasteiger charge," + "crippen logP and crippen MRReturns the 2D BCUT2D descriptors vector as " + "described in\n" + "returns [mass eigen value high, mass eigen value low,\n" + " gasteiger charge eigenvalue high, gasteiger charge low,\n" + " crippen lowgp eigenvalue high, crippen lowgp low,\n" + " crippen mr eigenvalue high, crippen mr low]\n" + ""; + + python::def("BCUT2D", BCUT, (python::arg("mol")), docString.c_str()); std_pair_to_python_converter(); - docString = \ - "Returns a 2D BCUT (eigen value hi, eigenvalue low) given the molecule and the specified atom props\n"\ - " there length ot atom_props must a list or tuple of floats equal in size to the number of atoms in mol"; - + docString = + "Returns a 2D BCUT (eigen value hi, eigenvalue low) given the molecule " + "and the specified atom props\n" + " atom_props must be a list or tuple of floats equal in " + "size to the number of atoms in mol"; + python::def("BCUT2D", BCUT2D_list, (python::arg("mol"), python::arg("atom_props")), docString.c_str()); @@ -1720,12 +1738,12 @@ BOOST_PYTHON_MODULE(rdMolDescriptors) { (python::arg("mol"), python::arg("atom_props")), docString.c_str()); - docString = \ - "Returns a 2D BCUT (eigen value high, eigen value low) given the molecule and the specified atom prop name\n" \ - "atom_propname must exist on each aton and be convertable to a float"; + docString = + "Returns a 2D BCUT (eigen value high, eigen value low) given the " + "molecule and the specified atom prop name\n" + "atom_propname must exist on each atom and be convertible to a float"; python::def("BCUT2D", BCUT_atomprops, (python::arg("mol"), python::arg("atom_propname")), docString.c_str()); #endif - -} \ No newline at end of file +} diff --git a/Code/GraphMol/Descriptors/testCoulombMat.cpp b/Code/GraphMol/Descriptors/testCoulombMat.cpp new file mode 100644 index 000000000..66487ca58 --- /dev/null +++ b/Code/GraphMol/Descriptors/testCoulombMat.cpp @@ -0,0 +1,86 @@ +// Created by Guillaume GODIN +// Copyright (C) 2012-2018 Greg Landrum +// @@ 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. + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +std::vector tokenize(const std::string &s) { + boost::char_separator sep(", \n\r\t"); + boost::tokenizer> tok(s, sep); + std::vector tokens; + std::copy(tok.begin(), tok.end(), std::back_inserter >(tokens)); + return tokens; +} + + + +void testCoulombMat1(){ + std::cerr << "===================== Testing CoulombMat =======================\n"; + + std::string pathName = getenv("RDBASE"); + + // CM test 1 + std::string fNameCM = + pathName + "/Code/GraphMol/Descriptors/test_data/CM1.out"; + + std::string mol_file = + pathName + "/Code/GraphMol/Descriptors/test_data/bobmol.sdf"; + + std::ifstream instrmCM(fNameCM.c_str()); + std::string line; + std::vector tokens; + + RDKit::ROMOL_SPTR mol( RDKit::MolFileToMol( mol_file , true, false) ); + + + std::vector> Mres; + int confId = -1; + RDKit::Descriptors::CoulombMat(*mol, Mres, confId); + + for ( const auto &v : Mres ) { + std::getline(instrmCM, line); + tokens = tokenize(line); + + unsigned int ti = 0; + for ( double x : v ) { + // std::cout << x << ","; + TEST_ASSERT(std::fabs(std::stof(tokens[ti]) - x)< 0.001); + ti++; + } + // std::cout << "\n"; + } + // std::cout << "\n"; + + std::cerr << "CM test 1 mat Done\n"; +} + + +int main() { + RDLog::InitLogs(); + testCoulombMat1(); +} \ No newline at end of file diff --git a/Code/GraphMol/Descriptors/test_data/CM1.out b/Code/GraphMol/Descriptors/test_data/CM1.out new file mode 100644 index 000000000..f7d631e96 --- /dev/null +++ b/Code/GraphMol/Descriptors/test_data/CM1.out @@ -0,0 +1,20 @@ +53.3587,28.8825,17.6339,11.2255,8.72554,6.8079,8.65217,6.58329,6.60033,3.37463,3.27199,2.63155,2.64092,1.74263,1.72194,1.44489,1.42311,1.00482,1.09884,1.07798, +28.8825,36.8581,23.6974,14.4435,9.49153,7.21172,9.20791,2.98292,2.9352,5.42428,5.50323,2.73963,2.77078,2.22406,2.18813,1.49618,1.46838,1.01732,1.17144,1.16919, +17.6339,23.6974,36.8581,23.6352,14.7987,9.45274,11.3184,2.45891,2.31393,2.82976,2.77261,5.45955,5.43053,2.81015,2.76517,2.31103,2.2502,1.30863,1.46182,1.36136, +11.2255,14.4435,23.6352,36.8581,23.8652,14.3849,16.7549,1.54952,1.5025,2.19477,2.21189,2.81032,2.82206,5.44131,5.4622,2.83714,2.76358,1.73626,2.2004,2.03985, +8.72554,9.49153,14.7987,23.8652,36.8581,23.7572,20.3322,1.2663,1.22377,1.46916,1.46352,2.30132,2.2869,2.75577,2.78271,5.41015,5.41041,2.78637,2.80358,2.07202, +6.8079,7.21172,9.45274,14.3849,23.7572,36.8581,34.6815,0.974949,0.9523,1.16749,1.17062,1.48822,1.48602,2.14083,2.17601,2.74887,2.79197,5.41823,5.35397,3.00643, +8.65217,9.20791,11.3184,16.7549,20.3322,34.6815,73.5167,1.25019,1.1738,1.5804,1.47696,1.67679,1.85156,3.07811,2.46567,2.40036,3.0276,3.9189,3.87097,8.10652, +6.58329,2.98292,2.45891,1.54952,1.2663,0.974949,1.25019,0.5,0.582079,0.428958,0.337355,0.35725,0.447337,0.245789,0.2247,0.206959,0.217354,0.1458,0.153006,0.15277, +6.60033,2.9352,2.31393,1.5025,1.22377,0.9523,1.1738,0.582079,0.5,0.338807,0.412752,0.414909,0.33205,0.221307,0.23854,0.213156,0.196437,0.142353,0.155748,0.145129, +3.37463,5.42428,2.82976,2.19477,1.46916,1.16749,1.5804,0.428958,0.338807,0.5,0.551678,0.326312,0.416222,0.404512,0.319664,0.219255,0.235493,0.16426,0.186241,0.203829, +3.27199,5.50323,2.77261,2.21189,1.46352,1.17062,1.47696,0.337355,0.412752,0.551678,0.5,0.388841,0.325132,0.331923,0.39925,0.236454,0.21549,0.16411,0.199393,0.192683, +2.63155,2.73963,5.45955,2.81032,2.30132,1.48822,1.67679,0.35725,0.414909,0.326312,0.388841,0.5,0.557344,0.329364,0.399923,0.438276,0.340401,0.214469,0.236977,0.199229, +2.64092,2.77078,5.43053,2.82206,2.2869,1.48602,1.85156,0.447337,0.33205,0.416222,0.325132,0.557344,0.5,0.412667,0.327377,0.341573,0.421277,0.213861,0.217612,0.215306, +1.74263,2.22406,2.81015,5.44131,2.75577,2.14083,3.07811,0.245789,0.221307,0.404512,0.331923,0.329364,0.412667,0.5,0.558072,0.326758,0.393978,0.262764,0.317627,0.398743, +1.72194,2.18813,2.76517,5.4622,2.78271,2.17601,2.46567,0.2247,0.23854,0.319664,0.39925,0.399923,0.327377,0.558072,0.5,0.405057,0.324055,0.265435,0.397947,0.324765, +1.44489,1.49618,2.31103,2.83714,5.41015,2.74887,2.40036,0.206959,0.213156,0.219255,0.236454,0.438276,0.341573,0.326758,0.405057,0.5,0.562862,0.396722,0.401525,0.260011, +1.42311,1.46838,2.2502,2.76358,5.41041,2.79197,3.0276,0.217354,0.196437,0.235493,0.21549,0.340401,0.421277,0.393978,0.324055,0.562862,0.5,0.407749,0.326456,0.297044, +1.00482,1.01732,1.30863,1.73626,2.78637,5.41823,3.9189,0.1458,0.142353,0.16426,0.16411,0.214469,0.213861,0.262764,0.265435,0.396722,0.407749,0.5,0.546132,0.353014, +1.09884,1.17144,1.46182,2.2004,2.80358,5.35397,3.87097,0.153006,0.155748,0.186241,0.199393,0.236977,0.217612,0.317627,0.397947,0.401525,0.326456,0.546132,0.5,0.452921, +1.07798,1.16919,1.36136,2.03985,2.07202,3.00643,8.10652,0.15277,0.145129,0.203829,0.192683,0.199229,0.215306,0.398743,0.324765,0.260011,0.297044,0.353014,0.452921,0.5, diff --git a/Code/GraphMol/Descriptors/test_data/bobmol.sdf b/Code/GraphMol/Descriptors/test_data/bobmol.sdf new file mode 100644 index 000000000..ae4f24c7d --- /dev/null +++ b/Code/GraphMol/Descriptors/test_data/bobmol.sdf @@ -0,0 +1,44 @@ + + RDKit 3D + + 20 19 0 0 0 0 0 0 0 0999 V2000 + 3.3709 -0.4775 0.3706 N 0 0 0 0 0 0 0 0 0 0 0 0 + 2.0231 -1.0195 0.3052 C 0 0 0 0 0 0 0 0 0 0 0 0 + 1.0808 0.1720 0.2902 C 0 0 0 0 0 0 0 0 0 0 0 0 + -0.3637 -0.3063 0.2221 C 0 0 0 0 0 0 0 0 0 0 0 0 + -1.2298 0.9287 0.2108 C 0 0 0 0 0 0 0 0 0 0 0 0 + -2.7029 0.5798 0.1440 C 0 0 0 0 0 0 0 0 0 0 0 0 + -2.9424 -0.1384 -1.0146 O 0 0 0 0 0 0 0 0 0 0 0 0 + 3.4168 0.2414 -0.4115 H 0 0 0 0 0 0 0 0 0 0 0 0 + 3.4680 0.0361 1.2934 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1.9149 -1.5221 -0.6742 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1.7832 -1.6968 1.1252 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1.1861 0.8040 1.1831 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1.2583 0.7881 -0.6096 H 0 0 0 0 0 0 0 0 0 0 0 0 + -0.4776 -0.9003 -0.6999 H 0 0 0 0 0 0 0 0 0 0 0 0 + -0.6039 -0.9393 1.0871 H 0 0 0 0 0 0 0 0 0 0 0 0 + -0.9864 1.4993 1.1301 H 0 0 0 0 0 0 0 0 0 0 0 0 + -0.9764 1.5878 -0.6443 H 0 0 0 0 0 0 0 0 0 0 0 0 + -3.3013 1.5115 0.1327 H 0 0 0 0 0 0 0 0 0 0 0 0 + -2.9476 -0.0336 1.0494 H 0 0 0 0 0 0 0 0 0 0 0 0 + -2.9701 -1.1150 -0.8754 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 1 0 + 2 3 1 0 + 3 4 1 0 + 4 5 1 0 + 5 6 1 0 + 6 7 1 0 + 1 8 1 0 + 1 9 1 0 + 2 10 1 0 + 2 11 1 0 + 3 12 1 0 + 3 13 1 0 + 4 14 1 0 + 4 15 1 0 + 5 16 1 0 + 5 17 1 0 + 6 18 1 0 + 6 19 1 0 + 7 20 1 0 +M END