wip, gaussian surface

src/apps/canvas/structure-view.ts

@@ -68,10 +68,10 @@ export async function StructureView(viewer: Viewer, models: ReadonlyArray<Model>
     const active: { [k: string]: boolean } = {
         cartoon: true,
         point: false,
-        surface: false,
+        surface: true,
         ballAndStick: false,
         carbohydrate: false,
-        spacefill: false,
+        spacefill: true,
         distanceRestraint: false,
         symmetryAxes: false,
         // polymerSphere: false,

src/mol-geo/representation/structure/visual/gaussian-surface-mesh.ts

@@ -13,8 +13,8 @@ import { StructureElementIterator, getElementLoci, markElement } from './util/el
 import { computeMarchingCubes } from '../../../util/marching-cubes/algorithm';
 import { Tensor, Vec3, Mat4 } from 'mol-math/linear-algebra';
 import { Box3D } from 'mol-math/geometry';
-import { smoothstep } from 'mol-math/interpolate';
 import { SizeThemeProps, SizeTheme } from 'mol-view/theme/size';
+import { Color } from 'mol-util/color';
 
 export interface GaussianSurfaceMeshProps {
     sizeTheme: SizeThemeProps
@@ -79,50 +79,58 @@ async function createGaussianSurfaceMesh(ctx: RuntimeContext, unit: Unit, struct
     const data = space.create()
     const field = Tensor.create(space, data)
 
+    const c = Vec3.zero()
+
+    const radius2 = Vec3.scale(Vec3.zero(), delta, 2)
+    const alpha = 1.5 // smoothness
+
+    console.time('surface density')
     for (let i = 0; i < elementCount; i++) {
         l.element = elements[i]
         pos(elements[i], v)
 
-        Vec3.mul(v, Vec3.sub(v, v, min), delta)
+        Vec3.sub(v, v, min)
+        Vec3.mul(c, v, delta)
 
-        const size = sizeTheme.size(l) + probeRadius
-        const radius = size * delta[0]
+        const radius = sizeTheme.size(l) // + probeRadius
 
-        const minX = Math.floor(v[0] - radius)
-        const minY = Math.floor(v[1] - radius)
-        const minZ = Math.floor(v[2] - radius)
-        const maxX = Math.ceil(v[0] + radius)
-        const maxY = Math.ceil(v[1] + radius)
-        const maxZ = Math.ceil(v[2] + radius)
+        const minX = Math.floor(c[0] - radius2[0]) - 1
+        const minY = Math.floor(c[1] - radius2[1]) - 1
+        const minZ = Math.floor(c[2] - radius2[2]) - 1
+        const maxX = Math.ceil(c[0] + radius2[0]) + 1
+        const maxY = Math.ceil(c[1] + radius2[1]) + 1
+        const maxZ = Math.ceil(c[2] + radius2[2]) + 1
 
         for (let x = minX; x < maxX; ++x) {
             for (let y = minY; y < maxY; ++y) {
                 for (let z = minZ; z < maxZ; ++z) {
-                    const dist = Vec3.distance(Vec3.set(p, x, y, z), v)
-                    if (dist <= radius) {
-                        // TODO use actual gaussian
-                        const density = 1.0 - smoothstep(0.0, radius * 1.0, dist)
-                        space.set(data, x, y, z, space.get(data, x, y, z) + density)
-                    }
+                    Vec3.set(p, x, y, z)
+                    Vec3.div(p, p, delta)
+                    const dist = Vec3.distance(p, v)
+                    const density = Math.exp(-alpha * Math.pow((dist) / radius, 2))
+                    space.set(data, x, y, z, space.get(data, x, y, z) + density)
                 }
             }
         }
 
-        if (i % 10000 === 0 && ctx.shouldUpdate) {
+        if (i % 100 === 0 && ctx.shouldUpdate) {
             await ctx.update({ message: 'Gaussian surface', current: i, max: elementCount });
         }
     }
+    console.timeEnd('surface density')
 
     // console.log('data', data)
 
+    console.time('surface mc')
     const surface = await computeMarchingCubes({
-        isoLevel: isoValue,
+        isoLevel: isoValue * Math.exp(-alpha),
         scalarField: field,
         oldSurface: mesh
     }).runAsChild(ctx)
+    console.timeEnd('surface mc')
 
     const t = Mat4.identity()
-    Mat4.fromUniformScaling(t, 1 / delta[0])
+    Mat4.fromScaling(t, Vec3.inverse(Vec3.zero(), delta))
     Mat4.setTranslation(t, expandedBox.min)
 
     Mesh.transformImmediate(surface, t)
@@ -140,10 +148,12 @@ export const DefaultGaussianSurfaceProps = {
 
     flipSided: true,
     // flatShaded: true,
+    alpha: 0.7,
+    colorTheme: { name: 'uniform' as 'uniform', value: Color(0xDDDDDD) },
 
     resolutionFactor: 7,
     probeRadius: 1.4,
-    isoValue: 0.1,
+    isoValue: 1,
 }
 export type GaussianSurfaceProps = typeof DefaultGaussianSurfaceProps