source: GTP/trunk/Lib/Vis/Preprocessing/src/VssRay.cpp @ 1942

#include "VssRay.h"
#include "AxisAlignedBox3.h"

namespace GtpVisibilityPreprocessor {

// Static variables
int
VssRay::mailID = 0;

#define OLD_PARAM 0

VssRay::VssRay(
                 const Vector3 &origin,
                 const Vector3 &termination,
                 Intersectable *originObject,
                 Intersectable *terminationObject,
                 const int pass,
                 const float pdf
                 ):
        mMailbox(-1),
    mOrigin(origin),
    mTermination(termination),
        mOriginObject(originObject),
        mTerminationObject(terminationObject),
        mRefCount(0),
    mFlags(0),
        mPass(pass),
        mViewCells(0),
        mWeightedPvsContribution(0),
        mPdf(pdf),
        mTerminationNode(NULL),
        mOriginNode(NULL),
        mPvsContribution(0)
{
        Precompute();
}


VssRay::VssRay(const Ray &ray):
        mRefCount(0),
        mFlags(0),
        mMailbox(-1),
        mOriginObject(ray.sourceObject.mObject),
        mPass(0),
        mViewCells(0),
        mPdf(1.0f),
        mTerminationNode(NULL),
        mOriginNode(NULL),
        mPvsContribution(0)
{
        if (ray.sourceObject.mObject)
                mOrigin = ray.Extrap(ray.sourceObject.mT);
        else
                mOrigin = ray.GetLoc(); 

        //Debug << "origin: " << mOrigin << endl;
        if (!ray.intersections.empty())
        {
                mTermination = ray.Extrap(ray.intersections[0].mT);
                mTerminationObject = ray.intersections[0].mObject;
        }
        else
        {
                mTermination = ray.Extrap(1e6);//TODO: should be Limits::Infinity
                mTerminationObject = NULL;
        }

        Precompute();
}


void VssRay::Precompute() 
{
        mFlags = 0;
        Vector3 dir = GetDir();
//mFlags |= BorderSample;
#define BIDIRECTIONAL_RAY 0
#if BIDIRECTIONAL_RAY
        if (dir.y < 0) {
                // swap objects and poits       
                swap(mOriginObject, mTerminationObject);
                swap(mOrigin, mTermination);
                dir = -dir;
        }
#endif
        if (dir.x > 0.0f) mFlags |= FPosDirX;
        if (dir.y > 0.0f) mFlags |= FPosDirY;
        if (dir.z > 0.0f) mFlags |= FPosDirZ;

        mInvSize = 1.0f/Magnitude(GetDir());
}


bool
VssRay::ComputeMinMaxT(const AxisAlignedBox3 &box,
                                           float &tmin,
                                           float &tmax) const
{
        
  float minx, maxx, miny, maxy, minz, maxz;
  Vector3 dir = GetDir();
  
  if (fabs(dir.x) < 0.001) {
    if (box.Min().x < GetOrigin().x && box.Max().x > GetOrigin().x) {
      minx = -MAXFLOAT;
      maxx = MAXFLOAT;
    }
    else
      return false;
  }
  else {
    float t1 = (box.Min().x - GetOrigin().x) / dir.x;
    float t2 = (box.Max().x - GetOrigin().x) / dir.x;
    if (t1 < t2) {
      minx = t1;
      maxx = t2;
    }
    else {
      minx = t2;
      maxx = t1;
    }
    if (maxx < 0)
      return false;
  }

  if (fabs(dir.y) < 0.001) {
    if (box.Min().y < GetOrigin().y && box.Max().y > GetOrigin().y) {
      miny = -MAXFLOAT;
      maxy = MAXFLOAT;
    }
    else
      return false;
  }
  else {
    float t1 = (box.Min().y - GetOrigin().y) / dir.y;
    float t2 = (box.Max().y - GetOrigin().y) / dir.y;
    if (t1 < t2) {
      miny = t1;
      maxy = t2;
    }
    else {
      miny = t2;
      maxy = t1;
    }
    if (maxy < 0.0)
      return false;
  }
        
  if (fabs(dir.z) < 0.001) {
    if (box.Min().z < GetOrigin().z && box.Max().z > GetOrigin().z) {
      minz = -MAXFLOAT;
      maxz = MAXFLOAT;
    }
    else
      return false;
  }
  else {
    float t1 = (box.Min().z - GetOrigin().z) / dir.z;
    float t2 = (box.Max().z - GetOrigin().z) / dir.z;
    if (t1 < t2) {
      minz = t1;
      maxz = t2;
    }
    else {
      minz = t2;
      maxz = t1;
    }
    if (maxz < 0.0)
      return false;
  }

  tmin = minx;
  if (miny > tmin)
    tmin = miny;
  if (minz > tmin)
    tmin = minz;

  tmax = maxx;
  if (maxy < tmax)
    tmax = maxy;
  if (maxz < tmax)
    tmax = maxz;

  return 1; // yes, intersection was found
}



bool
VssRay::Intersects(const AxisAlignedBox3 &box,
                                   float &tmin,
                                   float &tmax) const
{
  if (!ComputeMinMaxT(box, tmin, tmax))
    return false;
        
  if ( tmax < tmin)
    return false; // the ray passes outside the box
  
  if ( tmax < 0.0f)
    return false; // the intersection is not on the positive halfline
  
  return true; // ray hits the box .. origin can be outside or inside the box
}


bool
VssRay::IntersectsSphere(const Vector3 &center,
                                                 const float sqrRadius,
                                                 Vector3 &point,
                                                 float &t) const
{
  // compute ray/plane intersection
  Vector3 dir = GetDir();
        
  t = -DotProd(GetOrigin() - center, dir)*sqr(GetInvSize());
        
  if (t < 0.0f)
    t = 0.0f;
  else
    if (t > 1.0f)
      t = 1.0f;
  
  point = GetOrigin() + t*dir;
  
  return GtpVisibilityPreprocessor::SqrDistance(point, center) < sqrRadius;
}


float
VssRay::GetDirParam(const int axis, const Vector3 dir)
{
  Vector3 d = Normalize(dir);
#if OLD_PARAM
  return (axis == 0) ? atan2(dir.x, dir.z) : asin(dir.y);
#else
  //  x = cos(p0)sin(p1)
  //  y = cos(p1)
  //  z = sin(p0)sin(p1)
  return (axis == 0) ? atan2(dir.z, dir.x) : acos(dir.y);
#endif
}

Vector3 VssRay::GetInvDirParam(const float alpha, const float beta)
{
#if OLD_PARAM
  return Normalize(Vector3(sin(alpha), sin(beta), cos(alpha)));

#else
  return Normalize(Vector3(cos(alpha)*sin(beta),
                                                   cos(beta),
                                                   sin(alpha)*sin(beta)));
#endif
}

float
VssRay::GetDirParametrization(const int axis) const
{
  Vector3 dir = GetNormalizedDir();
  return GetDirParam(axis, dir);
}

float
VssRay::GetOpositeDirParametrization(const int axis) const
{
  Vector3 dir = -GetNormalizedDir();
  return GetDirParam(axis, dir);
}


void
GenerateExtendedConvexCombinationWeights2(float &w1,
                                                                                  float &w2,
                                                                                  const float overlap
                                                                                  )
{
  w1 = RandomValue(-overlap, 1.0f + overlap);
  w2 = 1.0f - w1;
}

void
GenerateExtendedConvexCombinationWeights(float &w1,
                                                                                 float &w2,
                                                                                 float &w3,
                                                                                 const float overlap
                                                                                 )
{

  while (1) {
        w1 = RandomValue(-overlap, 1.0f + overlap);
        w2 = RandomValue(-overlap, 1.0f + overlap);
        w3 = 1.0f + overlap - w1 - w2;
        if (w3 >= -overlap && w3 <= 1.0f + overlap)
          break;
        /*
        w3 = RandomValue(0.0f, 1.0f + overlap);
        
        float sum = w1 + w2 + w3;
        if (sum>0.0f) {
          w1/=sum;
          w2/=sum;
          w3/=sum;
          break;
        }
        */
  }
}


void
VssRayContainer::PrintStatistics(ostream &s)
{

  VssRayContainer::const_iterator it = begin(), it_end = end();
  int sumContributions = 0;
  float sumRelContributions = 0.0f;
  
  for (; it != it_end; ++it) {
        VssRay *ray = *it;
        sumContributions += ray->mPvsContribution;
        sumRelContributions += ray->mRelativePvsContribution;
  }
  
  s<<"##### VSS RAY STAT ##########\n";
  s<<"#RAYS\n"<<(int)size()<<endl;
  s<<"#AVG_RAY_PVS_CONTRIBUTION\n"<<sumContributions/(float)size()<<endl;
  s<<"#AVG_RAY_RELATIVE_PVS_CONTRIBUTION\n"<<sumRelContributions/
        (float)size()<<endl;
  
}


int
VssRayContainer::SelectRays(const int number,
                                                        VssRayContainer &selected,
                                                        const bool copy) const
{
  float p = number/(float)size();
  VssRayContainer::const_iterator it = begin();
  
  for (; it != end(); it++)
        if (Random(1.0f) < p) {
          if (!copy)
                selected.push_back(*it);
          else
                selected.push_back(new VssRay(**it));
        }
  
  return (int)selected.size();
}


int
VssRayContainer::GetContributingRays(VssRayContainer &selected,
                                                                         const int minPass
                                                                         ) const
{
  VssRayContainer::const_iterator it = begin();
  //  ofstream s("contrib.log");
  for (; it != end(); it++) {
        VssRay *ray = *it;
        //      cout<<"(pass="<<ray->mPass<<", c="<<ray->mPvsContribution<<")"<<endl;
        if (ray->mPass >= minPass && ray->mPvsContribution > 0)
          selected.push_back(ray);
  }
  
  return (int)selected.size();
}



}