1 | // Copyright 2009-2021 Intel Corporation |
2 | // SPDX-License-Identifier: Apache-2.0 |
3 | |
4 | #pragma once |
5 | |
6 | #include "../common/ray.h" |
7 | #include "curve_intersector_precalculations.h" |
8 | |
9 | namespace embree |
10 | { |
11 | namespace isa |
12 | { |
13 | namespace __coneline_internal |
14 | { |
15 | template<int M, typename Epilog, typename ray_tfar_func> |
16 | static __forceinline bool intersectCone(const vbool<M>& valid_i, |
17 | const Vec3vf<M>& ray_org_in, const Vec3vf<M>& ray_dir, |
18 | const vfloat<M>& ray_tnear, const ray_tfar_func& ray_tfar, |
19 | const Vec4vf<M>& v0, const Vec4vf<M>& v1, |
20 | const vbool<M>& cL, const vbool<M>& cR, |
21 | const Epilog& epilog) |
22 | { |
23 | vbool<M> valid = valid_i; |
24 | |
25 | /* move ray origin closer to make calculations numerically stable */ |
26 | const vfloat<M> dOdO = sqr(ray_dir); |
27 | const vfloat<M> rcp_dOdO = rcp(dOdO); |
28 | const Vec3vf<M> center = vfloat<M>(0.5f)*(v0.xyz()+v1.xyz()); |
29 | const vfloat<M> dt = dot(center-ray_org_in,ray_dir)*rcp_dOdO; |
30 | const Vec3vf<M> ray_org = ray_org_in + dt*ray_dir; |
31 | |
32 | const Vec3vf<M> dP = v1.xyz() - v0.xyz(); |
33 | const Vec3vf<M> p0 = ray_org - v0.xyz(); |
34 | const Vec3vf<M> p1 = ray_org - v1.xyz(); |
35 | |
36 | const vfloat<M> dPdP = sqr(dP); |
37 | const vfloat<M> dP0 = dot(p0,dP); |
38 | const vfloat<M> dP1 = dot(p1,dP); |
39 | const vfloat<M> dOdP = dot(ray_dir,dP); |
40 | |
41 | // intersect cone body |
42 | const vfloat<M> dr = v0.w - v1.w; |
43 | const vfloat<M> hy = dPdP + sqr(dr); |
44 | const vfloat<M> dO0 = dot(ray_dir,p0); |
45 | const vfloat<M> OO = sqr(p0); |
46 | const vfloat<M> dPdP2 = sqr(dPdP); |
47 | const vfloat<M> dPdPr0 = dPdP*v0.w; |
48 | |
49 | const vfloat<M> A = dPdP2 - sqr(dOdP)*hy; |
50 | const vfloat<M> B = dPdP2*dO0 - dP0*dOdP*hy + dPdPr0*(dr*dOdP); |
51 | const vfloat<M> C = dPdP2*OO - sqr(dP0)*hy + dPdPr0*(2.0f*dr*dP0 - dPdPr0); |
52 | |
53 | const vfloat<M> D = B*B - A*C; |
54 | valid &= D >= 0.0f; |
55 | if (unlikely(none(valid))) { |
56 | return false; |
57 | } |
58 | |
59 | /* standard case for "non-parallel" rays */ |
60 | const vfloat<M> Q = sqrt(D); |
61 | const vfloat<M> rcp_A = rcp(A); |
62 | /* special case for rays that are "parallel" to the cone - assume miss */ |
63 | const vbool<M> isParallel = abs(A) <= min_rcp_input; |
64 | |
65 | vfloat<M> t_cone_lower = select (isParallel, neg_inf, (-B-Q)*rcp_A); |
66 | vfloat<M> t_cone_upper = select (isParallel, pos_inf, (-B+Q)*rcp_A); |
67 | const vfloat<M> y_lower = dP0 + t_cone_lower*dOdP; |
68 | const vfloat<M> y_upper = dP0 + t_cone_upper*dOdP; |
69 | t_cone_lower = select(valid & y_lower > 0.0f & y_lower < dPdP, t_cone_lower, pos_inf); |
70 | t_cone_upper = select(valid & y_upper > 0.0f & y_upper < dPdP, t_cone_upper, neg_inf); |
71 | |
72 | const vbool<M> hitDisk0 = valid & cL; |
73 | const vbool<M> hitDisk1 = valid & cR; |
74 | const vfloat<M> rcp_dOdP = rcp(dOdP); |
75 | const vfloat<M> t_disk0 = select (hitDisk0, select (sqr(p0*dOdP-ray_dir*dP0)<(sqr(v0.w)*sqr(dOdP)), -dP0*rcp_dOdP, pos_inf), pos_inf); |
76 | const vfloat<M> t_disk1 = select (hitDisk1, select (sqr(p1*dOdP-ray_dir*dP1)<(sqr(v1.w)*sqr(dOdP)), -dP1*rcp_dOdP, pos_inf), pos_inf); |
77 | const vfloat<M> t_disk_lower = min(t_disk0, t_disk1); |
78 | const vfloat<M> t_disk_upper = max(t_disk0, t_disk1); |
79 | |
80 | const vfloat<M> t_lower = min(t_cone_lower, t_disk_lower); |
81 | const vfloat<M> t_upper = max(t_cone_upper, select(t_lower==t_disk_lower, |
82 | select(t_disk_upper==vfloat<M>(pos_inf),neg_inf,t_disk_upper), |
83 | select(t_disk_lower==vfloat<M>(pos_inf),neg_inf,t_disk_lower))); |
84 | |
85 | const vbool<M> valid_lower = valid & ray_tnear <= dt+t_lower & dt+t_lower <= ray_tfar() & t_lower != vfloat<M>(pos_inf); |
86 | const vbool<M> valid_upper = valid & ray_tnear <= dt+t_upper & dt+t_upper <= ray_tfar() & t_upper != vfloat<M>(neg_inf); |
87 | |
88 | const vbool<M> valid_first = valid_lower | valid_upper; |
89 | if (unlikely(none(valid_first))) |
90 | return false; |
91 | |
92 | const vfloat<M> t_first = select(valid_lower, t_lower, t_upper); |
93 | const vfloat<M> y_first = select(valid_lower, y_lower, y_upper); |
94 | |
95 | const vfloat<M> rcp_dPdP = rcp(dPdP); |
96 | const Vec3vf<M> dP2drr0dP = dPdP*dr*v0.w*dP; |
97 | const Vec3vf<M> dPhy = dP*hy; |
98 | const vbool<M> cone_hit_first = valid & (t_first == t_cone_lower | t_first == t_cone_upper); |
99 | const vbool<M> disk0_hit_first = valid & (t_first == t_disk0); |
100 | const Vec3vf<M> Ng_first = select(cone_hit_first, dPdP2*(p0+t_first*ray_dir)+dP2drr0dP-dPhy*y_first, select(disk0_hit_first, -dP, dP)); |
101 | const vfloat<M> u_first = select(cone_hit_first, y_first*rcp_dPdP, select(disk0_hit_first, vfloat<M>(zero), vfloat<M>(one))); |
102 | |
103 | /* invoke intersection filter for first hit */ |
104 | RoundLineIntersectorHitM<M> hit(u_first,zero,dt+t_first,Ng_first); |
105 | const bool is_hit_first = epilog(valid_first, hit); |
106 | |
107 | /* check for possible second hits before potentially accepted hit */ |
108 | const vfloat<M> t_second = t_upper; |
109 | const vfloat<M> y_second = y_upper; |
110 | const vbool<M> valid_second = valid_lower & valid_upper & (dt+t_upper <= ray_tfar()); |
111 | if (unlikely(none(valid_second))) |
112 | return is_hit_first; |
113 | |
114 | /* invoke intersection filter for second hit */ |
115 | const vbool<M> cone_hit_second = t_second == t_cone_lower | t_second == t_cone_upper; |
116 | const vbool<M> disk0_hit_second = t_second == t_disk0; |
117 | const Vec3vf<M> Ng_second = select(cone_hit_second, dPdP2*(p0+t_second*ray_dir)+dP2drr0dP-dPhy*y_second, select(disk0_hit_second, -dP, dP)); |
118 | const vfloat<M> u_second = select(cone_hit_second, y_second*rcp_dPdP, select(disk0_hit_first, vfloat<M>(zero), vfloat<M>(one))); |
119 | |
120 | hit = RoundLineIntersectorHitM<M>(u_second,zero,dt+t_second,Ng_second); |
121 | const bool is_hit_second = epilog(valid_second, hit); |
122 | |
123 | return is_hit_first | is_hit_second; |
124 | } |
125 | } |
126 | |
127 | template<int M> |
128 | struct ConeLineIntersectorHitM |
129 | { |
130 | __forceinline ConeLineIntersectorHitM() {} |
131 | |
132 | __forceinline ConeLineIntersectorHitM(const vfloat<M>& u, const vfloat<M>& v, const vfloat<M>& t, const Vec3vf<M>& Ng) |
133 | : vu(u), vv(v), vt(t), vNg(Ng) {} |
134 | |
135 | __forceinline void finalize() {} |
136 | |
137 | __forceinline Vec2f uv (const size_t i) const { return Vec2f(vu[i],vv[i]); } |
138 | __forceinline float t (const size_t i) const { return vt[i]; } |
139 | __forceinline Vec3fa Ng(const size_t i) const { return Vec3fa(vNg.x[i],vNg.y[i],vNg.z[i]); } |
140 | |
141 | public: |
142 | vfloat<M> vu; |
143 | vfloat<M> vv; |
144 | vfloat<M> vt; |
145 | Vec3vf<M> vNg; |
146 | }; |
147 | |
148 | template<int M> |
149 | struct ConeCurveIntersector1 |
150 | { |
151 | typedef CurvePrecalculations1 Precalculations; |
152 | |
153 | struct ray_tfar { |
154 | Ray& ray; |
155 | __forceinline ray_tfar(Ray& ray) : ray(ray) {} |
156 | __forceinline vfloat<M> operator() () const { return ray.tfar; }; |
157 | }; |
158 | |
159 | template<typename Epilog> |
160 | static __forceinline bool intersect(const vbool<M>& valid_i, |
161 | Ray& ray, |
162 | IntersectContext* context, |
163 | const LineSegments* geom, |
164 | const Precalculations& pre, |
165 | const Vec4vf<M>& v0i, const Vec4vf<M>& v1i, |
166 | const vbool<M>& cL, const vbool<M>& cR, |
167 | const Epilog& epilog) |
168 | { |
169 | const Vec3vf<M> ray_org(ray.org.x, ray.org.y, ray.org.z); |
170 | const Vec3vf<M> ray_dir(ray.dir.x, ray.dir.y, ray.dir.z); |
171 | const vfloat<M> ray_tnear(ray.tnear()); |
172 | const Vec4vf<M> v0 = enlargeRadiusToMinWidth<M>(context,geom,ray_org,v0i); |
173 | const Vec4vf<M> v1 = enlargeRadiusToMinWidth<M>(context,geom,ray_org,v1i); |
174 | return __coneline_internal::intersectCone<M>(valid_i,ray_org,ray_dir,ray_tnear,ray_tfar(ray),v0,v1,cL,cR,epilog); |
175 | } |
176 | }; |
177 | |
178 | template<int M, int K> |
179 | struct ConeCurveIntersectorK |
180 | { |
181 | typedef CurvePrecalculationsK<K> Precalculations; |
182 | |
183 | struct ray_tfar { |
184 | RayK<K>& ray; |
185 | size_t k; |
186 | __forceinline ray_tfar(RayK<K>& ray, size_t k) : ray(ray), k(k) {} |
187 | __forceinline vfloat<M> operator() () const { return ray.tfar[k]; }; |
188 | }; |
189 | |
190 | template<typename Epilog> |
191 | static __forceinline bool intersect(const vbool<M>& valid_i, |
192 | RayK<K>& ray, size_t k, |
193 | IntersectContext* context, |
194 | const LineSegments* geom, |
195 | const Precalculations& pre, |
196 | const Vec4vf<M>& v0i, const Vec4vf<M>& v1i, |
197 | const vbool<M>& cL, const vbool<M>& cR, |
198 | const Epilog& epilog) |
199 | { |
200 | const Vec3vf<M> ray_org(ray.org.x[k], ray.org.y[k], ray.org.z[k]); |
201 | const Vec3vf<M> ray_dir(ray.dir.x[k], ray.dir.y[k], ray.dir.z[k]); |
202 | const vfloat<M> ray_tnear = ray.tnear()[k]; |
203 | const Vec4vf<M> v0 = enlargeRadiusToMinWidth<M>(context,geom,ray_org,v0i); |
204 | const Vec4vf<M> v1 = enlargeRadiusToMinWidth<M>(context,geom,ray_org,v1i); |
205 | return __coneline_internal::intersectCone<M>(valid_i,ray_org,ray_dir,ray_tnear,ray_tfar(ray,k),v0,v1,cL,cR,epilog); |
206 | } |
207 | }; |
208 | } |
209 | } |
210 | |