From a8b2d8a10539a31ea5feb5a12db97bc194e57bb7 Mon Sep 17 00:00:00 2001
From: David Sehnal <david.sehnal@gmail.com>
Date: Wed, 29 May 2019 11:51:08 +0200
Subject: [PATCH] wip structure superposition
---
.../linear-algebra/3d/minimize-rmsd.ts | 9 ++-
src/mol-model/structure/structure/element.ts | 6 ++
.../structure/structure/util/superposition.ts | 55 +++++++++++++++++++
3 files changed, 68 insertions(+), 2 deletions(-)
create mode 100644 src/mol-model/structure/structure/util/superposition.ts
diff --git a/src/mol-math/linear-algebra/3d/minimize-rmsd.ts b/src/mol-math/linear-algebra/3d/minimize-rmsd.ts
index 84444b667..1296e10f3 100644
--- a/src/mol-math/linear-algebra/3d/minimize-rmsd.ts
+++ b/src/mol-math/linear-algebra/3d/minimize-rmsd.ts
@@ -19,6 +19,11 @@ namespace MinimizeRmsd {
}
export interface Positions { x: ArrayLike<number>, y: ArrayLike<number>, z: ArrayLike<number> }
+ export namespace Positions {
+ export function empty(n: number) {
+ return { x: new Float64Array(n), y: new Float64Array(n), z: new Float64Array(n) };
+ }
+ }
export interface Input {
a: Positions,
@@ -54,10 +59,10 @@ class RmsdTransformState {
this.b = data.b;
if (data.centerA) this.centerA = data.centerA;
- else this.centerA = CentroidHelper.compute(data.a, Vec3.zero());
+ else this.centerA = data.centerA = CentroidHelper.compute(data.a, Vec3.zero());
if (data.centerB) this.centerB = data.centerB;
- else this.centerB = CentroidHelper.compute(data.b, Vec3.zero());
+ else this.centerB = data.centerB = CentroidHelper.compute(data.b, Vec3.zero());
this.result = into;
}
diff --git a/src/mol-model/structure/structure/element.ts b/src/mol-model/structure/structure/element.ts
index 3e0ec4ec8..da1e5a307 100644
--- a/src/mol-model/structure/structure/element.ts
+++ b/src/mol-model/structure/structure/element.ts
@@ -118,6 +118,12 @@ namespace StructureElement {
}
export namespace Loci {
+ export function size(loci: Loci) {
+ let s = 0;
+ for (const u of loci.elements) s += OrderedSet.size(u.indices);
+ return s;
+ }
+
export function all(structure: Structure): Loci {
return Loci(structure, structure.units.map(unit => ({
unit,
diff --git a/src/mol-model/structure/structure/util/superposition.ts b/src/mol-model/structure/structure/util/superposition.ts
new file mode 100644
index 000000000..8c41ba996
--- /dev/null
+++ b/src/mol-model/structure/structure/util/superposition.ts
@@ -0,0 +1,55 @@
+/**
+ * Copyright (c) 2019 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author David Sehnal <david.sehnal@gmail.com>
+ */
+
+import { MinimizeRmsd } from 'mol-math/linear-algebra';
+import StructureElement from '../element';
+import { OrderedSet } from 'mol-data/int';
+
+export function superposeStructures(xs: StructureElement.Loci[]): MinimizeRmsd.Result[] {
+ const ret: MinimizeRmsd.Result[] = [];
+ if (xs.length <= 0) return ret;
+
+ const n = getMinSize(xs);
+ const input: MinimizeRmsd.Input = { a: getPositionTable(xs[0], n), b: getPositionTable(xs[1], n) };
+ ret[0] = MinimizeRmsd.compute(input);
+ for (let i = 2; i < xs.length; i++) {
+ input.b = getPositionTable(xs[i], n);
+ input.centerB = void 0;
+ ret.push(MinimizeRmsd.compute(input));
+ }
+
+ return ret;
+}
+
+function getPositionTable(xs: StructureElement.Loci, n: number): MinimizeRmsd.Positions {
+ const ret = MinimizeRmsd.Positions.empty(n);
+ let o = 0;
+ for (const u of xs.elements) {
+ const { unit, indices } = u;
+ const { elements } = unit;
+ const { x, y, z } = unit.conformation;
+ for (let i = 0, _i = OrderedSet.size(indices); i < _i; i++) {
+ const e = elements[OrderedSet.getAt(indices, i)];
+ ret.x[o] = x(e);
+ ret.y[o] = y(e);
+ ret.z[o] = z(e);
+ o++;
+ if (o >= n) break;
+ }
+ if (o >= n) break;
+ }
+ return ret;
+}
+
+function getMinSize(xs: StructureElement.Loci[]) {
+ if (xs.length === 0) return 0;
+ let s = StructureElement.Loci.size(xs[0]);
+ for (let i = 1; i < xs.length; i++) {
+ const t = StructureElement.Loci.size(xs[i]);
+ if (t < s) s = t;
+ }
+ return s;
+}
\ No newline at end of file
--
GitLab