Lomap MCS throws ValueError if rdkit.Mol created from SMILES

[dwhswenson]

from rdkit import Chem
from lomap import mcs

test_dir = "/Users/dwhs/omsf/src/openfe/openfe/tests/data/lomap_basic/"
toluene = Chem.MolFromMol2File(test_dir + "toluene.mol2")
methylcyclohexane = Chem.MolFromMol2File(test_dir + "methylcyclohexane.mol2")
mapping = mcs.MCS(toluene, methylcyclohexane)
# behaves as desired

toluene2 = Chem.MolFromSmiles(Chem.MolToSmiles(toluene))
methylcyclohexane2 = Chem.MolFromSmiles(Chem.MolToSmiles(methylcyclohexane))
mapping = mcs.MCS(toluene2, methylcyclohexane2)

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
~/mambaforge/envs/openfe/lib/python3.9/site-packages/lomap/mcs.py in __init__(self, moli, molj, time, verbose, max3d, threed)
    586             try:
--> 587                 trim_mcs_mol(max_deviation=max3d)
    588             except Exception as e:

~/mambaforge/envs/openfe/lib/python3.9/site-packages/lomap/mcs.py in trim_mcs_mol(max_deviation)
    156             while True:
--> 157                 (mapi,mapj) = best_substruct_match_to_mcs(self._moli_noh, self._molj_noh, by_rmsd=True)
    158                 # Compute the translation to bring molj's centre over moli

~/mambaforge/envs/openfe/lib/python3.9/site-packages/lomap/mcs.py in best_substruct_match_to_mcs(moli, molj, by_rmsd)
    121                     if by_rmsd:
--> 122                         coord_delta = (substructure_centre(moli,mapi)
    123                                  - substructure_centre(molj,mapj))

~/mambaforge/envs/openfe/lib/python3.9/site-packages/lomap/mcs.py in substructure_centre(mol, mol_sub)
     90             for i in mol_sub:
---> 91                 sum += mol.GetConformer().GetAtomPosition(i)
     92             return sum / len(mol_sub)

ValueError: Bad Conformer Id

During handling of the above exception, another exception occurred:

ValueError                                Traceback (most recent call last)
/var/folders/vj/28c107496sq5z4y10rjz5gdm0000gn/T/ipykernel_30741/2330220248.py in <module>
      1 toluene2 = Chem.MolFromSmiles(Chem.MolToSmiles(toluene))
      2 methylcyclohexane2 = Chem.MolFromSmiles(Chem.MolToSmiles(methylcyclohexane))
----> 3 mapping = mcs.MCS(toluene2, methylcyclohexane2)

~/mambaforge/envs/openfe/lib/python3.9/site-packages/lomap/mcs.py in __init__(self, moli, molj, time, verbose, max3d, threed)
    587                 trim_mcs_mol(max_deviation=max3d)
    588             except Exception as e:
--> 589                 raise ValueError(str(e))
    590 
    591         # Trim the MCS further to remove chirality mismatches

ValueError: Bad Conformer Id

[IAlibay]

What happens if you create 2D coordinates?

edit: ok yeah so it seems fine if you do Chem.rdCoordGen.AddCoords, so it's expecting "any conformer info" I guess?

[dwhswenson]

What happens if you create 2D coordinates?

That fixes:

Chem.rdDepictor.Compute2DCoords(toluene2)
Chem.rdDepictor.Compute2DCoords(methylcyclohexane2)
mapping = mcs.MCS(toluene2, methylcyclohexane2)
# no error

Should this be done internally when needed?

[IAlibay]

(follow up from a quick slack discusion)

My main question is - I have no idea if max3d should even be applied if there are non-3D coordinates, does it make sense to trim the MCS in 2D space?

The two options I could see are:

1. At init MCS checks for conformers, if there are non then max3d gets set to 0 (I'm not 100% sure what the impact of this is, but it would reflect not having any positional data)
2. At init if there are no conformers then either gracefully fail and ask for coords, or attempt to get some coords via some user-defined flag?