diff --git a/src/mol-math/linear-algebra/3d/vec3.ts b/src/mol-math/linear-algebra/3d/vec3.ts index 3e8fd5d5ee7765fcb88215a1455655b75357b10e..bb16cc3530375060824f061d2a8176b4112b2325 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); @@ -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 ===) */ 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 0000000000000000000000000000000000000000..bb98228acaf7b051186c8dd54cc985d81159ce52 --- /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 49618f0071fc1974657ffd2f01cd1ffaf5c18d31..9434507ebc9ae1a11b2fc02f4cb9e00d0d4c4d0c 100644 --- a/src/mol-model/structure/model/properties/seconday-structure.ts +++ b/src/mol-model/structure/model/properties/seconday-structure.ts @@ -17,12 +17,17 @@ interface SecondaryStructure { } namespace SecondaryStructure { - export type Element = None | Helix | Sheet + export type Element = None | Turn | Helix | Sheet export interface None { kind: 'none' } + export interface Turn { + kind: 'turn', + flags: SecondaryStructureType + } + export interface Helix { kind: 'helix', flags: SecondaryStructureType, 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 afa8199475f42a616da01e9c8ec056a416c0fd53..68d8c96acbe8615d4f9f4987269799b281acecae 100644 --- a/src/mol-model/structure/model/properties/utils/secondary-structure.ts +++ b/src/mol-model/structure/model/properties/utils/secondary-structure.ts @@ -2,6 +2,7 @@ * 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 { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure'; @@ -14,13 +15,120 @@ import { IntAdjacencyGraph } from 'mol-math/graph'; import { BitFlags } from 'mol-util'; import { ElementIndex } from 'mol-model/structure/model/indexing'; 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 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 backboneIndices = calcBackboneAtomIndices(hierarchy, proteinResidues) const hbonds = calcBackboneHbonds(hierarchy, conformation, proteinResidues, backboneIndices, lookup3d) @@ -28,69 +136,103 @@ export function computeModelDSSP(hierarchy: AtomicHierarchy, conformation: Atomi const residueCount = proteinResidues.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 ladders: Ladder[] = [] + const bridges: Bridge[] = [] + + const getResidueFlag = params.oldDefinition ? getOriginalResidueFlag : getUpdatedResidueFlag + const getFlagName = params.oldOrdering ? getOriginalFlagName : getUpdatedFlagName + const ctx: DSSPContext = { + params, + 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]) + // 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 = { type, - key: [], // TODO - elements: [] // TODO + key: keys, + elements: elements } + return secondaryStructure } -interface DSSPContext { - hierarchy: AtomicHierarchy - proteinResidues: SortedArray<ResidueIndex> - /** flags for each residue */ - flags: Uint32Array - - hbonds: DsspHbonds +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 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, +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 if (assignment === SecondaryStructureType.SecondaryStructureDssp.T) { + return 'turn' + } else if (assignment === SecondaryStructureType.SecondaryStructureDssp.S) { + return 'bend' + } else { + return 'none' } } -/** max distance between two C-alpha atoms to check for hbond */ -const caMaxDist = 7.0; - -/** min distance between two C-alpha atoms to check for hbond */ -const caMinDist = 4.0; - function calcAtomicTraceLookup3D(hierarchy: AtomicHierarchy, conformation: AtomicConformation) { const { x, y, z } = conformation; const { moleculeType, traceElementIndex } = hierarchy.derived.residue @@ -141,6 +283,104 @@ 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() + + 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 { const { cIndices, hIndices, nIndices, oIndices } = backboneIndices const { index } = hierarchy @@ -163,8 +403,6 @@ function calcBackboneHbonds(hierarchy: AtomicHierarchy, conformation: AtomicConf const cPosPrev = Vec3.zero() const oPosPrev = Vec3.zero() - const caMinDistSq = caMinDist * caMinDist - for (let i = 0, il = proteinResidues.length; i < il; ++i) { const oPI = i const oRI = proteinResidues[i] @@ -183,11 +421,9 @@ 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 - const nPI = indices[j] // 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 /** 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,55 +490,50 @@ 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 getDSSPAssignment(flags: Uint32Array, useOriginal = false) { - const getResidueFlag = useOriginal ? getOriginalResidueFlag : getUpdatedResidueFlag +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, 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> } -/** - * 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)) */ @@ -314,7 +545,13 @@ function calcHbondEnergy(oPos: Vec3, cPos: Vec3, nPos: Vec3, hPos: Vec3) { 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 } /** @@ -329,24 +566,31 @@ function assignTurns(ctx: DSSPContext) { const { chains, residueAtomSegments, chainAtomSegments } = hierarchy 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) { - const rI = proteinResidues[i] - const cI = chainAtomSegments.index[residueAtomSegments.offsets[rI]] + 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 - for (let k = 3; k <= 5; ++k) { - if (i + k >= proteinResidues.length) continue - - 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 + if (hbonds.getDirectedEdgeIndex(i, i + idx + 3) !== -1) { + 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) { * 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 +629,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] = new Bridge(i, j, BridgeType.PARALLEL) } // Parallel Bridge(i, j) =: [Hbond(j - 1, i) and Hbond(i, j + 1)] @@ -393,6 +639,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] = new Bridge(j, i, BridgeType.PARALLEL) } // Antiparallel Bridge(i, j) =: [Hbond(i, j) and Hbond(j, i)] @@ -401,6 +648,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] = new Bridge(j, i, BridgeType.ANTI_PARALLEL) } // Antiparallel Bridge(i, j) =: [Hbond(i - 1, j + 1) and Hbond(j - 1, i + l)] @@ -409,9 +657,12 @@ 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] = 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) { 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.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 fI1 = DSSPType.create(flags[i - 1]) + const fI2 = DSSPType.create(flags[i + 1]) - 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] + // 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 + } + } + + 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 + if (ctx.params.oldDefinition) { + 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,23 +726,137 @@ function assignHelices(ctx: DSSPContext) { * Type: E */ 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 + } + } + } } /** - * 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 + + // 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) } /** - * 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 0000000000000000000000000000000000000000..4f1a1f1bba2be14065663ee39ef36e21fb2d575a --- /dev/null +++ b/src/mol-model/structure/model/properties/utils/secondary-structure.validation @@ -0,0 +1,75 @@ +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 diff --git a/src/mol-model/structure/model/types.ts b/src/mol-model/structure/model/types.ts index 89293da9605b464d3eec6efe7ec3f9d136258b0e..76c0982a1f65792f5597106ad90b011c01420d79 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 384c5fb0b7d3c37031ef670b0d17a9e7b9f79ade..cb8878f28b77889d830573b5609a3b1c32feed6e 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) {