diff --git a/src/mol-math/geometry/spacegroup/construction.ts b/src/mol-math/geometry/spacegroup/construction.ts index dee2c694c6f9ab60647d6db90ac8097c0f2505de..2a49c7843fdefb24128d51425d5f5615980fada8 100644 --- a/src/mol-math/geometry/spacegroup/construction.ts +++ b/src/mol-math/geometry/spacegroup/construction.ts @@ -1,7 +1,8 @@ /** - * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. + * Copyright (c) 2018-2019 mol* contributors, licensed under MIT, See LICENSE file for more info. * * @author David Sehnal <david.sehnal@gmail.com> + * @author Alexander Rose <alexander.rose@weirdbyte.de> */ import { Vec3, Mat4 } from '../../linear-algebra' @@ -73,7 +74,7 @@ namespace SpacegroupCell { namespace Spacegroup { - // P1 with [1, 1, 1] cell. + /** P1 with [1, 1, 1] cell */ export const ZeroP1 = create(SpacegroupCell.Zero); export function create(cell: SpacegroupCell): Spacegroup { @@ -81,18 +82,56 @@ namespace Spacegroup { return { name: SpacegroupNames[cell.index], cell, operators }; } - const _tempVec = Vec3.zero(), _tempMat = Mat4.zero(); - export function updateOperatorMatrix(spacegroup: Spacegroup, index: number, i: number, j: number, k: number, target: Mat4) { - _tempVec[0] = i; - _tempVec[1] = j; - _tempVec[2] = k; - - Mat4.fromTranslation(_tempMat, _tempVec); - return Mat4.mul(target, Mat4.mul(target, Mat4.mul(target, spacegroup.cell.fromFractional, _tempMat), spacegroup.operators[index]), spacegroup.cell.toFractional); + const _ijkVec = Vec3(); + const _tempMat = Mat4(); + export function setOperatorMatrix(spacegroup: Spacegroup, index: number, i: number, j: number, k: number, target: Mat4) { + Vec3.set(_ijkVec, i, j, k); + + Mat4.fromTranslation(_tempMat, _ijkVec); + return Mat4.mul( + target, + Mat4.mul( + target, + Mat4.mul(target, spacegroup.cell.fromFractional, _tempMat), + spacegroup.operators[index] + ), + spacegroup.cell.toFractional + ); } export function getSymmetryOperator(spacegroup: Spacegroup, index: number, i: number, j: number, k: number): SymmetryOperator { - const operator = updateOperatorMatrix(spacegroup, index, i, j, k, Mat4.zero()); + const operator = setOperatorMatrix(spacegroup, index, i, j, k, Mat4.zero()); + return SymmetryOperator.create(`${index + 1}_${5 + i}${5 + j}${5 + k}`, operator, { id: '', operList: [] }, '', Vec3.create(i, j, k), index); + } + + const _translationRef = Vec3() + const _translationRefSymop = Vec3() + const _translationSymop = Vec3() + export function setOperatorMatrixRef(spacegroup: Spacegroup, index: number, i: number, j: number, k: number, ref: Vec3, target: Mat4) { + Vec3.set(_ijkVec, i, j, k); + Vec3.floor(_translationRef, ref) + + Mat4.copy(target, spacegroup.operators[index]) + + Vec3.floor(_translationRefSymop, Vec3.transformMat4(_translationRefSymop, ref, target)) + + Mat4.getTranslation(_translationSymop, target) + Vec3.sub(_translationSymop, _translationSymop, _translationRefSymop) + Vec3.add(_translationSymop, _translationSymop, _translationRef) + Vec3.add(_translationSymop, _translationSymop, _ijkVec) + + Mat4.setTranslation(target, _translationSymop) + Mat4.mul(target, spacegroup.cell.fromFractional, target) + Mat4.mul(target, target, spacegroup.cell.toFractional) + return target + } + + /** + * Get Symmetry operator for transformation around the given + * reference point `ref` in fractional coordinates + */ + export function getSymmetryOperatorRef(spacegroup: Spacegroup, index: number, i: number, j: number, k: number, ref: Vec3) { + const operator = setOperatorMatrixRef(spacegroup, index, i, j, k, ref, Mat4.zero()); return SymmetryOperator.create(`${index + 1}_${5 + i}${5 + j}${5 + k}`, operator, { id: '', operList: [] }, '', Vec3.create(i, j, k), index); } diff --git a/src/mol-model/structure/model/properties/symmetry.ts b/src/mol-model/structure/model/properties/symmetry.ts index ee3a9c5e28103a0df624f6ad7513291c6389b937..38b9e5e12642ac80ad7c8e5ae1daa661bd6334a9 100644 --- a/src/mol-model/structure/model/properties/symmetry.ts +++ b/src/mol-model/structure/model/properties/symmetry.ts @@ -1,5 +1,5 @@ /** - * Copyright (c) 2017 mol* contributors, licensed under MIT, See LICENSE file for more info. + * Copyright (c) 2017-2019 mol* contributors, licensed under MIT, See LICENSE file for more info. * * @author David Sehnal <david.sehnal@gmail.com> */ @@ -9,6 +9,7 @@ import { arrayFind } from '../../../../mol-data/util' import { StructureQuery } from '../../query' import { Model } from '../../model' import { Spacegroup } from '../../../../mol-math/geometry'; +import { Vec3 } from '../../../../mol-math/linear-algebra'; /** Determine an atom set and a list of operators that should be applied to that set */ export interface OperatorGroup { @@ -47,8 +48,14 @@ interface ModelSymmetry { readonly isNonStandardCrytalFrame: boolean, readonly ncsOperators?: ReadonlyArray<SymmetryOperator>, - // optionally cached operators from [-3, -3, -3] to [3, 3, 3] - _operators_333?: SymmetryOperator[] + /** + * optionally cached operators from [-3, -3, -3] to [3, 3, 3] + * around reference point `ref` in fractional coordinates + */ + _operators_333?: { + ref: Vec3, + operators: SymmetryOperator[] + } } namespace ModelSymmetry { diff --git a/src/mol-model/structure/structure/symmetry.ts b/src/mol-model/structure/structure/symmetry.ts index b472051b1e9cd5473d486353441768dfa06cc34e..6dbd5de05659c58bfaf2d7a46346d2baf49f8249 100644 --- a/src/mol-model/structure/structure/symmetry.ts +++ b/src/mol-model/structure/structure/symmetry.ts @@ -10,7 +10,7 @@ import { EquivalenceClasses } from '../../../mol-data/util'; import { Spacegroup, SpacegroupCell, SymmetryOperator } from '../../../mol-math/geometry'; import { Vec3, Mat4 } from '../../../mol-math/linear-algebra'; import { RuntimeContext, Task } from '../../../mol-task'; -import { ModelSymmetry } from '../model'; +import { ModelSymmetry, Model } from '../model'; import { QueryContext, StructureSelection } from '../query'; import Structure from './structure'; import Unit from './unit'; @@ -95,7 +95,7 @@ namespace StructureSymmetry { } } -function getOperators(symmetry: ModelSymmetry, ijkMin: Vec3, ijkMax: Vec3) { +function getOperators(symmetry: ModelSymmetry, ijkMin: Vec3, ijkMax: Vec3, modelCenter: Vec3) { const { spacegroup, ncsOperators } = symmetry; const ncsCount = (ncsOperators && ncsOperators.length) || 0 const operators: SymmetryOperator[] = []; @@ -107,13 +107,17 @@ function getOperators(symmetry: ModelSymmetry, ijkMin: Vec3, ijkMax: Vec3) { operators[0] = Spacegroup.getSymmetryOperator(spacegroup, 0, 0, 0, 0) } + const { toFractional } = spacegroup.cell + const ref = Vec3.transformMat4(Vec3(), modelCenter, toFractional) + 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++) { // check if we have added identity as the 1st operator. if (!ncsCount && op === 0 && i === 0 && j === 0 && k === 0) continue; - const symOp = Spacegroup.getSymmetryOperator(spacegroup, op, i, j, k); + + const symOp = Spacegroup.getSymmetryOperatorRef(spacegroup, op, i, j, k, ref) if (ncsCount) { for (let u = 0; u < ncsCount; ++u) { const ncsOp = ncsOperators![u] @@ -131,10 +135,15 @@ function getOperators(symmetry: ModelSymmetry, ijkMin: Vec3, ijkMax: Vec3) { return operators; } -function getOperatorsCached333(symmetry: ModelSymmetry) { - if (typeof symmetry._operators_333 !== 'undefined') return symmetry._operators_333; - symmetry._operators_333 = getOperators(symmetry, Vec3.create(-3, -3, -3), Vec3.create(3, 3, 3)); - return symmetry._operators_333; +function getOperatorsCached333(symmetry: ModelSymmetry, ref: Vec3) { + if (symmetry._operators_333 && Vec3.equals(ref, symmetry._operators_333.ref)) { + return symmetry._operators_333.operators; + } + symmetry._operators_333 = { + ref: Vec3.clone(ref), + operators: getOperators(symmetry, Vec3.create(-3, -3, -3), Vec3.create(3, 3, 3), ref) + }; + return symmetry._operators_333.operators; } function assembleOperators(structure: Structure, operators: ReadonlyArray<SymmetryOperator>) { @@ -164,7 +173,8 @@ async function findSymmetryRange(ctx: RuntimeContext, structure: Structure, ijkM const { spacegroup } = models[0].symmetry; if (SpacegroupCell.isZero(spacegroup.cell)) return structure; - const operators = getOperators(models[0].symmetry, ijkMin, ijkMax); + const modelCenter = Model.getCenter(models[0]) + const operators = getOperators(models[0].symmetry, ijkMin, ijkMax, modelCenter); return assembleOperators(structure, operators); } @@ -177,7 +187,8 @@ async function findMatesRadius(ctx: RuntimeContext, structure: Structure, radius if (SpacegroupCell.isZero(spacegroup.cell)) return structure; if (ctx.shouldUpdate) await ctx.update('Initialing...'); - const operators = getOperatorsCached333(symmetry); + const modelCenter = Model.getCenter(models[0]) + const operators = getOperatorsCached333(symmetry, modelCenter); const lookup = structure.lookup3d; const assembler = Structure.Builder();