DoorSubmergence.h

#pragma once
 
#include <algorithm>
#include <limits>
 
namespace door_submergence
{
 
inline double PolygonSubmergedFraction(const double u[], const double v[], const double s[], int n)
{
   bool anyWet = false, allWet = true;
   for (int i = 0; i < n; ++i) {
      anyWet = anyWet || s[i] > 0.0;
      allWet = allWet && s[i] > 0.0;
   }
   if (allWet)
      return 1.0;
   if (!anyWet)
      return 0.0;
 
   double cu[16], cv[16];
   int    m = 0;
   for (int i = 0; i < n; ++i) {
      const int j = (i + 1) % n;
      if (s[i] > 0.0) {
         cu[m]   = u[i];
         cv[m++] = v[i];
      }
      if ((s[i] > 0.0) != (s[j] > 0.0)) {
         const double t = s[i] / (s[i] - s[j]);
         cu[m]          = u[i] + t * (u[j] - u[i]);
         cv[m++]        = v[i] + t * (v[j] - v[i]);
      }
   }
   double wetArea = 0.0, area = 0.0;
   for (int i = 0; i < m; ++i)
      wetArea += cu[i] * cv[(i + 1) % m] - cu[(i + 1) % m] * cv[i];
   for (int i = 0; i < n; ++i)
      area += u[i] * v[(i + 1) % n] - u[(i + 1) % n] * v[i];
   return std::clamp(wetArea / area, 0.0, 1.0);
}
 
template<class Environment>
double DoorWetFraction(Environment& environment, double time, const double origin[3], const double cornerBody[4][3],
      const double cornerNed[4][3], double waveBound)
{
   double top = std::numeric_limits<double>::infinity();
   for (int i = 0; i < 4; ++i)
      top = std::min(top, cornerNed[i][2]);
   if (top > waveBound)
      return 1.0;
   double elevation = 0.0;
   environment.GetSurfaceElevation(time, origin, elevation);
   double u[4], v[4], s[4];
   for (int i = 0; i < 4; ++i) {
      u[i] = cornerBody[i][0];
      v[i] = cornerBody[i][2];
      s[i] = cornerNed[i][2] + elevation;
   }
   return PolygonSubmergedFraction(u, v, s, 4);
}
 
} // namespace door_submergence