Reflectance at normal incidence from the complex refractive index¶
A surface's reflectance at normal incidence from vacuum, F0 = ((n - 1)^2 + k^2) / ((n + 1)^2 + k^2), from its phase's complex refractive index n + ik at one wavelength (lusine.materials.
refractive_index, at 450, 550, 589.3 and 650 nm): a metal's colour, its F0 at the three, and a dielectric's with k = 0, ((n - 1) / (n + 1))^2. lusine.properties.fresnel at cos_theta = 1. Source: Born & Wolf, Principles of Optics, 7th ed., 14.2 (eq. 14.2.17); Hecht, Optics, 5th ed., 4.8. Checked against a reference in numpy, which has its own tests.
- Source
- Born & Wolf, Principles of Optics, 7th ed., 14.2 (eq. 14.2.17); Hecht, Optics, 5th ed., 4.8.
- Proof
- reference, limit.
fresnel_f0.trait
namespace lusine.properties;
float fresnel_f0(float n, float k)
{
let less = n - 1.0;
let more = n + 1.0;
let k2 = k * k;
return (less * less + k2) / (more * more + k2);
}
A surface's reflectance at normal incidence from vacuum, F0 = ((n − 1)² + k²)/((n + 1)² + k²), from its phase's complex index n + ik at one wavelength (lusine.materials.refractive_index, at 450, 550, 589.3 and 650 nm); a dielectric is k = 0, ((n − 1)/(n + 1))². A metal's colour is its F0 at the three wavelengths a renderer's channels sample, and the library stores none: this kernel derives it.
import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
import matplotlib.pyplot as plt
from lusine import trait_run
from lusine.library_sources import parse_rows
EPS = float(np.finfo(np.float32).eps) # Trait's float is binary32
The kernel¶
The law as the library states it.
print(trait_run.source('properties/fresnel_f0.trait'))
// lusine.properties.fresnel_f0: a surface's reflectance at normal incidence from vacuum, F0 = ((n - 1)^2 + k^2) /
// ((n + 1)^2 + k^2), from its phase's complex refractive index n + ik at one wavelength (lusine.materials.
// refractive_index, at 450, 550, 589.3 and 650 nm): a metal's colour, its F0 at the three, and a dielectric's with
// k = 0, ((n - 1) / (n + 1))^2. lusine.properties.fresnel at cos_theta = 1.
// Source: Born & Wolf, Principles of Optics, 7th ed., 14.2 (eq. 14.2.17); Hecht, Optics, 5th ed., 4.8.
// Checked against a reference in numpy, which has its own tests.
namespace lusine.properties;
float fresnel_f0(float n, float k)
{
let less = n - 1.0;
let more = n + 1.0;
let k2 = k * k;
return (less * less + k2) / (more * more + k2);
}
The program¶
Each section runs the kernel over a lattice of items, one a case.
print(trait_run.law_program('properties/fresnel_f0.trait', 'fresnel_f0'))
namespace probe;
lusine.version version = { { number = 1 } };
struct item {
host float n;
host float k;
float result_0;
};
lattice<item, 1> items(spacing = 1 m, boundary = wall);
void main() { items.result_0 = lusine.properties.fresnel_f0(items.n, items.k); }
Reference: every index in the library, against |(N − 1)/(N + 1)|²¶
Each row of refractive_index, as committed (k absent is a transparent phase, k = 0), is held to |(N − 1)/(N + 1)|² of the complex index N = n + ik, computed in complex arithmetic, to binary32's rounding. Soda-lime glass at the sodium D line reflects 4.30%, water 2.03%.
rows = parse_rows(trait_run.source('materials/generated/optical.trait'), 'refractive_index')
n = np.array([r['n'] for r in rows])
k = np.array([r.get('k', 0.0) for r in rows])
f0 = trait_run.law('properties/fresnel_f0.trait', 'fresnel_f0', n=n, k=k)
N = np.float32(n).astype(float) + 1j * np.float32(k).astype(float)
want = np.abs((N - 1) / (N + 1)) ** 2
assert np.all(np.abs(f0 - want) <= 8 * EPS * want), np.abs(f0 / want - 1).max()
at = {(r['material'], r['state'], r['wavelength']): v for r, v in zip(rows, f0)}
assert abs(at[('soda_lime_glass', 'solid', 589.3)] - 0.0430) < 1e-4 and abs(at[('water', 'liquid', 589.3)] - 0.0203) < 1e-4
for (m, s, w), v in at.items():
if w == 589.3:
print(f'{m:15s} {s:7s} F0 at 589.3 nm: {v:.4f}')
boxwood solid F0 at 589.3 nm: 0.0439 copper solid F0 at 589.3 nm: 0.7649 soda_lime_glass solid F0 at 589.3 nm: 0.0430 grass solid F0 at 589.3 nm: 0.0290 hematite solid F0 at 589.3 nm: 0.2894 iron solid F0 at 589.3 nm: 0.5209 mercury liquid F0 at 589.3 nm: 0.7803 silicon solid F0 at 589.3 nm: 0.3566 silicon_dioxide solid F0 at 589.3 nm: 0.0458 tungsten solid F0 at 589.3 nm: 0.5057 water solid F0 at 589.3 nm: 0.0180 water liquid F0 at 589.3 nm: 0.0203
A metal's colour, against measured F0¶
Copper's and iron's F0 at 650, 550 and 450 nm (red, green, blue) from Johnson & Christy's n and k, against Hoffman's measured F0 in linear RGB ("Physics and Math of Shading", SIGGRAPH 2015 course notes). Copper's agree within 2.3%; iron's within 15%, the derived ones darker throughout (−14.2% at blue): Hoffman's F0 is another measurement's, integrated over each sRGB channel's response, where these sample one wavelength a channel, and Johnson & Christy's iron is an evaporated film.
# Hoffman's measured F0, linear RGB
METALS = {'copper': (0.955, 0.638, 0.538), 'iron': (0.562, 0.565, 0.578)}
index = {(r['material'], r['wavelength']): (r['n'], r['k']) for r in rows if r['state'] == 'solid' and 'k' in r}
for metal, bound in (('copper', 0.03), ('iron', 0.15)):
nk = np.array([index[(metal, w)] for w in (650.0, 550.0, 450.0)])
rgb = trait_run.law('properties/fresnel_f0.trait', 'fresnel_f0', n=nk[:, 0], k=nk[:, 1])
hoffman = np.array(METALS[metal])
change = rgb / hoffman - 1
print(f'{metal:7s} derived {np.round(rgb, 3)} METALS {hoffman} largest change {change[np.argmax(np.abs(change))]:+.2%}')
assert np.all(np.abs(change) < bound), change
fig, ax = plt.subplots(figsize=(7, 3.2))
for i, metal in enumerate(('copper', 'iron')):
nk = np.array([index[(metal, w)] for w in (650.0, 550.0, 450.0)])
rgb = trait_run.law('properties/fresnel_f0.trait', 'fresnel_f0', n=nk[:, 0], k=nk[:, 1])
ax.bar(np.arange(3) + 4 * i - 0.2, rgb, 0.4, color=['#c33', '#3a3', '#33c'], label=None if i else 'derived')
ax.bar(np.arange(3) + 4 * i + 0.2, METALS[metal], 0.4, color='grey', label=None if i else 'METALS (Hoffman)')
ax.set_xticks([1, 5], ['copper', 'iron'])
ax.set_ylabel('F0 at 650, 550, 450 nm')
ax.legend()
plt.show()
copper derived [0.935 0.624 0.538] METALS [0.955 0.638 0.538] largest change -2.27% iron derived [0.531 0.512 0.496] METALS [0.562 0.565 0.578] largest change -14.22%
Limit: a dielectric, vacuum and a perfect conductor¶
At k = 0 it is ((n − 1)/(n + 1))², and glass of n = 1.5 reflects exactly 4% (both exact in binary32 here); n = 1 with k = 0 reflects nothing; as k grows without bound everything is reflected.
f = trait_run.law('properties/fresnel_f0.trait', 'fresnel_f0', n=[1.5, 1.0, 1.0, 2.0], k=[0.0, 0.0, 1e6, 0.0])
assert abs(f[0] - 0.04) <= EPS * 0.04 and f[1] == 0.0 and f[2] > 1 - 1e-6 and abs(f[3] - 1 / 9) <= EPS / 9, f
print(f)
[0.04 0. 1. 0.11111111]