279 ->
typename std::remove_reference<
decltype(data[0])>::type {
280 typedef typename std::remove_reference<
decltype(data[0])>::type DataT;
282 size_t index0, index0_hi, index1, index1_hi, index2, index2_hi;
284 if (!prepareInterpolation(point, p, index0, index0_hi, index1, index1_hi, index2, index2_hi, flags))
285 return NaN<
decltype(data[0])>();
287 Vec<3> pa = fullMesh.at(index0, index1, index2);
289 size_t step0 = (p.c0 < pa.c0)?
290 (index0 == 0)? 0 : -1 :
291 (index0_hi == fullMesh.axis[0]->size())? 0 : 1;
292 size_t step1 = (p.c1 < pa.c1)?
293 (index1 == 0)? 0 : -1 :
294 (index1_hi == fullMesh.axis[1]->size())? 0 : 1;
295 size_t step2 = (p.c2 < pa.c2)?
296 (index2 == 0)? 0 : -1 :
297 (index2_hi == fullMesh.axis[2]->size())? 0 : 1;
299 size_t index_aaa = index(index0, index1, index2), index_aab, index_aba, index_abb,
300 index_baa, index_bab, index_bba, index_bbb;
302 typename std::remove_const<DataT>::type data_aaa = data[index_aaa], data_aab, data_aba, data_abb,
303 data_baa, data_bab, data_bba, data_bbb;
305 if (step0 == 0 && step1 == 0 && step2 == 0) {
306 index_aab = index_aba = index_abb = index_baa = index_bab = index_bba = index_bbb = index_aaa;
307 data_aab = data_aba = data_abb = data_baa = data_bab = data_bba = data_bbb = data_aaa;
309 index_aab = index(index0, index1, index2+step2);
310 index_aba = index(index0, index1+step1, index2);
311 index_abb = index(index0, index1+step1, index2+step2);
312 index_baa = index(index0+step0, index1, index2);
313 index_bab = index(index0+step0, index1, index2+step2);
314 index_bba = index(index0+step0, index1+step1, index2);
315 index_bbb = index(index0+step0, index1+step1, index2+step2);
316 data_aab = (index_aab != Element::UNKNOWN_ELEMENT_INDEX)? data[index_aab] : data_aaa;
317 data_aba = (index_aba != Element::UNKNOWN_ELEMENT_INDEX)? data[index_aba] : data_aaa;
318 data_baa = (index_baa != Element::UNKNOWN_ELEMENT_INDEX)? data[index_baa] : data_aaa;
319 data_abb = (index_abb != Element::UNKNOWN_ELEMENT_INDEX)? data[index_abb] : data_aab + data_aba - data_aaa;
320 data_bab = (index_bab != Element::UNKNOWN_ELEMENT_INDEX)? data[index_bab] : data_baa + data_aab - data_aaa;
321 data_bba = (index_bba != Element::UNKNOWN_ELEMENT_INDEX)? data[index_bba] : data_aba + data_baa - data_aaa;
322 data_bbb = (index_bbb != Element::UNKNOWN_ELEMENT_INDEX)? data[index_bbb] : data_aab + data_aba + data_baa - 2. * data_aaa;
325 Vec<3> pb = fullMesh.at(index0+step0, index1+step1, index2+step2);
326 if (step0 == 0) pb.c0 += 1.;
327 if (step1 == 0) pb.c1 += 1.;
328 if (step2 == 0) pb.c2 += 1.;
330 return flags.postprocess(point,
331 interpolation::trilinear(pa.c0, pb.c0, pa.c1, pb.c1, pa.c2, pb.c2,
332 data_aaa, data_baa, data_bba, data_aba,
333 data_aab, data_bab, data_bbb, data_abb,
345 ->
typename std::remove_reference<
decltype(data[0])>::type {
347 std::size_t index0_lo, index0_hi, index1_lo, index1_hi, index2_lo, index2_hi;
349 if (!originalMesh->prepareInterpolation(point, wrapped_point, index0_lo, index0_hi, index1_lo, index1_hi, index2_lo, index2_hi, flags))
350 return NaN<
decltype(data[0])>();
352 return flags.postprocess(point, data[this->index(index0_lo, index1_lo, index2_lo)]);