diff --git a/src/apps/render-test/utils/mcubes.ts b/src/apps/render-test/utils/mcubes.ts index 350cdc9e53d622051ade629587041d872ff8eee5..92905c7b6f5d111486ee7a1bd6fab7b95a5cd5e6 100644 --- a/src/apps/render-test/utils/mcubes.ts +++ b/src/apps/render-test/utils/mcubes.ts @@ -5,7 +5,7 @@ */ import { Run } from 'mol-task' -import { compute } from 'mol-geo/util/marching-cubes/algorithm' +import { computeMarchingCubes } from 'mol-geo/util/marching-cubes/algorithm' import { Mesh } from 'mol-geo/shape/mesh' import { Tensor, Mat4, Vec3 } from 'mol-math/linear-algebra' @@ -38,7 +38,7 @@ export default async function computeSurface(f: (x: number, y: number, z: number const min = Vec3.create(-1.1, -1.1, -1.1), max = Vec3.create(1.1, 1.1, 1.1); fillField(field, f, min, max); - const surface = await Run(compute({ + const surface = await Run(computeMarchingCubes({ scalarField: field, isoLevel: 0, oldSurface: data ? data.surface : void 0 diff --git a/src/apps/structure-info/volume.ts b/src/apps/structure-info/volume.ts index e3b4b2e1efdba2514a0f052feea69932ac21f2cd..bfe15e630c66fb2b3e41cea24969a50c99fb4852 100644 --- a/src/apps/structure-info/volume.ts +++ b/src/apps/structure-info/volume.ts @@ -1,35 +1,77 @@ - /** * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. * * @author David Sehnal <david.sehnal@gmail.com> */ +import * as fs from 'fs' import * as argparse from 'argparse' -import { VolumeData, parseDensityServerData } from 'mol-model/volume' +import * as util from 'util' + +import { VolumeData, parseDensityServerData, VolumeIsoValue } from 'mol-model/volume' import { downloadCif } from './helpers' import CIF from 'mol-io/reader/cif' import { DensityServer_Data_Database } from 'mol-io/reader/cif/schema/density-server'; import { Run } from 'mol-task'; import { Table } from 'mol-data/db'; +import { computeVolumeMesh } from 'mol-geo/representation/volume/Mesh'; +import { StringBuilder } from 'mol-util'; + +require('util.promisify').shim(); +const writeFileAsync = util.promisify(fs.writeFile); type Volume = { source: DensityServer_Data_Database, volume: VolumeData } async function getVolume(url: string): Promise<Volume> { - const cif = await downloadCif(url, false); + const cif = await downloadCif(url, true); const data = CIF.schema.densityServer(cif.blocks[1]); return { source: data, volume: await Run(parseDensityServerData(data)) }; } -function print(volume: Volume) { - const { volume_data_3d_info } = volume.source; +function print(data: Volume) { + const { volume_data_3d_info } = data.source; const row = Table.getRow(volume_data_3d_info, 0); console.log(row); + console.log(data.volume.cell); + console.log(data.volume.dataStats); + console.log(data.volume.fractionalBox); +} + +async function doMesh(data: Volume, filename: string) { + const mesh = await Run(computeVolumeMesh(data.volume, VolumeIsoValue.relative(data.volume.dataStats, 1.5))); + console.log({ vc: mesh.vertexCount, tc: mesh.triangleCount }); + + // Export the mesh in OBJ format. + const { vertexCount, triangleCount } = mesh; + + const vs = mesh.vertexBuffer.ref.value; + const ts = mesh.indexBuffer.ref.value; + + const obj = StringBuilder.create(); + for (let i = 0; i < vertexCount; i++) { + StringBuilder.write(obj, 'v '); + StringBuilder.writeFloat(obj, vs[3 * i + 0], 100); + StringBuilder.whitespace1(obj); + StringBuilder.writeFloat(obj, vs[3 * i + 1], 100); + StringBuilder.whitespace1(obj); + StringBuilder.writeFloat(obj, vs[3 * i + 2], 100); + StringBuilder.newline(obj); + } + for (let i = 0; i < triangleCount; i++) { + StringBuilder.write(obj, 'f '); + StringBuilder.writeIntegerAndSpace(obj, ts[3 * i + 0] + 1); + StringBuilder.writeIntegerAndSpace(obj, ts[3 * i + 1] + 1); + StringBuilder.writeInteger(obj, ts[3 * i + 2] + 1); + StringBuilder.newline(obj); + } + + await writeFileAsync(filename, StringBuilder.getString(obj)); } -async function run(url: string) { +async function run(url: string, meshFilename: string) { const volume = await getVolume(url); print(volume); + await doMesh(volume, meshFilename); } const parser = new argparse.ArgumentParser({ @@ -39,9 +81,14 @@ description: 'Info about VolumeData from mol-model module' parser.addArgument([ '--emdb', '-e' ], { help: 'EMDB id, for example 8116', }); +parser.addArgument([ '--mesh' ], { + help: 'Mesh filename', + required: true +}); interface Args { - emdb?: string + emdb?: string, + mesh: string } const args: Args = parser.parseArgs(); -run(`https://webchem.ncbr.muni.cz/DensityServer/em/emd-${args.emdb}/cell?detail=1&encoding=cif`); \ No newline at end of file +run(`https://webchem.ncbr.muni.cz/DensityServer/em/emd-${args.emdb}/cell?detail=4`, args.mesh); \ No newline at end of file diff --git a/src/mol-geo/representation/volume/Mesh.ts b/src/mol-geo/representation/volume/Mesh.ts new file mode 100644 index 0000000000000000000000000000000000000000..317b51b3ecb968e1daaa8a6c807bc1f9c6bb529a --- /dev/null +++ b/src/mol-geo/representation/volume/Mesh.ts @@ -0,0 +1,29 @@ +/** + * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. + * + * @author David Sehnal <david.sehnal@gmail.com> + */ + +import { VolumeData, VolumeIsoValue } from 'mol-model/volume' +import { Task } from 'mol-task' +import { Mesh } from '../../shape/mesh'; +import { computeMarchingCubes } from '../../util/marching-cubes/algorithm'; + +function computeVolumeMesh(volume: VolumeData, isoValue: VolumeIsoValue) { + return Task.create<Mesh>('Volume Surface', async ctx => { + ctx.update({ message: 'Marching cubes...' }); + + const mesh = await ctx.runChild(computeMarchingCubes({ + isoLevel: VolumeIsoValue.toAbsolute(isoValue).absoluteValue, + scalarField: volume.data + })); + + const transform = VolumeData.getGridToCartesianTransform(volume); + ctx.update({ message: 'Transforming mesh...' }); + Mesh.transformImmediate(mesh, transform); + + return mesh; + }); +} + +export { computeVolumeMesh } \ No newline at end of file diff --git a/src/mol-geo/util/marching-cubes/algorithm.ts b/src/mol-geo/util/marching-cubes/algorithm.ts index 25b266773ce26c0afd7a22332083f915f63aa913..974864858173b7d6e70dfa5ea9d22ac10ad77ea9 100644 --- a/src/mol-geo/util/marching-cubes/algorithm.ts +++ b/src/mol-geo/util/marching-cubes/algorithm.ts @@ -26,7 +26,7 @@ export interface MarchingCubesParameters { oldSurface?: Mesh } -export function compute(parameters: MarchingCubesParameters) { +export function computeMarchingCubes(parameters: MarchingCubesParameters) { return Task.create('Marching Cubes', async ctx => { let comp = new MarchingCubesComputation(parameters, ctx); return await comp.run(); diff --git a/src/mol-math/geometry/primitives/box3d.ts b/src/mol-math/geometry/primitives/box3d.ts index efceedb4f763741e1ff6ba3b7476f58a287ac18d..6d0b5bd5f7da032fed3942b44558ef6d4f279ecc 100644 --- a/src/mol-math/geometry/primitives/box3d.ts +++ b/src/mol-math/geometry/primitives/box3d.ts @@ -12,7 +12,7 @@ interface Box3D { min: Vec3, max: Vec3 } namespace Box3D { export function create(min: Vec3, max: Vec3): Box3D { return { min, max }; } - + export function computeBounding(data: PositionData): Box3D { const min = [Number.MAX_VALUE, Number.MAX_VALUE, Number.MAX_VALUE]; const max = [-Number.MAX_VALUE, -Number.MAX_VALUE, -Number.MAX_VALUE]; diff --git a/src/mol-math/geometry/spacegroup/construction.ts b/src/mol-math/geometry/spacegroup/construction.ts index 626782b73af8f1392da9dd353ed24f1d7397a4f3..65aa22e849283499982fbf6fefdd71ade328105b 100644 --- a/src/mol-math/geometry/spacegroup/construction.ts +++ b/src/mol-math/geometry/spacegroup/construction.ts @@ -15,8 +15,7 @@ interface SpacegroupCell { /** Transfrom cartesian -> fractional coordinates within the cell */ readonly toFractional: Mat4, /** Transfrom fractional coordinates within the cell -> cartesian */ - readonly fromFractional: Mat4, - readonly fractionalBasis: Vec3[] + readonly fromFractional: Mat4 } interface Spacegroup { @@ -59,14 +58,8 @@ namespace SpacegroupCell { ]); const toFractional = Mat4.invert(Mat4.zero(), fromFractional)!; - const fractionalBasis = [ - Vec3.transformMat4(Vec3.zero(), Vec3.create(1, 0, 0), toFractional), - Vec3.transformMat4(Vec3.zero(), Vec3.create(0, 1, 0), toFractional), - Vec3.transformMat4(Vec3.zero(), Vec3.create(0, 0, 1), toFractional) - ]; - const num = typeof nameOrNumber === 'number' ? nameOrNumber : SpacegroupNumbers[nameOrNumber]; - return { number: num, size, anglesInRadians, toFractional, fromFractional, fractionalBasis }; + return { number: num, size, anglesInRadians, toFractional, fromFractional }; } } diff --git a/src/mol-model/volume/_spec/data.spec.ts b/src/mol-model/volume/_spec/data.spec.ts deleted file mode 100644 index 879b4fb842e26b1581ef57647aa48c69ba9ab9f5..0000000000000000000000000000000000000000 --- a/src/mol-model/volume/_spec/data.spec.ts +++ /dev/null @@ -1,6 +0,0 @@ -/** - * Copyright (c) 2018 mol* contributors, licensed under MIT, See LICENSE file for more info. - * - * @author David Sehnal <david.sehnal@gmail.com> - */ - diff --git a/src/mol-model/volume/data.ts b/src/mol-model/volume/data.ts index 19826278bb33abda9b02a791c90066cb287cda9a..e0313d893b7e373774868005566cfbde0f02c1eb 100644 --- a/src/mol-model/volume/data.ts +++ b/src/mol-model/volume/data.ts @@ -30,4 +30,36 @@ namespace VolumeData { } } -export { VolumeData } \ No newline at end of file +type VolumeIsoValue = VolumeIsoValue.Absolute | VolumeIsoValue.Relative + +namespace VolumeIsoValue { + export type Relative = Readonly<{ kind: 'relative', stats: VolumeData['dataStats'], relativeValue: number }> + export type Absolute = Readonly<{ kind: 'absolute', stats: VolumeData['dataStats'], absoluteValue: number }> + + export function absolute(stats: VolumeData['dataStats'], value: number): Absolute { return { kind: 'absolute', stats, absoluteValue: value }; } + export function relative(stats: VolumeData['dataStats'], value: number): Relative { return { kind: 'relative', stats, relativeValue: value }; } + + export function toAbsolute(value: VolumeIsoValue): Absolute { + if (value.kind === 'absolute') return value; + + const { mean, sigma } = value.stats + return { + kind: 'absolute', + stats: value.stats, + absoluteValue: value.relativeValue * sigma + mean + } + } + + export function toRelative(value: VolumeIsoValue): Relative { + if (value.kind === 'relative') return value; + + const { mean, sigma } = value.stats + return { + kind: 'relative', + stats: value.stats, + relativeValue: (mean - value.absoluteValue) / sigma + } + } +} + +export { VolumeData, VolumeIsoValue } \ No newline at end of file diff --git a/src/mol-model/volume/formats/density-server.ts b/src/mol-model/volume/formats/density-server.ts index c7972e690ef4ae5f40ced2651f8c49d5c6cea6fc..300d7039c56e329b2eea986c2753b9499225d5c1 100644 --- a/src/mol-model/volume/formats/density-server.ts +++ b/src/mol-model/volume/formats/density-server.ts @@ -13,7 +13,9 @@ import { Tensor, Vec3 } from 'mol-math/linear-algebra'; function parseDensityServerData(source: DensityServer_Data_Database): Task<VolumeData> { return Task.create<VolumeData>('Parse Volume Data', async ctx => { const { volume_data_3d_info: info, volume_data_3d: values } = source; - const cell = SpacegroupCell.create(info.spacegroup_number.value(0), Vec3.ofArray(info.spacegroup_cell_size.value(0)), Vec3.ofArray(info.spacegroup_cell_angles.value(0))); + const cell = SpacegroupCell.create(info.spacegroup_number.value(0), + Vec3.ofArray(info.spacegroup_cell_size.value(0)), + Vec3.scale(Vec3.zero(), Vec3.ofArray(info.spacegroup_cell_angles.value(0)), Math.PI / 180)); const tensorSpace = Tensor.Space(info.sample_count.value(0), info.axis_order.value(0), Float32Array); const data = Tensor.create(tensorSpace, Tensor.Data1(values.values.toArray()));