#ifndef __BRDF_INC
#define __BRDF_INC

#include "pbr.cginc"
#include "lighting.cginc"
#include "lysenko.cginc"
#include "math.cginc"

// Schlick "An Inexpensive BRDF Model for Physically-based Rendering".
// Equation 24.
// f0: Reflectance at normal incidence. Typically around 0.04.
// f90: Reflectance at grazing incidence. Typically around 1.0.
float F_Schlick(float LoH, float f0, float f90) {
  float term = 1.0f - LoH;
  float term2 = term * term;
  float term5 = term2 * term2 * term;
  return f0 + (f90 - f0) * term5;
}

// Walter "Microfacet Models for Refraction through Rough Surfaces"
// Equation 33.
// In the paper:
//  - m = microsurface normal
//  - n = macrosurface normal
//  - theta_m = angle between micro- & macrosurface normals
//  - alpha = roughness
//  - cos(theta_m) = NoH
// Per sohcahtoa:
//  tan(theta) = sin(theta) / cos(theta)
//  tan^2(theta) = sin^2(theta) / cos^2(theta)
//               = (1 - cos^2(theta)) / cos^2(theta)
//               = -1 + 1 / cos^2(theta)
float D_GGX(float roughness, float NoH) {
  float r2 = roughness * roughness;
  float NoH2 = NoH * NoH;
  float NoH4 = NoH2 * NoH2;

  float k = rcp(NoH2) - 1;
  float r2_plus_k = r2 + k;
  // Not sure why, but not using the factor of PI here makes the specular match
  // the Unity standard much more closely. Maybe the author was just folding
  // the 4.0 (historically used to be a PI) into the GGX calculation?
  float denom = NoH4 * r2_plus_k * r2_plus_k;

  return r2 / denom;
}

// Hammon "PBR Diffuse Lighting for GGX+Smith Microsurfaces"
// Slide 84. Note that we remove the (4 * NoL * NoV) from the
// denominator of the specular lobe because of some cancellations.
// The original, un-optimized equation is:
//  2 * NoL * NoV / lerp(2 * NoL * NoV, NoL + NoV, roughness)
float G_GGXSmith(float roughness, float NoL, float NoV) {
  float denom = 2.0f * lerp(2.0f * NoL * NoV, NoL + NoV, roughness);
  return rcp(denom);
}

// Estevez "Production Friendly Microfacet Sheen BRDF"
// Equation 2.
// The original equation is:
//  (2 + 1/r) * sin^(1-r)(theta) / (2 pi)
// Recall that:
//  cos^2(theta) + sin^2(theta) = 1
// So:
//  sin^2(theta) = 1 - cos^2(theta)
//  sin(theta) = (1 - cos^2(theta)) ^ (1/2)
//  sin^k(theta) = (1 - cos^2(theta)) ^ (k/2)
float D_Cloth(float roughness, float NoH) {
  float r_rcp = rcp(roughness);
  return (2.0f + r_rcp) * pow(1.0f - NoH * NoH, r_rcp * 0.5f) / TAU;
}

float G_Cloth_L(float x, float a, float b, float c, float d, float e) {
  return a / (1.0f + b * pow(x, c)) + d * x + e;
}

// Estevez "Production Friendly Microfacet Sheen BRDF"
// Equations 3 and 4.
float G_Cloth(float roughness, float LoH) {
  // Table 1
  float a0 = 25.3245f;
  float a1 = 21.5473f;
  float b0 = 3.32435f;
  float b1 = 3.82987f;
  float c0 = 0.16801f;
  float c1 = 0.19823f;
  float d0 = -1.27393f;
  float d1 = -1.97760f;
  float e0 = -4.85967f;
  float e1 = -4.32054f;
  float one_minus_r = 1.0f - roughness;
  float one_minus_r_sq = one_minus_r * one_minus_r;
  float one_minus_LoH = 1.0f - LoH;

  float lambda;
  [branch]
  if (LoH < 0.5f) {
    float L0 = G_Cloth_L(LoH, a0, b0, c0, d0, e0);
    float L1 = G_Cloth_L(LoH, a1, b1, c1, d1, e1);
    float L = lerp(L0, L1, one_minus_r_sq);
    lambda = exp(L);
  } else {
    float L_05_0 = G_Cloth_L(0.5f, a0, b0, c0, d0, e0);
    float L_05_1 = G_Cloth_L(0.5f, a1, b1, c1, d1, e1);
    float L_05 = lerp(L_05_0, L_05_1, one_minus_r_sq);

    float L_LoH_0 = G_Cloth_L(one_minus_LoH, a0, b0, c0, d0, e0);
    float L_LoH_1 = G_Cloth_L(one_minus_LoH, a1, b1, c1, d1, e1);
    float L_LoH = lerp(L_LoH_0, L_LoH_1, one_minus_r_sq);

    lambda = exp(2.0f * L_05 - L_LoH);
  }
  // Apply terminator softening (equation 4).
  return pow(lambda, 1.0f + 2.0f * pow(one_minus_LoH, 8));
}

float4 brdf(Pbr pbr, LightData data) {
  float3 specular = 0;
  float3 diffuse = 0;

  float f0 = 0.04f;
  const float f90 = 1.0f;

//#define FURNACE_TEST_DIRECT
#if defined(FURNACE_TEST_DIRECT)
  // Create the conditions for the standard BRDF furnace test.
  // Only applies to the direct lighting stage. The only variable left over is
  // NoV.
  f0 = 1;
  data.direct.color = 1;
  data.direct.NoL = 1;
  data.direct.NoH = 1;
  data.direct.LoH = 1;
#endif

  // Direct
  if (true) {
    float remainder = 1.0f;

#if defined(_CLEARCOAT)
    float cc_f0 = 0.04f;
    float Fcc = F_Schlick(data.direct.LoH, cc_f0, f90);
    float Dcc = D_GGX(pbr.cc_roughness, data.direct.NoH_cc);
    float Gcc = G_GGXSmith(pbr.cc_roughness, data.direct.NoL_cc, data.common.NoV_cc);
    float DFGcc = Fcc * Dcc * Gcc;
    float3 direct_specular_cc = DFGcc * data.direct.color * data.direct.NoL_cc * pbr.cc_strength;
    direct_specular_cc = max(0, direct_specular_cc);
    specular += direct_specular_cc;
    remainder -= Fcc * pbr.cc_strength;
#endif

#if defined(_CLOTH_SHEEN)
    float cl_f0 = 0.04f;
    float Fcl = 1;
    float Dcl = D_Cloth(pbr.roughness, data.direct.NoH);
    float Gcl = G_Cloth(pbr.roughness, data.direct.LoH);
    float DFGcl = Fcl * Dcl * Gcl;
    float3 direct_specular_cl = DFGcl * data.direct.color * pbr.cl_strength * pbr.cl_color * data.direct.NoL;
    direct_specular_cl = max(0, direct_specular_cl);
    specular += direct_specular_cl;
    remainder -= direct_specular_cl;
#endif

    float F = F_Schlick(data.direct.LoH, f0, f90);
    float D = D_GGX(pbr.roughness, data.direct.NoH);
    float G = G_GGXSmith(pbr.roughness, data.direct.NoL, data.common.NoV);

    float FDG = F * D * G;

    float3 direct_specular = FDG * remainder * data.direct.color * data.direct.NoL * lerp(1.0f, pbr.albedo.xyz, pbr.metallic);
    direct_specular = max(0, direct_specular);
    specular += direct_specular;
    remainder -= F;

    float Fd = Fd_OrenNayar(pbr.roughness, data.common.NoV, data.direct.NoL, data.direct.LoV) / PI;
    float3 direct_diffuse = Fd * remainder * (1.0f - pbr.metallic) * pbr.albedo.xyz * data.direct.color;
    direct_diffuse = max(0, direct_diffuse);
    diffuse += direct_diffuse;
  }

  // Indirect
#if defined(FORWARD_BASE_PASS)
  if (true) {
    float remainder = 1.0f;
    float2 dfg_uv = float2(data.common.NoV, pbr.roughness);

#if defined(_CLEARCOAT)
    float cc_f0 = 0.04f;
    float Fcc = F_Schlick(data.common.NoV, cc_f0, 1.0f);
    float3 indirect_specular_cc = Fcc * data.indirect.specular_cc * pbr.cc_strength;
    specular += indirect_specular_cc;
    remainder *= (1.0f - Fcc * pbr.cc_strength);
#endif

#if defined(_CLOTH_SHEEN)
    float DFGcl = _Cloth_Sheen_DFG_LUT.Sample(bilinear_clamp_s, dfg_uv).r;
    float3 indirect_specular_cl = DFGcl * data.indirect.specular * pbr.cl_strength * pbr.cl_color;
    specular += indirect_specular_cl * remainder;
    // Energy conservation for cloth is tricky with IBL.
    // A simple approximation is to use the Fresnel of the sheen layer.
    float Fcl = F_Schlick(data.common.NoV, 0.04, 1.0);
    remainder *= (1.0f - Fcl * pbr.cl_strength);
#endif

    // Standard PBR IBL using split-sum approximation
    float2 dfg = _DFG_LUT.Sample(bilinear_clamp_s, dfg_uv).rg;
    float3 f0_spec = lerp(0.04f, pbr.albedo.xyz, pbr.metallic);
    float3 ibl_specular_reflectance = f0_spec * dfg.x + dfg.y;
    float3 indirect_specular = data.indirect.specular * ibl_specular_reflectance;
    specular += indirect_specular * remainder;

    // For energy conservation with the diffuse term, we use the view-dependent Fresnel.
    float3 F = F_Schlick(data.common.NoV, f0_spec, 1.0f);
    remainder *= (1.0f - F);

    // Diffuse is Lambertian, which is pre-integrated into the SH diffuse probe
    float3 indirect_diffuse = pbr.albedo.xyz * data.indirect.diffuse * remainder * (1.0 - pbr.metallic);
    diffuse  += indirect_diffuse;
  }
#endif

  return float4(diffuse + specular, 1);
}

#endif  // __BRDF_INC