fix stereo

mojaie · Nov 16, 2024 · af7e0a7 · af7e0a7
1 parent be6aeed
commit af7e0a7
Show file tree

Hide file tree

Showing 2 changed files with 63 additions and 35 deletions.
diff --git a/src/stereo.jl b/src/stereo.jl
@@ -174,33 +174,53 @@ function stereocenter_from_sdf2d(g::SimpleGraph{T}, e_order, e_notation, e_isord
         dwcnt = count(x -> x === :down, drs)
         (upcnt == 0 && dwcnt == 0) && (@debug "$(i): unspecified"; continue)
         sortorder = [1, map(x -> x + 1, anglesort(v_coords2d, i, nbrs[1], nbrs[2:end]))...]
-        ons = nbrs[sortorder]  # node indices of neighbors in order
-        ods = drs[sortorder]  # direction symbols of neighbors in order
+        ons = nbrs[sortorder]  # node indices of neighbors in angular order
+        ods = drs[sortorder]  # direction symbols of neighbors in angular order
         if length(nbrs) == 3
-            if upcnt + dwcnt == 1
-                # if there is im
8000
plicit hydrogen, an up-wedge -> clockwise, a down-edge -> ccw
+            if upcnt == 1 && dwcnt == 1 || upcnt + dwcnt == 3
+                @debug "$(i): a reference plane cannot be defined $(ons) $(ods)"
+                push!(comments, "$(i): a reference plane cannot be defined")
+            elseif upcnt == 2 || dwcnt == 2
+                # up-wedge -> ccw, a down-edge -> clockwise
+                centers[i] = (ons[1], ons[2], ons[3], upcnt == 0)
+            else  # upcnt == 1 || dwcnt == 1
+                # up-wedge -> clockwise, a down-edge -> ccw
                 centers[i] = (ons[1], ons[2], ons[3], dwcnt == 0)
-            else
-                @debug "$(i): up and down wedges pointing in deg=3 stereocenter $(ons) $(ods)"
-                # Ambiguous stereochemistry
-                push!(comments, "$(i): up and down wedges pointing in deg=3 stereocenter")
-                continue
             end
-        elseif (ods[1] !== :unspecified && ods[3] !== :unspecified && ods[1] !== ods[3] ||
-                ods[2] !== :unspecified && ods[4] !== :unspecified && ods[2] !== ods[4])
-            @debug "$(i): two non-adjacent wedges are on opposite sides of the plane $(ons) $(ods)"
-            push!(comments, "$(i): two non-adjacent wedges are on opposite sides of the plane")
-            continue
-        elseif (ods[1] !== :unspecified && ods[1] === ods[2] ||
-                ods[2] !== :unspecified && ods[2] === ods[3] ||
-                ods[3] !== :unspecified && ods[3] === ods[4] ||
-                ods[4] !== :unspecified && ods[4] === ods[1])
-            @debug "$(i): two adjacent wedges are on the same side of the plane $(ons) $(ods)"
-            push!(comments, "$(i): two adjacent wedges are on the same side of the plane")
             continue
-        elseif upcnt != 0
+        end
+        if upcnt == 4 || dwcnt == 4
+            @debug "$(i): a reference plane cannot be defined $(ons) $(ods)"
+            push!(comments, "$(i): a reference plane cannot be defined")
+        elseif upcnt == 3
+            centers[i] = (ons[1], ons[2], ons[3], ods[1] === :up && ods[3] === :up)
+        elseif dwcnt == 3
+            centers[i] = (ons[1], ons[2], ons[3], ods[2] === :down && ods[4] === :down)
+        elseif upcnt == 2
+            if (ods[1] === :up && ods[3] === :up) || (ods[2] === :up && ods[4] === :up)
+                centers[i] = (ons[1], ons[2], ons[3], ods[1] === :up)
+            else
+                @debug "$(i): adjacent wedges lie on the same side of the reference plane $(ons) $(ods)"
+                push!(comments, "$(i): adjacent wedges lie on the same side of the reference plane")
+            end

8000
+        elseif dwcnt == 2
+            if (ods[1] === :down && ods[3] === :down) || (ods[2] === :down && ods[4] === :down)
+                centers[i] = (ons[1], ons[2], ons[3], ods[1] !== :down)
+            else
+                @debug "$(i): adjacent wedges lie on the same side of the reference plane $(ons) $(ods)"
+                push!(comments, "$(i): adjacent wedges lie on the same side of the reference plane")
+            end
+        elseif upcnt == 1 && dwcnt == 1
+            if ((ods[1] === :unspecified && ods[3] === :unspecified)
+                    || (ods[2] === :unspecified && ods[4] === :unspecified))
+                @debug "$(i): non-adjacent wedges lie on the opposite side of the reference plane $(ons) $(ods)"
+                push!(comments, "$(i): non-adjacent wedges lie on the opposite side of the reference plane")
+            else
+                centers[i] = (ons[1], ons[2], ons[3], ods[1] === :up || ods[3] === :up)
+            end
+        elseif upcnt == 1
             centers[i] = (ons[1], ons[2], ons[3], ods[1] === :up || ods[3] === :up)
-        else  # down wedges only
+        else  # dwcnt == 1
             centers[i] = (ons[1], ons[2], ons[3], ods[2] === :down || ods[4] === :down)
         end
     end

diff --git a/test/stereo.jl b/test/stereo.jl
@@ -139,42 +139,50 @@ end
         SDFBond(1, 6),
         SDFBond(1, 0)
     ]
-    wrong1 = MolGraph(edges, atoms, bonds, config_map=default_config)
-    @test isempty(get_prop(wrong1, :stereocenter))
-    @test haskey(wrong1.gprops, :stereocenter_ignored)
-    # serialization check
-    @test nv(MolGraph(to_json(wrong1))) == 5
+    mol7 = MolGraph(edges, atoms, bonds, config_map=default_config)
+    @test get_prop(mol7, :stereocenter)[1] == (2, 3, 5, true)
 
     bonds = [
         SDFBond(1, 1),
         SDFBond(1, 0),
         SDFBond(1, 1),
         SDFBond(1, 1)
     ]
-    wrong2 = MolGraph(edges, atoms, bonds, config_map=default_config)
-    @test isempty(get_prop(wrong2, :stereocenter))
-    @test haskey(wrong2.gprops, :stereocenter_ignored)
+    mol8 = MolGraph(edges, atoms, bonds, config_map=default_config)
+    @test get_prop(mol8, :stereocenter)[1] == (2, 3, 5, true)
 
     bonds = [
         SDFBond(1, 1),
         SDFBond(1, 0),
         SDFBond(1, 0),
         SDFBond(1, 6)
     ]
-    wrong3 = MolGraph(edges, atoms, bonds, config_map=default_config)
-    @test isempty(get_prop(wrong3, :stereocenter))
-    @test haskey(wrong3.gprops, :stereocenter_ignored)
+    wrong1 = MolGraph(edges, atoms, bonds, config_map=default_config)
+    @test isempty(get_prop(wrong1, :stereocenter))
+    @test haskey(wrong1.gprops, :stereocenter_ignored)
+    # serialization check
+    @test nv(MolGraph(to_json(wrong1))) == 5
 
     bonds = [
         SDFBond(1, 1),
         SDFBond(1, 1),
         SDFBond(1, 0),
         SDFBond(1, 0)
     ]
-    wrong4 = MolGraph(edges, atoms, bonds, config_map=default_config)
-    @test isempty(get_prop(wrong4, :stereocenter))
+    wrong2 = MolGraph(edges, atoms, bonds, config_map=default_config)
+    @test isempty(get_prop(wrong2
     @test haskey(wrong4.gprops, :stereocenter_ignored)
 
+    bonds = [
+        SDFBond(1, 1),
+        SDFBond(1, 1),
+        SDFBond(1, 1),
+        SDFBond(1, 1)
+    ]
+    wrong3 = MolGraph(edges, atoms, bonds, config_map=default_config)
+    @test isempty(get_prop(wrong3, :stereocenter))
+    @test haskey(wrong_all.gprops, :stereocenter_ignored)
+
     # degree=3
     atoms = [
         SDFAtom(:C, 0, 1, nothing, [0.0, 0.0]),