TTK
Loading...
Searching...
No Matches
DepthImageBasedGeometryApproximation.h
Go to the documentation of this file.
1
21
22#pragma once
23
24// base code includes
25#include <Debug.h>
26
27// std includes
28#include <cmath>
29#include <limits>
30
31namespace ttk {
32
34
35 public:
37 this->setDebugMsgPrefix("DIBGA");
38 }
40
48 template <class dataType, class idType>
49 int execute(
50 // output
51 float *pointCoordinates,
52 double *triangleDistortions,
53 idType *connectivityList,
54 idType *offsetArray,
55
56 // input
57 const dataType *depthValues,
58 const double *camPos,
59 const double *camDir,
60 const double *camUp,
61 const double *camNearFar,
62 const double *camHeight,
63 const double *resolution) const;
64 };
65} // namespace ttk
66
67template <class dataType, class idType>
69 // output
70 float *pointCoordinates,
71 double *triangleDistortions,
72 idType *connectivityList,
73 idType *offsetArray,
74
75 // input
76 const dataType *depthValues,
77 const double *camPos,
78 const double *camDir,
79 const double *camUp,
80 const double *camNearFar,
81 const double *camHeight,
82 const double *resolution) const {
83
84 const double camDirMag = std::sqrt(
85 camDir[0] * camDir[0] + camDir[1] * camDir[1] + camDir[2] * camDir[2]);
86 const double camDirN[3]{
87 camDir[0] / camDirMag, camDir[1] / camDirMag, camDir[2] / camDirMag};
88
91 this->printMsg({{"Resolution", std::to_string((int)resolution[0]) + "x"
92 + std::to_string((int)resolution[1])},
93 {"CamPos", "[" + std::to_string(camPos[0]) + ","
94 + std::to_string(camPos[1]) + ","
95 + std::to_string(camPos[2]) + "]"},
96 {"CamDir", "[" + std::to_string(camDirN[0]) + ","
97 + std::to_string(camDirN[1]) + ","
98 + std::to_string(camDirN[2]) + "]"},
99 {"CamHeight", std::to_string(camHeight[0])},
100 {"CamNearFar", "[" + std::to_string(camNearFar[0]) + ","
101 + std::to_string(camNearFar[1]) + "]"}},
105
106 Timer timer;
107
108 const size_t resolutionST[2] = {(size_t)resolution[0], (size_t)resolution[1]};
109
110 this->printMsg("Processing image (" + std::to_string(resolutionST[0]) + "x"
111 + std::to_string(resolutionST[1]) + ")",
113 // -------------------------------------------------------------------------
114 // Compute Camera Vectors
115 // -------------------------------------------------------------------------
116
117 // Compute camera size
118 const double camSize[2]
119 = {resolution[0] / resolution[1] * camHeight[0], camHeight[0]};
120
121 // Compute camRight = camDirN x CamUp
122 double camRight[3] = {camDirN[1] * camUp[2] - camDirN[2] * camUp[1],
123 camDirN[2] * camUp[0] - camDirN[0] * camUp[2],
124 camDirN[0] * camUp[1] - camDirN[1] * camUp[0]};
125 double temp = sqrt(camRight[0] * camRight[0] + camRight[1] * camRight[1]
126 + camRight[2] * camRight[2]);
127 camRight[0] /= temp;
128 camRight[1] /= temp;
129 camRight[2] /= temp;
130
131 // Compute true up vector
132 double camUpTrue[3]
133 = {camDirN[1] * (-camRight[2]) - camDirN[2] * (-camRight[1]),
134 camDirN[2] * (-camRight[0]) - camDirN[0] * (-camRight[2]),
135 camDirN[0] * (-camRight[1]) - camDirN[1] * (-camRight[0])};
136 temp = sqrt(camUpTrue[0] * camUpTrue[0] + camUpTrue[1] * camUpTrue[1]
137 + camUpTrue[2] * camUpTrue[2]);
138 camUpTrue[0] /= temp;
139 camUpTrue[1] /= temp;
140 camUpTrue[2] /= temp;
141
142 // -------------------------------------------------------------------------
143 // Create Vertices
144 // -------------------------------------------------------------------------
145 {
146 // Compute pixel size in world coordinates
147 const double pixelWidthWorld = camSize[0] / resolution[0];
148 const double pixelHeightWorld = camSize[1] / resolution[1];
149
150 // Optimization: precompute half of the camera size to reduce the number of
151 // operations in the for loop Include a half pixel offset (-0.5) to center
152 // vertices at pixel centers
153 const double camWidthWorldHalf = 0.5 * camSize[0] - 0.5 * pixelWidthWorld;
154 const double camHeightWorldHalf = 0.5 * camSize[1] - 0.5 * pixelHeightWorld;
155
156 // Compute depth delta
157 const double delta = camNearFar[1] - camNearFar[0];
158
159 // Optimization: reorient camera model to bottom left corner to reduce
160 // operations in for loop
161 const double camPosCorner[3] = {camPos[0] - camRight[0] * camWidthWorldHalf
162 - camUpTrue[0] * camHeightWorldHalf,
163 camPos[1] - camRight[1] * camWidthWorldHalf
164 - camUpTrue[1] * camHeightWorldHalf,
165 camPos[2] - camRight[2] * camWidthWorldHalf
166 - camUpTrue[2] * camHeightWorldHalf};
167
168// Compute vertex coordinates while parallelizing over rows
169#ifdef TTK_ENABLE_OPENMP
170#pragma omp parallel for num_threads(threadNumber_)
171#endif
172 for(size_t y = 0; y < resolutionST[1]; y++) {
173 const double v = ((double)y) * pixelHeightWorld;
174 const double vTimesUp[3]
175 = {v * camUpTrue[0], v * camUpTrue[1], v * camUpTrue[2]};
176
177 const size_t yOffset = y * resolutionST[0];
178 for(size_t x = 0; x < resolutionST[0]; x++) {
179 const size_t pixelIndex = x + yOffset;
180
181 // double d = (double)(depthValues[ pixelIndex ])*delta+camNearFar[0];
182 const double depth = ((double)depthValues[pixelIndex]);
183 // double d = depth > 0.98 ? 0 : depth * delta + camNearFar[0];
184 const double d = depth * delta + camNearFar[0];
185 const double u = ((double)x) * pixelWidthWorld;
186
187 // compute vertex coordinate
188 const size_t pointCoordinateOffset = pixelIndex * 3;
189 pointCoordinates[pointCoordinateOffset]
190 = camPosCorner[0] + u * camRight[0] + vTimesUp[0] + d * camDirN[0];
191 pointCoordinates[pointCoordinateOffset + 1]
192 = camPosCorner[1] + u * camRight[1] + vTimesUp[1] + d * camDirN[1];
193 pointCoordinates[pointCoordinateOffset + 2]
194 = camPosCorner[2] + u * camRight[2] + vTimesUp[2] + d * camDirN[2];
195 }
196 }
197 }
198
199 // -------------------------------------------------------------------------
200 // Create Triangles
201 // -------------------------------------------------------------------------
202 {
203 auto absDiff = [](const dataType &a, const dataType &b) {
204 return a > b ? a - b : b - a;
205 };
206 auto isNaN = [](const double &a) { return std::isnan(a) || a >= 1.0; };
207
208 const size_t nTriangles = 2 * (resolution[0] - 1) * (resolution[1] - 1);
209
210 for(size_t t = 0; t < nTriangles; t++)
211 offsetArray[t] = t * 3;
212 offsetArray[nTriangles] = 3 * nTriangles;
213
214 /* Index Structure:
215 0 - 1
216 | / |
217 2 - 3
218 */
219 const size_t xl = resolutionST[0] - 1;
220 const size_t yl = resolutionST[1] - 1;
221 const size_t trianglesPerRow = xl * 2;
222
223 const double myNan = std::numeric_limits<double>::quiet_NaN();
224
225#ifdef TTK_ENABLE_OPENMP
226#pragma omp parallel for num_threads(this->threadNumber_)
227#endif
228 for(size_t y = 0; y < yl; y++) {
229 size_t yOffset = y * resolutionST[0];
230 size_t triangleIndexOffset = y * trianglesPerRow * 3;
231 size_t triangleDistortionOffset = y * trianglesPerRow;
232
233 for(size_t x = 0; x < xl; x++) {
234 size_t i0 = x + yOffset;
235 size_t i1 = i0 + 1;
236 size_t i2 = i0 + resolutionST[0];
237 size_t i3 = i2 + 1;
238
239 connectivityList[triangleIndexOffset++] = i0;
240 connectivityList[triangleIndexOffset++] = i2;
241 connectivityList[triangleIndexOffset++] = i1;
242
243 connectivityList[triangleIndexOffset++] = i1;
244 connectivityList[triangleIndexOffset++] = i2;
245 connectivityList[triangleIndexOffset++] = i3;
246
247 const double i0Depth = (double)depthValues[i0];
248 const double i1Depth = (double)depthValues[i1];
249 const double i2Depth = (double)depthValues[i2];
250 const double i3Depth = (double)depthValues[i3];
251
252 // wow
253 triangleDistortions[triangleDistortionOffset++]
254 = isNaN(i0Depth) || isNaN(i2Depth) || isNaN(i1Depth)
255 ? myNan
256 : std::max(
257 absDiff(i0Depth, i1Depth),
258 std::max(absDiff(i1Depth, i2Depth), absDiff(i0Depth, i2Depth)));
259
260 triangleDistortions[triangleDistortionOffset++]
261 = isNaN(i1Depth) || isNaN(i2Depth) || isNaN(i3Depth)
262 ? myNan
263 : std::max(
264 absDiff(i1Depth, i3Depth),
265 std::max(absDiff(i3Depth, i2Depth), absDiff(i2Depth, i1Depth)));
266 }
267 }
268 }
269
270 // Print performance
271 this->printMsg("Processing image (" + std::to_string(resolutionST[0]) + "x"
272 + std::to_string(resolutionST[1]) + ")",
273 1, timer.getElapsedTime(), this->threadNumber_);
274
275 return 1;
276}
long long int idType
int threadNumber_
Definition: BaseClass.h:95
Minimalist debugging class.
Definition: Debug.h:88
void setDebugMsgPrefix(const std::string &prefix)
Definition: Debug.h:364
int printMsg(const std::string &msg, const debug::Priority &priority=debug::Priority::INFO, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cout) const
Definition: Debug.h:118
This module approximates the depicted geometry of a depth image.
~DepthImageBasedGeometryApproximation() override=default
int execute(float *pointCoordinates, double *triangleDistortions, idType *connectivityList, idType *offsetArray, const dataType *depthValues, const double *camPos, const double *camDir, const double *camUp, const double *camNearFar, const double *camHeight, const double *resolution) const
double getElapsedTime()
Definition: Timer.h:15
The Topology ToolKit.