IMP logo
IMP Reference Guide  2.14.0
The Integrative Modeling Platform
rigid_fitting.h
Go to the documentation of this file.
1 /**
2  * \file IMP/em/rigid_fitting.h
3  * \brief preforms rigid fitting between a set of particles and a density map
4  *
5  * Copyright 2007-2020 IMP Inventors. All rights reserved.
6  */
7 
8 #ifndef IMPEM_RIGID_FITTING_H
9 #define IMPEM_RIGID_FITTING_H
10 
11 #include <IMP/core/MonteCarlo.h>
13 #include <IMP/RestraintSet.h>
15 #include <IMP/Particle.h>
16 #include <IMP/Model.h>
17 #include "DensityMap.h"
18 #include "FitRestraint.h"
19 #include <IMP/em/em_config.h>
20 #include <IMP/OptimizerState.h>
21 #include <IMP/core/rigid_bodies.h>
22 #include <IMP/ScoreState.h>
23 #include <algorithm>
24 IMPEM_BEGIN_NAMESPACE
25 
26 //! A simple list of fitting solutions.
27 /**
28  \see local_rigid_fitting_around_point
29  \see local_rigid_fitting_around_points
30  \see local_rigid_fitting_grid_search
31  \see compute_fitting_scores
32  */
33 class IMPEMEXPORT FittingSolutions {
34  typedef std::pair<algebra::Transformation3D, Float> FittingSolution;
35  struct sort_by_cc {
36  bool operator()(const FittingSolution &s1,
37  const FittingSolution &s2) const {
38  return s1.second < s2.second;
39  }
40  };
41 
42  public:
43  FittingSolutions() {}
44  //! Get the number of solutions in the set
45  inline int get_number_of_solutions() const { return fs_.size(); }
46  //! Get the score of the i'th solution
47  /**
48  \return the i'th transformation, or throw an exception
49  if the index is out of range
50  */
51  inline algebra::Transformation3D get_transformation(unsigned int i) const {
52  IMP_USAGE_CHECK(i < fs_.size(),
53  "The index requested ("
54  << i << ") in get_transformation is out of range ("
55  << fs_.size() << ")" << std::endl);
56  return fs_[i].first;
57  }
58  //! Get the score of the i'th solution
59  /**
60  \return the i'th score, or throw an exception
61  if the index is out of range
62  */
63  inline Float get_score(unsigned int i) const {
64  IMP_USAGE_CHECK(i < fs_.size(),
65  "The index requested ("
66  << i << ") in get_transformation is out of range ("
67  << fs_.size() << ")" << std::endl);
68  return fs_[i].second;
69  }
70  void set_score(unsigned int i, Float score) {
71  IMP_USAGE_CHECK(i < fs_.size(),
72  "The index requested ("
73  << i << ") in get_transformation is out of range ("
74  << fs_.size() << ")" << std::endl);
75  fs_[i].second = score;
76  }
77  //! Add a solution to the fitting solution set
78  void add_solution(const algebra::Transformation3D &t, Float score);
79  //! Sort solutions by cross-correlation scores
80  void sort(bool reverse = false);
81  //! Multiply each transformation (T) by t,
82  //!such that the new transformation are T*t
84  for (unsigned int i = 0; i < fs_.size(); i++)
85  fs_[i].first = fs_[i].first * t;
86  }
87  inline algebra::Transformation3Ds get_transformations() const {
89  for (unsigned int i = 0; i < fs_.size(); i++)
90  all_ts.push_back(fs_[i].first);
91  return all_ts;
92  }
93  void show(std::ostream &out = std::cout) const {
94  for (std::vector<FittingSolution>::const_iterator it = fs_.begin();
95  it != fs_.end(); it++) {
96  out << "(" << it->first << " , " << it->second << ")" << std::endl;
97  }
98  }
99 
100  protected:
101  std::vector<FittingSolution> fs_;
102 };
104 
105 //! Local rigid fitting of a rigid body around a center point
106 /**
107 \brief Fit a set of particles to a density map around an anchor point.
108  The fitting is assessed using the cross-correlation score.
109  The optimization is a standard MC/CG procedure.
110  The function returns a list of solutions sorted by the cross-correlation
111  score.
112 \note The returned cross-correlation score is 1-cc, as we usually want to
113  minimize a scoring function. Thus a score of 1 means no-correlation
114  and a score of 0. is perfect correlation.
115 \note The input rigid body should also be IMP::atom::Hierarchy
116 \param[in] p The rigid body to fit
117 \param[in] refiner Refiner to yield rigid body members
118 \param[in] weight_key The weight key of the particles in the rigid body
119 \param[in] dmap The density map to fit to
120 \param[in] anchor_centroid The point to fit the particles around
121 \param[in] display_log If provided, then intermediate states
122  in during the optimization procedure are printed
123 \param[in] number_of_optimization_runs number of Monte Carlo optimizations
124 \param[in] number_of_mc_steps number of steps in a Monte Carlo optimization
125 \param[in] number_of_cg_steps number of Conjugate Gradients steps in
126  a Monte Carlo step
127 \param[in] max_translation maximum translation step in a MC optimization step
128 \param[in] max_rotation maximum rotation step in a single MC optimization step
129 \param[in] fast if true the density map of the rigid body is not resampled
130  but transformed at each iteration of the optimization
131 \return the refined fitting solutions
132 */
134  Particle *p, Refiner *refiner, const FloatKey &weight_key,
135  DensityMap *dmap, const algebra::Vector3D &anchor_centroid,
136  OptimizerStates display_log, Int number_of_optimization_runs = 5,
137  Int number_of_mc_steps = 10, Int number_of_cg_steps = 100,
138  Float max_translation = 2., Float max_rotation = .3, bool fast = false);
139 
140 //! Local rigid fitting of a rigid body
141 /** Fit a set of particles to a density map around their centroid.
142  The fitting is assessed using the cross-correlation score.
143  The optimization is a standard MC/CG procedure.
144  The function returns a list of solutions sorted by the cross-correlation
145  score.
146 \note The returned cross-correlation score is 1-cc, as we usually want to
147  minimize a scoring function. Thus a score of 1 means no correlation
148  and a score of 0 is perfect correlation.
149 \note The input rigid body should also be IMP::atom::Hierarchy
150 \param[in] p The root of the hierarchy to fit
151 \param[in] refiner The refiner to get the leaves of the particle
152 \param[in] weight_key The weight key of the particles in the rigid body
153 \param[in] dmap The density map to fit to
154 \param[in] display_log If provided, then intermediate states
155  during the optimization procedure are printed
156 \param[in] number_of_optimization_runs number of Monte Carlo optimizations
157 \param[in] number_of_mc_steps number of steps in a Monte Carlo optimization
158 \param[in] number_of_cg_steps number of Conjugate Gradients steps in
159  a Monte Carlo step
160 \param[in] max_translation maximum translation step in a MC optimization step
161 \param[in] max_rotation maximum rotation step in radians in a single
162  MC optimization step
163 \param[in] fast if true the density map of the rigid body is not resampled
164  but transformed at each iteration of the optimization
165 \return the refined fitting solutions
166 */
167 
169  Particle *p, Refiner *refiner, const FloatKey &weight_key,
170  DensityMap *dmap, OptimizerStates display_log,
171  Int number_of_optimization_runs = 5, Int number_of_mc_steps = 10,
172  Int number_of_cg_steps = 100, Float max_translation = 2.,
173  Float max_rotation = .3, bool fast = true) {
174  IMP_LOG_VERBOSE("Start: local_rigid_fitting\n");
175  algebra::Vector3D rb_cen =
177  IMP_LOG_VERBOSE("centroid is:" << rb_cen << "\n");
179  p, refiner, weight_key, dmap, rb_cen, display_log,
180  number_of_optimization_runs, number_of_mc_steps, number_of_cg_steps,
181  max_translation, max_rotation, fast);
182 }
183 
184 //! Local rigid fitting of a rigid body around a set of center points
185 /**
186 \brief Fit a set of particles to a density map around each of the input points.
187  The function apply local_rigid_fitting_around_point around each center.
188 \note The input rigid body should also be IMP::atom::Hierarchy
189 \param[in] p The rigid body to fit
190 \param[in] refiner Refiner to yield rigid body members
191 \param[in] wei_key The weight key of the particles in the rigid body
192 \param[in] dmap The density map to fit to
193 \param[in] anchor_centroids The points to fit the particles around
194 \param[in] display_log If provided, then intermediate states
195  in during the optimization procedure are printed
196 \param[in] number_of_optimization_runs number of Monte Carlo optimizations
197 \param[in] number_of_mc_steps number of steps in a Monte Carlo optimization
198 \param[in] number_of_cg_steps number of Conjugate Gradients steps in
199  a Monte Carlo step
200 \param[in] max_translation maximum translation step in a MC optimization step
201 \param[in] max_rotation maximum rotation step in a single MC optimization step
202 \return the refined fitting solutions
203 */
204 IMPEMEXPORT FittingSolutions local_rigid_fitting_around_points(
205  Particle *p, Refiner *refiner, const FloatKey &wei_key,
206  DensityMap *dmap, const algebra::Vector3Ds &anchor_centroids,
207  OptimizerStates display_log, Int number_of_optimization_runs = 5,
208  Int number_of_mc_steps = 10, Int number_of_cg_steps = 100,
209  Float max_translation = 2., Float max_rotation = .3);
210 
211 //! Local grid search rigid fitting
212 /**
213 \brief Fit a set of particles to a density map around their centroid.
214  The fitting is assessed using the cross-correlation score.
215  The optimization is a grid search
216 \note The transformations are not clustered.
217 \note The returned cross-correlation score is 1-cc, as we usually want to
218  minimize a scoring function. Thus a score of 1 means no-correlation
219  and a score of 0. is perfect correlation.
220 \param[in] ps The particles to be fitted (treated rigid)
221 \param[in] wei_key The weight key of the particles in the rigid body
222 \param[in] dmap The density map to fit to
223 \param[in] max_voxels_translation Sample translations within
224  -max_voxel_translation to max_voxel_translation
225 \param[in] translation_step The translation sampling step
226 \param[in] max_angle_in_radians Sample rotations with +- max_angle_in_radians
227  around the current rotation
228 \param[in] number_of_rotations The number of rotations to sample
229 \return the refined fitting solutions
230 */
231 IMPEMEXPORT FittingSolutions local_rigid_fitting_grid_search(
232  const ParticlesTemp &ps, const FloatKey &wei_key, DensityMap *dmap,
233  Int max_voxels_translation = 2, Int translation_step = 1,
234  Float max_angle_in_radians = 0.174, Int number_of_rotations = 100);
235 
236 //! Compute fitting scores for a given set of rigid transformations
237 /**
238 \brief Score how well a set of particles fit to the map in various
239  rigid transformations.
240 \param[in] ps The particles to be fitted (treated rigid)
241 \param[in] em_map The density map to fit to
242 \param[in] wei_key The weight key of the particles in the rigid body
243 \param[in] fast_version If fast_version is used the sampled density map of the
244  input particles (ps) is not resampled for each
245  transformation but instead the corresponding grid
246  is rotated. This option significantly improves the
247  running times but the returned scores are less accurate
248 \param[in] transformations A set of rigid transformations
249 \param[in] fast_version if true the density map of each transformation
250  is interpolated
251 \param[in] local_score if true a local cross correlation score is used
252 \return The scored fitting solutions
253 \note the function assumes the density map holds its density
254  */
255 IMPEMEXPORT FittingSolutions
256  compute_fitting_scores(const ParticlesTemp &ps, DensityMap *em_map,
257  const algebra::Transformation3Ds &transformations,
258  bool fast_version = false, bool local_score = false,
259  const FloatKey &wei_key =
260  atom::Mass::get_mass_key());
261 
262 //! Compute fitting scores for a given set of rigid transformations
263 /**
264 \brief Score how well a rigid body fits to the map
265 \param[in] em_map The density map to fit to
266 \param[in] rb The rigid body
267 \param[in] refiner The rigid body refiner
268 \param[in] transformations Transformations of the rigid body
269 \return The scored fitting solutions
270 \note the function assumes the density map holds its density
271  */
273  DensityMap *em_map, core::RigidBody rb, Refiner *refiner,
274  const algebra::Transformation3Ds &transformations) {
275  return compute_fitting_scores(refiner->get_refined(rb), em_map,
276  transformations, true);
277 }
278 
279 //! Compute fitting scores for a given set of rigid transformations
280 /**
281 \brief Scores how well a set of particles fit a map
282 \param[in] ps The particles to be fitted
283 \param[in] em_map The density map to fit to
284 \param[in] wei_key The weight key of the particles in the rigid body
285 \note the function assumes the density map holds its density
286  */
287 IMPEMEXPORT Float
288  compute_fitting_score(const ParticlesTemp &ps, DensityMap *em_map,
289  FloatKey wei_key = atom::Mass::get_mass_key());
290 
291 IMPEM_END_NAMESPACE
292 #endif /* IMPEM_RIGID_FITTING_H */
Simple 3D transformation class.
Simple conjugate gradients optimizer.
Calculate score based on fit to EM map.
Used to hold a set of related restraints.
algebra::Transformation3D get_transformation(unsigned int i) const
Get the score of the i'th solution.
Definition: rigid_fitting.h:51
Simple Monte Carlo optimizer.
Storage of a model, its restraints, constraints and particles.
A simple list of fitting solutions.
Definition: rigid_fitting.h:33
#define IMP_LOG_VERBOSE(expr)
Definition: log_macros.h:83
algebra::Vector3D get_centroid(const XYZs &ps)
Get the centroid.
Class for handling density maps.
Class for handling density maps.
Definition: DensityMap.h:93
#define IMP_VALUES(Name, PluralName)
Define the type for storing sets of values.
Definition: value_macros.h:23
functionality for defining rigid bodies
Float compute_fitting_score(const ParticlesTemp &ps, DensityMap *em_map, FloatKey wei_key=atom::Mass::get_mass_key())
Compute fitting scores for a given set of rigid transformations.
FittingSolutions compute_fitting_scores(DensityMap *em_map, core::RigidBody rb, Refiner *refiner, const algebra::Transformation3Ds &transformations)
Compute fitting scores for a given set of rigid transformations.
FittingSolutions local_rigid_fitting_around_point(Particle *p, Refiner *refiner, const FloatKey &weight_key, DensityMap *dmap, const algebra::Vector3D &anchor_centroid, OptimizerStates display_log, Int number_of_optimization_runs=5, Int number_of_mc_steps=10, Int number_of_cg_steps=100, Float max_translation=2., Float max_rotation=.3, bool fast=false)
Local rigid fitting of a rigid body around a center point.
Shared score state.
virtual const ParticlesTemp get_refined(Particle *a) const =0
Refine the passed particle into a set of particles.
Simple 3D transformation class.
Classes to handle individual model particles. (Note that implementation of inline functions is in int...
std::ostream & show(Hierarchy h, std::ostream &out=std::cout)
Print the hierarchy using a given decorator to display each node.
FittingSolutions local_rigid_fitting(Particle *p, Refiner *refiner, const FloatKey &weight_key, DensityMap *dmap, OptimizerStates display_log, Int number_of_optimization_runs=5, Int number_of_mc_steps=10, Int number_of_cg_steps=100, Float max_translation=2., Float max_rotation=.3, bool fast=true)
Local rigid fitting of a rigid body.
VectorD< 3 > Vector3D
Definition: VectorD.h:421
double Float
Basic floating-point value (could be float, double...)
Definition: types.h:20
Abstract class to implement hierarchical methods.
Definition: Refiner.h:34
void multiply(const algebra::Transformation3D &t)
Definition: rigid_fitting.h:83
Class to handle individual particles of a Model object.
Definition: Particle.h:41
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
Definition: check_macros.h:168
int Int
Basic integer value.
Definition: types.h:35
Shared optimizer state.
FittingSolutions local_rigid_fitting_around_points(Particle *p, Refiner *refiner, const FloatKey &wei_key, DensityMap *dmap, const algebra::Vector3Ds &anchor_centroids, OptimizerStates display_log, Int number_of_optimization_runs=5, Int number_of_mc_steps=10, Int number_of_cg_steps=100, Float max_translation=2., Float max_rotation=.3)
Local rigid fitting of a rigid body around a set of center points.
FittingSolutions local_rigid_fitting_grid_search(const ParticlesTemp &ps, const FloatKey &wei_key, DensityMap *dmap, Int max_voxels_translation=2, Int translation_step=1, Float max_angle_in_radians=0.174, Int number_of_rotations=100)
Local grid search rigid fitting.
A decorator for a rigid body.
Definition: rigid_bodies.h:82
int get_number_of_solutions() const
Get the number of solutions in the set.
Definition: rigid_fitting.h:45
Float get_score(unsigned int i) const
Get the score of the i'th solution.
Definition: rigid_fitting.h:63