Openholo  v4.2
Open Source Digital Holographic Library
vec.cpp
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 #include "vec.h"
47 
48 namespace oph {
49 
50 
51 const int vec2::n = 2;
52 
53 bool vec2::unit()
54 {
55  Real val = norm(*this);
56 
57  if (val < epsilon) return false;
58  (*this) = (*this)/val;
59  return true;
60 }
61 
63 {
64  return norm(*this);
65 }
66 
68  // returns 1: this and other vectors are and parallel
69  // -1: this and other vectors are anti-parallel
70  // 0: this and other vectors are not parallel
71  // or at least one of the vectors is zero
72  const vec2& vv,
73  Real angle_tolerance // (default=ON_DEFAULT_ANGLE_TOLERANCE) radians
74  ) const
75 {
76  int rc = 0;
77  const Real ll = norm(*this) * norm(vv);
78  if ( ll > 0.0 ) {
79  const Real cos_angle = (inner(*this, vv))/ll;
80  const Real cos_tol = cos(angle_tolerance);
81  if ( cos_angle >= cos_tol )
82  rc = 1;
83  else if ( cos_angle <= -cos_tol )
84  rc = -1;
85  }
86  return rc;
87 }
88 
90  // returns true: this and other vectors are perpendicular
91  // false: this and other vectors are not perpendicular
92  // or at least one of the vectors is zero
93  const vec2& vv,
94  Real angle_tolerance // (default=ON_DEFAULT_ANGLE_TOLERANCE) radians
95  ) const
96 {
97  bool rc = false;
98  const Real ll = norm(*this)*norm(vv);
99  if ( ll > 0.0 ) {
100  if ( fabs(inner(*this, vv)/ll) <= sin(angle_tolerance) )
101  rc = true;
102  }
103  return rc;
104 }
105 
106 // set this vector to be perpendicular to another vector
107 bool vec2::perpendicular( // Result is not unitized.
108  // returns false if input vector is zero
109  const vec2& vv
110  )
111 {
112  v[1] = vv[0];
113  v[0] = -vv[1];
114  return (v[0] != 0.0 || v[1] != 0.0) ? true : false;
115 }
116 
117 // set this vector to be perpendicular to a line defined by 2 points
119  const vec2& p,
120  const vec2& q
121  )
122 {
123  return perpendicular(q-p);
124 }
125 
126 
127 void store(FILE* fp, const vec2& v)
128 {
129  fprintf(fp, "(%lg", v[0]);
130  for(int i = 1; i < 2;++i){
131  fprintf(fp, " %lg", v[i]);
132  }
133  fprintf(fp, ")\n");
134 }
135 
136 //int scan(FILE* fp, const vec2& v)
137 //{
138 // int a = fscanf(fp, " (");
139 // for(int i = 0; i < 2;++i){
140 // a += fscanf(fp, " %lg", const_cast<Real*>(&v[i]));
141 // }
142 // a += fscanf(fp, " )");
143 // return a;
144 //}
145 
146 int apx_equal(const vec2& a, const vec2& b)
147 {
148  int c = 1;
149 
150  for (int i = 0 ; i < 2 ;++i){
151  c = c && apx_equal(a[i], b[i]);
152  }
153 
154  return c;
155 }
156 
157 int apx_equal(const vec2& a, const vec2& b, Real eps)
158 {
159  int c = 1;
160 
161  for (int i = 0 ; i < 2 ;++i){
162  c = c && apx_equal(a[i], b[i], eps);
163  }
164 
165  return c;
166 }
167 
168 
169 
170 
171 const int vec3::n = 3;
172 
173 bool vec3::unit()
174 {
175  Real val = norm(*this);
176 
177  if (val < epsilon) return false;
178  (*this) = (*this)/val;
179  return true;
180 }
181 
183 {
184  return norm(*this);
185 }
186 
188  // returns 1: this and other vectors are and parallel
189  // -1: this and other vectors are anti-parallel
190  // 0: this and other vectors are not parallel
191  // or at least one of the vectors is zero
192  const vec3& vv,
193  Real angle_tolerance // (default=ON_DEFAULT_ANGLE_TOLERANCE) radians
194  ) const
195 {
196  int rc = 0;
197  const Real ll = norm(*this) * norm(vv);
198  if ( ll > 0.0 ) {
199  const Real cos_angle = (inner(*this, vv))/ll;
200  const Real cos_tol = cos(angle_tolerance);
201  if ( cos_angle >= cos_tol )
202  rc = 1;
203  else if ( cos_angle <= -cos_tol )
204  rc = -1;
205  }
206  return rc;
207 }
208 
210  // returns true: this and other vectors are perpendicular
211  // false: this and other vectors are not perpendicular
212  // or at least one of the vectors is zero
213  const vec3& vv,
214  Real angle_tolerance // (default=ON_DEFAULT_ANGLE_TOLERANCE) radians
215  ) const
216 {
217  bool rc = false;
218  const Real ll = norm(*this) * norm(vv);
219  if ( ll > 0.0 ) {
220  if ( fabs(inner(oph::unit(*this), oph::unit(vv))/ll) <= sin(angle_tolerance) )
221  rc = true;
222  }
223  return rc;
224 }
225 
226 
227 bool vec3::perpendicular( const vec3& vv )
228 {
229  //bool rc = false;
230  int i, j, k;
231  Real a, b;
232  k = 2;
233  if ( fabs(vv[1]) > fabs(vv[0]) ) {
234  if ( fabs(vv[2]) > fabs(vv[1]) ) {
235  i = 2;
236  j = 1;
237  k = 0;
238  a = vv[2];
239  b = -vv[1];
240  }
241  else if ( fabs(vv[2]) >= fabs(vv[0]) ){
242  i = 1;
243  j = 2;
244  k = 0;
245  a = vv[1];
246  b = -vv[2];
247  }
248  else {
249  // |vv[1]| > |vv[0]| > |vv[2]|
250  i = 1;
251  j = 0;
252  k = 2;
253  a = vv[1];
254  b = -vv[0];
255  }
256  }
257  else if ( fabs(vv[2]) > fabs(vv[0]) ) {
258  // |vv[2]| > |vv[0]| >= |vv[1]|
259  i = 2;
260  j = 0;
261  k = 1;
262  a = vv[2];
263  b = -vv[0];
264  }
265  else if ( fabs(vv[2]) > fabs(vv[1]) ) {
266  // |vv[0]| >= |vv[2]| > |vv[1]|
267  i = 0;
268  j = 2;
269  k = 1;
270  a = vv[0];
271  b = -vv[2];
272  }
273  else {
274  // |vv[0]| >= |vv[1]| >= |vv[2]|
275  i = 0;
276  j = 1;
277  k = 2;
278  a = vv[0];
279  b = -vv[1];
280  }
281 
282  v[i] = b;
283  v[j] = a;
284  v[k] = 0.0;
285  return (a != 0.0) ? true : false;
286 }
287 
288 bool
290  const vec3& P0, const vec3& P1, const vec3& P2
291  )
292 {
293  // Find a the unit normal to a triangle defined by 3 points
294  vec3 V0, V1, V2, N0, N1, N2;
295 
296  v[0] = v[1] = v[2] = 0.0;
297 
298  V0 = P2 - P1;
299  V1 = P0 - P2;
300  V2 = P1 - P0;
301 
302  N0 = cross( V1, V2 );
303 
304  if (!N0.unit())
305  return false;
306 
307  N1 = cross( V2, V0 );
308 
309  if (!N1.unit())
310  return false;
311 
312  N2 = cross( V0, V1 );
313 
314  if (!N2.unit())
315  return false;
316 
317  const Real s0 = 1.0/V0.length();
318  const Real s1 = 1.0/V1.length();
319  const Real s2 = 1.0/V2.length();
320 
321  // choose normal with smallest total error
322  const Real e0 = s0*fabs(inner(N0,V0)) + s1*fabs(inner(N0,V1)) + s2*fabs(inner(N0,V2));
323  const Real e1 = s0*fabs(inner(N1,V0)) + s1*fabs(inner(N1,V1)) + s2*fabs(inner(N1,V2));
324  const Real e2 = s0*fabs(inner(N2,V0)) + s1*fabs(inner(N2,V1)) + s2*fabs(inner(N2,V2));
325 
326  if ( e0 <= e1 ) {
327  if ( e0 <= e2 ) {
328  *this = N0;
329  }
330  else {
331  *this = N2;
332  }
333  }
334  else if (e1 <= e2) {
335  *this = N1;
336  }
337  else {
338  *this = N2;
339  }
340 
341  return true;
342 }
343 
344 void store(FILE* fp, const vec3& v)
345 {
346  fprintf(fp, "(%lg", v[0]);
347  for(int i = 1; i < 3;++i){
348  fprintf(fp, " %lg", v[i]);
349  }
350  fprintf(fp, ")\n");
351 }
352 
353 //int scan(FILE* fp, const vec3& v)
354 //{
355 // int a = fscanf(fp, " (");
356 // for(int i = 0; i < 3;++i){
357 // a += fscanf(fp, " %lg", const_cast<Real*>(&v[i]));
358 // }
359 // a += fscanf(fp, " )");
360 // return a;
361 //}
362 
363 int apx_equal(const vec3& a, const vec3& b)
364 {
365  int c = 1;
366 
367  for (int i = 0 ; i < 3 ;++i){
368  c = c && apx_equal(a[i], b[i]);
369  }
370 
371  return c;
372 }
373 
374 int apx_equal(const vec3& a, const vec3& b, Real eps)
375 {
376  int c = 1;
377 
378  for (int i = 0 ; i < 3 ;++i){
379  c = c && apx_equal(a[i], b[i], eps);
380  }
381 
382  return c;
383 }
384 
385 
386 
387 const int vec4::n = 4;
388 
389 bool vec4::unit()
390 {
391  Real val = norm(*this);
392 
393  if (val < epsilon) return false;
394  (*this) = (*this)/val;
395  return true;
396 }
397 
399 {
400  return norm(*this);
401 }
402 
403 void store(FILE* fp, const vec4& v)
404 {
405  fprintf(fp, "(%lg", v[0]);
406  for(int i = 1; i < 4;++i){
407  fprintf(fp, " %lg", v[i]);
408  }
409  fprintf(fp, ")\n");
410 }
411 
412 //int scan(FILE* fp, const vec4& v)
413 //{
414 // int a = fscanf(fp, " (");
415 // for(int i = 0; i < 4;++i){
416 // a += fscanf(fp, " %lg", const_cast<Real*>(&v[i]));
417 // }
418 // a += fscanf(fp, " )");
419 // return a;
420 //}
421 
422 int apx_equal(const vec4& a, const vec4& b)
423 {
424  int c = 1;
425 
426  for (int i = 0 ; i < 4 ;++i){
427  c = c && oph::apx_equal(a[i], b[i]);
428  }
429 
430  return c;
431 }
432 
433 int apx_equal(const vec4& a, const vec4& b, Real eps)
434 {
435  int c = 1;
436 
437  for (int i = 0 ; i < 4 ;++i){
438  c = c && oph::apx_equal(a[i], b[i], eps);
439  }
440 
441  return c;
442 }
443 
444 vec3 cross(const vec3& a, const vec3& b)
445 {
446  vec3 c;
447 
448  c(0) = a(0 + 1) * b(0 + 2) - a(0 + 2) * b(0 + 1);
449 
450  c(1) = a(1 + 1) * b(1 + 2) - a(1 + 2) * b(1 + 1);
451 
452  c(2) = a(2 + 1) * b(2 + 2) - a(2 + 2) * b(2 + 1);
453 
454 
455  return c;
456 }
457 
458 }; // namespace graphics
Real inner(const vec2 &a, const vec2 &b)
Definition: vec.h:411
Real length() const
Definition: vec.cpp:182
bool is_perpendicular(const vec3 &, Real=angle_tolerance) const
Definition: vec.cpp:209
Real length() const
Definition: vec.cpp:398
Real norm(const vec2 &a)
Definition: vec.h:417
return true
Definition: Openholo.cpp:434
float Real
Definition: typedef.h:55
int is_parallel(const vec3 &, Real=angle_tolerance) const
Definition: vec.cpp:187
def k(wvl)
Definition: Depthmap.py:16
bool unit()
Definition: vec.cpp:173
int is_parallel(const vec2 &, Real=angle_tolerance) const
Definition: vec.cpp:67
void store(FILE *fp, const vec2 &v)
Definition: vec.cpp:127
Real length() const
Definition: vec.cpp:62
Real epsilon
Definition: epsilon.cpp:53
int apx_equal(Real x, Real y)
Definition: epsilon.cpp:91
vec3 cross(const vec3 &a, const vec3 &b)
Definition: vec.cpp:444
bool perpendicular(const vec3 &)
Definition: vec.cpp:227
bool unit()
Definition: vec.cpp:389
structure for 2-dimensional Real type vector and its arithmetic.
Definition: vec.h:66
bool perpendicular(const vec2 &)
Definition: vec.cpp:107
structure for 3-dimensional Real type vector and its arithmetic.
Definition: vec.h:466
static const int n
Definition: vec.h:68
bool unit()
Definition: vec.cpp:53
Real v[3]
Definition: vec.h:467
Real angle_tolerance
Definition: epsilon.cpp:60
static const int n
Definition: vec.h:468
Real v[2]
Definition: vec.h:67
Definition: Bitmap.h:49
structure for 4-dimensional Real type vector and its arithmetic.
Definition: vec.h:864
static const int n
Definition: vec.h:866
vec2 unit(const vec2 &a)
Definition: vec.h:426
bool is_perpendicular(const vec2 &, Real=angle_tolerance) const
Definition: vec.cpp:89