diff --git a/src/apps/structure-info/model.ts b/src/apps/structure-info/model.ts index 0728dde5d6a95f069cc3fabc9d38037694f421ab..9c7f0d9a9876f39c78cd5887f9696ef377c1ea4b 100644 --- a/src/apps/structure-info/model.ts +++ b/src/apps/structure-info/model.ts @@ -18,6 +18,7 @@ import { mmCIF_Database } from 'mol-io/reader/cif/schema/mmcif'; import { openCif, downloadCif } from './helpers'; import { BitFlags } from 'mol-util'; import { SecondaryStructureType } from 'mol-model/structure/model/types'; +import { UnitRings } from 'mol-model/structure/structure/unit/rings'; async function downloadFromPdb(pdb: string) { @@ -100,6 +101,21 @@ export function printSequence(model: Model) { console.log(); } +export function printRings(structure: Structure) { + console.log('\nRings\n============='); + for (const unit of structure.units) { + if (!Unit.isAtomic(unit)) continue; + const { all, byFingerprint } = unit.rings; + const fps: string[] = []; + for (let i = 0, _i = Math.min(5, all.length); i < _i; i++) { + fps[fps.length] = UnitRings.getRingFingerprint(unit, all[i]); + } + if (all.length > 5) fps.push('...') + console.log(`Unit ${unit.id}, ${all.length} ring(s), ${byFingerprint.size} different fingerprint(s).\n ${fps.join(', ')}`); + } + console.log(); +} + export function printUnits(structure: Structure) { console.log('\nUnits\n============='); const l = Element.Location(); @@ -143,6 +159,7 @@ async function run(mmcif: mmCIF_Database) { printSequence(models[0]); //printIHMModels(models[0]); printUnits(structure); + printRings(structure); //printBonds(structure); //printSecStructure(models[0]); } diff --git a/src/mol-io/writer/cif/encoder.ts b/src/mol-io/writer/cif/encoder.ts index 3e51663ac42e82effffa731e7b38529c42623d3a..05dcf0052bb191bb92370f30a54493c38ff76a79 100644 --- a/src/mol-io/writer/cif/encoder.ts +++ b/src/mol-io/writer/cif/encoder.ts @@ -21,9 +21,10 @@ import { ArrayEncoder, ArrayEncoding } from '../../common/binary-cif'; export interface Field<Key = any, Data = any> { name: string, type: Field.Type, + value(key: Key, data: Data): string | number valueKind?: (key: Key, data: Data) => Column.ValueKind, defaultFormat?: Field.Format, - value(key: Key, data: Data): string | number + shouldInclude?: (data: Data) => boolean } export namespace Field { @@ -35,31 +36,31 @@ export namespace Field { typedArray?: ArrayEncoding.TypedArrayCtor } - export function getDigitCount(field: Field) { - if (field.defaultFormat && typeof field.defaultFormat.digitCount !== 'undefined') return Math.max(0, Math.min(field.defaultFormat.digitCount, 16)); - return 6; - } + export type ParamsBase<K, D> = { valueKind?: (k: K, d: D) => Column.ValueKind, encoder?: ArrayEncoder, shouldInclude?: (data: D) => boolean } - export function str<K, D = any>(name: string, value: (k: K, d: D) => string, params?: { valueKind?: (k: K, d: D) => Column.ValueKind, encoder?: ArrayEncoder }): Field<K, D> { - return { name, type: Type.Str, value, valueKind: params && params.valueKind, defaultFormat: params && params.encoder ? { encoder: params.encoder } : void 0 }; + export function str<K, D = any>(name: string, value: (k: K, d: D) => string, params?: ParamsBase<K, D>): Field<K, D> { + return { name, type: Type.Str, value, valueKind: params && params.valueKind, defaultFormat: params && params.encoder ? { encoder: params.encoder } : void 0, shouldInclude: params && params.shouldInclude }; } - export function int<K, D = any>(name: string, value: (k: K, d: D) => number, params?: { valueKind?: (k: K, d: D) => Column.ValueKind, encoder?: ArrayEncoder, typedArray?: ArrayEncoding.TypedArrayCtor }): Field<K, D> { + export function int<K, D = any>(name: string, value: (k: K, d: D) => number, params?: ParamsBase<K, D> & { typedArray?: ArrayEncoding.TypedArrayCtor }): Field<K, D> { return { name, type: Type.Int, value, valueKind: params && params.valueKind, - defaultFormat: params ? { encoder: params.encoder, typedArray: params.typedArray } : void 0 }; + defaultFormat: params ? { encoder: params.encoder, typedArray: params.typedArray } : void 0, + shouldInclude: params && params.shouldInclude + }; } - export function float<K, D = any>(name: string, value: (k: K, d: D) => number, params?: { valueKind?: (k: K, d: D) => Column.ValueKind, encoder?: ArrayEncoder, typedArray?: ArrayEncoding.TypedArrayCtor, digitCount?: number }): Field<K, D> { + export function float<K, D = any>(name: string, value: (k: K, d: D) => number, params?: ParamsBase<K, D> & { typedArray?: ArrayEncoding.TypedArrayCtor, digitCount?: number }): Field<K, D> { return { name, type: Type.Float, value, valueKind: params && params.valueKind, - defaultFormat: params ? { encoder: params.encoder, typedArray: params.typedArray, digitCount: typeof params.digitCount !== 'undefined' ? params.digitCount : void 0 } : void 0 + defaultFormat: params ? { encoder: params.encoder, typedArray: params.typedArray, digitCount: typeof params.digitCount !== 'undefined' ? params.digitCount : void 0 } : void 0, + shouldInclude: params && params.shouldInclude }; } } diff --git a/src/mol-io/writer/cif/encoder/binary.ts b/src/mol-io/writer/cif/encoder/binary.ts index e10b86940a15efd63c6b0d00d1c046a185300fe0..564adbfd8059a37d877d6617be08e82ae48991e6 100644 --- a/src/mol-io/writer/cif/encoder/binary.ts +++ b/src/mol-io/writer/cif/encoder/binary.ts @@ -14,6 +14,7 @@ import { } from '../../../common/binary-cif' import { Field, Category, Encoder } from '../encoder' import Writer from '../../writer' +import { getIncludedFields } from './util'; export default class BinaryEncoder implements Encoder<Uint8Array> { private data: EncodedFile; @@ -57,12 +58,17 @@ export default class BinaryEncoder implements Encoder<Uint8Array> { const first = categories[0]!; const cat: EncodedCategory = { name: '_' + first.name, columns: [], rowCount: count }; const data = categories.map(c => ({ data: c.data, keys: () => c.keys ? c.keys() : Iterator.Range(0, c.rowCount - 1) })); - for (const f of first.fields) { + const fields = getIncludedFields(first); + + for (const f of fields) { if (!this.filter.includeField(first.name, f.name)) continue; const format = this.formatter.getFormat(first.name, f.name); cat.columns.push(encodeField(f, data, count, getArrayCtor(f, format), getEncoder(f, format))); } + // no columns included. + if (!cat.columns.length) return; + this.dataBlocks[this.dataBlocks.length - 1].categories.push(cat); } diff --git a/src/mol-io/writer/cif/encoder/text.ts b/src/mol-io/writer/cif/encoder/text.ts index 7352dbace857d189afe873ffd3677c8bf5c35e22..55ec8ec85db359d964f72526ba0813de77fe9be5 100644 --- a/src/mol-io/writer/cif/encoder/text.ts +++ b/src/mol-io/writer/cif/encoder/text.ts @@ -11,6 +11,7 @@ import { Column } from 'mol-data/db' import StringBuilder from 'mol-util/string-builder' import { Category, Field, Encoder } from '../encoder' import Writer from '../../writer' +import { getFieldDigitCount, getIncludedFields } from './util'; export default class TextEncoder implements Encoder<string> { private builder = StringBuilder.create(); @@ -102,13 +103,13 @@ function getFloatPrecisions(categoryName: string, fields: Field[], formatter: Ca for (const f of fields) { const format = formatter.getFormat(categoryName, f.name); if (format && typeof format.digitCount !== 'undefined') ret[ret.length] = f.type === Field.Type.Float ? Math.pow(10, Math.max(0, Math.min(format.digitCount, 15))) : 0; - else ret[ret.length] = f.type === Field.Type.Float ? Math.pow(10, Field.getDigitCount(f)) : 0; + else ret[ret.length] = f.type === Field.Type.Float ? Math.pow(10, getFieldDigitCount(f)) : 0; } return ret; } function writeCifSingleRecord(category: Category<any>, builder: StringBuilder, filter: Category.Filter, formatter: Category.Formatter) { - const fields = category.fields; + const fields = getIncludedFields(category); const data = category.data; let width = fields.reduce((w, f) => filter.includeField(category.name, f.name) ? Math.max(w, f.name.length) : 0, 0); @@ -134,8 +135,11 @@ function writeCifSingleRecord(category: Category<any>, builder: StringBuilder, f function writeCifLoop(categories: Category[], builder: StringBuilder, filter: Category.Filter, formatter: Category.Formatter) { const first = categories[0]; - const fields = filter === Category.DefaultFilter ? first.fields : first.fields.filter(f => filter.includeField(first.name, f.name)); + const fieldSource = getIncludedFields(first); + const fields = filter === Category.DefaultFilter ? fieldSource : fieldSource.filter(f => filter.includeField(first.name, f.name)); const fieldCount = fields.length; + if (fieldCount === 0) return; + const precisions = getFloatPrecisions(first.name, fields, formatter); writeLine(builder, 'loop_'); diff --git a/src/mol-io/writer/cif/encoder/util.ts b/src/mol-io/writer/cif/encoder/util.ts new file mode 100644 index 0000000000000000000000000000000000000000..af8b3369f61d11939cca26f24b3bae02a0de4a6a --- /dev/null +++ b/src/mol-io/writer/cif/encoder/util.ts @@ -0,0 +1,18 @@ +/** + * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. + * + * @author David Sehnal <david.sehnal@gmail.com> + */ + +import { Field, Category } from '../encoder'; + +export function getFieldDigitCount(field: Field) { + if (field.defaultFormat && typeof field.defaultFormat.digitCount !== 'undefined') return Math.max(0, Math.min(field.defaultFormat.digitCount, 16)); + return 6; +} + +export function getIncludedFields(category: Category<any, any>) { + return category.fields.some(f => !!f.shouldInclude) + ? category.fields.filter(f => !f.shouldInclude || f.shouldInclude(category.data)) + : category.fields; +} \ No newline at end of file diff --git a/src/mol-math/geometry/symmetry-operator.ts b/src/mol-math/geometry/symmetry-operator.ts index 8f9cac24028275ce2a8882e3859e1bd8c8e63c87..2e9776d8dd50ad3a642d1a08266ae4ced996bed2 100644 --- a/src/mol-math/geometry/symmetry-operator.ts +++ b/src/mol-math/geometry/symmetry-operator.ts @@ -18,7 +18,7 @@ interface SymmetryOperator { } namespace SymmetryOperator { - export const DefaultName = 'identity' + export const DefaultName = '1_555' export const Default: SymmetryOperator = create(DefaultName, Mat4.identity()); const RotationEpsilon = 0.0001; diff --git a/src/mol-model/structure/export/mmcif.ts b/src/mol-model/structure/export/mmcif.ts index 0a522303c91d4d124a9260cbb0c2031fed8d6e75..ea699f649fe5e900813d315af5cad2717b39bbb5 100644 --- a/src/mol-model/structure/export/mmcif.ts +++ b/src/mol-model/structure/export/mmcif.ts @@ -47,7 +47,9 @@ const atom_site_fields: CifField<Element.Location>[] = [ CifField.str('auth_asym_id', P.chain.auth_asym_id), CifField.int('pdbx_PDB_model_num', P.unit.model_num, { encoder: E.deltaRLE }), - CifField.str('operator_name', P.unit.operator_name) + CifField.str<Element.Location, Structure>('operator_name', P.unit.operator_name, { + shouldInclude: structure => structure.units.some(u => !u.conformation.operator.isIdentity) + }) ]; function copy_mmCif_cat(name: keyof mmCIF_Schema) { @@ -66,7 +68,7 @@ function _entity({ model, structure }: Context): CifCategory { function _atom_site({ structure }: Context): CifCategory { return { - data: void 0, + data: structure, name: 'atom_site', fields: atom_site_fields, rowCount: structure.elementCount, diff --git a/src/mol-model/structure/structure/unit.ts b/src/mol-model/structure/structure/unit.ts index 753d73eb211feef7b1b8921352a9daf9e2d5f7c4..e5ef8d9080ebf11d6b2b5b6e6f2a76353aee9a61 100644 --- a/src/mol-model/structure/structure/unit.ts +++ b/src/mol-model/structure/structure/unit.ts @@ -12,6 +12,7 @@ import { idFactory } from 'mol-util/id-factory'; import { IntraUnitBonds, computeIntraUnitBonds } from './unit/bonds' import { CoarseElements, CoarseSphereConformation, CoarseGaussianConformation } from '../model/properties/coarse'; import { ValueRef } from 'mol-util'; +import { UnitRings } from './unit/rings'; // A building block of a structure that corresponds to an atomic or a coarse grained representation // 'conveniently grouped together'. @@ -105,6 +106,12 @@ namespace Unit { return this.props.bonds.ref; } + get rings() { + if (this.props.rings.ref) return this.props.rings.ref; + this.props.rings.ref = UnitRings.create(this); + return this.props.rings.ref; + } + constructor(id: number, invariantId: number, model: Model, elements: SortedArray, conformation: SymmetryOperator.ArrayMapping, props: AtomicProperties) { this.id = id; this.invariantId = invariantId; @@ -121,10 +128,11 @@ namespace Unit { interface AtomicProperties { lookup3d: ValueRef<Lookup3D | undefined>, bonds: ValueRef<IntraUnitBonds | undefined>, + rings: ValueRef<UnitRings | undefined> } function AtomicProperties() { - return { lookup3d: ValueRef.create(void 0), bonds: ValueRef.create(void 0) }; + return { lookup3d: ValueRef.create(void 0), bonds: ValueRef.create(void 0), rings: ValueRef.create(void 0) }; } class Coarse<K extends Kind.Gaussians | Kind.Spheres, C extends CoarseSphereConformation | CoarseGaussianConformation> implements Base { diff --git a/src/mol-model/structure/structure/unit/rings.ts b/src/mol-model/structure/structure/unit/rings.ts new file mode 100644 index 0000000000000000000000000000000000000000..ed22f5418a9944c69b76099e3070d406ec6d4a9d --- /dev/null +++ b/src/mol-model/structure/structure/unit/rings.ts @@ -0,0 +1,107 @@ +/** + * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. + * + * @author David Sehnal <david.sehnal@gmail.com> + */ + +import computeRings from './rings/compute' +import Unit from '../unit'; + +interface UnitRings { + /** Each ring is specified as an array of indices in Unit.elements. */ + readonly all: ReadonlyArray<ReadonlyArray<number>>, + readonly byFingerprint: Map<string, ReadonlyArray<number>> +} + +namespace UnitRings { + export function getRingFingerprint(unit: Unit.Atomic, ring: ArrayLike<number>) { + const { elements } = unit; + const { type_symbol } = unit.model.atomicHierarchy.atoms; + + const symbols: string[] = []; + for (let i = 0, _i = ring.length; i < _i; i++) symbols[symbols.length] = type_symbol.value(elements[ring[i]]) as String as string; + return getFingerprint(symbols); + } + + export function create(unit: Unit.Atomic): UnitRings { + const rings = computeRings(unit); + const byFingerprint = new Map<string, number[]>(); + + let idx = 0; + for (const r of rings) { + const fp = getRingFingerprint(unit, r); + if (byFingerprint.has(fp)) byFingerprint.get(fp)!.push(idx); + else byFingerprint.set(fp, [idx]); + idx++; + } + + return { all: rings, byFingerprint }; + } +} + +export { UnitRings } + +function getFingerprint(elements: string[]) { + const len = elements.length; + const reversed: string[] = new Array(len); + + for (let i = 0; i < len; i++) reversed[i] = elements[len - i - 1]; + + const rotNormal = getMinimalRotation(elements); + const rotReversed = getMinimalRotation(reversed); + + let isNormalSmaller = false; + + for (let i = 0; i < len; i++) { + const u = elements[(i + rotNormal) % len], v = reversed[(i + rotReversed) % len]; + if (u !== v) { + isNormalSmaller = u < v; + break; + } + } + + if (isNormalSmaller) return buildFinderprint(elements, rotNormal); + return buildFinderprint(reversed, rotReversed); +} + +function getMinimalRotation(elements: string[]) { + // adapted from http://en.wikipedia.org/wiki/Lexicographically_minimal_string_rotation + + const len = elements.length; + const f = new Int32Array(len * 2); + for (let i = 0; i < f.length; i++) f[i] = -1; + + let u = '', v = '', k = 0; + + for (let j = 1; j < f.length; j++) { + let i = f[j - k - 1]; + while (i !== -1) { + u = elements[j % len]; v = elements[(k + i + 1) % len]; + if (u === v) break; + if (u < v) k = j - i - 1; + i = f[i]; + } + + if (i === -1) { + u = elements[j % len]; v = elements[(k + i + 1) % len]; + if (u !== v) { + if (u < v) k = j; + f[j - k] = -1; + } else f[j - k] = i + 1; + } else f[j - k] = i + 1; + } + + return k; +} + +function buildFinderprint(elements: string[], offset: number) { + const len = elements.length; + const ret: string[] = []; + let i; + for (i = 0; i < len - 1; i++) { + ret.push(elements[(i + offset) % len]); + ret.push('-'); + } + ret.push(elements[(i + offset) % len]); + return ret.join(''); +} \ No newline at end of file diff --git a/src/mol-model/structure/structure/unit/rings/compute.ts b/src/mol-model/structure/structure/unit/rings/compute.ts new file mode 100644 index 0000000000000000000000000000000000000000..a76cd6e47c29b67c4d57d2a7565f91a6cd3f2b76 --- /dev/null +++ b/src/mol-model/structure/structure/unit/rings/compute.ts @@ -0,0 +1,177 @@ +/** + * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. + * + * @author David Sehnal <david.sehnal@gmail.com> + */ + +import Unit from '../../unit'; +import { IntraUnitBonds } from '../bonds/intra-data'; +import { Segmentation } from 'mol-data/int'; + +export default function computeRings(unit: Unit.Atomic) { + const size = largestResidue(unit); + const state = State(unit, size); + + const residuesIt = Segmentation.transientSegments(unit.model.atomicHierarchy.residueSegments, unit.elements); + while (residuesIt.hasNext) { + const seg = residuesIt.move(); + processResidue(state, seg.start, seg.end); + } + + return state.rings; +} + +const enum Constants { + MaxDepth = 4 +} + +interface State { + startVertex: number, + endVertex: number, + count: number, + visited: Int32Array, + queue: Int32Array, + color: Int32Array, + pred: Int32Array, + + left: Int32Array, + right: Int32Array, + + currentColor: number, + + rings: number[][], + bonds: IntraUnitBonds, + unit: Unit.Atomic +} + +function State(unit: Unit.Atomic, capacity: number): State { + return { + startVertex: 0, + endVertex: 0, + count: 0, + visited: new Int32Array(capacity), + queue: new Int32Array(capacity), + pred: new Int32Array(capacity), + left: new Int32Array(Constants.MaxDepth), + right: new Int32Array(Constants.MaxDepth), + color: new Int32Array(capacity), + currentColor: 0, + rings: [], + unit, + bonds: unit.bonds + }; +} + +function resetState(state: State) { + state.count = state.endVertex - state.startVertex; + const { visited, pred, color } = state; + for (let i = 0; i < state.count; i++) { + visited[i] = -1; + pred[i] = -1; + color[i] = 0; + } + state.currentColor = 0; +} + +function largestResidue(unit: Unit.Atomic) { + const residuesIt = Segmentation.transientSegments(unit.model.atomicHierarchy.residueSegments, unit.elements); + let size = 0; + while (residuesIt.hasNext) { + const seg = residuesIt.move(); + size = Math.max(size, seg.end - seg.start); + } + return size; +} + +function processResidue(state: State, start: number, end: number) { + const { visited } = state; + state.startVertex = start; + state.endVertex = end; + + // no two atom rings + if (state.endVertex - state.startVertex < 3) return; + + resetState(state); + + for (let i = 0; i < state.count; i++) { + if (visited[i] >= 0) continue; + findRings(state, i); + } +} + +function addRing(state: State, a: number, b: number) { + // only "monotonous" rings + if (b < a) return; + + const { pred, color, left, right } = state; + const nc = ++state.currentColor; + + let current = a; + + for (let t = 0; t < Constants.MaxDepth; t++) { + color[current] = nc; + current = pred[current]; + if (current < 0) break; + } + + let leftOffset = 0, rightOffset = 0; + + let found = false, target = 0; + current = b; + for (let t = 0; t < Constants.MaxDepth; t++) { + if (color[current] === nc) { + target = current; + found = true; + break; + } + right[rightOffset++] = current; + current = pred[current]; + if (current < 0) break; + } + if (!found) return; + + current = a; + for (let t = 0; t < Constants.MaxDepth; t++) { + left[leftOffset++] = current; + if (target === current) break; + current = pred[current]; + if (current < 0) break; + } + + const ring = new Int32Array(leftOffset + rightOffset); + let ringOffset = 0; + for (let t = 0; t < leftOffset; t++) ring[ringOffset++] = state.startVertex + left[t]; + for (let t = rightOffset - 1; t >= 0; t--) ring[ringOffset++] = state.startVertex + right[t]; + + state.rings.push(ring as any as number[]); +} + +function findRings(state: State, from: number) { + const { bonds, startVertex, endVertex, visited, queue, pred } = state; + const { b: neighbor, flags: bondFlags, offset } = bonds; + visited[from] = 1; + queue[0] = from; + let head = 0, size = 1; + + while (head < size) { + const top = queue[head++]; + const a = startVertex + top; + const start = offset[a], end = offset[a + 1]; + + for (let i = start; i < end; i++) { + const b = neighbor[i]; + if (b < startVertex || b >= endVertex || !IntraUnitBonds.isCovalent(bondFlags[i])) continue; + + const other = b - startVertex; + + if (visited[other] > 0) { + if (pred[other] !== top && pred[top] !== other) addRing(state, top, other); + continue; + } + + visited[other] = 1; + queue[size++] = other; + pred[other] = top; + } + } +} \ No newline at end of file