<p>The polar coordinate system commonly used in BRDF definitions is set up relative to the surface being shaded, i.e. in tangent space. The <span class="math-container">$\theta$</span> angle measures how far you are from the surface normal while <span class="math-container">$\phi$</span> is the azimuth around the plane of the surface relative to some reference direction (which doesn't matter unless the BRDF is anisotropic). So if you want to convert to these coordinates, you have to first get your vector in tangent space (with the origin at the intersection point and two of the axes aligned with the surface), then apply the usual Cartesian-to-spherical transformation.</p>
<p>However, normally you can and should evaluate BRDFs <em>without</em> using polar coordinates, trigonometric functions, or angles at all, but just using vector math primitives such as dot products. This is usually more efficient <em>and</em> it's more robust, as you don't have to deal with angle wraparound, factors of pi, out-of-range arguments to inverse trig functions and so on. For instance, you probably know that the cosine of the angle between vectors can be obtained by just dotting the (normalized) vectors. The sine and tangent can be obtained through trig identities from the cosine (i.e. from the dot product).</p>
<p>Fabian Giesen wrote an article on this topic, <a href="https://fgiesen.wordpress.com/2010/10/21/finish-your-derivations-please/" rel="noreferrer">Finish your derivations, please</a>, that refers to the exact Oren-Nayar article you linked and gives an alternate, trig-free form:</p>
<blockquote>
<p>Time to brush up your trigonometric identities. A particularly bad offender can be found here – the relevant section from the simplified shader is this:</p>
<pre><code>float alpha = max( acos( dot( v, n ) ), acos( dot( l, n ) ) );
float beta  = min( acos( dot( v, n ) ), acos( dot( l, n ) ) );
C = sin(alpha) * tan(beta);
</code></pre>
<p>Ouch! If you use some trig identities and the fact that acos is monotonically decreasing over its domain, this reduces to:</p>
<pre><code>float vdotn = dot(v, n);
float ldotn = dot(l, n);
C = sqrt((1.0 - vdotn*vdotn) * (1.0 - ldotn*ldotn))
  / max(vdotn, ldotn);
</code></pre>
<p>..and suddenly there’s no need to use a lookup texture anymore (and by the way, this has way higher accuracy too).</p>
</blockquote>
<p>And note that this formulation only has one transcendental function (<code>sqrt</code>) in place of four (two <code>acos</code>, a <code>sin</code>, and a <code>tan</code>).</p>

CC BY-SA 3.0
2015-12-11T12:19:49.197Z
Thank you so much for your help @NathanReed. One last question: do you have any good reference material on how to calculate the tangent space and convert my vectors wi and wo to the coordinate space for ray tracing/BRDF? At the moment I didn't find any useful one. In this way I would be able to do some comparison between the way with tangent space and the way using vectors math.
Fabrizio Duroni
CC BY-SA 3.0
2015-12-11T19:11:04.707Z
@FabrizioDuroni I assume you're familiar with how to convert between coordinate systems in general? For tangent space you just have to set up coordinates using the surface normal plus some two vectors perpendicular to it as the axes. For normal mapping, the two vectors are often chosen to match the texture space U and V axes (as mapped to the particular surface). For isotropic BRDFs without normal mapping, it doesn't really matter.
Nathan Reed
CC BY-SA 3.0
2015-12-11T06:18:50.267Z
2015-12-11T06:27:07.310Z
2020-06-18T08:32:12.800Z
Nathan Reed
<p>I developed a ray tracer that use standard phong/blinn phong lighting model. Now I'm modifying it to support physically based rendering, so I'm implementing various BRDF models. At the moment I'm focused on Oren-Nayar and Torrance-Sparrow model. Each one of these is based on spherical coordinates used to express incident wi and outgoing wo light direction. </p>

<p><strong>My question is: which way is the right one convert wi and wo from cartesian coordinate to spherical coordinate?</strong> </p>

<p>I'm applying the standard formula reported here <a href="https://en.wikipedia.org/wiki/Spherical_coordinate_system#Coordinate_system_conversions" rel="nofollow">https://en.wikipedia.org/wiki/Spherical_coordinate_system#Coordinate_system_conversions</a> but I'm not sure I'm doing the right thing, <strong>because my vector are not with tail at the origin of the cartesian coordinate system, but are centered on the intersection point of the ray with the object.</strong> </p>

<p>Here you can find my current implementation:</p>

<ul>
<li><p><a href="https://github.com/chicio/Multispectral-Ray-tracing/tree/brdf/RayTracing/RayTracer/Objects/BRDF" rel="nofollow">https://github.com/chicio/Multispectral-Ray-tracing/tree/brdf/RayTracing/RayTracer/Objects/BRDF</a></p></li>
<li><p><a href="https://github.com/chicio/Multispectral-Ray-tracing/blob/brdf/RayTracing/RayTracer/Math/Vector3D.cpp" rel="nofollow">https://github.com/chicio/Multispectral-Ray-tracing/blob/brdf/RayTracing/RayTracer/Math/Vector3D.cpp</a></p></li>
</ul>

<p>Can anyone help me giving an explanation of the correct way to convert the wi and wo vector from cartesian to spherical coordinate?</p>

<p><strong>UPDATE</strong></p>

<p>I copy here the relevant part of code:</p>

<p>spherical coordinate calculation</p>

<pre><code>float Vector3D::sphericalTheta() const {

    float sphericalTheta = acosf(Utils::clamp(y, -1.f, 1.f));

    return sphericalTheta;
}

float Vector3D::sphericalPhi() const {

    float phi = atan2f(z, x);

    return (phi &lt; 0.f) ? phi + 2.f * M_PI : phi;
}
</code></pre>

<p>Oren Nayar</p>

<pre><code>OrenNayar::OrenNayar(Spectrum&lt;constant::spectrumSamples&gt; reflectanceSpectrum, float degree) : reflectanceSpectrum{reflectanceSpectrum} {

    float sigma = Utils::degreeToRadian(degree);
    float sigmaPowerTwo = sigma * sigma;

    A = 1.0f - (sigmaPowerTwo / 2.0f * (sigmaPowerTwo + 0.33f));
    B = 0.45f * sigmaPowerTwo / (sigmaPowerTwo + 0.09f);
};

Spectrum&lt;constant::spectrumSamples&gt; OrenNayar::f(const Vector3D&amp; wi, const Vector3D&amp; wo, const Intersection* intersection) const {

    float thetaI = wi.sphericalTheta();
    float phiI = wi.sphericalPhi();

    float thetaO = wo.sphericalTheta();
    float phiO = wo.sphericalPhi();

    float alpha = std::fmaxf(thetaI, thetaO);
    float beta = std::fminf(thetaI, thetaO);

    Spectrum&lt;constant::spectrumSamples&gt; orenNayar = reflectanceSpectrum * constant::inversePi * (A + B * std::fmaxf(0, cosf(phiI - phiO) * sinf(alpha) * tanf(beta)));

    return orenNayar;
}
</code></pre>

<p>Torrance-Sparrow</p>

<pre><code>float TorranceSparrow::G(const Vector3D&amp; wi, const Vector3D&amp; wo, const Vector3D&amp; wh, const Intersection* intersection) const {

    Vector3D normal = intersection-&gt;normal;
    normal.normalize();

    float normalDotWh = fabsf(normal.dot(wh));
    float normalDotWo = fabsf(normal.dot(wo));
    float normalDotWi = fabsf(normal.dot(wi));
    float woDotWh = fabsf(wo.dot(wh));

    float G = fminf(1.0f, std::fminf((2.0f * normalDotWh * normalDotWo)/woDotWh, (2.0f * normalDotWh * normalDotWi)/woDotWh));

    return G;
}

float TorranceSparrow::D(const Vector3D&amp; wh, const Intersection* intersection) const {

    Vector3D normal = intersection-&gt;normal;
    normal.normalize();

    float cosThetaH = fabsf(wh.dot(normal));

    float Dd = (exponent + 2) * constant::inverseTwoPi * powf(cosThetaH, exponent);

    return Dd;
}

Spectrum&lt;constant::spectrumSamples&gt; TorranceSparrow::f(const Vector3D&amp; wi, const Vector3D&amp; wo, const Intersection* intersection) const {

    Vector3D normal = intersection-&gt;normal;
    normal.normalize();

    float thetaI = wi.sphericalTheta();
    float thetaO = wo.sphericalTheta();

    float cosThetaO = fabsf(cosf(thetaO));
    float cosThetaI = fabsf(cosf(thetaI));

    if(cosThetaI == 0 || cosThetaO == 0) {

        return reflectanceSpectrum * 0.0f;
    }

    Vector3D wh = (wi + wo);
    wh.normalize();

    float cosThetaH = wi.dot(wh);

    float F = Fresnel::dieletricFresnel(cosThetaH, refractiveIndex);
    float g = G(wi, wo, wh, intersection);
    float d = D(wh, intersection);

    printf("f %f g %f d %f \n", F, g, d);
    printf("result %f \n", ((d * g * F) / (4.0f * cosThetaI * cosThetaO)));

    Spectrum&lt;constant::spectrumSamples&gt; torranceSparrow = reflectanceSpectrum * ((d * g * F) / (4.0f * cosThetaI * cosThetaO));

    return torranceSparrow;
}
</code></pre>

<p><strong>UPDATE 2</strong></p>

<p>After some search I found this implementation of Oren-Nayar BRDF.</p>

<p><a href="http://content.gpwiki.org/index.php/D3DBook:(Lighting)_Oren-Nayar" rel="nofollow">http://content.gpwiki.org/index.php/D3DBook:(Lighting)_Oren-Nayar</a></p>

<p>In the implementation above theta for wi and wo is obtained simply doing arccos(wo.dotProduct(Normal)) and arccos(wi.dotProduct(Normal)). This seems reasonable to me, as we can use the normal of the intersection point as the zenith direction for our spherical coordinate system and do the calculation. The calculation of gamma = cos(phi_wi - phi_wo) do some sort of projection of wi and wo on what it calls "tangent space". Assuming everything is correct in this implementation, can i just use the formulas |View - Normal x (View.dotProduct(Normal))| and |Light - Normal x (Light.dotProduct(Normal))| to obtain the phi coordinate (instead of using arctan("something"))? </p>

<hr>

<p>N.B.: I'm also reading "Physically based rendering: from theory to implementation", that is the book released to pbrt (<a href="http://www.pbrt.org" rel="nofollow">http://www.pbrt.org</a>). In this implementation there's some sort of change of coordinate system for the point of intersection (using partial derivates and the parametric coordinate of the surface) (I'm reading it now so what I'm saying could not be accurate). I want to find a straight way (maybe the one above in UPDATE 2 is what I'm searching for, if anyone could confirm it).</p>

CC BY-SA 3.0
2015-12-10T13:10:30.087Z
I already asked this question here http://gamedev.stackexchange.com/questions/112165/brdf-and-spherical-coordinate-in-ray-tracing. Nobody seems to have an answer. I think that the arguments in my "UPDATE 2" could be the way to follow (because it seems that the project of vector in this way described could be the correct way to calculate the azimuth angle). Anyone could help me with an answer and maybe a canonical reference to be used as study reference?
Fabrizio Duroni
<raytracing><brdf>
BRDF and Spherical coordinate in ray tracing
Fabrizio Duroni
# BRDF and Spherical coordinate in ray tracing

(Asked by: username_0 on 2015-12-10)

I developed a ray tracer that use standard phong/blinn phong lighting model. Now I'm modifying it to support physically based rendering, so I'm implementing various BRDF models. At the moment I'm focused on Oren-Nayar and Torrance-Sparrow model. Each one of these is based on spherical coordinates used to express incident wi and outgoing wo light direction.

**My question is: which way is the right one convert wi and wo from cartesian coordinate to spherical coordinate?**

I'm applying the standard formula reported here <https://en.wikipedia.org/wiki/Spherical_coordinate_system#Coordinate_system_conversions> but I'm not sure I'm doing the right thing, **because my vector are not with tail at the origin of the cartesian coordinate system, but are centered on the intersection point of the ray with the object.**

Here you can find my current implementation:

* <https://github.com/chicio/Multispectral-Ray-tracing/tree/brdf/RayTracing/RayTracer/Objects/BRDF>
* <https://github.com/chicio/Multispectral-Ray-tracing/blob/brdf/RayTracing/RayTracer/Math/Vector3D.cpp>

Can anyone help me giving an explanation of the correct way to convert the wi and wo vector from cartesian to spherical coordinate?

**UPDATE**

I copy here the relevant part of code:

spherical coordinate calculation

```
float Vector3D::sphericalTheta() const {

    float sphericalTheta = acosf(Utils::clamp(y, -1.f, 1.f));

    return sphericalTheta;
}

float Vector3D::sphericalPhi() const {

    float phi = atan2f(z, x);

    return (phi < 0.f) ? phi + 2.f * M_PI : phi;
}

```

Oren Nayar

```
OrenNayar::OrenNayar(Spectrum<constant::spectrumSamples> reflectanceSpectrum, float degree) : reflectanceSpectrum{reflectanceSpectrum} {

    float sigma = Utils::degreeToRadian(degree);
    float sigmaPowerTwo = sigma * sigma;

    A = 1.0f - (sigmaPowerTwo / 2.0f * (sigmaPowerTwo + 0.33f));
    B = 0.45f * sigmaPowerTwo / (sigmaPowerTwo + 0.09f);
};

Spectrum<constant::spectrumSamples> OrenNayar::f(const Vector3D& wi, const Vector3D& wo, const Intersection* intersection) const {

    float thetaI = wi.sphericalTheta();
    float phiI = wi.sphericalPhi();

    float thetaO = wo.sphericalTheta();
    float phiO = wo.sphericalPhi();

    float alpha = std::fmaxf(thetaI, thetaO);
    float beta = std::fminf(thetaI, thetaO);

    Spectrum<constant::spectrumSamples> orenNayar = reflectanceSpectrum * constant::inversePi * (A + B * std::fmaxf(0, cosf(phiI - phiO) * sinf(alpha) * tanf(beta)));

    return orenNayar;
}

```

Torrance-Sparrow

```
float TorranceSparrow::G(const Vector3D& wi, const Vector3D& wo, const Vector3D& wh, const Intersection* intersection) const {

    Vector3D normal = intersection->normal;
    normal.normalize();

    float normalDotWh = fabsf(normal.dot(wh));
    float normalDotWo = fabsf(normal.dot(wo));
    float normalDotWi = fabsf(normal.dot(wi));
    float woDotWh = fabsf(wo.dot(wh));

    float G = fminf(1.0f, std::fminf((2.0f * normalDotWh * normalDotWo)/woDotWh, (2.0f * normalDotWh * normalDotWi)/woDotWh));

    return G;
}

float TorranceSparrow::D(const Vector3D& wh, const Intersection* intersection) const {

    Vector3D normal = intersection->normal;
    normal.normalize();

    float cosThetaH = fabsf(wh.dot(normal));

    float Dd = (exponent + 2) * constant::inverseTwoPi * powf(cosThetaH, exponent);

    return Dd;
}

Spectrum<constant::spectrumSamples> TorranceSparrow::f(const Vector3D& wi, const Vector3D& wo, const Intersection* intersection) const {

    Vector3D normal = intersection->normal;
    normal.normalize();

    float thetaI = wi.sphericalTheta();
    float thetaO = wo.sphericalTheta();

    float cosThetaO = fabsf(cosf(thetaO));
    float cosThetaI = fabsf(cosf(thetaI));

    if(cosThetaI == 0 || cosThetaO == 0) {

        return reflectanceSpectrum * 0.0f;
    }

    Vector3D wh = (wi + wo);
    wh.normalize();

    float cosThetaH = wi.dot(wh);

    float F = Fresnel::dieletricFresnel(cosThetaH, refractiveIndex);
    float g = G(wi, wo, wh, intersection);
    float d = D(wh, intersection);

    printf("f %f g %f d %f \n", F, g, d);
    printf("result %f \n", ((d * g * F) / (4.0f * cosThetaI * cosThetaO)));

    Spectrum<constant::spectrumSamples> torranceSparrow = reflectanceSpectrum * ((d * g * F) / (4.0f * cosThetaI * cosThetaO));

    return torranceSparrow;
}

```

**UPDATE 2**

After some search I found this implementation of Oren-Nayar BRDF.

<http://content.gpwiki.org/index.php/D3DBook:(Lighting)_Oren-Nayar>

In the implementation above theta for wi and wo is obtained simply doing arccos(wo.dotProduct(Normal)) and arccos(wi.dotProduct(Normal)). This seems reasonable to me, as we can use the normal of the intersection point as the zenith direction for our spherical coordinate system and do the calculation. The calculation of gamma = cos(phi\_wi - phi\_wo) do some sort of projection of wi and wo on what it calls "tangent space". Assuming everything is correct in this implementation, can i just use the formulas |View - Normal x (View.dotProduct(Normal))| and |Light - Normal x (Light.dotProduct(Normal))| to obtain the phi coordinate (instead of using arctan("something"))?

---

N.B.: I'm also reading "Physically based rendering: from theory to implementation", that is the book released to pbrt (<http://www.pbrt.org>). In this implementation there's some sort of change of coordinate system for the point of intersection (using partial derivates and the parametric coordinate of the surface) (I'm reading it now so what I'm saying could not be accurate). I want to find a straight way (maybe the one above in UPDATE 2 is what I'm searching for, if anyone could confirm it).

**Tags:** raytracing, brdf

> **Comments:**
> **username_0**: I already asked this question here http://gamedev.stackexchange.com/questions/112165/brdf-and-spherical-coordinate-in-ray-tracing. Nobody seems to have an answer. I think that the arguments in my "UPDATE 2" could be the way to follow (because it seems that the project of vector in this way described could be the correct way to calculate the azimuth angle). Anyone could help me with an answer and maybe a canonical reference to be used as study reference?  

---

## 1 answers

### Answer by username_1 on 2015-12-11. Score: 9 | Accepted answer

The polar coordinate system commonly used in BRDF definitions is set up relative to the surface being shaded, i.e. in tangent space. The $\theta$ angle measures how far you are from the surface normal while $\phi$ is the azimuth around the plane of the surface relative to some reference direction (which doesn't matter unless the BRDF is anisotropic). So if you want to convert to these coordinates, you have to first get your vector in tangent space (with the origin at the intersection point and two of the axes aligned with the surface), then apply the usual Cartesian-to-spherical transformation.

However, normally you can and should evaluate BRDFs *without* using polar coordinates, trigonometric functions, or angles at all, but just using vector math primitives such as dot products. This is usually more efficient *and* it's more robust, as you don't have to deal with angle wraparound, factors of pi, out-of-range arguments to inverse trig functions and so on. For instance, you probably know that the cosine of the angle between vectors can be obtained by just dotting the (normalized) vectors. The sine and tangent can be obtained through trig identities from the cosine (i.e. from the dot product).

Fabian Giesen wrote an article on this topic, [Finish your derivations, please](https://fgiesen.wordpress.com/2010/10/21/finish-your-derivations-please/), that refers to the exact Oren-Nayar article you linked and gives an alternate, trig-free form:

> Time to brush up your trigonometric identities. A particularly bad offender can be found here – the relevant section from the simplified shader is this:
>
> ```
> float alpha = max( acos( dot( v, n ) ), acos( dot( l, n ) ) );
> float beta  = min( acos( dot( v, n ) ), acos( dot( l, n ) ) );
> C = sin(alpha) * tan(beta);
>
> ```
>
> Ouch! If you use some trig identities and the fact that acos is monotonically decreasing over its domain, this reduces to:
>
> ```
> float vdotn = dot(v, n);
> float ldotn = dot(l, n);
> C = sqrt((1.0 - vdotn*vdotn) * (1.0 - ldotn*ldotn))
>   / max(vdotn, ldotn);
>
> ```
>
> ..and suddenly there’s no need to use a lookup texture anymore (and by the way, this has way higher accuracy too).

And note that this formulation only has one transcendental function (`sqrt`) in place of four (two `acos`, a `sin`, and a `tan`).

> **Comments:**
> **username_0**: Thank you so much for your help @username_1. One last question: do you have any good reference material on how to calculate the tangent space and convert my vectors wi and wo to the coordinate space for ray tracing/BRDF? At the moment I didn't find any useful one. In this way I would be able to do some comparison between the way with tangent space and the way using vectors math.  
> **username_1**: @username_0 I assume you're familiar with how to convert between coordinate systems in general? For tangent space you just have to set up coordinates using the surface normal plus some two vectors perpendicular to it as the axes. For normal mapping, the two vectors are often chosen to match the texture space U and V axes (as mapped to the particular surface). For isotropic BRDFs without normal mapping, it doesn't really matter.  

---

