#include "config.h"
#include "oren_nayar.h"
#include "texture.h"
#include "functions.h"

__BEGIN_BRDFLIB_SOURCE

// describes the surface to a given stream
//========================================================================
void
CSurfaceOrenNayar::Describe(ostream &app, int ind) const
{
	indent(app, ind);
	app << GetName() << ":\n";

	if (_diffuse) {
		indent(app, ind+1);
		app << "  diffuse  = \t";
		_diffuse->Describe(app, ind); app << "\n";
	}


	if (_slope) {
		indent(app, ind+1);
		app << "  slope    = \t";
		_slope->Describe(app, ind); app << "\n";
	}

  app << flush;

  return;
}

/*
  Integral of this function cannot be evaluated analitically. However, it can
  be evaluated using Taylor's series. Unfortunately, it requres a lot of
  computations. So, for a now we leave this albedo function in this way.

  Possibly later on, we'll provide this integral evaluation in the way
  similar to the Arvo's procedures (Arvo'1995 "Applications of Irradiance
  Tensors to the Simulation of Non-Lambertian Phenomena").
*/
//========================================================================
bool
CSurfaceOrenNayar::GetAlbedo(CValues &albedo,
                          const CHitPointInfo &hitInfo,
                          EDistributionType   distType) const
{
  //CVector3D N       = hitInfo->GetNormal();
  //CVector3D dirInit = hitInfo->GetRay()->GetDir();

  // As default for all the Microfacet-based models. Actually, for each
  // model albedo is a function of the angle bewtween the incoming
  // direction and normal.
  CValues c(CValues::EE_Weigth3); // result - zero by default
  CValues d(CValues::EE_Weigth3);
  bool retval = false;

  if ((distType & EE_REFLECTED_DIFFUSE) && _diffuse) {
    hitInfo.GetTextureValue(d, *_diffuse, CMappingBase::EE_KD);
    c += d;
    retval = true;
  }

  if (retval) {
    albedo = c;
    return true;
  }

  return false;
}


//========================================================================

bool
CSurfaceOrenNayar::EvaluateBSDF(CValues &result,
			 const CHitPointInfo &hitInfo,
			 const CVector3D     &dirIn,
			 EDistributionType   distType) const
{
	bool retval = false;
	float d = Magnitude(dirIn);
	CVector3D N          =  hitInfo.GetNormal();
	CVector3D anisAxis   =  hitInfo.GetAnisAxis();

	assert(fabs(DotProd(N, anisAxis)) < 1e-06); // pi/2

	float lightCos       =  DotProd(dirIn, N);
	float incEyeCos  =  -(hitInfo.GetDotProd());

	if(lightCos >= -1e-06 && incEyeCos >= -1e-06)
	{
		if(lightCos < 0) lightCos = 0;
		else
		{
			if(lightCos < WTINY) lightCos       = WTINY;
		}

		if(incEyeCos < 0) incEyeCos = 0;
		else
		{
			if(incEyeCos < WTINY) incEyeCos  = WTINY;
		}

		if ( (d < WTINY1L) || (d > WTINY1G) )
		{
			Normalize(dirIn);
		}

		if(_diffuse && (distType & EE_REFLECTED_DIFFUSE) && _slope) 
		{
			CValues rho_d;
			hitInfo.GetTextureValue(rho_d, *_diffuse, CMappingBase::EE_KD);

			// -------------------------------------------------------------------------------
			// in practice, simplified version of computationally expensive Oren-Nayar BRDF
			// is used (according to Schwenk, K. "A Survey of Shading Models for Real-time 
			// Rendering", 2008
			// note : in Neumann phd thesis there is a mistake in computation of A term, where
			// in denominator is 0.5+sigma^2, but it has to be 0.33+sigma^2
			//
			// fON = kd/Pi * (A + B * max(0,(cos(phiV-phiL))) * sin(alfa) * tan(beta))
			//
			//                         max(0,(LV + (NL) * (NV)))
			//     = kd/Pi * (A + B * ---------------------------
			//								max(NL,NV)
			//
			// where
			//
			//		    0.5*sigma^2			0.45*sigma^2
			// A = 1 - ------------- , B = --------------
			//		    0.33+sigma^2		0.09+sigma^2
			//
			// where sigma = standard deviation of surface slope, 
			//		alfa = max(thetaV,thetaL)
			//		beta = min(thetaV,thetaL)
			//		  LV = dotproduct between light and viewer direction
			//		  NL = dotproduct between normal and light direction
			//		  NV = dotproduct betweem normal and viewer direction
			// --------------------------------------------------------------------------------

			float sigma, sigma_2;
			CTexture *slopeX = _slope->GetSlopeX();
			if(slopeX)
			{
				CValues slp;
				hitInfo.GetTextureValue(slp, *slopeX, CMappingBase::EE_KSI);
				sigma = slp.Get(0);
			}
			sigma_2 = sigma*sigma;

			float A = 1 - (K05*sigma_2) / (K033+sigma_2);
			float B = (K045*sigma_2) / (K009+sigma_2);

			float LV = DotProd(dirIn, -hitInfo.GetRay()->GetDir());
			float num, denom; // numerator and denominator of the fraction
			
			num = Max(0.0f,(LV + lightCos*incEyeCos));
			denom = Max(lightCos,incEyeCos);
			
			result = rho_d * M_1_PI * (A + B * num/denom);
			retval = true;
		}

		if(retval)
		{
			result *= lightCos;
		}
	} // not under the surface
	else // under the surface
	{
		retval = false;
	}
	return retval;
}


//========================================================================


int
CSurfaceOrenNayar::SpawnStochasticRays(CContributionArray  &carray,
							const CHitPointInfo &hitInfo,
                            CSamplingInfo       &samplingInfo,
							EDistributionType	distType) const
{
	int maxRays = samplingInfo.GetMaxSamples();
	assert(_diffuse);
	assert(_slope);

	if ( !(distType & EE_REFLECTED_DIFFUSE) ||
	   (maxRays < 1) )
	return 0; // no rays

	// store the seed to remain the same if possible, depending
	// on the random generator.
	samplingInfo.StoreSeed();

	// build coordinate system
	CVector3D N = hitInfo.GetNormal();
	if (hitInfo.GetDotProd() > 0.0)
	{
		N = -N;
	}

	// spawn the rays using stratification
	int spawned =
	SpawnStochasticRaysStratify(carray, hitInfo, EE_REFLECTED_DIFFUSE,
				samplingInfo, N);

	samplingInfo.IncDimension(2);

	samplingInfo.RecoverSeed();

#ifdef __TRAVERSAL_STATISTICS
	world->_statistics->reflectedRays += spawned;
#endif

	return spawned;
}

//========================================================================

int
CSurfaceOrenNayar::SpawnStochasticRaysStratify(CContributionArray  &carray,
			    const CHitPointInfo &hitInfo,
			    EDistributionType   distType,
			    CSamplingInfo       &samplingInfo,
			    const CVector3D &N) const
{
    ERayType rayType = REFLECTED_RAY;


	// build the coordinate system
	CVector3D U = ArbitraryNormal(N);     // U is perpendicular to N
	CVector3D V = CrossProd(N,U);         // orthonormal base with N and U
	CVector3D dirN;

	int maxRays = samplingInfo.GetMaxSamples();

	CValues kd_local( kd(hitInfo) );

	for(int i = 0; i < maxRays ; i++) {
	CContributionItem *contrib = carray.AllocateNext();

	if (i)
	  // increase the seed not to generate the same random sample again
	  samplingInfo.IncSeed();

	// sample the hemisphere in canonical coordinates with density 
	// cos(Phi)/M_PI
	// 
	// Sets the probability for generated direction
	CVector2D vv(samplingInfo.Sample(0), samplingInfo.Sample(1));

	contrib->prob = RND::hemisph_cos(dirN, vv);

	// transform from local coordinates
	CVector3D dirRes = dirN.x * U + dirN.y * V + dirN.z * N ;
	// set the reflected ray, direction is normalized here
	contrib->ray.Init(hitInfo.GetPoint(), dirRes, rayType);
	// set the flags concerning ray type
	contrib->SetFlags(CContributionItem::EE_ReflectedRay);
	// type of interaction
	contrib->interaction = distType;
	// diffuseness of the bounce
	contrib->diffuseness = 1.0;
	// Evaluate the value of BSDF for given input and generated
	// output direction
	/* weight = BSDF * cos(Light, N) / PDF =
			  = pd / M_PI * cos(Light, N) / (cos(Light, N)/M_PI) = pd    */
	contrib->weigth = kd_local ;
	} // for i

	return maxRays;
}



int 
CSurfaceOrenNayar::GetNumParameters() const
{
	return 2; // diffuse, slope
}
//========================================================================

CValues
CSurfaceOrenNayar::GetParameter(int parameterIndex) const
{
	CValues c(0,0,0);
	switch (parameterIndex)
	{
		case 0 :
			_diffuse->GetValue(c,CHitPointInfo::constMapInfo);
			break;
		case 1 :
			_slope->GetSlope(0)->GetValue(c,CHitPointInfo::constMapInfo);
			break;
	}
	return c;
}
//========================================================================

CValues
CSurfaceOrenNayar::GetParameterConstraints(int parameterIndex) const
{
	CValues c(CValues::EE_WeigthN);
	switch (parameterIndex)
	{
		case 0 :
			_setBasicConstraints(EE_ColorRGB,c);
			break;
		case 1 :
			c.SetRepresentation(CValues::EE_Weigth3);
			c.Set(0,0); // min slope - 0 radians
			c.Set(1,M_PI_2); // max slope - Pi/2 radians
			c.Set(2,0.5); // default
			break;
	}
	return c;
}
//========================================================================

const char* 
CSurfaceOrenNayar::GetParameterDescription(int parameterIndex) const
{
	switch (parameterIndex)
	{
	case 0 :
		return BrdfStrings::KDiffuseDesc;
	case 1 :
		return BrdfStrings::KSlopeONDesc;
	default :
		return BrdfStrings::KInvalidParameter;
	}
}
//========================================================================

CSurface::EParameterUIType
CSurfaceOrenNayar::GetParameterSyntax(int parameterIndex) const
{
	switch (parameterIndex)
	{
	case 0 : // diffuse 
		return EE_ColorRGB;
	case 1 : // slope[0]
		return EE_Float;
	default :
		return EE_Unknown;
	}
}
//========================================================================

const char*
CSurfaceOrenNayar::GetParameterTitle(int parameterIndex) const
{
	switch (parameterIndex)
	{
	case 0 :
		return BrdfStrings::KDiffuseTitle;
	case 1 :
		return BrdfStrings::KSlopeTitle;
	default :
		return BrdfStrings::KInvalidParameter;
	}
}
//========================================================================

void
CSurfaceOrenNayar::SetParameter(int parameterIndex, const CValues& value)
{

	switch (parameterIndex)
	{
	case 0 :
		{
			CConstantTexture* newTex = new CConstantTexture(value);
			SetKd(newTex);
			break;
		}
	case 1 :
		{
			float w = 1.f;
			CValues slope_value(value[0],value[0],value[0]);
			CConstantTexture* newTex = new CConstantTexture(slope_value);
			if (!_slope)
			{
				_slope = new CSlopeMicrofacets();
			}
			_slope->Set(newTex,0);
			_slope->SetWeight(w,0);
			break;
		}
	}
}
//========================================================================

const char* 
CSurfaceOrenNayar::GetHelp() const
{
	return BrdfStrings::KOrenNayarDesc;
}

__END_BRDFLIB_SOURCE