calculate gausian density as unit property

src/apps/canvas/component/structure-representation.tsx

@@ -36,6 +36,8 @@ export interface StructureRepresentationComponentState {
     pointSizeAttenuation?: boolean
     pointFilledCircle?: boolean
     pointEdgeBleach?: number
+
+    visuals?: { [k: string]: boolean }
 }
 
 export class StructureRepresentationComponent extends React.Component<StructureRepresentationComponentProps, StructureRepresentationComponentState> {
@@ -58,6 +60,8 @@ export class StructureRepresentationComponent extends React.Component<StructureR
             pointSizeAttenuation: (repr.props as any).pointSizeAttenuation,
             pointFilledCircle: (repr.props as any).pointFilledCircle,
             pointEdgeBleach: (repr.props as any).pointEdgeBleach,
+
+            visuals: (repr.props as any).visuals,
         }
     }
 
@@ -84,12 +88,13 @@ export class StructureRepresentationComponent extends React.Component<StructureR
         if (state.pointFilledCircle !== undefined) (props as any).pointFilledCircle = state.pointFilledCircle
         if (state.pointEdgeBleach !== undefined) (props as any).pointEdgeBleach = state.pointEdgeBleach
 
+        if (state.visuals !== undefined) (props as any).visuals = state.visuals
+
         await this.props.app.runTask(repr.createOrUpdate(props).run(
             progress => console.log(Progress.format(progress))
         ), 'Create/update representation')
         this.props.viewer.add(repr)
         this.props.viewer.draw(true)
-        console.log(this.props.viewer.stats)
 
         this.setState(this.stateFromRepresentation(repr))
     }
@@ -109,6 +114,18 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         {visible ? 'Hide' : 'Show'}
                     </button>
                 </div>
+                { this.state.visuals !== undefined ? <div>
+                    <span>Visuals: </span>
+                    { Object.keys(this.state.visuals).map(k => {
+                        return <span key={k}>{k} <input
+                            type='checkbox'
+                            checked={this.state.visuals[k]}
+                            onChange={e => {
+                                this.update({ visuals: { ...this.state.visuals, [k]: !!e.target.checked } })
+                            }}
+                        ></input> </span>
+                    }) }
+                </div> : '' }
                 <div>
                     <span>Depth Mask </span>
                     <button onClick={(e) => this.update({ depthMask: !depthMask }) }>
@@ -123,7 +140,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                 </div> : '' }
                 <div>
                     <span>Quality </span>
-                    <select value={quality} onChange={(e) => this.update({ quality: e.target.value as VisualQuality }) }>
+                    <select value={quality} onChange={e => this.update({ quality: e.target.value as VisualQuality }) }>
                         {VisualQualityNames.map(name => <option key={name} value={name}>{name}</option>)}
                     </select>
                 </div>
@@ -134,7 +151,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='0'
                         max='1'
                         step='0.05'
-                        onInput={(e) => this.update({ alpha: parseFloat(e.currentTarget.value) })}
+                        onInput={e => this.update({ alpha: parseFloat(e.currentTarget.value) })}
                     >
                     </input>
                 </div>
@@ -145,7 +162,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='4'
                         max='9'
                         step='1'
-                        onInput={(e) => this.update({ resolutionFactor: parseInt(e.currentTarget.value) })}
+                        onChange={(e) => this.update({ resolutionFactor: parseInt(e.currentTarget.value) })}
                     >
                     </input>
                 </div> : '' }
@@ -156,7 +173,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='1'
                         max='3'
                         step='0.1'
-                        onInput={(e) => this.update({ smoothness: parseFloat(e.currentTarget.value) })}
+                        onChange={e => this.update({ smoothness: parseFloat(e.currentTarget.value) })}
                     >
                     </input>
                 </div> : '' }
@@ -167,19 +184,19 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='0'
                         max='4'
                         step='0.1'
-                        onInput={(e) => this.update({ radiusOffset: parseFloat(e.currentTarget.value) })}
+                        onChange={e => this.update({ radiusOffset: parseFloat(e.currentTarget.value) })}
                     >
                     </input>
                 </div> : '' }
                 { this.state.pointSizeAttenuation !== undefined ? <div>
                     <span>Size Attenuation </span>
-                    <button onClick={(e) => this.update({ pointSizeAttenuation: !this.state.pointSizeAttenuation }) }>
+                    <button onClick={e => this.update({ pointSizeAttenuation: !this.state.pointSizeAttenuation }) }>
                         {this.state.pointSizeAttenuation ? 'Deactivate' : 'Activate'}
                     </button>
                 </div> : '' }
                 { this.state.pointFilledCircle !== undefined ? <div>
                     <span>Filled Circle </span>
-                    <button onClick={(e) => this.update({ pointFilledCircle: !this.state.pointFilledCircle }) }>
+                    <button onClick={e => this.update({ pointFilledCircle: !this.state.pointFilledCircle }) }>
                         {this.state.pointFilledCircle ? 'Deactivate' : 'Activate'}
                     </button>
                 </div> : '' }
@@ -190,7 +207,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='0'
                         max='1'
                         step='0.05'
-                        onInput={(e) => this.update({ pointEdgeBleach: parseFloat(e.currentTarget.value) })}
+                        onInput={e => this.update({ pointEdgeBleach: parseFloat(e.currentTarget.value) })}
                     >
                     </input>
                 </div> : '' }
@@ -201,7 +218,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                         min='0'
                         max='10'
                         step='0.1'
-                        onInput={(e) => this.update({
+                        onInput={e => this.update({
                             sizeTheme: { name: 'uniform', value: parseFloat(e.currentTarget.value) }
                         })}
                     >
@@ -209,7 +226,7 @@ export class StructureRepresentationComponent extends React.Component<StructureR
                 </div> : '' }
                 <div>
                     <span>Color Theme </span>
-                    <select value={colorTheme.name} onChange={(e) => this.update({ colorTheme: { name: e.target.value as ColorThemeName } }) }>
+                    <select value={colorTheme.name} onChange={e => this.update({ colorTheme: { name: e.target.value as ColorThemeName } }) }>
                         {ColorThemeNames.map(name => <option key={name} value={name}>{name}</option>)}
                     </select>
                     {ct.description ? <div><i>{ct.description}</i></div> : ''}

src/mol-geo/representation/structure/representation/surface.ts

@@ -12,13 +12,13 @@ import { MarkerAction } from '../../../geometry/marker-data';
 import { Loci } from 'mol-model/loci';
 import { PickingId } from '../../../geometry/picking';
 import { Task } from 'mol-task';
-import { GaussianDensityPointVisual } from '../visual/gaussian-density-point';
 import { DefaultGaussianWireframeProps, GaussianWireframeVisual } from '../visual/gaussian-surface-wireframe';
 
 export const DefaultSurfaceProps = {
     ...DefaultGaussianSurfaceProps,
-    // ...DefaultGaussianDensityPointProps,
     ...DefaultGaussianWireframeProps,
+
+    visuals: { surface: true, wireframe: false }
 }
 export type SurfaceProps = typeof DefaultSurfaceProps
 
@@ -26,23 +26,26 @@ export type SurfaceRepresentation = StructureRepresentation<SurfaceProps>
 
 export function SurfaceRepresentation(): SurfaceRepresentation {
     let currentProps: SurfaceProps
+    let currentStructure: Structure
     const gaussianSurfaceRepr = UnitsRepresentation('Gaussian surface', GaussianSurfaceVisual)
-    const gaussianPointRepr = UnitsRepresentation('Gaussian point grid', GaussianDensityPointVisual)
     const gaussianWireframeRepr = UnitsRepresentation('Gaussian wireframe', GaussianWireframeVisual)
     return {
         label: 'Surface',
         get renderObjects() {
-            return [ ...gaussianSurfaceRepr.renderObjects, ...gaussianPointRepr.renderObjects, ...gaussianWireframeRepr.renderObjects ]
+            const renderObjects = []
+            if (currentProps.visuals.surface) renderObjects.push(...gaussianSurfaceRepr.renderObjects)
+            if (currentProps.visuals.wireframe) renderObjects.push(...gaussianWireframeRepr.renderObjects)
+            return renderObjects
         },
         get props() {
-            return { ...gaussianSurfaceRepr.props, ...gaussianPointRepr.props, ...gaussianWireframeRepr.props }
+            return { ...gaussianSurfaceRepr.props, ...gaussianWireframeRepr.props, visuals: currentProps.visuals }
         },
         createOrUpdate: (props: Partial<SurfaceProps> = {}, structure?: Structure) => {
             currentProps = Object.assign({}, DefaultSurfaceProps, currentProps, props)
+            if (structure) currentStructure = structure
             return Task.create('Creating SurfaceRepresentation', async ctx => {
-                await gaussianSurfaceRepr.createOrUpdate(currentProps, structure).runInContext(ctx)
-                // await gaussianPointRepr.createOrUpdate(currentProps, structure).runInContext(ctx)
-                await gaussianWireframeRepr.createOrUpdate(currentProps, structure).runInContext(ctx)
+                if (currentProps.visuals.surface) await gaussianSurfaceRepr.createOrUpdate(currentProps, currentStructure).runInContext(ctx)
+                if (currentProps.visuals.wireframe) await gaussianWireframeRepr.createOrUpdate(currentProps, currentStructure).runInContext(ctx)
             })
         },
         getLoci: (pickingId: PickingId) => {
@@ -53,7 +56,6 @@ export function SurfaceRepresentation(): SurfaceRepresentation {
         },
         destroy() {
             gaussianSurfaceRepr.destroy()
-            gaussianPointRepr.destroy()
             gaussianWireframeRepr.destroy()
         }
     }

src/mol-geo/representation/structure/visual/gaussian-density-point.ts

@@ -11,10 +11,10 @@ import { StructureElementIterator } from './util/element';
 import { EmptyLoci } from 'mol-model/loci';
 import { Vec3 } from 'mol-math/linear-algebra';
 import { UnitsPointsVisual, DefaultUnitsPointsProps } from '../units-visual';
-import { computeGaussianDensity, DefaultGaussianDensityProps } from './util/gaussian';
 import { Points } from '../../../geometry/points/points';
 import { PointsBuilder } from '../../../geometry/points/points-builder';
 import { SizeThemeProps } from 'mol-view/theme/size';
+import { DefaultGaussianDensityProps, GaussianDensityProps } from 'mol-model/structure/structure/unit/gaussian-density';
 
 export const DefaultGaussianDensityPointProps = {
     ...DefaultUnitsPointsProps,
@@ -25,8 +25,8 @@ export const DefaultGaussianDensityPointProps = {
 }
 export type GaussianDensityPointProps = typeof DefaultGaussianDensityPointProps
 
-export async function createGaussianDensityPoint(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianDensityPointProps, points?: Points) {
-    const { transform, field: { space, data } } = await computeGaussianDensity(unit, structure, props).runAsChild(ctx)
+export async function createGaussianDensityPoint(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianDensityProps, points?: Points) {
+    const { transform, field: { space, data } } = await unit.computeGaussianDensity(props, ctx)
 
     const { dimensions, get } = space
     const [ xn, yn, zn ] = dimensions

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

@@ -11,11 +11,11 @@ import { Mesh } from '../../../geometry/mesh/mesh';
 import { UnitsMeshVisual, DefaultUnitsMeshProps } from '../units-visual';
 import { StructureElementIterator, getElementLoci, markElement } from './util/element';
 import { computeMarchingCubes } from '../../../util/marching-cubes/algorithm';
-import { computeGaussianDensity, DefaultGaussianDensityProps } from './util/gaussian';
+import { DefaultGaussianDensityProps, GaussianDensityProps } from 'mol-model/structure/structure/unit/gaussian-density';
 
-async function createGaussianSurfaceMesh(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianSurfaceProps, mesh?: Mesh): Promise<Mesh> {
+async function createGaussianSurfaceMesh(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianDensityProps, mesh?: Mesh): Promise<Mesh> {
     const { smoothness } = props
-    const { transform, field, idField } = await computeGaussianDensity(unit, structure, props).runAsChild(ctx)
+    const { transform, field, idField } = await unit.computeGaussianDensity(props, ctx)
 
     const surface = await computeMarchingCubes({
         isoLevel: Math.exp(-smoothness),

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

@@ -11,13 +11,13 @@ import { Mesh } from '../../../geometry/mesh/mesh';
 import { UnitsLinesVisual, DefaultUnitsLinesProps } from '../units-visual';
 import { StructureElementIterator, getElementLoci, markElement } from './util/element';
 import { computeMarchingCubes } from '../../../util/marching-cubes/algorithm';
-import { computeGaussianDensity, DefaultGaussianDensityProps, GaussianDensityProps } from './util/gaussian';
 import { Lines } from '../../../geometry/lines/lines';
 import { LinesBuilder } from '../../../geometry/lines/lines-builder';
+import { GaussianDensityProps, DefaultGaussianDensityProps } from 'mol-model/structure/structure/unit/gaussian-density';
 
 async function createGaussianWireframe(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianDensityProps, lines?: Lines): Promise<Lines> {
     const { smoothness } = props
-    const { transform, field, idField } = await computeGaussianDensity(unit, structure, props).runAsChild(ctx)
+    const { transform, field, idField } = await unit.computeGaussianDensity(props, ctx)
 
     const surface = await computeMarchingCubes({
         isoLevel: Math.exp(-smoothness),

src/mol-geo/representation/structure/visual/util/gaussian.ts

@@ -3,128 +3,3 @@
  *
  * @author Alexander Rose <alexander.rose@weirdbyte.de>
  */
-
-import { Unit, Structure, StructureElement } from 'mol-model/structure';
-import { RuntimeContext, Task } from 'mol-task'
-import { Tensor, Vec3, Mat4 } from 'mol-math/linear-algebra';
-import { Box3D } from 'mol-math/geometry';
-import { SizeTheme } from 'mol-view/theme/size';
-
-export const DefaultGaussianDensityProps = {
-    resolutionFactor: 6,
-    radiusOffset: 0,
-    smoothness: 1.5,
-}
-export type GaussianDensityProps = typeof DefaultGaussianDensityProps
-
-function getDelta(box: Box3D, resolutionFactor: number) {
-    const extent = Vec3.sub(Vec3.zero(), box.max, box.min)
-    const n = Math.pow(Math.pow(2, resolutionFactor), 3)
-    const f = (extent[0] * extent[1] * extent[2]) / n
-    const s = Math.pow(f, 1 / 3)
-    const size = Vec3.zero()
-    Vec3.ceil(size, Vec3.scale(size, extent, s))
-    const delta = Vec3.div(Vec3.zero(), extent, size)
-    return delta
-}
-
-type Density = { transform: Mat4, field: Tensor, idField: Tensor }
-
-export function computeGaussianDensity(unit: Unit, structure: Structure, props: GaussianDensityProps) {
-    return Task.create('Gaussian Density', async ctx => {
-        return await GaussianDensity(ctx, unit, structure, props);
-    });
-}
-
-export async function GaussianDensity(ctx: RuntimeContext, unit: Unit, structure: Structure, props: GaussianDensityProps): Promise<Density> {
-    const { resolutionFactor, radiusOffset, smoothness } = props
-
-    const { elements } = unit;
-    const elementCount = elements.length;
-    const sizeTheme = SizeTheme({ name: 'physical' })
-
-    const v = Vec3.zero()
-    const p = Vec3.zero()
-    const pos = unit.conformation.invariantPosition
-    const l = StructureElement.create(unit)
-
-    const pad = (radiusOffset + 3) * 3 // TODO calculate max radius
-    const box = unit.lookup3d.boundary.box
-    const expandedBox = Box3D.expand(Box3D.empty(), box, Vec3.create(pad, pad, pad));
-    const extent = Vec3.sub(Vec3.zero(), expandedBox.max, expandedBox.min)
-    const min = expandedBox.min
-
-    const delta = getDelta(Box3D.expand(Box3D.empty(), structure.boundary.box, Vec3.create(pad, pad, pad)), resolutionFactor)
-    const dim = Vec3.zero()
-    Vec3.ceil(dim, Vec3.mul(dim, extent, delta))
-
-    const space = Tensor.Space(dim, [0, 1, 2], Float32Array)
-    const data = space.create()
-    const field = Tensor.create(space, data)
-
-    const idData = space.create()
-    const idField = Tensor.create(space, idData)
-
-    const densData = space.create()
-
-    const c = Vec3.zero()
-
-    const alpha = smoothness
-
-    const _r2 = (radiusOffset + 1.4 * 2)
-    const _radius2 = Vec3.create(_r2, _r2, _r2)
-    Vec3.mul(_radius2, _radius2, delta)
-    const updateChunk = Math.ceil(10000 / (_radius2[0] * _radius2[1] * _radius2[2]))
-
-    const beg = Vec3.zero()
-    const end = Vec3.zero()
-
-    const gridPad = 1 / Math.max(...delta)
-
-    for (let i = 0; i < elementCount; i++) {
-        l.element = elements[i]
-        pos(elements[i], v)
-
-        Vec3.sub(v, v, min)
-        Vec3.mul(c, v, delta)
-
-        const radius = sizeTheme.size(l) + radiusOffset
-        const rSq = radius * radius
-
-        const r2 = radiusOffset + radius * 2 + gridPad
-        const radius2 = Vec3.create(r2, r2, r2)
-        Vec3.mul(radius2, radius2, delta)
-        const r2sq = r2 * r2
-
-        const [ begX, begY, begZ ] = Vec3.floor(beg, Vec3.sub(beg, c, radius2))
-        const [ endX, endY, endZ ] = Vec3.ceil(end, Vec3.add(end, c, radius2))
-
-        for (let x = begX; x < endX; ++x) {
-            for (let y = begY; y < endY; ++y) {
-                for (let z = begZ; z < endZ; ++z) {
-                    Vec3.set(p, x, y, z)
-                    Vec3.div(p, p, delta)
-                    const distSq = Vec3.squaredDistance(p, v)
-                    if (distSq <= r2sq) {
-                        const dens = Math.exp(-alpha * (distSq / rSq))
-                        space.add(data, x, y, z, dens)
-                        if (dens > space.get(densData, x, y, z)) {
-                            space.set(densData, x, y, z, dens)
-                            space.set(idData, x, y, z, i)
-                        }
-                    }
-                }
-            }
-        }
-
-        if (i % updateChunk === 0 && ctx.shouldUpdate) {
-            await ctx.update({ message: 'filling density grid', current: i, max: elementCount });
-        }
-    }
-
-    const transform = Mat4.identity()
-    Mat4.fromScaling(transform, Vec3.inverse(Vec3.zero(), delta))
-    Mat4.setTranslation(transform, expandedBox.min)
-
-    return { field, idField, transform }
-}

src/mol-math/geometry/common.ts

@@ -5,6 +5,7 @@
  */
 
 import { OrderedSet } from 'mol-data/int'
+import { Mat4, Tensor } from '../linear-algebra';
 
 export interface PositionData {
     x: ArrayLike<number>,
@@ -15,3 +16,5 @@ export interface PositionData {
     // optional element radius
     radius?: ArrayLike<number>
 }
+
+export type DensityData = { transform: Mat4, field: Tensor, idField: Tensor }

src/mol-math/geometry/gaussian-density.ts

@@ -0,0 +1,124 @@
+/**
+ * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author Alexander Rose <alexander.rose@weirdbyte.de>
+ */
+
+import { Box3D } from '../geometry';
+import { Vec3, Mat4, Tensor } from '../linear-algebra';
+import { RuntimeContext, Task } from 'mol-task';
+import { PositionData, DensityData } from './common';
+import { OrderedSet } from 'mol-data/int';
+
+export const DefaultGaussianDensityProps = {
+    resolutionFactor: 6,
+    radiusOffset: 0,
+    smoothness: 1.5,
+    box: Box3D.empty() // TODO remove
+}
+export type GaussianDensityProps = typeof DefaultGaussianDensityProps
+
+function getDelta(box: Box3D, resolutionFactor: number) {
+    const extent = Vec3.sub(Vec3.zero(), box.max, box.min)
+    const n = Math.pow(Math.pow(2, resolutionFactor), 3)
+    const f = (extent[0] * extent[1] * extent[2]) / n
+    const s = Math.pow(f, 1 / 3)
+    const size = Vec3.zero()
+    Vec3.ceil(size, Vec3.scale(size, extent, s))
+    const delta = Vec3.div(Vec3.zero(), extent, size)
+    return delta
+}
+
+export function computeGaussianDensity(position: PositionData, box: Box3D, radius: (index: number) => number,  props: GaussianDensityProps) {
+    return Task.create('Gaussian Density', async ctx => await GaussianDensity(ctx, position, box, radius, props));
+}
+
+export async function GaussianDensity(ctx: RuntimeContext, position: PositionData, box: Box3D, radius: (index: number) => number,  props: GaussianDensityProps): Promise<DensityData> {
+    const { resolutionFactor, radiusOffset, smoothness } = props
+
+    const { indices, x, y, z } = position
+    const n = OrderedSet.size(indices)
+
+    const v = Vec3.zero()
+    const p = Vec3.zero()
+
+    const pad = (radiusOffset + 3) * 3 // TODO calculate max radius
+    const expandedBox = Box3D.expand(Box3D.empty(), box, Vec3.create(pad, pad, pad));
+    const extent = Vec3.sub(Vec3.zero(), expandedBox.max, expandedBox.min)
+    const min = expandedBox.min
+
+    const delta = getDelta(Box3D.expand(Box3D.empty(), props.box, Vec3.create(pad, pad, pad)), resolutionFactor)
+    const dim = Vec3.zero()
+    Vec3.ceil(dim, Vec3.mul(dim, extent, delta))
+
+    const space = Tensor.Space(dim, [0, 1, 2], Float32Array)
+    const data = space.create()
+    const field = Tensor.create(space, data)
+
+    const idData = space.create()
+    const idField = Tensor.create(space, idData)
+
+    const densData = space.create()
+
+    const c = Vec3.zero()
+
+    const alpha = smoothness
+
+    const _r2 = (radiusOffset + 1.4 * 2)
+    const _radius2 = Vec3.create(_r2, _r2, _r2)
+    Vec3.mul(_radius2, _radius2, delta)
+    const updateChunk = Math.ceil(10000 / (_radius2[0] * _radius2[1] * _radius2[2]))
+
+    const beg = Vec3.zero()
+    const end = Vec3.zero()
+
+    const gridPad = 1 / Math.max(...delta)
+
+    for (let i = 0; i < n; ++i) {
+        const j = OrderedSet.getAt(indices, i);
+
+        Vec3.set(v, x[j], y[j], z[j])
+
+        Vec3.sub(v, v, min)
+        Vec3.mul(c, v, delta)
+
+        const rad = radius(j) + radiusOffset
+        const rSq = rad * rad
+
+        const r2 = radiusOffset + rad * 2 + gridPad
+        const rad2 = Vec3.create(r2, r2, r2)
+        Vec3.mul(rad2, rad2, delta)
+        const r2sq = r2 * r2
+
+        const [ begX, begY, begZ ] = Vec3.floor(beg, Vec3.sub(beg, c, rad2))
+        const [ endX, endY, endZ ] = Vec3.ceil(end, Vec3.add(end, c, rad2))
+
+        for (let xi = begX; xi < endX; ++xi) {
+            for (let yi = begY; yi < endY; ++yi) {
+                for (let zi = begZ; zi < endZ; ++zi) {
+                    Vec3.set(p, xi, yi, zi)
+                    Vec3.div(p, p, delta)
+                    const distSq = Vec3.squaredDistance(p, v)
+                    if (distSq <= r2sq) {
+                        const dens = Math.exp(-alpha * (distSq / rSq))
+                        space.add(data, xi, yi, zi, dens)
+                        if (dens > space.get(densData, xi, yi, zi)) {
+                            space.set(densData, xi, yi, zi, dens)
+                            space.set(idData, xi, yi, zi, i)
+                        }
+                    }
+                }
+            }
+        }
+
+        if (i % updateChunk === 0 && ctx.shouldUpdate) {
+            await ctx.update({ message: 'filling density grid', current: i, max: n });
+        }
+    }
+
+    const transform = Mat4.identity()
+    Mat4.fromScaling(transform, Vec3.inverse(Vec3.zero(), delta))
+    Mat4.setTranslation(transform, expandedBox.min)
+
+    return { field, idField, transform }
+}

src/mol-model/structure/structure/unit.ts

@@ -7,7 +7,7 @@
 
 import { SymmetryOperator } from 'mol-math/geometry/symmetry-operator'
 import { Model } from '../model'
-import { GridLookup3D, Lookup3D } from 'mol-math/geometry'
+import { GridLookup3D, Lookup3D, DensityData } from 'mol-math/geometry'
 import { IntraUnitLinks, computeIntraUnitBonds } from './unit/links'
 import { CoarseElements, CoarseSphereConformation, CoarseGaussianConformation } from '../model/properties/coarse';
 import { ValueRef } from 'mol-util';
@@ -18,6 +18,8 @@ import { IntMap, SortedArray } from 'mol-data/int';
 import { hash2 } from 'mol-data/util';
 import { getAtomicPolymerElements, getCoarsePolymerElements, getAtomicGapElements, getCoarseGapElements } from './util/polymer';
 import { getNucleotideElements } from './util/nucleotide';
+import { GaussianDensityProps, computeUnitGaussianDensityCached } from './unit/gaussian-density';
+import { RuntimeContext } from 'mol-task';
 
 // A building block of a structure that corresponds to an atomic or a coarse grained representation
 // 'conveniently grouped together'.
@@ -180,6 +182,10 @@ namespace Unit {
             return this.model.atomicHierarchy.residueAtomSegments.index[this.elements[elementIndex]];
         }
 
+        async computeGaussianDensity(props: GaussianDensityProps, ctx?: RuntimeContext) {
+            return computeUnitGaussianDensityCached(this, props, this.props.gaussianDensities, ctx);
+        }
+
         constructor(id: number, invariantId: number, model: Model, elements: StructureElement.Set, conformation: SymmetryOperator.ArrayMapping, props: AtomicProperties) {
             this.id = id;
             this.invariantId = invariantId;
@@ -200,6 +206,7 @@ namespace Unit {
         polymerElements: ValueRef<SortedArray<ElementIndex> | undefined>
         gapElements: ValueRef<SortedArray<ElementIndex> | undefined>
         nucleotideElements: ValueRef<SortedArray<ElementIndex> | undefined>
+        gaussianDensities: Map<string, DensityData>
     }
 
     function AtomicProperties(): AtomicProperties {
@@ -210,6 +217,7 @@ namespace Unit {
             polymerElements: ValueRef.create(void 0),
             gapElements: ValueRef.create(void 0),
             nucleotideElements: ValueRef.create(void 0),
+            gaussianDensities: new Map()
         };
     }
 
@@ -234,7 +242,7 @@ namespace Unit {
 
         applyOperator(id: number, operator: SymmetryOperator, dontCompose = false): Unit {
             const op = dontCompose ? operator : SymmetryOperator.compose(this.conformation.operator, operator);
-            const ret = createCoarse(id, this.invariantId, this.model, this.kind, this.elements, SymmetryOperator.createMapping(op, this.getCoarseElements(), this.conformation.r), this.props);
+            const ret = createCoarse(id, this.invariantId, this.model, this.kind, this.elements, SymmetryOperator.createMapping(op, this.getCoarseConformation(), this.conformation.r), this.props);
             // (ret as Coarse<K, C>)._lookup3d = this._lookup3d;
             return ret;
         }
@@ -242,27 +250,31 @@ namespace Unit {
         get lookup3d() {
             if (this.props.lookup3d.ref) return this.props.lookup3d.ref;
             // TODO: support sphere radius?
-            const { x, y, z } = this.getCoarseElements();
+            const { x, y, z } = this.getCoarseConformation();
             this.props.lookup3d.ref = GridLookup3D({ x, y, z, indices: this.elements });
             return this.props.lookup3d.ref;
         }
 
         get polymerElements() {
             if (this.props.polymerElements.ref) return this.props.polymerElements.ref;
-            this.props.polymerElements.ref = getCoarsePolymerElements(this as Unit.Spheres | Unit.Gaussians); // TODO
+            this.props.polymerElements.ref = getCoarsePolymerElements(this as Unit.Spheres | Unit.Gaussians); // TODO get rid of casting
             return this.props.polymerElements.ref;
         }
 
         get gapElements() {
             if (this.props.gapElements.ref) return this.props.gapElements.ref;
-            this.props.gapElements.ref = getCoarseGapElements(this as Unit.Spheres | Unit.Gaussians); // TODO
+            this.props.gapElements.ref = getCoarseGapElements(this as Unit.Spheres | Unit.Gaussians); // TODO get rid of casting
             return this.props.gapElements.ref;
         }
 
-        private getCoarseElements() {
+        private getCoarseConformation() {
             return this.kind === Kind.Spheres ? this.model.coarseConformation.spheres : this.model.coarseConformation.gaussians;
         }
 
+        async computeGaussianDensity(props: GaussianDensityProps, ctx?: RuntimeContext): Promise<DensityData> {
+            return computeUnitGaussianDensityCached(this as Unit.Spheres | Unit.Gaussians, props, this.props.gaussianDensities, ctx); // TODO get rid of casting
+        }
+
         constructor(id: number, invariantId: number, model: Model, kind: K, elements: StructureElement.Set, conformation: SymmetryOperator.ArrayMapping, props: CoarseProperties) {
             this.kind = kind;
             this.id = id;
@@ -278,6 +290,7 @@ namespace Unit {
 
     interface CoarseProperties {
         lookup3d: ValueRef<Lookup3D | undefined>,
+        gaussianDensities: Map<string, DensityData>
         polymerElements: ValueRef<SortedArray<ElementIndex> | undefined>
         gapElements: ValueRef<SortedArray<ElementIndex> | undefined>
     }
@@ -285,6 +298,7 @@ namespace Unit {
     function CoarseProperties(): CoarseProperties {
         return {
             lookup3d: ValueRef.create(void 0),
+            gaussianDensities: new Map(),
             polymerElements: ValueRef.create(void 0),
             gapElements: ValueRef.create(void 0),
         };

src/mol-model/structure/structure/unit/gaussian-density.ts

@@ -0,0 +1,57 @@
+/**
+ * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info.
+ *
+ * @author Alexander Rose <alexander.rose@weirdbyte.de>
+ */
+
+import { Unit, StructureElement, ElementIndex } from 'mol-model/structure';
+import { SizeTheme } from 'mol-view/theme/size';
+import { GaussianDensity } from 'mol-math/geometry/gaussian-density';
+import { Task, RuntimeContext } from 'mol-task';
+import { DensityData } from 'mol-math/geometry';
+
+export const DefaultGaussianDensityProps = {
+    resolutionFactor: 6,
+    radiusOffset: 0,
+    smoothness: 1.5,
+}
+export type GaussianDensityProps = typeof DefaultGaussianDensityProps
+
+function getConformation(unit: Unit) {
+    switch (unit.kind) {
+        case Unit.Kind.Atomic: return unit.model.atomicConformation
+        case Unit.Kind.Spheres: return unit.model.coarseConformation.spheres
+        case Unit.Kind.Gaussians: return unit.model.coarseConformation.gaussians
+    }
+}
+
+export function computeUnitGaussianDensity(unit: Unit, props: GaussianDensityProps) {
+    const conformation = getConformation(unit)
+    const { elements } = unit
+    const position = {
+        indices: elements,
+        x: conformation.x,
+        y: conformation.y,
+        z: conformation.z
+    }
+
+    const l = StructureElement.create(unit)
+    const sizeTheme = SizeTheme({ name: 'physical' })
+    const radius = (index: number) => {
+        l.element = index as ElementIndex
+        return sizeTheme.size(l)
+    }
+
+    return Task.create('Gaussian Density', async ctx => {
+        return await GaussianDensity(ctx, position, unit.lookup3d.boundary.box, radius, { ...props, box: unit.lookup3d.boundary.box });
+    });
+}
+
+export async function computeUnitGaussianDensityCached(unit: Unit, props: GaussianDensityProps, cache: Map<string, DensityData>, ctx?: RuntimeContext) {
+    const key = `${props.radiusOffset}|${props.resolutionFactor}|${props.smoothness}`
+    let density = cache.get(key)
+    if (density) return density
+    density = ctx ? await computeUnitGaussianDensity(unit, props).runInContext(ctx) : await computeUnitGaussianDensity(unit, props).run()
+    cache.set(key, density)
+    return density
+}