Add option to avoid bonding cations. Issue #776

pymol has a similar option.
3dmol · Mar 11, 2024 · 1d5de21 · 1d5de21
1 parent 5e70a63
commit 1d5de21
Show file tree

Hide file tree

Showing 14 changed files with 12,241 additions and 17 deletions.
diff --git a/src/GLModel.ts b/src/GLModel.ts
@@ -2146,7 +2146,7 @@ export class GLModel {
      * consider valence constraints and will only create single bonds.
      */
     public assignBonds() {
-        assignBonds(this.atoms);
+        assignBonds(this.atoms, {assignBonds: true});
     }
 
     /** Remove specified atoms from model

diff --git a/src/parsers/CIF.ts b/src/parsers/CIF.ts
@@ -407,7 +407,7 @@ export function CIF(str: string, options: ParserOptionsSpec = {}) {
   }
   for (let i = 0; i < atoms.length; i++) {
     if (assignbonds && !(options.duplicateAssemblyAtoms && !options.dontConnectDuplicatedAtoms)) {
-      assignPDBBonds(atoms[i]);
+      assignPDBBonds(atoms[i], options);
     }
     computeSecondaryStructure(atoms[i],options.hbondCutoff);
     processSymmetries(
@@ -421,7 +421,7 @@ export function CIF(str: string, options: ParserOptionsSpec = {}) {
       !options.dontConnectDuplicatedAtoms &&
       assignbonds
     )
-    assignPDBBonds(atoms[i]);
+    assignPDBBonds(atoms[i],options);
   }
 
   return atoms;

diff --git a/src/parsers/CUBE.ts b/src/parsers/CUBE.ts
@@ -103,7 +103,7 @@ export function CUBE(str: string, options: ParserOptionsSpec) {
 
   if(assignbonds) {
       for (let i = 0; i < atoms.length; i++)
-          assignBonds(atoms[i]);
+          assignBonds(atoms[i], options);
   }
   return atoms;
 };
diff --git a/src/parsers/GRO.ts b/src/parsers/GRO.ts
@@ -1,3 +1,4 @@
+import { ParserOptionsSpec } from "./ParserOptionsSpec";
 import { assignPDBBonds } from "./utils/assignPDBBonds";
 import { atomNameToElem } from "./utils/atomNameToElem";
 
@@ -11,7 +12,7 @@ import { atomNameToElem } from "./utils/atomNameToElem";
  * @category Parsers
 */
 
-export function GRO(str: string /*, options*/) {
+export function GRO(str: string,  options: ParserOptionsSpec) {
   var allatoms: any[][] & Record<string, any> = [];
   var lines = str.split(/\r?\n|\r/);
   while (lines.length > 0) {
@@ -62,7 +63,7 @@ export function GRO(str: string /*, options*/) {
   }
 
   for (let i = 0; i < allatoms.length; i++) {
-    assignPDBBonds(allatoms[i]);
+    assignPDBBonds(allatoms[i], options);
   }
   return allatoms;
 }
diff --git a/src/parsers/LAMMPSTRJ.ts b/src/parsers/LAMMPSTRJ.ts
@@ -66,7 +66,7 @@ export function LAMMPSTRJ(str: string, options: ParserOptionsSpec) {
     start = offset + atomCount - 1;
   }
   if (options.assignBonds) {
-    for (var i = 0; i < atoms.length; i++) assignBonds(atoms[i]);
+    for (var i = 0; i < atoms.length; i++) assignBonds(atoms[i], options);
   }
   return atoms;
 }
diff --git a/src/parsers/PQR.ts b/src/parsers/PQR.ts
@@ -106,7 +106,7 @@ export function PQR(str: string, options: ParserOptionsSpec) {
 
       // assign bonds - yuck, can't count on connect records
       for (let i = 0; i < atoms.length; i++) {
-          assignPDBBonds(atoms[i]);
+          assignPDBBonds(atoms[i],options);
           if (computeStruct)
               computeSecondaryStructure(atoms[i],options.hbondCutoff);
       }

diff --git a/src/parsers/ParserOptionsSpec.ts b/src/parsers/ParserOptionsSpec.ts
@@ -63,6 +63,8 @@ import { AtomStyleSpec } from "GLModel";
   assemblyIndex?: number;
   /** for formats without explicit bonds (e.g. PDB, xyz) infer bonding (default true).  */
   assignBonds?: boolean;
+  /** for formats without explicit bonds, if assigning bonds to not assign bonds to common cations */
+  unboundCations?: boolean;
   /** set model to this style after parsing */
   style?: AtomStyleSpec;
  };
diff --git a/src/parsers/XYZ.ts b/src/parsers/XYZ.ts
@@ -77,7 +77,7 @@ export function XYZ(str: string, options: ParserOptionsSpec) {
 
   if (assignbonds) {
     for (let i = 0; i < atoms.length; i++) {
-      assignBonds(atoms[i]);
+      assignBonds(atoms[i], options);
     }
   }
 

diff --git a/src/parsers/utils/areConnected.ts b/src/parsers/utils/areConnected.ts
@@ -1,10 +1,16 @@
 import { AtomSpec } from "specs";
 import { bondLength } from "./bondLength";
+import { ParserOptionsSpec } from "parsers/ParserOptionsSpec";
+
+const cations = new Set(["Na","K","Ca","Mg","Mn","Sr"]);
 
 /*
  * Return true if atom1 and atom2 are probably bonded to each other based on distance alone
  */
-export function areConnected(atom1: AtomSpec, atom2: AtomSpec) {
+export function areConnected(atom1: AtomSpec, atom2: AtomSpec, options: ParserOptionsSpec) {
+  if(options && options.unboundCations && (cations.has(atom1.elem) || cations.has(atom2.elem))) {
+    return false;
+  }
   let maxsq = bondLength(atom1.elem) + bondLength(atom2.elem);
   maxsq += 0.25; // fudge factor, especially important for md frames, also see 1i3d
   maxsq *= maxsq;

diff --git a/src/parsers/utils/assignBonds.ts b/src/parsers/utils/assignBonds.ts
@@ -1,5 +1,6 @@
 import { AtomSpec } from "specs";
 import { areConnected } from "./areConnected";
+import { ParserOptionsSpec } from "parsers/ParserOptionsSpec";
 
 /**
  * @param {AtomSpec[]} atoms
@@ -21,7 +22,7 @@ const OFFSETS = [
 ];
 const MAX_BOND_LENGTH = 4.95; // (largest bond length, Cs) 2.25 * 2 * 1.1 (fudge factor)
 
-export function assignBonds(atoms: AtomSpec[]) {
+export function assignBonds(atoms: AtomSpec[], options: ParserOptionsSpec) {
   // Assign bonds - yuck, can't count on connect records
 
   for (let i = 0, n = atoms.length; i < n; i++) {
@@ -70,7 +71,7 @@ export function assignBonds(atoms: AtomSpec[]) {
       for (let j = 0; j < otherPoints.length; j++) {
         const atom2 = otherPoints[j];
 
-        if (areConnected(atom1, atom2)) {
+        if (areConnected(atom1, atom2, options)) {
           //gracefully handle one-sided bonds
           const a2i = atom1.bonds.indexOf(atom2.index);
           const a1i = atom2.bonds.indexOf(atom1.index);
@@ -103,7 +104,7 @@ export function assignBonds(atoms: AtomSpec[]) {
           const atom1 = points[i];
           for (let j = i + 1; j < points.length; j++) {
             const atom2 = points[j];
-            if (areConnected(atom1, atom2)) {
+            if (areConnected(atom1, atom2,options)) {
               if (atom1.bonds.indexOf(atom2.index) == -1) {
                 atom1.bonds.push(atom2.index);
                 atom1.bondOrder.push(1);

diff --git a/src/parsers/utils/assignPDBBonds.ts b/src/parsers/utils/assignPDBBonds.ts
@@ -4,13 +4,15 @@ import { AtomSpec } from "specs";
 import { areConnected } from "./areConnected";
 import { assignBonds } from "./assignBonds";
 import { standardResidues } from "./standardResidues";
+import { ParserOptionsSpec } from "parsers/ParserOptionsSpec";
+
 
 /**
  * @param {AtomSpec[]}
  *            atomsarray
  */
 
-export function assignPDBBonds(atomsarray: AtomSpec[]) {
+export function assignPDBBonds(atomsarray: AtomSpec[], options: ParserOptionsSpec) {
   // assign bonds - yuck, can't count on connect records
   const protatoms: Array<AtomSpec> = [];
   const hetatoms: Array<AtomSpec> = [];
@@ -21,7 +23,7 @@ export function assignPDBBonds(atomsarray: AtomSpec[]) {
     else protatoms.push(atom);
   }
 
-  assignBonds(hetatoms);
+  assignBonds(hetatoms, options);
 
   // sort by resid
   protatoms.sort(function (a, b) {
@@ -50,7 +52,7 @@ export function assignPDBBonds(atomsarray: AtomSpec[]) {
       const aj = protatoms[j];
       if (aj.chain !== ai.chain || aj.resi - ai.resi > 1) break;
 
-      if (areConnected(ai, aj)) {
+      if (areConnected(ai, aj, options)) {
         if (ai.bonds.indexOf(aj.index) === -1) {
           // only add if not already there
           ai.bonds.push(aj.index);

diff --git a/src/parsers/utils/getSinglePDB.ts b/src/parsers/utils/getSinglePDB.ts
@@ -254,7 +254,7 @@ export function getSinglePDB(
   //fix any "one-way" bonds in CONECT records
   validateBonds(atoms, serialToIndex);
   // assign bonds - yuck, can't count on connect records
-  if (assignbonds) assignPDBBonds(atoms);
+  if (assignbonds) assignPDBBonds(atoms, options);
 
   if (!noAssembly)
     processSymmetries(modelData.symmetries, atoms, options, modelData.cryst);