Split atoms by atomic mass in the canon algorithm, rather than do this by default which would lead to incorrect results in "unique smiles"

johnmay · johnmay · commit 7ec685d64b63 · 2021-01-17T22:29:22.000Z
diff --git a/base/standard/src/main/java/org/openscience/cdk/graph/invariant/Canon.java b/base/standard/src/main/java/org/openscience/cdk/graph/invariant/Canon.java
@@ -33,7 +33,6 @@
 import java.io.IOException;
 import java.io.InputStreamReader;
 import java.util.Arrays;
-import java.util.Collections;
 import java.util.Comparator;
 
 /**
@@ -102,6 +101,25 @@ private Canon(int[][] g, long[] partition, boolean[] hydrogens, boolean symOnly)
         symmetry = refine(labelling, hydrogens);
     }
 
+    /**
+     * Compute the canonical labels for the provided structure. The labelling
+     * does not consider isomer information or stereochemistry. The current
+     * implementation does not fully distinguish all structure topologies
+     * but in practise performs well in the majority of cases. A complete
+     * canonical labelling can be obtained using the {@link InChINumbersTools}
+     * but is computationally much more expensive.
+     *
+     * @param container structure
+     * @param g         adjacency list graph representation
+     * @param opts      canonical generation options see {@link CanonOpts}
+     * @return the canonical labelling
+     * @see EquivalentClassPartitioner
+     * @see InChINumbersTools
+     */
+    public static long[] label(IAtomContainer container, int[][] g, int opts) {
+        return label(container, g, basicInvariants(container, g, opts));
+    }
+
     /**
      * Compute the canonical labels for the provided structure. The labelling
      * does not consider isomer information or stereochemistry. The current
@@ -117,7 +135,7 @@ private Canon(int[][] g, long[] partition, boolean[] hydrogens, boolean symOnly)
      * @see InChINumbersTools
      */
     public static long[] label(IAtomContainer container, int[][] g) {
-        return label(container, g, basicInvariants(container, g));
+        return label(container, g, CanonOpts.Default);
     }
 
     /**
@@ -182,13 +200,32 @@ public static long[] label(IAtomContainer    container,
      *
      * @param container structure
      * @param g         adjacency list graph representation
+     * @param opts      canonical generation options see {@link CanonOpts}
      * @return symmetry classes
      * @see EquivalentClassPartitioner
      */
+    public static long[] symmetry(IAtomContainer container, int[][] g, int opts) {
+        return new Canon(g, basicInvariants(container, g, opts), terminalHydrogens(container, g), true).symmetry;
+    }
+
+    /**
+     * Compute the symmetry classes for the provided structure. There are known
+     * examples where symmetry is incorrectly found. The {@link
+     * EquivalentClassPartitioner} gives more accurate symmetry perception but
+     * this method is very quick and in practise successfully portions the
+     * majority of chemical structures.
+     *
+     * @param container structure
+     * @param g         adjacency list graph representation
+     * @return symmetry classes
+     * @see EquivalentClassPartitioner
+     * @see #basicInvariants(IAtomContainer, int[][], int)
+     */
     public static long[] symmetry(IAtomContainer container, int[][] g) {
-        return new Canon(g, basicInvariants(container, g), terminalHydrogens(container, g), true).symmetry;
+        return symmetry(container, g, CanonOpts.Default);
     }
 
+
     /**
      * Internal - refine invariants to a canonical labelling and
      * symmetry classes.
@@ -304,6 +341,18 @@ private long primeProduct(int[] ws, long[] ranks, boolean[] hydrogens) {
         return prod;
     }
 
+    /**
+     * See {@link #basicInvariants(IAtomContainer, int[][], int)}.
+     * @param container an atom container to generate labels for
+     * @param graph     graph representation (adjacency list)
+
+     * @return the initial invariants
+     * @see #basicInvariants(IAtomContainer, int[][], int)
+     */
+    public static long[] basicInvariants(IAtomContainer container, int[][] graph) {
+        return basicInvariants(container, graph, CanonOpts.Default);
+    }
+
     /**
      * Generate the initial invariants for each atom in the {@code container}.
      * The labels use the invariants described in {@cdk.cite WEI89}. 
@@ -328,11 +377,12 @@ private long primeProduct(int[] ws, long[] ranks, boolean[] hydrogens) {
      *
      * @param container an atom container to generate labels for
      * @param graph     graph representation (adjacency list)
+     * @param flav      bit mask canon flavor (see {@link CanonOpts})
      * @return initial invariants
      * @throws NullPointerException an atom had unset atomic number, hydrogen
      *                              count or formal charge
      */
-    public static long[] basicInvariants(IAtomContainer container, int[][] graph) {
+    public static long[] basicInvariants(IAtomContainer container, int[][] graph, int flav) {
 
         long[] labels = new long[graph.length];
 
@@ -362,6 +412,16 @@ public static long[] basicInvariants(IAtomContainer container, int[][] graph) {
             label <<= 4; // hydrogen count <= 15 (4 bits)
             label |= impH + expH & 0xf;
 
+            // atomic mass to split ties (if flavour requested), we can't do this
+            // by default because "unique" smiles doesn't include the isotopic mass
+            // so splitting on something that doesn't appear in the output would not
+            // function correctly
+            // n.b. the comparator based invariants are much more flexible still
+            if ((flav & CanonOpts.AtomicMass) != 0 && atom.getMassNumber() != null) {
+                label <<= 10;
+                label |= atom.getMassNumber();
+            }
+            
             labels[v] = label;
         }
         return labels;
diff --git a/base/standard/src/main/java/org/openscience/cdk/graph/invariant/CanonOpts.java b/base/standard/src/main/java/org/openscience/cdk/graph/invariant/CanonOpts.java
@@ -0,0 +1,39 @@
+/*
+ * Copyright (c) 2021 John Mayfield
+ *
+ * Contact: cdk-devel@lists.sourceforge.net
+ *
+ * This program is free software; you can redistribute it and/or modify it
+ * under the terms of the GNU Lesser General Public License as published by
+ * the Free Software Foundation; either version 2.1 of the License, or (at
+ * your option) any later version. All we ask is that proper credit is given
+ * for our work, which includes - but is not limited to - adding the above
+ * copyright notice to the beginning of your source code files, and to any
+ * copyright notice that you may distribute with programs based on this work.
+ *
+ * This program is distributed in the hope that it will be useful, but WITHOUT
+ * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ * FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
+ * License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 U
+ */
+
+package org.openscience.cdk.graph.invariant;
+
+/**
+ * Basic flavor options to tweak canonical invariants, note these deliberately mirror some fields
+ * from the {@link SmiFlavor} settings.
+ */
+public class CanonOpts {
+    /**
+     * Default canon flavour options.
+     */
+    public static final int Default    = 0;
+    /**
+     * Distinguish atoms based on atomic mass.
+     */
+    public static final int AtomicMass = 0x008;
+}
diff --git a/base/test-standard/src/test/java/org/openscience/cdk/graph/invariant/CanonTest.java b/base/test-standard/src/test/java/org/openscience/cdk/graph/invariant/CanonTest.java
@@ -24,7 +24,6 @@
 
 package org.openscience.cdk.graph.invariant;
 
-import org.junit.Assert;
 import org.junit.Test;
 import org.openscience.cdk.CDKConstants;
 import org.openscience.cdk.graph.GraphUtil;
@@ -34,7 +33,6 @@
 import org.openscience.cdk.smiles.SmiFlavor;
 import org.openscience.cdk.smiles.SmilesGenerator;
 import org.openscience.cdk.smiles.SmilesParser;
-import org.xmlcml.euclid.Int;
 
 import java.util.Comparator;
 
@@ -134,6 +132,15 @@ public void bridgingExplicitHydrogensAreIncluded() throws Exception {
         assertThat(mask, is(new boolean[]{false, false, false, false}));
     }
 
+    @Test
+    public void isotopeFlavor() throws Exception {
+        IAtomContainer m = smi("CC[13CH3]");
+        long[] symmetry = Canon.symmetry(m, GraphUtil.toAdjList(m), CanonOpts.Default);
+        assertThat(symmetry, is(new long[]{1,3,1}));
+        long[] symmetry2 = Canon.symmetry(m, GraphUtil.toAdjList(m), CanonOpts.AtomicMass);
+        assertThat(symmetry2, is(new long[]{1,2,3}));
+    }
+
     @Test
     public void explicitHydrogensIonsAreIncluded() throws Exception {
         IAtomContainer m = smi("[H+]");
