Dear Guillaume,
please find attached a C++ snippet which exemplifies how to do what you
need.
For instance, if you issue the following:
$ ./conjGrp c1ccccc1Cc1ccccc1CCc1ccccc1
6 0,1
13 1
14 2
you get a printout of a list where the first column is the index of the
heavy atom which is a neighbour of a conjugated group, and the second
column represents a comma-separated list of the indices of the
conjugated groups it is a neighbour of.
Cheers,
p.
On 30/10/2015 04:03, Guillaume GODIN wrote:
Dear All,
I looking for a c++ method to enumate all conjugated pi system of a
molecule.
The final goal is to list all neighbors of each conjugated pi systems.
I found the SetConjugation and GetisConjugated bonds method in rdkit api.
Is there any over functions I could use ?
related to this question can we use a SMART pattern to do this to ?
best regards
guillaume
**********************************************************************
DISCLAIMER
This email and any files transmitted with it, including replies and
forwarded copies (which may contain alterations) subsequently
transmitted from Firmenich, are confidential and solely for the use of
the intended recipient. The contents do not represent the opinion of
Firmenich except to the extent that it relates to their official
business.
**********************************************************************
------------------------------------------------------------------------------
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
#include <iostream>
#include <vector>
#include <map>
#include <GraphMol/RDKitBase.h>
#include <GraphMol/SmilesParse/SmilesParse.h>
#include <GraphMol/FileParsers/MolWriters.h>
/*
g++ -g -Wall -I$RDBASE/Code conjGrp.cpp -o conjGrp \
-L$RDBASE/lib -lRDGeneral -lGraphMol -lSmilesParse \
-lFileParsers
cl /MD /O2 /I. /I%RDBASE%\Code conjGrp.cpp ^
/I C:\32\boost_1_59_0_py34 ^
%RDBASE%\lib\RDGeneral.lib ^
%RDBASE%\lib\GraphMol.lib ^
%RDBASE%\lib\SmilesParse.lib ^
%RDBASE%\lib\RDGeometryLib.lib ^
%RDBASE%\lib\FileParsers.lib ^
C:\32\boost_1_59_0_py34\lib32-msvc-12.0\boost_system-vc120-mt-1_59.lib ^
/FeconjGrp.exe
*/
using namespace RDKit;
unsigned int enumNbrConjGrp(ROMol &mol,
std::vector<unsigned int> &atomConjGrpIdx,
std::vector<unsigned int> &bondConjGrpIdx,
std::map<unsigned int, std::vector<unsigned int> > &m)
{
unsigned int nConjGrp = 0;
unsigned int nb = mol.getNumBonds();
bondConjGrpIdx.resize(nb, nb);
unsigned int na = mol.getNumAtoms();
atomConjGrpIdx.resize(na, na);
for (unsigned int i = 0; i < nb; ++i) {
const Bond *bi = mol.getBondWithIdx(i);
unsigned int biBeginIdx = bi->getBeginAtomIdx();
unsigned int biEndIdx = bi->getEndAtomIdx();
if (bi->getIsConjugated() && (bondConjGrpIdx[i] == nb)) {
// assign this conjugate bond to the matching group, if any
for (unsigned int j = 0;
(bondConjGrpIdx[i] == nb) && (j < nb); ++j) {
if ((i == j) || (bondConjGrpIdx[j] == nb))
continue;
const Bond *bj = mol.getBondWithIdx(j);
if ((bj->getBeginAtomIdx() == biBeginIdx)
|| (bj->getBeginAtomIdx() == biEndIdx)
|| (bj->getEndAtomIdx() == biBeginIdx)
|| (bj->getEndAtomIdx() == biEndIdx))
bondConjGrpIdx[i] = bondConjGrpIdx[j];
}
// no existing group matches: create a new group
if (bondConjGrpIdx[i] == nb)
bondConjGrpIdx[i] = nConjGrp++;
}
}
for (unsigned int i = 0; i < nb; ++i) {
if (bondConjGrpIdx[i] != nb) {
const Bond *bi = mol.getBondWithIdx(i);
unsigned int biBeginIdx = bi->getBeginAtomIdx();
unsigned int biEndIdx = bi->getEndAtomIdx();
atomConjGrpIdx[biBeginIdx] = bondConjGrpIdx[i];
atomConjGrpIdx[biEndIdx] = bondConjGrpIdx[i];
}
}
for (unsigned int ai = 0; ai < na; ++ai) {
if (atomConjGrpIdx[ai] == na)
continue;
ROMol::ADJ_ITER nbrIdx, endNbrs;
boost::tie(nbrIdx, endNbrs) =
mol.getAtomNeighbors(mol.getAtomWithIdx(ai));
std::vector <unsigned int> v;
for (; nbrIdx != endNbrs; ++nbrIdx) {
const Atom *atom = mol[*nbrIdx].get();
if (atom->getAtomicNum() == 1)
continue;
unsigned int aiNbr = atom->getIdx();
if (atomConjGrpIdx[aiNbr] == na) {
if (m.find(aiNbr) == m.end())
m[aiNbr] = v;
m[aiNbr].push_back(atomConjGrpIdx[ai]);
}
}
}
return nConjGrp;
}
int main (int argc, char **argv)
{
bool showHelp = (argc < 2);
std::string smi;
if (!showHelp) {
smi = argv[1];
showHelp = ((smi == "-h") || (smi == "-H") || (smi == "--help"));
}
if (showHelp) {
std::cerr
<< "./conjGrp <mol SMILES>" << std::endl;
exit(0);
}
ROMol *mol = SmilesToMol(smi);
std::vector<unsigned int> atomConjGrpIdx;
std::vector<unsigned int> bondConjGrpIdx;
std::map<unsigned int, std::vector<unsigned int> > m;
enumNbrConjGrp(*mol, atomConjGrpIdx, bondConjGrpIdx, m);
for (std::map<unsigned int, std::vector<unsigned int> >::const_iterator
mit = m.begin(); mit != m.end(); ++mit) {
std::cout << mit->first << "\t";
for (std::vector<unsigned int>::const_iterator
vit = mit->second.begin(); vit != mit->second.end(); ++vit)
std::cout << *vit << (vit == (mit->second.end() - 1) ? "\n" : ",");
}
return 0;
}
------------------------------------------------------------------------------
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
