diff --git a/src/apps/canvas/component/structure-representation.tsx b/src/apps/canvas/component/structure-representation.tsx index 2cdb76e8f2b496b6850b50633ec3d3fc670f2aa5..84ea120c157d99754dbc4f24c307b0a30cdeb6ff 100644 --- a/src/apps/canvas/component/structure-representation.tsx +++ b/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> : ''} diff --git a/src/mol-geo/representation/structure/representation/surface.ts b/src/mol-geo/representation/structure/representation/surface.ts index 98c50c3926b1ba8d54f2715370f1afbc2a78bf2a..dc05389aa9530dea438d6b5265480a59ad45ee3e 100644 --- a/src/mol-geo/representation/structure/representation/surface.ts +++ b/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() } } diff --git a/src/mol-geo/representation/structure/visual/gaussian-density-point.ts b/src/mol-geo/representation/structure/visual/gaussian-density-point.ts index 838df00b7ce741689dafd5fd3463fadea6517ed8..86d32e2818dfeb76d333f91ef4b9f47973cf3b46 100644 --- a/src/mol-geo/representation/structure/visual/gaussian-density-point.ts +++ b/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 diff --git a/src/mol-geo/representation/structure/visual/gaussian-surface-mesh.ts b/src/mol-geo/representation/structure/visual/gaussian-surface-mesh.ts index 8c99efa9716391433384bde15728e10951173277..42a281055660c184a203a94f017b5698317def6e 100644 --- a/src/mol-geo/representation/structure/visual/gaussian-surface-mesh.ts +++ b/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), diff --git a/src/mol-geo/representation/structure/visual/gaussian-surface-wireframe.ts b/src/mol-geo/representation/structure/visual/gaussian-surface-wireframe.ts index b1ddba32c48ceaef6ec83745e41e98beca931722..5c671537c9c07ce38e6fe1c2965ffb0be6a9466c 100644 --- a/src/mol-geo/representation/structure/visual/gaussian-surface-wireframe.ts +++ b/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), diff --git a/src/mol-geo/representation/structure/visual/util/gaussian.ts b/src/mol-geo/representation/structure/visual/util/gaussian.ts index ce06f079db8ad584c7b3814cf4c301aab98e6aa4..54cd72a4b8ba4dcf581d87842432de530588078b 100644 --- a/src/mol-geo/representation/structure/visual/util/gaussian.ts +++ b/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 } -} \ No newline at end of file diff --git a/src/mol-math/geometry/common.ts b/src/mol-math/geometry/common.ts index e203afdf228197b4b18a08bfd33ea0b50b4a455c..ae32e5438e3ccd75dbbd339ac1e02bf524f3e7f1 100644 --- a/src/mol-math/geometry/common.ts +++ b/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 } \ No newline at end of file diff --git a/src/mol-math/geometry/gaussian-density.ts b/src/mol-math/geometry/gaussian-density.ts new file mode 100644 index 0000000000000000000000000000000000000000..d6d2d8cd1249949b3239aa096bed3c906e859d8f --- /dev/null +++ b/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 } +} \ No newline at end of file diff --git a/src/mol-model/structure/structure/unit.ts b/src/mol-model/structure/structure/unit.ts index 4092d9d3cbb5521fe300c6c8cf45e35e9f845c0b..bce4b286b1aa335ef0e8f6ec589d37161593cf70 100644 --- a/src/mol-model/structure/structure/unit.ts +++ b/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), }; diff --git a/src/mol-model/structure/structure/unit/gaussian-density.ts b/src/mol-model/structure/structure/unit/gaussian-density.ts new file mode 100644 index 0000000000000000000000000000000000000000..df7bf0a89595b15a9b82d014902dd5c969e97fd1 --- /dev/null +++ b/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 +} \ No newline at end of file