Dear Paul,
I'm afraid currently there is no public method to get r0 and kb from
ForceField contribs. I have attached a C++ sample program which shows
how to get the information you need from a ROMol. Please feel free to
contact me if you need further help.
Cheers,
p.
On 09/29/2014 07:36 PM, Paul Emsley wrote:
Dear Paolo (I presume :-),
I'd like to get the ro values for the bonds in an MMFF. I have the
following so far:
void mmff_bonds(RDKit::ROMol &mol_in) {
RDKit::MMFF::MMFFMolProperties *mmffMolProperties = new
RDKit::MMFF::MMFFMolProperties(mol_in);
ForceFields::ForceField *field = new ForceFields::ForceField();
RDKit::MMFF::Tools::addBonds( mol_in, mmffMolProperties, field);
ForceFields::ContribPtrVect::const_iterator contrib;
for(contrib=field->contribs().begin();
contrib!=field->contribs().end();
contrib++) {
ForceFields::MMFF::MMFFBond *mmffBondParams;
// something xxx!!!
double d = something(contrib???)->calcBondRestLength(mmffBondParams);
}
}
As you can see, I got lost between contrib and bond parameters. Please
can you tell me what I am missing? I hope it's only a few lines - but I
can't determine which ones :)
Thanks,
Paul.
------------------------------------------------------------------------------
Slashdot TV. Videos for Nerds. Stuff that Matters.
http://pubads.g.doubleclick.net/gampad/clk?id=160591471&iu=/4140/ostg.clktrk
_______________________________________________
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
--
==========================================================
Paolo Tosco, Ph.D.
Department of Drug Science and Technology
Via Pietro Giuria, 9 - 10125 Torino (Italy)
Tel: +39 011 670 7680 | Mob: +39 348 5537206
Fax: +39 011 670 7687 | E-mail: paolo.to...@unito.it
http://open3dqsar.org | http://open3dalign.org
==========================================================
/*
Build with:
g++ -Wall -o get_r0_kb get_r0_kb.cpp \
-I$BOOST_INCLUDE -I$RDBASE/Code -L$BOOST_LIB -L$RDBASE/lib \
-lFileParsers -lRDGeneral -lForceFieldHelpers
*/
#include <GraphMol/RDKitBase.h>
#include <GraphMol/FileParsers/FileParsers.h>
#include <ForceField/MMFF/Params.h>
#include <ForceField/MMFF/Contribs.h>
#include <ForceField/MMFF/BondStretch.h>
#include <GraphMol/ForceFieldHelpers/MMFF/AtomTyper.h>
#include <GraphMol/ForceFieldHelpers/MMFF/Builder.h>
#include <ForceField/ForceField.h>
#include <string>
#include <iostream>
using namespace RDKit;
int main(int argc, char **argv)
{
std::stringstream ss;
if (argc < 2) {
std::cerr << "Usage: get_r0_kb <molFile>" << std::endl;
return 0;
}
std::string sdf(argv[1]);
ROMol *mol = static_cast<ROMol *>(MolFileToMol(sdf, true, false));
MMFF::MMFFMolProperties *mp = new MMFF::MMFFMolProperties(*mol);
if (!(mp->isValid())) {
std::cerr << "Error setting up force-field" << std::endl;
}
else {
ForceFields::MMFF::MMFFBondCollection *mmffBond =
ForceFields::MMFF::MMFFBondCollection::getMMFFBond();
std::cout << std::setw(8) << "idx1" << std::setw(8) << "idx2"
<< std::setw(8) << "r0" << std::setw(8) << "kb" << std::endl;
for (ROMol::BondIterator bi = mol->beginBonds();
bi != mol->endBonds(); ++bi) {
unsigned int idx1 = (*bi)->getBeginAtomIdx();
unsigned int idx2 = (*bi)->getEndAtomIdx();
unsigned int iAtomType = mp->getMMFFAtomType(idx1);
unsigned int jAtomType = mp->getMMFFAtomType(idx2);
unsigned int bondType = mp->getMMFFBondType(*bi);
bool areMMFFBondParamsEmpirical = false;
const ForceFields::MMFF::MMFFBond *mmffBondParams =
(*mmffBond)(bondType, iAtomType, jAtomType);
if (!mmffBondParams) {
mmffBondParams = mp->getMMFFBondStretchEmpiricalRuleParams(*mol, *bi);
areMMFFBondParamsEmpirical = true;
}
double r0 = ForceFields::MMFF::Utils::calcBondRestLength(mmffBondParams);
double kb = ForceFields::MMFF::Utils::calcBondForceConstant(mmffBondParams);
if (areMMFFBondParamsEmpirical) {
delete mmffBondParams;
}
std::cout << std::setw(8) << idx1 + 1 << std::setw(8) << idx2 + 1
<< std::setw(8) << std::fixed << std::setprecision(4) << r0
<< std::setw(8) << std::fixed << std::setprecision(4) << kb
<< std::endl;
}
}
if (mp) delete mp;
if (mol) delete mol;
return 0;
}
------------------------------------------------------------------------------
Slashdot TV. Videos for Nerds. Stuff that Matters.
http://pubads.g.doubleclick.net/gampad/clk?id=160591471&iu=/4140/ostg.clktrk
_______________________________________________
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
