[Rdkit-discuss] 3D alignment in Python: align conformers of 2 molecules

[Thomas S. <beg...@ho...>]

I'm trying to align all conformers of 2 molecules (and keep the best ones) using the python api by following some of the tutorials:

http://nbviewer.ipython.org/gist/greglandrum/4316435/Working%20in%203D.ipynb

http://nbviewer.ipython.org/github/greglandrum/rdkit_blog/blob/master/notebooks/Using%20ConstrainedEmbed.ipynb

However whatever I try the alignment is considerably different than the one created using the Open 3D alignment node in Knime.

Currently the issue is that the molecules do not seem to be aligned properly. There is some sort of small shift. See code below.

Best Regards,

Thomas

        refCids = generateConformers(refMol, numConformers)
        mcs = MCS.FindMCS([refMol,mol], ringMatchesRingOnly=matchesRingOnly)        
        if mcs.completed == 1 and mcs.numAtoms > 0:        
            core = Chem.MolFromSmarts(mcs.smarts)
            logger.info('MCS: %s', Chem.MolToSmiles(core))
            
            refMatch = refMol.GetSubstructMatch(core) 
            match = mol.GetSubstructMatch(core)           

            # conformers for current target
            cids = generateConformers(mol, numConformers, coordMap=coordMap)            
            
            minRmsd = 1000;              
            for refCid in refCids:
                for cid in cids:
                    rmsd = AllChem.AlignMol(refMol, mol, prbCid=refCid, refCid=cid, atomMap=zip(refMatch,match))    
                    logger.debug('RMSD: %.2f', rmsd)
                    if rmsd < minRmsd:
                        logger.debug('New min RMSD: %.2f', rmsd)
                        minRmsd = rmsd
                        refConformerId = refCid


def generateConformers(mol, numConformers):
    AllChem.EmbedMolecule(mol)
    AllChem.MMFFOptimizeMolecule(mol)    
    cids=AllChem.EmbedMultipleConfs(mol, numConfs=numConformers, maxAttempts=50, pruneRmsThresh=0.5, coordMap=coordMap)
    for cid in cids: AllChem.MMFFOptimizeMolecule(mol,confId=cid)
    return cids