wip, improved bounding sphere calculation

- Sphere3D.Hierarchy (one level)
- use hierarchy when merging spheres
- split sphere into two when data unevenly spread

src/mol-canvas3d/helper/bounding-sphere-helper.ts

@@ -53,19 +53,17 @@ export class BoundingSphereHelper {
             const newObjectData = updateBoundingSphereData(this.scene, r.values.boundingSphere.ref.value, objectData, ColorNames.tomato)
             if (newObjectData) this.objectsData.set(ro, newObjectData)
 
-            if (ro.type === 'mesh' || ro.type === 'lines' || ro.type === 'points') {
-                const instanceData = this.instancesData.get(ro)
-                const newInstanceData = updateBoundingSphereData(this.scene, r.values.invariantBoundingSphere.ref.value, instanceData, ColorNames.skyblue, {
-                    aTransform: ro.values.aTransform,
-                    matrix: ro.values.matrix,
-                    transform: ro.values.transform,
-                    extraTransform: ro.values.extraTransform,
-                    uInstanceCount: ro.values.uInstanceCount,
-                    instanceCount: ro.values.instanceCount,
-                    aInstance: ro.values.aInstance,
-                })
-                if (newInstanceData) this.instancesData.set(ro, newInstanceData)
-            }
+            const instanceData = this.instancesData.get(ro)
+            const newInstanceData = updateBoundingSphereData(this.scene, r.values.invariantBoundingSphere.ref.value, instanceData, ColorNames.skyblue, {
+                aTransform: ro.values.aTransform,
+                matrix: ro.values.matrix,
+                transform: ro.values.transform,
+                extraTransform: ro.values.extraTransform,
+                uInstanceCount: ro.values.uInstanceCount,
+                instanceCount: ro.values.instanceCount,
+                aInstance: ro.values.aInstance,
+            })
+            if (newInstanceData) this.instancesData.set(ro, newInstanceData)
         })
 
         this.objectsData.forEach((objectData, ro) => {
@@ -133,7 +131,11 @@ function createBoundingSphereMesh(boundingSphere: Sphere3D, mesh?: Mesh) {
     const detail = 2
     const vertexCount = sphereVertexCount(detail)
     const builderState = MeshBuilder.createState(vertexCount, vertexCount / 2, mesh)
-    if (boundingSphere.radius) addSphere(builderState, boundingSphere.center, boundingSphere.radius, detail)
+    if (boundingSphere.radius) {
+        for (const b of Sphere3D.getList(boundingSphere)) {
+            addSphere(builderState, b.center, b.radius, detail)
+        }
+    }
     return MeshBuilder.getMesh(builderState)
 }
 

src/mol-gl/renderable/util.ts

@@ -6,7 +6,7 @@
 
 import { Sphere3D } from '../../mol-math/geometry'
 import { Vec3 } from '../../mol-math/linear-algebra'
-import { Epos14, Epos98 } from '../../mol-math/geometry/epos-helper';
+import { EposHelper, HierarchyHelper } from '../../mol-math/geometry/epos-helper';
 
 export function calculateTextureInfo (n: number, itemSize: number) {
     const sqN = Math.sqrt(n)
@@ -78,37 +78,79 @@ export function printImageData(imageData: ImageData, scale = 1, pixelated = fals
 //
 
 const v = Vec3.zero()
-const eposHelper14 = Epos14()
-const eposHelper98 = Epos98()
+const eposHelperCoarse = new EposHelper('14')
+const eposHelperFine = new EposHelper('98')
+const hierarchyHelperCoarse = new HierarchyHelper('14')
+const hierarchyHelperFine = new HierarchyHelper('98')
+
+function getHelper(count: number) {
+    return count > 500_000 ? {
+        eposHelper: eposHelperCoarse,
+        hierarchyHelper: hierarchyHelperCoarse
+    } : {
+        eposHelper: eposHelperFine,
+        hierarchyHelper: hierarchyHelperFine
+    }
+}
 
 export function calculateInvariantBoundingSphere(position: Float32Array, positionCount: number, stepFactor: number): Sphere3D {
     const step = stepFactor * 3
-    eposHelper14.reset()
+    const { eposHelper, hierarchyHelper } = getHelper(positionCount)
+
+    eposHelper.reset()
     for (let i = 0, _i = positionCount * 3; i < _i; i += step) {
         Vec3.fromArray(v, position, i)
-        eposHelper14.includeStep(v)
+        eposHelper.includeStep(v)
     }
-    eposHelper14.finishedIncludeStep()
+    eposHelper.finishedIncludeStep()
     for (let i = 0, _i = positionCount * 3; i < _i; i += step) {
         Vec3.fromArray(v, position, i)
-        eposHelper14.radiusStep(v)
+        eposHelper.radiusStep(v)
+    }
+
+    const hierarchyInput = eposHelper.getHierarchyInput()
+    if (hierarchyInput) {
+        hierarchyHelper.reset(hierarchyInput.sphere, hierarchyInput.normal)
+        for (let i = 0, _i = positionCount * 3; i < _i; i += step) {
+            Vec3.fromArray(v, position, i)
+            hierarchyHelper.includeStep(v)
+        }
+        hierarchyHelper.finishedIncludeStep()
+        for (let i = 0, _i = positionCount * 3; i < _i; i += step) {
+            Vec3.fromArray(v, position, i)
+            hierarchyHelper.radiusStep(v)
+        }
+        return hierarchyHelper.getSphere()
+    } else {
+        return eposHelper.getSphere()
     }
-    return eposHelper14.getSphere()
 }
 
 export function calculateTransformBoundingSphere(invariantBoundingSphere: Sphere3D, transform: Float32Array, transformCount: number): Sphere3D {
-    const { center, radius } = invariantBoundingSphere
-    eposHelper98.reset()
-    for (let i = 0, _i = transformCount; i < _i; ++i) {
-        Vec3.transformMat4Offset(v, center, transform, 0, 0, i * 16)
-        eposHelper98.includeSphereStep(v, radius)
+    const { eposHelper } = getHelper(transformCount)
+    eposHelper.reset()
+
+    const transformedSpheres: Sphere3D[] = []
+    for (const b of Sphere3D.getList(invariantBoundingSphere)) {
+        for (let i = 0, _i = transformCount; i < _i; ++i) {
+            const c = Vec3.transformMat4Offset(Vec3(), b.center, transform, 0, 0, i * 16)
+            transformedSpheres.push(Sphere3D.create(c as Vec3, b.radius))
+        }
     }
-    eposHelper98.finishedIncludeStep()
-    for (let i = 0, _i = transformCount; i < _i; ++i) {
-        Vec3.transformMat4Offset(v, center, transform, 0, 0, i * 16)
-        eposHelper98.radiusSphereStep(v, radius)
+
+    for (const b of transformedSpheres) {
+        eposHelper.includeSphereStep(b.center, b.radius)
+    }
+    eposHelper.finishedIncludeStep()
+    for (const b of transformedSpheres) {
+        eposHelper.radiusSphereStep(b.center, b.radius)
+    }
+
+    const sphere = eposHelper.getSphere()
+    if (transformedSpheres.length > 1) {
+        (sphere as Sphere3D.Hierarchy).hierarchy = transformedSpheres
     }
-    return eposHelper98.getSphere()
+    return sphere
 }
 
 export function calculateBoundingSphere(position: Float32Array, positionCount: number, transform: Float32Array, transformCount: number, padding = 0, stepFactor = 1): { boundingSphere: Sphere3D, invariantBoundingSphere: Sphere3D } {

src/mol-gl/scene.ts

@@ -14,26 +14,32 @@ import { Sphere3D } from '../mol-math/geometry';
 import { CommitQueue } from './commit-queue';
 import { now } from '../mol-util/now';
 import { arraySetRemove } from '../mol-util/array';
-import { Epos98 } from '../mol-math/geometry/epos-helper';
+import { EposHelper } from '../mol-math/geometry/epos-helper';
 
-const eposHelper98 = Epos98()
+const eposHelper = new EposHelper('98')
 
 function calculateBoundingSphere(renderables: Renderable<RenderableValues & BaseValues>[], boundingSphere: Sphere3D): Sphere3D {
-    eposHelper98.reset();
+    eposHelper.reset();
 
     for (let i = 0, il = renderables.length; i < il; ++i) {
         const boundingSphere = renderables[i].values.boundingSphere.ref.value
         if (!boundingSphere.radius) continue;
-        eposHelper98.includeSphereStep(boundingSphere.center, boundingSphere.radius);
+
+        for (const b of Sphere3D.getList(boundingSphere)) {
+            eposHelper.includeSphereStep(b.center, b.radius);
+        }
     }
-    eposHelper98.finishedIncludeStep();
+    eposHelper.finishedIncludeStep();
     for (let i = 0, il = renderables.length; i < il; ++i) {
         const boundingSphere = renderables[i].values.boundingSphere.ref.value
         if (!boundingSphere.radius) continue;
-        eposHelper98.radiusSphereStep(boundingSphere.center, boundingSphere.radius);
+
+        for (const b of Sphere3D.getList(boundingSphere)) {
+            eposHelper.radiusSphereStep(b.center, b.radius);
+        }
     }
 
-    return eposHelper98.getSphere(boundingSphere);
+    return eposHelper.getSphere(boundingSphere);
 }
 
 function renderableSort(a: Renderable<RenderableValues & BaseValues>, b: Renderable<RenderableValues & BaseValues>) {

src/mol-math/geometry/epos-helper.ts

@@ -71,8 +71,28 @@ export class EposHelper {
         this.centroidHelper.radiusSphereStep(center, radius);
     }
 
+    getHierarchyInput() {
+        const sphere = this.centroidHelper.getSphere();
+        const normal = Vec3()
+        const t = sphere.radius * this.hierarchyThresholdFactor
+
+        let maxDist = -Infinity
+        let belowThreshold = false
+
+        for (let i = 0; i < this.extrema.length; i += 2) {
+            const halfDist = Vec3.distance(this.extrema[i], this.extrema[i + 1]) / 2
+            if (halfDist > maxDist) {
+                maxDist = halfDist
+                Vec3.normalize(normal, Vec3.sub(normal, this.extrema[i], this.extrema[i + 1]))
+            }
+            if (halfDist < t) belowThreshold = true
+        }
+
+        return belowThreshold ? { sphere, normal } : false
+    }
+
     getSphere(sphere?: Sphere3D) {
-        return this.centroidHelper.getSphere(sphere);
+        return this.centroidHelper.getSphere(sphere)
     }
 
     reset() {
@@ -85,15 +105,79 @@ export class EposHelper {
         this.centroidHelper.reset()
     }
 
-    constructor(dir: number[][]) {
-        this.dir = dir.map(a => {
-            const v = Vec3.create(a[0], a[1], a[2])
-            return Vec3.normalize(v, v)
-        })
+    constructor(quality: EposQuality, private hierarchyThresholdFactor = 0.66) {
+        this.dir = getEposDir(quality)
         this.reset()
     }
 }
 
+const tmpV = Vec3()
+
+export class HierarchyHelper {
+    private sphere = Sphere3D()
+    private normal = Vec3()
+    private helperA = new EposHelper(this.quality)
+    private helperB = new EposHelper(this.quality)
+
+    private checkSide(p: Vec3) {
+        return Vec3.dot(this.normal, Vec3.sub(tmpV, this.sphere.center, p)) > 0
+    }
+
+    includeStep(p: Vec3) {
+        if (this.checkSide(p)) {
+            this.helperA.includeStep(p)
+        } else {
+            this.helperB.includeStep(p)
+        }
+    }
+
+    finishedIncludeStep() {
+        this.helperA.finishedIncludeStep();
+        this.helperB.finishedIncludeStep();
+    }
+
+    radiusStep(p: Vec3) {
+        if (this.checkSide(p)) {
+            this.helperA.radiusStep(p)
+        } else {
+            this.helperB.radiusStep(p)
+        }
+    }
+
+    getSphere(): Sphere3D.Hierarchy {
+        return {
+            center: this.sphere.center,
+            radius: this.sphere.radius,
+            hierarchy: [this.helperA.getSphere(), this.helperB.getSphere()]
+        }
+    }
+
+    reset(sphere: Sphere3D, normal: Vec3) {
+        Sphere3D.copy(this.sphere, sphere)
+        Vec3.copy(this.normal, normal)
+        this.helperA.reset()
+        this.helperB.reset()
+    }
+
+    constructor(private quality: EposQuality) { }
+}
+
+type EposQuality = '6' | '14' | '26' | '98'
+
+function getEposDir(quality: EposQuality) {
+    let dir: number[][]
+    switch (quality) {
+        case '6': dir = [ ...Type001 ]; break
+        case '14': dir = [ ...Type001, ...Type111 ]; break
+        case '26': dir = [ ...Type001, ...Type111, ...Type011 ]; break
+        case '98': dir = [ ...Type001, ...Type111, ...Type011, ...Type012, ...Type112, ...Type122 ]; break
+    }
+    return dir.map(a => {
+        const v = Vec3.create(a[0], a[1], a[2])
+        return Vec3.normalize(v, v)
+    })
+}
+
 const Type001 = [
     [1, 0, 0], [0, 1, 0], [0, 0, 1]
 ]
@@ -119,24 +203,4 @@ const Type112 = [
 const Type122 = [
     [2, 2, 1], [1, 2, 2], [2, 1, 2], [2, -2, 1], [2, 2, -1], [2, -2, -1],
     [1, -2, 2], [1, 2, -2], [1, -2, -2], [2, -1, 2], [2, 1, -2], [2, -1, -2]
-]
-
-const DirEpos6 = [ ...Type001 ]
-export function Epos6() {
-    return new EposHelper(DirEpos6)
-}
-
-const DirEpos14 = [ ...Type001, ...Type111 ]
-export function Epos14() {
-    return new EposHelper(DirEpos14)
-}
-
-const DirEpos26 = [ ...Type001, ...Type111, ...Type011 ]
-export function Epos26() {
-    return new EposHelper(DirEpos26)
-}
-
-const DirEpos98 = [ ...Type001, ...Type111, ...Type011, ...Type012, ...Type112, ...Type122 ]
-export function Epos98() {
-    return new EposHelper(DirEpos98)
-}
+]

src/mol-math/geometry/primitives/sphere3d.ts

@@ -13,26 +13,42 @@ import { Box3D } from './box3d';
 import { Axes3D } from './axes3d';
 import { BoundaryHelper } from '../boundary-helper';
 
-interface Sphere3D { center: Vec3, radius: number }
+type Sphere3D = Sphere3D.Data | Sphere3D.Hierarchy
 
 function Sphere3D() {
     return Sphere3D.zero();
 }
 
 namespace Sphere3D {
+    export interface Data { center: Vec3, radius: number }
+    export interface Hierarchy extends Data { hierarchy: Sphere3D[] }
+    export function isHierarchy(x: Sphere3D | Hierarchy): x is Hierarchy {
+        return 'hierarchy' in x
+    }
+    export function getList(sphere: Sphere3D) {
+        return Sphere3D.isHierarchy(sphere) ? sphere.hierarchy : [sphere]
+    }
+
     export function create(center: Vec3, radius: number): Sphere3D { return { center, radius }; }
     export function zero(): Sphere3D { return { center: Vec3.zero(), radius: 0 }; }
 
     export function clone(a: Sphere3D): Sphere3D {
-        const out = zero();
-        Vec3.copy(out.center, a.center);
-        out.radius = a.radius
+        const out = create(Vec3.clone(a.center), a.radius)
+        if (isHierarchy(a)) (out as Hierarchy).hierarchy = a.hierarchy
         return out;
     }
 
     export function copy(out: Sphere3D, a: Sphere3D) {
         Vec3.copy(out.center, a.center)
         out.radius = a.radius
+        if (isHierarchy(a)) {
+            if (isHierarchy(out)) {
+                out.hierarchy.length = 0
+                out.hierarchy.push(...a.hierarchy)
+            } else {
+                (out as Hierarchy).hierarchy = a.hierarchy
+            }
+        }
         return out;
     }
 
@@ -127,6 +143,15 @@ namespace Sphere3D {
     export function expand(out: Sphere3D, sphere: Sphere3D, delta: number): Sphere3D {
         Vec3.copy(out.center, sphere.center)
         out.radius = sphere.radius + delta
+        if (isHierarchy(sphere)) {
+            const hierarchy = sphere.hierarchy.map(s => expand(Sphere3D(), s, delta))
+            if (isHierarchy(out)) {
+                out.hierarchy.length = 0
+                out.hierarchy.push(...hierarchy)
+            } else {
+                (out as Hierarchy).hierarchy = hierarchy
+            }
+        }
         return out
     }