[Rdkit-discuss] bug in rdMolTransforms.SetDihedralDeg?

[Michal Krompiec]

Hello, I am trying to use the new functionality for manipulation of
dihedral angles in a function similar to OpenBabel's obrotate tool.
But it doesn't work: I get a ValueError exception. Here is an example
that replicates the error:

from rdkit import Chem
from rdkit.Chem import AllChem
from rdkit.Chem import rdMolTransforms

mol=Chem.MolFromSmiles('c2ccsc2c1sccc1') #2,2'-bithiophene
mol=Chem.AddHs(mol)
AllChem.EmbedMolecule(mol)
sp=Chem.MolFromSmarts("c2([H])ccsc2c1sccc1")
atoms=(5,6,7,8)
newangle=180.0
maplist = mol.GetSubstructMatches(sp)
if (len(maplist)>0):
for match in maplist :
a=[]
for i in range (0,4) :
a.append(match[atoms[i]])
angle=rdMolTransforms.GetDihedralDeg(mol.GetConformer(), a[0],
a[1], a[2], a[3])
print("angle between atoms {}".format(a)+" is
{}".format(angle))
print("trying to set to the same angle")
rdMolTransforms.SetDihedralDeg(mol.GetConformer(), a[0], a[1],
a[2], a[3], angle) #if you comment this line out
print("trying to set to another angle")
rdMolTransforms.SetDihedralDeg(mol.GetConformer(), a[0], a[1],
a[2], a[3], newangle) #this one will crash as well
angle=rdMolTransforms.GetDihedralDeg(mol.GetConformer(), a[0],
a[1], a[2], a[3])
print("angle between atoms {}".format(a)+" is
{}".format(angle))

Best wishes,
Michal Krompiec

[Paolo Tosco]

Dear Michal,

I have the impression that you wanted to select atoms 5,6,7,8, which do
indicate the bond between the two rings, but the RDKit counts atoms
starting from 0, so the tuple

atoms=(5,6,7,8)

should probably read

atoms=(4,5,6,7)

and this actually works; instead, the first tuple indicates a bond
inside a ring.

I hope that I understood the problem correctly; if this is not the case
please let me know.

Kind regards,
Paolo

Post by Michal Krompiec
Hello, I am trying to use the new functionality for manipulation of
dihedral angles in a function similar to OpenBabel's obrotate tool.
But it doesn't work: I get a ValueError exception. Here is an example
from rdkit import Chem
from rdkit.Chem import AllChem
from rdkit.Chem import rdMolTransforms
mol=Chem.MolFromSmiles('c2ccsc2c1sccc1') #2,2'-bithiophene
mol=Chem.AddHs(mol)
AllChem.EmbedMolecule(mol)
sp=Chem.MolFromSmarts("c2([H])ccsc2c1sccc1")
atoms=(5,6,7,8)
newangle=180.0
maplist = mol.GetSubstructMatches(sp)
a=[]
a.append(match[atoms[i]])
angle=rdMolTransforms.GetDihedralDeg(mol.GetConformer(), a[0],
a[1], a[2], a[3])
print("angle between atoms {}".format(a)+" is
{}".format(angle))
print("trying to set to the same angle")
rdMolTransforms.SetDihedralDeg(mol.GetConformer(), a[0], a[1],
a[2], a[3], angle) #if you comment this line out
print("trying to set to another angle")
rdMolTransforms.SetDihedralDeg(mol.GetConformer(), a[0], a[1],
a[2], a[3], newangle) #this one will crash as well
angle=rdMolTransforms.GetDihedralDeg(mol.GetConformer(), a[0],
a[1], a[2], a[3])
print("angle between atoms {}".format(a)+" is
{}".format(angle))
Best wishes,
Michal Krompiec
------------------------------------------------------------------------------
October Webinars: Code for Performance
Free Intel webinars can help you accelerate application performance.
Explore tips for MPI, OpenMP, advanced profiling, and more. Get the most from
the latest Intel processors and coprocessors. See abstracts and register >
http://pubads.g.doubleclick.net/gampad/clk?id=60135991&iu=/4140/ostg.clktrk
_______________________________________________
Rdkit-discuss mailing list
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss