OpenMS
BilinearInterpolation.h
Go to the documentation of this file.
1 // Copyright (c) 2002-2023, The OpenMS Team -- EKU Tuebingen, ETH Zurich, and FU Berlin
2 // SPDX-License-Identifier: BSD-3-Clause
3 //
4 // --------------------------------------------------------------------------
5 // $Maintainer: Timo Sachsenberg $
6 // $Authors: $
7 // --------------------------------------------------------------------------
8 
9 #pragma once
10 
12 
13 namespace OpenMS
14 {
15 
16  namespace Math
17  {
18 
45  template <typename Key = double, typename Value = Key>
47  {
48 
49 public:
50 
52 
53  typedef Value value_type;
54 
55  typedef Key key_type;
57 
59  typedef key_type KeyType;
62 
63 public:
64 
68 
71  scale_0_(1),
72  offset_0_(0),
73  scale_1_(1),
74  offset_1_(0),
75  inside_0_(0),
76  outside_0_(0),
77  inside_1_(0),
78  outside_1_(0),
79  data_()
80  {}
81 
84  scale_0_(arg.scale_0_),
85  offset_0_(arg.offset_0_),
86  scale_1_(arg.scale_1_),
87  offset_1_(arg.offset_1_),
88  inside_0_(arg.inside_0_),
90  inside_1_(arg.inside_1_),
92  data_(arg.data_)
93  {}
94 
97  {
98  if (&arg == this)
99  return *this;
100 
101  scale_0_ = arg.scale_0_;
102  offset_0_ = arg.offset_0_;
103  scale_1_ = arg.scale_1_;
104  offset_1_ = arg.offset_1_;
105  inside_0_ = arg.inside_0_;
106  outside_1_ = arg.outside_1_;
107  inside_1_ = arg.inside_1_;
108  outside_0_ = arg.outside_0_;
109  data_ = arg.data_;
110  return *this;
111  }
112 
115  {}
116 
118 
119  // ----------------------------------------------------------------------
120 
122 
123 
125  ValueType value(KeyType arg_pos_0, KeyType arg_pos_1) const
126  {
127  // apply the key transformations
128  KeyType const pos_0 = key2index_0(arg_pos_0);
129  KeyType const pos_1 = key2index_1(arg_pos_1);
130 
131  // ???? should use modf() here!
132 
133  SignedSize const size_0 = data_.rows();
134  SignedSize const lower_0 = SignedSize(pos_0); // this rounds towards zero
135  SignedSize const size_1 = data_.cols();
136  SignedSize const lower_1 = SignedSize(pos_1); // this rounds towards zero
137 
138  // small pos_0
139  if (pos_0 <= 0)
140  {
141  if (lower_0 != 0)
142  {
143  return 0;
144  }
145  else // that is: -1 < pos_0 <= 0
146  { // small pos_1
147  if (pos_1 <= 0)
148  {
149  if (lower_1 != 0)
150  {
151  return 0;
152  }
153  else // that is: -1 < pos_1 <= 0
154  {
155  return data_(0, 0) * (1. + pos_0) * (1. + pos_1);
156  }
157  }
158 
159  // big pos_1
160  if (lower_1 >= size_1 - 1)
161  {
162  if (lower_1 != size_1 - 1)
163  {
164  return 0;
165  }
166  else
167  {
168  return data_(0, lower_1) * (1. + pos_0) * (size_1 - pos_1);
169  }
170  }
171 
172  // medium pos_1
173  KeyType const factor_1 = pos_1 - KeyType(lower_1);
174  KeyType const factor_1_complement = KeyType(1.) - factor_1;
175  return (
176  data_(0, lower_1 + 1) * factor_1 +
177  data_(0, lower_1) * factor_1_complement
178  ) * (1. + pos_0);
179  }
180  }
181 
182  // big pos_0
183  if (lower_0 >= size_0 - 1)
184  {
185  if (lower_0 != size_0 - 1)
186  {
187  return 0;
188  }
189  else // that is: size_0 - 1 <= pos_0 < size_0
190  { // small pos_1
191  if (pos_1 <= 0)
192  {
193  if (lower_1 != 0)
194  {
195  return 0;
196  }
197  else // that is: -1 < pos_1 <= 0
198  {
199  return data_(lower_0, 0) * (size_0 - pos_0) * (1. + pos_1);
200  }
201  }
202 
203  // big pos_1
204  if (lower_1 >= size_1 - 1)
205  {
206  if (lower_1 != size_1 - 1)
207  {
208  return 0;
209  }
210  else
211  {
212  return data_(lower_0, lower_1) * (size_0 - pos_0) * (size_1 - pos_1);
213  }
214  }
215 
216  // medium pos_1
217  KeyType const factor_1 = pos_1 - KeyType(lower_1);
218  KeyType const factor_1_complement = KeyType(1.) - factor_1;
219  return (
220  data_(lower_0, lower_1 + 1) * factor_1 +
221  data_(lower_0, lower_1) * factor_1_complement
222  )
223  * (size_0 - pos_0);
224  }
225  }
226 
227  // medium pos_0
228  {
229  KeyType const factor_0 = pos_0 - KeyType(lower_0);
230  KeyType const factor_0_complement = KeyType(1.) - factor_0;
231 
232  // small pos_1
233  if (pos_1 <= 0)
234  {
235  if (lower_1 != 0)
236  {
237  return 0;
238  }
239  else // that is: -1 < pos_1 <= 0
240  {
241  return (
242  data_(lower_0 + 1, 0) * factor_0
243  +
244  data_(lower_0, 0) * factor_0_complement
245  )
246  * (1. + pos_1);
247  }
248  }
249 
250  // big pos_1
251  if (lower_1 >= size_1 - 1)
252  {
253  if (lower_1 != size_1 - 1)
254  {
255  return 0;
256  }
257  else
258  {
259  return (
260  data_(lower_0 + 1, lower_1) * factor_0
261  +
262  data_(lower_0, lower_1) * factor_0_complement
263  )
264  * (size_1 - pos_1);
265  }
266  }
267  KeyType const factor_1 = pos_1 - KeyType(lower_1);
268  KeyType const factor_1_complement = KeyType(1.) - factor_1;
269 
270  // medium pos_0 and medium pos_1 --> "within" the matrix
271  return (
272  data_(lower_0 + 1, lower_1 + 1) * factor_0
273  +
274  data_(lower_0, lower_1 + 1) * factor_0_complement
275  )
276  * factor_1
277  +
278  (
279  data_(lower_0 + 1, lower_1) * factor_0
280  +
281  data_(lower_0, lower_1) * factor_0_complement
282  )
283  * factor_1_complement;
284  }
285  }
286 
290  void addValue(KeyType arg_pos_0, KeyType arg_pos_1, ValueType arg_value)
291  {
292 
293  typedef typename container_type::difference_type DiffType;
294 
295  // apply key transformation _0
296  KeyType const pos_0 = key2index_0(arg_pos_0);
297  KeyType lower_0_key;
298  KeyType const frac_0 = std::modf(pos_0, &lower_0_key);
299  DiffType const lower_0 = DiffType(lower_0_key);
300 
301  // Small pos_0 ?
302  if (pos_0 < 0)
303  {
304  if (lower_0)
305  {
306  return;
307  }
308  else // lower_0 == 0
309  { // apply key transformation _1
310  KeyType const pos_1 = key2index_1(arg_pos_1);
311  KeyType lower_1_key;
312  KeyType const frac_1 = std::modf(pos_1, &lower_1_key);
313  DiffType const lower_1 = DiffType(lower_1_key);
314 
315  // Small pos_1 ?
316  if (pos_1 < 0)
317  {
318  if (lower_1)
319  {
320  return;
321  }
322  else // lower_1 == 0
323  {
324  data_(0, 0) += arg_value * (1 + frac_0) * (1 + frac_1);
325  return;
326  }
327  }
328  else // pos_1 >= 0
329  {
330  DiffType const back_1 = data_.cols() - 1;
331  // big pos_1
332  if (lower_1 >= back_1)
333  {
334  if (lower_1 != back_1)
335  {
336  return;
337  }
338  else // lower_1 == back_1
339  {
340  data_(0, lower_1) += arg_value * (1 + frac_0) * (1 - frac_1);
341  return;
342  }
343  }
344  else
345  {
346  // medium pos_1
347  KeyType const tmp_prod = KeyType(arg_value * (1. + frac_0));
348  data_(0, lower_1 + 1) += tmp_prod * frac_1;
349  data_(0, lower_1) += tmp_prod * (1. - frac_1);
350  return;
351  }
352  }
353  }
354  }
355  else // pos_0 >= 0
356  {
357  DiffType const back_0 = data_.rows() - 1;
358  if (lower_0 >= back_0)
359  {
360  if (lower_0 != back_0)
361  {
362  return;
363  }
364  else // lower_0 == back_0
365  {
366 
367  KeyType const tmp_prod = KeyType(arg_value * (1. - frac_0));
368 
369  // apply key transformation _1
370  KeyType const pos_1 = key2index_1(arg_pos_1);
371  KeyType lower_1_key;
372  KeyType const frac_1 = std::modf(pos_1, &lower_1_key);
373  DiffType const lower_1 = DiffType(lower_1_key);
374 
375  // Small pos_1 ?
376  if (pos_1 < 0)
377  {
378  if (lower_1)
379  {
380  return;
381  }
382  else // lower_1 == 0
383  {
384  data_(lower_0, 0) += tmp_prod * (1 + frac_1);
385  return;
386  }
387  }
388  else // pos_1 >= 0
389  {
390  DiffType const back_1 = data_.cols() - 1;
391  // big pos_1
392  if (lower_1 >= back_1)
393  {
394  if (lower_1 != back_1)
395  {
396  return;
397  }
398  else // lower_1 == back_1
399  {
400  data_(lower_0, lower_1) += tmp_prod * (1 - frac_1);
401  return;
402  }
403  }
404  else
405  {
406  // medium pos_1
407  data_(lower_0, lower_1 + 1) += tmp_prod * frac_1;
408  data_(lower_0, lower_1) += tmp_prod * (1 - frac_1);
409  return;
410  }
411  }
412  }
413  }
414  else // lower_0 < back_0
415  {
416 
417  // Medium pos_0 !
418 
419  // apply key transformation _1
420  KeyType const pos_1 = key2index_1(arg_pos_1);
421  KeyType lower_1_key;
422  KeyType const frac_1 = std::modf(pos_1, &lower_1_key);
423  DiffType const lower_1 = DiffType(lower_1_key);
424 
425  // Small pos_1 ?
426  if (pos_1 < 0)
427  {
428  if (lower_1)
429  {
430  return;
431  }
432  else // lower_1 == 0
433  {
434  KeyType const tmp_prod = KeyType(arg_value * (1 + frac_1));
435  data_(lower_0 + 1, 0) += tmp_prod * frac_0;
436  data_(lower_0, 0) += tmp_prod * (1 - frac_0);
437  return;
438  }
439  }
440  else // pos_1 >= 0
441  {
442  DiffType const back_1 = data_.cols() - 1;
443  // big pos_1
444  if (lower_1 >= back_1)
445  {
446  if (lower_1 != back_1)
447  {
448  return;
449  }
450  else // lower_1 == back_1
451  {
452  KeyType const tmp_prod = KeyType(arg_value * (1 - frac_1));
453  data_(lower_0 + 1, lower_1) += tmp_prod * frac_0;
454  data_(lower_0, lower_1) += tmp_prod * (1 - frac_0);
455  return;
456  }
457  }
458  else
459  {
460  // Medium pos_1 !
461 
462  // medium pos_0 and medium pos_1 --> "within" the matrix
463  KeyType tmp_prod = KeyType(arg_value * frac_0);
464  data_(lower_0 + 1, lower_1 + 1) += tmp_prod * frac_1;
465  data_(lower_0 + 1, lower_1) += tmp_prod * (1 - frac_1);
466  tmp_prod = KeyType(arg_value * (1 - frac_0));
467  data_(lower_0, lower_1 + 1) += tmp_prod * frac_1;
468  data_(lower_0, lower_1) += tmp_prod * (1 - frac_1);
469  return;
470  }
471  }
472  }
473  }
474  }
475 
477 
478  // ----------------------------------------------------------------------
479 
481 
482 
485  {
486  return data_;
487  }
488 
490  ContainerType const & getData() const
491  {
492  return data_;
493  }
494 
500  template <typename SourceContainer>
501  void setData(SourceContainer const & data)
502  {
503  data_ = data;
504  }
505 
507  bool empty() const
508  {
509  return data_.empty();
510  }
511 
513 
514  // ----------------------------------------------------------------------
515 
517 
518 
521  {
522  if (scale_0_)
523  {
524  pos -= offset_0_;
525  pos /= scale_0_;
526  return pos;
527  }
528  else
529  {
530  return 0;
531  }
532  }
533 
536  {
537  pos *= scale_0_;
538  pos += offset_0_;
539  return pos;
540  }
541 
544  {
545  if (scale_1_)
546  {
547  pos -= offset_1_;
548  pos /= scale_1_;
549  return pos;
550  }
551  else
552  {
553  return 0;
554  }
555  }
556 
559  {
560  pos *= scale_1_;
561  pos += offset_1_;
562  return pos;
563  }
564 
566  KeyType const & getScale_0() const
567  {
568  return scale_0_;
569  }
570 
572  KeyType const & getScale_1() const
573  {
574  return scale_1_;
575  }
576 
582  void setScale_0(KeyType const & scale)
583  {
584  scale_0_ = scale;
585  }
586 
592  void setScale_1(KeyType const & scale)
593  {
594  scale_1_ = scale;
595  }
596 
598  KeyType const & getOffset_0() const
599  {
600  return offset_0_;
601  }
602 
604  KeyType const & getOffset_1() const
605  {
606  return offset_1_;
607  }
608 
615  void setOffset_0(KeyType const & offset)
616  {
617  offset_0_ = offset;
618  }
619 
626  void setOffset_1(KeyType const & offset)
627  {
628  offset_1_ = offset;
629  }
630 
644  void setMapping_0(KeyType const & scale, KeyType const & inside_low, KeyType const & outside_low)
645  {
646  scale_0_ = scale;
647  inside_0_ = inside_low;
648  outside_0_ = outside_low;
649  offset_0_ = outside_low - scale * inside_low;
650  return;
651  }
652 
659  void setMapping_0(KeyType const & inside_low, KeyType const & outside_low,
660  KeyType const & inside_high, KeyType const & outside_high)
661  {
662  if (inside_high != inside_low)
663  {
664  setMapping_0((outside_high - outside_low) / (inside_high - inside_low),
665  inside_low, outside_low);
666  }
667  else
668  {
669  setMapping_0(0, inside_low, outside_low);
670  }
671  return;
672  }
673 
687  void setMapping_1(KeyType const & scale, KeyType const & inside_low, KeyType const & outside_low)
688  {
689  scale_1_ = scale;
690  inside_1_ = inside_low;
691  outside_1_ = outside_low;
692  offset_1_ = outside_low - scale * inside_low;
693  return;
694  }
695 
702  void setMapping_1(KeyType const & inside_low, KeyType const & outside_low,
703  KeyType const & inside_high, KeyType const & outside_high)
704  {
705  if (inside_high != inside_low)
706  {
707  setMapping_1((outside_high - outside_low) / (inside_high - inside_low),
708  inside_low, outside_low);
709  }
710  else
711  {
712  setMapping_1(0, inside_low, outside_low);
713  }
714  return;
715  }
716 
719  {
720  return inside_0_;
721  }
722 
725  {
726  return inside_1_;
727  }
728 
731  {
732  return outside_0_;
733  }
734 
737  {
738  return outside_1_;
739  }
740 
743  {
744  return index2key_0(empty() ? KeyType(0.) : KeyType(-1.));
745  }
746 
749  {
750  return index2key_1(empty() ? KeyType(0.) : KeyType(-1.));
751  }
752 
755  {
756  return index2key_0(KeyType(data_.rows()));
757  }
758 
761  {
762  return index2key_1(KeyType(data_.cols()));
763  }
764 
766 
767 protected:
768 
781  };
782 
783  } // namespace Math
784 
785 } // namespace OpenMS
786 
Provides access to bilinearly interpolated values (and derivatives) from discrete data points....
Definition: BilinearInterpolation.h:47
KeyType const & getOffset_1() const
Accessor. "Offset" is the point (in "outside" units) which corresponds to "Data(0,...
Definition: BilinearInterpolation.h:604
ContainerType const & getData() const
Returns the internal random access container storing the data.
Definition: BilinearInterpolation.h:490
Value value_type
Definition: BilinearInterpolation.h:53
KeyType supportMin_1() const
Lower boundary of the support, in "outside" coordinates.
Definition: BilinearInterpolation.h:748
ValueType value(KeyType arg_pos_0, KeyType arg_pos_1) const
Returns the interpolated value ("backward resampling")
Definition: BilinearInterpolation.h:125
KeyType const & getOffset_0() const
Accessor. "Offset" is the point (in "outside" units) which corresponds to "Data(0,...
Definition: BilinearInterpolation.h:598
Key key_type
Definition: BilinearInterpolation.h:55
Matrix< value_type > container_type
Definition: BilinearInterpolation.h:56
ContainerType data_
Definition: BilinearInterpolation.h:779
container_type ContainerType
Definition: BilinearInterpolation.h:60
void setData(SourceContainer const &data)
Assigns data to the internal random access container storing the data.
Definition: BilinearInterpolation.h:501
KeyType outside_0_
Definition: BilinearInterpolation.h:776
void setMapping_1(KeyType const &inside_low, KeyType const &outside_low, KeyType const &inside_high, KeyType const &outside_high)
Specifies the mapping from "outside" to "inside" coordinates by the following data:
Definition: BilinearInterpolation.h:702
KeyType supportMax_1() const
Upper boundary of the support, in "outside" coordinates.
Definition: BilinearInterpolation.h:760
BilinearInterpolation & operator=(BilinearInterpolation const &arg)
Assignment operator.
Definition: BilinearInterpolation.h:96
KeyType inside_1_
Definition: BilinearInterpolation.h:777
BilinearInterpolation(BilinearInterpolation const &arg)
Copy constructor.
Definition: BilinearInterpolation.h:83
KeyType const & getScale_1() const
Accessor. "Scale" is the difference (in "outside" units) between consecutive entries in "Data".
Definition: BilinearInterpolation.h:572
void setMapping_1(KeyType const &scale, KeyType const &inside_low, KeyType const &outside_low)
Specifies the mapping from "outside" to "inside" coordinates by the following data:
Definition: BilinearInterpolation.h:687
bool empty() const
Returns true if getData() is empty.
Definition: BilinearInterpolation.h:507
KeyType supportMax_0() const
Upper boundary of the support, in "outside" coordinates.
Definition: BilinearInterpolation.h:754
KeyType const & getOutsideReferencePoint_1() const
Accessor. See setMapping().
Definition: BilinearInterpolation.h:736
void setOffset_1(KeyType const &offset)
Accessor. "Offset" is the point (in "outside" units) which corresponds to "Data(0,...
Definition: BilinearInterpolation.h:626
KeyType const & getOutsideReferencePoint_0() const
Accessor. See setMapping().
Definition: BilinearInterpolation.h:730
void setMapping_0(KeyType const &scale, KeyType const &inside_low, KeyType const &outside_low)
Specifies the mapping from "outside" to "inside" coordinates by the following data:
Definition: BilinearInterpolation.h:644
KeyType const & getInsideReferencePoint_1() const
Accessor. See setMapping().
Definition: BilinearInterpolation.h:724
value_type ValueType
Definition: BilinearInterpolation.h:58
KeyType scale_1_
Definition: BilinearInterpolation.h:773
BilinearInterpolation()
Constructors and destructor.
Definition: BilinearInterpolation.h:70
KeyType const & getInsideReferencePoint_0() const
Accessor. See setMapping().
Definition: BilinearInterpolation.h:718
KeyType index2key_0(KeyType pos) const
The transformation from "inside" to "outside" coordinates.
Definition: BilinearInterpolation.h:535
KeyType offset_0_
Definition: BilinearInterpolation.h:772
KeyType scale_0_
Data members.
Definition: BilinearInterpolation.h:771
KeyType inside_0_
Definition: BilinearInterpolation.h:775
ContainerType & getData()
Returns the internal random access container storing the data.
Definition: BilinearInterpolation.h:484
void setMapping_0(KeyType const &inside_low, KeyType const &outside_low, KeyType const &inside_high, KeyType const &outside_high)
Specifies the mapping from "outside" to "inside" coordinates by the following data:
Definition: BilinearInterpolation.h:659
void addValue(KeyType arg_pos_0, KeyType arg_pos_1, ValueType arg_value)
Performs bilinear resampling. The arg_value is split up and added to the data points around arg_pos....
Definition: BilinearInterpolation.h:290
key_type KeyType
Definition: BilinearInterpolation.h:59
KeyType offset_1_
Definition: BilinearInterpolation.h:774
~BilinearInterpolation()
Destructor.
Definition: BilinearInterpolation.h:114
KeyType outside_1_
Definition: BilinearInterpolation.h:778
void setScale_1(KeyType const &scale)
Accessor. "Scale" is the difference (in "outside" units) between consecutive entries in "Data".
Definition: BilinearInterpolation.h:592
KeyType key2index_1(KeyType pos) const
The transformation from "outside" to "inside" coordinates.
Definition: BilinearInterpolation.h:543
void setOffset_0(KeyType const &offset)
Accessor. "Offset" is the point (in "outside" units) which corresponds to "Data(0,...
Definition: BilinearInterpolation.h:615
KeyType index2key_1(KeyType pos) const
The transformation from "inside" to "outside" coordinates.
Definition: BilinearInterpolation.h:558
KeyType const & getScale_0() const
Accessor. "Scale" is the difference (in "outside" units) between consecutive entries in "Data".
Definition: BilinearInterpolation.h:566
KeyType supportMin_0() const
Lower boundary of the support, in "outside" coordinates.
Definition: BilinearInterpolation.h:742
KeyType key2index_0(KeyType pos) const
The transformation from "outside" to "inside" coordinates.
Definition: BilinearInterpolation.h:520
void setScale_0(KeyType const &scale)
Accessor. "Scale" is the difference (in "outside" units) between consecutive entries in "Data".
Definition: BilinearInterpolation.h:582
A two-dimensional matrix. Similar to std::vector, but uses a binary operator(,) for element access.
Definition: Matrix.h:51
SizeType rows() const
Number of rows.
Definition: Matrix.h:232
Base::difference_type difference_type
Definition: Matrix.h:61
SizeType cols() const
Number of columns.
Definition: Matrix.h:238
ptrdiff_t SignedSize
Signed Size type e.g. used as pointer difference.
Definition: Types.h:108
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:22