49 #include <cuda_runtime_api.h> 73 void cudaFFT(CUstream_st* stream,
int nx,
int ny, cufftDoubleComplex* in_filed, cufftDoubleComplex* output_field,
int direction,
bool bNormailized =
false);
89 void cudaCropFringe(CUstream_st* stream,
int nx,
int ny, cufftDoubleComplex* in_field, cufftDoubleComplex* out_field,
int cropx1,
int cropx2,
int cropy1,
int cropy2);
108 void cudaGetFringe(CUstream_st* stream,
int pnx,
int pny, cufftDoubleComplex* in_field, cufftDoubleComplex* out_field,
int sig_locationx,
int sig_locationy,
118 , m_lpEncoded(nullptr)
119 , m_lpNormalized(nullptr)
122 , m_dFieldLength(0.0)
126 , normalized(nullptr)
135 , m_bRandomPhase(false)
153 for (
uint i = 0; i < m_nOldChannel; i++) {
164 for (
uint i = 0; i < nChannel; i++) {
170 for (
uint i = 0; i < m_nOldChannel; i++) {
180 for (
uint i = 0; i < nChannel; i++) {
186 for (
uint i = 0; i < m_nOldChannel; i++) {
196 for (
uint i = 0; i < nChannel; i++) {
201 m_nOldChannel = nChannel;
229 LOG(
"<FAILED> Wrong file ext.\n");
235 LOG(
"<FAILED> Loading file (%d)\n", ret);
241 char szNodeName[32] = { 0, };
243 sprintf(szNodeName,
"SLM_WaveNum");
245 if (!next ||
XML_SUCCESS != next->QueryIntText(&nWave))
247 LOG(
"<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
254 for (
int i = 1; i <= nWave; i++) {
255 sprintf(szNodeName,
"SLM_WaveLength_%d", i);
259 LOG(
"<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
264 sprintf(szNodeName,
"SLM_PixelNumX");
268 LOG(
"<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
272 sprintf(szNodeName,
"SLM_PixelNumY");
276 LOG(
"<FAILED> Not found node : \'%s\' (Integer) \n", szNodeName);
280 sprintf(szNodeName,
"SLM_PixelPitchX");
284 LOG(
"<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
288 sprintf(szNodeName,
"SLM_PixelPitchY");
292 LOG(
"<FAILED> Not found node : \'%s\' (Double) \n", szNodeName);
335 for (
int i = 0; i < nWave; i++)
338 LOG(
"**************************************************\n");
339 LOG(
" Read Config (Common) \n");
345 LOG(
"4) Image Rotate : %s\n",
imgCfg.
rotate ?
"Y" :
"N");
349 LOG(
"6) Image Merge : %s\n",
imgCfg.
merge ?
"Y" :
"N");
350 LOG(
"**************************************************\n");
361 const int pnXY = pnX * pnY;
364 const Real ssX = pConfig->
ss[
_X] = pnX * ppX;
365 const Real ssY = pConfig->
ss[
_Y] = pnY * ppY;
366 const int offsetX = pConfig->
offset[
_X];
367 const int offsetY = pConfig->
offset[
_Y];
369 const Real tx = lambda / (2 * ppX);
370 const Real ty = lambda / (2 * ppY);
371 const Real sqrtX = sqrt(1 - (tx * tx));
372 const Real sqrtY = sqrt(1 - (ty * ty));
373 const Real x = -ssX / 2;
374 const Real y = -ssY / 2;
378 Real ampZ = amplitude * z;
391 floor((_xbound[
_X] - x) / ppX) + 1,
392 floor((_xbound[
_Y] - x) / ppX) + 1
396 pnY - floor((_ybound[
_Y] - y) / ppY),
397 pnY - floor((_ybound[
_X] - y) / ppY)
400 if (Xbound[
_X] > pnX) Xbound[
_X] = pnX;
401 if (Xbound[
_Y] < 0) Xbound[
_Y] = 0;
402 if (Ybound[
_X] > pnY) Ybound[
_X] = pnY;
403 if (Ybound[
_Y] < 0) Ybound[
_Y] = 0;
406 for (
int yytr = Ybound[
_Y]; yytr < Ybound[
_X]; ++yytr)
408 int offset = yytr * pnX;
409 Real yyy = y + ((pnY - yytr + offsetY) * ppY);
416 for (
int xxtr = Xbound[
_Y]; xxtr < Xbound[
_X]; ++xxtr)
418 Real xxx = x + ((xxtr - 1 + offsetX) * ppX);
425 if (((xxx < range_x[
_X]) && (xxx > range_x[
_Y])) && ((yyy < range_y[
_X]) && (yyy > range_y[
_Y]))) {
427 Real operand = lambda * r * r;
428 Real res_real = (ampZ * sin(kr)) / operand;
429 Real res_imag = (-ampZ * cos(kr)) / operand;
432 dst[offset + xxtr][
_RE] += res_real;
436 dst[offset + xxtr][
_IM] += res_imag;
448 const int pnXY = pnX * pnY;
451 const Real ssX = pConfig->
ss[
_X] = pnX * ppX;
452 const Real ssY = pConfig->
ss[
_Y] = pnY * ppY;
453 const int offsetX = pConfig->
offset[
_X];
454 const int offsetY = pConfig->
offset[
_Y];
462 Real operand = lambda * z;
465 src.
pos[
_X] +
abs(operand / (2 * ppX)),
466 src.
pos[
_X] -
abs(operand / (2 * ppX))
470 src.
pos[
_Y] +
abs(operand / (2 * ppY)),
471 src.
pos[
_Y] -
abs(operand / (2 * ppY))
475 floor((_xbound[
_X] - x) / ppX) + 1,
476 floor((_xbound[
_Y] - x) / ppX) + 1
480 pnY - floor((_ybound[
_Y] - y) / ppY),
481 pnY - floor((_ybound[
_X] - y) / ppY)
484 if (Xbound[
_X] > pnX) Xbound[
_X] = pnX;
485 if (Xbound[
_Y] < 0) Xbound[
_Y] = 0;
486 if (Ybound[
_X] > pnY) Ybound[
_X] = pnY;
487 if (Ybound[
_Y] < 0) Ybound[
_Y] = 0;
489 for (
int yytr = Ybound[
_Y]; yytr < Ybound[
_X]; ++yytr)
491 Real yyy = (y + (pnY - yytr + offsetY) * ppY) - src.
pos[
_Y];
492 int offset = yytr * pnX;
493 for (
int xxtr = Xbound[
_Y]; xxtr < Xbound[
_X]; ++xxtr)
495 Real xxx = (x + (xxtr - 1 + offsetX) * ppX) - src.
pos[
_X];
496 Real p = k * (xxx * xxx + yyy * yyy + 2 * zz) / (2 * z);
498 Real res_real = amplitude * sin(p) / operand;
499 Real res_imag = amplitude * (-cos(p)) / operand;
504 dst[offset + xxtr][
_RE] += res_real;
508 dst[offset + xxtr][
_IM] += res_imag;
518 const int pnXY = pnX * pnY;
521 const Real ssX = pConfig->
ss[
_X] = pnX * ppX;
522 const Real ssY = pConfig->
ss[
_Y] = pnY * ppY;
523 const Real ssX2 = ssX * 2;
524 const Real ssY2 = ssY * 2;
527 int newSize = pnXY * 4;
532 int half_pnX = pnX >> 1;
533 int half_pnY = pnY >> 1;
537 for (
int i = 0; i < pnY; i++) {
538 for (
int j = 0; j < pnX; j++) {
539 in2x[(i + half_pnY) * pnX2 + (j + half_pnX)] = src[idxIn++];
553 for (
int idxFy = -pnY; idxFy < pnY; idxFy++) {
554 for (
int idxFx = -pnX; idxFx < pnX; idxFx++) {
555 fx[i] = idxFx / ssX2;
556 fy[i] = idxFy / ssY2;
567 Real lambda_square = lambda * lambda;
570 for (
int i = 0; i < newSize; i++) {
571 sqrtPart = sqrt(1 / lambda_square - fx[i] * fx[i] - fy[i] * fy[i]);
573 prop[i][
_IM] *= sqrtPart;
574 temp2[i] = temp1[i] * exp(prop[i]);
582 for (
int i = 0; i < pnY; i++) {
583 for (
int j = 0; j < pnX; j++) {
584 dst[idxOut++] = temp3[(i + half_pnY) * pnX2 + (j + half_pnX)];
601 const int N = pnX * pnY;
612 Real kd = k * distance;
613 Real fx = -1 / (ppX * 2);
614 Real fy = 1 / (ppY * 2);
617 #pragma omp parallel for firstprivate(pnX, dfx, dfy, lambda, kd, kk) 619 for (
int i = 0; i < N; i++)
624 Real fxx = fx + dfx * x;
625 Real fyy = fy - dfy - dfy * y;
627 Real fxxx = lambda * fxx;
628 Real fyyy = lambda * fyy;
630 Real sval = sqrt(1 - (fxxx * fxxx) - (fyyy * fyyy));
635 bool prop_mask = ((fxx * fxx + fyy * fyy) < kk) ?
true :
false;
639 u_frequency = kernel * src[i];
640 dst[i][
_RE] += u_frequency[
_RE];
641 dst[i][
_IM] += u_frequency[
_IM];
648 int N = size[
_X] * size[
_Y];
658 for (
int i = 0; i < N; i++)
659 dstFT[i] = src1FT[i] * src2FT[i];
680 const long long int N = pnX * pnY;
683 for (
uint ch = 0; ch < nWave; ch++)
695 for (
uint ch = 0; ch < nWave; ch++)
698 #pragma omp parallel for firstprivate(min, gap, pnX) 700 for (
int j = 0; j < pnY; j++) {
701 for (
int i = 0; i < pnX; i++) {
702 int idx = j * pnX + i;
713 if (fname ==
nullptr)
return bOK;
720 if (px == 0 && py == 0)
724 std::string file = fname;
725 std::replace(file.begin(), file.end(),
'\\',
'/');
728 std::vector<std::string> components;
729 std::stringstream ss(file);
733 while (std::getline(ss, item, token)) {
734 components.push_back(item);
739 for (
size_t i = 0; i < components.size() - 1; i++)
741 dir += components[i];
745 std::string filename = components[components.size() - 1];
749 size_t ext_pos = file.rfind(
".");
750 hasExt = (ext_pos == string::npos) ?
false :
true;
753 filename.append(
".bmp");
755 std::string fullpath = dir + filename;
757 if (src ==
nullptr) {
760 saveAsImg(fullpath.c_str(), bitsperpixel, source, p[
_X], p[
_Y]);
762 else if (nChannel == 3) {
764 uint nSize = (((p[
_X] * bitsperpixel >> 3) + 3) & ~3) * p[
_Y];
765 source =
new uchar[nSize];
767 for (
uint i = 0; i < nChannel; i++) {
771 saveAsImg(fullpath.c_str(), bitsperpixel, source, p[
_X], p[
_Y]);
772 if (bAlloc)
delete[] source;
775 for (
uint i = 0; i < nChannel; i++) {
776 char path[FILENAME_MAX] = { 0, };
777 sprintf(path,
"%s%d_%s", dir.c_str(), i, filename.c_str());
779 saveAsImg(path, bitsperpixel / nChannel, source, p[
_X], p[
_Y]);
786 saveAsImg(fullpath.c_str(), bitsperpixel, source, p[
_X], p[
_Y]);
812 for (
uint ch = 0; ch < nChannel; ch++) {
844 Real *encoded =
nullptr;
858 LOG(
"<FAILED> WRONG PARAMETERS.\n");
891 default: LOG(
"<FAILED> WRONG PARAMETERS.\n");
return;
894 if (holo ==
nullptr) holo =
complex_H[0];
904 holo ==
nullptr ? holo = *
complex_H : holo;
905 encoded ==
nullptr ? encoded = *
m_lpEncoded : encoded;
911 for (
uint ch = 0; ch < nChannel; ch++) {
930 LOG(
"ENCODE_OFFSSB");
939 LOG(
"<FAILED> WRONG PARAMETERS.\n");
951 switch (BIN_ENCODE_FLAG) {
953 LOG(
"ENCODE_SIMPLEBINARY\n");
954 if (holo ==
nullptr || encoded ==
nullptr)
955 for (
uint ch = 0; ch < nChannel; ch++) {
962 LOG(
"ENCODE_EDBINARY\n");
964 LOG(
"<FAILED> WRONG PARAMETERS : %d\n",
ENCODE_FLAG);
967 if (holo ==
nullptr || encoded ==
nullptr)
968 for (
uint ch = 0; ch < nChannel; ch++) {
975 LOG(
"<FAILED> WRONG PARAMETERS.\n");
990 for (
uint ch = 0; ch < nChannel; ch++) {
1004 LOG(
"ENCODE_SIMPLENI\n");
1008 LOG(
"ENCODE_REAL\n");
1012 LOG(
"ENCODE_BURCKHARDT\n");
1016 LOG(
"ENCODE_TWOPHASE\n");
1020 LOG(
"ENCODE_PHASE\n");
1024 LOG(
"ENCODE_AMPLITUDE\n");
1028 LOG(
"ENCODE_SSB\n");
1032 LOG(
"ENCODE_OFFSSB\n");
1037 LOG(
"<FAILED> WRONG PARAMETERS.\n");
1045 const int nX = holosize[
_X];
1046 const int nY = holosize[
_Y];
1047 const int half_nX = nX >> 1;
1048 const int half_nY = nY >> 1;
1050 const long long int N = nX * nY;
1051 const long long int half_N = N >> 1;
1062 for (
int i = 0; i < nY; i++)
1065 for (
int j = half_nX; j < nX; j++)
1072 for (
int i = 0; i < nY; i++)
1075 for (
int j = 0; j < half_nX; j++)
1098 oph::realPart<Real>(filtered, realFiltered, N);
1104 delete[] realFiltered;
1109 int N = holosize[
_X] * holosize[
_Y];
1117 circShift<Complex<Real>>(
AS, shifted, shift_x, shift_y, holosize.v[
_X], holosize.v[
_Y]);
1128 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)
1131 int nx = holosize[
_X], ny = holosize[
_Y];
1139 cout <<
"W saved" << endl;
1141 oph::absCplxArr<Real>(
weightC, temp1, ss);
1144 cout <<
"WC saved" << endl;
1146 oph::absCplxArr<Real>(
AS, temp1, ss);
1149 cout <<
"AS saved" << endl;
1153 cout <<
"SSB saved" << endl;
1157 cout <<
"HP saved" << endl;
1159 oph::absCplxArr<Real>(
freqW, temp1, ss);
1162 cout <<
"Freq saved" << endl;
1166 cout <<
"Real saved" << endl;
1170 cout <<
"Bin saved" << endl;
1174 for (
int i = 0; i < ss; i++) {
1180 for (
int i = 0; i < ss; i++) {
1191 oph::absCplxArr<Real>(reconBin, temp1, ss);
1192 for (
int i = 0; i < ss; i++) {
1193 temp1[i] = temp1[i] * temp1[i];
1196 saveAsImg(fnameReconBin, 8, temp2, nx, ny);
1197 cout <<
"recon bin saved" << endl;
1201 for (
int i = 0; i < ss; i++) {
1207 for (
int i = 0; i < ss; i++) {
1217 oph::absCplxArr<Real>(reconBin, temp1, ss);
1218 for (
int i = 0; i < ss; i++) {
1219 temp1[i] = temp1[i] * temp1[i];
1222 saveAsImg(fnameReconErr, 8, temp2, nx, ny);
1223 cout <<
"recon error saved" << endl;
1229 for (
int i = 0; i < ss; i++) {
1235 for (
int i = 0; i < ss; i++) {
1245 oph::absCplxArr<Real>(reconBin, temp1, ss);
1246 for (
int i = 0; i < ss; i++) {
1247 temp1[i] = temp1[i] * temp1[i];
1250 saveAsImg(fnameReconNo, 8, temp2, nx, ny);
1260 int ss = holosize[
_X] * holosize[
_Y];
1266 memsetArr<Real>(
weight, 0.0, 0, ss - 1);
1277 for (
int i = 0; i < ss; i++) {
1292 LOG(
"normalize finishied..\n");
1293 if (encoded ==
nullptr)
1294 encoded =
new Real[ss];
1297 LOG(
"shiftW finishied..\n");
1302 for (
int i = 0; i < ss; i++) {
1303 oph::absCplx<Real>(
freqW[i], absFW);
1304 if (((
Real)i / (
Real)holosize[
_X]) < ((
Real)holosize[
_Y] / 2.0) && absFW < 0.6)
1314 for (
int i = 0; i < ss; i++) {
1323 for (
int i = 0; i < ss; i++) {
1327 int cx = (holosize[
_X] + 1) / 2, cy = (holosize[
_Y] + 1) / 2;
1329 for (
int iy = 0; iy < holosize[
_Y] - Nw[
_Y]; iy++) {
1330 for (
int ix = Nw[
_X]; ix < holosize[
_X] - Nw[
_X]; ix++) {
1332 ii = ix + iy*holosize[
_X];
1333 if (ix >= Nw[
_X] && ix < (holosize[
_X] - Nw[
_X]) && iy < (holosize[
_Y] - Nw[
_Y])) {
1335 if (toBeBin[ii][
_RE] > threshold)
1340 error = toBeBin[ii][
_RE] - encoded[ii];
1342 for (
int iwy = 0; iwy < Nw[
_Y] + 1; iwy++) {
1343 for (
int iwx = -Nw[
_X]; iwx < Nw[
_X] + 1; iwx++) {
1344 iii = (ix + iwx) + (iy + iwy)*holosize[
_X];
1345 jjj = (cx + iwx) + (cy + iwy)*holosize[
_X];
1347 toBeBin[iii] +=
weightC[jjj] * error;
1356 LOG(
"binary finishied..\n");
1365 int cx = (holosize[
_X] + 1) / 2;
1366 int cy = (holosize[
_Y] + 1) / 2;
1372 LOG(
"ERROR DIFFUSION : FLOYD_STEINBERG\n");
1373 weight[(cx + 1) + cy*holosize[
_X]] = 7.0 / 16.0;
1374 weight[(cx - 1) + (cy + 1)*holosize[
_X]] = 3.0 / 16.0;
1375 weight[(cx)+(cy + 1)*holosize[
_X]] = 5.0 / 16.0;
1376 weight[(cx + 1) + (cy + 1)*holosize[
_X]] = 1.0 / 16.0;
1377 Nw[
_X] = 1; Nw[
_Y] = 1;
1380 LOG(
"ERROR DIFFUSION : SINGLE_RIGHT\n");
1381 weight[(cx + 1) + cy*holosize[
_X]] = 1.0;
1382 Nw[
_X] = 1; Nw[
_Y] = 0;
1385 LOG(
"ERROR DIFFUSION : SINGLE_DOWN\n");
1386 weight[cx + (cy + 1)*holosize[
_X]] = 1.0;
1387 Nw[
_X] = 0; Nw[
_Y] = 1;
1390 LOG(
"<FAILED> WRONG PARAMETERS.\n");
1401 int ss = holosize[
_X] * holosize[
_Y];
1406 for (
int i = 0; i < ss; i++) {
1408 x = i % holosize[
_X] - holosize[
_X] / 2; y = i / holosize[
_X] - holosize[
_Y];
1418 for (
int i = 0; i < ss; i++) {
1431 #pragma omp parallel for firstprivate(threshold) 1433 for (
int i = 0; i < size; i++) {
1434 if (src[i][
_RE] > threshold)
1445 const long long int pnXY = pnX * pnY;
1452 int beginY = pnY >> 1;
1453 int beginX = pnX >> 1;
1454 int endY = pnY + beginY;
1455 int endX = pnX + beginX;
1457 for (
int idxnY = beginY; idxnY < endY; idxnY++) {
1458 for (
int idxnX = beginX; idxnX < endX; idxnX++) {
1459 in2x[idxnY * pnX * 2 + idxnX] = in[idxIn++];
1473 for (
int idxFy = -pnY; idxFy < pnY; idxFy++) {
1474 for (
int idxFx = -pnX; idxFx < pnX; idxFx++) {
1482 memsetArr<Complex<Real>>(prop, zero, 0, pnXY * 4 - 1);
1488 for (
int i = 0; i < pnXY * 4; i++) {
1490 prop[i][
_IM] = 2 *
M_PI * distance;
1491 prop[i][
_IM] *= sqrtPart;
1492 temp2[i] = temp1[i] * exp(prop[i]);
1502 for (
int idxNy = pnY / 2; idxNy < pnY + (pnY / 2); idxNy++) {
1503 for (
int idxNx = pnX / 2; idxNx < pnX + (pnX / 2); idxNx++) {
1504 out[idxOut] = temp3[idxNy * pnX * 2 + idxNx];
1524 const int pnXY = pnX * pnY;
1526 const Real ssX = pnX * ppX * 2;
1527 const Real ssY = pnY * ppY * 2;
1528 const Real z = 2 *
M_PI * distance;
1529 const Real v = 1 / (lambda * lambda);
1530 const int hpnX = pnX / 2;
1531 const int hpnY = pnY / 2;
1532 const int pnX2 = pnX * 2;
1533 const int pnY2 = pnY * 2;
1539 #pragma omp parallel for firstprivate(pnX, pnX2, hpnX, hpnY) 1541 for (
int i = 0; i < pnY; i++)
1544 int dst = pnX2 * (i + hpnY) + hpnX;
1551 #pragma omp parallel for firstprivate(ssX, ssY, z, v) 1553 for (
int j = 0; j < pnY2; j++)
1555 Real fy = (-pnY + j) / ssY;
1557 int iWidth = j * pnX2;
1558 for (
int i = 0; i < pnX2; i++)
1560 Real fx = (-pnX + i) / ssX;
1563 Real sqrtPart = sqrt(v - fxx - fyy);
1565 temp[iWidth + i] *= prop.
exp();
1572 #pragma omp parallel for firstprivate(pnX, pnX2, hpnX, hpnY) 1574 for (
int i = 0; i < pnY; i++)
1576 int src = pnX2 * (i + hpnY) + hpnX;
1586 if (x == 0.0 && y == 0.0)
return false;
1588 bool bAxisX = (x == 0.0) ? false :
true;
1589 bool bAxisY = (y == 0.0) ? false :
true;
1596 Real ppY2 = ppY * 2;
1597 Real ppX2 = ppX * 2;
1601 Real *waveRatio =
new Real[nChannel];
1605 for (
uint i = 0; i < nChannel; i++) {
1608 Real ratioX = x * waveRatio[i];
1609 Real ratioY = y * waveRatio[i];
1612 #pragma omp parallel for firstprivate(ppX, ppY, ppX2, ppY2, ssX, ssY, pi2) 1614 for (
int y = 0; y < pnY; y++) {
1617 Real startY = ssY + (ppY * y);
1618 Real shiftY = startY / ppY2 * ratioY;
1619 yy = (pi2 * shiftY).exp();
1621 int offset = y * pnX;
1623 for (
int x = 0; x < pnX; x++) {
1628 Real startX = ssX + (ppX * x);
1629 Real shiftX = startX / ppX2 * ratioX;
1645 const long long int pnXY = pnX * pnY;
1647 Real dfx = 1 / ppX / pnX;
1648 Real dfy = 1 / ppY / pnY;
1650 int beginX = -pnX >> 1;
1651 int endX = pnX >> 1;
1652 int beginY = pnY >> 1;
1653 int endY = -pnY >> 1;
1655 Real carryX = distance * tan(carryingAngleX);
1656 Real carryY = distance * tan(carryingAngleY);
1658 for (
uint ch = 0; ch < nChannel; ch++) {
1660 for (
int i = beginY; i > endY; i--)
1662 Real _carryY = carryY * i * dfy;
1664 for (
int j = beginX; j < endX; j++)
1682 const long long int pnXY = pnX * pnY;
1683 Real dfx = 1 / ppX / pnX;
1684 Real dfy = 1 / ppY / pnY;
1686 int beginX = -pnX >> 1;
1687 int endX = pnX >> 1;
1688 int beginY = pnY >> 1;
1689 int endY = -pnY >> 1;
1693 for (
int i = beginY; i > endY; i--)
1695 Real y = i * dfy * ppY * carryIdxY;
1697 for (
int j = beginX; j < endX; j++)
1699 Real x = j * dfx * ppX * carryIdxX;
1702 dst[k] = src[k] * exp(carrier);
1712 LOG(
"Not found diffracted data.");
1719 int cropx1, cropx2, cropx, cropy1, cropy2, cropy;
1720 if (sig_location[1] == 0) {
1725 cropy = (int)floor(((
Real)pnY) / 2);
1726 cropy1 = cropy - (int)floor(((
Real)cropy) / 2);
1727 cropy2 = cropy1 + cropy - 1;
1730 if (sig_location[0] == 0) {
1735 cropx = (int)floor(((
Real)pnX) / 2);
1736 cropx1 = cropx - (int)floor(((
Real)cropx) / 2);
1737 cropx2 = cropx1 + cropx - 1;
1755 const long long int pnXY = pnX * pnY;
1760 for (
uint ch = 0; ch < nChannel; ch++) {
1764 #pragma omp parallel for firstprivate(cropx1, cropx2, cropy1, cropy2, pnX) 1766 for (
int p = 0; p < pnXY; p++)
1770 if (x >= cropx1 && x <= cropx2 && y >= cropy1 && y <= cropy2)
1781 #pragma omp parallel for firstprivate(sig_location) 1783 for (
int i = 0; i < pnXY; i++) {
1787 m_lpEncoded[ch][i] = (h_crop[i] * shift_phase).real();
1797 const long long int pnXY = pnX * pnY;
1804 cufftDoubleComplex *k_temp_d_, *u_complex_gpu_;
1805 cudaStream_t stream_;
1806 cudaStreamCreate(&stream_);
1808 cudaMalloc((
void**)&u_complex_gpu_,
sizeof(cufftDoubleComplex) * pnXY);
1809 cudaMalloc((
void**)&k_temp_d_,
sizeof(cufftDoubleComplex) * pnXY);
1811 for (
uint ch = 0; ch < nChannel; ch++) {
1812 cudaMemcpy(u_complex_gpu_,
complex_H[ch],
sizeof(cufftDoubleComplex) * pnXY, cudaMemcpyHostToDevice);
1814 cudaMemsetAsync(k_temp_d_, 0,
sizeof(cufftDoubleComplex) * pnXY, stream_);
1815 cudaCropFringe(stream_, pnX, pnY, u_complex_gpu_, k_temp_d_, cropx1, cropx2, cropy1, cropy2);
1817 cudaMemsetAsync(u_complex_gpu_, 0,
sizeof(cufftDoubleComplex) * pnXY, stream_);
1818 cudaFFT(stream_, pnX, pnY, k_temp_d_, u_complex_gpu_, 1,
true);
1820 cudaMemsetAsync(k_temp_d_, 0,
sizeof(cufftDoubleComplex) * pnXY, stream_);
1821 cudaGetFringe(stream_, pnX, pnY, u_complex_gpu_, k_temp_d_, sig_location[0], sig_location[1], ssX, ssY, ppX, ppY,
M_PI);
1823 cufftDoubleComplex* sample_fd = (cufftDoubleComplex*)malloc(
sizeof(cufftDoubleComplex) * pnXY);
1824 memset(sample_fd, 0.0,
sizeof(cufftDoubleComplex) * pnXY);
1826 cudaMemcpyAsync(sample_fd, k_temp_d_,
sizeof(cufftDoubleComplex) * pnXY, cudaMemcpyDeviceToHost, stream_);
1829 for (
int i = 0; i < pnX * pnY; i++)
1832 cudaFree(sample_fd);
1834 cudaStreamDestroy(stream_);
1846 if (sig_location[1] != 0)
1850 Real yy = (ssY / 2.0) - (ppY)*r - ppY;
1853 if (sig_location[1] == 1)
1854 val[
_IM] = 2 *
M_PI * (yy / (4 * ppY));
1856 val[
_IM] = 2 *
M_PI * (-yy / (4 * ppY));
1859 shift_phase_val *= val;
1862 if (sig_location[0] != 0)
1866 Real xx = (-ssX / 2.0) - (ppX)*c - ppX;
1869 if (sig_location[0] == -1)
1870 val[
_IM] = 2 *
M_PI * (-xx / (4 * ppX));
1872 val[
_IM] = 2 *
M_PI * (xx / (4 * ppX));
1875 shift_phase_val *= val;
1883 rand_phase_val[
_RE] = 0.0;
1885 #if REAL_IS_DOUBLE & true 1893 rand_phase_val.
exp();
1897 rand_phase_val[
_RE] = 1.0;
1898 rand_phase_val[
_IM] = 0.0;
1912 template <
typename T>
1916 #pragma omp parallel for 1918 for (
int i = 0; i < size; i++) {
1919 encoded[i] = real(holo[i]);
1923 template <
typename T>
1927 #pragma omp parallel for 1929 for (
int i = 0; i < size; i++) {
1930 encoded[i] = imag(holo[i]);
1934 template <
typename T>
1937 double pi2 =
M_PI * 2;
1939 #pragma omp parallel for firstprivate(pi2) 1941 for (
int i = 0; i < size; i++) {
1942 encoded[i] = (holo[i].
angle() +
M_PI) / pi2;
1946 template <
typename T>
1950 #pragma omp parallel for 1952 for (
int i = 0; i < size; i++) {
1953 encoded[i] = holo[i].
mag();
1957 template <
typename T>
1960 int resize = size >> 1;
1964 #pragma omp parallel for 1966 for (
int i = 0; i < resize; i++) {
1967 normCplx[i] = holo[i * 2];
1970 oph::normalize<T>(normCplx, normCplx, resize);
1972 T* ampl =
new T[resize];
1975 T* phase =
new T[resize];
1976 Phase(normCplx, phase, resize);
1979 #pragma omp parallel for 1981 for (
int i = 0; i < resize; i++) {
1982 T delPhase = acos(ampl[i]);
1983 encoded[i * 2] = (phase[i] +
M_PI) + delPhase;
1984 encoded[i * 2 + 1] = (phase[i] +
M_PI) - delPhase;
1992 template <
typename T>
1995 int resize = size / 3;
1998 #pragma omp parallel for 2000 for (
int i = 0; i < resize; i++) {
2001 norm[i] = holo[i * 3];
2006 T* phase =
new T[resize];
2009 T* ampl =
new T[resize];
2017 #pragma omp parallel for firstprivate(pi2, pi4, sqrt3) 2019 for (
int i = 0; i < resize; i++) {
2021 if (phase[i] >= 0 && phase[i] < (pi2 / 3))
2023 encoded[idx] = ampl[i] * (cos(phase[i]) + sin(phase[i]) / sqrt3);
2024 encoded[idx + 1] = 2 * sin(phase[i]) / sqrt3;
2026 else if (phase[i] >= (pi2 / 3) && phase[i] < (pi4 / 3))
2028 encoded[idx + 1] = ampl[i] * (cos(phase[i] - (pi2 / 3)) + sin(phase[i] - (pi2 / 3)) / sqrt3);
2029 encoded[idx + 2] = 2 * sin(phase[i] - (pi2 / 3)) / sqrt3;
2031 else if (phase[i] >= (pi4 / 3) && phase[i] < (pi2))
2033 encoded[idx + 2] = ampl[i] * (cos(phase[i] - (pi4 / 3)) + sin(phase[i] - (pi4 / 3)) / sqrt3);
2034 encoded[idx] = 2 * sin(phase[i] - (pi4 / 3)) / sqrt3;
2044 template <
typename T>
2047 T* tmp1 =
new T[size];
2049 #pragma omp parallel for 2051 for (
int i = 0; i < size; i++) {
2052 tmp1[i] = holo[i].
mag();
2059 #pragma omp parallel for firstprivate(max) 2061 for (
int i = 0; i < size; i++) {
2062 T tmp = (holo[i] + max).mag();
2063 encoded[i] = tmp * tmp;
2070 for (
int i = 0; i < nVertex; i++)
2081 for (
int i = 0; i < nVertex; i++) {
2082 dst[i] = -fieldLens * src[i] / (src[i] - fieldLens);
2093 #pragma omp parallel for firstprivate(x, y, z) 2095 for (
int i = 0; i < nSize; i++) {
2096 dst[i + 0] = src[i + 0] * x;
2097 dst[i + 1] = src[i + 1] * y;
2098 dst[i + 2] = src[i + 2] * z;
2107 for (
int i = 0; i < len; i++) {
2108 if (src[i] > maxTmp)
2110 if (src[i] < minTmp)
2123 const int pnXY = width * height;
2127 int scaleX = round(width * 4 * waveRatio);
2128 int scaleY = round(height * 4 * waveRatio);
2131 int hh = height * 4;
2132 int nSize = ww * hh;
2133 int nScaleSize = scaleX * scaleY;
2138 memset(img_tmp, 0,
sizeof(
uchar) * nSize);
2140 int h1 = round((hh - scaleY) / 2);
2141 int w1 = round((ww - scaleX) / 2);
2143 #pragma omp parallel for firstprivate(w1, h1, scaleX, ww) 2145 for (
int y = 0; y < scaleY; y++) {
2146 for (
int x = 0; x < scaleX; x++) {
2147 img_tmp[(y + h1)*ww + x + w1] = imgScaled[y*scaleX + x];
2165 for (
uint i = 0; i < nChannel; i++) {
2176 for (
uint i = 0; i < nChannel; i++) {
bool loadPLY(const std::string &fileName, ulonglong &n_points, Vertex **vertices)
void abs(const oph::Complex< T > &src, oph::Complex< T > &dst)
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.
virtual bool saveAsImg(const char *fname, uint8_t bitsperpixel, uchar *src, int width, int height)
Function for creating image files.
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.
void setResolution(ivec2 resolution)
Function for setting buffer size.
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)
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.
void binarization(Complex< Real > *src, Real *dst, const int size, int ENCODE_FLAG, Real threshold)
void GetRandomPhaseValue(Complex< Real > &rand_phase_val, bool rand_phase)
Assign random phase value if random_phase == 1.
ulonglong n_points
Number of points.
void setPixelNumberOHC(const ivec2 pixel_number)
getter/setter for OHC file read and write
Complex< Real > * weightC
Real ** m_lpEncoded
buffer to encoded.
void initialize(void)
Initialize variables for Hologram complex field, encoded data, normalized data.
void RealPart(Complex< T > *holo, T *encoded, const int size)
Encoding method.
void setPixelPitchOHC(const vec2 pixel_pitch)
bool shiftW(ivec2 holosize)
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 TwoPhase(Complex< T > *holo, T *encoded, const int size)
bool checkExtension(const char *fname, const char *ext)
Functions for extension checking.
structure for 2-dimensional integer vector and its arithmetic.
const XMLNode * FirstChild() const
Get the first child node, or null if none exists.
void normalize(const Complex< T > *src, Complex< T > *dst, const int &size)
Normalize all elements of Complex<T>* src from 0 to 1.
void Burckhardt(Complex< T > *holo, T *encoded, const int size)
void setWavelengthOHC(const Real wavelength, const LenUnit wavelength_unit)
ImgEncoderOhc * OHC_encoder
OHC file format Variables for read and write.
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)
void normalize(void)
Normalization function to save as image file after hologram creation.
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.
void singleSideBand(Complex< Real > *holo, Real *encoded, const ivec2 holosize, int passband)
Encoding method.
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.
void getShiftPhaseValue(Complex< Real > &shift_phase_val, int idx, ivec2 sig_location)
Calculate the shift phase value.
void imgScaleBilinear(uchar *src, uchar *dst, int w, int h, int neww, int newh, int channels=1)
Function for change image size.
void * load(const char *fname)
Function for loading image files.
void SimpleNI(Complex< T > *holo, T *encoded, const int size)
void CorrectionChromaticAberration(uchar *src, uchar *dst, int width, int height, int ch)
void waveCarry(Real carryingAngleX, Real carryingAngleY, Real distance)
Wave carry.
void fft2(ivec2 n, Complex< Real > *in, int sign=OPH_FORWARD, uint flag=OPH_ESTIMATE)
Functions for performing fftw 2-dimension operations inside Openholo.
structure for 2-dimensional Real type vector and its arithmetic.
void fresnelPropagation(OphConfig context, Complex< Real > *in, Complex< Real > *out, Real distance)
Fresnel propagation.
#define ELAPSED_TIME(x, y)
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 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)
bool saveRefImages(char *fnameW, char *fnameWC, char *fnameAS, char *fnameSSB, char *fnameHP, char *fnameFreq, char *fnameReal, char *fnameBin, char *fnameReconBin, char *fnameReconErr, char *fnameReconNo)
Vertex * vertices
Data of point clouds.
virtual bool loadAsOhc(const char *fname)
Function to read OHC file.
bool readConfig(const char *fname)
load to configuration file.
XMLError LoadFile(const char *filename)
virtual uchar * loadAsImg(const char *fname)
Function for loading image files.
Real rand(const Real min, const Real max, oph::ulong _SEED_VALUE=0)
Get random Real value from min to max.
void getAmplitude(oph::Complex< Real > *src, Real *dst, const int &size)
void ScaleChange(Real *src, Real *dst, int nSize, Real scaleX, Real scaleY, Real scaleZ)
uchar ** m_lpNormalized
buffer to normalized.
Complex< Real > * AS
Binary Encoding - Error diffusion.
Real minOfArr(const std::vector< Real > &arr)
Real maxOfArr(const std::vector< Real > &arr)
void resetBuffer()
reset buffer
void setPixelNumber(ivec2 n)
Function for setting the output resolution.
Complex< Real > * normalized
void Fresnel_Diffraction(Point src, Complex< Real > *dst, Real lambda, Real distance, Real amplitude)
Fresnel-diffraction method.
Complex< Real > ** complex_H
virtual bool loadAsOhc(const char *fname)
Function to read OHC file.
void transVW(int nVertex, Vertex *dst, Vertex *src)
int ENCODE_METHOD
Encoding method flag.
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.
const XMLElement * FirstChildElement(const char *name=0) const
bool mergeColor(int idx, int width, int height, uchar *src, uchar *dst)
Function for generate RGB image from each grayscale image.
virtual void ophFree(void)
Pure virtual function for override in child classes.
virtual void ophFree(void)
Pure virtual function for override in child classes.
virtual ~ophGen(void)=0
Destructor.