Browse Source

structure boundary

David Sehnal 7 years ago
parent
commit
0b038559eb

+ 1 - 1
src/mol-math/geometry/lookup3d/common.ts

@@ -32,5 +32,5 @@ export interface Lookup3D<T = number> {
     // The result is mutated with each call to find.
     find(x: number, y: number, z: number, radius: number): Result<T>,
     check(x: number, y: number, z: number, radius: number): boolean,
-    boundary: { box: Box3D, sphere: Sphere3D }
+    readonly boundary: { readonly box: Box3D, readonly sphere: Sphere3D }
 }

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

@@ -65,6 +65,14 @@ namespace Vec3 {
         return out;
     }
 
+    export function ofArray(array: ArrayLike<number>) {
+        const out = zero();
+        out[0] = array[0];
+        out[1] = array[1];
+        out[2] = array[2];
+        return out;
+    }
+
     export function set(out: Vec3, x: number, y: number, z: number): Vec3 {
         out[0] = x;
         out[1] = y;

+ 4 - 0
src/mol-model/structure/structure/structure.ts

@@ -91,6 +91,10 @@ namespace Structure {
         return s.__lookup3d__;
     }
 
+    export function getBoundary(s: Structure) {
+        return getLookup3d(s).boundary;
+    }
+
     // TODO: "lift" atom set operators?
     // TODO: "diff"
 }

+ 70 - 0
src/mol-model/structure/structure/util/boundary.ts

@@ -0,0 +1,70 @@
+/**
+ * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author David Sehnal <david.sehnal@gmail.com>
+ */
+
+import Structure from '../structure'
+import ElementSet from '../element/set'
+import { ElementGroup } from '../../structure'
+import { Box3D, Sphere3D } from 'mol-math/geometry';
+import { Vec3 } from 'mol-math/linear-algebra';
+
+function computeStructureBoundary(s: Structure): { box: Box3D, sphere: Sphere3D } {
+    const min = [Number.MAX_VALUE, Number.MAX_VALUE, Number.MAX_VALUE];
+    const max = [-Number.MAX_VALUE, -Number.MAX_VALUE, -Number.MAX_VALUE];
+
+    const { units, elements } = s;
+
+    let cx = 0, cy = 0, cz = 0;
+    let radiusSq = 0;
+    let size = 0;
+
+    for (let i = 0, _i = ElementSet.unitCount(elements); i < _i; i++) {
+        const group = ElementSet.unitGetByIndex(elements, i);
+        const { x, y, z } = units[ElementSet.unitGetId(elements, i)];
+
+        size += ElementGroup.size(group);
+        for (let j = 0, _j = ElementGroup.size(group); j < _j; j++) {
+            const e = ElementGroup.getAt(group, j);
+            const xx = x(e), yy = y(e), zz = z(e);
+
+            min[0] = Math.min(xx, min[0]);
+            min[1] = Math.min(yy, min[1]);
+            min[2] = Math.min(zz, min[2]);
+            max[0] = Math.max(xx, max[0]);
+            max[1] = Math.max(yy, max[1]);
+            max[2] = Math.max(zz, max[2]);
+
+            cx += xx;
+            cy += yy;
+            cz += zz;
+        }
+    }
+
+    if (size > 0) {
+        cx /= size;
+        cy /= size;
+        cz /= size;
+    }
+
+    for (let i = 0, _i = ElementSet.unitCount(elements); i < _i; i++) {
+        const group = ElementSet.unitGetByIndex(elements, i);
+        const { x, y, z } = units[ElementSet.unitGetId(elements, i)];
+
+        size += ElementGroup.size(group);
+        for (let j = 0, _j = ElementGroup.size(group); j < _j; j++) {
+            const e = ElementGroup.getAt(group, j);
+            const dx = x(e) - cx, dy = y(e) - cy, dz = z(e) - cz;
+            const d = dx * dx + dy * dy + dz * dz;
+            if (d > radiusSq) radiusSq = d;
+        }
+    }
+
+    return {
+        box: { min: Vec3.ofArray(min), max: Vec3.ofArray(max) },
+        sphere: { center: Vec3.create(cx, cy, cz), radius: Math.sqrt(radiusSq) }
+    };
+}
+
+export { computeStructureBoundary }

+ 23 - 8
src/mol-model/structure/structure/util/lookup3d.ts

@@ -10,6 +10,7 @@ import { Lookup3D, GridLookup3D, Result, Box3D, Sphere3D } from 'mol-math/geomet
 import { ElementGroup, ElementSet } from '../../structure';
 import { Vec3 } from 'mol-math/linear-algebra';
 import { OrderedSet } from 'mol-data/int';
+import { computeStructureBoundary } from './boundary';
 
 interface StructureLookup3D extends Lookup3D<Element> {}
 
@@ -25,8 +26,6 @@ namespace StructureLookup3D {
             const closeUnits = this.unitLookup.find(x, y, z, radius);
             if (closeUnits.count === 0) return this.result;
 
-            console.log({ closeUnits });
-
             for (let t = 0, _t = closeUnits.count; t < _t; t++) {
                 const i = closeUnits.indices[t];
                 const unitId = ElementSet.unitGetId(elements, i);
@@ -38,8 +37,6 @@ namespace StructureLookup3D {
                 }
                 const groupLookup = ElementGroup.getLookup3d(unit, group);
                 const groupResult = groupLookup.find(this.pivot[0], this.pivot[1], this.pivot[2], radius);
-                //console.log(groupLookup);
-                //console.log({ groupCount: groupResult.count });
                 for (let j = 0, _j = groupResult.count; j < _j; j++) {
                     Result.add(this.result, Element.create(unitId, groupResult.indices[j]), groupResult.squaredDistances[j]);
                 }
@@ -47,10 +44,29 @@ namespace StructureLookup3D {
 
             return this.result;
         }
+
         check(x: number, y: number, z: number, radius: number): boolean {
-            throw new Error("Method not implemented.");
+            const { units, elements } = this.structure;
+            const closeUnits = this.unitLookup.find(x, y, z, radius);
+            if (closeUnits.count === 0) return false;
+
+            for (let t = 0, _t = closeUnits.count; t < _t; t++) {
+                const i = closeUnits.indices[t];
+                const unitId = ElementSet.unitGetId(elements, i);
+                const group = ElementSet.unitGetByIndex(elements, i);
+                const unit = units[unitId];
+                Vec3.set(this.pivot, x, y, z);
+                if (!unit.operator.isIdentity) {
+                    Vec3.transformMat4(this.pivot, this.pivot, unit.operator.inverse);
+                }
+                const groupLookup = ElementGroup.getLookup3d(unit, group);
+                if (groupLookup.check(this.pivot[0], this.pivot[1], this.pivot[2], radius)) return true;
+            }
+
+            return false;
         }
-        boundary: { box: Box3D; sphere: Sphere3D; } = 0 as any;
+
+        boundary: { box: Box3D; sphere: Sphere3D; };
 
         constructor(private structure: Structure) {
             const { units, elements } = structure;
@@ -75,9 +91,8 @@ namespace StructureLookup3D {
                 radius[i] = s.radius;
             }
 
-            console.log({ xs, ys, zs, radius, unitCount })
-
             this.unitLookup = GridLookup3D({ x: xs, y: ys, z: zs, radius, indices: OrderedSet.ofBounds(0, unitCount) });
+            this.boundary = computeStructureBoundary(structure);
         }
     }
 

+ 3 - 0
src/perf-tests/lookup3d.ts

@@ -56,6 +56,9 @@ export async function test() {
     const sl = Structure.getLookup3d(structures[0]);
     const result1 = sl.find(-30.07, 8.178, -13.897, 10);
     console.log(result1.count);//, result1.indices);
+
+    console.log(Structure.getBoundary(structures[0]));
+    console.log(lookup.boundary);
 }
 
 test();