From 8e31f702789ea9f26dbf3e16d5ea7c1fe2f5d622 Mon Sep 17 00:00:00 2001
From: Sebastian Bittrich <bittrich@hs-mittweida.de>
Date: Wed, 20 Mar 2019 11:41:41 -0700
Subject: [PATCH] fixes to DSSP impl, detection of bends and sheets
---
src/mol-math/linear-algebra/3d/vec3.ts | 1 +
.../linear-algebra/_spec/vec3.spec.ts | 14 +
.../model/properties/seconday-structure.ts | 3 +-
.../properties/utils/secondary-structure.ts | 615 ++++++++++++++++--
.../utils/secondary-structure.validation | 73 +++
src/mol-model/structure/model/types.ts | 63 +-
src/mol-theme/color/secondary-structure.ts | 6 +-
src/tests/browser/render-structure.ts | 9 +-
8 files changed, 688 insertions(+), 96 deletions(-)
create mode 100644 src/mol-math/linear-algebra/_spec/vec3.spec.ts
create mode 100644 src/mol-model/structure/model/properties/utils/secondary-structure.validation
diff --git a/src/mol-math/linear-algebra/3d/vec3.ts b/src/mol-math/linear-algebra/3d/vec3.ts
index 3e8fd5d5e..59308ddf0 100644
--- a/src/mol-math/linear-algebra/3d/vec3.ts
+++ b/src/mol-math/linear-algebra/3d/vec3.ts
@@ -414,6 +414,7 @@ namespace Vec3 {
}
const angleTempA = zero(), angleTempB = zero();
+ /** Computes the angle between 2 vectors, reports in rad. */
export function angle(a: Vec3, b: Vec3) {
copy(angleTempA, a);
copy(angleTempB, b);
diff --git a/src/mol-math/linear-algebra/_spec/vec3.spec.ts b/src/mol-math/linear-algebra/_spec/vec3.spec.ts
new file mode 100644
index 000000000..bb98228ac
--- /dev/null
+++ b/src/mol-math/linear-algebra/_spec/vec3.spec.ts
@@ -0,0 +1,14 @@
+import { Vec3 } from '../3d'
+
+describe('vec3', () => {
+ const vec1 = [ 1, 2, 3 ] as Vec3;
+ const vec2 = [ 2, 3, 1 ] as Vec3;
+ const orthVec1 = [ 0, 1, 0 ] as Vec3;
+ const orthVec2 = [ 1, 0, 0 ] as Vec3;
+
+ it('angle calculation', () => {
+ expect(Vec3.angle(vec1, vec1) * 360 / (2 * Math.PI)).toBe(0.0);
+ expect(Vec3.angle(orthVec1, orthVec2) * 360 / (2 * Math.PI)).toBe(90.0);
+ expect(Vec3.angle(vec1, vec2)).toBeCloseTo(0.666946);
+ });
+})
\ No newline at end of file
diff --git a/src/mol-model/structure/model/properties/seconday-structure.ts b/src/mol-model/structure/model/properties/seconday-structure.ts
index 49618f007..45f9a41cf 100644
--- a/src/mol-model/structure/model/properties/seconday-structure.ts
+++ b/src/mol-model/structure/model/properties/seconday-structure.ts
@@ -20,7 +20,8 @@ namespace SecondaryStructure {
export type Element = None | Helix | Sheet
export interface None {
- kind: 'none'
+ kind: 'none',
+ flags?: SecondaryStructureType // TODO should this level of detail be added to non-defined secondary structure elements?
}
export interface Helix {
diff --git a/src/mol-model/structure/model/properties/utils/secondary-structure.ts b/src/mol-model/structure/model/properties/utils/secondary-structure.ts
index afa819947..e21e4e24b 100644
--- a/src/mol-model/structure/model/properties/utils/secondary-structure.ts
+++ b/src/mol-model/structure/model/properties/utils/secondary-structure.ts
@@ -15,12 +15,25 @@ import { BitFlags } from 'mol-util';
import { ElementIndex } from 'mol-model/structure/model/indexing';
import { AtomicHierarchy, AtomicConformation } from '../atomic';
-export function computeSecondaryStructure(hierarchy: AtomicHierarchy, conformation: AtomicConformation): SecondaryStructure {
+/**
+ * TODO bugs to fix:
+ * - some turns are not detected correctly: see e.g. pdb:1acj - maybe more than 2 hbonds require some residue to donate electrons
+ * - some sheets are not extended correctly: see e.g. pdb:1acj
+ * - validate new helix definition
+ * - validate new ordering of secondary structure elements
+ */
+
+export function computeSecondaryStructure(hierarchy: AtomicHierarchy,
+ conformation: AtomicConformation): SecondaryStructure {
// TODO use Zhang-Skolnik for CA alpha only parts or for coarse parts with per-residue elements
return computeModelDSSP(hierarchy, conformation)
}
-export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: AtomicConformation) {
+export function computeModelDSSP(hierarchy: AtomicHierarchy,
+ conformation: AtomicConformation,
+ oldDefinition = true,
+ oldOrdering = true) {
+ // console.log(`calculating secondary structure elements using ${ oldDefinition ? 'old' : 'revised'} definition and ${ oldOrdering ? 'old' : 'revised'} ordering of secondary structure elements`)
const { lookup3d, proteinResidues } = calcAtomicTraceLookup3D(hierarchy, conformation)
const backboneIndices = calcBackboneAtomIndices(hierarchy, proteinResidues)
const hbonds = calcBackboneHbonds(hierarchy, conformation, proteinResidues, backboneIndices, lookup3d)
@@ -28,42 +41,153 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi
const residueCount = proteinResidues.length
const flags = new Uint32Array(residueCount)
+ const torsionAngles = calculateDihedralAngles(hierarchy, conformation, proteinResidues, backboneIndices)
+
+ const ladders: Ladder[] = []
+ const bridges: Bridge[] = []
+
+ const getResidueFlag = oldOrdering ? getOriginalResidueFlag : getUpdatedResidueFlag
+ const getFlagName = oldOrdering ? getOriginalFlagName : getUpdatedFlagName
+
const ctx: DSSPContext = {
+ oldDefinition,
+ oldOrdering,
+ getResidueFlag,
+ getFlagName,
+
hierarchy,
proteinResidues,
flags,
- hbonds
+ hbonds,
+
+ torsionAngles,
+ backboneIndices,
+ conformation,
+ ladders,
+ bridges
}
- assignBends(ctx)
assignTurns(ctx)
assignHelices(ctx)
+ assignBends(ctx)
assignBridges(ctx)
assignLadders(ctx)
assignSheets(ctx)
- const assignment = getDSSPAssignment(flags)
+ const assignment = getDSSPAssignment(flags, getResidueFlag)
const type = new Uint32Array(hierarchy.residues._rowCount) as unknown as SecondaryStructureType[]
+
+ const keys: number[] = []
+ const elements: SecondaryStructure.Element[] = []
+
for (let i = 0, il = proteinResidues.length; i < il; ++i) {
- type[proteinResidues[i]] = assignment[i]
+ const assign = assignment[i]
+ type[proteinResidues[i]] = assign
+ const flag = getResidueFlag(flags[i])
+ // const kind = mapToKind(assign)
+ // TODO is this expected behavior? elements will be strictly split depending on 'winning' flag
+ if (elements.length === 0 || // check ought to fail at very start
+ // elements[elements.length - 1].kind !== kind || // new element if overall kind changed
+ flag !== (elements[elements.length - 1] as SecondaryStructure.Helix | SecondaryStructure.Sheet).flags) { // exact flag changed
+ elements[elements.length] = createElement(mapToKind(assign), flags[i], getResidueFlag)
+ }
+ keys[i] = elements.length - 1
}
const secondaryStructure: SecondaryStructure = {
type,
- key: [], // TODO
- elements: [] // TODO
+ key: keys,
+ elements: elements
}
+
+ // console.log(keys)
+ // console.log(elements)
+ // printDSSPString(flags, getFlagName)
+
return secondaryStructure
}
+// function printDSSPString(flags: Uint32Array, getFlagName: (f: DSSPType) => String) {
+// let out = ''
+// for (let i = 0, il = flags.length; i < il; ++i) {
+// const f = DSSPType.create(flags[i])
+// out += getFlagName(f)
+// }
+// console.log(out)
+// }
+
+function createElement(kind: string, flag: DSSPType.Flag, getResidueFlag: (f: DSSPType) => SecondaryStructureType): SecondaryStructure.Element {
+ // TODO would be nice to add more detailed information
+ if (kind === 'helix') {
+ return {
+ kind: 'helix',
+ flags: getResidueFlag(flag)
+ } as SecondaryStructure.Helix
+ } else if (kind === 'sheet') {
+ return {
+ kind: 'sheet',
+ flags: getResidueFlag(flag)
+ } as SecondaryStructure.Sheet
+ } else {
+ return {
+ kind: 'none',
+ flags: getResidueFlag(flag)
+ }
+ }
+}
+
+function mapToKind(assignment: SecondaryStructureType.Flag) {
+ if (assignment === SecondaryStructureType.SecondaryStructureDssp.H || assignment === SecondaryStructureType.SecondaryStructureDssp.G ||
+ assignment === SecondaryStructureType.SecondaryStructureDssp.I) {
+ return 'helix'
+ } else if (assignment === SecondaryStructureType.SecondaryStructureDssp.B || assignment === SecondaryStructureType.SecondaryStructureDssp.E) {
+ return 'sheet'
+ } else {
+ return 'none'
+ }
+}
+
interface DSSPContext {
+ /** Whether to use the old DSSP convention for the annotation of turns and helices, causes them to be two residues shorter */
+ oldDefinition: boolean,
+ /** alpha-helices are preferred over 3-10 helix */
+ oldOrdering: boolean,
+ getResidueFlag: (f: DSSPType) => SecondaryStructureType,
+ getFlagName: (f: DSSPType) => String,
+
hierarchy: AtomicHierarchy
proteinResidues: SortedArray<ResidueIndex>
/** flags for each residue */
flags: Uint32Array
+ hbonds: DsspHbonds,
- hbonds: DsspHbonds
+ torsionAngles: { phi: number[], psi: number[], omega: number[] },
+ backboneIndices: BackboneAtomIndices,
+ conformation: AtomicConformation,
+ ladders: Ladder[],
+ bridges: Bridge[]
+}
+
+interface Ladder {
+ previousLadder: number,
+ nextLadder: number,
+ firstStart: number,
+ secondStart: number,
+ secondEnd: number,
+ firstEnd: number,
+ type: BridgeType
+}
+
+const enum BridgeType {
+ PARALLEL = 0x0,
+ ANTI_PARALLEL = 0x1
+}
+
+interface Bridge {
+ partner1: number,
+ partner2: number,
+ type: BridgeType
}
type DSSPType = BitFlags<DSSPType.Flag>
@@ -82,6 +206,9 @@ namespace DSSPType {
T3 = 0x80,
T4 = 0x100,
T5 = 0x200,
+ T3S = 0x400, // marks 3-turn start
+ T4S = 0x800,
+ T5S = 0x1000
}
}
@@ -89,7 +216,7 @@ namespace DSSPType {
const caMaxDist = 7.0;
/** min distance between two C-alpha atoms to check for hbond */
-const caMinDist = 4.0;
+const caMinDist = 0.5;
function calcAtomicTraceLookup3D(hierarchy: AtomicHierarchy, conformation: AtomicConformation) {
const { x, y, z } = conformation;
@@ -141,6 +268,150 @@ function calcBackboneAtomIndices(hierarchy: AtomicHierarchy, proteinResidues: So
type DsspHbonds = IntAdjacencyGraph<{ readonly energies: ArrayLike<number> }>
+/**
+ * Bend(i) =: [angle ((CW - Ca(i - 2)),(C"(i + 2) - C"(i))) > 70"]
+ *
+ * Type: S
+ */
+function assignBends(ctx: DSSPContext) {
+ const flags = ctx.flags
+ const { x, y, z } = ctx.conformation
+ const { traceElementIndex } = ctx.hierarchy.derived.residue
+
+ const proteinResidues = ctx.proteinResidues
+ const residueCount = proteinResidues.length
+
+ const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
+
+ const caPosPrev2 = Vec3.zero()
+ const caPos = Vec3.zero()
+ const caPosNext2 = Vec3.zero()
+ // let bends = 0
+
+ const nIndices = ctx.backboneIndices.nIndices
+ const cPos = Vec3.zero()
+ const nPosNext = Vec3.zero()
+
+ f1: for (let i = 2; i < residueCount - 2; i++) {
+ // 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)
+ 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 caAtomPrev2 = traceElementIndex[oRIprev2]
+ const caAtom = traceElementIndex[oRI]
+ const caAtomNext2 = traceElementIndex[oRInext2]
+
+ position(caAtomPrev2, caPosPrev2)
+ position(caAtom, caPos)
+ position(caAtomNext2, caPosNext2)
+
+ const caMinus2 = Vec3.zero()
+ const caPlus2 = Vec3.zero()
+
+ Vec3.sub(caMinus2, caPosPrev2, caPos)
+ Vec3.sub(caPlus2, caPos, caPosNext2)
+
+ const angle = Vec3.angle(caMinus2, caPlus2) * 360 / (2 * Math.PI)
+ if (angle && angle > 70.00) {
+ // console.log(`found bend at ${ i } with kappa ${ angle }`)
+ flags[i] |= DSSPType.Flag.S
+ // bends++
+ }
+ }
+ // console.log(bends + ' bends')
+}
+
+function calculateDihedralAngles(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices): { phi: number[], psi: number[], omega: number[] } {
+ const { cIndices, nIndices } = backboneIndices
+ const { x, y, z } = conformation
+ const { traceElementIndex } = hierarchy.derived.residue
+
+ const residueCount = proteinResidues.length
+ const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
+
+ const cPosPrev = Vec3.zero()
+ const caPosPrev = Vec3.zero()
+ const nPosPrev = Vec3.zero()
+ const cPos = Vec3.zero()
+ const caPos = Vec3.zero()
+ const nPos = Vec3.zero()
+ const cPosNext = Vec3.zero()
+ const caPosNext = Vec3.zero()
+ const nPosNext = Vec3.zero()
+
+ const phi: number[] = []
+ const psi: number[] = []
+ const omega: number[] = []
+
+ for (let i = 0; i < residueCount - 1; ++i) {
+ const oPIprev = i - 1
+ const oRIprev = proteinResidues[i - 1]
+ const oPI = i
+ const oRI = proteinResidues[i]
+ const oPInext = i + 1
+ const oRInext = proteinResidues[i + 1]
+
+ const cAtomPrev = cIndices[oPIprev]
+ const caAtomPrev = traceElementIndex[oRIprev]
+ const nAtomPrev = nIndices[oPIprev]
+ const cAtom = cIndices[oPI]
+ const caAtom = traceElementIndex[oRI]
+ const nAtom = nIndices[oPI]
+ const cAtomNext = cIndices[oPInext]
+ const caAtomNext = traceElementIndex[oRInext]
+ const nAtomNext = nIndices[oRInext]
+
+ // ignore C-terminal residue as acceptor
+ // if (index.findAtomOnResidue(oRI, 'OXT') !== -1) continue
+
+ position(cAtomPrev, cPosPrev)
+ position(caAtomPrev, caPosPrev)
+ position(nAtomPrev, nPosPrev)
+ position(cAtom, cPos)
+ position(caAtom, caPos)
+ position(nAtom, nPos)
+ position(cAtomNext, cPosNext)
+ position(caAtomNext, caPosNext)
+ position(nAtomNext, nPosNext)
+
+ const tPhi = calculatePhi(cPosPrev, nPos, caPos, cPos)
+ const tPsi = calculatePsi(nPos, caPos, cPos, nPosNext)
+ const tOmega = calculateOmega(caPos, cPos, nPosNext, caPosNext)
+
+ // console.log(`phi: ${ tPhi }, psi: ${ tPsi }, omega: ${ tOmega }`)
+
+ phi[phi.length] = tPhi
+ psi[psi.length] = tPsi
+ omega[omega.length] = tOmega
+ }
+
+ return { phi, psi, omega };
+}
+
+// angle calculations return NaN upon missing atoms
+
+function calculatePhi(c: Vec3, nNext: Vec3, caNext: Vec3, cNext: Vec3) {
+ return angle(c, nNext, caNext, cNext);
+}
+
+function calculatePsi(n: Vec3, ca: Vec3, c: Vec3, nNext: Vec3) {
+ return angle(n, ca, c, nNext);
+}
+
+function calculateOmega(ca: Vec3, c: Vec3, nNext: Vec3, caNext: Vec3) {
+ return angle(ca, c, nNext, caNext);
+}
+
function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices, lookup3d: GridLookup3D): DsspHbonds {
const { cIndices, hIndices, nIndices, oIndices } = backboneIndices
const { index } = hierarchy
@@ -163,7 +434,7 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf
const cPosPrev = Vec3.zero()
const oPosPrev = Vec3.zero()
- const caMinDistSq = caMinDist * caMinDist
+ // const caMinDistSq = caMinDist * caMinDist
for (let i = 0, il = proteinResidues.length; i < il; ++i) {
const oPI = i
@@ -183,10 +454,10 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf
position(cAtom, cPos)
position(caAtom, caPos)
- const { indices, count, squaredDistances } = lookup3d.find(caPos[0], caPos[1], caPos[2], caMaxDist)
+ const { indices, count } = lookup3d.find(caPos[0], caPos[1], caPos[2], caMaxDist)
for (let j = 0; j < count; ++j) {
- if (squaredDistances[j] < caMinDistSq) continue
+ // if (squaredDistances[j] < caMinDistSq) continue
const nPI = indices[j]
@@ -221,6 +492,7 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf
const e = calcHbondEnergy(oPos, cPos, nPos, hPos)
if (e > hbondEnergyCutoff) continue
+ // console.log(`detected hbond between ${ oPI } and ${ nPI } with energy ${ e }`)
oAtomResidues[oAtomResidues.length] = oPI
nAtomResidues[nAtomResidues.length] = nPI
energies[energies.length] = e
@@ -245,8 +517,8 @@ function buildHbondGraph(residueCount: number, oAtomResidues: number[], nAtomRes
/** Original priority: H,B,E,G,I,T,S */
function getOriginalResidueFlag(f: DSSPType) {
if (DSSPType.is(f, DSSPType.Flag.H)) return SecondaryStructureType.SecondaryStructureDssp.H
- if (DSSPType.is(f, DSSPType.Flag.B)) return SecondaryStructureType.SecondaryStructureDssp.B
if (DSSPType.is(f, DSSPType.Flag.E)) return SecondaryStructureType.SecondaryStructureDssp.E
+ if (DSSPType.is(f, DSSPType.Flag.B)) return SecondaryStructureType.SecondaryStructureDssp.B
if (DSSPType.is(f, DSSPType.Flag.G)) return SecondaryStructureType.SecondaryStructureDssp.G
if (DSSPType.is(f, DSSPType.Flag.I)) return SecondaryStructureType.SecondaryStructureDssp.I
if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T
@@ -254,37 +526,47 @@ function getOriginalResidueFlag(f: DSSPType) {
return SecondaryStructureType.Flag.None
}
+function getOriginalFlagName(f: DSSPType) {
+ if (DSSPType.is(f, DSSPType.Flag.H)) return 'H'
+ if (DSSPType.is(f, DSSPType.Flag.E)) return 'E'
+ if (DSSPType.is(f, DSSPType.Flag.B)) return 'B'
+ if (DSSPType.is(f, DSSPType.Flag.G)) return 'G'
+ if (DSSPType.is(f, DSSPType.Flag.I)) return 'I'
+ if (DSSPType.is(f, DSSPType.Flag.T)) return 'T'
+ if (DSSPType.is(f, DSSPType.Flag.S)) return 'S'
+ return '-'
+}
+
/** Version 2.1.0 priority: I,H,B,E,G,T,S */
function getUpdatedResidueFlag(f: DSSPType) {
if (DSSPType.is(f, DSSPType.Flag.I)) return SecondaryStructureType.SecondaryStructureDssp.I
if (DSSPType.is(f, DSSPType.Flag.H)) return SecondaryStructureType.SecondaryStructureDssp.H
- if (DSSPType.is(f, DSSPType.Flag.B)) return SecondaryStructureType.SecondaryStructureDssp.B
if (DSSPType.is(f, DSSPType.Flag.E)) return SecondaryStructureType.SecondaryStructureDssp.E
+ if (DSSPType.is(f, DSSPType.Flag.B)) return SecondaryStructureType.SecondaryStructureDssp.B
if (DSSPType.is(f, DSSPType.Flag.G)) return SecondaryStructureType.SecondaryStructureDssp.G
if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T
if (DSSPType.is(f, DSSPType.Flag.S)) return SecondaryStructureType.SecondaryStructureDssp.S
return SecondaryStructureType.Flag.None
}
-// function geFlagName(f: DSSPType) {
-// if (DSSPType.is(f, DSSPType.Flag.I)) return 'I'
-// if (DSSPType.is(f, DSSPType.Flag.H)) return 'H'
-// if (DSSPType.is(f, DSSPType.Flag.B)) return 'B'
-// if (DSSPType.is(f, DSSPType.Flag.E)) return 'E'
-// if (DSSPType.is(f, DSSPType.Flag.G)) return 'G'
-// if (DSSPType.is(f, DSSPType.Flag.T)) return 'T'
-// if (DSSPType.is(f, DSSPType.Flag.S)) return 'S'
-// return '-'
-// }
+function getUpdatedFlagName(f: DSSPType) {
+ if (DSSPType.is(f, DSSPType.Flag.I)) return 'I'
+ if (DSSPType.is(f, DSSPType.Flag.H)) return 'H'
+ if (DSSPType.is(f, DSSPType.Flag.E)) return 'E'
+ if (DSSPType.is(f, DSSPType.Flag.B)) return 'B'
+ if (DSSPType.is(f, DSSPType.Flag.G)) return 'G'
+ if (DSSPType.is(f, DSSPType.Flag.T)) return 'T'
+ if (DSSPType.is(f, DSSPType.Flag.S)) return 'S'
+ return '-'
+}
-function getDSSPAssignment(flags: Uint32Array, useOriginal = false) {
- const getResidueFlag = useOriginal ? getOriginalResidueFlag : getUpdatedResidueFlag
+function getDSSPAssignment(flags: Uint32Array, getResidueFlag: (f: DSSPType) => SecondaryStructureType) {
const type = new Uint32Array(flags.length)
for (let i = 0, il = flags.length; i < il; ++i) {
const f = DSSPType.create(flags[i])
- // console.log(i, geFlagName(f))
type[i] = getResidueFlag(f)
}
+
return type as unknown as ArrayLike<SecondaryStructureType>
}
@@ -303,6 +585,8 @@ const Q = -27.888
/** cutoff for hbonds in kcal/mol, must be lower to be consider as an hbond */
const hbondEnergyCutoff = -0.5
+const hbondEnergyMinimal = -9.9
+
/**
* E = Q * (1/r(ON) + l/r(CH) - l/r(OH) - l/r(CN))
*/
@@ -312,9 +596,19 @@ function calcHbondEnergy(oPos: Vec3, cPos: Vec3, nPos: Vec3, hPos: Vec3) {
const distCN = Vec3.distance(cPos, nPos)
const distON = Vec3.distance(oPos, nPos)
+ if (distOH < caMinDist || distCH < caMinDist || distCN < caMinDist || distON < caMinDist) {
+ return hbondEnergyMinimal
+ }
+
const e1 = Q / distOH - Q / distCH
const e2 = Q / distCN - Q / distON
- return e1 + e2
+ const e = e1 + e2
+
+ // cap lowest possible energy
+ if (e < hbondEnergyMinimal)
+ return hbondEnergyMinimal
+
+ return e
}
/**
@@ -328,28 +622,43 @@ function assignTurns(ctx: DSSPContext) {
const { proteinResidues, hbonds, flags, hierarchy } = ctx
const { chains, residueAtomSegments, chainAtomSegments } = hierarchy
const { label_asym_id } = chains
+ // let turns = 0, turnOpenings = 0
- const turnFlag = [0, 0, 0, DSSPType.Flag.T3, DSSPType.Flag.T4, DSSPType.Flag.T5]
-
- for (let i = 0, il = proteinResidues.length; i < il; ++i) {
- const rI = proteinResidues[i]
- const cI = chainAtomSegments.index[residueAtomSegments.offsets[rI]]
+ const turnFlag = [DSSPType.Flag.T3S, DSSPType.Flag.T4S, DSSPType.Flag.T5S, DSSPType.Flag.T3, DSSPType.Flag.T4, DSSPType.Flag.T5]
- // TODO should take sequence gaps into account
- for (let k = 3; k <= 5; ++k) {
- if (i + k >= proteinResidues.length) continue
+ 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]]
- const rN = proteinResidues[i + k]
+ // 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
// check if hbond exists
- if (hbonds.getDirectedEdgeIndex(i, i + k) !== -1) {
- flags[i] |= turnFlag[k] | DSSPType.Flag.T
+ // console.log(`${ i } has ${ idx + 3}-turn opening: ${ hbonds.getDirectedEdgeIndex(i, i + idx + 3) !== -1 }`)
+ if (hbonds.getDirectedEdgeIndex(i, i + idx + 3) !== -1) {
+ // console.log(`${ idx + 3 }-turn opening at ${ i } to ${ i + idx + 3 }`)
+ flags[i] |= turnFlag[idx + 3] | turnFlag[idx]
+ // turnOpenings++
+ if (ctx.oldDefinition) {
+ for (let k = 1; k < idx + 3; ++k) {
+ flags[i + k] |= turnFlag[idx + 3] | DSSPType.Flag.T
+ // turns++
+ }
+ } else {
+ for (let k = 0; k <= idx + 3; ++k) {
+ flags[i + k] |= turnFlag[idx + 3] | DSSPType.Flag.T
+ // turns++
+ }
+ }
}
}
}
+
+ // console.log(`${ turns } turns, ${ turnOpenings } turn openings`)
}
/**
@@ -369,7 +678,7 @@ function assignTurns(ctx: DSSPContext) {
* Type: B
*/
function assignBridges(ctx: DSSPContext) {
- const { proteinResidues, hbonds, flags } = ctx
+ const { proteinResidues, hbonds, flags, bridges } = ctx
const { offset, b } = hbonds
let i: number, j: number
@@ -385,6 +694,8 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j, i + 1) !== -1) {
flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B
+ // TODO move to constructor, actually omit object all together
+ bridges[bridges.length] = { partner1: Math.min(i, j), partner2: Math.max(i, j), type: BridgeType.PARALLEL }
}
// Parallel Bridge(i, j) =: [Hbond(j - 1, i) and Hbond(i, j + 1)]
@@ -393,6 +704,7 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i) !== -1) {
flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B
+ bridges[bridges.length] = { partner1: Math.min(j, i), partner2: Math.max(j, i), type: BridgeType.PARALLEL }
}
// Antiparallel Bridge(i, j) =: [Hbond(i, j) and Hbond(j, i)]
@@ -401,6 +713,7 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j, i) !== -1) {
flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B
+ bridges[bridges.length] = { partner1: Math.min(j, i), partner2: Math.max(j, i), type: BridgeType.ANTI_PARALLEL }
}
// Antiparallel Bridge(i, j) =: [Hbond(i - 1, j + 1) and Hbond(j - 1, i + l)]
@@ -409,9 +722,14 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i + 1) !== -1) {
flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B
+ bridges[bridges.length] = { partner1: Math.min(j, i), partner2: Math.max(j, i), type: BridgeType.ANTI_PARALLEL }
}
}
}
+
+ bridges.sort((a, b) => a.partner1 > b.partner1 ? 1 : a.partner1 < b.partner1 ? -1 : 0)
+ // console.log(`${ bridges.length } bridges`)
+ // bridges.forEach(b => console.log(b))
}
/**
@@ -428,21 +746,77 @@ function assignBridges(ctx: DSSPContext) {
function assignHelices(ctx: DSSPContext) {
const { proteinResidues, flags } = ctx
- const turnFlag = [0, 0, 0, DSSPType.Flag.T3, DSSPType.Flag.T4, DSSPType.Flag.T5]
+ 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]
- for (let i = 1, il = proteinResidues.length; i < il; ++i) {
- const fI = DSSPType.create(flags[i])
- const fI1 = DSSPType.create(flags[i - 1])
+ const helixCheckOrder = ctx.oldOrdering ? [4, 3, 5] : [3, 4, 5]
+ // const count = [0, 0, 0, 0, 0, 0]
+ for (let ni = 0; ni < helixCheckOrder.length; ni++) {
+ const n = helixCheckOrder[ni]
+
+ for (let i = 1, il = proteinResidues.length - n; i < il; i++) {
+ const fI = DSSPType.create(flags[i])
+ const fI1 = DSSPType.create(flags[i - 1])
+ const fI2 = DSSPType.create(flags[i + 1])
+
+ // TODO rework to elegant solution which will not break instantly
+ if (ctx.oldOrdering) {
+ if ((n === 3 && (DSSPType.is(fI, DSSPType.Flag.H) || DSSPType.is(fI2, DSSPType.Flag.H)) || // for 3-10 yield to alpha helix
+ (n === 5 && ((DSSPType.is(fI, DSSPType.Flag.H) || DSSPType.is(fI, DSSPType.Flag.G)) || (DSSPType.is(fI2, DSSPType.Flag.H) || DSSPType.is(fI2, DSSPType.Flag.G)))))) { // for pi yield to all other helices
+ // console.log('skipping in favor of more preferable helix already present')
+ continue
+ }
+ } else {
+ if ((n === 4 && (DSSPType.is(fI, DSSPType.Flag.G) || DSSPType.is(fI2, DSSPType.Flag.G)) || // for alpha helix yield to 3-10
+ (n === 5 && ((DSSPType.is(fI, DSSPType.Flag.H) || DSSPType.is(fI, DSSPType.Flag.G)) || (DSSPType.is(fI2, DSSPType.Flag.H) || DSSPType.is(fI2, DSSPType.Flag.G)))))) { // for pi yield to all other helices
+ // console.log('skipping in favor of more preferable helix already present')
+ continue
+ }
+ }
- for (let k = 3; k <= 5; ++k) {
- if (DSSPType.is(fI, turnFlag[k]) && DSSPType.is(fI1, turnFlag[k])) {
- for (let l = 0; l < k; ++l) {
- flags[i + l] |= helixFlag[k]
+ if (DSSPType.is(fI, turnFlag[n]) && DSSPType.is(fI, turnFlag[n - 3]) && // check fI for turn start of proper type
+ DSSPType.is(fI1, turnFlag[n]) && DSSPType.is(fI1, turnFlag[n - 3])) { // check fI1 accordingly
+ // console.log(`valid hit for ${n} from ${i} to ${i+n}`)
+ if (ctx.oldDefinition) {
+ for (let k = 0; k < n; k++) {
+ flags[i + k] |= helixFlag[n]
+ // count[n]++
+ }
+ } else {
+ for (let k = -1; k <= n; k++) {
+ flags[i + k] |= helixFlag[n]
+ // count[n]++
+ }
}
}
}
}
+ // console.log(`${ count[3] + count[4] + count[5] } helix elements - ${ count[3] } 3-10, ${ count[4] } alpha, ${ count[5] } pi`)
+}
+
+function angle(a: Vec3, b: Vec3, c: Vec3, d: Vec3) {
+ const ab = Vec3.zero();
+ const cb = Vec3.zero();
+ const bc = Vec3.zero();
+ const dc = Vec3.zero();
+
+ Vec3.sub(ab, a, b);
+ Vec3.sub(cb, c, b);
+ Vec3.sub(bc, b, c);
+ Vec3.sub(dc, d, c);
+
+ const abc = Vec3.zero();
+ const bcd = Vec3.zero();
+
+ Vec3.cross(abc, ab, cb);
+ Vec3.cross(bcd, bc, dc);
+
+ const angle = Vec3.angle(abc, bcd) * 360.0 / (2 * Math.PI);
+ const cross = Vec3.zero();
+ Vec3.cross(cross, abc, bcd);
+ const value = Vec3.dot(cb, cross);
+
+ return value > 0 ? angle : -angle;
}
/**
@@ -451,23 +825,146 @@ function assignHelices(ctx: DSSPContext) {
* Type: E
*/
function assignLadders(ctx: DSSPContext) {
- // TODO
+ const { bridges, ladders } = ctx
+
+ // let ext = 0
+ // create ladders
+ for (let bridgeIndex = 0; bridgeIndex < bridges.length; bridgeIndex++) {
+ const bridge = bridges[bridgeIndex]
+ let found = false
+ for (let ladderIndex = 0; ladderIndex < ladders.length; ladderIndex++) {
+ const ladder = ladders[ladderIndex]
+ if (shouldExtendLadder(ladder, bridge)) {
+ // ext++
+ found = true
+ ladder.firstEnd++
+ if (bridge.type === BridgeType.PARALLEL) {
+ ladder.secondEnd++
+ } else {
+ ladder.secondStart--
+ }
+ }
+ }
+
+ // no suitable assignment: create new ladder with single bridge as content
+ if (!found) {
+ ladders[ladders.length] = {
+ previousLadder: 0,
+ nextLadder: 0,
+ firstStart: bridge.partner1,
+ firstEnd: bridge.partner1,
+ secondStart: bridge.partner2,
+ secondEnd: bridge.partner2,
+ type: bridge.type
+ }
+ }
+ }
+ // console.log(`${ ext } extension operations`)
+
+ // connect ladders
+ for (let ladderIndex1 = 0; ladderIndex1 < ladders.length; ladderIndex1++) {
+ const ladder1 = ladders[ladderIndex1]
+ for (let ladderIndex2 = ladderIndex1; ladderIndex2 < ladders.length; ladderIndex2++) {
+ const ladder2 = ladders[ladderIndex2]
+ if (resemblesBulge(ladder1, ladder2)) {
+ // console.log(`connecting ladders ${ ladderIndex1 } and ${ ladderIndex2 } with inbetween bulge`)
+ ladder1.nextLadder = ladderIndex2
+ ladder2.previousLadder = ladderIndex1
+ }
+ }
+ }
+
+ // console.log(`${ ladders.length } ladders`)
}
/**
- * sheet=: set of one or more ladders connected by shared residues
- *
- * Type: E
+ * For beta structures, we define: a bulge-linked ladder consists of two ladders or bridges of the same type
+ * connected by at most one extra residue of one strand and at most four extra residues on the other strand,
+ * all residues in bulge-linked ladders are marked E, including any extra residues.
*/
-function assignSheets(ctx: DSSPContext) {
- // TODO
+function resemblesBulge(ladder1: Ladder, ladder2: Ladder) {
+ if (!(ladder1.type === ladder2.type && ladder2.firstStart - ladder1.firstEnd < 6 &&
+ ladder1.firstStart < ladder2.firstStart && ladder2.nextLadder === 0)) return false
+
+ if (ladder1.type === BridgeType.PARALLEL) {
+ return bulgeCriterion2(ladder1, ladder2)
+ } else {
+ return bulgeCriterion2(ladder2, ladder1)
+ }
+}
+
+function bulgeCriterion2(ladder1: Ladder, ladder2: Ladder) {
+ return ladder2.secondStart - ladder1.secondEnd > 0 && ((ladder2.secondStart - ladder1.secondEnd < 6 &&
+ ladder2.firstStart - ladder1.firstEnd < 3) || ladder2.secondStart - ladder1.secondEnd < 3)
+}
+
+function shouldExtendLadder(ladder: Ladder, bridge: Bridge): boolean {
+ // in order to extend ladders, same type must be present
+ if (bridge.type !== ladder.type) return false
+
+ // console.log(`${ bridge.partner1 }-${ bridge.partner2 } ${ ladder.firstEnd }`)
+ // only extend if residue 1 is sequence neighbor with regard to ladder
+ if (bridge.partner1 !== ladder.firstEnd + 1) return false
+
+ if (bridge.type === BridgeType.PARALLEL) {
+ if (bridge.partner2 === ladder.secondEnd + 1) {
+ // console.log(ladder)
+ return true
+ }
+ } else {
+ if (bridge.partner2 === ladder.secondStart - 1) {
+ // console.log(ladder)
+ return true
+ }
+ }
+
+ return false
+}
+
+function isHelixType(f: DSSPType) {
+ return DSSPType.is(f, DSSPType.Flag.G) || DSSPType.is(f, DSSPType.Flag.H) || DSSPType.is(f, DSSPType.Flag.I)
}
/**
- * Bend(i) =: [angle ((CW - Ca(i - 2)),(C"(i + 2) - C"(i))) > 70"]
+ * sheet=: set of one or more ladders connected by shared residues
*
- * Type: S
+ * Type: E
*/
-function assignBends(ctx: DSSPContext) {
- // TODO
+function assignSheets(ctx: DSSPContext) {
+ const { ladders, flags } = ctx
+ for (let ladderIndex = 0; ladderIndex < ladders.length; ladderIndex++) {
+ const ladder = ladders[ladderIndex]
+ for (let lcount = ladder.firstStart; lcount <= ladder.firstEnd; lcount++) {
+ const diff = ladder.firstStart - lcount
+ const l2count = ladder.secondStart - diff
+
+ if (ladder.firstStart !== ladder.firstEnd) {
+ flags[lcount] |= DSSPType.Flag.E
+ flags[l2count] |= DSSPType.Flag.E
+ } else {
+ if (!isHelixType(flags[lcount]) && DSSPType.is(flags[lcount], DSSPType.Flag.E)) {
+ flags[lcount] |= DSSPType.Flag.B
+ }
+ if (!isHelixType(flags[l2count]) && DSSPType.is(flags[l2count], DSSPType.Flag.E)) {
+ flags[l2count] |= DSSPType.Flag.B
+ }
+ }
+ }
+
+ if (ladder.nextLadder === 0) continue
+
+ const conladder = ladders[ladder.nextLadder]
+ for (let lcount = ladder.firstStart; lcount <= conladder.firstEnd; lcount++) {
+ flags[lcount] |= DSSPType.Flag.E
+ }
+ if (ladder.type === BridgeType.PARALLEL) {
+ for (let lcount = ladder.secondStart; lcount <= conladder.secondEnd; lcount++) {
+ flags[lcount] |= DSSPType.Flag.E
+ }
+ } else {
+ for (let lcount = conladder.secondEnd; lcount <= ladder.secondStart; lcount++) {
+ flags[lcount] |= DSSPType.Flag.E
+ }
+ }
+ }
}
\ No newline at end of file
diff --git a/src/mol-model/structure/model/properties/utils/secondary-structure.validation b/src/mol-model/structure/model/properties/utils/secondary-structure.validation
new file mode 100644
index 000000000..9b76268f3
--- /dev/null
+++ b/src/mol-model/structure/model/properties/utils/secondary-structure.validation
@@ -0,0 +1,73 @@
+compares Mol* port of DSSP (with default parameters) to the BioJava implementation
+
+### pdb:1pga ###
+# turns #
+Mol*: ----------------------TTTTTTTTTTTTTTTT--------TTTT------
+53 turns, 18 openings
+DSSP: ----------------------TTTTTTTTTTTTTTTT--------TTTT------
+53 turns, 18 openings
+
+# bends #
+Mol*: ---------SS---------SSSSSSSSSSSSSSSSSS--------SSS-------
+23 bends
+DSSP: ---------SS---------SSSSSSSSSSSSSSSSSS--------SSS-------
+23 bends
+
+# helices #
+Mol*: ----------------------HHHHHHHHHHHHHHTT--------TTTT------
+44 helix elements - 0 3-10, 44 alpha, 0 pi
+DSSP: ----------------------HHHHHHHHHHHHHHTT--------TTTT------
+44 helix elements - 0 3-10, 44 alpha, 0 pi
+
+# all #
+Mol*: -EEEEEEE-SS-EEEEEEE-SSHHHHHHHHHHHHHHTT---EEEEETTTTEEEEE-
+DSSP: -EEEEEEE-SS-EEEEEEE-SSHHHHHHHHHHHHHHTT---EEEEETTTTEEEEE-
+
+
+### pdb:1bta ###
+# turns #
+Mol*: ------TTT---TTTTTTTTTTTTT--TT----TTTTTTTTTTT-----------TTTTTTTTT--TTTTTTTTTTTTTTT--------
+127 turns, 44 openings
+DSSP: ------TTT---TTTTTTTTTTTTT--TT----TTTTTTTTTTT-----------TTTTTTTTT--TTTTTTTTTTTTTTT--------
+127 turns, 44 openings
+
+# bends #
+Mol*: ------SSS--SSSSSSSSSSSSS---SS---SSSSSSSSSSSSS-SS------SSSSSSSSSSSSSSSSSSSSSSSSSSS--------
+60 bends
+DSSP: ------SSS--SSSSSSSSSSSSS---SS---SSSSSSSSSSSSS-SS------SSSSSSSSSSSSSSSSSSSSSSSSSSS--------
+60 bends
+
+# helices #
+Mol*: ------TTT---HHHHHHHHHHHHT--TT----HHHHHHHHTTT-----------TTHHHHTTT--HHHHHHHHHHHHHTT--------
+100 helix elements - 0 3-10, 100 alpha, 0 pi
+DSSP: ------TTT---HHHHHHHHHHHHT--TT----HHHHHHHHTTT-----------TTHHHHTTT--HHHHHHHHHHHHHTT--------
+100 helix elements - 0 3-10, 100 alpha, 0 pi
+
+# all #
+Mol*: -EEEEETTT--SHHHHHHHHHHHHT--TT---SHHHHHHHHTTTS-SSEEEEEESTTHHHHTTTSSHHHHHHHHHHHHHTT--EEEEE-
+DSSP: -EEEEETTT--SHHHHHHHHHHHHT--TT---SHHHHHHHHTTTS-SSEEEEEESTTHHHHTTTSSHHHHHHHHHHHHHTT--EEEEE-
+
+
+### pdb:1acj ###
+# turns #
+Mol*: -------TT----------TT----------------TTTTT----------------------------TTTT-TTTTTT----------------------------------TTT------TTT-TTTTTTTTT-----------TTTT---TT-------TTTTTTTTTTTTTTTTTTTTT--TT-------TTTTTTTTTTTT-TTTTTT----------TT-TT----TTTTTTTTTTTTTTTT-----TTTTTTTTTT--TTTTTTTTTTT-----------------------TTTTTTTT----------------TTTTTTT-TT--TT------TTTTTTTTTTTTTT--TTTTTTTTTTT--TTTTT-TTTTTTTTTTTTTTTTTTTTTTTTTTTTTTT-------------TT----TTT---TTTTTTTTTTTTT-TTT---TTTTTTTTTTTTTTTTTTTT--------------------------------TTTTTTTTTTTTTTTTTTT-
+614 turns, 223 openings
+DSSP: -------TT----------TT----------------TTTTT------------------------------TT-TTTTTT----------------------------------TTT------TTT-TTTTTTTTT-----------TTTT---TT-------TTTTTTTTTTTTTTTTTTTTT--TT-------TTTTTTTTTTTT-TTTTTT----------TT-TT----TTTTTTTTTTTTTTTT-----TTTTTTTTTT--TTTTTTTTTTT-----------------------TTTTTTTT----------------TTTTTTT-TT--TT------TTTTTTTTTTT-TT--TTTTTTTTTTT--TTTTT-TTTTTTTTTTTTTTTTTTTTTTTTTTTTTTT-------------TT----TTT---TTTTTTTTTTTTT-TTT---TTTTTTTTTTTTTTTTTTTT--------------------------------TTTTTTTTTTTTTTTTTTT-
+606 turns, 220 openings
+
+# bends #
+Mol*: --S----SS----------SSS----------S----SS-SSS--------SS-----S-----------SSSS-SSSSSSS--S---S----------SS--SS---------SSSS---S-SSSS--SSSSSSS-----------SSSS----SS-SSS-S-SSSSSSSSSSSSSSSSSSSSS--SSS------SSSSSSSSSSSS-SSSSSS-S----SS--SSSSSS---SSSSSSSSSSSSSSS----S-SSSSSSSSSSS-SSSSSSSSSS--SS--SS--S------SSSSSS-SSSSSSS--S--S-------S--SSSSSSSSSSS--SSS-----SSSSSSSSSSS-SS--SSSSSSSSSSS--SSSSS-SSSSSSSSSSSSSSSSS----SSSSSSSSSSSS-----------SS--S-SSS-S-SS-SSSSSS--SS-SSS---SSSSSSSSSSSSSSSSSSSSSSSS---------SSS------SSSS----SSSS-S----SSSSSSSSSS--
+305 bends
+DSSP: --S----SS----------SSS----------S----SS-SSS--------SS-----S-----------SSSS-SSSSSSS--S---S----------SS--SS---------SSSS---S-SSSS--SSSSSSS-----------SSSS----SS-SSS-S-SSSSSSSSSSSSSSSSSSSSS--SSS------SSSSSSSSSSSS-SSSSSS-S----SS--SSSSSS---SSSSSSSSSSSSSSS----S-SSSSSSSSSSS-SSSSSSSSSS--SS--SS--S------SSSSSS-SSSSSSS--S--S-------S--SSSSSSSSSSS--SSS-----SSSSSSSSSSS-SS--SSSSSSSSSSS--SSSSS-SSSSSSSSSSSSSSSSS----SSSSSSSSSSSS-----------SS--S-SSS-S-SS-SSSSSS--SS-SSS---SSSSSSSSSSSSSSSSSSSSSSSS---------SSS------SSSS----SSSS-S----SSSSSSSSSS--
+305 bends
+
+# helices #
+Mol*: -------TT----------TT----------------GGGTT----------------------------TTTT-HHHHTT----------------------------------TTT------GGG-THHHHHHHT-----------HHHH---TT-------HHHHHHHHHHHHHHHHGGGGT--TT-------THHHHHHHHHHH-HHHHTT----------TT-TT----HHHHHHHHHHHHHHTT-----HHHHHHHHHH--HHHHHHHHGGG-----------------------HHHHHHHT----------------HHHHHHH-TT--TT------HHHHHHHHHHHTTT--HHHHHHHHHHH--GGGTT-HHHHHHHHHHHHHHHHTHHHHHHHHHHHHTT-------------TT----GGG---TTTTHHHHTTGGG-GGG---HHHHHHHHHHHHHHHHHHHT--------------------------------TTHHHHHHHHTHHHHHHHH-
+523 helix elements - 27 3-10, 496 alpha, 0 pi
+DSSP: -------TT----------TT----------------GGGTT------------------------------TT-HHHHTT----------------------------------TTT------GGG-THHHHHHHT-----------HHHH---TT-------HHHHHHHHHHHHHHHHGGGGT--TT-------THHHHHHHHHHH-HHHHTT----------TT-TT----HHHHHHHHHHHHHHTT-----HHHHHHHHHH--HHHHHHHHGGG-----------------------HHHHHHHT----------------HHHHHHH-TT--TT------HHHHHHHHHHH-TT--HHHHHHHHHHH--GGGTT-HHHHHHHHHHHHHHHHTHHHHHHHHHHHHTT-------------TT----GGG---TTTTHHHHTTGGG-GGG---HHHHHHHHHHHHHHHHHHHT--------------------------------TTHHHHHHHHTHHHHHHHH-
+523 helix elements - 27 3-10, 496 alpha, 0 pi
+
+# all #
+Mol*: --SEEEETTEEEE-EEEEETTEEEEEEEEEE-EE---GGGTTS--EE----SSEEE--S---B-------TTTT-HHHHTTS--S-B-S---EEEEEE-SS--SSEEEEEEE--STTT---S-SGGG-THHHHHHHT-EEEE-----SHHHH---TT-SSS-S-HHHHHHHHHHHHHHHHGGGGTEEEEEEEEEEETHHHHHHHHHHH-HHHHTT-SEEEEES--TTSTTSEEEHHHHHHHHHHHHHHTT---S-HHHHHHHHHHS-HHHHHHHHGGG-SS--SS--S--EEE-SSSSSS-HHHHHHHT-S--S-EEEEEESB-SHHHHHHHSTT--TTS-----HHHHHHHHHHHTTT--HHHHHHHHHHH--GGGTT-HHHHHHHHHHHHHHHHTHHHHHHHHHHHHTTSS-EEEEEE----TT--S-GGG-SBTTTTHHHHTTGGG-GGG---HHHHHHHHHHHHHHHHHHHTSSSS---------SSS-EEEEESSSS--EEESTTHHHHHHHHTHHHHHHHH-
+DSSP: --SEEEETTEEEE-EEEEETTEEEEEEEEEE-EE---GGGTTS--EE----SSEEE--S---B-------SSTT-HHHHTTS--S-B-S---EEEEEE-SS--SSEEEEEEE--STTT---S-SGGG-THHHHHHHT-EEEE-----SHHHH---TT-SSS-S-HHHHHHHHHHHHHHHHGGGGTEEEEEEEEEEETHHHHHHHHHHH-HHHHTT-SEEEEES--TTSTTS-EEHHHHHHHHHHHHHHTT---S-HHHHHHHHHHS-HHHHHHHHGGG-SS--SS--S---EE-SSSSSS-HHHHHHHT-S--S-EEEEEESB-SHHHHHHHSTT--TTS-----HHHHHHHHHHH-TT--HHHHHHHHHHH--GGGTT-HHHHHHHHHHHHHHHHTHHHHHHHHHHHHTTSS-EEEEEE----TT--S-GGG-SBTTTTHHHHTTGGG-GGG---HHHHHHHHHHHHHHHHHHHTSSSS---------SSS-EEEEESSSS--EEESTTHHHHHHHHTHHHHHHHH-
+TODO fix mismatches e.g. here
\ No newline at end of file
diff --git a/src/mol-model/structure/model/types.ts b/src/mol-model/structure/model/types.ts
index 89293da96..76c0982a1 100644
--- a/src/mol-model/structure/model/types.ts
+++ b/src/mol-model/structure/model/types.ts
@@ -275,41 +275,42 @@ export namespace SecondaryStructureType {
DoubleHelix = 0x1,
Helix = 0x2,
Beta = 0x4,
- Turn = 0x8,
+ Bend = 0x8,
+ Turn = 0x10,
// category variant
- LeftHanded = 0x10, // helix
- RightHanded = 0x20,
+ LeftHanded = 0x20, // helix
+ RightHanded = 0x40,
- ClassicTurn = 0x40, // turn
- InverseTurn = 0x80,
+ ClassicTurn = 0x80, // turn
+ InverseTurn = 0x100,
// sub-category
- HelixOther = 0x100, // protein
- Helix27 = 0x200,
- Helix3Ten = 0x400,
- HelixAlpha = 0x800,
- HelixGamma = 0x1000,
- HelixOmega = 0x2000,
- HelixPi = 0x4000,
- HelixPolyproline = 0x8000,
-
- DoubleHelixOther = 0x10000, // nucleic
- DoubleHelixZ = 0x20000,
- DoubleHelixA = 0x40000,
- DoubleHelixB = 0x80000,
-
- BetaOther = 0x100000, // protein
- BetaStrand = 0x200000, // single strand
- BetaSheet = 0x400000, // multiple hydrogen bonded strands
- BetaBarell = 0x800000, // closed series of sheets
-
- TurnOther = 0x1000000, // protein
- Turn1 = 0x2000000,
- Turn2 = 0x4000000,
- Turn3 = 0x8000000,
-
- NA = 0x10000000, // not applicable/available
+ HelixOther = 0x200, // protein
+ Helix27 = 0x400,
+ Helix3Ten = 0x800,
+ HelixAlpha = 0x1000,
+ HelixGamma = 0x2000,
+ HelixOmega = 0x4000,
+ HelixPi = 0x8000,
+ HelixPolyproline = 0x10000,
+
+ DoubleHelixOther = 0x20000, // nucleic
+ DoubleHelixZ = 0x40000,
+ DoubleHelixA = 0x80000,
+ DoubleHelixB = 0x100000,
+
+ BetaOther = 0x200000, // protein
+ BetaStrand = 0x400000, // single strand
+ BetaSheet = 0x800000, // multiple hydrogen bonded strands
+ BetaBarell = 0x1000000, // closed series of sheets
+
+ TurnOther = 0x2000000, // protein
+ Turn1 = 0x4000000,
+ Turn2 = 0x8000000,
+ Turn3 = 0x10000000,
+
+ NA = 0x20000000, // not applicable/available
}
export const SecondaryStructureMmcif: { [value: string]: number } = {
@@ -386,7 +387,7 @@ export namespace SecondaryStructureType {
G: Flag.Helix | Flag.Helix3Ten, // 3-helix (310 helix)
I: Flag.Helix | Flag.HelixPi, // 5 helix (pi-helix)
T: Flag.Turn, // hydrogen bonded turn
- S: Flag.Turn, // bend
+ S: Flag.Bend, // bend
}
}
diff --git a/src/mol-theme/color/secondary-structure.ts b/src/mol-theme/color/secondary-structure.ts
index 384c5fb0b..cb8878f28 100644
--- a/src/mol-theme/color/secondary-structure.ts
+++ b/src/mol-theme/color/secondary-structure.ts
@@ -22,6 +22,8 @@ const SecondaryStructureColors = ColorMap({
'betaTurn': 0x6080FF,
'betaStrand': 0xFFC800,
'coil': 0xFFFFFF,
+ 'bend': 0x66D8C9 /* biting original color used 0x00FF00 */,
+ 'turn': 0x00B266,
'dna': 0xAE00FE,
'rna': 0xFD0162,
@@ -53,8 +55,10 @@ export function secondaryStructureColor(unit: Unit, element: ElementIndex): Colo
return SecondaryStructureColors.alphaHelix
} else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Beta)) {
return SecondaryStructureColors.betaStrand
+ } else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Bend)) {
+ return SecondaryStructureColors.bend
} else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Turn)) {
- return SecondaryStructureColors.coil
+ return SecondaryStructureColors.turn
} else {
const moleculeType = getElementMoleculeType(unit, element)
if (moleculeType === MoleculeType.DNA) {
diff --git a/src/tests/browser/render-structure.ts b/src/tests/browser/render-structure.ts
index 9dae8a78c..e3732508f 100644
--- a/src/tests/browser/render-structure.ts
+++ b/src/tests/browser/render-structure.ts
@@ -65,9 +65,10 @@ async function init() {
const cif = await downloadFromPdb('3j3q')
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 secondaryStructure = computeModelDSSP(models[0].atomicHierarchy,
+ models[0].atomicConformation)
+ console.timeEnd('computeModelDSSP');
+ (models[0].properties as any).secondaryStructure = secondaryStructure
const structure = await getStructure(models[0])
const cartoonRepr = getCartoonRepr()
@@ -76,7 +77,7 @@ async function init() {
size: reprCtx.sizeThemeRegistry.create('uniform', { structure })
})
await cartoonRepr.createOrUpdate({ ...CartoonRepresentationProvider.defaultValues, quality: 'auto' }, structure).run()
-
+
canvas3d.add(cartoonRepr)
canvas3d.resetCamera()
}
--
GitLab