assembly parsing from mmCIF

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

@@ -24,7 +24,7 @@ namespace SymmetryOperator {
     export function create(name: string, matrix: Mat4, hkl?: number[]): SymmetryOperator {
         const _hkl = hkl ? Vec3.create(hkl[0], hkl[1], hkl[2]) : Vec3.zero();
         if (Mat4.isIdentity(matrix)) return { name, matrix, inverse: Mat4.identity(), isIdentity: true, hkl: _hkl };
-        if (!Mat4.isRotationAndTranslation(matrix)) throw new Error('Symmetry operator must be a composition of rotation and translation.');
+        if (!Mat4.isRotationAndTranslation(matrix)) throw new Error(`Symmetry operator (${name}) must be a composition of rotation and translation.`);
         return { name, matrix, inverse: Mat4.invert(Mat4.zero(), matrix), isIdentity: false, hkl: _hkl };
     }
 

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

@@ -214,7 +214,7 @@ export namespace Mat4 {
         return mul(out, mul(out, a, b), c);
     }
 
-    export function translate(out: Mat4, a: Mat4, v: Mat4) {
+    export function translate(out: Mat4, a: Mat4, v: Vec3) {
         const x = v[0], y = v[1], z = v[2];
         let a00: number, a01: number, a02: number, a03: number,
             a10: number, a11: number, a12: number, a13: number,
@@ -243,7 +243,7 @@ export namespace Mat4 {
         return out;
     }
 
-    export function fromTranslation(out: Mat4, v: Mat4) {
+    export function fromTranslation(out: Mat4, v: Vec3) {
         out[0] = 1;
         out[1] = 0;
         out[2] = 0;
@@ -263,6 +263,13 @@ export namespace Mat4 {
         return out;
     }
 
+    export function setTranslation(out: Mat4, v: Vec3) {
+        out[12] = v[0];
+        out[13] = v[1];
+        out[14] = v[2];
+        return out;
+    }
+
     export function rotate(out: Mat4, a: Mat4, rad: number, axis: Mat4) {
         let x = axis[0], y = axis[1], z = axis[2],
             len = Math.sqrt(x * x + y * y + z * z),
@@ -432,15 +439,23 @@ export namespace Mat4 {
         return b00 * b11 - b01 * b10 + b02 * b09 + b03 * b08 - b04 * b07 + b05 * b06;
     }
 
-    export function isRotationAndTranslation(a: Mat4) {
+    export function isRotationAndTranslation(a: Mat4, eps?: number) {
+        return _isRotationAndTranslation(a, typeof eps !== 'undefined' ? eps : EPSILON.Value)
+    }
+
+    function _isRotationAndTranslation(a: Mat4, eps: number) {
         const a00 = a[0], a01 = a[1], a02 = a[2], a03 = a[3],
             a10 = a[4], a11 = a[5], a12 = a[6], a13 = a[7],
             a20 = a[8], a21 = a[9], a22 = a[10], a23 = a[11],
             /* a30 = a[12], a31 = a[13], a32 = a[14],*/ a33 = a[15];
 
-        if (a33 !== 1 || a03 !== 0 || a13 !== 0 || a23 !== 0) return false;
-        const det3x3 = a00 * (a11 * a22 - a21 * a21) - a01 * (a10 * a22 - a12 * a20) + a02 * (a10 * a21 - a11 * a20);
-        if (det3x3 < 1 - EPSILON.Value || det3x3 > 1 + EPSILON.Value) return false;
+        if (a33 !== 1 || a03 !== 0 || a13 !== 0 || a23 !== 0) {
+            return false;
+        }
+        const det3x3 = a00 * (a11 * a22 - a12 * a21) - a01 * (a10 * a22 - a12 * a20) + a02 * (a10 * a21 - a11 * a20);
+        if (det3x3 < 1 - eps || det3x3 > 1 + eps) {
+            return false;
+        }
         return true;
     }
 }

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

@@ -14,6 +14,7 @@ import Conformation from '../properties/conformation'
 import Symmetry from '../properties/symmetry'
 import findHierarchyKeys from '../utils/hierarchy-keys'
 import { ElementSymbol} from '../types'
+import createAssemblies from './mmcif/assembly'
 
 import mmCIF_Format = Format.mmCIF
 
@@ -77,8 +78,8 @@ function getConformation({ data }: mmCIF_Format, bounds: Interval): Conformation
     }
 }
 
-function getSymmetry({ data }: mmCIF_Format): Symmetry {
-    return Symmetry.Empty;
+function getSymmetry(format: mmCIF_Format): Symmetry {
+    return { assemblies: createAssemblies(format) };
 }
 
 function isHierarchyDataEqual(a: Hierarchy.Hierarchy, b: Hierarchy.Data) {

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

@@ -0,0 +1,149 @@
+/**
+ * Copyright (c) 2017 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author David Sehnal <david.sehnal@gmail.com>
+ */
+
+import { Mat4, Tensor } from 'mol-math/linear-algebra'
+import SymmetryOperator from 'mol-math/geometry/symmetry-operator'
+import Format from '../../format'
+import { Assembly, OperatorGroup, OperatorGroups } from '../../properties/symmetry'
+import { Queries as Q } from '../../../query'
+
+import mmCIF_Format = Format.mmCIF
+
+export default function create(format: mmCIF_Format): ReadonlyArray<Assembly> {
+    const { pdbx_struct_assembly } = format.data;
+    if (!pdbx_struct_assembly._rowCount) return [];
+
+    const matrices = getMatrices(format);
+    const assemblies: Assembly[] = [];
+    for (let i = 0; i < pdbx_struct_assembly._rowCount; i++) {
+        assemblies[assemblies.length] = createAssembly(format, i, matrices);
+    }
+    return assemblies;
+}
+
+type Matrices = Map<string, Mat4>
+type Generator = { expression: string, asymIds: string[] }
+
+function createAssembly(format: mmCIF_Format, index: number, matrices: Matrices): Assembly {
+    const { pdbx_struct_assembly, pdbx_struct_assembly_gen } = format.data;
+
+    const id = pdbx_struct_assembly.id.value(index);
+    const details = pdbx_struct_assembly.details.value(index);
+    const generators: Generator[] = [];
+
+    const { assembly_id, oper_expression, asym_id_list } = pdbx_struct_assembly_gen;
+
+    for (let i = 0, _i = pdbx_struct_assembly_gen._rowCount; i < _i; i++) {
+        if (assembly_id.value(i) !== id) continue;
+        generators[generators.length] = {
+            expression: oper_expression.value(i),
+            asymIds: asym_id_list.value(i).split(',').map(x => x.trim()).filter(x => !!x)
+        };
+    }
+
+    return Assembly.create(id, details, operatorGroupsProvider(generators, matrices));
+}
+
+function operatorGroupsProvider(generators: Generator[], matrices: Matrices): () => OperatorGroups {
+    return () => {
+        const groups: OperatorGroup[] = [];
+
+        let operatorOffset = 0;
+        for (let i = 0; i < generators.length; i++) {
+            const gen = generators[i];
+            const operatorList = parseOperatorList(gen.expression);
+            const operatorNames = expandOperators(operatorList);
+            const operators = getAssemblyOperators(matrices, operatorNames, operatorOffset);
+            const selector = Q.generators.atoms({ chainTest: Q.pred.and(
+                Q.pred.eq(Q.props.unit.operator_name, SymmetryOperator.DefaultName),
+                Q.pred.inSet(Q.props.chain.label_asym_id, gen.asymIds)
+            )});
+            groups[groups.length] = { selector, operators };
+            operatorOffset += operators.length;
+        }
+
+        return groups;
+    }
+}
+
+function getMatrices({ data }: mmCIF_Format): Matrices {
+    const { pdbx_struct_oper_list } = data;
+    const { id, matrix, vector, _schema } = pdbx_struct_oper_list;
+    const matrices = new Map<string, Mat4>();
+
+    for (let i = 0, _i = pdbx_struct_oper_list._rowCount; i < _i; i++) {
+        const m = Tensor.toMat4(_schema.matrix.space, matrix.value(i));
+        const t = Tensor.toVec3(_schema.vector.space, vector.value(i));
+        Mat4.setTranslation(m, t);
+        Mat4.setValue(m, 3, 3, 1);
+        matrices.set(id.value(i), m);
+    }
+
+    return matrices;
+}
+
+function expandOperators(operatorList: string[][]) {
+    const ops: string[][] = [];
+    const currentOp: string[] = [];
+    for (let i = 0; i < operatorList.length; i++) currentOp[i] = '';
+    expandOperators1(operatorList, ops, operatorList.length - 1, currentOp);
+    return ops;
+}
+
+function expandOperators1(operatorNames: string[][], list: string[][], i: number, current: string[]) {
+    if (i < 0) {
+        list[list.length] = current.slice(0);
+        return;
+    }
+
+    let ops = operatorNames[i], len = ops.length;
+    for (let j = 0; j < len; j++) {
+        current[i] = ops[j];
+        expandOperators1(operatorNames, list, i - 1, current);
+    }
+}
+
+function getAssemblyOperators(matrices: Matrices, operatorNames: string[][], startIndex: number) {
+    const operators: SymmetryOperator[] = [];
+
+    let index = startIndex;
+    for (let op of operatorNames) {
+        let m = Mat4.identity();
+        for (let i = 0; i < op.length; i++) {
+            Mat4.mul(m, m, matrices.get(op[i])!);
+        }
+        index++;
+        operators[operators.length] = SymmetryOperator.create(`A-${index}`, m);
+    }
+
+    return operators;
+}
+
+function parseOperatorList(value: string): string[][] {
+    // '(X0)(1-5)' becomes [['X0'], ['1', '2', '3', '4', '5']]
+    // kudos to Glen van Ginkel.
+
+    const oeRegex = /\(?([^\(\)]+)\)?]*/g, groups: string[] = [], ret: string[][] = [];
+
+    let g: any;
+    while (g = oeRegex.exec(value)) groups[groups.length] = g[1];
+
+    groups.forEach(g => {
+        const group: string[] = [];
+        g.split(',').forEach(e => {
+            const dashIndex = e.indexOf('-');
+            if (dashIndex > 0) {
+                const from = parseInt(e.substring(0, dashIndex)), to = parseInt(e.substr(dashIndex + 1));
+                for (let i = from; i <= to; i++) group[group.length] = i.toString();
+            } else {
+                group[group.length] = e.trim();
+            }
+        });
+        ret[ret.length] = group;
+    });
+
+    return ret;
+}

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

@@ -8,12 +8,13 @@ import SymmetryOperator from 'mol-math/geometry/symmetry-operator'
 import { Query } from '../../query'
 import { Model } from '../../model'
 
+/** Determine an atom set and a list of operators that should be applied to that set  */
 export interface OperatorGroup {
-    readonly query: Query,
+    readonly selector: Query,
     readonly operators: ReadonlyArray<SymmetryOperator>
 }
 
-export type OperatorGroups = ReadonlyArray<ReadonlyArray<OperatorGroup>>
+export type OperatorGroups = ReadonlyArray<OperatorGroup>
 
 export class Assembly {
     readonly id: string;
@@ -33,7 +34,7 @@ export class Assembly {
 }
 
 export namespace Assembly {
-    export function create(id: string, details: string, operatorsProvider: () => OperatorGroups) {
+    export function create(id: string, details: string, operatorsProvider: () => OperatorGroups): Assembly {
         return new Assembly(id, details, operatorsProvider);
     }
 }
@@ -45,7 +46,12 @@ interface Symmetry {
 namespace Symmetry {
     export const Empty: Symmetry = { assemblies: [] };
 
-    export function findAssembly(model: Model): Assembly | undefined {
+    export function findAssembly(model: Model, id: string): Assembly | undefined {
+        const { assemblies } = model.symmetry;
+        const _id = id.toLocaleLowerCase();
+        for (let i = 0; i < assemblies.length; i++) {
+            if (assemblies[i].id.toLowerCase() === _id) return assemblies[i];
+        }
         return void 0;
     }
 }

src/perf-tests/structure.ts

@@ -242,12 +242,16 @@ export namespace PropertyAccess {
     }
 
     export async function run() {
-        const { structures, models, mmcif } = await readCIF('./examples/1cbs_full.bcif');
+        //const { structures, models, mmcif } = await readCIF('./examples/1cbs_full.bcif');
         //const { structures, models, mmcif } = await readCIF('e:/test/quick/3j3q_full.bcif');
         //const { structures, models, mmcif } = await readCIF('e:/test/quick/1cbs_updated.cif');
+        const { structures, models/*, mmcif*/ } = await readCIF('e:/test/quick/1hrv_updated.cif');
+
+        //console.log(mmcif.pdbx_struct_oper_list.matrix.toArray());
+        // console.log(mmcif.pdbx_struct_oper_list.vector.toArray());
+
+        console.log(models[0].symmetry.assemblies);
 
-        console.log(mmcif.pdbx_struct_oper_list.matrix.toArray());
-        console.log(mmcif.pdbx_struct_oper_list.vector.toArray());
 
         //const { structures, models } = await readCIF('e:/test/molstar/3j3q.bcif');