Ver Fonte

wip, per-unit dssp

Alexander Rose há 6 anos atrás
pai
commit
992d0f10ac

+ 2 - 2
src/mol-model-formats/structure/mmcif/parser.ts

@@ -200,7 +200,7 @@ function createStandardModel(format: mmCIF_Format, atom_site: AtomSite, sourceIn
         coarseHierarchy: coarse.hierarchy,
         coarseConformation: coarse.conformation,
         properties: {
-            secondaryStructure: getSecondaryStructure(format.data, atomic.hierarchy, atomic.conformation),
+            secondaryStructure: getSecondaryStructure(format.data, atomic.hierarchy),
             ...formatData
         },
         customProperties: new CustomProperties(),
@@ -226,7 +226,7 @@ function createModelIHM(format: mmCIF_Format, data: IHMData, formatData: FormatD
         coarseHierarchy: coarse.hierarchy,
         coarseConformation: coarse.conformation,
         properties: {
-            secondaryStructure: getSecondaryStructure(format.data, atomic.hierarchy, atomic.conformation),
+            secondaryStructure: getSecondaryStructure(format.data, atomic.hierarchy),
             ...formatData
         },
         customProperties: new CustomProperties(),

+ 2 - 12
src/mol-model-formats/structure/mmcif/secondary-structure.ts

@@ -7,22 +7,12 @@
 
 import { mmCIF_Database as mmCIF, mmCIF_Database } from 'mol-io/reader/cif/schema/mmcif'
 import { SecondaryStructureType } from 'mol-model/structure/model/types';
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
+import { AtomicHierarchy } from 'mol-model/structure/model/properties/atomic';
 import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
 import { Column } from 'mol-data/db';
 import { ChainIndex, ResidueIndex } from 'mol-model/structure/model/indexing';
-import { computeModelDSSP } from 'mol-model-props/computed/secondary-structure/dssp';
-
-// TODO add parameter to allow forcing computation
-export function getSecondaryStructure(data: mmCIF_Database, hierarchy: AtomicHierarchy, conformation: AtomicConformation): SecondaryStructure {
-    if (!data.struct_conf._rowCount && !data.struct_sheet_range._rowCount) {
-        return computeModelDSSP(hierarchy, conformation)
-    } else {
-        return getSecondaryStructureMmCif(data, hierarchy)
-    }
-}
 
-export function getSecondaryStructureMmCif(data: mmCIF_Database, hierarchy: AtomicHierarchy): SecondaryStructure {
+export function getSecondaryStructure(data: mmCIF_Database, hierarchy: AtomicHierarchy): SecondaryStructure {
     const map: SecondaryStructureMap = new Map();
     const elements: SecondaryStructure.Element[] = [{ kind: 'none' }];
     addHelices(data.struct_conf, map, elements);

+ 22 - 7
src/mol-model-props/computed/secondary-structure.ts

@@ -6,12 +6,19 @@
 
 import { CustomPropertyDescriptor, Structure } from 'mol-model/structure';
 import { Task } from 'mol-task';
-import { DSSPComputationParams, computeModelDSSP } from './secondary-structure/dssp';
+import { DSSPComputationParams, computeUnitDSSP } from './secondary-structure/dssp';
 import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
 import { ParamDefinition as PD } from 'mol-util/param-definition';
+import { Unit } from 'mol-model/structure/structure';
+import { idFactory } from 'mol-util/id-factory';
+
+const nextSecondaryStructureId = idFactory()
 
 export namespace ComputedSecondaryStructure {
-    export type Property = SecondaryStructure
+    export type Property = {
+        id: number
+        map: Map<number, SecondaryStructure>
+    }
 
     export function get(structure: Structure): Property | undefined {
         return structure.inheritedPropertyData.__ComputedSecondaryStructure__;
@@ -36,7 +43,7 @@ export namespace ComputedSecondaryStructure {
     export async function attachFromCifOrCompute(structure: Structure, params: Partial<SecondaryStructureComputationProps> = {}) {
         if (structure.customPropertyDescriptors.has(Descriptor)) return true;
 
-        const compSecStruc = computeSecondaryStructure(structure, params)
+        const compSecStruc = await computeSecondaryStructure(structure, params)
 
         structure.customPropertyDescriptors.add(Descriptor);
         set(structure, compSecStruc);
@@ -50,9 +57,17 @@ export const SecondaryStructureComputationParams = {
 export type SecondaryStructureComputationParams = typeof SecondaryStructureComputationParams
 export type SecondaryStructureComputationProps = PD.Values<SecondaryStructureComputationParams>
 
-function computeSecondaryStructure(structure: Structure, params: Partial<SecondaryStructureComputationProps>): ComputedSecondaryStructure.Property {
-    // TODO compute from structure not from model
+async function computeSecondaryStructure(structure: Structure, params: Partial<SecondaryStructureComputationProps>): Promise<ComputedSecondaryStructure.Property> {
+    const p = { ...PD.getDefaultValues(SecondaryStructureComputationParams), params }
+    // TODO take inter-unit hbonds into account
     // TODO use Zhang-Skolnik for CA alpha only parts or for coarse parts with per-residue elements
-    const { atomicHierarchy, atomicConformation } = structure.model
-    return computeModelDSSP(atomicHierarchy, atomicConformation, params)
+    const map = new Map<number, SecondaryStructure>()
+    for (let i = 0, il = structure.unitSymmetryGroups.length; i < il; ++i) {
+        const u = structure.unitSymmetryGroups[i].units[0]
+        if (Unit.isAtomic(u)) {
+            const secondaryStructure = await computeUnitDSSP(u, p)
+            map.set(u.invariantId, secondaryStructure)
+        }
+    }
+    return { id: nextSecondaryStructureId(), map }
 }

+ 20 - 20
src/mol-model-props/computed/secondary-structure/dssp.ts

@@ -8,18 +8,18 @@
 import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
 import { SecondaryStructureType } from 'mol-model/structure/model/types';
 import { ParamDefinition as PD } from 'mol-util/param-definition'
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
 import { assignBends } from './dssp/bends';
-import { calcBackboneHbonds } from './dssp/backbone-hbonds';
+import { calcUnitBackboneHbonds } from './dssp/backbone-hbonds';
 import { Ladder, Bridge, DSSPContext, DSSPType } from './dssp/common';
 import { assignTurns } from './dssp/turns';
 import { assignHelices } from './dssp/helices';
 import { assignLadders } from './dssp/ladders';
 import { assignBridges } from './dssp/bridges';
 import { assignSheets } from './dssp/sheets';
-import { calculateDihedralAngles } from './dssp/dihedral-angles';
-import { calcBackboneAtomIndices } from './dssp/backbone-indices';
-import { calcAtomicTraceLookup3D } from './dssp/trace-lookup';
+import { calculateUnitDihedralAngles } from './dssp/dihedral-angles';
+import { calcUnitProteinTraceLookup3D } from './dssp/trace-lookup';
+import { Unit } from 'mol-model/structure';
+import { getUnitProteinInfo } from './dssp/protein-info';
 
 /**
  * TODO bugs to fix:
@@ -36,19 +36,20 @@ export const DSSPComputationParams = {
 export type DSSPComputationParams = typeof DSSPComputationParams
 export type DSSPComputationProps = PD.Values<DSSPComputationParams>
 
-export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: AtomicConformation, params: Partial<DSSPComputationProps> = {}): SecondaryStructure {
-    const p = { ...PD.getDefaultValues(DSSPComputationParams), ...params };
+export async function computeUnitDSSP(unit: Unit.Atomic, params: DSSPComputationProps): Promise<SecondaryStructure> {
 
-    const { lookup3d, proteinResidues } = calcAtomicTraceLookup3D(hierarchy, conformation)
-    const backboneIndices = calcBackboneAtomIndices(hierarchy, proteinResidues)
-    const hbonds = calcBackboneHbonds(hierarchy, conformation, proteinResidues, backboneIndices, lookup3d)
+    const proteinInfo = getUnitProteinInfo(unit)
+    const { residueIndices } = proteinInfo
+    const lookup3d = calcUnitProteinTraceLookup3D(unit, residueIndices)
 
-    const residueCount = proteinResidues.length
+    const hbonds = calcUnitBackboneHbonds(unit, proteinInfo, lookup3d)
+
+    const residueCount = residueIndices.length
     const flags = new Uint32Array(residueCount)
 
     // console.log(`calculating secondary structure elements using ${ params.oldDefinition ? 'old' : 'revised'} definition and ${ params.oldOrdering ? 'old' : 'revised'} ordering of secondary structure elements`)
 
-    const torsionAngles = calculateDihedralAngles(hierarchy, conformation, proteinResidues, backboneIndices)
+    const torsionAngles = calculateUnitDihedralAngles(unit, proteinInfo)
 
     const ladders: Ladder[] = []
     const bridges: Bridge[] = []
@@ -57,18 +58,16 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi
     const getFlagName = params.oldOrdering ? getOriginalFlagName : getUpdatedFlagName
 
     const ctx: DSSPContext = {
-        params: p,
+        params,
         getResidueFlag,
         getFlagName,
 
-        hierarchy,
-        proteinResidues,
+        unit,
+        proteinInfo,
         flags,
         hbonds,
 
         torsionAngles,
-        backboneIndices,
-        conformation,
         ladders,
         bridges
     }
@@ -81,13 +80,13 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi
     assignSheets(ctx)
 
     const assignment = getDSSPAssignment(flags, getResidueFlag)
-    const type = new Uint32Array(hierarchy.residues._rowCount) as unknown as SecondaryStructureType[]
+    const type = new Uint32Array(residueCount) as unknown as SecondaryStructureType[]
     const keys: number[] = []
     const elements: SecondaryStructure.Element[] = []
 
-    for (let i = 0, il = proteinResidues.length; i < il; ++i) {
+    for (let i = 0, il = residueIndices.length; i < il; ++i) {
         const assign = assignment[i]
-        type[proteinResidues[i]] = assign
+        type[residueIndices[i]] = assign
         const flag = getResidueFlag(flags[i])
         // TODO is this expected behavior? elements will be strictly split depending on 'winning' flag
         if (elements.length === 0 /* would fail at very start */ || flag !== (elements[elements.length - 1] as SecondaryStructure.Helix | SecondaryStructure.Sheet | SecondaryStructure.Turn).flags /* flag changed */) {
@@ -95,6 +94,7 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi
         }
         keys[i] = elements.length - 1
     }
+    // TODO expose model to unit residueIndex mapping
 
     return SecondaryStructure(type, keys, elements)
 }

+ 21 - 21
src/mol-model-props/computed/secondary-structure/dssp/backbone-hbonds.ts

@@ -5,13 +5,13 @@
  * @author Sebastian Bittrich <sebastian.bittrich@rcsb.org>
  */
 
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
-import { SortedArray } from 'mol-data/int';
 import { IntAdjacencyGraph } from 'mol-math/graph';
-import { ResidueIndex } from 'mol-model/structure';
+import { Unit } from 'mol-model/structure';
 import { GridLookup3D } from 'mol-math/geometry';
 import { Vec3 } from 'mol-math/linear-algebra';
-import { DsspHbonds, BackboneAtomIndices } from './common';
+import { DsspHbonds } from './common';
+import { ProteinInfo } from './protein-info';
+import { ElementIndex } from 'mol-model/structure/model';
 
 /** max distance between two C-alpha atoms to check for hbond */
 const caMaxDist = 9.0;
@@ -53,14 +53,14 @@ function calcHbondEnergy(oPos: Vec3, cPos: Vec3, nPos: Vec3, hPos: Vec3) {
     return e
 }
 
-export function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices, lookup3d: GridLookup3D): DsspHbonds {
-    const { cIndices, hIndices, nIndices, oIndices } = backboneIndices
-    const { index } = hierarchy
-    const { x, y, z } = conformation
-    const { traceElementIndex } = hierarchy.derived.residue
+export function calcUnitBackboneHbonds(unit: Unit.Atomic, proteinInfo: ProteinInfo, lookup3d: GridLookup3D): DsspHbonds {
+    const { residueIndices, cIndices, hIndices, nIndices, oIndices } = proteinInfo
 
-    const residueCount = proteinResidues.length
-    const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
+    const { index } = unit.model.atomicHierarchy
+    const { invariantPosition } = unit.conformation
+    const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue
+
+    const residueCount = residueIndices.length
 
     const oAtomResidues: number[] = [];
     const nAtomResidues: number[] = [];
@@ -75,9 +75,9 @@ export function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: Ato
     const cPosPrev = Vec3()
     const oPosPrev = Vec3()
 
-    for (let i = 0, il = proteinResidues.length; i < il; ++i) {
+    for (let i = 0, il = residueIndices.length; i < il; ++i) {
         const oPI = i
-        const oRI = proteinResidues[i]
+        const oRI = residueIndices[i]
 
         const oAtom = oIndices[oPI]
         const cAtom = cIndices[oPI]
@@ -89,22 +89,22 @@ export function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: Ato
         // ignore C-terminal residue as acceptor
         if (index.findAtomOnResidue(oRI, 'OXT') !== -1) continue
 
-        position(oAtom, oPos)
-        position(cAtom, cPos)
-        position(caAtom, caPos)
+        invariantPosition(oAtom, oPos)
+        invariantPosition(cAtom, cPos)
+        invariantPosition(caAtom as ElementIndex, caPos)
 
         const { indices, count } = lookup3d.find(caPos[0], caPos[1], caPos[2], caMaxDist)
 
         for (let j = 0; j < count; ++j) {
             const nPI = indices[j]
 
-            // ignore bonds within a residue or to prev or next residue, TODO take chain border into account
+            // ignore bonds within a residue or to prev or next residue
             if (nPI === oPI || nPI - 1 === oPI || nPI + 1 === oPI) continue
 
             const nAtom = nIndices[nPI]
             if (nAtom === -1) continue
 
-            position(nAtom, nPos)
+            invariantPosition(nAtom, nPos)
 
             const hAtom = hIndices[nPI]
             if (hAtom === -1) {
@@ -116,14 +116,14 @@ export function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: Ato
                 const cAtomPrev = cIndices[nPIprev]
                 if (oAtomPrev === -1 || cAtomPrev === -1) continue
 
-                position(oAtomPrev, oPosPrev)
-                position(cAtomPrev, cPosPrev)
+                invariantPosition(oAtomPrev, oPosPrev)
+                invariantPosition(cAtomPrev, cPosPrev)
 
                 Vec3.sub(hPos, cPosPrev, oPosPrev)
                 const dist = Vec3.distance(oPosPrev, cPosPrev)
                 Vec3.scaleAndAdd(hPos, nPos, hPos, 1 / dist)
             } else {
-                position(hAtom, hPos)
+                invariantPosition(hAtom, hPos)
             }
 
             const e = calcHbondEnergy(oPos, cPos, nPos, hPos)

+ 16 - 15
src/mol-model-props/computed/secondary-structure/dssp/bends.ts

@@ -8,6 +8,7 @@
 import { Vec3 } from 'mol-math/linear-algebra';
 import { radToDeg } from 'mol-math/misc';
 import { DSSPContext, DSSPType } from './common';
+import { ElementIndex } from 'mol-model/structure';
 
 /**
  * Bend(i) =: [angle ((CW - Ca(i - 2)),(C"(i + 2) - C"(i))) > 70"]
@@ -15,20 +16,20 @@ import { DSSPContext, DSSPType } from './common';
  * Type: S
  */
 export function assignBends(ctx: DSSPContext) {
-    const flags = ctx.flags
-    const { x, y, z } = ctx.conformation
-    const { traceElementIndex } = ctx.hierarchy.derived.residue
+    const { unit, flags, proteinInfo } = ctx
+    const { position } = unit.conformation
+    const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue
 
-    const proteinResidues = ctx.proteinResidues
-    const residueCount = proteinResidues.length
+    const { residueIndices, nIndices } = proteinInfo
+    const residueCount = residueIndices.length
 
-    const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
+    // const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
+    const p = (i: ElementIndex | -1, v: Vec3) => i === -1 ? Vec3.setNaN(v) : position(i, v)
 
     const caPosPrev2 = Vec3()
     const caPos = Vec3()
     const caPosNext2 = Vec3()
 
-    const nIndices = ctx.backboneIndices.nIndices
     const cPos = Vec3()
     const nPosNext = Vec3()
 
@@ -39,24 +40,24 @@ export function assignBends(ctx: DSSPContext) {
         // check for peptide bond
         for (let k = 0; k < 4; k++) {
             let index = i + k - 2
-            position(traceElementIndex[index], cPos)
-            position(nIndices[index + 1], nPosNext)
+            p(traceElementIndex[index], cPos)
+            p(nIndices[index + 1], nPosNext)
             if (Vec3.squaredDistance(cPos, nPosNext) > 6.25 /* max squared peptide bond distance allowed */) {
                 continue f1
             }
         }
 
-        const oRIprev2 = proteinResidues[i - 2]
-        const oRI = proteinResidues[i]
-        const oRInext2 = proteinResidues[i + 2]
+        const oRIprev2 = residueIndices[i - 2]
+        const oRI = residueIndices[i]
+        const oRInext2 = residueIndices[i + 2]
 
         const caAtomPrev2 = traceElementIndex[oRIprev2]
         const caAtom = traceElementIndex[oRI]
         const caAtomNext2 = traceElementIndex[oRInext2]
 
-        position(caAtomPrev2, caPosPrev2)
-        position(caAtom, caPos)
-        position(caAtomNext2, caPosNext2)
+        p(caAtomPrev2, caPosPrev2)
+        p(caAtom, caPos)
+        p(caAtomNext2, caPosNext2)
 
         Vec3.sub(caMinus2, caPosPrev2, caPos)
         Vec3.sub(caPlus2, caPos, caPosNext2)

+ 2 - 2
src/mol-model-props/computed/secondary-structure/dssp/bridges.ts

@@ -24,12 +24,12 @@ import { DSSPContext, DSSPType, BridgeType, Bridge } from './common';
  * Type: B
  */
 export function assignBridges(ctx: DSSPContext) {
-    const { proteinResidues, hbonds, flags, bridges } = ctx
+    const { proteinInfo, hbonds, flags, bridges } = ctx
 
     const { offset, b } = hbonds
     let i: number, j: number
 
-    for (let k = 0, kl = proteinResidues.length; k < kl; ++k) {
+    for (let k = 0, kl = proteinInfo.residueIndices.length; k < kl; ++k) {
         for (let t = offset[k], _t = offset[k + 1]; t < _t; t++) {
             const l = b[t]
             if (k > l) continue

+ 4 - 14
src/mol-model-props/computed/secondary-structure/dssp/common.ts

@@ -7,10 +7,9 @@
 
 import { BitFlags } from 'mol-util';
 import { SecondaryStructureType } from 'mol-model/structure/model/types';
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
-import { SortedArray } from 'mol-data/int';
-import { ResidueIndex, ElementIndex } from 'mol-model/structure';
 import { IntAdjacencyGraph } from 'mol-math/graph';
+import { Unit } from 'mol-model/structure/structure';
+import { ProteinInfo } from './protein-info';
 
 export interface DSSPContext {
     params: {
@@ -20,15 +19,13 @@ export interface DSSPContext {
     getResidueFlag: (f: DSSPType) => SecondaryStructureType,
     getFlagName: (f: DSSPType) => String,
 
-    hierarchy: AtomicHierarchy
-    proteinResidues: SortedArray<ResidueIndex>
+    unit: Unit.Atomic
+    proteinInfo: ProteinInfo
     /** flags for each residue */
     flags: Uint32Array
     hbonds: DsspHbonds,
 
     torsionAngles: { phi: Float32Array, psi: Float32Array },
-    backboneIndices: BackboneAtomIndices,
-    conformation: AtomicConformation,
     ladders: Ladder[],
     bridges: Bridge[]
 }
@@ -56,13 +53,6 @@ namespace DSSPType {
     }
 }
 
-export interface BackboneAtomIndices {
-    cIndices: ArrayLike<ElementIndex | -1>
-    hIndices: ArrayLike<ElementIndex | -1>
-    oIndices: ArrayLike<ElementIndex | -1>
-    nIndices: ArrayLike<ElementIndex | -1>
-}
-
 export type DsspHbonds = IntAdjacencyGraph<{ readonly energies: ArrayLike<number> }>
 
 export interface Ladder {

+ 28 - 24
src/mol-model-props/computed/secondary-structure/dssp/dihedral-angles.ts

@@ -5,20 +5,24 @@
  * @author Sebastian Bittrich <sebastian.bittrich@rcsb.org>
  */
 
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
-import { BackboneAtomIndices } from './common';
-import { ResidueIndex } from 'mol-model/structure';
-import { SortedArray } from 'mol-data/int';
+import { Unit } from 'mol-model/structure';
 import { Vec3 } from 'mol-math/linear-algebra';
+import { ProteinInfo } from './protein-info';
+import { ElementIndex } from 'mol-model/structure/model';
 
-export function calculateDihedralAngles(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices): { phi: Float32Array, psi: Float32Array } {
-    const { cIndices, nIndices } = backboneIndices
-    const { index } = hierarchy
-    const { x, y, z } = conformation
-    const { traceElementIndex } = hierarchy.derived.residue
+export interface DihedralAngles {
+    phi: Float32Array
+    psi: Float32Array
+}
 
-    const residueCount = proteinResidues.length
-    const position = (i: number, v: Vec3) => i === -1 ? Vec3.setNaN(v) : Vec3.set(v, x[i], y[i], z[i])
+export function calculateUnitDihedralAngles(unit: Unit.Atomic, proteinInfo: ProteinInfo): DihedralAngles {
+    const { cIndices, nIndices, residueIndices } = proteinInfo
+    const { position } = unit.conformation
+    const { index } = unit.model.atomicHierarchy
+    const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue
+
+    const residueCount = residueIndices.length
+    const p = (i: ElementIndex | -1, v: Vec3) => i === -1 ? Vec3.setNaN(v) : position(i, v)
 
     let cPosPrev = Vec3(), caPosPrev = Vec3(), nPosPrev = Vec3()
     let cPos = Vec3(), caPos = Vec3(), nPos = Vec3()
@@ -29,21 +33,21 @@ export function calculateDihedralAngles(hierarchy: AtomicHierarchy, conformation
     const phi: Float32Array = new Float32Array(residueCount - 1)
     const psi: Float32Array = new Float32Array(residueCount - 1)
 
-    position(-1, cPosPrev)
-    position(-1, caPosPrev)
-    position(-1, nPosPrev)
+    p(-1, cPosPrev)
+    p(-1, caPosPrev)
+    p(-1, nPosPrev)
 
-    position(cIndices[0], cPos)
-    position(traceElementIndex[proteinResidues[0]], caPos)
-    position(nIndices[0], nPos)
+    p(cIndices[0], cPos)
+    p(traceElementIndex[residueIndices[0]], caPos)
+    p(nIndices[0], nPos)
 
-    position(cIndices[1], cPosNext)
-    position(traceElementIndex[proteinResidues[1]], caPosNext)
-    position(nIndices[1], nPosNext)
+    p(cIndices[1], cPosNext)
+    p(traceElementIndex[residueIndices[1]], caPosNext)
+    p(nIndices[1], nPosNext)
 
     for (let i = 0; i < residueCount - 1; ++i) {
         // ignore C-terminal residue as acceptor
-        if (index.findAtomOnResidue(proteinResidues[i], 'OXT') !== -1) continue
+        if (index.findAtomOnResidue(residueIndices[i], 'OXT') !== -1) continue
 
         // returns NaN for missing atoms
         phi[i] = Vec3.dihedralAngle(cPosPrev, nPos, caPos, cPos)
@@ -52,9 +56,9 @@ export function calculateDihedralAngles(hierarchy: AtomicHierarchy, conformation
         cPosPrev = cPos, caPosPrev = caPos, nPosPrev = nPos
         cPos = cPosNext, caPos = caPosNext, nPos = nPosNext
 
-        position(cIndices[i + 1], cPosNext)
-        position(traceElementIndex[proteinResidues[i + 1]], caPosNext)
-        position(nIndices[i + 1], nPosNext)
+        p(cIndices[i + 1], cPosNext)
+        p(traceElementIndex[residueIndices[i + 1]], caPosNext)
+        p(nIndices[i + 1], nPosNext)
     }
 
     return { phi, psi };

+ 3 - 2
src/mol-model-props/computed/secondary-structure/dssp/helices.ts

@@ -19,7 +19,8 @@ import { DSSPContext, DSSPType } from './common';
  * Type: G (n=3), H (n=4), I (n=5)
  */
 export function assignHelices(ctx: DSSPContext) {
-    const { proteinResidues, flags } = ctx
+    const { proteinInfo, flags } = ctx
+    const residueCount = proteinInfo.residueIndices.length
 
     const turnFlag = [DSSPType.Flag.T3S, DSSPType.Flag.T4S, DSSPType.Flag.T5S, DSSPType.Flag.T3, DSSPType.Flag.T4, DSSPType.Flag.T5]
     const helixFlag = [0, 0, 0, DSSPType.Flag.G, DSSPType.Flag.H, DSSPType.Flag.I]
@@ -28,7 +29,7 @@ export function assignHelices(ctx: DSSPContext) {
     for (let ni = 0; ni < helixCheckOrder.length; ni++) {
         const n = helixCheckOrder[ni]
 
-        for (let i = 1, il = proteinResidues.length - n; i < il; i++) {
+        for (let i = 1, il = residueCount - n; i < il; i++) {
             const fI = DSSPType.create(flags[i])
             const fI1 = DSSPType.create(flags[i - 1])
             const fI2 = DSSPType.create(flags[i + 1])

+ 18 - 8
src/mol-model-props/computed/secondary-structure/dssp/backbone-indices.ts → src/mol-model-props/computed/secondary-structure/dssp/protein-info.ts

@@ -5,22 +5,31 @@
  * @author Sebastian Bittrich <sebastian.bittrich@rcsb.org>
  */
 
-import { AtomicHierarchy } from 'mol-model/structure/model/properties/atomic';
+import { Unit, ResidueIndex, ElementIndex } from 'mol-model/structure';
 import { SortedArray } from 'mol-data/int';
-import { ResidueIndex, ElementIndex } from 'mol-model/structure';
-import { BackboneAtomIndices } from './common';
 
-export function calcBackboneAtomIndices(hierarchy: AtomicHierarchy, proteinResidues: SortedArray<ResidueIndex>): BackboneAtomIndices {
-    const residueCount = proteinResidues.length
-    const { index } = hierarchy
+export interface ProteinInfo {
+    readonly residueIndices: SortedArray<ResidueIndex>
+    readonly cIndices: ArrayLike<ElementIndex | -1>
+    readonly hIndices: ArrayLike<ElementIndex | -1>
+    readonly oIndices: ArrayLike<ElementIndex | -1>
+    readonly nIndices: ArrayLike<ElementIndex | -1>
+}
 
+export function getUnitProteinInfo(unit: Unit.Atomic): ProteinInfo {
+    const { index } = unit.model.atomicHierarchy
+    const { proteinElements, residueIndex } = unit
+    const residueCount = proteinElements.length
+
+    const unitProteinResidues = new Uint32Array(residueCount)
     const c = new Int32Array(residueCount)
     const h = new Int32Array(residueCount)
     const o = new Int32Array(residueCount)
     const n = new Int32Array(residueCount)
 
-    for (let i = 0, il = residueCount; i < il; ++i) {
-        const rI = proteinResidues[i]
+    for (let i = 0; i < residueCount; ++i) {
+        const rI = residueIndex[proteinElements[i]]
+        unitProteinResidues[i] = rI
         c[i] = index.findAtomOnResidue(rI, 'C')
         h[i] = index.findAtomOnResidue(rI, 'H')
         o[i] = index.findAtomOnResidue(rI, 'O')
@@ -28,6 +37,7 @@ export function calcBackboneAtomIndices(hierarchy: AtomicHierarchy, proteinResid
     }
 
     return {
+        residueIndices: SortedArray.ofSortedArray(unitProteinResidues),
         cIndices: c as unknown as ArrayLike<ElementIndex | -1>,
         hIndices: h as unknown as ArrayLike<ElementIndex | -1>,
         oIndices: o as unknown as ArrayLike<ElementIndex | -1>,

+ 8 - 16
src/mol-model-props/computed/secondary-structure/dssp/trace-lookup.ts

@@ -2,27 +2,19 @@
  * Copyright (c) 2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
- * @author Sebastian Bittrich <sebastian.bittrich@rcsb.org>
  */
 
-import { AtomicHierarchy, AtomicConformation } from 'mol-model/structure/model/properties/atomic';
-import { MoleculeType } from 'mol-model/structure/model/types';
 import { GridLookup3D } from 'mol-math/geometry';
 import { SortedArray } from 'mol-data/int';
+import { Unit } from 'mol-model/structure/structure';
 import { ResidueIndex } from 'mol-model/structure';
 
-export function calcAtomicTraceLookup3D(hierarchy: AtomicHierarchy, conformation: AtomicConformation) {
-    const { x, y, z } = conformation;
-    const { moleculeType, traceElementIndex } = hierarchy.derived.residue
-    const indices: number[] = []
-    const _proteinResidues: number[] = []
-    for (let i = 0, il = moleculeType.length; i < il; ++i) {
-        if (moleculeType[i] === MoleculeType.protein) {
-            indices[indices.length] = traceElementIndex[i]
-            _proteinResidues[_proteinResidues.length] = i
-        }
+export function calcUnitProteinTraceLookup3D(unit: Unit.Atomic, unitProteinResidues: SortedArray<ResidueIndex>): GridLookup3D {
+    const { x, y, z } = unit.model.atomicConformation;
+    const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue
+    const indices = new Uint32Array(unitProteinResidues.length)
+    for (let i = 0, il = unitProteinResidues.length; i < il; ++i) {
+        indices[indices.length] = traceElementIndex[unitProteinResidues[i]]
     }
-    const lookup3d = GridLookup3D({ x, y, z, indices: SortedArray.ofSortedArray(indices) }, 4);
-    const proteinResidues = SortedArray.ofSortedArray<ResidueIndex>(_proteinResidues)
-    return { lookup3d, proteinResidues }
+    return GridLookup3D({ x, y, z, indices: SortedArray.ofSortedArray(indices) });
 }

+ 2 - 12
src/mol-model-props/computed/secondary-structure/dssp/turns.ts

@@ -16,22 +16,12 @@ import { DSSPContext, DSSPType } from './common';
  * Type: T
  */
 export function assignTurns(ctx: DSSPContext) {
-    const { proteinResidues, hbonds, flags, hierarchy } = ctx
-    const { chains, residueAtomSegments, chainAtomSegments } = hierarchy
-    const { label_asym_id } = chains
+    const { proteinInfo, hbonds, flags } = ctx
 
     const turnFlag = [DSSPType.Flag.T3S, DSSPType.Flag.T4S, DSSPType.Flag.T5S, DSSPType.Flag.T3, DSSPType.Flag.T4, DSSPType.Flag.T5]
 
     for (let idx = 0; idx < 3; idx++) {
-        for (let i = 0, il = proteinResidues.length - 1; i < il; ++i) {
-            const rI = proteinResidues[i]
-            const cI = chainAtomSegments.index[residueAtomSegments.offsets[rI]]
-
-            // TODO should take sequence gaps into account
-            const rN = proteinResidues[i + idx + 3]
-            const cN = chainAtomSegments.index[residueAtomSegments.offsets[rN]]
-            // check if on same chain
-            if (!label_asym_id.areValuesEqual(cI, cN)) continue
+        for (let i = 0, il = proteinInfo.residueIndices.length - 1; i < il; ++i) {
 
             // check if hbond exists
             if (hbonds.getDirectedEdgeIndex(i, i + idx + 3) !== -1) {

+ 2 - 6
src/mol-model/structure/model/properties/seconday-structure.ts

@@ -1,17 +1,13 @@
 /**
- * Copyright (c) 2017 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2017-2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author David Sehnal <david.sehnal@gmail.com>
  */
 
 import { SecondaryStructureType } from '../types';
-import { idFactory } from 'mol-util/id-factory';
-
-const getNextSecondaryStructureId = idFactory()
 
 /** Secondary structure "indexed" by residues. */
 interface SecondaryStructure {
-    readonly id: number
     readonly type: ArrayLike<SecondaryStructureType>,
     /** index into the elements array */
     readonly key: ArrayLike<number>,
@@ -20,7 +16,7 @@ interface SecondaryStructure {
 }
 
 function SecondaryStructure(type: SecondaryStructure['type'], key: SecondaryStructure['key'], elements: SecondaryStructure['elements']) {
-    return { id: getNextSecondaryStructureId(), type, key, elements }
+    return { type, key, elements }
 }
 
 namespace SecondaryStructure {

+ 1 - 2
src/mol-repr/structure/visual/polymer-trace-mesh.ts

@@ -18,7 +18,6 @@ import { addTube } from 'mol-geo/geometry/mesh/builder/tube';
 import { VisualContext } from 'mol-repr/visual';
 import { Theme } from 'mol-theme/theme';
 import { ComputedSecondaryStructure } from 'mol-model-props/computed/secondary-structure';
-import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
 
 export const PolymerTraceMeshParams = {
     sizeFactor: PD.Numeric(0.2, { min: 0, max: 10, step: 0.01 }),
@@ -123,7 +122,7 @@ export function PolymerTraceVisual(materialId: number): UnitsVisual<PolymerTrace
             )
 
             const computedSecondaryStructure = ComputedSecondaryStructure.get(newStructureGroup.structure)
-            if ((state.info.computedSecondaryStructure as SecondaryStructure) !== computedSecondaryStructure) {
+            if ((state.info.computedSecondaryStructure as ComputedSecondaryStructure.Property) !== computedSecondaryStructure) {
                 state.createGeometry = true;
                 state.info.computedSecondaryStructure = computedSecondaryStructure
             }

+ 6 - 2
src/mol-repr/structure/visual/util/polymer/trace-iterator.ts

@@ -1,5 +1,5 @@
 /**
- * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2018-2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
  */
@@ -249,8 +249,12 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
         this.value = createPolymerTraceElement(unit)
         this.hasNext = this.residueIt.hasNext && this.polymerIt.hasNext
 
+        this.secondaryStructureType = unit.model.properties.secondaryStructure.type
         const computedSecondaryStructure = ComputedSecondaryStructure.get(structure)
-        this.secondaryStructureType = (computedSecondaryStructure && computedSecondaryStructure.type) || unit.model.properties.secondaryStructure.type
+        if (computedSecondaryStructure) {
+            const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId)
+            if (secondaryStructure) this.secondaryStructureType = secondaryStructure.type
+        }
     }
 }
 

+ 7 - 6
src/mol-theme/color/secondary-structure.ts

@@ -1,5 +1,5 @@
 /**
- * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ * Copyright (c) 2018-2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
  */
@@ -14,7 +14,6 @@ import { ParamDefinition as PD } from 'mol-util/param-definition'
 import { ThemeDataContext } from '../theme';
 import { TableLegend } from 'mol-util/color/tables';
 import { ComputedSecondaryStructure } from 'mol-model-props/computed/secondary-structure';
-import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
 
 // from Jmol http://jmol.sourceforge.net/jscolors/ (shapely)
 const SecondaryStructureColors = ColorMap({
@@ -42,12 +41,14 @@ export function getSecondaryStructureColorThemeParams(ctx: ThemeDataContext) {
     return SecondaryStructureColorThemeParams // TODO return copy
 }
 
-export function secondaryStructureColor(unit: Unit, element: ElementIndex, computedSecondaryStructure: SecondaryStructure | undefined): Color {
+export function secondaryStructureColor(unit: Unit, element: ElementIndex, computedSecondaryStructure: ComputedSecondaryStructure.Property | undefined): Color {
     let secStrucType = SecondaryStructureType.create(SecondaryStructureType.Flag.None)
     if (Unit.isAtomic(unit)) {
-        secStrucType = computedSecondaryStructure ?
-            computedSecondaryStructure.type[unit.residueIndex[element]] :
-            unit.model.properties.secondaryStructure.type[unit.residueIndex[element]]
+        secStrucType = unit.model.properties.secondaryStructure.type[unit.residueIndex[element]]
+        if (computedSecondaryStructure) {
+            const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId)
+            if (secondaryStructure) secStrucType = secondaryStructure.type[unit.residueIndex[element]]
+        }
     }
 
     if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Helix)) {

+ 4 - 5
src/tests/browser/render-structure.ts

@@ -15,7 +15,7 @@ import { trajectoryFromMmCIF } from 'mol-model-formats/structure/mmcif';
 import { MolecularSurfaceRepresentationProvider } from 'mol-repr/structure/representation/molecular-surface';
 import { BallAndStickRepresentationProvider } from 'mol-repr/structure/representation/ball-and-stick';
 import { GaussianSurfaceRepresentationProvider } from 'mol-repr/structure/representation/gaussian-surface';
-import { computeModelDSSP } from 'mol-model-props/computed/secondary-structure/dssp';
+import { ComputedSecondaryStructure } from 'mol-model-props/computed/secondary-structure';
 
 const parent = document.getElementById('app')!
 parent.style.width = '100%'
@@ -79,11 +79,10 @@ function getGaussianSurfaceRepr() {
 async function init() {
     const cif = await downloadFromPdb('1crn')
     const models = await getModels(cif)
-    console.time('computeModelDSSP')
-    const secondaryStructure = computeModelDSSP(models[0].atomicHierarchy, models[0].atomicConformation)
-    console.timeEnd('computeModelDSSP');
-    (models[0].properties as any).secondaryStructure = secondaryStructure
     const structure = await getStructure(models[0])
+    console.time('computeDSSP')
+    await ComputedSecondaryStructure.attachFromCifOrCompute(structure)
+    console.timeEnd('computeDSSP');
 
     const show = {
         cartoon: false,