Field3D
MIPUtil.h
Go to the documentation of this file.
1//----------------------------------------------------------------------------//
2
3/*
4 * Copyright (c) 2013 Sony Pictures Imageworks Inc
5 *
6 * All rights reserved.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions
10 * are met:
11 *
12 * Redistributions of source code must retain the above copyright
13 * notice, this list of conditions and the following disclaimer.
14 * Redistributions in binary form must reproduce the above copyright
15 * notice, this list of conditions and the following disclaimer in the
16 * documentation and/or other materials provided with the
17 * distribution. Neither the name of Sony Pictures Imageworks nor the
18 * names of its contributors may be used to endorse or promote
19 * products derived from this software without specific prior written
20 * permission.
21 *
22 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
23 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
24 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
25 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
26 * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
27 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
28 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
29 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
31 * STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
32 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
33 * OF THE POSSIBILITY OF SUCH DAMAGE.
34 */
35
36//----------------------------------------------------------------------------//
37
42//----------------------------------------------------------------------------//
43
44#ifndef _INCLUDED_Field3D_MIPUtil_H_
45#define _INCLUDED_Field3D_MIPUtil_H_
46
47//----------------------------------------------------------------------------//
48
49#include <vector>
50
51#include <boost/thread/thread.hpp>
52#include <boost/thread/mutex.hpp>
53
54#include "Resample.h"
55#include "SparseField.h"
56#include "Types.h"
57
58//----------------------------------------------------------------------------//
59
60#include "ns.h"
61
63
64//----------------------------------------------------------------------------//
65// Functions
66//----------------------------------------------------------------------------//
67
69V3i computeOffset(const FieldRes &f);
70
75template <typename MIPField_T, typename Filter_T>
76typename MIPField_T::Ptr
77makeMIP(const typename MIPField_T::NestedType &base, const int minSize,
78 const V3i &offset, const size_t numThreads);
79
81template <typename MIPField_T, typename Filter_T>
82typename MIPField_T::Ptr
83makeMIP(const typename MIPField_T::NestedType &base, const int minSize,
84 const size_t numThreads);
85
86//----------------------------------------------------------------------------//
87// Implementation details
88//----------------------------------------------------------------------------//
89
90namespace detail {
91
92 //--------------------------------------------------------------------------//
93
94 extern const std::string k_mipOffsetStr;
95
96 //--------------------------------------------------------------------------//
97
99 template <typename T>
101 {
102 typedef T type;
103 };
104
106 template <>
107 struct ComputationType<Field3D::half>
108 {
109 typedef float type;
110 };
111
112 //--------------------------------------------------------------------------//
113
114 FIELD3D_API V3i mipResolution(const V3i &baseRes, const size_t level,
115 const V3i &add);
116
117 //--------------------------------------------------------------------------//
118
120 template <typename Data_T>
122 {
123 return 16;
124 }
125
127 template <typename Data_T>
129 {
130 return f.blockSize();
131 }
132
133 //--------------------------------------------------------------------------//
134
135 template <typename Data_T>
137 const SparseField<Data_T> &/*tgt*/,
138 const Box3i &tgtBox, const float support,
139 const size_t dim)
140 {
141 const int intSupport = static_cast<int>(std::ceil(support * 0.5));
142 const int pad = std::max(0, intSupport);
143 Box3i tgtBoxPad = tgtBox;
144 tgtBoxPad.min[dim] -= pad;
145 tgtBoxPad.max[dim] += pad;
146 Box3i srcBoxPad = tgtBoxPad;
147 srcBoxPad.min[dim] *= 2;
148 srcBoxPad.max[dim] *= 2;
149
150 // Get the block coordinates
151 const Box3i dbsBounds = blockCoords(clipBounds(srcBoxPad, src.dataWindow()),
152 &src);
153
154 static boost::mutex mutex;
155 boost::mutex::scoped_lock lock(mutex);
156
157 // Check all blocks
158 for (int k = dbsBounds.min.z; k <= dbsBounds.max.z; ++k) {
159 for (int j = dbsBounds.min.y; j <= dbsBounds.max.y; ++j) {
160 for (int i = dbsBounds.min.x; i <= dbsBounds.max.x; ++i) {
161 if (src.blockIsAllocated(i, j, k) ||
162 src.getBlockEmptyValue(i, j, k) != static_cast<Data_T>(0)) {
163 return false;
164 }
165 }
166 }
167 }
168
169 // No hits. Empty
170 return true;
171 }
172
173 //--------------------------------------------------------------------------//
174
176 template <typename Field_T>
177 bool checkInputEmpty(const Field_T &/*src*/, const Field_T &/*tgt*/,
178 const Box3i &/*tgtBox*/, const float /*support*/,
179 const size_t /*dim*/)
180 {
181 return false;
182 }
183
184 //--------------------------------------------------------------------------//
185
186 template <typename Field_T, typename FilterOp_T, bool IsAnalytic_T>
188 {
189 typedef typename Field_T::value_type T;
190
191 MIPSeparableThreadOp(const Field_T &src, Field_T &tgt,
192 const size_t level, const V3i &add,
193 const FilterOp_T &filterOp,
194 const size_t dim,
195 const std::vector<Box3i> &blocks,
196 size_t &nextIdx, boost::mutex &mutex)
197 : m_src(src),
198 m_tgt(tgt),
199 m_filterOp(filterOp),
200 m_level(level),
201 m_add(add),
202 m_dim(dim),
203 m_blocks(blocks),
204 m_nextIdx(nextIdx),
205 m_mutex(mutex),
206 m_numBlocks(blocks.size())
207 {
208 // Empty
209 }
210
212 {
213 using namespace std;
214
215 // Defer to ComputationType to determine the processing data type
216 typedef typename Field_T::value_type Data_T;
217 typedef typename ComputationType<Data_T>::type Value_T;
218
219 // To ensure we don't sample outside source data
220 Box3i srcDw = m_src.dataWindow();
221
222 // Coordinate frame conversion constants
223 const float tgtToSrcMult = 2.0;
224 const float filterCoordMult = 1.0f / (tgtToSrcMult);
225
226 // Filter info, support size in target space
227 const float support = m_filterOp.support();
228
229 // Get next index to process
230 size_t idx;
231 {
232 boost::mutex::scoped_lock lock(m_mutex);
233 idx = m_nextIdx;
234 m_nextIdx++;
235 }
236 // Keep going while there is data to process
237 while (idx < m_numBlocks) {
238 // Grab the bounds
239 const Box3i box = m_blocks[idx];
240 // Early exit if input blocks are all empty
241 if (!detail::checkInputEmpty(m_src, m_tgt, box, support, m_dim)) {
242 // For each output voxel
243 for (int k = box.min.z; k <= box.max.z; ++k) {
244 for (int j = box.min.y; j <= box.max.y; ++j) {
245 for (int i = box.min.x; i <= box.max.x; ++i) {
246 Value_T accumValue(m_filterOp.initialValue());
247 if (IsAnalytic_T) {
248 // Transform from current point in target frame to source frame
249 const int curTgt = V3i(i, j, k)[m_dim];
250 const float curSrc = discToCont(curTgt) * tgtToSrcMult - m_add[m_dim];
251 // Find interval
252 int startSrc =
253 static_cast<int>(std::floor(curSrc - support * tgtToSrcMult));
254 int endSrc =
255 static_cast<int>(std::ceil(curSrc + support *
256 tgtToSrcMult)) - 1;
257 // Clamp coordinates
258 startSrc = std::max(startSrc, srcDw.min[m_dim]);
259 endSrc = std::min(endSrc, srcDw.max[m_dim]);
260 // Loop over source voxels
261 for (int s = startSrc; s <= endSrc; ++s) {
262 // Source index
263 const int xIdx = m_dim == 0 ? s : i;
264 const int yIdx = m_dim == 1 ? s : j;
265 const int zIdx = m_dim == 2 ? s : k;
266 // Source voxel in continuous coords
267 const float srcP = discToCont(s);
268 // Compute filter weight in source space (twice as wide)
269 const float weight = m_filterOp.eval(std::abs(srcP - curSrc) *
270 filterCoordMult);
271 // Value
272 const Value_T value = m_src.fastValue(xIdx, yIdx, zIdx);
273 // Update
274 if (weight > 0.0f) {
275 FilterOp_T::op(accumValue, value);
276 }
277 }
278 // Update final value
279 if (accumValue !=
280 static_cast<Value_T>(m_filterOp.initialValue())) {
281 m_tgt.fastLValue(i, j, k) = accumValue;
282 }
283 } else {
284 float accumWeight = 0.0f;
285 // Transform from current point in target frame to source frame
286 const int curTgt = V3i(i, j, k)[m_dim];
287 const float curSrc = discToCont(curTgt) * tgtToSrcMult - m_add[m_dim];
288 // Find interval
289 int startSrc =
290 static_cast<int>(std::floor(curSrc - support * tgtToSrcMult));
291 int endSrc =
292 static_cast<int>(std::ceil(curSrc + support *
293 tgtToSrcMult)) - 1;
294 // Clamp coordinates
295 startSrc = std::max(startSrc, srcDw.min[m_dim]);
296 endSrc = std::min(endSrc, srcDw.max[m_dim]);
297 // Loop over source voxels
298 for (int s = startSrc; s <= endSrc; ++s) {
299 // Source index
300 const int xIdx = m_dim == 0 ? s : i;
301 const int yIdx = m_dim == 1 ? s : j;
302 const int zIdx = m_dim == 2 ? s : k;
303 // Source voxel in continuous coords
304 const float srcP = discToCont(s);
305 // Compute filter weight in source space (twice as wide)
306 const float weight = m_filterOp.eval(std::abs(srcP - curSrc) *
307 filterCoordMult);
308 // Value
309 const Value_T value = m_src.fastValue(xIdx, yIdx, zIdx);
310 // Update
311 accumWeight += weight;
312 accumValue += value * weight;
313 }
314 // Update final value
315 if (accumWeight > 0.0f &&
316 accumValue != static_cast<Value_T>(0.0)) {
317 m_tgt.fastLValue(i, j, k) = accumValue / accumWeight;
318 }
319 } // if (IsAnalytic_T)
320 }
321 }
322 }
323 } // Empty input
324 // Get next index
325 {
326 boost::mutex::scoped_lock lock(m_mutex);
327 idx = m_nextIdx;
328 m_nextIdx++;
329 }
330 }
331 }
332
333 private:
334
335 // Data members ---
336
337 const Field_T &m_src;
338 Field_T &m_tgt;
339 const FilterOp_T &m_filterOp;
340 const size_t m_level;
341 const V3i &m_add;
342 const size_t m_dim;
343 const std::vector<Box3i> &m_blocks;
344 size_t &m_nextIdx;
345 boost::mutex &m_mutex;
346 const size_t m_numBlocks;
347
348 };
349
350 //--------------------------------------------------------------------------//
351
353 template <typename Field_T, typename FilterOp_T>
354 void mipSeparable(const Field_T &src, Field_T &tgt,
355 const V3i &oldRes, const V3i &newRes, const size_t level,
356 const V3i &add, const FilterOp_T &filterOp,
357 const size_t dim, const size_t numThreads)
358 {
359 using namespace std;
360
361 // To ensure we don't sample outside source data
362 Box3i srcDw = src.dataWindow();
363
364 // Compute new res
365 V3i res;
366 if (dim == 2) {
367 res = newRes;
368 } else if (dim == 1) {
369 res = V3i(newRes.x, newRes.y, oldRes.z);
370 } else {
371 res = V3i(newRes.x, oldRes.y, oldRes.z);
372 }
373
374 // Resize new field
375 tgt.setSize(res);
376
377 // Determine granularity
378 const size_t blockSize = threadingBlockSize(src);
379
380 // Build block list
381 std::vector<Box3i> blocks;
382 for (int k = 0; k < res.z; k += blockSize) {
383 for (int j = 0; j < res.y; j += blockSize) {
384 for (int i = 0; i < res.x; i += blockSize) {
385 Box3i box;
386 // Initialize block size
387 box.min = V3i(i, j, k);
388 box.max = box.min + V3i(blockSize - 1);
389 // Clip against resolution
390 box.max.x = std::min(box.max.x, res.x - 1);
391 box.max.y = std::min(box.max.y, res.y - 1);
392 box.max.z = std::min(box.max.z, res.z - 1);
393 // Add to list
394 blocks.push_back(box);
395 }
396 }
397 }
398
399 // Next index counter and mutex
400 size_t nextIdx = 0;
401 boost::mutex mutex;
402
403 // Launch threads ---
404
405 boost::thread_group threads;
406
407 for (size_t i = 0; i < numThreads; ++i) {
408 threads.create_thread(
410 (src, tgt, level, add, filterOp,
411 dim, blocks, nextIdx, mutex));
412 }
413
414 // Join
415 threads.join_all();
416 }
417
418 //--------------------------------------------------------------------------//
419
420 template <typename Field_T, typename FilterOp_T>
421 void mipResample(const Field_T &base, const Field_T &src, Field_T &tgt,
422 const size_t level, const V3i &offset,
423 const FilterOp_T &filterOp,
424 const size_t numThreads)
425 {
426 using std::ceil;
427
428 // Odd-numbered offsets need a pad of one in the negative directions
429 const V3i add((offset.x % 2 == 0) ? 0 : 1,
430 (offset.y % 2 == 0) ? 0 : 1,
431 (offset.z % 2 == 0) ? 0 : 1);
432
433 // Compute new res
434 const Box3i baseDw = base.dataWindow();
435 const V3i baseRes = baseDw.size() + V3i(1);
436 const V3i newRes = mipResolution(baseRes, level, add);
437
438 // Source res
439 const Box3i srcDw = src.dataWindow();
440 const V3i srcRes = srcDw.size() + V3i(1);
441
442 // Temporary field for y component
443 Field_T tmp;
444
445 // X axis (src into tgt)
446 mipSeparable(src, tgt, srcRes, newRes, level, add, filterOp, 0, numThreads);
447 // Y axis (tgt into temp)
448 mipSeparable(tgt, tmp, srcRes, newRes, level, add, filterOp, 1, numThreads);
449 // Z axis (temp into tgt)
450 mipSeparable(tmp, tgt, srcRes, newRes, level, add, filterOp, 2, numThreads);
451
452 // Update final target with mapping and metadata
453 tgt.name = base.name;
454 tgt.attribute = base.attribute;
455 tgt.setMapping(base.mapping());
456 tgt.copyMetadata(base);
457 }
458
459 //--------------------------------------------------------------------------//
460
463 const V3i &baseRes,
464 const Box3i &extents,
465 const size_t level);
466
467 //--------------------------------------------------------------------------//
468
469} // namespace detail
470
471//----------------------------------------------------------------------------//
472// Function implementations
473//----------------------------------------------------------------------------//
474
475template <typename MIPField_T, typename Filter_T>
476typename MIPField_T::Ptr
477makeMIP(const typename MIPField_T::NestedType &base, const int minSize,
478 const size_t numThreads)
479{
480 // By default, there is no offset
481 const V3i zero(0);
482 // Call out to perform actual work
483 return makeMIP<MIPField_T, Filter_T>(base, minSize, zero, numThreads);
484}
485
486//----------------------------------------------------------------------------//
487
488template <typename MIPField_T, typename Filter_T>
489typename MIPField_T::Ptr
490makeMIP(const typename MIPField_T::NestedType &base, const int minSize,
491 const V3i &baseOffset, const size_t numThreads)
492{
493 using namespace Field3D::detail;
494
495 typedef typename MIPField_T::value_type Data_T;
496 typedef typename MIPField_T::NestedType Src_T;
497 typedef typename Src_T::Ptr SrcPtr;
498 typedef typename MIPField_T::Ptr MIPPtr;
499 typedef std::vector<typename Src_T::Ptr> SrcVec;
500
501 if (base.extents() != base.dataWindow()) {
502 return MIPPtr();
503 }
504
505 // Initialize output vector with base resolution
506 SrcVec result;
507 result.push_back(field_dynamic_cast<Src_T>(base.clone()));
508
509 // Iteration variables
510 V3i res = base.extents().size() + V3i(1);
511 V3i offset = baseOffset;
512
513 // Loop until minimum size is found
514 size_t level = 1;
515 while ((res.x > minSize || res.y > minSize || res.z > minSize) &&
516 (res.x > 2 && res.y > 2 && res.z > 2)) {
517 // Perform filtering
518 SrcPtr nextField(new Src_T);
519 mipResample(base, *result.back(), *nextField, level, offset,
520 Filter_T(), numThreads);
521 // Add to vector of filtered fields
522 result.push_back(nextField);
523 // Set up for next iteration
524 res = nextField->dataWindow().size() + V3i(1);
525 // ... offset needs to be rounded towards negative inf, not towards zero
526 for (int i = 0; i < 3; ++i) {
527 if (offset[i] < 0) {
528 offset[i] = (offset[i] - 1) / 2;
529 } else {
530 offset[i] /= 2;
531 }
532 }
533 level++;
534 }
535
536 MIPPtr mipField(new MIPField_T);
537 mipField->name = base.name;
538 mipField->attribute = base.attribute;
539 mipField->copyMetadata(base);
540 mipField->setMIPOffset(baseOffset);
541 mipField->setup(result);
542
543 return mipField;
544}
545
546//----------------------------------------------------------------------------//
547
549
550//----------------------------------------------------------------------------//
551
552#endif // Include guard
Box3i clipBounds(const Box3i &bbox, const Box3i &bounds)
Definition Field.h:1145
double discToCont(int discCoord)
Goes from discrete coordinates to continuous coordinates See Graphics Gems - What is a pixel.
Definition Field.h:1070
MIPField_T::Ptr makeMIP(const typename MIPField_T::NestedType &base, const int minSize, const V3i &offset, const size_t numThreads)
Constructs a MIP representation of the given field, with optional offset vector. The offset vector in...
Definition MIPUtil.h:490
FIELD3D_NAMESPACE_OPEN V3i computeOffset(const FieldRes &f)
Computes the origin/offset of a field.
Definition MIPUtil.cpp:157
Contains functions for resampling fields.
Contains the SparseField class.
Box3i blockCoords(const Box3i &dvsBounds, const SparseField< Data_T > *f)
Imath::V3i V3i
Definition SpiMathLib.h:71
FIELD3D_NAMESPACE_OPENtypedef ::half half
Definition SpiMathLib.h:64
Imath::Box3i Box3i
Definition SpiMathLib.h:77
Contains typedefs for the commonly used types in Field3D.
This subclass of Field stores data in a contiguous std::vector.
Definition DenseField.h:87
boost::intrusive_ptr< FieldMapping > Ptr
const Box3i & dataWindow() const
Returns the data window. Any coordinate inside this window is safe to pass to value() in the Field su...
Definition Field.h:253
This Field subclass stores voxel data in block-allocated arrays.
const Data_T getBlockEmptyValue(int bi, int bj, int bk) const
Returns the constant value of an block, whether it's allocated already or not..
bool blockIsAllocated(int bi, int bj, int bk) const
Checks if a block is allocated.
int blockSize() const
Returns the block size.
Field_T::Ptr field_dynamic_cast(RefBase::Ptr field)
Dynamic cast that uses string-comparison in order to be safe even after an object crosses a shared li...
Definition RefCount.h:256
const std::string k_mipOffsetStr
Definition MIPUtil.cpp:66
void mipSeparable(const Field_T &src, Field_T &tgt, const V3i &oldRes, const V3i &newRes, const size_t level, const V3i &add, const FilterOp_T &filterOp, const size_t dim, const size_t numThreads)
Threaded implementation of separable MIP filtering.
Definition MIPUtil.h:354
size_t threadingBlockSize(const DenseField< Data_T > &)
Constant size for all dense fields.
Definition MIPUtil.h:121
FIELD3D_API V3i mipResolution(const V3i &baseRes, const size_t level, const V3i &add)
Definition MIPUtil.cpp:70
void mipResample(const Field_T &base, const Field_T &src, Field_T &tgt, const size_t level, const V3i &offset, const FilterOp_T &filterOp, const size_t numThreads)
Definition MIPUtil.h:421
FIELD3D_API FieldMapping::Ptr adjustedMIPFieldMapping(const FieldRes *base, const V3i &baseRes, const Box3i &extents, const size_t level)
Definition MIPUtil.cpp:82
bool checkInputEmpty(const SparseField< Data_T > &src, const SparseField< Data_T > &, const Box3i &tgtBox, const float support, const size_t dim)
Definition MIPUtil.h:136
#define FIELD3D_API
Definition ns.h:77
#define FIELD3D_NAMESPACE_HEADER_CLOSE
Definition ns.h:58
Used to delegate the choice of bit depth to process at.
Definition MIPUtil.h:101
MIPSeparableThreadOp(const Field_T &src, Field_T &tgt, const size_t level, const V3i &add, const FilterOp_T &filterOp, const size_t dim, const std::vector< Box3i > &blocks, size_t &nextIdx, boost::mutex &mutex)
Definition MIPUtil.h:191
const std::vector< Box3i > & m_blocks
Definition MIPUtil.h:343
Field_T::value_type T
Definition MIPUtil.h:189
const FilterOp_T & m_filterOp
Definition MIPUtil.h:339
const Field_T & m_src
Definition MIPUtil.h:337