diff --git a/CHANGELOG.md b/CHANGELOG.md index fe9b484a01af77ba4f4c46198be80a9f98dbac7c..8ceb4f160067d48d49d734c0ba0acf86c44afdcf 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,6 +9,10 @@ Note that since we don't clearly distinguish between a public and private interf - Add PDBj as a ``pdb-provider`` option - Move Viewer APP to a separate file to allow use without importing light theme & index.html - Add symmetry support for mol2 files (only spacegroup setting 1) +- Fix mol2 files element symbol assignment +- Improve bond assignment from ``IndexPairBonds`` + - Add ``key`` field for mapping to source data + - Fix assignment of bonds with unphysical length - Fix label/stats of single atom selection in multi-chain units - [Breaking] Add rock animation to trackball controls - Add ``animate`` to ``TrackballControlsParams``, remove ``spin`` and ``spinSpeed`` diff --git a/src/mol-model-formats/structure/mol2.ts b/src/mol-model-formats/structure/mol2.ts index dcf060801f39185f946ba016066e5ed449aaa30e..40708bef5f43a66a868c8673639c67be3a65d48b 100644 --- a/src/mol-model-formats/structure/mol2.ts +++ b/src/mol-model-formats/structure/mol2.ts @@ -31,8 +31,17 @@ async function getModels(mol2: Mol2File, ctx: RuntimeContext) { const A = Column.ofConst('A', atoms.count, Column.Schema.str); const type_symbol = new Array<string>(atoms.count); + let hasAtomType = false; for (let i = 0; i < atoms.count; ++i) { - type_symbol[i] = guessElementSymbolString(atoms.atom_name.value(i)); + if (atoms.atom_type.value(i).includes('.')) { + hasAtomType = true; + break; + } + } + for (let i = 0; i < atoms.count; ++i) { + type_symbol[i] = hasAtomType + ? atoms.atom_type.value(i).split('.')[0].toUpperCase() + : guessElementSymbolString(atoms.atom_name.value(i)); } const atom_site = Table.ofPartialColumns(BasicSchema.atom_site, { @@ -77,6 +86,7 @@ async function getModels(mol2: Mol2File, ctx: RuntimeContext) { if (_models.frameCount > 0) { const indexA = Column.ofIntArray(Column.mapToArray(bonds.origin_atom_id, x => x - 1, Int32Array)); const indexB = Column.ofIntArray(Column.mapToArray(bonds.target_atom_id, x => x - 1, Int32Array)); + const key = bonds.bond_id; const order = Column.ofIntArray(Column.mapToArray(bonds.bond_type, x => { switch (x) { case 'ar': // aromatic @@ -103,7 +113,7 @@ async function getModels(mol2: Mol2File, ctx: RuntimeContext) { return BondType.Flag.Covalent; } }, Int8Array)); - const pairBonds = IndexPairBonds.fromData({ pairs: { indexA, indexB, order, flag }, count: atoms.count }); + const pairBonds = IndexPairBonds.fromData({ pairs: { key, indexA, indexB, order, flag }, count: atoms.count }); const first = _models.representative; IndexPairBonds.Provider.set(first, pairBonds); diff --git a/src/mol-model-formats/structure/property/bonds/index-pair.ts b/src/mol-model-formats/structure/property/bonds/index-pair.ts index 518581a03f9b78eeab3c89f9c6ce492e97983197..8bc506fd6993fd987391ac261756a7d0aa1062d7 100644 --- a/src/mol-model-formats/structure/property/bonds/index-pair.ts +++ b/src/mol-model-formats/structure/property/bonds/index-pair.ts @@ -1,5 +1,5 @@ /** - * Copyright (c) 2019-2021 Mol* contributors, licensed under MIT, See LICENSE file for more info. + * Copyright (c) 2019-2022 Mol* contributors, licensed under MIT, See LICENSE file for more info. * * @author Alexander Rose <alexander.rose@weirdbyte.de> */ @@ -10,9 +10,9 @@ import { Column } from '../../../../mol-data/db'; import { FormatPropertyProvider } from '../../common/property'; import { BondType } from '../../../../mol-model/structure/model/types'; import { ElementIndex } from '../../../../mol-model/structure'; -import { DefaultBondMaxRadius } from '../../../../mol-model/structure/structure/unit/bonds/common'; export type IndexPairsProps = { + readonly key: ArrayLike<number> readonly order: ArrayLike<number> readonly distance: ArrayLike<number> readonly flag: ArrayLike<BondType.Flag> @@ -22,17 +22,19 @@ export type IndexPairBonds = { bonds: IndexPairs, maxDistance: number } function getGraph(indexA: ArrayLike<ElementIndex>, indexB: ArrayLike<ElementIndex>, props: Partial<IndexPairsProps>, count: number): IndexPairs { const builder = new IntAdjacencyGraph.EdgeBuilder(count, indexA, indexB); + const key = new Int32Array(builder.slotCount); const order = new Int8Array(builder.slotCount); const distance = new Array(builder.slotCount); const flag = new Array(builder.slotCount); for (let i = 0, _i = builder.edgeCount; i < _i; i++) { builder.addNextEdge(); + builder.assignProperty(key, props.key ? props.key[i] : -1); builder.assignProperty(order, props.order ? props.order[i] : 1); builder.assignProperty(distance, props.distance ? props.distance[i] : -1); builder.assignProperty(flag, props.flag ? props.flag[i] : BondType.Flag.Covalent); } - return builder.createGraph({ order, distance, flag }); + return builder.createGraph({ key, order, distance, flag }); } export namespace IndexPairBonds { @@ -45,15 +47,33 @@ export namespace IndexPairBonds { export type Data = { pairs: { indexA: Column<number>, - indexB: Column<number> + indexB: Column<number>, + key?: Column<number>, order?: Column<number>, + /** + * Useful for bonds in periodic cells. That is, only bonds within the given + * distance are added. This allows for bond between periodic image but + * avoids unwanted bonds with wrong distances. If negative, test using the + * `maxDistance` option from `Props`. + */ distance?: Column<number>, flag?: Column<BondType.Flag>, }, count: number } - export const DefaultProps = { maxDistance: DefaultBondMaxRadius }; + export const DefaultProps = { + /** + * If negative, test using element-based threshold, otherwise distance in Angstrom. + * + * This option exists to handle bonds in periodic cells. For systems that are + * made from beads (as opposed to atomic elements), set to a specific distance. + * + * Note that `Data` has a `distance` field which allows specifying a distance + * for each bond individually which takes precedence over this option. + */ + maxDistance: -1 + }; export type Props = typeof DefaultProps export function fromData(data: Data, props: Partial<Props> = {}): IndexPairBonds { @@ -61,11 +81,12 @@ export namespace IndexPairBonds { const { pairs, count } = data; const indexA = pairs.indexA.toArray() as ArrayLike<ElementIndex>; const indexB = pairs.indexB.toArray() as ArrayLike<ElementIndex>; + const key = pairs.key && pairs.key.toArray(); const order = pairs.order && pairs.order.toArray(); const distance = pairs.distance && pairs.distance.toArray(); const flag = pairs.flag && pairs.flag.toArray(); return { - bonds: getGraph(indexA, indexB, { order, distance, flag }, count), + bonds: getGraph(indexA, indexB, { key, order, distance, flag }, count), maxDistance: p.maxDistance }; } diff --git a/src/mol-model/structure/structure/unit/bonds/common.ts b/src/mol-model/structure/structure/unit/bonds/common.ts index 061b6abc6ca90739b3a40c3199d19baa78ea3a2c..c242d9eb51eadc261d1f3139b2376fb62926f7a8 100644 --- a/src/mol-model/structure/structure/unit/bonds/common.ts +++ b/src/mol-model/structure/structure/unit/bonds/common.ts @@ -1,5 +1,5 @@ /** - * Copyright (c) 2017-2021 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> @@ -110,6 +110,15 @@ export function getElementThreshold(i: number) { return r; } +export function getPairingThreshold(elementIndexA: number, elementIndexB: number, thresholdA: number, thresholdB: number) { + const thresholdAB = getElementPairThreshold(elementIndexA, elementIndexB); + return thresholdAB > 0 + ? thresholdAB + : elementIndexB < 0 + ? thresholdA + : (thresholdA + thresholdB) / 1.95; // not sure if avg or min but max is too big +} + const H_ID = __ElementIndex['H']!; export function isHydrogen(i: number) { return i === H_ID; diff --git a/src/mol-model/structure/structure/unit/bonds/inter-compute.ts b/src/mol-model/structure/structure/unit/bonds/inter-compute.ts index ebb81e46a92c5fff863679de485388bf24bd17a3..d19aa14d8b87b3e2baa83e3ad402a2f8e645c1f6 100644 --- a/src/mol-model/structure/structure/unit/bonds/inter-compute.ts +++ b/src/mol-model/structure/structure/unit/bonds/inter-compute.ts @@ -1,5 +1,5 @@ /** - * Copyright (c) 2017-2021 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> @@ -8,7 +8,7 @@ import { BondType, MoleculeType } from '../../../model/types'; import { Structure } from '../../structure'; import { Unit } from '../../unit'; -import { getElementIdx, getElementPairThreshold, getElementThreshold, isHydrogen, BondComputationProps, MetalsSet, DefaultBondComputationProps } from './common'; +import { getElementIdx, getElementThreshold, isHydrogen, BondComputationProps, MetalsSet, DefaultBondComputationProps, getPairingThreshold } from './common'; import { InterUnitBonds, InterUnitEdgeProps } from './data'; import { SortedArray } from '../../../../../mol-data/int'; import { Vec3, Mat4 } from '../../../../../mol-math/linear-algebra'; @@ -82,11 +82,31 @@ function findPairBonds(unitA: Unit.Atomic, unitB: Unit.Atomic, props: BondComput const _bI = SortedArray.indexOf(unitB.elements, bI) as StructureElement.UnitIndex; if (_bI < 0) continue; - if (type_symbolA.value(aI) === 'H' && type_symbolB.value(bI) === 'H') continue; + + const aeI = getElementIdx(type_symbolA.value(aI)); + const beI = getElementIdx(type_symbolA.value(bI)); const d = distance[i]; const dist = getDistance(unitA, aI, unitB, bI); - if ((d !== -1 && equalEps(dist, d, 0.5)) || dist < maxDistance) { + + let add = false; + if (d >= 0) { + add = equalEps(dist, d, 0.3); + } else if (maxDistance >= 0) { + add = dist < maxDistance; + } else { + const pairingThreshold = getPairingThreshold( + aeI, beI, getElementThreshold(aeI), getElementThreshold(beI) + ); + add = dist < pairingThreshold; + + if (isHydrogen(aeI) && isHydrogen(beI)) { + // TODO handle molecular hydrogen + add = false; + } + } + + if (add) { builder.add(_aI, _bI, { order: order[i], flag: flag[i] }); } } @@ -155,13 +175,7 @@ function findPairBonds(unitA: Unit.Atomic, unitB: Unit.Atomic, props: BondComput const dist = Math.sqrt(squaredDistances[ni]); if (dist === 0) continue; - const thresholdAB = getElementPairThreshold(aeI, beI); - const pairingThreshold = thresholdAB > 0 - ? thresholdAB - : beI < 0 - ? thresholdA - : (thresholdA + getElementThreshold(beI)) / 1.95; // not sure if avg or min but max is too big - + const pairingThreshold = getPairingThreshold(aeI, beI, thresholdA, getElementThreshold(beI)); if (dist <= pairingThreshold) { const atomIdB = label_atom_idB.value(bI); const compIdB = label_comp_idB.value(residueIndexB[bI]); diff --git a/src/mol-model/structure/structure/unit/bonds/intra-compute.ts b/src/mol-model/structure/structure/unit/bonds/intra-compute.ts index 9a77f179656b166a9dfdfc4a7465714ffc2522a3..f9eae7c0568aadfbfb37498123cd57b8ac0ebaa6 100644 --- a/src/mol-model/structure/structure/unit/bonds/intra-compute.ts +++ b/src/mol-model/structure/structure/unit/bonds/intra-compute.ts @@ -1,5 +1,5 @@ /** - * Copyright (c) 2017-2021 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> @@ -9,7 +9,7 @@ import { BondType } from '../../../model/types'; import { IntraUnitBonds } from './data'; import { Unit } from '../../unit'; import { IntAdjacencyGraph } from '../../../../../mol-math/graph'; -import { BondComputationProps, getElementIdx, MetalsSet, getElementThreshold, isHydrogen, getElementPairThreshold, DefaultBondComputationProps } from './common'; +import { BondComputationProps, getElementIdx, MetalsSet, getElementThreshold, isHydrogen, DefaultBondComputationProps, getPairingThreshold } from './common'; import { SortedArray } from '../../../../../mol-data/int'; import { getIntraBondOrderFromTable } from '../../../model/properties/atomic/bonds'; import { StructureElement } from '../../element'; @@ -62,7 +62,8 @@ function findIndexPairBonds(unit: Unit.Atomic) { for (let _aI = 0 as StructureElement.UnitIndex; _aI < atomCount; _aI++) { const aI = atoms[_aI]; - const isHa = type_symbol.value(aI) === 'H'; + const aeI = getElementIdx(type_symbol.value(aI)); + const isHa = isHydrogen(aeI); const srcA = sourceIndex.value(aI); @@ -72,11 +73,30 @@ function findIndexPairBonds(unit: Unit.Atomic) { const _bI = SortedArray.indexOf(unit.elements, bI) as StructureElement.UnitIndex; if (_bI < 0) continue; - if (isHa && type_symbol.value(bI) === 'H') continue; + + const beI = getElementIdx(type_symbol.value(bI)); const d = distance[i]; const dist = getDistance(unit, aI, bI); - if ((d !== -1 && equalEps(dist, d, 0.5)) || dist < maxDistance) { + + let add = false; + if (d >= 0) { + add = equalEps(dist, d, 0.3); + } else if (maxDistance >= 0) { + add = dist < maxDistance; + } else { + const pairingThreshold = getPairingThreshold( + aeI, beI, getElementThreshold(aeI), getElementThreshold(beI) + ); + add = dist < pairingThreshold; + + if (isHa && isHydrogen(beI)) { + // TODO handle molecular hydrogen + add = false; + } + } + + if (add) { atomA[atomA.length] = _aI; atomB[atomB.length] = _bI; orders[orders.length] = order[i]; @@ -214,13 +234,7 @@ function findBonds(unit: Unit.Atomic, props: BondComputationProps): IntraUnitBon const dist = Math.sqrt(squaredDistances[ni]); if (dist === 0) continue; - const thresholdAB = getElementPairThreshold(aeI, beI); - const pairingThreshold = thresholdAB > 0 - ? thresholdAB - : beI < 0 - ? thresholdA - : (thresholdA + getElementThreshold(beI)) / 1.95; // not sure if avg or min but max is too big - + const pairingThreshold = getPairingThreshold(aeI, beI, thresholdA, getElementThreshold(beI)); if (dist <= pairingThreshold) { atomA[atomA.length] = _aI; atomB[atomB.length] = _bI;