IMP logo
IMP Reference Guide  develop.27926d84dc,2024/04/15
The Integrative Modeling Platform
IMP::em::SampledDensityMap Class Reference

Class for sampling a density map from particles. More...

#include <IMP/em/SampledDensityMap.h>

+ Inheritance diagram for IMP::em::SampledDensityMap:

Detailed Description

Class for sampling a density map from particles.

Definition at line 31 of file SampledDensityMap.h.

Public Member Functions

 SampledDensityMap (KernelType kt=GAUSSIAN)
 Creates a new density map for sampled map. More...
 
 SampledDensityMap (const DensityHeader &header, KernelType kt=GAUSSIAN)
 The size of the map is determined by the header and the data is allocated. More...
 
 SampledDensityMap (const ParticlesTemp &ps, double resolution, double voxel_size, IMP::FloatKey mass_key=IMP::atom::Mass::get_mass_key(), int sig_cutoff=3, KernelType kt=GAUSSIAN)
 Generate a sampled density map from the particles. More...
 
KernelParametersget_kernel_params ()
 
float get_minimum_resampled_value ()
 Get minimum density value between voxels that correspond to particles. More...
 
const Particlesget_sampled_particles () const
 
FloatKey get_weight_key () const
 
const core::XYZRsget_xyzr_particles () const
 
void project (const ParticlesTemp &ps, int x_margin, int y_margin, int z_margin, algebra::Vector3D shift=algebra::Vector3D(0., 0., 0.), FloatKey mass_key=atom::Mass::get_mass_key())
 Project particles on the grid by their mass value. More...
 
virtual void resample ()
 Resample beads on an EM grid. More...
 
void set_particles (const ParticlesTemp &ps, IMP::FloatKey mass_key=IMP::atom::Mass::get_mass_key())
 setting particles in case they were not set by the constructor More...
 
void update_resolution (Float res)
 Update the simulation resolution. More...
 
- Public Member Functions inherited from IMP::em::DensityMap
 DensityMap (std::string name="DensityMap%1%")
 
 DensityMap (const DensityHeader &header, std::string name="DensityMap%1%")
 Construct a density map as instructed in the input header. More...
 
void add (const DensityMap *other)
 Sums two grids; the result is kept in the map. More...
 
void add (Float d)
 Add a number to every voxel in the map. More...
 
void calc_all_voxel2loc ()
 Calculates the coordinates that correspond to all voxels. More...
 
double calcRMS ()
 Calculates RMSD and mean of a map values as stored in the header. More...
 
void copy_map (const DensityMap *other)
 Copy map into this map. More...
 
algebra::Vector3D get_centroid (double threshold=0.0) const
 Get the centroid of all voxels with density above a given threshold. More...
 
DensityMapget_cropped (float threshold)
 Create and return a new cropped map. More...
 
DensityMapget_cropped (Particles ps, double distance, bool inverse=false, bool keep_map_dimensions=false)
 Create and return a new cropped map based on particles. More...
 
DensityMapget_cropped (const algebra::BoundingBox3D &bb)
 Create and return a new cropped map with the given bounding box extent. More...
 
int get_dim_index_by_location (const algebra::Vector3D &v, int ind) const
 Calculate dimension index of a given location. More...
 
int get_dim_index_by_location (float loc_val, int ind) const
 
const DensityHeaderget_header () const
 Returns a read-only pointer to the header of the map. More...
 
DensityHeaderget_header_writable ()
 Returns a pointer to the header of the map in a writable version. More...
 
algebra::Vector3D get_location_by_voxel (long index) const
 Calculate the location of a given voxel. More...
 
float get_location_in_dim_by_voxel (long index, int dim) const
 Calculate the location of a given voxel in a given dimension. More...
 
std::string get_locations_string (float t)
 Print the locations of all voxels with value above a given threshold. More...
 
double get_max_value () const
 Returns the value of the voxel with the highest density. More...
 
float get_maximum_value_in_xy_plane (int z_ind)
 Get the maximum value in a XY plane indicated by a Z index. More...
 
float get_maximum_value_in_xz_plane (int y_ind)
 Get the maximum value in a XZ plane indicated by a Y index. More...
 
float get_maximum_value_in_yz_plane (int x_ind)
 Get the maximum value in a YZ plane indicated by a X index. More...
 
double get_min_value () const
 Returns the value of the voxel with the lowest density. More...
 
long get_number_of_voxels () const
 Get the number of map voxels. More...
 
algebra::Vector3D get_origin () const
 
bool get_rms_calculated () const
 
Float get_spacing () const
 Get the voxel size of the map. More...
 
algebra::Vector3D get_top () const
 
virtual std::string get_type_name () const override
 
double get_value (float x, float y, float z) const
 Gets the value of the voxel located at (x,y,z) More...
 
double get_value (const algebra::Vector3D &point) const
 
double get_value (long index) const
 Gets the value of the voxel at a given index. More...
 
virtual ::IMP::VersionInfo get_version_info () const override
 Get information about the module and version of the object. More...
 
long get_voxel_by_location (float x, float y, float z) const
 Calculate the voxel of a given location. More...
 
long get_voxel_by_location (const algebra::Vector3D &v) const
 Calculate the voxel of a given location. More...
 
bool is_normalized () const
 
bool is_part_of_volume (float x, float y, float z) const
 Checks whether a given point is in the grid. More...
 
bool is_part_of_volume (const algebra::Vector3D &v) const
 Checks whether a given point is in the grid. More...
 
bool is_xyz_ind_part_of_volume (int ix, int iy, int iz) const
 
int lower_voxel_shift (double loc, double kdist, double orig, int ndim) const
 
void multiply (float factor)
 Multiply in place each voxel in the map by the input factor. More...
 
void pad (int nx, int ny, int nz, float val=0.0)
 Increase the dimension of the map. More...
 
DensityMappad_margin (int mrg_x, int mrg_y, int mrg_z, float val=0.0)
 Create a new padded map. More...
 
void pick_max (const DensityMap *other)
 Pick the max value between two corresponding voxels between two maps. More...
 
void reset_data (float value=0.0)
 Set the density voxels to some value and reset the management flags. More...
 
bool same_dimensions (const DensityMap *other) const
 Checks if two maps have the same dimensions. More...
 
bool same_origin (const DensityMap *other) const
 Checks if two maps have the same origin. More...
 
bool same_voxel_size (const DensityMap *other) const
 Checks if two maps have the same voxel size. More...
 
void set_origin (float x, float y, float z)
 Sets the origin of the header. More...
 
void set_origin (const IMP::algebra::Vector3D &v)
 
void set_value (long index, double value)
 Set the value of the voxel at a given index. More...
 
void set_value (float x, float y, float z, double value)
 Set the value of the voxel at given coordinates. More...
 
void set_void_map (int nx, int ny, int nz)
 Set the map dimension and reset all voxels to 0. More...
 
void std_normalize ()
 Normalize the density voxels according to standard deviation (stdv). More...
 
void update_voxel_size (float new_apix)
 Update the voxel size of the map. More...
 
int upper_voxel_shift (double loc, double kdist, double orig, int ndim) const
 
long xyz_ind2voxel (int x, int y, int z) const
 Calculate the voxel of a given xyz index. More...
 
- Public Member Functions inherited from IMP::Object
virtual void clear_caches ()
 
CheckLevel get_check_level () const
 
LogLevel get_log_level () const
 
void set_check_level (CheckLevel l)
 
void set_log_level (LogLevel l)
 Set the logging level used in this object. More...
 
void set_was_used (bool tf) const
 
void show (std::ostream &out=std::cout) const
 
const std::string & get_name () const
 
void set_name (std::string name)
 

Protected Member Functions

IMP::algebra::BoundingBoxD< 3 > calculate_particles_bounding_box (const Particles &ps)
 Calculate the parameters of the particles bounding box. More...
 
void determine_grid_size (double resolution, double voxel_size, int sig_cutoff)
 
void set_header (const algebra::Vector3D &lower_bound, const algebra::Vector3D &upper_bound, double maxradius, double resolution, double voxel_size, int sig_offset)
 
void set_neighbor_mask (float radius)
 
- Protected Member Functions inherited from IMP::em::DensityMap
void allocated_data ()
 
void float2real (float *f_data, boost::scoped_array< double > &r_data)
 
void real2float (double *r_data, boost::scoped_array< float > &f_data)
 
void reset_all_voxel2loc ()
 
void update_header ()
 update the header values – still in work More...
 
- Protected Member Functions inherited from IMP::Object
 Object (std::string name)
 Construct an object with the given name. More...
 
virtual void do_destroy ()
 

Protected Attributes

KernelParameters kernel_params_
 kernel handling More...
 
KernelType kt_
 
Particles ps_
 
FloatKey weight_key_
 
FloatKey x_key_
 
core::XYZRs xyzr_
 
FloatKey y_key_
 
FloatKey z_key_
 
- Protected Attributes inherited from IMP::em::DensityMap
boost::scoped_array< double > data_
 
bool data_allocated_
 
DensityHeader header_
 
bool loc_calculated_
 true if the locations have already been computed More...
 
bool normalized_
 
bool rms_calculated_
 
boost::scoped_array< float > x_loc_
 
boost::scoped_array< float > y_loc_
 
boost::scoped_array< float > z_loc_
 

Additional Inherited Members

Constructor & Destructor Documentation

IMP::em::SampledDensityMap::SampledDensityMap ( KernelType  kt = GAUSSIAN)

Creates a new density map for sampled map.

The header of the map is not determined and no data is being allocated

Definition at line 37 of file SampledDensityMap.h.

IMP::em::SampledDensityMap::SampledDensityMap ( const DensityHeader header,
KernelType  kt = GAUSSIAN 
)

The size of the map is determined by the header and the data is allocated.

IMP::em::SampledDensityMap::SampledDensityMap ( const ParticlesTemp ps,
double  resolution,
double  voxel_size,
IMP::FloatKey  mass_key = IMP::atom::Mass::get_mass_key(),
int  sig_cutoff = 3,
KernelType  kt = GAUSSIAN 
)

Generate a sampled density map from the particles.

Parameters
[in]psparticles with XYZ, radius and weight attributes
[in]resolutionhalf width the Gaussian
[in]voxel_sizevoxel size in angstroms
[in]mass_keykey to use to weight particles
[in]sig_cutoffChoose what should be the sigma cutoff for accurate sampling. It is used in two functions; (i) to determine the size of the grid dimensions (ii) to determine the voxels around the coords participating in the sampling procedure.
[in]ktType of kernel to use.

Member Function Documentation

IMP::algebra::BoundingBoxD<3> IMP::em::SampledDensityMap::calculate_particles_bounding_box ( const Particles ps)
protected

Calculate the parameters of the particles bounding box.

Parameters
[in]psparticles with XYZ, radius and weight attributes
Returns
the particles bounding box
void IMP::em::SampledDensityMap::determine_grid_size ( double  resolution,
double  voxel_size,
int  sig_cutoff 
)
protected

Determine the size of the grid as a function of the particles and the resolution.

float IMP::em::SampledDensityMap::get_minimum_resampled_value ( )

Get minimum density value between voxels that correspond to particles.

void IMP::em::SampledDensityMap::project ( const ParticlesTemp ps,
int  x_margin,
int  y_margin,
int  z_margin,
algebra::Vector3D  shift = algebra::Vector3D(0., 0., 0.),
FloatKey  mass_key = atom::Mass::get_mass_key() 
)

Project particles on the grid by their mass value.

Parameters
psthe particles to project
[in]x_marginsampling is restricted to [x_margin,nx-x_margin]
[in]y_marginsampling is restricted to [y_margin,ny-y_margin]
[in]z_marginsampling is restricted to [z_margin,nz-z_margin]
[in]shiftthe positions of all particles are shifted by this value before projection
[in]mass_keykey to obtain particle mass
virtual void IMP::em::SampledDensityMap::resample ( )
virtual

Resample beads on an EM grid.

Note
The density of a particle p centered at pl at position gl is: \(\frac{{Z}e^{\frac{{-0.5}({p_l}-{g_l})}{\sigma}}}{\sqrt{{2}{\pi}\sigma}}\) , such that \({Z}\) is the weight of the particle and \({\sigma}\) is defined to be \({0.425}\) times the resolution, to follow the 'full width at half maxima' criterion. For more details please refer to Topf et al, Structure, 2008.

Reimplemented in IMP::em::SurfaceShellDensityMap.

void IMP::em::SampledDensityMap::set_particles ( const ParticlesTemp ps,
IMP::FloatKey  mass_key = IMP::atom::Mass::get_mass_key() 
)

setting particles in case they were not set by the constructor

void IMP::em::SampledDensityMap::update_resolution ( Float  res)

Update the simulation resolution.

Member Data Documentation

KernelParameters IMP::em::SampledDensityMap::kernel_params_
protected

kernel handling

Definition at line 128 of file SampledDensityMap.h.


The documentation for this class was generated from the following file: