Skip to content
Snippets Groups Projects
Commit 606ef255 authored by Alexander Rose's avatar Alexander Rose
Browse files

per-unit dssp as structure property

parent 992d0f10
No related branches found
No related tags found
No related merge requests found
...@@ -19,14 +19,17 @@ export function getSecondaryStructure(data: mmCIF_Database, hierarchy: AtomicHie ...@@ -19,14 +19,17 @@ export function getSecondaryStructure(data: mmCIF_Database, hierarchy: AtomicHie
// must add Helices 1st because of 'key' value assignment. // must add Helices 1st because of 'key' value assignment.
addSheets(data.struct_sheet_range, map, data.struct_conf._rowCount, elements); addSheets(data.struct_sheet_range, map, data.struct_conf._rowCount, elements);
const n = hierarchy.residues._rowCount
const getIndex = (rI: ResidueIndex) => rI
const secStruct: SecondaryStructureData = { const secStruct: SecondaryStructureData = {
type: new Int32Array(hierarchy.residues._rowCount) as any, type: new Int32Array(n) as any,
key: new Int32Array(hierarchy.residues._rowCount) as any, key: new Int32Array(n) as any,
elements elements
}; };
if (map.size > 0) assignSecondaryStructureRanges(hierarchy, map, secStruct); if (map.size > 0) assignSecondaryStructureRanges(hierarchy, map, secStruct);
return SecondaryStructure(secStruct.type, secStruct.key, secStruct.elements); return SecondaryStructure(secStruct.type, secStruct.key, secStruct.elements, getIndex);
} }
type SecondaryStructureEntry = { type SecondaryStructureEntry = {
......
...@@ -59,7 +59,8 @@ export type SecondaryStructureComputationProps = PD.Values<SecondaryStructureCom ...@@ -59,7 +59,8 @@ export type SecondaryStructureComputationProps = PD.Values<SecondaryStructureCom
async function computeSecondaryStructure(structure: Structure, params: Partial<SecondaryStructureComputationProps>): Promise<ComputedSecondaryStructure.Property> { async function computeSecondaryStructure(structure: Structure, params: Partial<SecondaryStructureComputationProps>): Promise<ComputedSecondaryStructure.Property> {
const p = { ...PD.getDefaultValues(SecondaryStructureComputationParams), params } const p = { ...PD.getDefaultValues(SecondaryStructureComputationParams), params }
// TODO take inter-unit hbonds into account // TODO take inter-unit hbonds into account for bridge, ladder, sheet assignment
// TODO store unit-only secStruc as custom unit property???
// 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
const map = new Map<number, SecondaryStructure>() const map = new Map<number, SecondaryStructure>()
for (let i = 0, il = structure.unitSymmetryGroups.length; i < il; ++i) { for (let i = 0, il = structure.unitSymmetryGroups.length; i < il; ++i) {
......
...@@ -20,6 +20,8 @@ import { calculateUnitDihedralAngles } from './dssp/dihedral-angles'; ...@@ -20,6 +20,8 @@ import { calculateUnitDihedralAngles } from './dssp/dihedral-angles';
import { calcUnitProteinTraceLookup3D } from './dssp/trace-lookup'; import { calcUnitProteinTraceLookup3D } from './dssp/trace-lookup';
import { Unit } from 'mol-model/structure'; import { Unit } from 'mol-model/structure';
import { getUnitProteinInfo } from './dssp/protein-info'; import { getUnitProteinInfo } from './dssp/protein-info';
import { ResidueIndex } from 'mol-model/structure/model';
import { SortedArray } from 'mol-data/int';
/** /**
* TODO bugs to fix: * TODO bugs to fix:
...@@ -37,7 +39,6 @@ export type DSSPComputationParams = typeof DSSPComputationParams ...@@ -37,7 +39,6 @@ export type DSSPComputationParams = typeof DSSPComputationParams
export type DSSPComputationProps = PD.Values<DSSPComputationParams> export type DSSPComputationProps = PD.Values<DSSPComputationParams>
export async function computeUnitDSSP(unit: Unit.Atomic, params: DSSPComputationProps): Promise<SecondaryStructure> { export async function computeUnitDSSP(unit: Unit.Atomic, params: DSSPComputationProps): Promise<SecondaryStructure> {
const proteinInfo = getUnitProteinInfo(unit) const proteinInfo = getUnitProteinInfo(unit)
const { residueIndices } = proteinInfo const { residueIndices } = proteinInfo
const lookup3d = calcUnitProteinTraceLookup3D(unit, residueIndices) const lookup3d = calcUnitProteinTraceLookup3D(unit, residueIndices)
...@@ -83,20 +84,20 @@ export async function computeUnitDSSP(unit: Unit.Atomic, params: DSSPComputation ...@@ -83,20 +84,20 @@ export async function computeUnitDSSP(unit: Unit.Atomic, params: DSSPComputation
const type = new Uint32Array(residueCount) as unknown as SecondaryStructureType[] const type = new Uint32Array(residueCount) as unknown as SecondaryStructureType[]
const keys: number[] = [] const keys: number[] = []
const elements: SecondaryStructure.Element[] = [] const elements: SecondaryStructure.Element[] = []
const getIndex = (rI: ResidueIndex) => SortedArray.indexOf(residueIndices, rI)
for (let i = 0, il = residueIndices.length; i < il; ++i) { for (let i = 0, il = residueIndices.length; i < il; ++i) {
const assign = assignment[i] const assign = assignment[i]
type[residueIndices[i]] = assign type[i] = assign
const flag = getResidueFlag(flags[i]) const flag = getResidueFlag(flags[i])
// console.log(i, SortedArray.indexOf(residueIndices, i), getFlagName(flags[i]))
// TODO is this expected behavior? elements will be strictly split depending on 'winning' flag // 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 */) { 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) elements[elements.length] = createElement(mapToKind(assign), flags[i], getResidueFlag)
} }
keys[i] = elements.length - 1 keys[i] = elements.length - 1
} }
// TODO expose model to unit residueIndex mapping return SecondaryStructure(type, keys, elements, getIndex)
return SecondaryStructure(type, keys, elements)
} }
function createElement(kind: string, flag: DSSPType.Flag, getResidueFlag: (f: DSSPType) => SecondaryStructureType): SecondaryStructure.Element { function createElement(kind: string, flag: DSSPType.Flag, getResidueFlag: (f: DSSPType) => SecondaryStructureType): SecondaryStructure.Element {
......
...@@ -14,7 +14,7 @@ export function calcUnitProteinTraceLookup3D(unit: Unit.Atomic, unitProteinResid ...@@ -14,7 +14,7 @@ export function calcUnitProteinTraceLookup3D(unit: Unit.Atomic, unitProteinResid
const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue const { traceElementIndex } = unit.model.atomicHierarchy.derived.residue
const indices = new Uint32Array(unitProteinResidues.length) const indices = new Uint32Array(unitProteinResidues.length)
for (let i = 0, il = unitProteinResidues.length; i < il; ++i) { for (let i = 0, il = unitProteinResidues.length; i < il; ++i) {
indices[indices.length] = traceElementIndex[unitProteinResidues[i]] indices[i] = traceElementIndex[unitProteinResidues[i]]
} }
return GridLookup3D({ x, y, z, indices: SortedArray.ofSortedArray(indices) }); return GridLookup3D({ x, y, z, indices: SortedArray.ofSortedArray(indices) });
} }
\ No newline at end of file
...@@ -5,6 +5,7 @@ ...@@ -5,6 +5,7 @@
*/ */
import { SecondaryStructureType } from '../types'; import { SecondaryStructureType } from '../types';
import { ResidueIndex } from '../indexing';
/** Secondary structure "indexed" by residues. */ /** Secondary structure "indexed" by residues. */
interface SecondaryStructure { interface SecondaryStructure {
...@@ -13,10 +14,12 @@ interface SecondaryStructure { ...@@ -13,10 +14,12 @@ interface SecondaryStructure {
readonly key: ArrayLike<number>, readonly key: ArrayLike<number>,
/** indexed by key */ /** indexed by key */
readonly elements: ReadonlyArray<SecondaryStructure.Element> readonly elements: ReadonlyArray<SecondaryStructure.Element>
/** mapping from residue index */
readonly getIndex: (rI: ResidueIndex) => number,
} }
function SecondaryStructure(type: SecondaryStructure['type'], key: SecondaryStructure['key'], elements: SecondaryStructure['elements']) { function SecondaryStructure(type: SecondaryStructure['type'], key: SecondaryStructure['key'], elements: SecondaryStructure['elements'], getIndex: SecondaryStructure['getIndex']) {
return { type, key, elements } return { type, key, elements, getIndex }
} }
namespace SecondaryStructure { namespace SecondaryStructure {
......
...@@ -14,6 +14,7 @@ import { CoarseSphereConformation, CoarseGaussianConformation } from 'mol-model/ ...@@ -14,6 +14,7 @@ import { CoarseSphereConformation, CoarseGaussianConformation } from 'mol-model/
import { getPolymerRanges } from '../polymer'; import { getPolymerRanges } from '../polymer';
import { AtomicConformation } from 'mol-model/structure/model/properties/atomic'; import { AtomicConformation } from 'mol-model/structure/model/properties/atomic';
import { ComputedSecondaryStructure } from 'mol-model-props/computed/secondary-structure'; import { ComputedSecondaryStructure } from 'mol-model-props/computed/secondary-structure';
import { SecondaryStructure } from 'mol-model/structure/model/properties/seconday-structure';
/** /**
* Iterates over individual residues/coarse elements in polymers of a unit while * Iterates over individual residues/coarse elements in polymers of a unit while
...@@ -69,7 +70,8 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement> ...@@ -69,7 +70,8 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
private residueIt: Segmentation.SegmentIterator<ResidueIndex> private residueIt: Segmentation.SegmentIterator<ResidueIndex>
private polymerSegment: Segmentation.Segment<ResidueIndex> private polymerSegment: Segmentation.Segment<ResidueIndex>
private cyclicPolymerMap: Map<ResidueIndex, ResidueIndex> private cyclicPolymerMap: Map<ResidueIndex, ResidueIndex>
private secondaryStructureType: ArrayLike<SecondaryStructureType> private secondaryStructureType: SecondaryStructure['type']
private secondaryStructureGetIndex: SecondaryStructure['getIndex']
private residueSegmentMin: ResidueIndex private residueSegmentMin: ResidueIndex
private residueSegmentMax: ResidueIndex private residueSegmentMax: ResidueIndex
private prevSecStrucType: SecondaryStructureType private prevSecStrucType: SecondaryStructureType
...@@ -133,8 +135,12 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement> ...@@ -133,8 +135,12 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
return residueIndex as ResidueIndex return residueIndex as ResidueIndex
} }
private getSecStruc(residueIndex: number) {
return this.secondaryStructureType[this.secondaryStructureGetIndex(residueIndex as ResidueIndex)]
}
private setControlPoint(out: Vec3, p1: Vec3, p2: Vec3, p3: Vec3, residueIndex: ResidueIndex) { private setControlPoint(out: Vec3, p1: Vec3, p2: Vec3, p3: Vec3, residueIndex: ResidueIndex) {
const ss = this.secondaryStructureType[residueIndex] const ss = this.getSecStruc(residueIndex)
if (SecondaryStructureType.is(ss, SecondaryStructureType.Flag.Beta)) { if (SecondaryStructureType.is(ss, SecondaryStructureType.Flag.Beta)) {
Vec3.scale(out, Vec3.add(out, p1, Vec3.add(out, p3, Vec3.add(out, p2, p2))), 1/4) Vec3.scale(out, Vec3.add(out, p1, Vec3.add(out, p3, Vec3.add(out, p2, p2))), 1/4)
} else { } else {
...@@ -153,7 +159,7 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement> ...@@ -153,7 +159,7 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
if (residueIt.hasNext) { if (residueIt.hasNext) {
this.state = AtomicPolymerTraceIteratorState.nextResidue this.state = AtomicPolymerTraceIteratorState.nextResidue
this.currSecStrucType = SecStrucTypeNA this.currSecStrucType = SecStrucTypeNA
this.nextSecStrucType = this.secondaryStructureType[this.residueSegmentMin] this.nextSecStrucType = this.getSecStruc(this.residueSegmentMin)
this.currCoarseBackbone = false this.currCoarseBackbone = false
this.nextCoarseBackbone = this.directionElementIndex[this.residueSegmentMin] === -1 this.nextCoarseBackbone = this.directionElementIndex[this.residueSegmentMin] === -1
break break
...@@ -165,7 +171,7 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement> ...@@ -165,7 +171,7 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
const { index: residueIndex } = residueIt.move(); const { index: residueIndex } = residueIt.move();
this.prevSecStrucType = this.currSecStrucType this.prevSecStrucType = this.currSecStrucType
this.currSecStrucType = this.nextSecStrucType this.currSecStrucType = this.nextSecStrucType
this.nextSecStrucType = residueIt.hasNext ? this.secondaryStructureType[residueIndex + 1] : SecStrucTypeNA this.nextSecStrucType = residueIt.hasNext ? this.getSecStruc(residueIndex + 1) : SecStrucTypeNA
this.prevCoarseBackbone = this.currCoarseBackbone this.prevCoarseBackbone = this.currCoarseBackbone
this.currCoarseBackbone = this.nextCoarseBackbone this.currCoarseBackbone = this.nextCoarseBackbone
...@@ -250,10 +256,14 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement> ...@@ -250,10 +256,14 @@ export class AtomicPolymerTraceIterator implements Iterator<PolymerTraceElement>
this.hasNext = this.residueIt.hasNext && this.polymerIt.hasNext this.hasNext = this.residueIt.hasNext && this.polymerIt.hasNext
this.secondaryStructureType = unit.model.properties.secondaryStructure.type this.secondaryStructureType = unit.model.properties.secondaryStructure.type
this.secondaryStructureGetIndex = unit.model.properties.secondaryStructure.getIndex
const computedSecondaryStructure = ComputedSecondaryStructure.get(structure) const computedSecondaryStructure = ComputedSecondaryStructure.get(structure)
if (computedSecondaryStructure) { if (computedSecondaryStructure) {
const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId) const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId)
if (secondaryStructure) this.secondaryStructureType = secondaryStructure.type if (secondaryStructure) {
this.secondaryStructureType = secondaryStructure.type
this.secondaryStructureGetIndex = secondaryStructure.getIndex
}
} }
} }
} }
......
...@@ -47,7 +47,7 @@ export function secondaryStructureColor(unit: Unit, element: ElementIndex, compu ...@@ -47,7 +47,7 @@ export function secondaryStructureColor(unit: Unit, element: ElementIndex, compu
secStrucType = unit.model.properties.secondaryStructure.type[unit.residueIndex[element]] secStrucType = unit.model.properties.secondaryStructure.type[unit.residueIndex[element]]
if (computedSecondaryStructure) { if (computedSecondaryStructure) {
const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId) const secondaryStructure = computedSecondaryStructure.map.get(unit.invariantId)
if (secondaryStructure) secStrucType = secondaryStructure.type[unit.residueIndex[element]] if (secondaryStructure) secStrucType = secondaryStructure.type[secondaryStructure.getIndex(unit.residueIndex[element])]
} }
} }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment