ccl_device_inline void barycentric const float3 tri_a const float3 tri

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
ccl_device_inline void barycentric(const float3 tri_a, const float3 tri_b, const float3 tri_c, const float3 p, float &u, float &v)
{
const float3 v0 = tri_a - tri_c, v1 = tri_b - tri_c, v2 = p - tri_c;
const float d00 = dot(v0, v0);
const float d01 = dot(v0, v1);
const float d11 = dot(v1, v1);
const float d20 = dot(v2, v0);
const float d21 = dot(v2, v1);
const float denom = d00 * d11 - d01 * d01;
u = (d11 * d20 - d01 * d21) / denom;
v = (d00 * d21 - d01 * d20) / denom;
}
ccl_device_inline float calculate_dot_and_distance(const float3 tri_a, const float3 tri_b, const float3 tri_c,
const float3 tri_n_a, const float3 tri_n_b, const float3 tri_n_c,
const float3 point, const float w) {
const float3 T_a = tri_a + w * tri_n_a;
const float3 T_b = tri_b + w * tri_n_b;
const float3 T_c = tri_c + w * tri_n_c;
const float3 n = normalize(cross(T_b - T_a, T_c - T_a));
return dot(n, point - T_a);
}
ccl_device_inline float3 calculate_uwv(const float3 tri_a, const float3 tri_b, const float3 tri_c,
const float3 tri_n_a, const float3 tri_n_b, const float3 tri_n_c,
const float3 point, const float w) {
const float3 T_a = tri_a + w * tri_n_a;
const float3 T_b = tri_b + w * tri_n_b;
const float3 T_c = tri_c + w * tri_n_c;
const float3 n = normalize(cross(T_b - T_a, T_c - T_a));
float3 v0 = point - T_a;
float3 proj = point - dot(n, v0) * n;
float u, v;
barycentric(T_a, T_b, T_c, proj, u, v);
return make_float3(u, v , w);
}
ccl_device_inline float3 prism_point_uvw_coords(const float3 tri_a, const float3 tri_b, const float3 tri_c,
const float3 tri_n_a, const float3 tri_n_b, const float3 tri_n_c,
const float3 point, const float EPS) {
float w = 0.0f;
float min = 0.0f;
float max = 1.0f;
float f_min_dot = calculate_dot_and_distance(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, min);
if ((std::abs(f_min_dot) < EPS)) {
return calculate_uwv(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, min);
}
float f_max_dot = calculate_dot_and_distance(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, max);
if ((std::abs(f_max_dot) < EPS)) {
return calculate_uwv(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, max);
}
// usually 12 iterations are enough
float f_mid_dot = 0.0f;
for (int i = 0; i < 40; ++i) {
w = 0.5f * min + 0.5f * max;
f_mid_dot = calculate_dot_and_distance(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, w);
if (std::abs(f_mid_dot) < EPS) {
break;
}
if ((f_min_dot < 0.0f) == (f_mid_dot < 0.0f)) {
f_min_dot = f_mid_dot;
min = w;
} else {
max = w;
}
}
if (f_mid_dot > 0.01f) {
printf("prism_point_uvw_coords: %f %f %f dot %f\n", (double)point.x, (double)point.y, (double)point.z, (double)f_mid_dot);
}
return calculate_uwv(tri_a, tri_b, tri_c, tri_n_a, tri_n_b, tri_n_c, point, w);
}