Skip to content
Snippets Groups Projects
Commit 4e70c822 authored by David Sehnal's avatar David Sehnal
Browse files

Merge branch 'master' of https://github.com/molstar/molstar-proto into plugin

parents e13df4dc 702fc483
Branches
Tags
No related merge requests found
...@@ -414,6 +414,7 @@ namespace Vec3 { ...@@ -414,6 +414,7 @@ namespace Vec3 {
} }
const angleTempA = zero(), angleTempB = zero(); const angleTempA = zero(), angleTempB = zero();
/** Computes the angle between 2 vectors, reports in rad. */
export function angle(a: Vec3, b: Vec3) { export function angle(a: Vec3, b: Vec3) {
copy(angleTempA, a); copy(angleTempA, a);
copy(angleTempB, b); copy(angleTempB, b);
...@@ -433,6 +434,28 @@ namespace Vec3 { ...@@ -433,6 +434,28 @@ namespace Vec3 {
} }
} }
const tmp_dh_ab = Vec3.zero();
const tmp_dh_cb = Vec3.zero();
const tmp_dh_bc = Vec3.zero();
const tmp_dh_dc = Vec3.zero();
const tmp_dh_abc = Vec3.zero();
const tmp_dh_bcd = Vec3.zero();
const tmp_dh_cross = Vec3.zero();
/** Computes the dihedral angles of 4 related atoms. phi: C, N+1, CA+1, C+1 - psi: N, CA, C, N+1 - omega: CA, C, N+1, CA+1 */
export function dihedralAngle(a: Vec3, b: Vec3, c: Vec3, d: Vec3) {
Vec3.sub(tmp_dh_ab, a, b);
Vec3.sub(tmp_dh_cb, c, b);
Vec3.sub(tmp_dh_bc, b, c);
Vec3.sub(tmp_dh_dc, d, c);
Vec3.cross(tmp_dh_abc, tmp_dh_ab, tmp_dh_cb);
Vec3.cross(tmp_dh_bcd, tmp_dh_bc, tmp_dh_dc);
const angle = Vec3.angle(tmp_dh_abc, tmp_dh_bcd) * 360.0 / (2 * Math.PI);
Vec3.cross(tmp_dh_cross, tmp_dh_abc, tmp_dh_bcd);
return Vec3.dot(tmp_dh_cb, tmp_dh_cross) > 0 ? angle : -angle;
}
/** /**
* Returns whether or not the vectors have exactly the same elements in the same position (when compared with ===) * Returns whether or not the vectors have exactly the same elements in the same position (when compared with ===)
*/ */
......
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
...@@ -17,12 +17,17 @@ interface SecondaryStructure { ...@@ -17,12 +17,17 @@ interface SecondaryStructure {
} }
namespace SecondaryStructure { namespace SecondaryStructure {
export type Element = None | Helix | Sheet export type Element = None | Turn | Helix | Sheet
export interface None { export interface None {
kind: 'none' kind: 'none'
} }
export interface Turn {
kind: 'turn',
flags: SecondaryStructureType
}
export interface Helix { export interface Helix {
kind: 'helix', kind: 'helix',
flags: SecondaryStructureType, flags: SecondaryStructureType,
......
...@@ -2,6 +2,7 @@ ...@@ -2,6 +2,7 @@
* Copyright (c) 2019 mol* contributors, licensed under MIT, See LICENSE file for more info. * Copyright (c) 2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
* *
* @author Alexander Rose <alexander.rose@weirdbyte.de> * @author Alexander Rose <alexander.rose@weirdbyte.de>
* @author Sebastian Bittrich <sebastian.bittrich@rcsb.org>
*/ */
import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure'; import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
...@@ -14,13 +15,120 @@ import { IntAdjacencyGraph } from 'mol-math/graph'; ...@@ -14,13 +15,120 @@ import { IntAdjacencyGraph } from 'mol-math/graph';
import { BitFlags } from 'mol-util'; import { BitFlags } from 'mol-util';
import { ElementIndex } from 'mol-model/structure/model/indexing'; import { ElementIndex } from 'mol-model/structure/model/indexing';
import { AtomicHierarchy, AtomicConformation } from '../atomic'; import { AtomicHierarchy, AtomicConformation } from '../atomic';
import { ParamDefinition as PD } from 'mol-util/param-definition'
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
*/
/** max distance between two C-alpha atoms to check for hbond */
const caMaxDist = 9.0;
/**
* Constant for electrostatic energy in kcal/mol
* f * q1 * q2
* Q = -332 * 0.42 * 0.20
*
* f is the dimensional factor
*
* q1 and q2 are partial charges which are placed on the C,O
* (+q1,-q1) and N,H (-q2,+q2)
*/
const Q = -27.888
/** cutoff for hbonds in kcal/mol, must be lower to be consider as an hbond */
const hbondEnergyCutoff = -0.5
/** prevent extremely low hbond energies */
const hbondEnergyMinimal = -9.9
interface DSSPContext {
params: Partial<PD.Values<SecondaryStructureComputationParams>>,
getResidueFlag: (f: DSSPType) => SecondaryStructureType,
getFlagName: (f: DSSPType) => String,
hierarchy: AtomicHierarchy
proteinResidues: SortedArray<ResidueIndex>
/** flags for each residue */
flags: Uint32Array
hbonds: DsspHbonds,
torsionAngles: { phi: Float32Array, psi: Float32Array },
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
}
class Bridge {
partner1: number;
partner2: number;
type: BridgeType;
constructor(p1: number, p2: number, type: BridgeType) {
this.partner1 = Math.min(p1, p2)
this.partner2 = Math.max(p1, p2)
this.type = type
}
}
type DSSPType = BitFlags<DSSPType.Flag>
namespace DSSPType {
export const is: (t: DSSPType, f: Flag) => boolean = BitFlags.has
export const create: (f: Flag) => DSSPType = BitFlags.create
export const enum Flag {
_ = 0x0,
H = 0x1,
B = 0x2,
E = 0x4,
G = 0x8,
I = 0x10,
S = 0x20,
T = 0x40,
T3 = 0x80,
T4 = 0x100,
T5 = 0x200,
T3S = 0x400, // marks 3-turn start
T4S = 0x800,
T5S = 0x1000
}
}
export const SecondaryStructureComputationParams = {
oldDefinition: PD.Boolean(true, { description: 'Whether to use the old DSSP convention for the annotation of turns and helices, causes them to be two residues shorter' }),
oldOrdering: PD.Boolean(true, { description: 'Alpha-helices are preferred over 3-10 helices' })
}
export type SecondaryStructureComputationParams = typeof SecondaryStructureComputationParams
export function computeSecondaryStructure(hierarchy: AtomicHierarchy,
conformation: AtomicConformation) {
// TODO use Zhang-Skolnik for CA alpha only parts or for coarse parts with per-residue elements // TODO use Zhang-Skolnik for CA alpha only parts or for coarse parts with per-residue elements
return computeModelDSSP(hierarchy, conformation) return computeModelDSSP(hierarchy, conformation)
} }
export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: AtomicConformation) { export function computeModelDSSP(hierarchy: AtomicHierarchy,
conformation: AtomicConformation,
params: Partial<PD.Values<SecondaryStructureComputationParams>> = {}): SecondaryStructure {
params = { ...PD.getDefaultValues(SecondaryStructureComputationParams), ...params };
const { lookup3d, proteinResidues } = calcAtomicTraceLookup3D(hierarchy, conformation) const { lookup3d, proteinResidues } = calcAtomicTraceLookup3D(hierarchy, conformation)
const backboneIndices = calcBackboneAtomIndices(hierarchy, proteinResidues) const backboneIndices = calcBackboneAtomIndices(hierarchy, proteinResidues)
const hbonds = calcBackboneHbonds(hierarchy, conformation, proteinResidues, backboneIndices, lookup3d) const hbonds = calcBackboneHbonds(hierarchy, conformation, proteinResidues, backboneIndices, lookup3d)
...@@ -28,68 +136,102 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi ...@@ -28,68 +136,102 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi
const residueCount = proteinResidues.length const residueCount = proteinResidues.length
const flags = new Uint32Array(residueCount) 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 ladders: Ladder[] = []
const bridges: Bridge[] = []
const getResidueFlag = params.oldDefinition ? getOriginalResidueFlag : getUpdatedResidueFlag
const getFlagName = params.oldOrdering ? getOriginalFlagName : getUpdatedFlagName
const ctx: DSSPContext = { const ctx: DSSPContext = {
params,
getResidueFlag,
getFlagName,
hierarchy, hierarchy,
proteinResidues, proteinResidues,
flags, flags,
hbonds hbonds,
torsionAngles,
backboneIndices,
conformation,
ladders,
bridges
} }
assignBends(ctx)
assignTurns(ctx) assignTurns(ctx)
assignHelices(ctx) assignHelices(ctx)
assignBends(ctx)
assignBridges(ctx) assignBridges(ctx)
assignLadders(ctx) assignLadders(ctx)
assignSheets(ctx) assignSheets(ctx)
const assignment = getDSSPAssignment(flags) const assignment = getDSSPAssignment(flags, getResidueFlag)
const type = new Uint32Array(hierarchy.residues._rowCount) as unknown as SecondaryStructureType[] 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) { 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])
// 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 */) {
elements[elements.length] = createElement(mapToKind(assign), flags[i], getResidueFlag)
}
keys[i] = elements.length - 1
} }
const secondaryStructure: SecondaryStructure = { const secondaryStructure: SecondaryStructure = {
type, type,
key: [], // TODO key: keys,
elements: [] // TODO elements: elements
} }
return secondaryStructure return secondaryStructure
} }
interface DSSPContext { function createElement(kind: string, flag: DSSPType.Flag, getResidueFlag: (f: DSSPType) => SecondaryStructureType): SecondaryStructure.Element {
hierarchy: AtomicHierarchy // TODO would be nice to add more detailed information
proteinResidues: SortedArray<ResidueIndex> if (kind === 'helix') {
/** flags for each residue */ return {
flags: Uint32Array kind: 'helix',
flags: getResidueFlag(flag)
hbonds: DsspHbonds } as SecondaryStructure.Helix
} else if (kind === 'sheet') {
return {
kind: 'sheet',
flags: getResidueFlag(flag)
} as SecondaryStructure.Sheet
} else if (kind === 'turn' || kind === 'bend') {
return {
kind: 'turn',
flags: getResidueFlag(flag)
}
} else {
return {
kind: 'none'
} }
type DSSPType = BitFlags<DSSPType.Flag>
namespace DSSPType {
export const is: (t: DSSPType, f: Flag) => boolean = BitFlags.has
export const create: (f: Flag) => DSSPType = BitFlags.create
export const enum Flag {
_ = 0x0,
H = 0x1,
B = 0x2,
E = 0x4,
G = 0x8,
I = 0x10,
S = 0x20,
T = 0x40,
T3 = 0x80,
T4 = 0x100,
T5 = 0x200,
} }
} }
/** max distance between two C-alpha atoms to check for hbond */ function mapToKind(assignment: SecondaryStructureType.Flag) {
const caMaxDist = 7.0; if (assignment === SecondaryStructureType.SecondaryStructureDssp.H || assignment === SecondaryStructureType.SecondaryStructureDssp.G || assignment === SecondaryStructureType.SecondaryStructureDssp.I) {
return 'helix'
/** min distance between two C-alpha atoms to check for hbond */ } else if (assignment === SecondaryStructureType.SecondaryStructureDssp.B || assignment === SecondaryStructureType.SecondaryStructureDssp.E) {
const caMinDist = 4.0; return 'sheet'
} else if (assignment === SecondaryStructureType.SecondaryStructureDssp.T) {
return 'turn'
} else if (assignment === SecondaryStructureType.SecondaryStructureDssp.S) {
return 'bend'
} else {
return 'none'
}
}
function calcAtomicTraceLookup3D(hierarchy: AtomicHierarchy, conformation: AtomicConformation) { function calcAtomicTraceLookup3D(hierarchy: AtomicHierarchy, conformation: AtomicConformation) {
const { x, y, z } = conformation; const { x, y, z } = conformation;
...@@ -141,6 +283,104 @@ function calcBackboneAtomIndices(hierarchy: AtomicHierarchy, proteinResidues: So ...@@ -141,6 +283,104 @@ function calcBackboneAtomIndices(hierarchy: AtomicHierarchy, proteinResidues: So
type DsspHbonds = IntAdjacencyGraph<{ readonly energies: ArrayLike<number> }> 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()
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) {
flags[i] |= DSSPType.Flag.S
}
}
}
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
const residueCount = proteinResidues.length
const position = (i: number, v: Vec3) => Vec3.set(v, x[i], y[i], z[i])
let cPosPrev = Vec3.zero(), caPosPrev = Vec3.zero(), nPosPrev = Vec3.zero()
let cPos = Vec3.zero(), caPos = Vec3.zero(), nPos = Vec3.zero()
let cPosNext = Vec3.zero(), caPosNext = Vec3.zero(), nPosNext = Vec3.zero()
const phi: Float32Array = new Float32Array(residueCount - 1)
const psi: Float32Array = new Float32Array(residueCount - 1)
const cAtomPrev = cIndices[-1], caAtomPrev = traceElementIndex[proteinResidues[-1]], nAtomPrev = nIndices[-1]
position(cAtomPrev, cPosPrev), position(caAtomPrev, caPosPrev), position(nAtomPrev, nPosPrev)
const cAtom = cIndices[0], caAtom = traceElementIndex[proteinResidues[0]], nAtom = nIndices[0]
position(cAtom, cPos), position(caAtom, caPos), position(nAtom, nPos)
const cAtomNext = cIndices[1], caAtomNext = traceElementIndex[proteinResidues[1]], nAtomNext = nIndices[1]
position(cAtomNext, cPosNext), position(caAtomNext, caPosNext), position(nAtomNext, nPosNext)
for (let i = 0; i < residueCount - 1; ++i) {
// ignore C-terminal residue as acceptor
if (index.findAtomOnResidue(proteinResidues[i], 'OXT') !== -1) continue
// returns NaN for missing atoms
phi[i] = Vec3.dihedralAngle(cPosPrev, nPos, caPos, cPos)
psi[i] = Vec3.dihedralAngle(nPos, caPos, cPos, nPosNext)
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)
}
return { phi, psi };
}
function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices, lookup3d: GridLookup3D): DsspHbonds { function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConformation, proteinResidues: SortedArray<ResidueIndex>, backboneIndices: BackboneAtomIndices, lookup3d: GridLookup3D): DsspHbonds {
const { cIndices, hIndices, nIndices, oIndices } = backboneIndices const { cIndices, hIndices, nIndices, oIndices } = backboneIndices
const { index } = hierarchy const { index } = hierarchy
...@@ -163,8 +403,6 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf ...@@ -163,8 +403,6 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf
const cPosPrev = Vec3.zero() const cPosPrev = Vec3.zero()
const oPosPrev = Vec3.zero() const oPosPrev = Vec3.zero()
const caMinDistSq = caMinDist * caMinDist
for (let i = 0, il = proteinResidues.length; i < il; ++i) { for (let i = 0, il = proteinResidues.length; i < il; ++i) {
const oPI = i const oPI = i
const oRI = proteinResidues[i] const oRI = proteinResidues[i]
...@@ -183,11 +421,9 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf ...@@ -183,11 +421,9 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf
position(cAtom, cPos) position(cAtom, cPos)
position(caAtom, caPos) 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) { for (let j = 0; j < count; ++j) {
if (squaredDistances[j] < caMinDistSq) continue
const nPI = indices[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, TODO take chain border into account
...@@ -245,8 +481,8 @@ function buildHbondGraph(residueCount: number, oAtomResidues: number[], nAtomRes ...@@ -245,8 +481,8 @@ function buildHbondGraph(residueCount: number, oAtomResidues: number[], nAtomRes
/** Original priority: H,B,E,G,I,T,S */ /** Original priority: H,B,E,G,I,T,S */
function getOriginalResidueFlag(f: DSSPType) { function getOriginalResidueFlag(f: DSSPType) {
if (DSSPType.is(f, DSSPType.Flag.H)) return SecondaryStructureType.SecondaryStructureDssp.H 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.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.G)) return SecondaryStructureType.SecondaryStructureDssp.G
if (DSSPType.is(f, DSSPType.Flag.I)) return SecondaryStructureType.SecondaryStructureDssp.I if (DSSPType.is(f, DSSPType.Flag.I)) return SecondaryStructureType.SecondaryStructureDssp.I
if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T
...@@ -254,55 +490,50 @@ function getOriginalResidueFlag(f: DSSPType) { ...@@ -254,55 +490,50 @@ function getOriginalResidueFlag(f: DSSPType) {
return SecondaryStructureType.Flag.None 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 */ /** Version 2.1.0 priority: I,H,B,E,G,T,S */
function getUpdatedResidueFlag(f: DSSPType) { function getUpdatedResidueFlag(f: DSSPType) {
if (DSSPType.is(f, DSSPType.Flag.I)) return SecondaryStructureType.SecondaryStructureDssp.I 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.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.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.G)) return SecondaryStructureType.SecondaryStructureDssp.G
if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T if (DSSPType.is(f, DSSPType.Flag.T)) return SecondaryStructureType.SecondaryStructureDssp.T
if (DSSPType.is(f, DSSPType.Flag.S)) return SecondaryStructureType.SecondaryStructureDssp.S if (DSSPType.is(f, DSSPType.Flag.S)) return SecondaryStructureType.SecondaryStructureDssp.S
return SecondaryStructureType.Flag.None return SecondaryStructureType.Flag.None
} }
// function geFlagName(f: DSSPType) { function getUpdatedFlagName(f: DSSPType) {
// if (DSSPType.is(f, DSSPType.Flag.I)) return 'I' if (DSSPType.is(f, DSSPType.Flag.I)) return 'I'
// if (DSSPType.is(f, DSSPType.Flag.H)) return 'H' 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.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.G)) return 'G'
// if (DSSPType.is(f, DSSPType.Flag.T)) return 'T' if (DSSPType.is(f, DSSPType.Flag.T)) return 'T'
// if (DSSPType.is(f, DSSPType.Flag.S)) return 'S' if (DSSPType.is(f, DSSPType.Flag.S)) return 'S'
// return '-' return '-'
// } }
function getDSSPAssignment(flags: Uint32Array, useOriginal = false) { function getDSSPAssignment(flags: Uint32Array, getResidueFlag: (f: DSSPType) => SecondaryStructureType) {
const getResidueFlag = useOriginal ? getOriginalResidueFlag : getUpdatedResidueFlag
const type = new Uint32Array(flags.length) const type = new Uint32Array(flags.length)
for (let i = 0, il = flags.length; i < il; ++i) { for (let i = 0, il = flags.length; i < il; ++i) {
const f = DSSPType.create(flags[i]) const f = DSSPType.create(flags[i])
// console.log(i, geFlagName(f))
type[i] = getResidueFlag(f) type[i] = getResidueFlag(f)
} }
return type as unknown as ArrayLike<SecondaryStructureType> return type as unknown as ArrayLike<SecondaryStructureType>
} }
/**
* Constant for electrostatic energy in kcal/mol
* f * q1 * q2
* Q = -332 * 0.42 * 0.20
*
* f is the dimensional factor
*
* q1 and q2 are partial charges which are placed on the C,O
* (+q1,-q1) and N,H (-q2,+q2)
*/
const Q = -27.888
/** cutoff for hbonds in kcal/mol, must be lower to be consider as an hbond */
const hbondEnergyCutoff = -0.5
/** /**
* E = Q * (1/r(ON) + l/r(CH) - l/r(OH) - l/r(CN)) * E = Q * (1/r(ON) + l/r(CH) - l/r(OH) - l/r(CN))
*/ */
...@@ -314,7 +545,13 @@ function calcHbondEnergy(oPos: Vec3, cPos: Vec3, nPos: Vec3, hPos: Vec3) { ...@@ -314,7 +545,13 @@ function calcHbondEnergy(oPos: Vec3, cPos: Vec3, nPos: Vec3, hPos: Vec3) {
const e1 = Q / distOH - Q / distCH const e1 = Q / distOH - Q / distCH
const e2 = Q / distCN - Q / distON const e2 = Q / distCN - Q / distON
return e1 + e2 const e = e1 + e2
// cap lowest possible energy
if (e < hbondEnergyMinimal)
return hbondEnergyMinimal
return e
} }
/** /**
...@@ -329,24 +566,31 @@ function assignTurns(ctx: DSSPContext) { ...@@ -329,24 +566,31 @@ function assignTurns(ctx: DSSPContext) {
const { chains, residueAtomSegments, chainAtomSegments } = hierarchy const { chains, residueAtomSegments, chainAtomSegments } = hierarchy
const { label_asym_id } = chains const { label_asym_id } = chains
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]
for (let i = 0, il = proteinResidues.length; i < il; ++i) { for (let idx = 0; idx < 3; idx++) {
for (let i = 0, il = proteinResidues.length - 1; i < il; ++i) {
const rI = proteinResidues[i] const rI = proteinResidues[i]
const cI = chainAtomSegments.index[residueAtomSegments.offsets[rI]] const cI = chainAtomSegments.index[residueAtomSegments.offsets[rI]]
// TODO should take sequence gaps into account // TODO should take sequence gaps into account
for (let k = 3; k <= 5; ++k) { const rN = proteinResidues[i + idx + 3]
if (i + k >= proteinResidues.length) continue
const rN = proteinResidues[i + k]
const cN = chainAtomSegments.index[residueAtomSegments.offsets[rN]] const cN = chainAtomSegments.index[residueAtomSegments.offsets[rN]]
// check if on same chain // check if on same chain
if (!label_asym_id.areValuesEqual(cI, cN)) continue if (!label_asym_id.areValuesEqual(cI, cN)) continue
// check if hbond exists // check if hbond exists
if (hbonds.getDirectedEdgeIndex(i, i + k) !== -1) { if (hbonds.getDirectedEdgeIndex(i, i + idx + 3) !== -1) {
flags[i] |= turnFlag[k] | DSSPType.Flag.T flags[i] |= turnFlag[idx + 3] | turnFlag[idx]
if (ctx.params.oldDefinition) {
for (let k = 1; k < idx + 3; ++k) {
flags[i + k] |= turnFlag[idx + 3] | DSSPType.Flag.T
}
} else {
for (let k = 0; k <= idx + 3; ++k) {
flags[i + k] |= turnFlag[idx + 3] | DSSPType.Flag.T
}
}
} }
} }
} }
...@@ -369,7 +613,7 @@ function assignTurns(ctx: DSSPContext) { ...@@ -369,7 +613,7 @@ function assignTurns(ctx: DSSPContext) {
* Type: B * Type: B
*/ */
function assignBridges(ctx: DSSPContext) { function assignBridges(ctx: DSSPContext) {
const { proteinResidues, hbonds, flags } = ctx const { proteinResidues, hbonds, flags, bridges } = ctx
const { offset, b } = hbonds const { offset, b } = hbonds
let i: number, j: number let i: number, j: number
...@@ -385,6 +629,8 @@ function assignBridges(ctx: DSSPContext) { ...@@ -385,6 +629,8 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j, i + 1) !== -1) { if (i !== j && hbonds.getDirectedEdgeIndex(j, i + 1) !== -1) {
flags[i] |= DSSPType.Flag.B flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B flags[j] |= DSSPType.Flag.B
// TODO move to constructor, actually omit object all together
bridges[bridges.length] = new Bridge(i, j, BridgeType.PARALLEL)
} }
// Parallel Bridge(i, j) =: [Hbond(j - 1, i) and Hbond(i, j + 1)] // Parallel Bridge(i, j) =: [Hbond(j - 1, i) and Hbond(i, j + 1)]
...@@ -393,6 +639,7 @@ function assignBridges(ctx: DSSPContext) { ...@@ -393,6 +639,7 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i) !== -1) { if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i) !== -1) {
flags[i] |= DSSPType.Flag.B flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B flags[j] |= DSSPType.Flag.B
bridges[bridges.length] = new Bridge(j, i, BridgeType.PARALLEL)
} }
// Antiparallel Bridge(i, j) =: [Hbond(i, j) and Hbond(j, i)] // Antiparallel Bridge(i, j) =: [Hbond(i, j) and Hbond(j, i)]
...@@ -401,6 +648,7 @@ function assignBridges(ctx: DSSPContext) { ...@@ -401,6 +648,7 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j, i) !== -1) { if (i !== j && hbonds.getDirectedEdgeIndex(j, i) !== -1) {
flags[i] |= DSSPType.Flag.B flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B flags[j] |= DSSPType.Flag.B
bridges[bridges.length] = new Bridge(j, i, BridgeType.ANTI_PARALLEL)
} }
// Antiparallel Bridge(i, j) =: [Hbond(i - 1, j + 1) and Hbond(j - 1, i + l)] // Antiparallel Bridge(i, j) =: [Hbond(i - 1, j + 1) and Hbond(j - 1, i + l)]
...@@ -409,9 +657,12 @@ function assignBridges(ctx: DSSPContext) { ...@@ -409,9 +657,12 @@ function assignBridges(ctx: DSSPContext) {
if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i + 1) !== -1) { if (i !== j && hbonds.getDirectedEdgeIndex(j - 1, i + 1) !== -1) {
flags[i] |= DSSPType.Flag.B flags[i] |= DSSPType.Flag.B
flags[j] |= DSSPType.Flag.B flags[j] |= DSSPType.Flag.B
bridges[bridges.length] = new Bridge(j, i, BridgeType.ANTI_PARALLEL)
} }
} }
} }
bridges.sort((a, b) => a.partner1 > b.partner1 ? 1 : a.partner1 < b.partner1 ? -1 : 0)
} }
/** /**
...@@ -428,17 +679,41 @@ function assignBridges(ctx: DSSPContext) { ...@@ -428,17 +679,41 @@ function assignBridges(ctx: DSSPContext) {
function assignHelices(ctx: DSSPContext) { function assignHelices(ctx: DSSPContext) {
const { proteinResidues, flags } = ctx 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] 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 helixCheckOrder = ctx.params.oldOrdering ? [4, 3, 5] : [3, 4, 5]
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 fI = DSSPType.create(flags[i])
const fI1 = DSSPType.create(flags[i - 1]) 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.params.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
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
continue
}
}
for (let k = 3; k <= 5; ++k) { if (DSSPType.is(fI, turnFlag[n]) && DSSPType.is(fI, turnFlag[n - 3]) && // check fI for turn start of proper type
if (DSSPType.is(fI, turnFlag[k]) && DSSPType.is(fI1, turnFlag[k])) { DSSPType.is(fI1, turnFlag[n]) && DSSPType.is(fI1, turnFlag[n - 3])) { // check fI1 accordingly
for (let l = 0; l < k; ++l) { if (ctx.params.oldDefinition) {
flags[i + l] |= helixFlag[k] for (let k = 0; k < n; k++) {
flags[i + k] |= helixFlag[n]
}
} else {
for (let k = -1; k <= n; k++) {
flags[i + k] |= helixFlag[n]
}
} }
} }
} }
...@@ -451,7 +726,95 @@ function assignHelices(ctx: DSSPContext) { ...@@ -451,7 +726,95 @@ function assignHelices(ctx: DSSPContext) {
* Type: E * Type: E
*/ */
function assignLadders(ctx: DSSPContext) { function assignLadders(ctx: DSSPContext) {
// TODO const { bridges, ladders } = ctx
// 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)) {
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
}
}
}
// 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)) {
ladder1.nextLadder = ladderIndex2
ladder2.previousLadder = ladderIndex1
}
}
}
}
/**
* 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 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
// 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) {
return true
}
} else {
if (bridge.partner2 === ladder.secondStart - 1) {
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)
} }
/** /**
...@@ -460,14 +823,40 @@ function assignLadders(ctx: DSSPContext) { ...@@ -460,14 +823,40 @@ function assignLadders(ctx: DSSPContext) {
* Type: E * Type: E
*/ */
function assignSheets(ctx: DSSPContext) { function assignSheets(ctx: DSSPContext) {
// TODO 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
* Bend(i) =: [angle ((CW - Ca(i - 2)),(C"(i + 2) - C"(i))) > 70"]
* const conladder = ladders[ladder.nextLadder]
* Type: S for (let lcount = ladder.firstStart; lcount <= conladder.firstEnd; lcount++) {
*/ flags[lcount] |= DSSPType.Flag.E
function assignBends(ctx: DSSPContext) { }
// TODO 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
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
TODO move to spec.ts once tests are running
\ No newline at end of file
...@@ -275,41 +275,42 @@ export namespace SecondaryStructureType { ...@@ -275,41 +275,42 @@ export namespace SecondaryStructureType {
DoubleHelix = 0x1, DoubleHelix = 0x1,
Helix = 0x2, Helix = 0x2,
Beta = 0x4, Beta = 0x4,
Turn = 0x8, Bend = 0x8,
Turn = 0x10,
// category variant // category variant
LeftHanded = 0x10, // helix LeftHanded = 0x20, // helix
RightHanded = 0x20, RightHanded = 0x40,
ClassicTurn = 0x40, // turn ClassicTurn = 0x80, // turn
InverseTurn = 0x80, InverseTurn = 0x100,
// sub-category // sub-category
HelixOther = 0x100, // protein HelixOther = 0x200, // protein
Helix27 = 0x200, Helix27 = 0x400,
Helix3Ten = 0x400, Helix3Ten = 0x800,
HelixAlpha = 0x800, HelixAlpha = 0x1000,
HelixGamma = 0x1000, HelixGamma = 0x2000,
HelixOmega = 0x2000, HelixOmega = 0x4000,
HelixPi = 0x4000, HelixPi = 0x8000,
HelixPolyproline = 0x8000, HelixPolyproline = 0x10000,
DoubleHelixOther = 0x10000, // nucleic DoubleHelixOther = 0x20000, // nucleic
DoubleHelixZ = 0x20000, DoubleHelixZ = 0x40000,
DoubleHelixA = 0x40000, DoubleHelixA = 0x80000,
DoubleHelixB = 0x80000, DoubleHelixB = 0x100000,
BetaOther = 0x100000, // protein BetaOther = 0x200000, // protein
BetaStrand = 0x200000, // single strand BetaStrand = 0x400000, // single strand
BetaSheet = 0x400000, // multiple hydrogen bonded strands BetaSheet = 0x800000, // multiple hydrogen bonded strands
BetaBarell = 0x800000, // closed series of sheets BetaBarell = 0x1000000, // closed series of sheets
TurnOther = 0x1000000, // protein TurnOther = 0x2000000, // protein
Turn1 = 0x2000000, Turn1 = 0x4000000,
Turn2 = 0x4000000, Turn2 = 0x8000000,
Turn3 = 0x8000000, Turn3 = 0x10000000,
NA = 0x10000000, // not applicable/available NA = 0x20000000, // not applicable/available
} }
export const SecondaryStructureMmcif: { [value: string]: number } = { export const SecondaryStructureMmcif: { [value: string]: number } = {
...@@ -386,7 +387,7 @@ export namespace SecondaryStructureType { ...@@ -386,7 +387,7 @@ export namespace SecondaryStructureType {
G: Flag.Helix | Flag.Helix3Ten, // 3-helix (310 helix) G: Flag.Helix | Flag.Helix3Ten, // 3-helix (310 helix)
I: Flag.Helix | Flag.HelixPi, // 5 helix (pi-helix) I: Flag.Helix | Flag.HelixPi, // 5 helix (pi-helix)
T: Flag.Turn, // hydrogen bonded turn T: Flag.Turn, // hydrogen bonded turn
S: Flag.Turn, // bend S: Flag.Bend, // bend
} }
} }
......
...@@ -22,6 +22,8 @@ const SecondaryStructureColors = ColorMap({ ...@@ -22,6 +22,8 @@ const SecondaryStructureColors = ColorMap({
'betaTurn': 0x6080FF, 'betaTurn': 0x6080FF,
'betaStrand': 0xFFC800, 'betaStrand': 0xFFC800,
'coil': 0xFFFFFF, 'coil': 0xFFFFFF,
'bend': 0x66D8C9 /* biting original color used 0x00FF00 */,
'turn': 0x00B266,
'dna': 0xAE00FE, 'dna': 0xAE00FE,
'rna': 0xFD0162, 'rna': 0xFD0162,
...@@ -53,8 +55,10 @@ export function secondaryStructureColor(unit: Unit, element: ElementIndex): Colo ...@@ -53,8 +55,10 @@ export function secondaryStructureColor(unit: Unit, element: ElementIndex): Colo
return SecondaryStructureColors.alphaHelix return SecondaryStructureColors.alphaHelix
} else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Beta)) { } else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Beta)) {
return SecondaryStructureColors.betaStrand return SecondaryStructureColors.betaStrand
} else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Bend)) {
return SecondaryStructureColors.bend
} else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Turn)) { } else if (SecondaryStructureType.is(secStrucType, SecondaryStructureType.Flag.Turn)) {
return SecondaryStructureColors.coil return SecondaryStructureColors.turn
} else { } else {
const moleculeType = getElementMoleculeType(unit, element) const moleculeType = getElementMoleculeType(unit, element)
if (moleculeType === MoleculeType.DNA) { if (moleculeType === MoleculeType.DNA) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment