Openholo  v4.2
Open Source Digital Holographic Library
ophGen.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 "ophGen.h"
47 #include "sys.h"
48 #include "function.h"
49 #include <cuda_runtime_api.h>
50 #include <cufft.h>
51 #include <omp.h>
52 #include <iostream>
53 #include <sstream>
54 #include <algorithm>
55 #include "ImgControl.h"
56 #include "tinyxml2.h"
57 #include "PLYparser.h"
58 
59 extern "C"
60 {
73  void cudaFFT(CUstream_st* stream, int nx, int ny, cufftDoubleComplex* in_filed, cufftDoubleComplex* output_field, int direction, bool bNormailized = false);
74 
89  void cudaCropFringe(CUstream_st* stream, int nx, int ny, cufftDoubleComplex* in_field, cufftDoubleComplex* out_field, int cropx1, int cropx2, int cropy1, int cropy2);
90 
108  void cudaGetFringe(CUstream_st* stream, int pnx, int pny, cufftDoubleComplex* in_field, cufftDoubleComplex* out_field, int sig_locationx, int sig_locationy,
109  Real ssx, Real ssy, Real ppx, Real ppy, Real PI);
110 }
111 
112 ophGen::ophGen(void)
113  : Openholo()
114  , m_vecEncodeSize()
115  , ENCODE_METHOD(0)
116  , SSB_PASSBAND(0)
117  , m_elapsedTime(0.0)
118  , m_lpEncoded(nullptr)
119  , m_lpNormalized(nullptr)
120  , m_nOldChannel(0)
121  , m_precision(PRECISION::DOUBLE)
122  , m_dFieldLength(0.0)
123  , m_nStream(1)
124  , m_mode(0)
125  , AS(nullptr)
126  , normalized(nullptr)
127  , fftTemp(nullptr)
128  , weight(nullptr)
129  , weightC(nullptr)
130  , freqW(nullptr)
131  , realEnc(nullptr)
132  , binary(nullptr)
133  , maskSSB(nullptr)
134  , maskHP(nullptr)
135  , m_bRandomPhase(false)
136 {
137 
138 }
139 
140 ophGen::~ophGen(void)
141 {
142 
143 }
144 
145 void ophGen::initialize(void)
146 {
147  // Output Image Size
149  const uint nChannel = context_.waveNum;
150 
151  // Memory Location for Result Image
152  if (complex_H != nullptr) {
153  for (uint i = 0; i < m_nOldChannel; i++) {
154  if (complex_H[i] != nullptr) {
155  delete[] complex_H[i];
156  complex_H[i] = nullptr;
157  }
158  }
159  delete[] complex_H;
160  complex_H = nullptr;
161  }
162 
163  complex_H = new Complex<Real>*[nChannel];
164  for (uint i = 0; i < nChannel; i++) {
165  complex_H[i] = new Complex<Real>[rc.size];
166  memset(complex_H[i], 0, sizeof(Complex<Real>) * rc.size);
167  }
168 
169  if (m_lpEncoded != nullptr) {
170  for (uint i = 0; i < m_nOldChannel; i++) {
171  if (m_lpEncoded[i] != nullptr) {
172  delete[] m_lpEncoded[i];
173  m_lpEncoded[i] = nullptr;
174  }
175  }
176  delete[] m_lpEncoded;
177  m_lpEncoded = nullptr;
178  }
179  m_lpEncoded = new Real*[nChannel];
180  for (uint i = 0; i < nChannel; i++) {
181  m_lpEncoded[i] = new Real[rc.size];
182  memset(m_lpEncoded[i], 0, sizeof(Real) * rc.size);
183  }
184 
185  if (m_lpNormalized != nullptr) {
186  for (uint i = 0; i < m_nOldChannel; i++) {
187  if (m_lpNormalized[i] != nullptr) {
188  delete[] m_lpNormalized[i];
189  m_lpNormalized[i] = nullptr;
190  }
191  }
192  delete[] m_lpNormalized;
193  m_lpNormalized = nullptr;
194  }
195  m_lpNormalized = new uchar*[nChannel];
196  for (uint i = 0; i < nChannel; i++) {
197  m_lpNormalized[i] = new uchar[rc.size];
198  memset(m_lpNormalized[i], 0, sizeof(uchar) * rc.size);
199  }
200 
201  m_nOldChannel = nChannel;
204 }
205 
206 int ophGen::loadPointCloud(const char* pc_file, OphPointCloudData *pc_data_)
207 {
208  int n_points = 0;
209  auto begin = CUR_TIME;
210 
211  PLYparser plyIO;
212  if (!plyIO.loadPLY(pc_file, pc_data_->n_points, &pc_data_->vertices))
213  n_points = -1;
214  else
215  n_points = pc_data_->n_points;
216  LOG("%s : %.5lf (sec)\n", __FUNCTION__, ELAPSED_TIME(begin, CUR_TIME));
217  return n_points;
218 }
219 
220 bool ophGen::readConfig(const char* fname)
221 {
222  bool bRet = true;
223  using namespace tinyxml2;
224  tinyxml2::XMLDocument xml_doc;
225  XMLNode *xml_node = nullptr;
226 
227  if (!checkExtension(fname, ".xml"))
228  {
229  LOG("<FAILED> Wrong file ext.\n");
230  return false;
231  }
232  auto ret = xml_doc.LoadFile(fname);
233  if (ret != XML_SUCCESS)
234  {
235  LOG("<FAILED> Loading file (%d)\n", ret);
236  return false;
237  }
238  xml_node = xml_doc.FirstChild();
239 
240  int nWave = 1;
241  char szNodeName[32] = { 0, };
242 
243  sprintf(szNodeName, "SLM_WaveNum");
244  auto next = xml_node->FirstChildElement(szNodeName); // OffsetInDepth
245  if (!next || XML_SUCCESS != next->QueryIntText(&nWave))
246  {
247  LOG("<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
248  bRet = false;
249  }
250  context_.waveNum = nWave;
252  context_.wave_length = new Real[nWave];
253 
254  for (int i = 1; i <= nWave; i++) {
255  sprintf(szNodeName, "SLM_WaveLength_%d", i);
256  next = xml_node->FirstChildElement(szNodeName);
257  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.wave_length[i - 1]))
258  {
259  LOG("<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
260  bRet = false;
261  }
262  }
263 
264  sprintf(szNodeName, "SLM_PixelNumX");
265  next = xml_node->FirstChildElement(szNodeName);
266  if (!next || XML_SUCCESS != next->QueryIntText(&context_.pixel_number[_X]))
267  {
268  LOG("<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
269  bRet = false;
270  }
271 
272  sprintf(szNodeName, "SLM_PixelNumY");
273  next = xml_node->FirstChildElement(szNodeName);
274  if (!next || XML_SUCCESS != next->QueryIntText(&context_.pixel_number[_Y]))
275  {
276  LOG("<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
277  bRet = false;
278  }
279 
280  sprintf(szNodeName, "SLM_PixelPitchX");
281  next = xml_node->FirstChildElement("SLM_PixelPitchX");
282  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.pixel_pitch[_X]))
283  {
284  LOG("<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
285  bRet = false;
286  }
287 
288  sprintf(szNodeName, "SLM_PixelPitchY");
289  next = xml_node->FirstChildElement("SLM_PixelPitchY");
290  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.pixel_pitch[_Y]))
291  {
292  LOG("<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
293  bRet = false;
294  }
295 
296  // option
297  next = xml_node->FirstChildElement("IMG_Rotation");
298  if (!next || XML_SUCCESS != next->QueryBoolText(&imgCfg.rotate))
299  imgCfg.rotate = false;
300  next = xml_node->FirstChildElement("IMG_Merge");
301  if (!next || XML_SUCCESS != next->QueryBoolText(&imgCfg.merge))
302  imgCfg.merge = false;
303  next = xml_node->FirstChildElement("IMG_Flip");
304  if (!next || XML_SUCCESS != next->QueryIntText(&imgCfg.flip))
305  imgCfg.flip = 0;
306  next = xml_node->FirstChildElement("DoublePrecision");
307  if (!next || XML_SUCCESS != next->QueryBoolText(&context_.bUseDP))
308  context_.bUseDP = true;
309  next = xml_node->FirstChildElement("ShiftX");
310  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.shift[_X]))
311  context_.shift[_X] = 0.0;
312  next = xml_node->FirstChildElement("ShiftY");
313  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.shift[_Y]))
314  context_.shift[_Y] = 0.0;
315  next = xml_node->FirstChildElement("ShiftZ");
316  if (!next || XML_SUCCESS != next->QueryDoubleText(&context_.shift[_Z]))
317  context_.shift[_Z] = 0.0;
318  next = xml_node->FirstChildElement("FieldLength");
319  if (!next || XML_SUCCESS != next->QueryDoubleText(&m_dFieldLength))
320  m_dFieldLength = 0.0;
321  next = xml_node->FirstChildElement("RandomPhase");
322  if (!next || XML_SUCCESS != next->QueryBoolText(&m_bRandomPhase))
323  m_bRandomPhase = false;
324  next = xml_node->FirstChildElement("NumOfStream");
325  if (!next || XML_SUCCESS != next->QueryIntText(&m_nStream))
326  m_nStream = 1;
327 
330 
333 
334  OHC_encoder->clearWavelength();
335  for (int i = 0; i < nWave; i++)
337 
338  // 2024.04.23. mwnam
339  // set variable for resolution
341 
342  LOG("**************************************************\n");
343  LOG(" Read Config (Common) \n");
344  LOG("1) SLM Number of Waves : %d\n", context_.waveNum);
345  for (uint i = 0; i < context_.waveNum; i++)
346  LOG(" 1-%d) SLM Wave length : %e\n", i + 1, context_.wave_length[i]);
347  LOG("2) SLM Resolution : %d x %d\n", context_.pixel_number[_X], context_.pixel_number[_Y]);
348  LOG("3) SLM Pixel Pitch : %e x %e\n", context_.pixel_pitch[_X], context_.pixel_pitch[_Y]);
349  LOG("4) Image Rotate : %s\n", imgCfg.rotate ? "Y" : "N");
350  LOG("5) Image Flip : %s\n", (imgCfg.flip == FLIP::NONE) ? "NONE" :
351  (imgCfg.flip == FLIP::VERTICAL) ? "VERTICAL" :
352  (imgCfg.flip == FLIP::HORIZONTAL) ? "HORIZONTAL" : "BOTH");
353  LOG("6) Image Merge : %s\n", imgCfg.merge ? "Y" : "N");
354  LOG("**************************************************\n");
355 
356  return bRet;
357 }
358 
359 
360 void ophGen::RS_Diffraction(Point src, Complex<Real> *dst, Real lambda, Real distance, Real amplitude)
361 {
362  OphConfig *pConfig = &context_;
363  const int pnX = pConfig->pixel_number[_X];
364  const int pnY = pConfig->pixel_number[_Y];
365  const Real ppX = pConfig->pixel_pitch[_X];
366  const Real ppY = pConfig->pixel_pitch[_Y];
367  const Real ssX = pConfig->ss[_X] = pnX * ppX;
368  const Real ssY = pConfig->ss[_Y] = pnY * ppY;
369  const int offsetX = pConfig->offset[_X];
370  const int offsetY = pConfig->offset[_Y];
371 
372  const Real tx = lambda / (2 * ppX);
373  const Real ty = lambda / (2 * ppY);
374  const Real sqrtX = sqrt(1 - (tx * tx));
375  const Real sqrtY = sqrt(1 - (ty * ty));
376  const Real x = -ssX / 2;
377  const Real y = -ssY / 2;
378  const Real k = (2 * M_PI) / lambda;
379  Real z = src.pos[_Z] + distance;
380  Real zz = z * z;
381  Real ampZ = amplitude * z;
382 
383  Real _xbound[2] = {
384  src.pos[_X] + abs(tx / sqrtX * z),
385  src.pos[_X] - abs(tx / sqrtX * z)
386  };
387 
388  Real _ybound[2] = {
389  src.pos[_Y] + abs(ty / sqrtY * z),
390  src.pos[_Y] - abs(ty / sqrtY * z)
391  };
392 
393  Real Xbound[2] = {
394  floor((_xbound[_X] - x) / ppX) + 1,
395  floor((_xbound[_Y] - x) / ppX) + 1
396  };
397 
398  Real Ybound[2] = {
399  pnY - floor((_ybound[_Y] - y) / ppY),
400  pnY - floor((_ybound[_X] - y) / ppY)
401  };
402 
403  if (Xbound[_X] > pnX) Xbound[_X] = pnX;
404  if (Xbound[_Y] < 0) Xbound[_Y] = 0;
405  if (Ybound[_X] > pnY) Ybound[_X] = pnY;
406  if (Ybound[_Y] < 0) Ybound[_Y] = 0;
407 
408 
409  for (int yytr = Ybound[_Y]; yytr < Ybound[_X]; ++yytr)
410  {
411  int offset = yytr * pnX;
412  Real yyy = y + ((pnY - yytr + offsetY) * ppY);
413 
414  Real range_x[2] = {
415  src.pos[_X] + abs(tx / sqrtX * sqrt((yyy - src.pos[_Y]) * (yyy - src.pos[_Y]) + zz)),
416  src.pos[_X] - abs(tx / sqrtX * sqrt((yyy - src.pos[_Y]) * (yyy - src.pos[_Y]) + zz))
417  };
418 
419  for (int xxtr = Xbound[_Y]; xxtr < Xbound[_X]; ++xxtr)
420  {
421  Real xxx = x + ((xxtr - 1 + offsetX) * ppX);
422  Real r = sqrt((xxx - src.pos[_X]) * (xxx - src.pos[_X]) + (yyy - src.pos[_Y]) * (yyy - src.pos[_Y]) + zz);
423  Real range_y[2] = {
424  src.pos[_Y] + abs(ty / sqrtY * sqrt((xxx - src.pos[_X]) * (xxx - src.pos[_X]) + zz)),
425  src.pos[_Y] - abs(ty / sqrtY * sqrt((xxx - src.pos[_X]) * (xxx - src.pos[_X]) + zz))
426  };
427 
428  if (((xxx < range_x[_X]) && (xxx > range_x[_Y])) && ((yyy < range_y[_X]) && (yyy > range_y[_Y]))) {
429  Real kr = k * r;
430  Real operand = lambda * r * r;
431  Real res_real = (ampZ * sin(kr)) / operand;
432  Real res_imag = (-ampZ * cos(kr)) / operand;
433 #ifdef _OPENMP
434 #pragma omp atomic
435  dst[offset + xxtr][_RE] += res_real;
436 #endif
437 #ifdef _OPENMP
438 #pragma omp atomic
439  dst[offset + xxtr][_IM] += res_imag;
440 #endif
441  }
442  }
443  }
444 }
445 
446 void ophGen::Fresnel_Diffraction(Point src, Complex<Real> *dst, Real lambda, Real distance, Real amplitude)
447 {
448  OphConfig *pConfig = &context_;
449  const int pnX = pConfig->pixel_number[_X];
450  const int pnY = pConfig->pixel_number[_Y];
451  const Real ppX = pConfig->pixel_pitch[_X];
452  const Real ppY = pConfig->pixel_pitch[_Y];
453  const Real ssX = pConfig->ss[_X] = pnX * ppX;
454  const Real ssY = pConfig->ss[_Y] = pnY * ppY;
455  const int offsetX = pConfig->offset[_X];
456  const int offsetY = pConfig->offset[_Y];
457  const Real k = (2 * M_PI) / lambda;
458 
459  // for performance
460  Real x = -ssX / 2;
461  Real y = -ssY / 2;
462  Real z = src.pos[_Z] + distance;
463  Real zz = z * z;
464  Real operand = lambda * z;
465 
466  Real _xbound[2] = {
467  src.pos[_X] + abs(operand / (2 * ppX)),
468  src.pos[_X] - abs(operand / (2 * ppX))
469  };
470 
471  Real _ybound[2] = {
472  src.pos[_Y] + abs(operand / (2 * ppY)),
473  src.pos[_Y] - abs(operand / (2 * ppY))
474  };
475 
476  Real Xbound[2] = {
477  floor((_xbound[_X] - x) / ppX) + 1,
478  floor((_xbound[_Y] - x) / ppX) + 1
479  };
480 
481  Real Ybound[2] = {
482  pnY - floor((_ybound[_Y] - y) / ppY),
483  pnY - floor((_ybound[_X] - y) / ppY)
484  };
485 
486  if (Xbound[_X] > pnX) Xbound[_X] = pnX;
487  if (Xbound[_Y] < 0) Xbound[_Y] = 0;
488  if (Ybound[_X] > pnY) Ybound[_X] = pnY;
489  if (Ybound[_Y] < 0) Ybound[_Y] = 0;
490 
491  for (int yytr = Ybound[_Y]; yytr < Ybound[_X]; ++yytr)
492  {
493  Real yyy = (y + (pnY - yytr + offsetY) * ppY) - src.pos[_Y];
494  int offset = yytr * pnX;
495  for (int xxtr = Xbound[_Y]; xxtr < Xbound[_X]; ++xxtr)
496  {
497  Real xxx = (x + (xxtr - 1 + offsetX) * ppX) - src.pos[_X];
498  Real p = k * (xxx * xxx + yyy * yyy + 2 * zz) / (2 * z);
499 
500  Real res_real = amplitude * sin(p) / operand;
501  Real res_imag = amplitude * (-cos(p)) / operand;
502 
503 #ifdef _OPENMP
504 #pragma omp atomic
505 #endif
506  dst[offset + xxtr][_RE] += res_real;
507 #ifdef _OPENMP
508 #pragma omp atomic
509 #endif
510  dst[offset + xxtr][_IM] += res_imag;
511  }
512  }
513 }
514 
515 void ophGen::Fresnel_FFT(Complex<Real> *src, Complex<Real> *dst, Real lambda, Real waveRatio, Real distance)
516 {
517  OphConfig *pConfig = &context_;
518  const ResolutionConfig* rc = &resCfg;
519  const int pnX = pConfig->pixel_number[_X];
520  const int pnY = pConfig->pixel_number[_Y];
521  const Real ppX = pConfig->pixel_pitch[_X];
522  const Real ppY = pConfig->pixel_pitch[_Y];
523  const Real ssX = pConfig->ss[_X] = pnX * ppX;
524  const Real ssY = pConfig->ss[_Y] = pnY * ppY;
525  const Real ssX2 = ssX * 2;
526  const Real ssY2 = ssY * 2;
527  const Real k = (2 * M_PI) / lambda;
528 
529  Complex<Real>* in2x = new Complex<Real>[rc->double_size];
530  memset(in2x, 0, sizeof(Complex<Real>) * rc->double_size);
531 
532  uint idxIn = 0;
533  for (int i = 0; i < pnY; i++) {
534  for (int j = 0; j < pnX; j++) {
535  in2x[(i + rc->half_pnY) * rc->double_pnX+ (j + rc->half_pnX)] = src[idxIn++];
536  }
537  }
538 
539  Complex<Real>* temp1 = new Complex<Real>[rc->double_size];
540 
541  fft2({ rc->double_pnX, rc->double_pnY}, in2x, OPH_FORWARD, OPH_ESTIMATE); // fft spatial domain
542  fft2(in2x, temp1, rc->double_pnX, rc->double_pnY, OPH_FORWARD, false); //
543 
544  Real* fx = new Real[rc->double_size];
545  Real* fy = new Real[rc->double_size];
546 
547  /*
548  for (int i = 0; i < rc->double_pnY; i++) {
549  uint k = i * rc->double_pnX;
550  Real val = (-rc->pnY + i) / ssY2;
551  memset(&fy[k], val, rc->double_pnX);
552  for (int j = 0; j < rc->double_pnX; j++) {
553  fx[k + j] = (-rc->pnX + j) / ssX2;
554  }
555  }
556  */
557 
558  uint i = 0;
559 
560  for (int idxFy = -pnY; idxFy < pnY; idxFy++) {
561  //fy[i] = idxFy / ssY2;
562  memset(&fy[i], idxFy / ssY2, rc->double_pnX);
563 
564  for (int idxFx = -pnX; idxFx < pnX; idxFx++) {
565  fx[i] = idxFx / ssX2;
566  i++;
567  }
568  }
569 
570  Complex<Real>* prop = new Complex<Real>[rc->double_size];
571  memset(prop, 0, sizeof(Complex<Real>) * rc->double_size);
572 
573  Real sqrtPart;
574 
575  Complex<Real>* temp2 = new Complex<Real>[rc->double_size];
576  Real lambda_square = lambda * lambda;
577  Real tmp = 2 * M_PI * distance;
578 
579  for (int i = 0; i < rc->double_size; i++) {
580  sqrtPart = sqrt(1 / lambda_square - fx[i] * fx[i] - fy[i] * fy[i]);
581  prop[i][_IM] = tmp;
582  prop[i][_IM] *= sqrtPart;
583  temp2[i] = temp1[i] * exp(prop[i]);
584  }
585 
586  Complex<Real>* temp3 = new Complex<Real>[rc->double_size];
587  fft2({ rc->double_pnX, rc->double_pnY}, temp2, OPH_BACKWARD, OPH_ESTIMATE);
588  fft2(temp2, temp3, rc->double_pnX, rc->double_pnY, OPH_BACKWARD, false);
589 
590  uint idxOut = 0;
591  for (int i = 0; i < pnY; i++) {
592  for (int j = 0; j < pnX; j++) {
593  dst[idxOut++] = temp3[(i + rc->half_pnY) * rc->double_pnX+ (j + rc->half_pnX)];
594  }
595  }
596 
597  delete[] in2x;
598  delete[] temp1;
599  delete[] fx;
600  delete[] fy;
601  delete[] prop;
602  delete[] temp2;
603  delete[] temp3;
604 }
605 
606 void ophGen::AngularSpectrumMethod(Complex<Real> *src, Complex<Real> *dst, Real lambda, Real distance)
607 {
608  const int pnX = context_.pixel_number[_X];
609  const int pnY = context_.pixel_number[_Y];
610  const int N = pnX * pnY;
611  const Real ppX = context_.pixel_pitch[_X];
612  const Real ppY = context_.pixel_pitch[_Y];
613  const Real ssX = context_.ss[_X] = pnX * ppX;
614  const Real ssY = context_.ss[_Y] = pnY * ppY;
615 
616  Real dfx = 1 / ssX;
617  Real dfy = 1 / ssY;
618 
619  Real k = context_.k = (2 * M_PI / lambda);
620  Real kk = k * k;
621  Real kd = k * distance;
622  Real fx = -1 / (ppX * 2);
623  Real fy = 1 / (ppY * 2);
624 
625 #ifdef _OPENMP
626 #pragma omp parallel for firstprivate(pnX, dfx, dfy, lambda, kd, kk)
627 #endif
628  for (int i = 0; i < N; i++)
629  {
630  Real x = i % pnX;
631  Real y = i / pnX;
632 
633  Real fxx = fx + dfx * x;
634  Real fyy = fy - dfy - dfy * y;
635 
636  Real fxxx = lambda * fxx;
637  Real fyyy = lambda * fyy;
638 
639  Real sval = sqrt(1 - (fxxx * fxxx) - (fyyy * fyyy));
640  sval = sval * kd;
641  Complex<Real> kernel(0, sval);
642  kernel.exp();
643 
644  bool prop_mask = ((fxx * fxx + fyy * fyy) < kk) ? true : false;
645 
646  Complex<Real> u_frequency;
647  if (prop_mask) {
648  u_frequency = kernel * src[i];
649  dst[i][_RE] += u_frequency[_RE];
650  dst[i][_IM] += u_frequency[_IM];
651  }
652  }
653 }
654 
655 void ophGen::conv_fft2(Complex<Real>* src1, Complex<Real>* src2, Complex<Real>* dst, ivec2 size)
656 {
657  int N = size[_X] * size[_Y];
658  if (N <= 0) return;
659 
660  Complex<Real>* src1FT = new Complex<Real>[N];
661  Complex<Real>* src2FT = new Complex<Real>[N];
662  Complex<Real>* dstFT = new Complex<Real>[N];
663 
664  fft2(src1, src1FT, size[_X], size[_Y], OPH_FORWARD, (bool)OPH_ESTIMATE);
665  fft2(src2, src2FT, size[_X], size[_Y], OPH_FORWARD, (bool)OPH_ESTIMATE);
666 
667  for (int i = 0; i < N; i++)
668  dstFT[i] = src1FT[i] * src2FT[i];
669 
670  fft2(dstFT, dst, size[_X], size[_Y], OPH_BACKWARD, (bool)OPH_ESTIMATE);
671 
672  delete[] src1FT;
673  delete[] src2FT;
674  delete[] dstFT;
675 }
676 
677 void ophGen::normalize(void)
678 {
679  const int pnX = context_.pixel_number[_X];
680  const int pnY = context_.pixel_number[_Y];
681  const uint nWave = context_.waveNum;
682  const long long int N = pnX * pnY;
683 
684  Real minVal = MAX_DOUBLE, maxVal = MIN_DOUBLE, gap = 0;
685  for (uint ch = 0; ch < nWave; ch++)
686  {
687  Real minTmp = minOfArr(m_lpEncoded[ch], N);
688  Real maxTmp = maxOfArr(m_lpEncoded[ch], N);
689 
690  if (minVal > minTmp)
691  minVal = minTmp;
692  if (maxVal < maxTmp)
693  maxVal = maxTmp;
694  }
695 
696  gap = maxVal - minVal;
697  for (uint ch = 0; ch < nWave; ch++)
698  {
699 #ifdef _OPENMP
700 #pragma omp parallel for firstprivate(minVal, gap, pnX)
701 #endif
702  for (int j = 0; j < pnY; j++) {
703  int start = j * pnX;
704  for (int i = 0; i < pnX; i++) {
705  int idx = start + i;
706  m_lpNormalized[ch][idx] = (((m_lpEncoded[ch][idx] - minVal) / gap) * 255 + 0.5);
707  }
708  }
709  }
710 }
711 
712 bool ophGen::save(const char* fname, uint8_t bitsperpixel, uchar* src, uint px, uint py)
713 {
714  bool bOK = false;
715 
716  if (fname == nullptr) return bOK;
717 
718  uchar* source = src;
719  bool bAlloc = false;
720  const uint nChannel = context_.waveNum;
721 
722  ivec2 p(px, py);
723  if (px == 0 && py == 0)
725 
726 
727  std::string file = fname;
728  std::replace(file.begin(), file.end(), '\\', '/');
729 
730  // split path
731  std::vector<std::string> components;
732  std::stringstream ss(file);
733  std::string item;
734  char token = '/';
735 
736  while (std::getline(ss, item, token)) {
737  components.push_back(item);
738  }
739 
740  std::string dir;
741 
742  for (size_t i = 0; i < components.size() - 1; i++)
743  {
744  dir += components[i];
745  dir += "/";
746  }
747 
748  std::string filename = components[components.size() - 1];
749 
750  // find extension
751  bool hasExt;
752  size_t ext_pos = file.rfind(".");
753  hasExt = (ext_pos == string::npos) ? false : true;
754 
755  if (!hasExt)
756  filename.append(".bmp");
757 
758  std::string fullpath = dir + filename;
759 
760  if (src == nullptr) {
761  if (nChannel == 1) {
762  source = m_lpNormalized[0];
763  saveAsImg(fullpath.c_str(), bitsperpixel, source, p[_X], p[_Y]);
764  }
765  else if (nChannel == 3) {
766  if (imgCfg.merge) {
767  uint nSize = (((p[_X] * bitsperpixel >> 3) + 3) & ~3) * p[_Y];
768  source = new uchar[nSize];
769  bAlloc = true;
770  for (uint i = 0; i < nChannel; i++) {
771  mergeColor(i, p[_X], p[_Y], m_lpNormalized[i], source);
772  }
773 
774  saveAsImg(fullpath.c_str(), bitsperpixel, source, p[_X], p[_Y]);
775  if (bAlloc) delete[] source;
776  }
777  else {
778  for (uint i = 0; i < nChannel; i++) {
779  char path[FILENAME_MAX] = { 0, };
780  sprintf(path, "%s%d_%s", dir.c_str(), i, filename.c_str());
781  source = m_lpNormalized[i];
782  saveAsImg(path, bitsperpixel / nChannel, source, p[_X], p[_Y]);
783  }
784  }
785  }
786  else return false;
787  }
788  else
789  saveAsImg(fullpath.c_str(), bitsperpixel, source, p[_X], p[_Y]);
790 
791  return true;
792 }
793 
794 void* ophGen::load(const char * fname)
795 {
796  if (checkExtension(fname, ".bmp")) {
797  return Openholo::loadAsImg(fname);
798  }
799  else { // when extension is not .bmp
800  return nullptr;
801  }
802 
803  return nullptr;
804 }
805 
806 bool ophGen::loadAsOhc(const char * fname)
807 {
808  if (!Openholo::loadAsOhc(fname)) return false;
809 
810  const uint nChannel = context_.waveNum;
811  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
812 
813  m_lpEncoded = new Real*[nChannel];
814  m_lpNormalized = new uchar*[nChannel];
815  for (uint ch = 0; ch < nChannel; ch++) {
816  m_lpEncoded[ch] = new Real[pnXY];
817  memset(m_lpEncoded[ch], 0, sizeof(Real) * pnXY);
818  m_lpNormalized[ch] = new uchar[pnXY];
819  memset(m_lpNormalized[ch], 0, sizeof(uchar) * pnXY);
820  }
821  return true;
822 }
823 
824 void ophGen::resetBuffer()
825 {
827  int N2 = m_vecEncodeSize[_X] * m_vecEncodeSize[_Y];
828 
829  for (uint ch = 0; ch < context_.waveNum; ch++) {
830  if (complex_H[ch])
831  memset(complex_H[ch], 0., sizeof(Complex<Real>) * N);
832  if (m_lpEncoded[ch])
833  memset(m_lpEncoded[ch], 0., sizeof(Real) * N2);
834  if (m_lpNormalized[ch])
835  memset(m_lpNormalized[ch], 0, sizeof(uchar) * N2);
836  }
837 }
838 
839 void ophGen::encoding(unsigned int ENCODE_FLAG)
840 {
841  auto begin = CUR_TIME;
842 
843  // func pointer
844  void (ophGen::*encodeFunc) (Complex<Real>*, Real*, const int) = nullptr;
845 
846  Complex<Real> *holo = nullptr;
847  Real *encoded = nullptr;
848  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
850 
851  switch (ENCODE_FLAG)
852  {
853  case ENCODE_PHASE: encodeFunc = &ophGen::Phase; LOG("ENCODE_PHASE\n"); break;
854  case ENCODE_AMPLITUDE: encodeFunc = &ophGen::Amplitude; LOG("ENCODE_AMPLITUDE\n"); break;
855  case ENCODE_REAL: encodeFunc = &ophGen::RealPart; LOG("ENCODE_REAL\n"); break;
856  case ENCODE_IMAGINARY: encodeFunc = &ophGen::ImaginaryPart; LOG("ENCODE_IMAGINARY\n"); break;
857  case ENCODE_SIMPLENI: encodeFunc = &ophGen::SimpleNI; LOG("ENCODE_SIMPLENI\n"); break;
858  case ENCODE_BURCKHARDT: encodeFunc = &ophGen::Burckhardt; LOG("ENCODE_BURCKHARDT\n"); break;
859  case ENCODE_TWOPHASE: encodeFunc = &ophGen::TwoPhase; LOG("ENCODE_TWOPHASE\n"); break;
860  default:
861  LOG("<FAILED> WRONG PARAMETERS.\n");
862  LOG("%s : %.5lf (sec)\n", __FUNCTION__, ELAPSED_TIME(begin, CUR_TIME));
863  return;
864  }
865 
866  for (uint ch = 0; ch < context_.waveNum; ch++) {
867  holo = complex_H[ch];
868  encoded = m_lpEncoded[ch];
869  (this->*encodeFunc)(holo, encoded, m_vecEncodeSize[_X] * m_vecEncodeSize[_Y]);
870  }
871  LOG("%s : %.5lf (sec)\n", __FUNCTION__, ELAPSED_TIME(begin, CUR_TIME));
872 }
873 
874 //template <typename T>
875 //void ophGen::encoding(unsigned int ENCODE_FLAG, Complex<T>* holo, T* encoded)
876 void ophGen::encoding(unsigned int ENCODE_FLAG, Complex<Real>* holo, Real* encoded)
877 {
878  auto begin = CUR_TIME;
879  // func pointer
880  void (ophGen::*encodeFunc) (Complex<Real>*, Real*, const int) = nullptr;
881 
882  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
884 
885  switch (ENCODE_FLAG)
886  {
887  case ENCODE_PHASE: encodeFunc = &ophGen::Phase; LOG("ENCODE_PHASE\n"); break;
888  case ENCODE_AMPLITUDE: encodeFunc = &ophGen::Amplitude; LOG("ENCODE_AMPLITUDE\n"); break;
889  case ENCODE_REAL: encodeFunc = &ophGen::RealPart; LOG("ENCODE_REAL\n"); break;
890  case ENCODE_IMAGINARY: encodeFunc = &ophGen::ImaginaryPart; LOG("ENCODE_IMAGINARY\n"); break;
891  case ENCODE_SIMPLENI: encodeFunc = &ophGen::SimpleNI; LOG("ENCODE_SIMPLENI\n"); break;
892  case ENCODE_BURCKHARDT: encodeFunc = &ophGen::Burckhardt; LOG("ENCODE_BURCKHARDT\n"); break;
893  case ENCODE_TWOPHASE: encodeFunc = &ophGen::TwoPhase; LOG("ENCODE_TWOPHASE\n"); break;
894  default: LOG("<FAILED> WRONG PARAMETERS.\n"); return;
895  }
896 
897  if (holo == nullptr) holo = complex_H[0];
898  if (encoded == nullptr) encoded = m_lpEncoded[0];
899 
900  (this->*encodeFunc)(holo, encoded, m_vecEncodeSize[_X] * m_vecEncodeSize[_Y]);
901 
902  LOG("%s : %.5lf (sec)\n", __FUNCTION__, ELAPSED_TIME(begin, CUR_TIME));
903 }
904 
905 void ophGen::encoding(unsigned int ENCODE_FLAG, unsigned int passband, Complex<Real>* holo, Real* encoded)
906 {
907  holo == nullptr ? holo = *complex_H : holo;
908  encoded == nullptr ? encoded = *m_lpEncoded : encoded;
909 
911  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
912  const uint nChannel = context_.waveNum;
913 
914  for (uint ch = 0; ch < nChannel; ch++) {
915  /* initialize */
916  long long int m_vecEncodeSize = pnXY;
917  if (m_lpEncoded[ch] != nullptr) delete[] m_lpEncoded[ch];
918  m_lpEncoded[ch] = new Real[m_vecEncodeSize];
919  memset(m_lpEncoded[ch], 0, sizeof(Real) * m_vecEncodeSize);
920 
921  if (m_lpNormalized[ch] != nullptr) delete[] m_lpNormalized[ch];
923  memset(m_lpNormalized[ch], 0, sizeof(uchar) * m_vecEncodeSize);
924 
925  switch (ENCODE_FLAG)
926  {
927  case ENCODE_SSB:
928  LOG("ENCODE_SSB");
929  singleSideBand((holo), m_lpEncoded[ch], context_.pixel_number, passband);
930  break;
931  case ENCODE_OFFSSB:
932  {
933  LOG("ENCODE_OFFSSB");
934  Complex<Real> *tmp = new Complex<Real>[pnXY];
935  memcpy(tmp, holo, sizeof(Complex<Real>) * pnXY);
936  freqShift(tmp, tmp, context_.pixel_number, 0, 100);
937  singleSideBand(tmp, m_lpEncoded[ch], context_.pixel_number, passband);
938  delete[] tmp;
939  break;
940  }
941  default:
942  LOG("<FAILED> WRONG PARAMETERS.\n");
943  return;
944  }
945  }
946 }
947 
948 void ophGen::encoding(unsigned int BIN_ENCODE_FLAG, unsigned int ENCODE_FLAG, Real threshold, Complex<Real>* holo, Real* encoded)
949 {
950  auto begin = CUR_TIME;
951  const uint nChannel = context_.waveNum;
952  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
953 
954  switch (BIN_ENCODE_FLAG) {
955  case ENCODE_SIMPLEBINARY:
956  LOG("ENCODE_SIMPLEBINARY\n");
957  if (holo == nullptr || encoded == nullptr)
958  for (uint ch = 0; ch < nChannel; ch++) {
959  binarization(complex_H[ch], m_lpEncoded[ch], pnXY, ENCODE_FLAG, threshold);
960  }
961  else
962  binarization(holo, encoded, pnXY, ENCODE_FLAG, threshold);
963  break;
964  case ENCODE_EDBINARY:
965  LOG("ENCODE_EDBINARY\n");
966  if (ENCODE_FLAG != ENCODE_REAL) {
967  LOG("<FAILED> WRONG PARAMETERS : %d\n", ENCODE_FLAG);
968  return;
969  }
970  if (holo == nullptr || encoded == nullptr)
971  for (uint ch = 0; ch < nChannel; ch++) {
973  }
974  else
975  binaryErrorDiffusion(holo, encoded, pnXY, FLOYD_STEINBERG, threshold);
976  break;
977  default:
978  LOG("<FAILED> WRONG PARAMETERS.\n");
979  return;
980  }
981 
982  LOG("%s : %.5lf (sec)\n", __FUNCTION__, ELAPSED_TIME(begin, CUR_TIME));
983 }
984 
985 void ophGen::encoding()
986 {
987  const long long int pnXY = context_.pixel_number[_X] * context_.pixel_number[_Y];
988  const uint nChannel = context_.waveNum;
989 
991  else if (ENCODE_METHOD == ENCODE_TWOPHASE) m_vecEncodeSize[_X] *= 2;
992 
993  for (uint ch = 0; ch < nChannel; ch++) {
994  /* initialize */
995  if (m_lpEncoded[ch] != nullptr) delete[] m_lpEncoded[ch];
997  memset(m_lpEncoded[ch], 0, sizeof(Real) * m_vecEncodeSize[_X] * m_vecEncodeSize[_Y]);
998 
999  if (m_lpNormalized[ch] != nullptr) delete[] m_lpNormalized[ch];
1001  memset(m_lpNormalized[ch], 0, sizeof(uchar) * m_vecEncodeSize[_X] * m_vecEncodeSize[_Y]);
1002 
1003 
1004  switch (ENCODE_METHOD)
1005  {
1006  case ENCODE_SIMPLENI:
1007  LOG("ENCODE_SIMPLENI\n");
1008  SimpleNI(complex_H[ch], m_lpEncoded[ch], pnXY);
1009  break;
1010  case ENCODE_REAL:
1011  LOG("ENCODE_REAL\n");
1012  realPart<Real>(complex_H[ch], m_lpEncoded[ch], pnXY);
1013  break;
1014  case ENCODE_BURCKHARDT:
1015  LOG("ENCODE_BURCKHARDT\n");
1016  Burckhardt(complex_H[ch], m_lpEncoded[ch], pnXY);
1017  break;
1018  case ENCODE_TWOPHASE:
1019  LOG("ENCODE_TWOPHASE\n");
1020  TwoPhase(complex_H[ch], m_lpEncoded[ch], pnXY);
1021  break;
1022  case ENCODE_PHASE:
1023  LOG("ENCODE_PHASE\n");
1024  Phase(complex_H[ch], m_lpEncoded[ch], pnXY);
1025  break;
1026  case ENCODE_AMPLITUDE:
1027  LOG("ENCODE_AMPLITUDE\n");
1028  getAmplitude(complex_H[ch], m_lpEncoded[ch], pnXY);
1029  break;
1030  case ENCODE_SSB:
1031  LOG("ENCODE_SSB\n");
1033  break;
1034  case ENCODE_OFFSSB:
1035  LOG("ENCODE_OFFSSB\n");
1036  freqShift(complex_H[ch], complex_H[ch], context_.pixel_number, 0, 100);
1038  break;
1039  default:
1040  LOG("<FAILED> WRONG PARAMETERS.\n");
1041  return;
1042  }
1043  }
1044 }
1045 
1046 void ophGen::singleSideBand(Complex<Real>* src, Real* dst, const ivec2 holosize, int SSB_PASSBAND)
1047 {
1048  const int nX = holosize[_X];
1049  const int nY = holosize[_Y];
1050  const int half_nX = nX >> 1;
1051  const int half_nY = nY >> 1;
1052 
1053  const long long int N = nX * nY;
1054  const long long int half_N = N >> 1;
1055 
1056  Complex<Real>* AS = new Complex<Real>[N];
1057  //fft2(holosize, holo, OPH_FORWARD, OPH_ESTIMATE);
1058  fft2(src, AS, nX, nY, OPH_FORWARD, false);
1059  //fftExecute(temp);
1060 
1061 
1062  switch (SSB_PASSBAND)
1063  {
1064  case SSB_LEFT:
1065  for (int i = 0; i < nY; i++)
1066  {
1067  int k = i * nX;
1068  for (int j = half_nX; j < nX; j++)
1069  {
1070  AS[k + j] = 0;
1071  }
1072  }
1073  break;
1074  case SSB_RIGHT:
1075  for (int i = 0; i < nY; i++)
1076  {
1077  int k = i * nX;
1078  for (int j = 0; j < half_nX; j++)
1079  {
1080  AS[k + j] = 0;
1081  }
1082  }
1083  break;
1084  case SSB_TOP:
1085  memset(&AS[half_N], 0, sizeof(Complex<Real>) * half_N);
1086  break;
1087 
1088  case SSB_BOTTOM:
1089  memset(&AS[0], 0, sizeof(Complex<Real>) * half_N);
1090  break;
1091  }
1092 
1093  Complex<Real>* filtered = new Complex<Real>[N];
1094  //fft2(holosize, AS, OPH_BACKWARD, OPH_ESTIMATE);
1095  fft2(AS, filtered, nX, nY, OPH_BACKWARD, false);
1096 
1097  //fftExecute(filtered);
1098 
1099 
1100  Real* realFiltered = new Real[N];
1101  oph::realPart<Real>(filtered, realFiltered, N);
1102 
1103  oph::normalize(realFiltered, dst, N);
1104 
1105  delete[] AS;
1106  delete[] filtered;
1107  delete[] realFiltered;
1108 }
1109 
1110 void ophGen::freqShift(Complex<Real>* src, Complex<Real>* dst, const ivec2 holosize, int shift_x, int shift_y)
1111 {
1112  int N = holosize[_X] * holosize[_Y];
1113 
1114  Complex<Real>* AS = new Complex<Real>[N];
1115  //fft2(holosize, src, OPH_FORWARD, OPH_ESTIMATE);
1116  fft2(src, AS, holosize[_X], holosize[_Y], OPH_FORWARD);
1117  //fftExecute(AS);
1118 
1119  Complex<Real>* shifted = new Complex<Real>[N];
1120  circShift<Complex<Real>>(AS, shifted, shift_x, shift_y, holosize.v[_X], holosize.v[_Y]);
1121 
1122  //fft2(holosize, shifted, OPH_BACKWARD, OPH_ESTIMATE);
1123  fft2(shifted, dst, holosize[_X], holosize[_Y], OPH_BACKWARD);
1124  //fftExecute(dst);
1125 
1126  delete[] AS;
1127  delete[] shifted;
1128 }
1129 
1130 
1131 bool ophGen::saveRefImages(char* fnameW, char* fnameWC, char* fnameAS, char* fnameSSB, char* fnameHP, char* fnameFreq, char* fnameReal, char* fnameBin, char* fnameReconBin, char* fnameReconErr, char* fnameReconNo)
1132 {
1133  ivec2 holosize = context_.pixel_number;
1134  int nx = holosize[_X], ny = holosize[_Y];
1135  int ss = nx*ny;
1136 
1137  Real* temp1 = new Real[ss];
1138  uchar* temp2 = new uchar[ss];
1139 
1140  oph::normalize(weight, temp2, nx, ny);
1141  saveAsImg(fnameW, 8, temp2, nx, ny);
1142  cout << "W saved" << endl;
1143 
1144  oph::absCplxArr<Real>(weightC, temp1, ss);
1145  oph::normalize(temp1, temp2, nx, ny);
1146  saveAsImg(fnameWC, 8, temp2, nx, ny);
1147  cout << "WC saved" << endl;
1148 
1149  oph::absCplxArr<Real>(AS, temp1, ss);
1150  oph::normalize(temp1, temp2, nx, ny);
1151  saveAsImg(fnameAS, 8, temp2, nx, ny);
1152  cout << "AS saved" << endl;
1153 
1154  oph::normalize(maskSSB, temp2, nx, ny);
1155  saveAsImg(fnameSSB, 8, temp2, nx, ny);
1156  cout << "SSB saved" << endl;
1157 
1158  oph::normalize(maskHP, temp2, nx, ny);
1159  saveAsImg(fnameHP, 8, temp2, nx, ny);
1160  cout << "HP saved" << endl;
1161 
1162  oph::absCplxArr<Real>(freqW, temp1, ss);
1163  oph::normalize(temp1, temp2, nx, ny);
1164  saveAsImg(fnameFreq, 8, temp2, nx, ny);
1165  cout << "Freq saved" << endl;
1166 
1167  oph::normalize(realEnc, temp2, nx, ny);
1168  saveAsImg(fnameReal, 8, temp2, nx, ny);
1169  cout << "Real saved" << endl;
1170 
1171  oph::normalize(binary, temp2, nx, ny);
1172  saveAsImg(fnameBin, 8, temp2, nx, ny);
1173  cout << "Bin saved" << endl;
1174 
1175 
1176  Complex<Real>* temp = new Complex<Real>[ss];
1177  for (int i = 0; i < ss; i++) {
1178  temp[i][_RE] = binary[i];
1179  temp[i][_IM] = 0;
1180  }
1181  fft2(ivec2(nx, ny), temp, OPH_FORWARD);
1182  fft2(temp, temp, nx, ny, OPH_FORWARD);
1183  for (int i = 0; i < ss; i++) {
1184  temp[i][_RE] *= maskSSB[i];
1185  temp[i][_IM] *= maskSSB[i];
1186  }
1187  fft2(ivec2(nx, ny), temp, OPH_BACKWARD);
1188  fft2(temp, temp, nx, ny, OPH_BACKWARD);
1189 
1190  Complex<Real>* reconBin = new Complex<Real>[ss];
1191  memset(reconBin, 0, sizeof(Complex<Real>) * ss);
1192  fresnelPropagation(temp, reconBin, 0.001, 0);
1193 
1194  oph::absCplxArr<Real>(reconBin, temp1, ss);
1195  for (int i = 0; i < ss; i++) {
1196  temp1[i] = temp1[i] * temp1[i];
1197  }
1198  oph::normalize(temp1, temp2, nx, ny);
1199  saveAsImg(fnameReconBin, 8, temp2, nx, ny);
1200  cout << "recon bin saved" << endl;
1201 
1202 
1203  temp = new Complex<Real>[ss];
1204  for (int i = 0; i < ss; i++) {
1205  temp[i][_RE] = m_lpEncoded[0][i];
1206  temp[i][_IM] = 0;
1207  }
1208  fft2(ivec2(nx, ny), temp, OPH_FORWARD);
1209  fft2(temp, temp, holosize[_X], holosize[_Y], OPH_FORWARD);
1210  for (int i = 0; i < ss; i++) {
1211  temp[i][_RE] *= maskHP[i];
1212  temp[i][_IM] *= maskHP[i];
1213  }
1214  fft2(ivec2(nx, ny), temp, OPH_BACKWARD);
1215  fft2(temp, temp, nx, ny, OPH_BACKWARD);
1216 
1217  reconBin = new Complex<Real>[ss];
1218  fresnelPropagation(temp, reconBin, 0.001, 0);
1219 
1220  oph::absCplxArr<Real>(reconBin, temp1, ss);
1221  for (int i = 0; i < ss; i++) {
1222  temp1[i] = temp1[i] * temp1[i];
1223  }
1224  oph::normalize(temp1, temp2, nx, ny);
1225  saveAsImg(fnameReconErr, 8, temp2, nx, ny);
1226  cout << "recon error saved" << endl;
1227 
1228 
1229 
1230 
1231  temp = new Complex<Real>[ss];
1232  for (int i = 0; i < ss; i++) {
1233  temp[i][_RE] = normalized[i][_RE];
1234  temp[i][_IM] = normalized[i][_IM];
1235  }
1236  fft2(ivec2(nx, ny), temp, OPH_FORWARD);
1237  fft2(temp, temp, holosize[_X], holosize[_Y], OPH_FORWARD);
1238  for (int i = 0; i < ss; i++) {
1239  temp[i][_RE] *= maskSSB[i];
1240  temp[i][_IM] *= maskSSB[i];
1241  }
1242  fft2(ivec2(nx, ny), temp, OPH_BACKWARD);
1243  fft2(temp, temp, nx, ny, OPH_BACKWARD);
1244 
1245  reconBin = new Complex<Real>[ss];
1246  fresnelPropagation(temp, reconBin, 0.001, 0);
1247 
1248  oph::absCplxArr<Real>(reconBin, temp1, ss);
1249  for (int i = 0; i < ss; i++) {
1250  temp1[i] = temp1[i] * temp1[i];
1251  }
1252  oph::normalize(temp1, temp2, nx, ny);
1253  saveAsImg(fnameReconNo, 8, temp2, nx, ny);
1254 
1255 
1256  return true;
1257 }
1258 
1259 bool ophGen::binaryErrorDiffusion(Complex<Real>* holo, Real* encoded, const ivec2 holosize, const int type, Real threshold)
1260 {
1261 
1262  //cout << "\nin?" << endl;
1263  int ss = holosize[_X] * holosize[_Y];
1264  weight = new Real[ss];
1265  //cout << "?" << endl;
1266  weightC = new Complex<Real>[ss];
1267  //cout << "??" << endl;
1268  ivec2 Nw;
1269  memsetArr<Real>(weight, 0.0, 0, ss - 1);
1270  //cout << "???" << endl;
1271  if (!getWeightED(holosize, type, &Nw))
1272  return false;
1273  //cout << "1?" << endl;
1274  AS = new Complex<Real>[ss];
1275  fft2(ivec2(holosize[_X], holosize[_Y]), holo, OPH_FORWARD);
1276  fft2(holo, AS, holosize[_X], holosize[_Y], OPH_FORWARD);
1277  //cout << "2?" << endl;
1278  // SSB mask generation
1279  maskSSB = new Real[ss];
1280  for (int i = 0; i < ss; i++) {
1281  if (((Real)i / (Real)holosize[_X]) < ((Real)holosize[_Y] / 2.0))
1282  maskSSB[i] = 1;
1283  else
1284  maskSSB[i] = 0;
1285  AS[i] *= maskSSB[i];
1286  }
1287 
1288  //cout << "3?" << endl;
1289  Complex<Real>* filtered = new Complex<Real>[ss];
1290  fft2(ivec2(holosize[_X], holosize[_Y]), AS, OPH_BACKWARD);
1291  fft2(AS, filtered, holosize[_X], holosize[_Y], OPH_BACKWARD);
1292  //cout << "4?" << endl;
1293  normalized = new Complex<Real>[ss];
1294  oph::normalize(filtered, normalized, ss);
1295  LOG("normalize finishied..\n");
1296  if (encoded == nullptr)
1297  encoded = new Real[ss];
1298 
1299  shiftW(holosize);
1300  LOG("shiftW finishied..\n");
1301 
1302  // HP mask generation
1303  maskHP = new Real[ss];
1304  Real absFW;
1305  for (int i = 0; i < ss; i++) {
1306  oph::absCplx<Real>(freqW[i], absFW);
1307  if (((Real)i / (Real)holosize[_X]) < ((Real)holosize[_Y] / 2.0) && absFW < 0.6)
1308  maskHP[i] = 1;
1309  else
1310  maskHP[i] = 0;
1311  }
1312  //cout << "5?" << endl;
1313 
1314  // For checking
1315  binary = new Real[ss];
1316  realEnc = new Real[ss];
1317  for (int i = 0; i < ss; i++) {
1318  realEnc[i] = normalized[i][_RE];
1319  if (normalized[i][_RE] > threshold)
1320  binary[i] = 1;
1321  else
1322  binary[i] = 0;
1323  }
1324 
1325  Complex<Real>* toBeBin = new Complex<Real>[ss];
1326  for (int i = 0; i < ss; i++) {
1327  toBeBin[i] = normalized[i];
1328  }
1329  int ii, iii, jjj;
1330  int cx = (holosize[_X] + 1) / 2, cy = (holosize[_Y] + 1) / 2;
1331  Real error;
1332  for (int iy = 0; iy < holosize[_Y] - Nw[_Y]; iy++) {
1333  for (int ix = Nw[_X]; ix < holosize[_X] - Nw[_X]; ix++) {
1334 
1335  ii = ix + iy*holosize[_X];
1336  if (ix >= Nw[_X] && ix < (holosize[_X] - Nw[_X]) && iy < (holosize[_Y] - Nw[_Y])) {
1337 
1338  if (toBeBin[ii][_RE] > threshold)
1339  encoded[ii] = 1;
1340  else
1341  encoded[ii] = 0;
1342 
1343  error = toBeBin[ii][_RE] - encoded[ii];
1344 
1345  for (int iwy = 0; iwy < Nw[_Y] + 1; iwy++) {
1346  for (int iwx = -Nw[_X]; iwx < Nw[_X] + 1; iwx++) {
1347  iii = (ix + iwx) + (iy + iwy)*holosize[_X];
1348  jjj = (cx + iwx) + (cy + iwy)*holosize[_X];
1349 
1350  toBeBin[iii] += weightC[jjj] * error;
1351  }
1352  }
1353  }
1354  else {
1355  encoded[ii] = 0;
1356  }
1357  }
1358  }
1359  LOG("binary finishied..\n");
1360 
1361  return true;
1362 }
1363 
1364 
1365 bool ophGen::getWeightED(const ivec2 holosize, const int type, ivec2* pNw)
1366 {
1367 
1368  int cx = (holosize[_X] + 1) / 2;
1369  int cy = (holosize[_Y] + 1) / 2;
1370 
1371  ivec2 Nw;
1372 
1373  switch (type) {
1374  case FLOYD_STEINBERG:
1375  LOG("ERROR DIFFUSION : FLOYD_STEINBERG\n");
1376  weight[(cx + 1) + cy*holosize[_X]] = 7.0 / 16.0;
1377  weight[(cx - 1) + (cy + 1)*holosize[_X]] = 3.0 / 16.0;
1378  weight[(cx)+(cy + 1)*holosize[_X]] = 5.0 / 16.0;
1379  weight[(cx + 1) + (cy + 1)*holosize[_X]] = 1.0 / 16.0;
1380  Nw[_X] = 1; Nw[_Y] = 1;
1381  break;
1382  case SINGLE_RIGHT:
1383  LOG("ERROR DIFFUSION : SINGLE_RIGHT\n");
1384  weight[(cx + 1) + cy*holosize[_X]] = 1.0;
1385  Nw[_X] = 1; Nw[_Y] = 0;
1386  break;
1387  case SINGLE_DOWN:
1388  LOG("ERROR DIFFUSION : SINGLE_DOWN\n");
1389  weight[cx + (cy + 1)*holosize[_X]] = 1.0;
1390  Nw[_X] = 0; Nw[_Y] = 1;
1391  break;
1392  default:
1393  LOG("<FAILED> WRONG PARAMETERS.\n");
1394  return false;
1395  }
1396 
1397  *pNw = Nw;
1398  return true;
1399 
1400 }
1401 
1402 bool ophGen::shiftW(ivec2 holosize)
1403 {
1404  int ss = holosize[_X] * holosize[_Y];
1405 
1406  Complex<Real> term(0, 0);
1407  Complex<Real> temp(0, 0);
1408  int x, y;
1409  for (int i = 0; i < ss; i++) {
1410 
1411  x = i % holosize[_X] - holosize[_X] / 2; y = i / holosize[_X] - holosize[_Y];
1412  term[_IM] = 2.0 * M_PI*((1.0 / 4.0)*(Real)x + (0.0)*(Real)y);
1413  temp[_RE] = weight[i];
1414  weightC[i] = temp *exp(term);
1415  }
1416 
1417  freqW = new Complex<Real>[ss];
1418 
1419  fft2(ivec2(holosize[_X], holosize[_Y]), weightC, OPH_FORWARD);
1420  fft2(weightC, freqW, holosize[_X], holosize[_Y], OPH_FORWARD);
1421  for (int i = 0; i < ss; i++) {
1422  freqW[i][_RE] -= 1.0;
1423  }
1424  return true;
1425 
1426 }
1427 
1428 void ophGen::binarization(Complex<Real>* src, Real* dst, const int size, int ENCODE_FLAG, Real threshold)
1429 {
1430  oph::normalize(src, src, size);
1432 
1433 #ifdef _OPENMP
1434 #pragma omp parallel for firstprivate(threshold)
1435 #endif
1436  for (int i = 0; i < size; i++) {
1437  if (src[i][_RE] > threshold)
1438  dst[i] = 1;
1439  else
1440  dst[i] = 0;
1441  }
1442 }
1443 
1444 void ophGen::fresnelPropagation(OphConfig context, Complex<Real>* in, Complex<Real>* out, Real distance)
1445 {
1446  const int pnX = context.pixel_number[_X];
1447  const int pnY = context.pixel_number[_Y];
1448  const int pnX2 = pnX << 1;
1449  const int pnY2 = pnY << 1;
1450  const long long int pnXY = pnX * pnY;
1451  const long long int size = pnXY << 2;
1452 
1453  Complex<Real>* in2x = new Complex<Real>[size];
1454  Complex<Real> zero(0, 0);
1455  memset(in2x, 0, sizeof(Complex<Real>) * size);
1456 
1457  uint idxIn = 0;
1458  int beginY = pnY >> 1;
1459  int beginX = pnX >> 1;
1460  int endY = pnY + beginY;
1461  int endX = pnX + beginX;
1462 
1463  for (int idxnY = beginY; idxnY < endY; idxnY++) {
1464  for (int idxnX = beginX; idxnX < endX; idxnX++) {
1465  in2x[idxnY * pnX2 + idxnX] = in[idxIn++];
1466  }
1467  }
1468 
1469 
1470  Complex<Real>* temp1 = new Complex<Real>[size];
1471 
1472  fft2({ pnX2, pnY2 }, in2x, OPH_FORWARD, OPH_ESTIMATE);
1473  fft2(in2x, temp1, pnX2, pnY2, OPH_FORWARD);
1474  //fftExecute(temp1);
1475  Real* fx = new Real[size];
1476  Real* fy = new Real[size];
1477 
1478  uint i = 0;
1479  for (int idxFy = -pnY; idxFy < pnY; idxFy++) {
1480  for (int idxFx = -pnX; idxFx < pnX; idxFx++) {
1481  fx[i] = idxFx / (pnX2 * context.pixel_pitch[_X]);
1482  fy[i] = idxFy / (pnY2 * context.pixel_pitch[_Y]);
1483  i++;
1484  }
1485  }
1486 
1487  Complex<Real>* prop = new Complex<Real>[size];
1488  memsetArr<Complex<Real>>(prop, zero, 0, size - 1);
1489 
1490  Real sqrtPart;
1491 
1492  Complex<Real>* temp2 = new Complex<Real>[size];
1493 
1494  for (int i = 0; i < size; i++) {
1495  sqrtPart = sqrt(1 / (context.wave_length[0] * context.wave_length[0]) - fx[i] * fx[i] - fy[i] * fy[i]);
1496  prop[i][_IM] = 2 * M_PI * distance;
1497  prop[i][_IM] *= sqrtPart;
1498  temp2[i] = temp1[i] * exp(prop[i]);
1499  }
1500 
1501  Complex<Real>* temp3 = new Complex<Real>[size];
1502  fft2({ pnX2, pnY2 }, temp2, OPH_BACKWARD, OPH_ESTIMATE);
1503  fft2(temp2, temp3, pnX2, pnY2, OPH_BACKWARD);
1504 
1505  uint idxOut = 0;
1506 
1507  for (int idxNy = pnY / 2; idxNy < pnY + (pnY / 2); idxNy++) {
1508  for (int idxNx = pnX / 2; idxNx < pnX + (pnX / 2); idxNx++) {
1509  out[idxOut] = temp3[idxNy * pnX * 2 + idxNx];
1510  idxOut++;
1511  }
1512  }
1513 
1514  delete[] in2x;
1515  delete[] temp1;
1516  delete[] fx;
1517  delete[] fy;
1518  delete[] prop;
1519  delete[] temp2;
1520  delete[] temp3;
1521 }
1522 
1523 void ophGen::fresnelPropagation(Complex<Real>* in, Complex<Real>* out, Real distance, uint channel)
1524 {
1525  const int pnX = context_.pixel_number[_X];
1526  const int pnY = context_.pixel_number[_Y];
1527  const Real ppX = context_.pixel_pitch[_X];
1528  const Real ppY = context_.pixel_pitch[_Y];
1529  const int pnXY = pnX * pnY;
1530  const Real lambda = context_.wave_length[channel];
1531  const Real ssX = pnX * ppX * 2;
1532  const Real ssY = pnY * ppY * 2;
1533  const Real z = 2 * M_PI * distance;
1534  const Real v = 1 / (lambda * lambda);
1535  const int hpnX = pnX / 2;
1536  const int hpnY = pnY / 2;
1537  const int pnX2 = pnX * 2;
1538  const int pnY2 = pnY * 2;
1539 
1540  Complex<Real>* temp = new Complex<Real>[pnXY * 4];
1541  memset(temp, 0, sizeof(Complex<Real>) * pnXY * 4);
1542 
1543 #ifdef _OPENMP
1544 #pragma omp parallel for firstprivate(pnX, pnX2, hpnX, hpnY)
1545 #endif
1546  for (int i = 0; i < pnY; i++)
1547  {
1548  int src = pnX * i;
1549  int dst = pnX2 * (i + hpnY) + hpnX;
1550  memcpy(&temp[dst], &in[src], sizeof(Complex<Real>) * pnX);
1551  }
1552 
1553  fft2(temp, temp, pnX2, pnY2, OPH_FORWARD, false);
1554 
1555 #ifdef _OPENMP
1556 #pragma omp parallel for firstprivate(ssX, ssY, z, v)
1557 #endif
1558  for (int j = 0; j < pnY2; j++)
1559  {
1560  Real fy = (-pnY + j) / ssY;
1561  Real fyy = fy * fy;
1562  int iWidth = j * pnX2;
1563  for (int i = 0; i < pnX2; i++)
1564  {
1565  Real fx = (-pnX + i) / ssX;
1566  Real fxx = fx * fx;
1567 
1568  Real sqrtPart = sqrt(v - fxx - fyy);
1569  Complex<Real> prop(0, z * sqrtPart);
1570  temp[iWidth + i] *= prop.exp();
1571  }
1572  }
1573 
1574  fft2(temp, temp, pnX2, pnY2, OPH_BACKWARD, true);
1575 
1576 #ifdef _OPENMP
1577 #pragma omp parallel for firstprivate(pnX, pnX2, hpnX, hpnY)
1578 #endif
1579  for (int i = 0; i < pnY; i++)
1580  {
1581  int src = pnX2 * (i + hpnY) + hpnX;
1582  int dst = pnX * i;
1583  memcpy(&out[dst], &temp[src], sizeof(Complex<Real>) * pnX);
1584  }
1585  delete[] temp;
1586 
1587 }
1588 
1589 bool ophGen::Shift(Real x, Real y)
1590 {
1591  if (x == 0.0 && y == 0.0) return false;
1592 
1593  bool bAxisX = (x == 0.0) ? false : true;
1594  bool bAxisY = (y == 0.0) ? false : true;
1595  const uint nChannel = context_.waveNum;
1596  const int pnX = context_.pixel_number[_X];
1597  const int pnY = context_.pixel_number[_Y];
1598  const Real ppX = context_.pixel_pitch[_X];
1599  const Real ppY = context_.pixel_pitch[_Y];
1600  const vec2 ss = context_.ss;
1601  Real ppY2 = ppY * 2;
1602  Real ppX2 = ppX * 2;
1603  Real ssX = -ss[_X] / 2;
1604  Real ssY = -ss[_Y] / 2;
1605 
1606  Real *waveRatio = new Real[nChannel];
1607 
1608  Complex<Real> pi2(0.0, -2 * M_PI);
1609 
1610  for (uint i = 0; i < nChannel; i++) {
1611  waveRatio[i] = context_.wave_length[nChannel - 1] / context_.wave_length[i];
1612 
1613  Real ratioX = x * waveRatio[i];
1614  Real ratioY = y * waveRatio[i];
1615 
1616 #ifdef _OPENMP
1617 #pragma omp parallel for firstprivate(ppX, ppY, ppX2, ppY2, ssX, ssY, pi2)
1618 #endif
1619  for (int y = 0; y < pnY; y++) {
1620  Complex<Real> yy(0, 0);
1621  if (bAxisY) {
1622  Real startY = ssY + (ppY * y);
1623  Real shiftY = startY / ppY2 * ratioY;
1624  yy = (pi2 * shiftY).exp();
1625  }
1626  int offset = y * pnX;
1627 
1628  for (int x = 0; x < pnX; x++) {
1629  if (bAxisY) {
1630  complex_H[i][offset + x] = complex_H[i][offset + x] * yy;
1631  }
1632  if (bAxisX) {
1633  Real startX = ssX + (ppX * x);
1634  Real shiftX = startX / ppX2 * ratioX;
1635  Complex<Real> xx = (pi2 * shiftX).exp();
1636  complex_H[i][offset + x] = complex_H[i][offset + x] * xx;
1637  }
1638  }
1639  }
1640  }
1641  return true;
1642 }
1643 
1644 void ophGen::waveCarry(Real carryingAngleX, Real carryingAngleY, Real distance)
1645 {
1646  const int pnX = context_.pixel_number[_X];
1647  const int pnY = context_.pixel_number[_Y];
1648  const Real ppX = context_.pixel_pitch[_X];
1649  const Real ppY = context_.pixel_pitch[_Y];
1650  const long long int pnXY = pnX * pnY;
1651  const uint nChannel = context_.waveNum;
1652  Real dfx = 1 / ppX / pnX;
1653  Real dfy = 1 / ppY / pnY;
1654 
1655  int beginX = -pnX >> 1;
1656  int endX = pnX >> 1;
1657  int beginY = pnY >> 1;
1658  int endY = -pnY >> 1;
1659 
1660  Real carryX = distance * tan(carryingAngleX);
1661  Real carryY = distance * tan(carryingAngleY);
1662 
1663  for (uint ch = 0; ch < nChannel; ch++) {
1664  int k = 0;
1665  for (int i = beginY; i > endY; i--)
1666  {
1667  Real _carryY = carryY * i * dfy;
1668 
1669  for (int j = beginX; j < endX; j++)
1670  {
1671  Real x = j * dfx;
1672  Complex<Real> carrier(0, carryX * x + _carryY);
1673  complex_H[ch][k] = complex_H[ch][k] * exp(carrier);
1674 
1675  k++;
1676  }
1677  }
1678  }
1679 }
1680 
1681 void ophGen::waveCarry(Complex<Real>* src, Complex<Real>* dst, Real wavelength, int carryIdxX, int carryIdxY)
1682 {
1683  const int pnX = context_.pixel_number[_X];
1684  const int pnY = context_.pixel_number[_Y];
1685  const Real ppX = context_.pixel_pitch[_X];
1686  const Real ppY = context_.pixel_pitch[_Y];
1687  const long long int pnXY = pnX * pnY;
1688  Real dfx = 1 / ppX / pnX;
1689  Real dfy = 1 / ppY / pnY;
1690 
1691  int beginX = -pnX >> 1;
1692  int endX = pnX >> 1;
1693  int beginY = pnY >> 1;
1694  int endY = -pnY >> 1;
1695  Real pi2 = M_PI * 2;
1696 
1697  int k = 0;
1698  for (int i = beginY; i > endY; i--)
1699  {
1700  Real y = i * dfy * ppY * carryIdxY;
1701 
1702  for (int j = beginX; j < endX; j++)
1703  {
1704  Real x = j * dfx * ppX * carryIdxX;
1705 
1706  Complex<Real> carrier(0, pi2 * (x + y));
1707  dst[k] = src[k] * exp(carrier);
1708 
1709  k++;
1710  }
1711  }
1712 }
1713 
1714 void ophGen::encodeSideBand(bool bCPU, ivec2 sig_location)
1715 {
1716  if (complex_H == nullptr) {
1717  LOG("Not found diffracted data.");
1718  return;
1719  }
1720 
1721  const uint pnX = context_.pixel_number[_X];
1722  const uint pnY = context_.pixel_number[_Y];
1723 
1724  int cropx1, cropx2, cropx, cropy1, cropy2, cropy;
1725  if (sig_location[1] == 0) { //Left or right half
1726  cropy1 = 1;
1727  cropy2 = pnY;
1728  }
1729  else {
1730  cropy = (int)floor(((Real)pnY) / 2);
1731  cropy1 = cropy - (int)floor(((Real)cropy) / 2);
1732  cropy2 = cropy1 + cropy - 1;
1733  }
1734 
1735  if (sig_location[0] == 0) { // Upper or lower half
1736  cropx1 = 1;
1737  cropx2 = pnX;
1738  }
1739  else {
1740  cropx = (int)floor(((Real)pnX) / 2);
1741  cropx1 = cropx - (int)floor(((Real)cropx) / 2);
1742  cropx2 = cropx1 + cropx - 1;
1743  }
1744 
1745  cropx1 -= 1;
1746  cropx2 -= 1;
1747  cropy1 -= 1;
1748  cropy2 -= 1;
1749 
1750  if (bCPU)
1751  encodeSideBand_CPU(cropx1, cropx2, cropy1, cropy2, sig_location);
1752  else
1753  encodeSideBand_GPU(cropx1, cropx2, cropy1, cropy2, sig_location);
1754 }
1755 
1756 void ophGen::encodeSideBand_CPU(int cropx1, int cropx2, int cropy1, int cropy2, oph::ivec2 sig_location)
1757 {
1758  const uint pnX = context_.pixel_number[_X];
1759  const uint pnY = context_.pixel_number[_Y];
1760  const long long int pnXY = pnX * pnY;
1761  const uint nChannel = context_.waveNum;
1762 
1763  Complex<Real>* h_crop = new Complex<Real>[pnXY];
1764 
1765  for (uint ch = 0; ch < nChannel; ch++) {
1766 
1767  memset(h_crop, 0.0, sizeof(Complex<Real>) * pnXY);
1768 #ifdef _OPENMP
1769 #pragma omp parallel for firstprivate(cropx1, cropx2, cropy1, cropy2, pnX)
1770 #endif
1771  for (int p = 0; p < pnXY; p++)
1772  {
1773  int x = p % pnX;
1774  int y = p / pnX;
1775  if (x >= cropx1 && x <= cropx2 && y >= cropy1 && y <= cropy2)
1776  h_crop[p] = complex_H[ch][p];
1777  }
1778 
1779  Complex<Real> *in = nullptr;
1780 
1781  fft2(ivec2(pnX, pnY), in, OPH_BACKWARD);
1782  fft2(h_crop, h_crop, pnX, pnY, OPH_BACKWARD, true);
1783 
1784  memset(m_lpEncoded[ch], 0.0, sizeof(Real) * pnXY);
1785 #ifdef _OPENMP
1786 #pragma omp parallel for firstprivate(sig_location)
1787 #endif
1788  for (int i = 0; i < pnXY; i++) {
1789  Complex<Real> shift_phase(1, 0);
1790  getShiftPhaseValue(shift_phase, i, sig_location);
1791 
1792  m_lpEncoded[ch][i] = (h_crop[i] * shift_phase).real();
1793  }
1794  }
1795  delete[] h_crop;
1796 }
1797 
1798 void ophGen::encodeSideBand_GPU(int cropx1, int cropx2, int cropy1, int cropy2, oph::ivec2 sig_location)
1799 {
1800  const int pnX = context_.pixel_number[_X];
1801  const int pnY = context_.pixel_number[_Y];
1802  const long long int pnXY = pnX * pnY;
1803  const Real ppX = context_.pixel_pitch[_X];
1804  const Real ppY = context_.pixel_pitch[_Y];
1805  const Real ssX = context_.ss[_X] = pnX * ppX;
1806  const Real ssY = context_.ss[_Y] = pnY * ppY;
1807  const uint nChannel = context_.waveNum;
1808 
1809  cufftDoubleComplex *k_temp_d_, *u_complex_gpu_;
1810  cudaStream_t stream_;
1811  cudaStreamCreate(&stream_);
1812 
1813  cudaMalloc((void**)&u_complex_gpu_, sizeof(cufftDoubleComplex) * pnXY);
1814  cudaMalloc((void**)&k_temp_d_, sizeof(cufftDoubleComplex) * pnXY);
1815 
1816  for (uint ch = 0; ch < nChannel; ch++) {
1817  cudaMemcpy(u_complex_gpu_, complex_H[ch], sizeof(cufftDoubleComplex) * pnXY, cudaMemcpyHostToDevice);
1818 
1819  cudaMemsetAsync(k_temp_d_, 0, sizeof(cufftDoubleComplex) * pnXY, stream_);
1820  cudaCropFringe(stream_, pnX, pnY, u_complex_gpu_, k_temp_d_, cropx1, cropx2, cropy1, cropy2);
1821 
1822  cudaMemsetAsync(u_complex_gpu_, 0, sizeof(cufftDoubleComplex) * pnXY, stream_);
1823  cudaFFT(stream_, pnX, pnY, k_temp_d_, u_complex_gpu_, 1, true);
1824 
1825  cudaMemsetAsync(k_temp_d_, 0, sizeof(cufftDoubleComplex) * pnXY, stream_);
1826  cudaGetFringe(stream_, pnX, pnY, u_complex_gpu_, k_temp_d_, sig_location[0], sig_location[1], ssX, ssY, ppX, ppY, M_PI);
1827 
1828  cufftDoubleComplex* sample_fd = (cufftDoubleComplex*)malloc(sizeof(cufftDoubleComplex) * pnXY);
1829  memset(sample_fd, 0.0, sizeof(cufftDoubleComplex) * pnXY);
1830 
1831  cudaMemcpyAsync(sample_fd, k_temp_d_, sizeof(cufftDoubleComplex) * pnXY, cudaMemcpyDeviceToHost, stream_);
1832  memset(m_lpEncoded[ch], 0.0, sizeof(Real) * pnXY);
1833 
1834  for (int i = 0; i < pnX * pnY; i++)
1835  m_lpEncoded[ch][i] = sample_fd[i].x;
1836 
1837  cudaFree(sample_fd);
1838  }
1839  cudaStreamDestroy(stream_);
1840 }
1841 
1842 void ophGen::getShiftPhaseValue(oph::Complex<Real>& shift_phase_val, int idx, oph::ivec2 sig_location)
1843 {
1844  const int pnX = context_.pixel_number[_X];
1845  const int pnY = context_.pixel_number[_Y];
1846  const Real ppX = context_.pixel_pitch[_X];
1847  const Real ppY = context_.pixel_pitch[_Y];
1848  const Real ssX = context_.ss[_X] = pnX * ppX;
1849  const Real ssY = context_.ss[_Y] = pnY * ppY;
1850 
1851  if (sig_location[1] != 0)
1852  {
1853  int r = idx / pnX;
1854  int c = idx % pnX;
1855  Real yy = (ssY / 2.0) - (ppY)*r - ppY;
1856 
1857  oph::Complex<Real> val;
1858  if (sig_location[1] == 1)
1859  val[_IM] = 2 * M_PI * (yy / (4 * ppY));
1860  else
1861  val[_IM] = 2 * M_PI * (-yy / (4 * ppY));
1862 
1863  val.exp();
1864  shift_phase_val *= val;
1865  }
1866 
1867  if (sig_location[0] != 0)
1868  {
1869  int r = idx / pnX;
1870  int c = idx % pnX;
1871  Real xx = (-ssX / 2.0) - (ppX)*c - ppX;
1872 
1873  oph::Complex<Real> val;
1874  if (sig_location[0] == -1)
1875  val[_IM] = 2 * M_PI * (-xx / (4 * ppX));
1876  else
1877  val[_IM] = 2 * M_PI * (xx / (4 * ppX));
1878 
1879  val.exp();
1880  shift_phase_val *= val;
1881  }
1882 }
1883 
1884 void ophGen::GetRandomPhaseValue(Complex<Real>& rand_phase_val, bool rand_phase)
1885 {
1886  if (rand_phase)
1887  {
1888  rand_phase_val[_RE] = 0.0;
1889  Real min, max;
1890 #if REAL_IS_DOUBLE & true
1891  min = 0.0;
1892  max = 1.0;
1893 #else
1894  min = 0.f;
1895  max = 1.f;
1896 #endif
1897  rand_phase_val[_IM] = 2 * M_PI * oph::rand(min, max);
1898  rand_phase_val.exp();
1899 
1900  }
1901  else {
1902  rand_phase_val[_RE] = 1.0;
1903  rand_phase_val[_IM] = 0.0;
1904  }
1905 }
1906 
1907 void ophGen::setResolution(ivec2 resolution)
1908 {
1909  // ±âÁ¸ ÇØ»óµµ¿Í ´Ù¸£¸é ¹öÆÛ¸¦ ´Ù½Ã »ý¼º.
1910  if (context_.pixel_number != resolution) {
1911  setPixelNumber(resolution);
1912  Openholo::setPixelNumberOHC(resolution);
1913  initialize();
1914  }
1915 }
1916 
1917 template <typename T>
1918 void ophGen::RealPart(Complex<T> *holo, T *encoded, const int size)
1919 {
1920 #ifdef _OPENMP
1921 #pragma omp parallel for
1922 #endif
1923  for (int i = 0; i < size; i++) {
1924  encoded[i] = real(holo[i]);
1925  }
1926 }
1927 
1928 template <typename T>
1929 void ophGen::ImaginaryPart(Complex<T> *holo, T *encoded, const int size)
1930 {
1931 #ifdef _OPENMP
1932 #pragma omp parallel for
1933 #endif
1934  for (int i = 0; i < size; i++) {
1935  encoded[i] = imag(holo[i]);
1936  }
1937 }
1938 
1939 template <typename T>
1940 void ophGen::Phase(Complex<T> *holo, T *encoded, const int size)
1941 {
1942  double pi2 = M_PI * 2;
1943 #ifdef _OPENMP
1944 #pragma omp parallel for firstprivate(pi2)
1945 #endif
1946  for (int i = 0; i < size; i++) {
1947  encoded[i] = (holo[i].angle() + M_PI) / pi2; // 0 ~ 1
1948  }
1949 }
1950 
1951 template <typename T>
1952 void ophGen::Amplitude(Complex<T> *holo, T *encoded, const int size)
1953 {
1954 #ifdef _OPENMP
1955 #pragma omp parallel for
1956 #endif
1957  for (int i = 0; i < size; i++) {
1958  encoded[i] = holo[i].mag();
1959  }
1960 }
1961 
1962 template <typename T>
1963 void ophGen::TwoPhase(Complex<T>* holo, T* encoded, const int size)
1964 {
1965  int resize = size >> 1;
1966  Complex<T>* normCplx = new Complex<T>[resize];
1967 
1968 #ifdef _OPENMP
1969 #pragma omp parallel for
1970 #endif
1971  for (int i = 0; i < resize; i++) {
1972  normCplx[i] = holo[i * 2];
1973  }
1974 
1975  oph::normalize<T>(normCplx, normCplx, resize);
1976 
1977  T* ampl = new T[resize];
1978  Amplitude(normCplx, ampl, resize);
1979 
1980  T* phase = new T[resize];
1981  Phase(normCplx, phase, resize);
1982 
1983 #ifdef _OPENMP
1984 #pragma omp parallel for
1985 #endif
1986  for (int i = 0; i < resize; i++) {
1987  T delPhase = acos(ampl[i]);
1988  encoded[i * 2] = (phase[i] + M_PI) + delPhase;
1989  encoded[i * 2 + 1] = (phase[i] + M_PI) - delPhase;
1990  }
1991 
1992  delete[] normCplx;
1993  delete[] ampl;
1994  delete[] phase;
1995 }
1996 
1997 template <typename T>
1998 void ophGen::Burckhardt(Complex<T>* holo, T* encoded, const int size)
1999 {
2000  int resize = size / 3;
2001  Complex<T>* norm = new Complex<T>[resize];
2002 #ifdef _OPENMP
2003 #pragma omp parallel for
2004 #endif
2005  for (int i = 0; i < resize; i++) {
2006  norm[i] = holo[i * 3];
2007  }
2008 
2009  oph::normalize(norm, norm, resize);
2010 
2011  T* phase = new T[resize];
2012  Phase(norm, phase, resize);
2013 
2014  T* ampl = new T[resize];
2015  Amplitude(norm, ampl, resize);
2016 
2017  T sqrt3 = sqrt(3);
2018  T pi2 = 2 * M_PI;
2019  T pi4 = 4 * M_PI;
2020 
2021 #ifdef _OPENMP
2022 #pragma omp parallel for firstprivate(pi2, pi4, sqrt3)
2023 #endif
2024  for (int i = 0; i < resize; i++) {
2025  int idx = 3 * i;
2026  if (phase[i] >= 0 && phase[i] < (pi2 / 3))
2027  {
2028  encoded[idx] = ampl[i] * (cos(phase[i]) + sin(phase[i]) / sqrt3);
2029  encoded[idx + 1] = 2 * sin(phase[i]) / sqrt3;
2030  }
2031  else if (phase[i] >= (pi2 / 3) && phase[i] < (pi4 / 3))
2032  {
2033  encoded[idx + 1] = ampl[i] * (cos(phase[i] - (pi2 / 3)) + sin(phase[i] - (pi2 / 3)) / sqrt3);
2034  encoded[idx + 2] = 2 * sin(phase[i] - (pi2 / 3)) / sqrt3;
2035  }
2036  else if (phase[i] >= (pi4 / 3) && phase[i] < (pi2))
2037  {
2038  encoded[idx + 2] = ampl[i] * (cos(phase[i] - (pi4 / 3)) + sin(phase[i] - (pi4 / 3)) / sqrt3);
2039  encoded[idx] = 2 * sin(phase[i] - (pi4 / 3)) / sqrt3;
2040  }
2041  }
2042 
2043  delete[] ampl;
2044  delete[] phase;
2045  delete[] norm;
2046 }
2047 
2048 
2049 template <typename T>
2050 void ophGen::SimpleNI(Complex<T>* holo, T* encoded, const int size)
2051 {
2052  T* tmp1 = new T[size];
2053 #ifdef _OPENMP
2054 #pragma omp parallel for
2055 #endif
2056  for (int i = 0; i < size; i++) {
2057  tmp1[i] = holo[i].mag();
2058  }
2059 
2060  T max = maxOfArr(tmp1, size);
2061  delete[] tmp1;
2062 
2063 #ifdef _OPENMP
2064 #pragma omp parallel for firstprivate(max)
2065 #endif
2066  for (int i = 0; i < size; i++) {
2067  T tmp = (holo[i] + max).mag();
2068  encoded[i] = tmp * tmp;
2069  }
2070 }
2071 
2072 void ophGen::transVW(int nVertex, Vertex* dst, Vertex* src)
2073 {
2074  Real fieldLens = m_dFieldLength;
2075  for (int i = 0; i < nVertex; i++)
2076  {
2077  dst[i].point.pos[_X] = -fieldLens * src[i].point.pos[_X] / (src[i].point.pos[_X] - fieldLens);
2078  dst[i].point.pos[_Y] = -fieldLens * src[i].point.pos[_Y] / (src[i].point.pos[_Y] - fieldLens);
2079  dst[i].point.pos[_Z] = -fieldLens * src[i].point.pos[_Z] / (src[i].point.pos[_Z] - fieldLens);
2080  }
2081 }
2082 
2083 void ophGen::transVW(int nVertex, Real *dst, Real *src)
2084 {
2085  Real fieldLens = m_dFieldLength;
2086  for (int i = 0; i < nVertex; i++) {
2087  dst[i] = -fieldLens * src[i] / (src[i] - fieldLens);
2088  }
2089 }
2090 
2091 void ophGen::ScaleChange(Real *src, Real *dst, int nSize, Real scaleX, Real scaleY, Real scaleZ)
2092 {
2093  Real x = scaleX;
2094  Real y = scaleY;
2095  Real z = scaleZ;
2096 
2097 #ifdef _OPENMP
2098 #pragma omp parallel for firstprivate(x, y, z)
2099 #endif
2100  for (int i = 0; i < nSize; i++) {
2101  dst[i + 0] = src[i + 0] * x;
2102  dst[i + 1] = src[i + 1] * y;
2103  dst[i + 2] = src[i + 2] * z;
2104  }
2105 }
2106 
2107 void ophGen::GetMaxMin(Real *src, int len, Real& max, Real& min)
2108 {
2109  Real maxTmp = MIN_DOUBLE;
2110  Real minTmp = MAX_DOUBLE;
2111 
2112  for (int i = 0; i < len; i++) {
2113  if (src[i] > maxTmp)
2114  maxTmp = src[i];
2115  if (src[i] < minTmp)
2116  minTmp = src[i];
2117  }
2118  max = maxTmp;
2119  min = minTmp;
2120 }
2121 
2122 
2123 void ophGen::CorrectionChromaticAberration(uchar* src, uchar* dst, int width, int height, int ch)
2124 {
2125  if (ch < 2)
2126  {
2127  const uint nWave = context_.waveNum;
2128  const int pnXY = width * height;
2129  const Real lambda = context_.wave_length[ch];
2130  const Real waveRatio = nWave > 1 ? context_.wave_length[nWave - 1] / lambda : 1.0;
2131 
2132  int scaleX = round(width * 4 * waveRatio);
2133  int scaleY = round(height * 4 * waveRatio);
2134 
2135  int ww = width * 4;
2136  int hh = height * 4;
2137  int nSize = ww * hh;
2138  int nScaleSize = scaleX * scaleY;
2139  uchar *img_tmp = new uchar[nSize];
2140  uchar *imgScaled = new uchar[nScaleSize];
2141  imgScaleBilinear(src, imgScaled, width, height, scaleX, scaleY);
2142 
2143  memset(img_tmp, 0, sizeof(uchar) * nSize);
2144 
2145  int h1 = round((hh - scaleY) / 2);
2146  int w1 = round((ww - scaleX) / 2);
2147 #ifdef _OPENMP
2148 #pragma omp parallel for firstprivate(w1, h1, scaleX, ww)
2149 #endif
2150  for (int y = 0; y < scaleY; y++) {
2151  for (int x = 0; x < scaleX; x++) {
2152  img_tmp[(y + h1)*ww + x + w1] = imgScaled[y*scaleX + x];
2153  }
2154  }
2155  imgScaleBilinear(img_tmp, dst, ww, hh, width, height);
2156 
2157  delete[] img_tmp;
2158  delete[] imgScaled;
2159  }
2160 }
2161 
2162 
2163 void ophGen::ophFree(void)
2164 {
2166 
2167  const uint nChannel = context_.waveNum;
2168 
2169  if (m_lpEncoded != nullptr) {
2170  for (uint i = 0; i < nChannel; i++) {
2171  if (m_lpEncoded[i] != nullptr) {
2172  delete[] m_lpEncoded[i];
2173  m_lpEncoded[i] = nullptr;
2174  }
2175  }
2176  delete[] m_lpEncoded;
2177  m_lpEncoded = nullptr;
2178  }
2179 
2180  if (m_lpNormalized != nullptr) {
2181  for (uint i = 0; i < nChannel; i++) {
2182  if (m_lpNormalized[i] != nullptr) {
2183  delete[] m_lpNormalized[i];
2184  m_lpNormalized[i] = nullptr;
2185  }
2186  }
2187  delete[] m_lpNormalized;
2188  m_lpNormalized = nullptr;
2189  }
2190 
2191  if (freqW != nullptr) delete[] freqW;
2192 }
bool loadPLY(const std::string &fileName, ulonglong &n_points, Vertex **vertices)
Definition: PLYparser.cpp:142
ENCODE_FLAG
Definition: ophGen.h:84
void abs(const oph::Complex< T > &src, oph::Complex< T > &dst)
Definition: function.h:113
void Phase(Complex< T > *holo, T *encoded, const int size)
void RS_Diffraction(Point src, Complex< Real > *dst, Real lambda, Real distance, Real amplitude)
RS-diffraction method.
Abstract class.
Definition: Openholo.h:145
#define M_PI
Definition: define.h:52
Point point
Definition: struct.h:103
#define OPH_FORWARD
Definition: define.h:66
virtual bool saveAsImg(const char *fname, uint8_t bitsperpixel, uchar *src, int width, int height)
Function for creating image files.
Definition: Openholo.cpp:135
Real k
Definition: Openholo.h:103
ResolutionConfig resCfg
Definition: Openholo.h:486
Complex< Real > * freqW
Definition: ophGen.h:431
void setResolution(ivec2 resolution)
Function for setting buffer size.
void cudaFFT(CUstream_st *stream, int nx, int ny, cufftDoubleComplex *in_filed, cufftDoubleComplex *output_field, int direction, bool bNormailized=false)
Convert data from the spatial domain to the frequency domain using 2D FFT on GPU. ...
void conv_fft2(Complex< Real > *src1, Complex< Real > *src2, Complex< Real > *dst, ivec2 size)
Convolution between Complex arrays which have same size.
void GetMaxMin(Real *src, int len, Real &max, Real &min)
Definition: ophGen.cpp:2107
Complex< Real > * normalized
Definition: ophGen.h:427
vec3 shift
Definition: Openholo.h:102
bool getWeightED(const ivec2 holosize, const int type, ivec2 *pNw)
void AngularSpectrumMethod(Complex< Real > *src, Complex< Real > *dst, Real lambda, Real distance)
Angular spectrum propagation method.
void encodeSideBand_CPU(int cropx1, int cropx2, int cropy1, int cropy2, ivec2 sig_location)
Encode the CGH according to a signal location parameter on the CPU.
Definition: ophGen.cpp:1756
void binarization(Complex< Real > *src, Real *dst, const int size, int ENCODE_FLAG, Real threshold)
#define MIN_DOUBLE
Definition: define.h:148
void GetRandomPhaseValue(Complex< Real > &rand_phase_val, bool rand_phase)
Assign random phase value if random_phase == 1.
ulonglong n_points
Number of points.
Definition: ophGen.h:584
void setPixelNumberOHC(const ivec2 pixel_number)
getter/setter for OHC file read and write
Definition: Openholo.h:501
unsigned char uchar
Definition: typedef.h:64
HORIZONTAL
Definition: ImgControl.h:25
Real * maskHP
Definition: ophGen.h:435
Real m_dFieldLength
Definition: ophGen.h:348
long long int double_size
Definition: Openholo.h:69
NONE
Definition: ImgControl.h:25
void initialize(void)
Initialize variables for Hologram complex field, encoded data, normalized data.
Definition: ophGen.cpp:145
long long int size
Definition: Openholo.h:67
Real norm(const vec2 &a)
Definition: vec.h:417
Definition: struct.h:86
void RealPart(Complex< T > *holo, T *encoded, const int size)
Encoding method.
void setPixelPitchOHC(const vec2 pixel_pitch)
Definition: Openholo.h:504
#define MAX_DOUBLE
Definition: define.h:140
bool shiftW(ivec2 holosize)
vec2 ss
Definition: Openholo.h:104
return true
Definition: Openholo.cpp:434
void TwoPhase(Complex< T > *holo, T *encoded, const int size)
float Real
Definition: typedef.h:55
bool checkExtension(const char *fname, const char *ext)
Functions for extension checking.
Definition: Openholo.cpp:86
structure for 2-dimensional integer vector and its arithmetic.
Definition: ivec.h:66
#define _IM
Definition: complex.h:58
def k(wvl)
Definition: Depthmap.py:16
const XMLNode * FirstChild() const
Get the first child node, or null if none exists.
Definition: tinyxml2.h:761
VERTICAL
Definition: ImgControl.h:25
ivec2 offset
Definition: Openholo.h:100
void normalize(const Complex< T > *src, Complex< T > *dst, const int &size)
Normalize all elements of Complex<T>* src from 0 to 1.
Definition: function.h:176
bool rotate
Definition: Openholo.h:124
void Burckhardt(Complex< T > *holo, T *encoded, const int size)
ImageConfig imgCfg
Definition: Openholo.h:487
void setWavelengthOHC(const Real wavelength, const LenUnit wavelength_unit)
Definition: Openholo.h:507
vec2 pixel_pitch
Definition: Openholo.h:101
void ImaginaryPart(Complex< T > *holo, T *encoded, const int size)
void encodeSideBand(bool bCPU, ivec2 sig_location)
Encode the CGH according to a signal location parameter.
bool Shift(Real x, Real y)
Definition: ophGen.cpp:1589
void normalize(void)
Normalization function to save as image file after hologram creation.
Definition: ophGen.cpp:677
bool binaryErrorDiffusion(Complex< Real > *holo, Real *encoded, const ivec2 holosize, const int type, Real threshold)
bool save(const char *fname, uint8_t bitsperpixel=8, uchar *src=nullptr, uint px=0, uint py=0)
Function for saving image files.
ivec2 m_vecEncodeSize
Encoded hologram size, varied from encoding type.
Definition: ophGen.h:330
#define _X
Definition: define.h:92
void singleSideBand(Complex< Real > *holo, Real *encoded, const ivec2 holosize, int passband)
Encoding method.
void cudaGetFringe(CUstream_st *stream, int pnx, int pny, cufftDoubleComplex *in_field, cufftDoubleComplex *out_field, int sig_locationx, int sig_locationy, Real ssx, Real ssy, Real ppx, Real ppy, Real PI)
Encode the CGH according to a signal location parameter on the GPU.
void encodeSideBand_GPU(int cropx1, int cropx2, int cropy1, int cropy2, ivec2 sig_location)
Encode the CGH according to a signal location parameter on the GPU.
Definition: ophGen.cpp:1798
void getShiftPhaseValue(Complex< Real > &shift_phase_val, int idx, ivec2 sig_location)
Calculate the shift phase value.
Definition: ophGen.cpp:1842
void imgScaleBilinear(uchar *src, uchar *dst, int w, int h, int neww, int newh, int channels=1)
Function for change image size.
Definition: Openholo.cpp:437
Real * weight
Definition: ophGen.h:429
void * load(const char *fname)
Function for loading image files.
Definition: ophGen.cpp:794
Real ** m_lpEncoded
buffer to encoded.
Definition: ophGen.h:338
void SimpleNI(Complex< T > *holo, T *encoded, const int size)
bool m_bRandomPhase
Definition: ophGen.h:436
void CorrectionChromaticAberration(uchar *src, uchar *dst, int width, int height, int ch)
#define _RE
Definition: complex.h:55
void waveCarry(Real carryingAngleX, Real carryingAngleY, Real distance)
Wave carry.
Definition: ophGen.cpp:1644
void fft2(ivec2 n, Complex< Real > *in, int sign=OPH_FORWARD, uint flag=OPH_ESTIMATE)
Functions for performing fftw 2-dimension operations inside Openholo.
Definition: Openholo.cpp:559
void fresnelPropagation(OphConfig context, Complex< Real > *in, Complex< Real > *out, Real distance)
Fresnel propagation.
Definition: struct.h:102
uint waveNum
Definition: Openholo.h:105
void Fresnel_FFT(Complex< Real > *src, Complex< Real > *dst, Real lambda, Real waveRatio, Real distance)
Fresnel-fft method.
void Amplitude(Complex< T > *holo, T *encoded, const int size)
Complex< Real > * weightC
Definition: ophGen.h:430
SSB_PASSBAND
Passband in Single-side band encoding.
Definition: ophGen.h:295
void cudaCropFringe(CUstream_st *stream, int nx, int ny, cufftDoubleComplex *in_field, cufftDoubleComplex *out_field, int cropx1, int cropx2, int cropy1, int cropy2)
Crop input data according to x, y coordinates on GPU.
ivec2 pixel_number
Definition: Openholo.h:99
bool saveRefImages(char *fnameW, char *fnameWC, char *fnameAS, char *fnameSSB, char *fnameHP, char *fnameFreq, char *fnameReal, char *fnameBin, char *fnameReconBin, char *fnameReconErr, char *fnameReconNo)
Definition: ophGen.cpp:1131
virtual bool loadAsOhc(const char *fname)
Function to read OHC file.
Definition: Openholo.cpp:280
bool readConfig(const char *fname)
load to configuration file.
Definition: ophGen.cpp:220
#define _Y
Definition: define.h:96
XMLError LoadFile(const char *filename)
Definition: tinyxml2.cpp:2150
Real * maskSSB
Definition: ophGen.h:434
bool bUseDP
Definition: Openholo.h:98
virtual uchar * loadAsImg(const char *fname)
Function for loading image files.
Definition: Openholo.cpp:321
Real rand(const Real min, const Real max, oph::ulong _SEED_VALUE=0)
Get random Real value from min to max.
Definition: function.h:294
void getAmplitude(oph::Complex< Real > *src, Real *dst, const int &size)
Definition: function.h:330
void ScaleChange(Real *src, Real *dst, int nSize, Real scaleX, Real scaleY, Real scaleZ)
Definition: ophGen.cpp:2091
Vertex * vertices
Data of point clouds.
Definition: ophGen.h:588
#define PI
Definition: ophACPAS.h:9
Complex< Real > ** complex_H
Definition: Openholo.h:489
Real minOfArr(const std::vector< Real > &arr)
Definition: function.h:73
Data for Point Cloud.
Definition: ophGen.h:582
#define OPH_ESTIMATE
Definition: define.h:76
Complex< Real > * AS
Binary Encoding - Error diffusion.
Definition: ophGen.h:426
Real maxOfArr(const std::vector< Real > &arr)
Definition: function.h:87
uchar ** m_lpNormalized
buffer to normalized.
Definition: ophGen.h:340
void encoding()
Definition: ophGen.cpp:985
Real * binary
Definition: ophGen.h:433
void resetBuffer()
reset buffer
Definition: ophGen.cpp:824
void setPixelNumber(ivec2 n)
Function for setting the output resolution.
Definition: Openholo.h:242
#define _Z
Definition: define.h:100
#define ELAPSED_TIME(x, y)
Definition: function.h:59
int m_nStream
Definition: ophGen.h:349
ImgEncoderOhc * OHC_encoder
OHC file format Variables for read and write.
Definition: Openholo.h:494
ophGen(void)
Constructor.
Definition: ophGen.cpp:112
OphConfig context_
Definition: Openholo.h:485
void Fresnel_Diffraction(Point src, Complex< Real > *dst, Real lambda, Real distance, Real amplitude)
Fresnel-diffraction method.
#define OPH_BACKWARD
Definition: define.h:67
virtual bool loadAsOhc(const char *fname)
Function to read OHC file.
Definition: ophGen.cpp:806
void transVW(int nVertex, Vertex *dst, Vertex *src)
Definition: ophGen.cpp:2072
int ENCODE_METHOD
Encoding method flag.
Definition: ophGen.h:332
Real * wave_length
Definition: Openholo.h:106
void freqShift(Complex< Real > *holo, Complex< Real > *encoded, const ivec2 holosize, int shift_x, int shift_y)
Frequency shift.
int loadPointCloud(const char *pc_file, OphPointCloudData *pc_data_)
load to point cloud data.
Definition: ophGen.cpp:206
const XMLElement * FirstChildElement(const char *name=0) const
Definition: tinyxml2.cpp:940
bool mergeColor(int idx, int width, int height, uchar *src, uchar *dst)
Function for generate RGB image from each grayscale image.
Definition: Openholo.cpp:103
unsigned int uint
Definition: typedef.h:62
#define CUR_TIME
Definition: function.h:58
virtual void ophFree(void)
Pure virtual function for override in child classes.
Definition: Openholo.cpp:805
bool merge
Definition: Openholo.h:125
Definition: ophGen.h:76
Real pos[3]
Definition: struct.h:87
virtual void ophFree(void)
Pure virtual function for override in child classes.
Definition: ophGen.cpp:2163
Real * realEnc
Definition: ophGen.h:432
virtual ~ophGen(void)=0
Destructor.
Definition: ophGen.cpp:140