IsoSpec  2.1.3
fixedEnvelopes.h
1 /*
2  * Copyright (C) 2015-2020 Mateusz Łącki and Michał Startek.
3  *
4  * This file is part of IsoSpec.
5  *
6  * IsoSpec is free software: you can redistribute it and/or modify
7  * it under the terms of the Simplified ("2-clause") BSD licence.
8  *
9  * IsoSpec is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
12  *
13  * You should have received a copy of the Simplified BSD Licence
14  * along with IsoSpec. If not, see <https://opensource.org/licenses/BSD-2-Clause>.
15  */
16 
17 #pragma once
18 
19 #include <cstdlib>
20 #include <algorithm>
21 #include <vector>
22 #include <utility>
23 
24 #include "isoSpec++.h"
25 
26 #ifdef DEBUG
27 #define ISOSPEC_INIT_TABLE_SIZE 16
28 #else
29 #define ISOSPEC_INIT_TABLE_SIZE 1024
30 #endif
31 
32 namespace IsoSpec
33 {
34 
35 class ISOSPEC_EXPORT_SYMBOL FixedEnvelope {
36  protected:
37  double* _masses;
38  double* _probs;
39  int* _confs;
40  size_t _confs_no;
41  int allDim;
42  bool sorted_by_mass;
43  bool sorted_by_prob;
44  double total_prob;
45  size_t current_size;
46  double* tmasses;
47  double* tprobs;
48  int* tconfs;
49  int allDimSizeofInt;
50 
51  public:
52  ISOSPEC_FORCE_INLINE FixedEnvelope() : _masses(nullptr),
53  _probs(nullptr),
54  _confs(nullptr),
55  _confs_no(0),
56  allDim(0),
57  sorted_by_mass(false),
58  sorted_by_prob(false),
59  total_prob(0.0),
60  current_size(0),
61  allDimSizeofInt(0)
62  // Deliberately not initializing tmasses, tprobs, tconfs
63  {};
64 
65  FixedEnvelope(const FixedEnvelope& other);
67 
68  FixedEnvelope(double* masses, double* probs, size_t confs_no, bool masses_sorted = false, bool probs_sorted = false, double _total_prob = NAN);
69 
70  virtual ~FixedEnvelope()
71  {
72  free(_masses);
73  free(_probs);
74  free(_confs);
75  };
76 
77  FixedEnvelope operator+(const FixedEnvelope& other) const;
78  FixedEnvelope operator*(const FixedEnvelope& other) const;
79 
80  inline size_t confs_no() const { return _confs_no; }
81  inline int getAllDim() const { return allDim; }
82 
83  inline const double* masses() const { return _masses; }
84  inline const double* probs() const { return _probs; }
85  inline const int* confs() const { return _confs; }
86 
87  inline double* release_masses() { double* ret = _masses; _masses = nullptr; return ret; }
88  inline double* release_probs() { double* ret = _probs; _probs = nullptr; return ret; }
89  inline int* release_confs() { int* ret = _confs; _confs = nullptr; return ret; }
90 
91 
92  inline double mass(size_t i) const { return _masses[i]; }
93  inline double prob(size_t i) const { return _probs[i]; }
94  inline const int* conf(size_t i) const { return _confs + i*allDim; }
95 
96  void sort_by_mass();
97  void sort_by_prob();
98 
99  double get_total_prob();
100  void scale(double factor);
101  void normalize();
102 
103  double empiric_average_mass();
104  double empiric_variance();
105  double empiric_stddev() { return sqrt(empiric_variance()); }
106 
107  double WassersteinDistance(FixedEnvelope& other);
108  double OrientedWassersteinDistance(FixedEnvelope& other);
109 
110  static FixedEnvelope LinearCombination(const std::vector<const FixedEnvelope*>& spectra, const std::vector<double>& intensities);
111  static FixedEnvelope LinearCombination(const FixedEnvelope* const * spectra, const double* intensities, size_t size);
112 
113 
114  FixedEnvelope bin(double bin_width = 1.0, double middle = 0.0);
115 
116  private:
117  void sort_by(double* order);
118 
119 
120  protected:
121  template<typename T, bool tgetConfs> ISOSPEC_FORCE_INLINE void store_conf(const T& generator)
122  {
123  *tmasses = generator.mass(); tmasses++;
124  *tprobs = generator.prob(); tprobs++;
125  constexpr_if(tgetConfs) { generator.get_conf_signature(tconfs); tconfs += allDim; }
126  }
127 
128  ISOSPEC_FORCE_INLINE void store_conf(double _mass, double _prob)
129  {
130  if(_confs_no == current_size)
131  reallocate_memory<false>(current_size*2);
132 
133  *tprobs = _prob;
134  *tmasses = _mass;
135  tprobs++;
136  tmasses++;
137 
138  _confs_no++;
139  }
140 
141  template<bool tgetConfs> ISOSPEC_FORCE_INLINE void swap(size_t idx1, size_t idx2, ISOSPEC_MAYBE_UNUSED int* conf_swapspace)
142  {
143  std::swap<double>(this->_probs[idx1], this->_probs[idx2]);
144  std::swap<double>(this->_masses[idx1], this->_masses[idx2]);
145  constexpr_if(tgetConfs)
146  {
147  int* c1 = this->_confs + (idx1*this->allDim);
148  int* c2 = this->_confs + (idx2*this->allDim);
149  memcpy(conf_swapspace, c1, this->allDimSizeofInt);
150  memcpy(c1, c2, this->allDimSizeofInt);
151  memcpy(c2, conf_swapspace, this->allDimSizeofInt);
152  }
153  }
154 
155  template<bool tgetConfs> void reallocate_memory(size_t new_size);
156  void slow_reallocate_memory(size_t new_size);
157 
158  public:
159  template<bool tgetConfs> void threshold_init(Iso&& iso, double threshold, bool absolute);
160 
161  template<bool tgetConfs, typename GenType = IsoLayeredGenerator> void addConfILG(const GenType& generator)
162  {
163  if(this->_confs_no == this->current_size)
164  this->template reallocate_memory<tgetConfs>(this->current_size*2);
165 
166  this->template store_conf<GenType, tgetConfs>(generator);
167  this->_confs_no++;
168  }
169 
170  template<bool tgetConfs> void total_prob_init(Iso&& iso, double target_prob, bool trim);
171 
172  static FixedEnvelope FromThreshold(Iso&& iso, double threshold, bool absolute, bool tgetConfs = false)
173  {
174  FixedEnvelope ret;
175 
176  if(tgetConfs)
177  ret.threshold_init<true>(std::move(iso), threshold, absolute);
178  else
179  ret.threshold_init<false>(std::move(iso), threshold, absolute);
180  return ret;
181  }
182 
183  inline static FixedEnvelope FromThreshold(const Iso& iso, double _threshold, bool _absolute, bool tgetConfs = false)
184  {
185  return FromThreshold(Iso(iso, false), _threshold, _absolute, tgetConfs);
186  }
187 
188  static FixedEnvelope FromTotalProb(Iso&& iso, double target_total_prob, bool optimize, bool tgetConfs = false)
189  {
190  FixedEnvelope ret;
191 
192  if(tgetConfs)
193  ret.total_prob_init<true>(std::move(iso), target_total_prob, optimize);
194  else
195  ret.total_prob_init<false>(std::move(iso), target_total_prob, optimize);
196 
197  return ret;
198  }
199 
200  inline static FixedEnvelope FromTotalProb(const Iso& iso, double _target_total_prob, bool _optimize, bool tgetConfs = false)
201  {
202  return FromTotalProb(Iso(iso, false), _target_total_prob, _optimize, tgetConfs);
203  }
204 
205  template<bool tgetConfs> void stochastic_init(Iso&& iso, size_t _no_molecules, double _precision, double _beta_bias);
206 
207  inline static FixedEnvelope FromStochastic(Iso&& iso, size_t _no_molecules, double _precision = 0.9999, double _beta_bias = 5.0, bool tgetConfs = false)
208  {
209  FixedEnvelope ret;
210 
211  if(tgetConfs)
212  ret.stochastic_init<true>(std::move(iso), _no_molecules, _precision, _beta_bias);
213  else
214  ret.stochastic_init<false>(std::move(iso), _no_molecules, _precision, _beta_bias);
215 
216  return ret;
217  }
218 
219  static FixedEnvelope FromStochastic(const Iso& iso, size_t _no_molecules, double _precision = 0.9999, double _beta_bias = 5.0, bool tgetConfs = false)
220  {
221  return FromStochastic(Iso(iso, false), _no_molecules, _precision, _beta_bias, tgetConfs);
222  }
223 
224  static FixedEnvelope Binned(Iso&& iso, double target_total_prob, double bin_width, double bin_middle = 0.0);
225  static FixedEnvelope Binned(const Iso& iso, double target_total_prob, double bin_width, double bin_middle = 0.0)
226  {
227  return Binned(Iso(iso, false), target_total_prob, bin_width, bin_middle);
228  }
229 };
230 
231 } // namespace IsoSpec
The Iso class for the calculation of the isotopic distribution.
Definition: isoSpec++.h:49