On Wed, Sep 14, 2016 at 4:16 AM, Guillaume GODIN < guillaume.go...@firmenich.com> wrote:
> Your solution is perfect! > > > glad it worked > I am currently implementing this article "http://www.mdpi.com/1420- > 3049/20/10/18279" using RDKit. > > Interesting, I hadn't seen that one. > It's now almost done, I need to check the results on Heat of Formation now > and understand why my code to get Angle in torsions is not working. > > def Angles(mol): > tors = mol.GetSubstructMatches(Chem.MolFromSmarts('[C]-[C;O;S]-[C; > O;S]-[C]')) > conf = mol.GetConformer() > A60=0 > A90=0 > A102=0 > for tor in tors: > Angle = abs(AllChem.GetDihedralDeg(conf,tor[0], tor[1], tor[2], > tor[3])) > if Angle<=60: > A60+=1 > elif Angle>60 and Angle<=90: > A90+=1 > elif Angle>90 and Angle<=102: > A102+=1 > return (A60,A90,A102) > > Comparing to the publication on Heat of Formation page 18294, I have to > much more occurences (ie. a factor of 50-100 more then expected). > I suspect that the problem is that the SMARTS you specify: '[C]-[C;O;S]-[C;O;S]-[C]' won't actually match anything. [C;O;S] means 'aliphatic C AND aliphatic O AND aliphatic S'. It think you like mean '[C,O,S]', which is 'aliphatic C OR aliphatic O OR aliphatic S'. You might also want to think about whether or not you should replace C, O, S (all of which only match aliphatic elements) with #6, #8, #16 (which can match aliphatic or aromatic) throughout the SMARTS > When it will be done I will share the code with you guys. > > Cool, I look forward to trying that out. Best, -greg > BR, > > Guillaume > ------------------------------ > *De :* Greg Landrum <greg.land...@gmail.com> > *Envoyé :* mercredi 14 septembre 2016 12:16 > *À :* Guillaume GODIN > *Cc :* RDKit Discuss > *Objet :* Re: [Rdkit-discuss] Angstroms Hydrogen bonding > > Hi Guillaume, > > On Tue, Sep 13, 2016 at 10:12 PM, Guillaume GODIN < > guillaume.go...@firmenich.com> wrote: > >> 1 Does 3D coordinates of a conformer is in Angstroms ? >> >> If you read the conformer from a file, for example a mol file, then the > 3D coordinates are in whatever units they were in in that file. This is > usually Angstrsom. > If you generate the conformation using one of the RDKit embedding > functions, then they are certainly in Angstroms. > >> 2 How to enumerate all HBonding to determine the bond length ? >> > > Interesting question. Here's a python function that might be a starting > point: > > def findHBonds(m,confId=-1,possiblePartners='[#8,#7]', > possibleHs='[#1][#8,#7]',distThresh=2.5): > conf = m.GetConformer(confId) > partners =[x[0] for x in m.GetSubstructMatches(Chem.MolFromSmarts( > possiblePartners))] > hs= [x[0] for x in m.GetSubstructMatches(Chem. > MolFromSmarts(possibleHs))] > res = [] > for h in hs: > ph = conf.GetAtomPosition(h) > for partner in partners: > if m.GetBondBetweenAtoms(h,partner) is not None: > continue > d = conf.GetAtomPosition(partner).Distance(ph) > if d<=distThresh: > res.append((h,partner,d)) > return res > > > In order to allow flexibility about what an H bond is, I left the > definitions of the acceptors (partners in the above code) and polar Hs > (just Hs in the above code) as SMARTS definitions so that they can be > customized. > > > > ********************************************************************** > 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 Rdkit-discuss@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/rdkit-discuss