Openholo  v4.2
Open Source Digital Holographic Library
complex.h
Go to the documentation of this file.
1 /*M///////////////////////////////////////////////////////////////////////////////////////
2 //
3 // IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4 //
5 // By downloading, copying, installing or using the software you agree to this license.
6 // If you do not agree to this license, do not download, install, copy or use the software.
7 //
8 //
9 // License Agreement
10 // For Open Source Digital Holographic Library
11 //
12 // Openholo library is free software;
13 // you can redistribute it and/or modify it under the terms of the BSD 2-Clause license.
14 //
15 // Copyright (C) 2017-2024, Korea Electronics Technology Institute. All rights reserved.
16 // E-mail : contact.openholo@gmail.com
17 // Web : http://www.openholo.org
18 //
19 // Redistribution and use in source and binary forms, with or without modification,
20 // are permitted provided that the following conditions are met:
21 //
22 // 1. Redistribution's of source code must retain the above copyright notice,
23 // this list of conditions and the following disclaimer.
24 //
25 // 2. Redistribution's in binary form must reproduce the above copyright notice,
26 // this list of conditions and the following disclaimer in the documentation
27 // and/or other materials provided with the distribution.
28 //
29 // This software is provided by the copyright holders and contributors "as is" and
30 // any express or implied warranties, including, but not limited to, the implied
31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
32 // In no event shall the copyright holder or contributors be liable for any direct,
33 // indirect, incidental, special, exemplary, or consequential damages
34 // (including, but not limited to, procurement of substitute goods or services;
35 // loss of use, data, or profits; or business interruption) however caused
36 // and on any theory of liability, whether in contract, strict liability,
37 // or tort (including negligence or otherwise) arising in any way out of
38 // the use of this software, even if advised of the possibility of such damage.
39 //
40 // This software contains opensource software released under GNU Generic Public License,
41 // NVDIA Software License Agreement, or CUDA supplement to Software License Agreement.
42 // Check whether software you use contains licensed software.
43 //
44 //M*/
45 
46 #ifndef __complex_h
47 #define __complex_h
48 
49 #include <iostream>
50 #include <cmath>
51 #include <complex> // c11 standard
52 #include "define.h"
53 
54 #ifndef _RE
55 #define _RE 0
56 #endif
57 #ifndef _IM
58 #define _IM 1
59 #endif
60 
61 #pragma once
62 #undef OPH_DLL
63 #ifdef OPH_EXPORT
64 #ifdef _WIN32
65 #define OPH_DLL __declspec(dllexport)
66 #else
67 #define OPH_DLL __attribute__((visibility("default")))
68 #endif
69 #else
70 #define OPH_DLL
71 #endif
72 
73 namespace oph {
82 #if _MSC_VER > 1900 || defined(__GNUC__)
83  template<typename T = double>
84  class OPH_DLL Complex : public std::complex<T>
85  {
86  public:
87  using cplx = typename std::enable_if<std::is_same<double, T>::value || std::is_same<float, T>::value || std::is_same<long double, T>::value, T>::type;
88 
89  public:
90  Complex() : std::complex<T>() {}
91  Complex(T p) : std::complex<T>(p) { this->real(p); this->imag((T)0.0); }
92  Complex(T tRe, T tIm) : std::complex<T>(tRe, tIm) {}
93  Complex(const Complex<T>& p) {
94  this->real(p.real());
95  this->imag(p.imag());
96  }
97  T mag2() const { return this->real() * this->real() + this->imag() * this->imag(); }
98  T mag() const { return sqrt(mag2()); }
99 
100  T arg() const
101  {
102  T r = mag();
103  T theta = acos(this->real() / r);
104 
105  if (sin(theta) - this->imag() / r < 10e-6)
106  return theta;
107  else
108  return 2.0 * M_PI - theta;
109  }
110 
111  void euler(T& r, T& theta)
112  {
113  r = mag();
114  theta = arg();
115  }
116 
117  T angle(void)
118  {
119  if (std::is_same<double, T>::value)
120  return atan2(this->imag(), this->real());
121  else if (std::is_same<float, T>::value)
122  return atan2f((float)this->imag(), (float)this->real());
123  }
124 
125  Complex<T>& exp() {
126  Complex<T> p(this->real(), this->imag());
127  if (std::is_same<double, T>::value) {
128  this->real(std::exp(p.real()) * cos(p.imag()));
129  this->imag(std::exp(p.real()) * sin(p.imag()));
130  }
131  else {
132 #ifdef _MSC_VER
133  this->real(std::expf(p.real()) * cos(p.imag()));
134  this->imag(std::expf(p.real()) * sin(p.imag()));
135 #else
136  this->real(std::exp(p.real()) * cos(p.imag()));
137  this->imag(std::exp(p.real()) * sin(p.imag()));
138 #endif
139  }
140  return *this;
141  }
142 
143  Complex<T> conj() const { return Complex<T>(this->real(), -this->imag()); }
144 
145  Complex<T>& operator()(T re, T im) {
146  this->real(re);
147  this->imag(im);
148 
149  return *this;
150  }
151 
152  // arithmetic
153  Complex<T>& operator= (const Complex<T>& p) {
154  this->real(p.real());
155  this->imag(p.imag());
156 
157  return *this;
158  }
159 
160  Complex<T>& operator = (const T& p) {
161  this->real(p);
162  this->imag(0.0);
163 
164  return *this;
165  }
166 
167  Complex<T> operator+ (const Complex<T>& p) {
168  Complex<T> n(this->real() + p.real(), this->imag() + p.imag());
169 
170  return n;
171  }
172 
173  Complex<T>& operator+= (const Complex<T>& p) {
174  this->real(this->real() + p.real());
175  this->imag(this->imag() + p.imag());
176 
177  return *this;
178  }
179 
180  Complex<T> operator+ (const T p) {
181  Complex<T> n(this->real() + p, this->imag());
182 
183  return n;
184  }
185 
186  Complex<T>& operator+= (const T p) {
187  this->real(this->real() + p);
188 
189  return *this;
190  }
191 
192  Complex<T> operator- (const Complex<T>& p) {
193  Complex<T> n(this->real() - p.real(), this->imag() - p.imag());
194 
195  return n;
196  }
197 
198  Complex<T>& operator-= (const Complex<T>& p) {
199  this->real(this->real() - p.real());
200  this->imag(this->imag() - p.imag());
201 
202  return *this;
203  }
204 
205  Complex<T> operator - (const T p) {
206  Complex<T> n(this->real() - p, this->imag());
207 
208  return n;
209  }
210 
211  Complex<T>& operator -= (const T p) {
212  this->real(this->real() - p);
213 
214  return *this;
215  }
216 
217  Complex<T> operator* (const T k) {
218  Complex<T> n(this->real() * k, this->imag() * k);
219 
220  return n;
221  }
222 
223  Complex<T>& operator*= (const T k) {
224  this->real(this->real() * k);
225  this->imag(this->imag() * k);
226 
227  return *this;
228  }
229 
230  Complex<T>& operator = (const std::complex<T>& p) {
231  this->real(p.real());
232  this->imag(p.imag());
233 
234  return *this;
235  }
236 
237  Complex<T> operator* (const Complex<T>& p) {
238  const T tRe = this->real();
239  const T tIm = this->imag();
240 
241  Complex<T> n(tRe * p.real() - tIm * p.imag(), tRe * p.imag() + tIm * p.real());
242 
243  return n;
244  }
245 
246  Complex<T>& operator*= (const Complex<T>& p) {
247  const T tRe = this->real();
248  const T tIm = this->imag();
249 
250  this->real(tRe * p.real() - tIm * p.imag());
251  this->imag(tRe * p.imag() + tIm * p.real());
252 
253  return *this;
254  }
255 
256  Complex<T> operator / (const T& p) {
257  Complex<T> n(this->real() / p, this->imag() / p);
258 
259  return n;
260  }
261 
262  Complex<T>& operator/= (const T k) {
263  this->real(this->real() / k);
264  this->imag(this->imag() / k);
265 
266  return *this;
267  }
268 
269  Complex<T> operator / (const Complex<T>& p) {
270  std::complex<T> a(this->real(), this->imag());
271  std::complex<T> b(p.real(), p.imag());
272 
273  std::complex<T> c = a / b;
274 
275  Complex<T> n(c.real(), c.imag());
276 
277  return n;
278  }
279 
280  Complex<T>& operator /= (const Complex<T>& p) {
281  std::complex<T> a(this->real(), this->imag());
282  std::complex<T> b(p.real(), p.imag());
283 
284  a /= b;
285  this->real(a.real());
286  this->imag(a.imag());
287 
288  return *this;
289  }
290 
291  T& operator [](const int idx) {
292  return reinterpret_cast<T*>(this)[idx];
293  }
294 
295  bool operator < (const Complex<T>& p) {
296  return (this->real() < p.real());
297  }
298 
299  bool operator > (const Complex<T>& p) {
300  return (this->real() > p.real());
301  }
302 
303  operator unsigned char() {
304  return uchar(this->real());
305  }
306 
307  operator int() {
308  return int(this->real());
309  }
310 
311  friend Complex<T> operator+ (const Complex<T>&p, const T q) {
312  return Complex<T>(p.real() + q, p.imag());
313  }
314 
315  friend Complex<T> operator- (const Complex<T>&p, const T q) {
316  return Complex<T>(p.real() - q, p.imag());
317  }
318 
319  friend Complex<T> operator* (const T k, const Complex<T>& p) {
320  return Complex<T>(p) *= k;
321  }
322 
323  friend Complex<T> operator* (const Complex<T>& p, const T k) {
324  return Complex<T>(p) *= k;
325  }
326 
327  friend Complex<T> operator* (const Complex<T>& p, const Complex<T>& q) {
328  return Complex<T>(p.real() * q.real() - p.imag() * q.imag(), p.real() * q.imag() + p.imag() * q.real());
329  }
330 
331  friend Complex<T> operator/ (const Complex<T>& p, const Complex<T>& q) {
332  return Complex<T>((1.0 / q.mag2())*(p*q.conj()));
333  }
334 
335  friend Complex<T> operator/ (const Complex<T>& p, const T& q) {
336  return Complex<T>(p.real() / q, p.imag() / q);
337  }
338 
339  friend Complex<T> operator/ (const T& p, const Complex<T>& q) {
340  return Complex<T>(p / q.real(), p / q.imag());
341  }
342 
343  friend bool operator< (const Complex<T>& p, const Complex<T>& q) {
344  return (p.real() < q.real());
345  }
346 
347  friend bool operator> (const Complex<T>& p, const Complex<T>& q) {
348  return (p.real() > q.real());
349  }
350  };
351 
352 #else
353 template<typename T = double>
354 class OPH_DLL Complex : public std::complex<T>
355 {
356 public:
357  using cplx = typename std::enable_if<std::is_same<double, T>::value || std::is_same<float, T>::value || std::is_same<long double, T>::value, T>::type;
358 
359 public:
360  Complex() : std::complex<T>() {}
361  Complex(T p) : std::complex<T>(p) { this->_Val[_RE] = p; this->_Val[_IM] = (T)0.0; }
362  Complex(T tRe, T tIm) : std::complex<T>(tRe, tIm) {}
363  Complex(const Complex<T>& p) {
364  this->_Val[_RE] = p._Val[_RE];
365  this->_Val[_IM] = p._Val[_IM];
366  }
367  T mag2() const { return this->_Val[_RE] * this->_Val[_RE] + this->_Val[_IM] * this->_Val[_IM]; }
368  T mag() const { return sqrt(mag2()); }
369 
370  T arg() const
371  {
372  T r = mag();
373  T theta = acos(this->_Val[_RE] / r);
374 
375  if (sin(theta) - this->_Val[_IM] / r < 10e-6)
376  return theta;
377  else
378  return 2.0 * M_PI - theta;
379  }
380 
381  void euler(T& r, T& theta)
382  {
383  r = mag();
384  theta = arg();
385  }
386 
387  T angle(void)
388  {
389  if (std::is_same<double, T>::value)
390  return atan2(this->_Val[_IM], this->_Val[_RE]);
391  else if (std::is_same<float, T>::value)
392  return atan2f((float)this->_Val[_IM], (float)this->_Val[_RE]);
393  }
394 
396  Complex<T> p(this->_Val[_RE], this->_Val[_IM]);
397  if (std::is_same<double, T>::value) {
398  this->_Val[_RE] = std::exp(p._Val[_RE]) * cos(p._Val[_IM]);
399  this->_Val[_IM] = std::exp(p._Val[_RE]) * sin(p._Val[_IM]);
400  }
401  else {
402 #ifdef _MSC_VER
403  this->_Val[_RE] = std::expf(p._Val[_RE]) * cos(p._Val[_IM]);
404  this->_Val[_IM] = std::expf(p._Val[_RE]) * sin(p._Val[_IM]);
405 #else
406  this->_Val[_RE] = std::exp(p._Val[_RE]) * cos(p._Val[_IM]);
407  this->_Val[_IM] = std::exp(p._Val[_RE]) * sin(p._Val[_IM]);
408 #endif
409  }
410  return *this;
411  }
412 
413  Complex<T> conj() const { return Complex<T>(this->_Val[_RE], -this->_Val[_IM]); }
414 
415  Complex<T>& operator()(T re, T im) {
416  this->_Val[_RE] = re;
417  this->_Val[_IM] = im;
418 
419  return *this;
420  }
421 
422  // arithmetic
423  Complex<T>& operator= (const Complex<T>& p) {
424  this->_Val[_RE] = p._Val[_RE];
425  this->_Val[_IM] = p._Val[_IM];
426 
427  return *this;
428  }
429 
430  Complex<T>& operator = (const T& p) {
431  this->_Val[_RE] = p;
432  this->_Val[_IM] = 0.0;
433 
434  return *this;
435  }
436 
438  Complex<T> n(this->_Val[_RE] + p._Val[_RE], this->_Val[_IM] + p._Val[_IM]);
439 
440  return n;
441  }
442 
444  this->_Val[_RE] = this->_Val[_RE] + p._Val[_RE];
445  this->_Val[_IM] = this->_Val[_IM] + p._Val[_IM];
446 
447  return *this;
448  }
449 
450  Complex<T> operator+ (const T p) {
451  Complex<T> n(this->_Val[_RE] + p, this->_Val[_IM]);
452 
453  return n;
454  }
455 
456  Complex<T>& operator+= (const T p) {
457  this->_Val[_RE] = this->_Val[_RE] + p;
458 
459  return *this;
460  }
461 
463  Complex<T> n(this->_Val[_RE] - p._Val[_RE], this->_Val[_IM] - p._Val[_IM]);
464 
465  return n;
466  }
467 
469  this->_Val[_RE] = this->_Val[_RE] - p._Val[_RE];
470  this->_Val[_IM] = this->_Val[_IM] - p._Val[_IM];
471 
472  return *this;
473  }
474 
475  Complex<T> operator - (const T p) {
476  Complex<T> n(this->_Val[_RE] - p, this->_Val[_IM]);
477 
478  return n;
479  }
480 
481  Complex<T>& operator -= (const T p) {
482  this->_Val[_RE] = this->_Val[_RE] - p;
483 
484  return *this;
485  }
486 
488  Complex<T> n(this->_Val[_RE] * k, this->_Val[_IM] * k);
489 
490  return n;
491  }
492 
494  this->_Val[_RE] = this->_Val[_RE] * k;
495  this->_Val[_IM] = this->_Val[_IM] * k;
496 
497  return *this;
498  }
499 
500  Complex<T>& operator = (const std::complex<T>& p) {
501  this->_Val[_RE] = p._Val[_RE];
502  this->_Val[_IM] = p._Val[_IM];
503 
504  return *this;
505  }
506 
508  const T tRe = this->_Val[_RE];
509  const T tIm = this->_Val[_IM];
510 
511  Complex<T> n(tRe * p._Val[_RE] - tIm * p._Val[_IM], tRe * p._Val[_IM] + tIm * p._Val[_RE]);
512 
513  return n;
514  }
515 
517  const T tRe = this->_Val[_RE];
518  const T tIm = this->_Val[_IM];
519 
520  this->_Val[_RE] = tRe * p._Val[_RE] - tIm * p._Val[_IM];
521  this->_Val[_IM] = tRe * p._Val[_IM] + tIm * p._Val[_RE];
522 
523  return *this;
524  }
525 
526  Complex<T> operator / (const T& p) {
527  Complex<T> n(this->_Val[_RE] / p, this->_Val[_IM] / p);
528 
529  return n;
530  }
531 
533  this->_Val[_RE] = this->_Val[_RE] / k;
534  this->_Val[_IM] = this->_Val[_IM] / k;
535 
536  return *this;
537  }
538 
540  std::complex<T> a(this->_Val[_RE], this->_Val[_IM]);
541  std::complex<T> b(p._Val[_RE], p._Val[_IM]);
542 
543  std::complex<T> c = a / b;
544 
545  Complex<T> n(c._Val[_RE], c._Val[_IM]);
546 
547  return n;
548  }
549 
551  std::complex<T> a(this->_Val[_RE], this->_Val[_IM]);
552  std::complex<T> b(p._Val[_RE], p._Val[_IM]);
553 
554  a /= b;
555  this->_Val[_RE] = a._Val[_RE];
556  this->_Val[_IM] = a._Val[_IM];
557 
558  return *this;
559  }
560 
561  T& operator [](const int idx) {
562  return reinterpret_cast<T*>(this)[idx];
563  }
564 
565  bool operator < (const Complex<T>& p) {
566  return (this->_Val[_RE] < p._Val[_RE]);
567  }
568 
569  bool operator > (const Complex<T>& p) {
570  return (this->_Val[_RE] > p._Val[_RE]);
571  }
572 
573  operator unsigned char() {
574  return uchar(this->_Val[_RE]);
575  }
576 
577  operator int() {
578  return int(this->_Val[_RE]);
579  }
580 
581  friend Complex<T> operator+ (const Complex<T>& p, const T q) {
582  return Complex<T>(p._Val[_RE] + q, p._Val[_IM]);
583  }
584 
585  friend Complex<T> operator- (const Complex<T>& p, const T q) {
586  return Complex<T>(p._Val[_RE] - q, p._Val[_IM]);
587  }
588 
589  friend Complex<T> operator* (const T k, const Complex<T>& p) {
590  return Complex<T>(p) *= k;
591  }
592 
593  friend Complex<T> operator* (const Complex<T>& p, const T k) {
594  return Complex<T>(p) *= k;
595  }
596 
597  friend Complex<T> operator* (const Complex<T>& p, const Complex<T>& q) {
598  return Complex<T>(p._Val[_RE] * q._Val[_RE] - p._Val[_IM] * q._Val[_IM], p._Val[_RE] * q._Val[_IM] + p._Val[_IM] * q._Val[_RE]);
599  }
600 
601  friend Complex<T> operator/ (const Complex<T>& p, const Complex<T>& q) {
602  return Complex<T>((1.0 / q.mag2()) * (p * q.conj()));
603  }
604 
605  friend Complex<T> operator/ (const Complex<T>& p, const T& q) {
606  return Complex<T>(p._Val[_RE] / q, p._Val[_IM] / q);
607  }
608 
609  friend Complex<T> operator/ (const T& p, const Complex<T>& q) {
610  return Complex<T>(p / q._Val[_RE], p / q._Val[_IM]);
611  }
612 
613  friend bool operator< (const Complex<T>& p, const Complex<T>& q) {
614  return (p._Val[_RE] < q._Val[_RE]);
615  }
616 
617  friend bool operator> (const Complex<T>& p, const Complex<T>& q) {
618  return (p._Val[_RE] > q._Val[_RE]);
619  }
620  };
621 
622 
623 
624 
625 #endif
626 }
627 
628 
629 #endif // !__complex_h_
T mag2() const
Definition: complex.h:367
ivec2 operator+=(ivec2 &a, const ivec2 &b)
Definition: ivec.h:204
#define M_PI
Definition: define.h:52
class for the complex number and its arithmetic. type T == type cplx type only float || double T real...
Definition: complex.h:354
T angle(void)
Definition: complex.h:387
Complex(T p)
Definition: complex.h:361
unsigned char uchar
Definition: typedef.h:64
#define _IM
Definition: complex.h:58
ivec2 operator-(const ivec2 &a, const ivec2 &b)
Definition: ivec.h:132
def k(wvl)
Definition: Depthmap.py:16
ivec2 operator+(const ivec2 &a, const ivec2 &b)
Definition: ivec.h:109
ivec2 operator*(const ivec2 &a, const ivec2 &b)
Definition: ivec.h:155
Complex< T > & operator()(T re, T im)
Definition: complex.h:415
vec2 operator/=(vec2 &a, const vec2 &b)
Definition: vec.h:262
T arg() const
Definition: complex.h:370
vec2 operator/(const vec2 &a, const vec2 &b)
Definition: vec.h:200
#define OPH_DLL
Definition: complex.h:70
ivec2 operator-=(ivec2 &a, const ivec2 &b)
Definition: ivec.h:216
void euler(T &r, T &theta)
Definition: complex.h:381
Complex< T > & exp()
Definition: complex.h:395
Complex< T > conj() const
Definition: complex.h:413
typename std::enable_if< std::is_same< double, Real >::value||std::is_same< float, Real >::value||std::is_same< long double, Real >::value, Real >::type cplx
Definition: complex.h:357
int operator>(const ivec2 &a, const ivec2 &b)
Definition: ivec.h:324
Complex(const Complex< T > &p)
Definition: complex.h:363
void angle(const std::vector< Complex< T >> &src, std::vector< T > &dst)
Definition: function.h:159
Definition: Bitmap.h:49
T mag() const
Definition: complex.h:368
Complex(T tRe, T tIm)
Definition: complex.h:362
ivec2 operator*=(ivec2 &a, const ivec2 &b)
Definition: ivec.h:228
#define _RE
Definition: complex.h:55