Openholo  v4.2
Open Source Digital Holographic Library
PointCloud.py
Go to the documentation of this file.
1 from numba import njit, types, typed
2 import numpy as np
3 import matplotlib.pyplot as plt
4 import plyfile
5 
6 # unit parameters
7 mm = 1e-3
8 um = mm * mm
9 nm = um * mm
10 
11 @njit
12 def k(wvl):
13  return (2 * np.pi) / wvl
14 
15 def get_pointcloud(path):
16  # converting point cloud data to numba type
17  ans = typed.List()
18  with open(path, 'rb') as f:
19  plydata = plyfile.PlyData.read(f)
20  plydata = np.asarray(plydata.elements[1].data)
21 
22  for i in range(len(plydata)):
23  raw = typed.List()
24  [raw.append(types.float64(n)) for n in plydata[i]]
25  ans.append(raw)
26  return ans
27 
28 
29 @njit
30 def RSIntegral(plydata, color, propagation=1, angleX=0, angleY=0, Red=639 * nm, Green=525 * nm, Blue=463 * nm,
31  SLM_width=3840, SLM_height=2160, scaleXY=0.03, scaleZ=0.25, pixel_pitch=3.6 * um):
32  """Fresnel approximation point propagation"""
33  # define shift
34  if color == 'green':
35  num_color = 4
36  wvl = Green
37  elif color == 'blue':
38  num_color = 5
39  wvl = Blue
40  else:
41  wvl = Red
42  num_color = 3
43  shiftx = np.tan(angleX * (np.pi / 180)) * propagation
44  shifty = np.tan(angleY * (np.pi / 180)) * propagation
45  CHR = np.zeros((SLM_height, SLM_width))
46  CHI = np.zeros((SLM_height, SLM_width))
47  num = 0
48  for n in plydata:
49  num += 1
50  ch_r = np.zeros((SLM_height, SLM_width))
51  ch_i = np.zeros((SLM_height, SLM_width))
52  x1 = (n[0] + shiftx) * scaleXY
53  y1 = (n[1] + shifty) * scaleXY
54  z1 = n[2] * scaleZ
55  z = propagation - z1
56  amp = n[num_color] * (z / wvl)
57  for i in np.arange(SLM_height):
58  for j in np.arange(SLM_width):
59  x2 = (j - SLM_width / 2) * pixel_pitch
60  y2 = -(i - SLM_height / 2) * pixel_pitch
61  r = np.sqrt((x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2) + z * z)
62  t = (wvl * r) / (2 * pixel_pitch)
63  if (x1 - t < x2 < x1 + t) and (y1 - t < y2 < y1 + t): # anti aliasing
64  ch_r[i, j] = np.sin(k(wvl) * r) / r ** 2
65  ch_i[i, j] = np.cos(k(wvl) * r) / r ** 2
66  CHR += ch_r * amp
67  CHI += ch_i * amp
68  print(num, 'th point', color, ' done')
69  return (CHR + 1j * CHI)
70 
71 
72 @njit
73 def FresnelIntegral(plydata, color, propagation=1, angleX=0, angleY=0, Red=639 * nm, Green=525 * nm, Blue=463 * nm,
74  SLM_width=3840, SLM_height=2160, scaleXY=0.03, scaleZ=0.25, pixel_pitch=3.6 * um):
75  """Fresnel approximation point propagation"""
76  # define shift
77  if color == 'green':
78  num_color = 4
79  wvl = Green
80  elif color == 'blue':
81  num_color = 5
82  wvl = Blue
83  else:
84  wvl = Red
85  num_color = 3
86  shiftx = np.tan(angleX * (np.pi / 180)) * propagation
87  shifty = np.tan(angleY * (np.pi / 180)) * propagation
88  CHR = np.zeros((SLM_height, SLM_width))
89  CHI = np.zeros((SLM_height, SLM_width))
90  num = 0
91  for n in plydata:
92  num += 1
93  ch_r = np.zeros((SLM_height, SLM_width))
94  ch_i = np.zeros((SLM_height, SLM_width))
95  x1 = (n[0] + shiftx) * scaleXY
96  y1 = (n[1] + shifty) * scaleXY
97  z1 = n[2] * scaleZ
98  z = propagation - z1
99  amp = n[num_color] * (z / wvl)
100  for i in np.arange(SLM_height):
101  for j in np.arange(SLM_width):
102  x2 = (j - SLM_width / 2) * pixel_pitch
103  y2 = -(i - SLM_height / 2) * pixel_pitch
104  r = ((x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2)) / (2 * z)
105  t = (wvl * z) / (2 * pixel_pitch)
106  if (x1 - t < x2 < x1 + t) and (y1 - t < y2 < y1 + t): # anti aliasing
107  ch_r[i, j] = np.cos(k(wvl) * r)
108  ch_i[i, j] = np.sin(k(wvl) * r)
109  CHR += ch_r * amp
110  CHI += ch_i * amp
111  print(num, 'th point done')
112  return (CHR + 1j * CHI)
113 
114 
def get_pointcloud(path)
Definition: PointCloud.py:15
def FresnelIntegral(plydata, color, propagation=1, angleX=0, angleY=0, Red=639 *nm, Green=525 *nm, Blue=463 *nm, SLM_width=3840, SLM_height=2160, scaleXY=0.03, scaleZ=0.25, pixel_pitch=3.6 *um)
Definition: PointCloud.py:74
def RSIntegral(plydata, color, propagation=1, angleX=0, angleY=0, Red=639 *nm, Green=525 *nm, Blue=463 *nm, SLM_width=3840, SLM_height=2160, scaleXY=0.03, scaleZ=0.25, pixel_pitch=3.6 *um)
Definition: PointCloud.py:31