Package rdkit :: Package Chem :: Module inchi
[hide private]
[frames] | no frames]

Source Code for Module rdkit.Chem.inchi

  1  # $Id$ 
  2  # 
  3  #  Copyright (c) 2011, Novartis Institutes for BioMedical Research Inc. 
  4  #  All rights reserved. 
  5  #  
  6  # Redistribution and use in source and binary forms, with or without 
  7  # modification, are permitted provided that the following conditions are 
  8  # met:  
  9  # 
 10  #     * Redistributions of source code must retain the above copyright  
 11  #       notice, this list of conditions and the following disclaimer. 
 12  #     * Redistributions in binary form must reproduce the above 
 13  #       copyright notice, this list of conditions and the following  
 14  #       disclaimer in the documentation and/or other materials provided  
 15  #       with the distribution. 
 16  #     * Neither the name of Novartis Institutes for BioMedical Research Inc.  
 17  #       nor the names of its contributors may be used to endorse or promote  
 18  #       products derived from this software without specific prior written permission. 
 19  # 
 20  # THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 
 21  # "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 
 22  # LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 
 23  # A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT 
 24  # OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 
 25  # SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT 
 26  # LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, 
 27  # DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 
 28  # THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 
 29  # (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 
 30  # OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
 31  # 
 32   
 33  INCHI_AVAILABLE = True 
 34   
 35  import logging 
 36   
 37  from rdkit.Chem import rdinchi 
 38  from rdkit import RDLogger 
 39  logger = RDLogger.logger() 
 40   
 41  logLevelToLogFunctionLookup = { 
 42          logging.INFO : logger.info, 
 43          logging.DEBUG : logger.debug, 
 44          logging.WARNING : logger.warning, 
 45          logging.CRITICAL : logger.critical, 
 46          logging.ERROR : logger.error 
 47          } 
 48   
49 -class InchiReadWriteError(Exception):
50 pass
51
52 -def MolFromInchi(inchi, sanitize=True, removeHs=True, logLevel=None, 53 treatWarningAsError=False):
54 """Construct a molecule from a InChI string 55 56 Keyword arguments: 57 sanitize -- set to True to enable sanitization of the molecule. Default is 58 True 59 removeHs -- set to True to remove Hydrogens from a molecule. This only 60 makes sense when sanitization is enabled 61 logLevel -- the log level used for logging logs and messages from InChI 62 API. set to None to diable the logging completely 63 treatWarningAsError -- set to True to raise an exception in case of a 64 molecule that generates warning in calling InChI API. The resultant 65 molecule and error message are part of the excpetion 66 67 Returns: 68 a rdkit.Chem.rdchem.Mol instance 69 """ 70 try: 71 mol, retcode, message, log = rdinchi.InchiToMol(inchi, sanitize, removeHs) 72 except ValueError as e : 73 logger.error(str(e)) 74 return None 75 76 if logLevel is not None: 77 if logLevel not in logLevelToLogFunctionLookup: 78 raise ValueError("Unsupported log level: %d" % logLevel) 79 log = logLevelToLogFunctionLookup[logLevel] 80 if retcode == 0: 81 log(message) 82 83 if retcode != 0: 84 if retcode == 1: 85 logger.warning(message) 86 else: 87 logger.error(message) 88 if treatWarningAsError and retcode != 0: 89 raise InchiReadWriteError(mol, message) 90 return mol
91
92 -def MolToInchiAndAuxInfo(mol, options="", logLevel=None, 93 treatWarningAsError=False):
94 """Returns the standard InChI string and InChI auxInfo for a molecule 95 96 Keyword arguments: 97 logLevel -- the log level used for logging logs and messages from InChI 98 API. set to None to diable the logging completely 99 treatWarningAsError -- set to True to raise an exception in case of a 100 molecule that generates warning in calling InChI API. The resultant InChI 101 string and AuxInfo string as well as the error message are encoded in the 102 exception. 103 104 Returns: 105 a tuple of the standard InChI string and the auxInfo string returned by 106 InChI API, in that order, for the input molecule 107 """ 108 inchi, retcode, message, logs, aux = rdinchi.MolToInchi(mol, options) 109 if logLevel is not None: 110 if logLevel not in logLevelToLogFunctionLookup: 111 raise ValueError("Unsupported log level: %d" % logLevel) 112 log = logLevelToLogFunctionLookup[logLevel] 113 if retcode == 0: 114 log(message) 115 if retcode != 0: 116 if retcode == 1: 117 logger.warning(message) 118 else: 119 logger.error(message) 120 121 if treatWarningAsError and retcode != 0: 122 raise InchiReadWriteError(inchi, aux, message) 123 return inchi, aux
124
125 -def MolToInchi(mol, options="", logLevel=None, treatWarningAsError=False):
126 """Returns the standard InChI string for a molecule 127 128 Keyword arguments: 129 logLevel -- the log level used for logging logs and messages from InChI 130 API. set to None to diable the logging completely 131 treatWarningAsError -- set to True to raise an exception in case of a 132 molecule that generates warning in calling InChI API. The resultant InChI 133 string and AuxInfo string as well as the error message are encoded in the 134 exception. 135 136 Returns: 137 the standard InChI string returned by InChI API for the input molecule 138 """ 139 if options.find('AuxNone')==-1: 140 if options: 141 options += " /AuxNone" 142 else: 143 options += "/AuxNone" 144 145 try: 146 inchi, aux = MolToInchiAndAuxInfo(mol, options, logLevel=logLevel, 147 treatWarningAsError=treatWarningAsError) 148 except InchiReadWriteError as inst: 149 inchi, aux, message = inst.args 150 raise InchiReadWriteError(inchi, message) 151 return inchi
152
153 -def InchiToInchiKey(inchi):
154 """Return the InChI key for the given InChI string. Return None on error""" 155 ret = rdinchi.InchiToInchiKey(inchi) 156 if ret: return ret 157 else: return None
158 159 __all__ = ['MolToInchiAndAuxInfo', 'MolToInchi', 'MolFromInchi', 160 'InchiReadWriteError', 'InchiToInchiKey', 'INCHI_AVAILABLE'] 161