/** * Copyright (c) 2018-2020 mol* contributors, licensed under MIT, See LICENSE file for more info. * * @author Alexander Rose */ import { Unit, Structure, ElementIndex, StructureElement, ResidueIndex } from '../../../../mol-model/structure'; import { Mat4, Vec3 } from '../../../../mol-math/linear-algebra'; import { TransformData, createTransform } from '../../../../mol-geo/geometry/transform-data'; import { OrderedSet, SortedArray } from '../../../../mol-data/int'; import { EmptyLoci, Loci } from '../../../../mol-model/loci'; import { PhysicalSizeTheme } from '../../../../mol-theme/size/physical'; import { AtomicNumbers } from '../../../../mol-model/structure/model/properties/atomic'; import { fillSerial } from '../../../../mol-util/array'; import { ParamDefinition as PD } from '../../../../mol-util/param-definition'; import { AssignableArrayLike } from '../../../../mol-util/type-helpers'; import { getBoundary } from '../../../../mol-math/geometry/boundary'; import { Box3D } from '../../../../mol-math/geometry'; /** Return a Loci for the elements of a whole residue the elementIndex belongs to. */ export function getResidueLoci(structure: Structure, unit: Unit.Atomic, elementIndex: ElementIndex): Loci { const { elements, model } = unit; if (OrderedSet.indexOf(elements, elementIndex) !== -1) { const { index, offsets } = model.atomicHierarchy.residueAtomSegments; const rI = index[elementIndex]; const _indices: number[] = []; for (let i = offsets[rI], il = offsets[rI + 1]; i < il; ++i) { const unitIndex = OrderedSet.indexOf(elements, i); if (unitIndex !== -1) _indices.push(unitIndex); } const indices = OrderedSet.ofSortedArray(SortedArray.ofSortedArray(_indices)); return StructureElement.Loci(structure, [{ unit, indices }]); } return EmptyLoci; } /** * Return a Loci for the elements of a whole residue the elementIndex belongs to but * restrict to elements that have the same label_alt_id or none */ export function getAltResidueLoci(structure: Structure, unit: Unit.Atomic, elementIndex: ElementIndex) { const { elements, model } = unit; const { label_alt_id } = model.atomicHierarchy.atoms; const elementAltId = label_alt_id.value(elementIndex); if (OrderedSet.indexOf(elements, elementIndex) !== -1) { const { index } = model.atomicHierarchy.residueAtomSegments; const rI = index[elementIndex]; return getAltResidueLociFromId(structure, unit, rI, elementAltId); } return StructureElement.Loci(structure, []); } export function getAltResidueLociFromId(structure: Structure, unit: Unit.Atomic, residueIndex: ResidueIndex, elementAltId: string) { const { elements, model } = unit; const { label_alt_id } = model.atomicHierarchy.atoms; const { offsets } = model.atomicHierarchy.residueAtomSegments; const _indices: number[] = []; for (let i = offsets[residueIndex], il = offsets[residueIndex + 1]; i < il; ++i) { const unitIndex = OrderedSet.indexOf(elements, i); if (unitIndex !== -1) { const altId = label_alt_id.value(i); if (elementAltId === altId || altId === '') { _indices.push(unitIndex); } } } const indices = OrderedSet.ofSortedArray(SortedArray.ofSortedArray(_indices)); return StructureElement.Loci(structure, [{ unit, indices }]); } // export function createUnitsTransform({ units }: Unit.SymmetryGroup, transformData?: TransformData) { const unitCount = units.length; const n = unitCount * 16; const array = transformData && transformData.aTransform.ref.value.length >= n ? transformData.aTransform.ref.value : new Float32Array(n); for (let i = 0; i < unitCount; i++) { Mat4.toArray(units[i].conformation.operator.matrix, array, i * 16); } return createTransform(array, unitCount, transformData); } export const UnitKindInfo = { 'atomic': {}, 'spheres': {}, 'gaussians': {}, }; export type UnitKind = keyof typeof UnitKindInfo export const UnitKindOptions = PD.objectToOptions(UnitKindInfo); export function includesUnitKind(unitKinds: UnitKind[], unit: Unit) { for (let i = 0, il = unitKinds.length; i < il; ++i) { if (Unit.isAtomic(unit) && unitKinds[i] === 'atomic') return true; if (Unit.isSpheres(unit) && unitKinds[i] === 'spheres') return true; if (Unit.isGaussians(unit) && unitKinds[i] === 'gaussians') return true; } return false; } // const MaxCells = 500_000_000; /** guard against overly high resolution for the given box size */ export function ensureReasonableResolution(box: Box3D, props: { resolution: number } & T) { const volume = Box3D.volume(box); const approxCells = volume / props.resolution; const resolution = approxCells > MaxCells ? volume / MaxCells : props.resolution; return { ...props, resolution }; } export 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 const CommonSurfaceParams = { ignoreHydrogens: PD.Boolean(false, { description: 'Whether or not to include hydrogen atoms in the surface calculation.' }), traceOnly: PD.Boolean(false, { description: 'Whether or not to only use trace atoms in the surface calculation.' }), includeParent: PD.Boolean(false, { description: 'Include elements of the parent structure in surface calculation to get a surface patch of the current structure.' }), }; export const DefaultCommonSurfaceProps = PD.getDefaultValues(CommonSurfaceParams); export type CommonSurfaceProps = typeof DefaultCommonSurfaceProps const v = Vec3(); function squaredDistance(x: number, y: number, z: number, center: Vec3) { return Vec3.squaredDistance(Vec3.set(v, x, y, z), center); } /** marks `indices` for filtering/ignoring in `id` when not in `elements` */ function filterId(id: AssignableArrayLike, elements: SortedArray, indices: SortedArray) { let start = 0; const end = elements.length; for (let i = 0, il = indices.length; i < il; ++i) { const idx = SortedArray.indexOfInRange(elements, indices[i], start, end); if (idx === -1) { id[i] = -2; } else { id[i] = idx; start = idx; } } } export function getUnitConformationAndRadius(structure: Structure, unit: Unit, props: CommonSurfaceProps) { const { ignoreHydrogens, traceOnly, includeParent } = props; const rootUnit = includeParent ? structure.root.unitMap.get(unit.id) : unit; const { x, y, z } = getConformation(rootUnit); const { elements } = rootUnit; const { center, radius: sphereRadius } = unit.boundary.sphere; const extraRadius = (2 + 1.5) * 2; // TODO should be twice (the max vdW/sphere radius plus the probe radius) const radiusSq = (sphereRadius + extraRadius) * (sphereRadius + extraRadius); let indices: SortedArray; let id: AssignableArrayLike; if (ignoreHydrogens || traceOnly || (includeParent && rootUnit !== unit)) { const _indices = []; const _id = []; for (let i = 0, il = elements.length; i < il; ++i) { const eI = elements[i]; if (ignoreHydrogens && isHydrogen(rootUnit, eI)) continue; if (traceOnly && !isTrace(rootUnit, eI)) continue; if (includeParent && squaredDistance(x[eI], y[eI], z[eI], center) > radiusSq) continue; _indices.push(eI); _id.push(i); } indices = SortedArray.ofSortedArray(_indices); id = _id; } else { indices = elements; id = fillSerial(new Int32Array(indices.length)); } if (includeParent && rootUnit !== unit) { filterId(id, unit.elements, indices); } const position = { indices, x, y, z, id }; const boundary = unit === rootUnit ? unit.boundary : getBoundary(position); const l = StructureElement.Location.create(structure, rootUnit); const sizeTheme = PhysicalSizeTheme({}, { scale: 1 }); const radius = (index: number) => { l.element = index as ElementIndex; return sizeTheme.size(l); }; return { position, boundary, radius }; } export function getStructureConformationAndRadius(structure: Structure, ignoreHydrogens: boolean, traceOnly: boolean) { const l = StructureElement.Location.create(structure); const sizeTheme = PhysicalSizeTheme({}, { scale: 1 }); let xs: ArrayLike; let ys: ArrayLike; let zs: ArrayLike; let rs: ArrayLike; let id: ArrayLike; if (ignoreHydrogens || traceOnly) { const _xs: number[] = []; const _ys: number[] = []; const _zs: number[] = []; const _rs: number[] = []; const _id: number[] = []; for (let i = 0, m = 0, il = structure.units.length; i < il; ++i) { const unit = structure.units[i]; const { elements } = unit; const { x, y, z } = unit.conformation; l.unit = unit; for (let j = 0, jl = elements.length; j < jl; ++j) { const eI = elements[j]; if (ignoreHydrogens && isHydrogen(unit, eI)) continue; if (traceOnly && !isTrace(unit, eI)) continue; _xs.push(x(eI)); _ys.push(y(eI)); _zs.push(z(eI)); l.element = eI; _rs.push(sizeTheme.size(l)); _id.push(m + j); } m += elements.length; } xs = _xs, ys = _ys, zs = _zs, rs = _rs; id = _id; } else { const { elementCount } = structure; const _xs = new Float32Array(elementCount); const _ys = new Float32Array(elementCount); const _zs = new Float32Array(elementCount); const _rs = new Float32Array(elementCount); for (let i = 0, m = 0, il = structure.units.length; i < il; ++i) { const unit = structure.units[i]; const { elements } = unit; const { x, y, z } = unit.conformation; l.unit = unit; for (let j = 0, jl = elements.length; j < jl; ++j) { const eI = elements[j]; const mj = m + j; _xs[mj] = x(eI); _ys[mj] = y(eI); _zs[mj] = z(eI); l.element = eI; _rs[mj] = sizeTheme.size(l); } m += elements.length; } xs = _xs, ys = _ys, zs = _zs, rs = _rs; id = fillSerial(new Uint32Array(elementCount)); } const position = { indices: OrderedSet.ofRange(0, id.length), x: xs, y: ys, z: zs, id }; const radius = (index: number) => rs[index]; return { position, radius }; } const _H = AtomicNumbers['H']; export function isHydrogen(unit: Unit, element: ElementIndex) { if (Unit.isCoarse(unit)) return false; return unit.model.atomicHierarchy.derived.atom.atomicNumber[element] === _H; } export function isH(atomicNumber: ArrayLike, element: ElementIndex) { return atomicNumber[element] === _H; } export function isTrace(unit: Unit, element: ElementIndex) { if (Unit.isCoarse(unit)) return true; const atomId = unit.model.atomicHierarchy.atoms.label_atom_id.value(element); if (atomId === 'CA' || atomId === 'P') return true; return false; } export function getUnitExtraRadius(unit: Unit) { if (Unit.isAtomic(unit)) return 4; let max = 0; const { elements } = unit; const { r } = unit.conformation; for (let i = 0, _i = elements.length; i < _i; i++) { const _r = r(elements[i]); if (_r > max) max = _r; } return max + 1; } export function getStructureExtraRadius(structure: Structure) { let max = 0; for (const ug of structure.unitSymmetryGroups) { const r = getUnitExtraRadius(ug.units[0]); if (r > max) max = r; } return max; }