IMP  2.3.0
The Integrative Modeling Platform
Profile.h
Go to the documentation of this file.
1 /**
2  * \file IMP/saxs/Profile.h
3  * \brief A class for profile storing and computation
4  *
5  * Copyright 2007-2014 IMP Inventors. All rights reserved.
6  *
7  */
8 
9 #ifndef IMPSAXS_PROFILE_H
10 #define IMPSAXS_PROFILE_H
11 
12 #include <IMP/saxs/saxs_config.h>
13 #include <IMP/base/Object.h>
15 
16 #include "FormFactorTable.h"
17 #include "Distribution.h"
18 
19 #include <iostream>
20 #include <vector>
21 
22 IMPSAXS_BEGIN_NAMESPACE
23 
24 class RadialDistributionFunction;
25 
26 /**
27  Basic profile class, can be initialized from the input file
28  (experimental or theoretical) or computed from a set of kernel::Model
29  kernel::Particles (theoretical)
30 */
31 class IMPSAXSEXPORT Profile : public base::Object {
32  public:
33  // Constructors
34 
35  //! init from file
36  Profile(const String& file_name, bool fit_file = false);
37 
38  //! init for theoretical profile
39  Profile(Float qmin = 0.0, Float qmax = 0.5, Float delta = 0.005);
40 
41  // Various ways to compute a profile
42 
43  //! computes theoretical profile
44  void calculate_profile(const kernel::Particles& particles,
45  FormFactorType ff_type = HEAVY_ATOMS,
46  bool reciprocal = false) {
47  if (!reciprocal)
48  calculate_profile_real(particles, ff_type);
49  else
50  calculate_profile_reciprocal(particles, ff_type);
51  }
52 
53  //! compute profile for fitting with hydration layer and excluded volume
54  /**
55  A partial profile is a pre-computed profile, where intensities are
56  split into 6 parts that can be summed up into a regular profile
57  given a pair of c1/c2 values by sum_partial_profiles function.
58  see FoXS paper for details.
59  */
60  void calculate_profile_partial(const kernel::Particles& particles,
61  const Floats& surface = Floats(),
62  FormFactorType ff_type = HEAVY_ATOMS);
63 
64  //! compute profile for fitting with hydration layer and excluded volume
65  void calculate_profile_partial(const kernel::Particles& particles1,
66  const kernel::Particles& particles2,
67  const Floats& surface1 = Floats(),
68  const Floats& surface2 = Floats(),
69  FormFactorType ff_type = HEAVY_ATOMS);
70 
71  void calculate_profile_reciprocal_partial(const kernel::Particles& particles,
72  const Floats& surface = Floats(),
73  FormFactorType ff_type =
74  HEAVY_ATOMS);
75 
76  //! computes theoretical profile contribution from inter-molecular
77  //! interactions between the particles
78  void calculate_profile(const kernel::Particles& particles1,
79  const kernel::Particles& particles2,
80  FormFactorType ff_type = HEAVY_ATOMS) {
81  calculate_profile_real(particles1, particles2, ff_type);
82  }
83 
84  //! calculate Intensity at zero (= squared number of electrons)
85  Float calculate_I0(const kernel::Particles& particles,
86  FormFactorType ff_type = HEAVY_ATOMS);
87 
88  //! calculate profile for any type of kernel::Particles that have coordinates
89  void calculate_profile_constant_form_factor(
90  const kernel::Particles& particles, Float form_factor = 1.0);
91 
92  // computes theoretical profile faster for cyclically symmetric particles
93  // assumes that the units particles are ordered one after another in the
94  // input particles vector (n - symmetry order)
95  void calculate_profile_symmetric(const kernel::Particles& particles,
96  unsigned int n,
97  FormFactorType ff_type = HEAVY_ATOMS);
98 
99  //! convert to real space P(r) function P(r) = 1/2PI^2 Sum(I(q)*qr*sin(qr))
100  void profile_2_distribution(RadialDistributionFunction& rd,
101  Float max_distance) const;
102 
103  //! convert to reciprocal space I(q) = Sum(P(r)*sin(qr)/qr)
104  void distribution_2_profile(const RadialDistributionFunction& r_dist);
105 
106  //! return a profile that is sampled on the q values of the exp_profile
107  void resample(const Profile* exp_profile, Profile* resampled_profile,
108  bool partial_profiles = false) const;
109 
110  //! downsample the profile to a given number of points
111  void downsample(Profile* downsampled_profile,
112  unsigned int point_number) const;
113 
114  //! compute radius of gyration with Guinier approximation
115  /** ln[I(q)]=ln[I(0)] - (q^2*rg^2)/3
116  \param[in] end_q_rg determines the range of profile used for approximation:
117  i.e. q*rg < end_q_rg. Use 1.3 for globular proteins, 0.8 for elongated
118  */
119  double radius_of_gyration(double end_q_rg = 1.3) const;
120 
121  // IO functions
122 
123  //! reads SAXS profile from file
124  /**
125  \param[in] file_name profile file name
126  \param[in] fit_file if true, intensities are read from column 3
127  */
128  void read_SAXS_file(const String& file_name, bool fit_file = false);
129 
130  //! print to file
131  /** \param[in] file_name output file name
132  \param[in] max_q output till maximal q value = max_q, or all if max_q<=0
133  */
134  void write_SAXS_file(const String& file_name, Float max_q = 0.0) const;
135 
136  //! read a partial profile from file (7 columns)
137  void read_partial_profiles(const String& file_name);
138 
139  //! write a partial profile to file
140  void write_partial_profiles(const String& file_name) const;
141 
142  // Access functions
143 
144  //! return sampling resolution
145  Float get_delta_q() const { return delta_q_; }
146 
147  //! return minimal sampling point
148  Float get_min_q() const { return min_q_; }
149 
150  //! return maximal sampling point
151  Float get_max_q() const { return max_q_; }
152 
153  Float get_intensity(unsigned int i) const { return intensity_[i]; }
154  Float get_q(unsigned int i) const { return q_[i]; }
155  Float get_error(unsigned int i) const { return error_[i]; }
156  Float get_weight(unsigned int i) const {
157  IMP_UNUSED(i);
158  return 1.0;
159  }
160  Float get_average_radius() const { return average_radius_; }
161 
162  //! return number of entries in SAXS profile
163  unsigned int size() const { return q_.size(); }
164 
165  //! checks the sampling of experimental profile
166  bool is_uniform_sampling() const;
167 
168  std::string get_name() const { return name_; }
169 
170  unsigned int get_id() const { return id_; }
171 
172  // Modifiers
173 
174  void set_intensity(unsigned int i, Float iq) { intensity_[i] = iq; }
175 
176  //! required for reciprocal space calculation
177  void set_ff_table(FormFactorTable* ff_table) { ff_table_ = ff_table; }
178 
179  void set_average_radius(Float r) { average_radius_ = r; }
180 
181  void set_average_volume(Float v) { average_volume_ = v; }
182 
183  void set_name(std::string name) { name_ = name; }
184 
185  void set_id(unsigned int id) { id_ = id; }
186 
187  void set_beam_profile(std::string beam_profile_file) {
188  beam_profile_ = new Profile(beam_profile_file);
189  }
190 
191  //! add intensity entry to profile
192  void add_entry(Float q, Float intensity, Float error = 1.0) {
193  q_.push_back(q);
194  intensity_.push_back(intensity);
195  error_.push_back(error);
196  }
197 
198  //! add simulated error
199  void add_errors();
200 
201  //! add simulated noise
202  void add_noise(Float percentage = 0.03);
203 
204  //! computes full profile for given fitting parameters
205  void sum_partial_profiles(Float c1, Float c2, bool check_cashed = true);
206 
207  //! add another profile - useful for rigid bodies
208  void add(const Profile* other_profile, Float weight = 1.0);
209 
210  //! add partial profiles
211  void add_partial_profiles(const Profile* other_profile, Float weight = 1.0);
212 
213  //! add other profiles - useful for weighted ensembles
214  void add(const std::vector<Profile*>& profiles,
215  const std::vector<Float>& weights = std::vector<Float>());
216 
217  //! add other partial profiles
218  void add_partial_profiles(const std::vector<Profile*>& profiles,
219  const std::vector<Float>& weights =
220  std::vector<Float>());
221 
222  //! background adjustment option
223  void background_adjust(double start_q);
224 
225  //! scale
226  void scale(Float c);
227 
228  //! offset profile by c, I(q) = I(q) - c
229  void offset(Float c);
230 
231  // copy error bars from the matching experimental profile
232  void copy_errors(const Profile* exp_profile);
233 
234  // parameter for E^2(q), used in faster calculation
235  static const Float modulation_function_parameter_;
236 
238 
239  protected:
240  void init();
241 
242  private:
243  void calculate_profile_reciprocal(const kernel::Particles& particles,
244  FormFactorType ff_type = HEAVY_ATOMS);
245 
246  void calculate_profile_reciprocal(const kernel::Particles& particles1,
247  const kernel::Particles& particles2,
248  FormFactorType ff_type = HEAVY_ATOMS);
249 
250  void calculate_profile_real(const kernel::Particles& particles,
251  FormFactorType ff_type = HEAVY_ATOMS);
252 
253  void calculate_profile_real(const kernel::Particles& particles1,
254  const kernel::Particles& particles2,
255  FormFactorType ff_type = HEAVY_ATOMS);
256 
257  void squared_distribution_2_profile(const RadialDistributionFunction& r_dist);
258 
259  void squared_distributions_2_partial_profiles(
260  const std::vector<RadialDistributionFunction>& r_dist);
261 
262  double radius_of_gyration_fixed_q(double end_q) const;
263 
264  protected:
265  std::vector<double> q_; // q sampling points
266  std::vector<double> intensity_;
267  std::vector<double> error_; // error bar of each point
268 
269  Float min_q_, max_q_; // minimal and maximal s values in the profile
270  Float delta_q_; // profile sampling resolution
271  FormFactorTable* ff_table_; // pointer to form factors table
272 
273  // stores the intensity split into 6 for c1/c2 enumeration
274  std::vector<std::vector<double> > partial_profiles_;
275  Float c1_, c2_;
276 
277  bool experimental_; // experimental profile read from file
278  Float average_radius_; // average radius of the particles
279  Float average_volume_; // average volume
280 
281  // mapping from q values to vector index for fast profile resampling
282  std::map<float, unsigned int> q_mapping_;
283 
284  std::string name_; // file name
285  unsigned int id_; // identifier
286 
287  Profile* beam_profile_;
288 };
289 
291 
292 IMPSAXS_END_NAMESPACE
293 
294 #endif /* IMPSAXS_PROFILE_H */
unsigned int size() const
return number of entries in SAXS profile
Definition: Profile.h:163
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Definition: object_macros.h:25
void calculate_profile(const kernel::Particles &particles1, const kernel::Particles &particles2, FormFactorType ff_type=HEAVY_ATOMS)
Definition: Profile.h:78
Float get_delta_q() const
return sampling resolution
Definition: Profile.h:145
Float get_min_q() const
return minimal sampling point
Definition: Profile.h:148
Various general useful macros for IMP.
#define IMP_UNUSED(variable)
A class for computation of atomic and residue level form factors for SAXS calculations.
void calculate_profile(const kernel::Particles &particles, FormFactorType ff_type=HEAVY_ATOMS, bool reciprocal=false)
computes theoretical profile
Definition: Profile.h:44
Float get_max_q() const
return maximal sampling point
Definition: Profile.h:151
computes distribution functions
Common base class for heavy weight IMP objects.
Definition: Object.h:106
#define IMP_OBJECTS(Name, PluralName)
Define the types for storing sets of objects.
Definition: object_macros.h:52
A shared base class to help in debugging and things.
Float radius_of_gyration(const kernel::Particles &particles)
compute radius_of_gyration
Definition: saxs/utility.h:73
IMP::base::Vector< Float > Floats
Standard way to pass a bunch of Float values.
Definition: types.h:47
double Float
Basic floating-point value (could be float, double...)
Definition: types.h:20
FormFactorType
type of the form factors for profile calculations
void set_ff_table(FormFactorTable *ff_table)
required for reciprocal space calculation
Definition: Profile.h:177
std::string String
Basic string value.
Definition: types.h:44
void add_entry(Float q, Float intensity, Float error=1.0)
add intensity entry to profile
Definition: Profile.h:192