Merge pull request #378 from molstar/ar-bonds

fix visuals for aromatic/delocalized bonds

CHANGELOG.md

@@ -9,6 +9,10 @@ Note that since we don't clearly distinguish between a public and private interf
 - Fix parsing contour-level from emdb v3 header files
 - Fix invalid CSS (#376)
 - Fix "texture not renderable" & "texture not bound" warnings (#319)
+- Fix visual for bonds between two aromatic rings
+- Fix visual for delocalized bonds (parsed from mmcif and mol2)
+- Fix ring computation algorithm
+- Add ``UnitResonance`` property with info about delocalized triplets
 - Resolve marking in main renderer loop to improve overall performance
 - Use ``throttleTime`` instead of ``debounceTime`` in sequence viewer for better responsiveness
 

docs/interesting-pdb-entries.md

@@ -34,6 +34,14 @@
     * ACE (many, e.g. 5AGU, 1E1X)
     * ACY in 7ABY
     * NH2 (many, e.g. 6Y13)
+* Ligands with many rings
+    * STU (e.g. 1U59) - many fused rings
+    * HT (e.g. 127D) - rings connected by a single bond
+    * J2C (e.g. 7EFJ) - rings connected by a single atom
+    * RBF (e.g. 7QF2) - three linearly fused rings
+    * TA1 (e.g. 1JFF) - many fused rings (incl. a 8-member rings)
+    * BPA (e.g. 1JDG) - many fused rings
+    * CLR (e.g. 3GKI) - four fused rings
 
 Assembly symmetries
 * 5M30 (Assembly 1, C3 local and pseudo)

src/mol-model-formats/structure/mol2.ts

@@ -103,11 +103,11 @@ async function getModels(mol2: Mol2File, ctx: RuntimeContext) {
             const flag = Column.ofIntArray(Column.mapToArray(bonds.bond_type, x => {
                 switch (x) {
                     case 'ar': // aromatic
+                    case 'am': // amide
                         return BondType.Flag.Aromatic | BondType.Flag.Covalent;
                     case 'du': // dummy
                     case 'nc': // not connected
                         return BondType.Flag.None;
-                    case 'am': // amide
                     case 'un': // unknown
                     default:
                         return BondType.Flag.Covalent;

src/mol-model-formats/structure/property/bonds/chem_comp.ts

@@ -1,5 +1,5 @@
 /**
- * Copyright (c) 2017-2020 Mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2017-2022 Mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author David Sehnal <david.sehnal@gmail.com>
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
@@ -56,10 +56,10 @@ export namespace ComponentBond {
         const entries: Map<string, Entry> = new Map();
 
         function addEntry(id: string) {
-            // weird behavior when 'PRO' is requested - will report a single bond between N and H because a later operation would override real content
-            if (entries.has(id)) {
-                return entries.get(id)!;
-            }
+            // weird behavior when 'PRO' is requested - will report a single bond
+            // between N and H because a later operation would override real content
+            if (entries.has(id)) return entries.get(id)!;
+
             const e = new Entry(id);
             entries.set(id, e);
             return e;
@@ -83,10 +83,8 @@ export namespace ComponentBond {
             let ord = 1;
             if (aromatic) flags |= BondType.Flag.Aromatic;
             switch (order.toLowerCase()) {
-                case 'doub':
-                case 'delo':
-                    ord = 2;
-                    break;
+                case 'delo': flags |= BondType.Flag.Aromatic; break;
+                case 'doub': ord = 2; break;
                 case 'trip': ord = 3; break;
                 case 'quad': ord = 4; break;
             }

src/mol-model/structure/structure/unit.ts

@@ -25,6 +25,7 @@ import { Mat4, Vec3 } from '../../../mol-math/linear-algebra';
 import { IndexPairBonds } from '../../../mol-model-formats/structure/property/bonds/index-pair';
 import { ElementSetIntraBondCache } from './unit/bonds/element-set-intra-bond-cache';
 import { ModelSymmetry } from '../../../mol-model-formats/structure/property/symmetry';
+import { getResonance, UnitResonance } from './unit/resonance';
 
 /**
  * A building block of a structure that corresponds to an atomic or
@@ -282,6 +283,12 @@ namespace Unit {
             return this.props.rings;
         }
 
+        get resonance() {
+            if (this.props.resonance) return this.props.resonance;
+            this.props.resonance = getResonance(this);
+            return this.props.resonance;
+        }
+
         get polymerElements() {
             if (this.props.polymerElements) return this.props.polymerElements;
             this.props.polymerElements = getAtomicPolymerElements(this);
@@ -342,6 +349,7 @@ namespace Unit {
     interface AtomicProperties extends BaseProperties {
         bonds?: IntraUnitBonds
         rings?: UnitRings
+        resonance?: UnitResonance
         nucleotideElements?: SortedArray<ElementIndex>
         proteinElements?: SortedArray<ElementIndex>
         residueCount?: number

src/mol-model/structure/structure/unit/resonance.ts

@@ -0,0 +1,83 @@
+/**
+ * Copyright (c) 2022 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author Alexander Rose <alexander.rose@weirdbyte.de>
+ */
+
+import { SortedArray } from '../../../../mol-data/int/sorted-array';
+import { sortedCantorPairing } from '../../../../mol-data/util';
+import { BondType } from '../../model/types';
+import { StructureElement } from '../element';
+import { Unit } from '../unit';
+
+export type UnitResonance = {
+    /**
+     * Lookup for triplets of atoms in delocalized bonds.
+     *
+     * Does not include triplets that are part of aromatic rings.
+     */
+    readonly delocalizedTriplets: {
+        /** Return 3rd element in triplet or undefined if `a` and `b` are not part of a triplet */
+        readonly getThirdElement: (a: StructureElement.UnitIndex, b: StructureElement.UnitIndex) => StructureElement.UnitIndex | undefined
+        /** Return index into `triplets` or undefined if `a` is not part of any triplet */
+        readonly getTripletIndices: (a: StructureElement.UnitIndex) => number[] | undefined
+        readonly triplets: SortedArray<StructureElement.UnitIndex>[]
+    }
+}
+
+export function getResonance(unit: Unit.Atomic): UnitResonance {
+    return {
+        delocalizedTriplets: getDelocalizedTriplets(unit)
+    };
+}
+
+function getDelocalizedTriplets(unit: Unit.Atomic) {
+    const bonds = unit.bonds;
+    const { b, edgeProps, offset } = bonds;
+    const { order: _order, flags: _flags } = edgeProps;
+    const { elementAromaticRingIndices } = unit.rings;
+
+    const triplets: SortedArray<StructureElement.UnitIndex>[] = [];
+    const thirdElementMap = new Map<number, StructureElement.UnitIndex>();
+    const indicesMap = new Map<number, number[]>();
+
+    const add = (a: StructureElement.UnitIndex, b: StructureElement.UnitIndex, c: StructureElement.UnitIndex) => {
+        const index = triplets.length;
+        triplets.push(SortedArray.ofUnsortedArray([a, b, c]));
+        thirdElementMap.set(sortedCantorPairing(a, b), c);
+        if (indicesMap.has(a)) indicesMap.get(a)!.push(index);
+        else indicesMap.set(a, [index]);
+    };
+
+    for (let i = 0 as StructureElement.UnitIndex; i < unit.elements.length; i++) {
+        if (elementAromaticRingIndices.has(i)) continue;
+
+        const count = offset[i + 1] - offset[i] + 1;
+        if (count < 2) continue;
+
+        const deloBonds: StructureElement.UnitIndex[] = [];
+        for (let t = offset[i], _t = offset[i + 1]; t < _t; t++) {
+            const f = _flags[t];
+            if (!BondType.is(f, BondType.Flag.Aromatic)) continue;
+
+            deloBonds.push(b[t]);
+        }
+
+        if (deloBonds.length >= 2) {
+            add(i, deloBonds[0], deloBonds[1]);
+            for (let j = 1, jl = deloBonds.length; j < jl; j++) {
+                add(i, deloBonds[j], deloBonds[0]);
+            }
+        }
+    }
+
+    return {
+        getThirdElement: (a: StructureElement.UnitIndex, b: StructureElement.UnitIndex) => {
+            return thirdElementMap.get(sortedCantorPairing(a, b));
+        },
+        getTripletIndices: (a: StructureElement.UnitIndex) => {
+            return indicesMap.get(a);
+        },
+        triplets,
+    };
+}

src/mol-model/structure/structure/unit/rings/compute.ts

@@ -28,17 +28,19 @@ export function computeRings(unit: Unit.Atomic) {
 }
 
 const enum Constants {
-    MaxDepth = 4
+    MaxDepth = 5
 }
 
 interface State {
     startVertex: number,
     endVertex: number,
     count: number,
-    visited: Int32Array,
+    isRingAtom: Int32Array,
+    marked: Int32Array,
     queue: Int32Array,
     color: Int32Array,
     pred: Int32Array,
+    depth: Int32Array,
 
     left: Int32Array,
     right: Int32Array,
@@ -59,9 +61,11 @@ function State(unit: Unit.Atomic, capacity: number): State {
         startVertex: 0,
         endVertex: 0,
         count: 0,
-        visited: new Int32Array(capacity),
+        isRingAtom: new Int32Array(capacity),
+        marked: new Int32Array(capacity),
         queue: new Int32Array(capacity),
         pred: new Int32Array(capacity),
+        depth: new Int32Array(capacity),
         left: new Int32Array(Constants.MaxDepth),
         right: new Int32Array(Constants.MaxDepth),
         color: new Int32Array(capacity),
@@ -78,17 +82,26 @@ function State(unit: Unit.Atomic, capacity: number): State {
 
 function resetState(state: State) {
     state.count = state.endVertex - state.startVertex;
-    const { visited, pred, color } = state;
+    const { isRingAtom, pred, color, depth, marked } = state;
     for (let i = 0; i < state.count; i++) {
-        visited[i] = -1;
+        isRingAtom[i] = 0;
         pred[i] = -1;
+        marked[i] = -1;
         color[i] = 0;
+        depth[i] = 0;
     }
     state.currentColor = 0;
     state.currentAltLoc = '';
     state.hasAltLoc = false;
 }
 
+function resetDepth(state: State) {
+    const { depth } = state;
+    for (let i = 0; i < state.count; i++) {
+        depth[i] = state.count + 1;
+    }
+}
+
 function largestResidue(unit: Unit.Atomic) {
     const residuesIt = Segmentation.transientSegments(unit.model.atomicHierarchy.residueAtomSegments, unit.elements);
     let size = 0;
@@ -99,8 +112,16 @@ function largestResidue(unit: Unit.Atomic) {
     return size;
 }
 
+function isStartIndex(state: State, i: number) {
+    const bondOffset = state.bonds.offset;
+    const a = state.startVertex + i;
+    const bStart = bondOffset[a], bEnd = bondOffset[a + 1];
+    const bondCount = bEnd - bStart;
+    if (bondCount <= 1 || (state.isRingAtom[i] && bondCount === 2)) return false;
+    return true;
+}
+
 function processResidue(state: State, start: number, end: number) {
-    const { visited } = state;
     state.startVertex = start;
     state.endVertex = end;
 
@@ -117,11 +138,13 @@ function processResidue(state: State, start: number, end: number) {
     }
     arraySetRemove(altLocs, '');
 
+    let mark = 1;
     if (altLocs.length === 0) {
         resetState(state);
         for (let i = 0; i < state.count; i++) {
-            if (visited[i] >= 0) continue;
-            findRings(state, i);
+            if (!isStartIndex(state, i)) continue;
+            resetDepth(state);
+            mark = findRings(state, i, mark);
         }
     } else {
         for (let aI = 0; aI < altLocs.length; aI++) {
@@ -129,12 +152,13 @@ function processResidue(state: State, start: number, end: number) {
             state.hasAltLoc = true;
             state.currentAltLoc = altLocs[aI];
             for (let i = 0; i < state.count; i++) {
-                if (visited[i] >= 0) continue;
+                if (!isStartIndex(state, i)) continue;
                 const altLoc = state.altLoc.value(elements[state.startVertex + i]);
                 if (altLoc && altLoc !== state.currentAltLoc) {
                     continue;
                 }
-                findRings(state, i);
+                resetDepth(state);
+                mark = findRings(state, i, mark);
             }
         }
     }
@@ -144,10 +168,10 @@ function processResidue(state: State, start: number, end: number) {
     }
 }
 
-function addRing(state: State, a: number, b: number) {
+function addRing(state: State, a: number, b: number, isRingAtom: Int32Array) {
     // only "monotonous" rings
     if (b < a) {
-        return;
+        return false;
     }
 
     const { pred, color, left, right } = state;
@@ -176,7 +200,7 @@ function addRing(state: State, a: number, b: number) {
         if (current < 0) break;
     }
     if (!found) {
-        return;
+        return false;
     }
 
     current = a;
@@ -190,50 +214,50 @@ function addRing(state: State, a: number, b: number) {
     const len = leftOffset + rightOffset;
     // rings must have at least three elements
     if (len < 3) {
-        return;
+        return false;
     }
 
     const ring = new Int32Array(len);
     let ringOffset = 0;
-    for (let t = 0; t < leftOffset; t++) ring[ringOffset++] = state.startVertex + left[t];
-    for (let t = rightOffset - 1; t >= 0; t--) ring[ringOffset++] = state.startVertex + right[t];
+    for (let t = 0; t < leftOffset; t++) {
+        ring[ringOffset++] = state.startVertex + left[t];
+        isRingAtom[left[t]] = 1;
+    }
+    for (let t = rightOffset - 1; t >= 0; t--) {
+        ring[ringOffset++] = state.startVertex + right[t];
+        isRingAtom[right[t]] = 1;
+    }
 
     sortArray(ring);
 
-    if (state.hasAltLoc) {
-        // we need to check if the ring was already added because alt locs are present.
-
-        for (let rI = 0, _rI = state.currentRings.length; rI < _rI; rI++) {
-            const r = state.currentRings[rI];
-            if (ring[0] !== r[0]) continue;
-            if (ring.length !== r.length) continue;
+    // Check if the ring is unique and another one is not it's subset
+    for (let rI = 0, _rI = state.currentRings.length; rI < _rI; rI++) {
+        const r = state.currentRings[rI];
 
-            let areSame = true;
-            for (let aI = 0, _aI = ring.length; aI < _aI; aI++) {
-                if (ring[aI] !== r[aI]) {
-                    areSame = false;
-                    break;
-                }
-            }
-            if (areSame) {
-                return;
-            }
+        if (ring.length === r.length) {
+            if (SortedArray.areEqual(ring as any, r)) return false;
+        } else if (ring.length > r.length) {
+            if (SortedArray.isSubset(ring as any, r)) return false;
         }
     }
 
     state.currentRings.push(SortedArray.ofSortedArray(ring));
+
+    return true;
 }
 
-function findRings(state: State, from: number) {
-    const { bonds, startVertex, endVertex, visited, queue, pred } = state;
+function findRings(state: State, from: number, mark: number) {
+    const { bonds, startVertex, endVertex, isRingAtom, marked, queue, pred, depth } = state;
     const { elements } = state.unit;
     const { b: neighbor, edgeProps: { flags: bondFlags }, offset } = bonds;
-    visited[from] = 1;
+    marked[from] = mark;
+    depth[from] = 0;
     queue[0] = from;
     let head = 0, size = 1;
 
     while (head < size) {
         const top = queue[head++];
+        const d = depth[top];
         const a = startVertex + top;
         const start = offset[a], end = offset[a + 1];
 
@@ -250,18 +274,25 @@ function findRings(state: State, from: number) {
 
             const other = b - startVertex;
 
-            if (visited[other] > 0) {
+            if (marked[other] === mark) {
                 if (pred[other] !== top && pred[top] !== other) {
-                    addRing(state, top, other);
+                    if (addRing(state, top, other, isRingAtom)) {
+                        return mark + 1;
+                    }
                 }
                 continue;
             }
 
-            visited[other] = 1;
+            const newDepth = Math.min(depth[other], d + 1);
+            if (newDepth > Constants.MaxDepth) continue;
+
+            depth[other] = newDepth;
+            marked[other] = mark;
             queue[size++] = other;
             pred[other] = top;
         }
     }
+    return mark + 1;
 }
 
 export function getFingerprint(elements: string[]) {

src/mol-repr/structure/visual/bond-intra-unit-cylinder.ts

@@ -1,5 +1,5 @@
 /**
- * Copyright (c) 2018-2021 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2018-2022 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
  * @author David Sehnal <david.sehnal@gmail.com>
@@ -79,12 +79,16 @@ function getIntraUnitBondCylinderBuilderProps(unit: Unit.Atomic, structure: Stru
     };
 
     const { elementRingIndices, elementAromaticRingIndices } = unit.rings;
+    const deloTriplets = aromaticBonds ? unit.resonance.delocalizedTriplets : undefined;
 
     return {
         linkCount: edgeCount * 2,
         referencePosition: (edgeIndex: number) => {
             let aI = a[edgeIndex], bI = b[edgeIndex];
 
+            const rI = deloTriplets?.getThirdElement(aI, bI);
+            if (rI !== undefined) return pos(elements[rI], vRef);
+
             if (aI > bI) [aI, bI] = [bI, aI];
             if (offset[aI + 1] - offset[aI] === 1) [aI, bI] = [bI, aI];
 
@@ -145,8 +149,10 @@ function getIntraUnitBondCylinderBuilderProps(unit: Unit.Atomic, structure: Stru
                 if (isBondType(f, BondType.Flag.Aromatic) || (arCount && !ignoreComputedAromatic)) {
                     if (arCount === 2) {
                         return LinkStyle.MirroredAromatic;
-                    } else {
+                    } else if (arCount === 1 || deloTriplets?.getThirdElement(aI, bI)) {
                         return LinkStyle.Aromatic;
+                    } else {
+                        // case for bonds between two aromatic rings
                     }
                 }
             }

src/mol-repr/structure/visual/bond-intra-unit-line.ts

@@ -1,5 +1,5 @@
 /**
- * Copyright (c) 2020-2021 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2020-2022 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
  */
@@ -52,12 +52,16 @@ function createIntraUnitBondLines(ctx: VisualContext, unit: Unit, structure: Str
     const pos = unit.conformation.invariantPosition;
 
     const { elementRingIndices, elementAromaticRingIndices } = unit.rings;
+    const deloTriplets = aromaticBonds ? unit.resonance.delocalizedTriplets : undefined;
 
     const builderProps: LinkBuilderProps = {
         linkCount: edgeCount * 2,
         referencePosition: (edgeIndex: number) => {
             let aI = a[edgeIndex], bI = b[edgeIndex];
 
+            const rI = deloTriplets?.getThirdElement(aI, bI);
+            if (rI !== undefined) return pos(elements[rI], vRef);
+
             if (aI > bI) [aI, bI] = [bI, aI];
             if (offset[aI + 1] - offset[aI] === 1) [aI, bI] = [bI, aI];
 
@@ -106,8 +110,10 @@ function createIntraUnitBondLines(ctx: VisualContext, unit: Unit, structure: Str
                 if (isBondType(f, BondType.Flag.Aromatic) || (arCount && !ignoreComputedAromatic)) {
                     if (arCount === 2) {
                         return LinkStyle.MirroredAromatic;
-                    } else {
+                    } else if (arCount === 1 || deloTriplets?.getThirdElement(aI, bI)) {
                         return LinkStyle.Aromatic;
+                    } else {
+                        // case for bonds between two aromatic rings
                     }
                 }
             }

src/mol-repr/structure/visual/util/bond.ts

@@ -264,4 +264,4 @@ export function eachInterBond(loci: Loci, structure: Structure, apply: (interval
         __unitMap.clear();
     }
     return changed;
-}
+}