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

data/mmcif-field-names.csv

@@ -144,6 +144,12 @@ struct_keywords.entry_id
 struct_keywords.pdbx_keywords
 struct_keywords.text
 
+struct_ncs_oper.id
+struct_ncs_oper.code
+struct_ncs_oper.matrix
+struct_ncs_oper.vector
+struct_ncs_oper.details
+
 struct_sheet_range.sheet_id
 struct_sheet_range.id
 struct_sheet_range.beg_label_comp_id

src/apps/structure-info/model.ts

@@ -18,6 +18,7 @@ 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';
+import { Vec3 } from 'mol-math/linear-algebra';
 
 
 async function downloadFromPdb(pdb: string) {
@@ -182,6 +183,16 @@ export function printUnits(structure: Structure) {
     }
 }
 
+export function printSymmetryInfo(model: Model) {
+    console.log('\nSymmetry Info\n=============');
+    const { symmetry } = model;
+    const { size, anglesInRadians } = symmetry.spacegroup.cell;
+    console.log(`Spacegroup: ${symmetry.spacegroup.name} size: ${Vec3.toString(size)} angles: ${Vec3.toString(anglesInRadians)}`);
+    console.log(`Assembly names: ${symmetry.assemblies.map(a => a.id).join(', ')}`);
+    // NCS example: 1auy
+    console.log(`NCS operators: ${symmetry.ncsOperators && symmetry.ncsOperators.map(a => a.name).join(', ')}`);
+}
+
 export function printIHMModels(model: Model) {
     if (!model.coarseHierarchy.isDefined) return false;
     console.log('\nIHM Models\n=============');
@@ -194,6 +205,7 @@ async function run(frame: CifFrame) {
     //printSequence(models[0]);
     //printIHMModels(models[0]);
     printUnits(structure);
+    printSymmetryInfo(models[0]);
     //printRings(structure);
     //printLinks(structure, true, true);
     //printModRes(models[0]);

src/mol-io/reader/cif/schema/mmcif.ts

@@ -409,6 +409,13 @@ export const mmCIF_Schema = {
         experiment_type: Aliased<'Fraction of bulk' | 'Single molecule'>(str),
         details: str,
     },
+    struct_ncs_oper: {
+        id: str,
+        code: str,
+        matrix: Matrix(3, 3),
+        vector: Vector(3),
+        details: str
+    },
     ihm_modeling_post_process: {
         id: int,
         protocol_id: int,

src/mol-math/geometry/symmetry-operator.ts

@@ -4,7 +4,7 @@
  * @author David Sehnal <david.sehnal@gmail.com>
  */
 
-import { Vec3, Mat4 } from '../linear-algebra/3d'
+import { Vec3, Mat4, Mat3 } from '../linear-algebra/3d'
 
 interface SymmetryOperator {
     readonly name: string,
@@ -30,6 +30,17 @@ namespace SymmetryOperator {
         return { name, matrix, inverse: Mat4.invert(Mat4.zero(), matrix), isIdentity: false, hkl: _hkl };
     }
 
+    export function ofRotationAndOffset(name: string, rot: Mat3, offset: Vec3) {
+        const t = Mat4.identity();
+        for (let i = 0; i < 3; i++) {
+            for (let j = 0; j < 3; j++) {
+                Mat4.setValue(t, i, j, Mat3.getValue(rot, i, j));
+            }
+        }
+        Mat4.setTranslation(t, offset);
+        return create(name, t);
+    }
+
     // Apply the 1st and then 2nd operator. ( = second.matrix * first.matrix)
     export function compose(first: SymmetryOperator, second: SymmetryOperator) {
         const matrix = Mat4.mul(Mat4.zero(), second.matrix, first.matrix);

src/mol-math/linear-algebra/3d/mat3.ts

@@ -108,6 +108,10 @@ namespace Mat3 {
         a[3 * j + i] = value;
     }
 
+    export function getValue(a: Mat3, i: number, j: number) {
+        return a[3 * j + i];
+    }
+
     /**
      * Copy the values from one Mat3 to another
      */

src/mol-math/linear-algebra/3d/vec3.ts

@@ -367,6 +367,10 @@ namespace Vec3 {
     export function isZero(v: Vec3) {
         return v[0] === 0 && v[1] === 0 && v[2] === 0
     }
+
+    export function toString(a: Vec3) {
+        return `[${a[0]} ${a[1]} ${a[2]}]`;
+    }
 }
 
 export default Vec3
\ No newline at end of file

src/mol-model/structure/model/formats/mmcif.ts

@@ -6,8 +6,8 @@
 
 import { Column, Table } from 'mol-data/db';
 import { Interval, Segmentation } from 'mol-data/int';
-import { Spacegroup, SpacegroupCell } from 'mol-math/geometry';
-import { Vec3 } from 'mol-math/linear-algebra';
+import { Spacegroup, SpacegroupCell, SymmetryOperator } from 'mol-math/geometry';
+import { Vec3, Tensor, Mat4 } from 'mol-math/linear-algebra';
 import { Task } from 'mol-task';
 import UUID from 'mol-util/uuid';
 import Format from '../format';
@@ -22,7 +22,7 @@ import { getIHMCoarse } from './mmcif/ihm';
 import { getSecondaryStructureMmCif } from './mmcif/secondary-structure';
 import { getSequence } from './mmcif/sequence';
 import { sortAtomSite } from './mmcif/sort';
-import { mmCIF_Database } from 'mol-io/reader/cif/schema/mmcif';
+import { mmCIF_Database, mmCIF_Schema } from 'mol-io/reader/cif/schema/mmcif';
 
 import mmCIF_Format = Format.mmCIF
 type AtomSite = mmCIF_Database['atom_site']
@@ -90,7 +90,7 @@ function getSymmetry(format: mmCIF_Format): ModelSymmetry {
     const assemblies = createAssemblies(format);
     const spacegroup = getSpacegroup(format);
     const isNonStandardCrytalFrame = checkNonStandardCrystalFrame(format, spacegroup);
-    return { assemblies, spacegroup, isNonStandardCrytalFrame };
+    return { assemblies, spacegroup, isNonStandardCrytalFrame, ncsOperators: getNcsOperators(format) };
 }
 
 function checkNonStandardCrystalFrame(format: mmCIF_Format, spacegroup: Spacegroup) {
@@ -111,6 +111,22 @@ function getSpacegroup(format: mmCIF_Format): Spacegroup {
     return Spacegroup.create(spaceCell);
 }
 
+function getNcsOperators(format: mmCIF_Format) {
+    const { struct_ncs_oper } = format.data;
+    if (struct_ncs_oper._rowCount === 0) return void 0;
+    const { id, matrix, vector } = struct_ncs_oper;
+
+    const matrixSpace = mmCIF_Schema.struct_ncs_oper.matrix.space, vectorSpace = mmCIF_Schema.struct_ncs_oper.vector.space;
+
+    const opers: SymmetryOperator[] = [];
+    for (let i = 0; i < struct_ncs_oper._rowCount; i++) {
+        const m = Tensor.toMat3(matrixSpace, matrix.value(i));
+        const v = Tensor.toVec3(vectorSpace, vector.value(i));
+        opers[i] = SymmetryOperator.ofRotationAndOffset(`ncs_${id.value(i)}`, m, v);
+    }
+    return opers;
+}
+
 function isHierarchyDataEqual(a: AtomicData, b: AtomicData) {
     // need to cast because of how TS handles type resolution for interfaces https://github.com/Microsoft/TypeScript/issues/15300
     return Table.areEqual(a.chains as Table<ChainsSchema>, b.chains as Table<ChainsSchema>)

src/mol-model/structure/model/properties/symmetry.ts

@@ -44,7 +44,11 @@ export namespace Assembly {
 interface ModelSymmetry {
     readonly assemblies: ReadonlyArray<Assembly>,
     readonly spacegroup: Spacegroup,
-    readonly isNonStandardCrytalFrame: boolean
+    readonly isNonStandardCrytalFrame: boolean,
+    readonly ncsOperators?: ReadonlyArray<SymmetryOperator>,
+
+    // optionally cached operators from [-3, -3, -3] to [3, 3, 3]
+    _operators_333?: SymmetryOperator[]
 }
 
 namespace ModelSymmetry {

src/mol-model/structure/structure/symmetry.ts

@@ -7,7 +7,7 @@
 import Structure from './structure'
 import { Selection } from '../query'
 import { ModelSymmetry } from '../model'
-import { Task } from 'mol-task';
+import { Task, RuntimeContext } from 'mol-task';
 import { SortedArray } from 'mol-data/int';
 import Unit from './unit';
 import { EquivalenceClasses, hash2 } from 'mol-data/util';
@@ -44,41 +44,16 @@ namespace StructureSymmetry {
     }
 
     export function builderSymmetryMates(structure: Structure, radius: number) {
-        // TODO: do it properly
-        return buildSymmetryRange(structure, Vec3.create(-3, -3, -3), Vec3.create(3, 3, 3));
+        return Task.create('Find Symmetry Mates', ctx => findMatesRadius(ctx, structure, radius));
     }
 
     export function buildSymmetryRange(structure: Structure, ijkMin: Vec3, ijkMax: Vec3) {
-        return Task.create('Build Assembly', async ctx => {
-            const models = Structure.getModels(structure);
-            if (models.length !== 1) throw new Error('Can only build symmetries from structures based on 1 model.');
-
-            const { spacegroup } = models[0].symmetry;
-            if (SpacegroupCell.isZero(spacegroup.cell)) return structure;
-
-            const operators: SymmetryOperator[] = [];
-            for (let op = 0; op < spacegroup.operators.length; op++) {
-                for (let i = ijkMin[0]; i < ijkMax[0]; i++) {
-                    for (let j = ijkMin[1]; j < ijkMax[1]; j++) {
-                        for (let k = ijkMin[2]; k < ijkMax[2]; k++) {
-                            operators[operators.length] = Spacegroup.getSymmetryOperator(spacegroup, op, i, j, k);
-                        }
-                    }
-                }
-            }
-
-            const assembler = Structure.Builder();
-
-            const { units } = structure;
-            for (const oper of operators) {
-                for (const unit of units) {
-                    assembler.addWithOperator(unit, oper);
-                }
-            }
-
+        return Task.create('Build Symmetry', ctx => findSymmetryRange(ctx, structure, ijkMin, ijkMax));
+    }
 
-            return assembler.getStructure();
-        });
+    /** Builds NCS structure, returns the original if NCS operators are not present. */
+    export function buildNcs(structure: Structure) {
+        return Task.create('Build NCS', ctx => _buildNCS(ctx, structure));
     }
 
     function hashUnit(u: Unit) {
@@ -107,4 +82,91 @@ namespace StructureSymmetry {
     }
 }
 
+function getOperators(symmetry: ModelSymmetry, ijkMin: Vec3, ijkMax: Vec3) {
+    const operators: SymmetryOperator[] = symmetry._operators_333 || [];
+    const { spacegroup } = symmetry;
+    if (operators.length === 0) {
+        operators[0] = Spacegroup.getSymmetryOperator(spacegroup, 0, 0, 0, 0)
+        for (let op = 0; op < spacegroup.operators.length; op++) {
+            for (let i = ijkMin[0]; i < ijkMax[0]; i++) {
+                for (let j = ijkMin[1]; j < ijkMax[1]; j++) {
+                    for (let k = ijkMin[2]; k < ijkMax[2]; k++) {
+                        // we have added identity as the 1st operator.
+                        if (op === 0 && i === 0 && j === 0 && k === 0) continue;
+                        operators[operators.length] = Spacegroup.getSymmetryOperator(spacegroup, op, i, j, k);
+                    }
+                }
+            }
+        }
+        symmetry._operators_333 = operators;
+    }
+    return operators;
+}
+
+function assembleOperators(structure: Structure, operators: ReadonlyArray<SymmetryOperator>) {
+    const assembler = Structure.Builder();
+    const { units } = structure;
+    for (const oper of operators) {
+        for (const unit of units) {
+            assembler.addWithOperator(unit, oper);
+        }
+    }
+    return assembler.getStructure();
+}
+
+async function _buildNCS(ctx: RuntimeContext, structure: Structure) {
+    const models = Structure.getModels(structure);
+    if (models.length !== 1) throw new Error('Can only build NCS from structures based on 1 model.');
+
+    const operators = models[0].symmetry.ncsOperators;
+    if (!operators || !operators.length) return structure;
+    return assembleOperators(structure, operators);
+}
+
+async function findSymmetryRange(ctx: RuntimeContext, structure: Structure, ijkMin: Vec3, ijkMax: Vec3) {
+    const models = Structure.getModels(structure);
+    if (models.length !== 1) throw new Error('Can only build symmetries from structures based on 1 model.');
+
+    const { spacegroup } = models[0].symmetry;
+    if (SpacegroupCell.isZero(spacegroup.cell)) return structure;
+
+    const operators = getOperators(models[0].symmetry, ijkMin, ijkMax);
+    return assembleOperators(structure, operators);
+}
+
+async function findMatesRadius(ctx: RuntimeContext, structure: Structure, radius: number) {
+    const models = Structure.getModels(structure);
+    if (models.length !== 1) throw new Error('Can only build symmetries from structures based on 1 model.');
+
+    const symmetry = models[0].symmetry;
+    const { spacegroup } = symmetry;
+    if (SpacegroupCell.isZero(spacegroup.cell)) return structure;
+
+    if (ctx.shouldUpdate) await ctx.update('Initialing...');
+    const operators = getOperators(symmetry, Vec3.create(-3, -3, -3), Vec3.create(3, 3, 3));
+    const lookup = structure.lookup3d;
+
+    const assembler = Structure.Builder();
+
+    const { units } = structure;
+    const center = Vec3.zero();
+    for (const oper of operators) {
+        for (const unit of units) {
+            const boundingSphere = unit.lookup3d.boundary.sphere;
+            Vec3.transformMat4(center, boundingSphere.center, oper.matrix);
+
+            const closeUnits = lookup.findUnitIndices(center[0], center[1], center[2], boundingSphere.radius + radius);
+            for (let uI = 0, _uI = closeUnits.count; uI < _uI; uI++) {
+                const closeUnit = units[closeUnits.indices[uI]];
+                if (!closeUnit.lookup3d.check(center[0], center[1], center[2], boundingSphere.radius + radius)) continue;
+                assembler.addWithOperator(unit, oper);
+            }
+        }
+        if (ctx.shouldUpdate) await ctx.update('Building symmetry...');
+    }
+
+
+    return assembler.getStructure();
+}
+
 export default StructureSymmetry;
\ No newline at end of file

src/mol-task/execution/observable.ts

@@ -133,9 +133,9 @@ function abortTree(root: Progress.Node) {
     for (const c of root.children) abortTree(c);
 }
 
-function shouldNotify(info: ProgressInfo, time: number) {
-    return time - info.lastNotified > info.updateRateMs;
-}
+// function shouldNotify(info: ProgressInfo, time: number) {
+//     return time - info.lastNotified > info.updateRateMs;
+// }
 
 function notifyObserver(info: ProgressInfo, time: number) {
     info.lastNotified = time;
@@ -194,6 +194,7 @@ class ObservableRuntimeContext implements RuntimeContext {
         this.lastUpdatedTime = now();
         this.updateProgress(progress);
 
+        // TODO: do the shouldNotify check here?
         if (!!dontNotify /*|| !shouldNotify(this.info, this.lastUpdatedTime)*/) return;
 
         notifyObserver(this.info, this.lastUpdatedTime);