SampleVolumeTask_Vanilla.cpp

#if 0
#include "SampleVolumeTask.h"
 
using namespace Cogs::Core;
 
using glm::uvec3;
using glm::ivec3;
using glm::vec3;
using glm::clamp;
using glm::max;
using glm::min;
using glm::floor;
using glm::ceil;
using glm::inverse;
 
void EchoSounder::SampleVolumeTask2::runVanilla(TileValue * dstValues,
                                                uint16_t * timeValues,
                                                const EchoSounder::PingInfo* pinf,
                                                const uint32_t upperFansToRemove,
                                                const uint32_t numPings,
                                                const uint32_t numSamples,
                                                const uint32_t gridSizeX,
                                                const uint32_t gridSizeY,
                                                const float depthOffset,
                                                const float depthStep,
                                                const float sampleSpacing,
                                                const float decay,
                                                const std::vector<float>& beamAngleAlongship,
                                                const std::vector<float>& beamAngleAthwartship,
                                                const glm::ivec3& this_minIndex,
                                                const glm::ivec3& this_maxIndex)
{
  const auto time = pinf->time;
  const auto weight = 1.f;// / data.numPings;
  const auto numBeamsMajor = static_cast<uint32_t>(beamAngleAlongship.size());
  const auto numBeamsMinor = static_cast<uint32_t>(beamAngleAthwartship.size());
  const auto s = 1.f / sampleSpacing;
  //const auto width = data.cellWidth;
 
  const vec3 minPolar(beamAngleAthwartship.front(),
                      beamAngleAlongship.front(),
                      depthOffset);
  const vec3 maxPolar(beamAngleAthwartship.back(),
                      beamAngleAlongship.back(),
                      depthOffset + (numSamples - 1)*depthStep);
 
 
  const uvec3 maxTau(max(2u, numBeamsMinor) - 2,
                     max(2u + upperFansToRemove, numBeamsMajor) - 2 - upperFansToRemove,
                     max(2u, numSamples) - 2);
 
  // map minPolar to 0 and maxPolar to N-1.
  const vec3 polarToIndex = vec3(max(2u, numBeamsMinor) - 2,
                                 max(2u, numBeamsMajor) - 2,
                                 max(2u, numSamples) - 2) / (maxPolar - minPolar);
 
  const auto * minorMu = beamAngleAthwartship.data();
  const auto * majorMu = beamAngleAlongship.data();
 
  for (uint32_t p = 0; p < numPings; p++) {
 
    const auto minIndex = this_minIndex;// max(this_minIndex, ivec3(floor(s * pinf[p].boundingBoxTile[0])));
    const auto maxIndex = this_maxIndex;// min(this_maxIndex, ivec3(ceil(s * pinf[p].boundingBoxTile[1])));
    const auto * srcValues = pinf[p].field;
 
    const auto rot = inverse(pinf[p].metaPing.arrayOrientationGlobal);
    const auto shift = pinf[p].arrayPositionTile;
    for (int k = minIndex.z; k < maxIndex.z; k++) {
      for (int j = minIndex.y; j < maxIndex.y; j++) {
        for (int i = minIndex.x; i < maxIndex.x; i++) {
 
          const auto samplePosGrid = sampleSpacing*vec3(i, j, k);
          if (samplePosGrid.z < pinf[p].depthRestriction) {
            continue;
          }
 
          const auto samplePosArray = (samplePosGrid - shift);
          const auto r2 = glm::dot(samplePosArray, samplePosArray);
 
          bool in = true;
 
          // Early exit
          bool in_0 = 0.f <= glm::dot(pinf[p].boundingFrustumNormals[0], samplePosArray);
          bool in_1 = 0.f <= glm::dot(pinf[p].boundingFrustumNormals[1], samplePosArray);
          bool in_2 = 0.f <= glm::dot(pinf[p].boundingFrustumNormals[2], samplePosArray);
          bool in_3 = 0.f <= glm::dot(pinf[p].boundingFrustumNormals[3], samplePosArray);
          bool in_4 = pinf[p].minDepthSquared <= r2;
          bool in_5 = r2 <= pinf[p].maxDepthSquared;
          if (!in_0 || !in_1 || !in_2 || !in_3 || !in_4 || !in_5) {
            continue;
          }
          const auto cartesianPos = rot * samplePosArray;
          float r = std::sqrt(r2);
 
          const float dirX = std::asin(cartesianPos.x / r);
          const float dirY = std::asin(cartesianPos.y / r);
          const auto polarPos = glm::vec3(dirY, dirX, r);
 
          auto xi = polarToIndex*(polarPos - minPolar);
          auto tau = uvec3(xi);
 
          bool range_1_0 = tau.x <= maxTau.x; // Negative is also handled here due to unsigned type.
          bool range_1_1 = tau.y <= maxTau.y;
          bool range_1_2 = tau.z <= maxTau.z;
          if (!range_1_0 || !range_1_1 || !range_1_2) {
            continue;
          }
 
          // Bi-linear interpolation, nearest nb along depth.
          float bx = xi.x - tau.x;
          float by = xi.y - tau.y;
 
          const auto val00 = srcValues[((tau.y + 0)*numBeamsMinor + (tau.x + 0))*numSamples + tau.z];
          const auto val01 = srcValues[((tau.y + 1)*numBeamsMinor + (tau.x + 0))*numSamples + tau.z];
          const auto val0 = (1.f - by)*PingValue_Decode(val00) + by*PingValue_Decode(val01);
 
          const auto val10 = srcValues[((tau.y + 0)*numBeamsMinor + (tau.x + 1))*numSamples + tau.z];
          const auto val11 = srcValues[((tau.y + 1)*numBeamsMinor + (tau.x + 1))*numSamples + tau.z];
          const auto val1 = (1.f - by)*PingValue_Decode(val10) + by*PingValue_Decode(val11);
 
          const auto val = (1.f - bx)*val0 + bx*val1;
 
          auto dstOffset = (k*gridSizeY + j)*gridSizeX + i;
 
          auto oldValue = TileValue_Decode(dstValues[dstOffset]);
          if (oldValue == 0.0) oldValue = val;
 
          dstValues[dstOffset] = TileValue_Encode(0.5f*oldValue + 0.5f*val);
          //dstValues[dstOffset] = max(oldValue, val);
          timeValues[dstOffset] = time;
        }
      }
    }
  }
}
#endif