CoCalc -- bspline

GitHub Repository: godotengine/godot
Path: blob/master/thirdparty/embree/kernels/subdiv/bspline_patch.h
⁹⁹¹³ views
1
// Copyright 2009-2021 Intel Corporation
2
// SPDX-License-Identifier: Apache-2.0
3

4
#pragma once
5

6
#include "catmullclark_patch.h"
7
#include "bspline_curve.h"
8

9
namespace embree
10
{
11
  template<typename Vertex, typename Vertex_t = Vertex>
12
    class __aligned(64) BSplinePatchT
13
    {
14
      typedef CatmullClark1RingT<Vertex,Vertex_t> CatmullClarkRing;
15
      typedef CatmullClarkPatchT<Vertex,Vertex_t> CatmullClarkPatch;
16
      
17
    public:
18
      
19
      __forceinline BSplinePatchT () {}
20

21
      __forceinline BSplinePatchT (const CatmullClarkPatch& patch) {
22
        init(patch);
23
      }
24

25
      __forceinline BSplinePatchT(const CatmullClarkPatch& patch,
26
                                  const BezierCurveT<Vertex>* border0,
27
                                  const BezierCurveT<Vertex>* border1,
28
                                  const BezierCurveT<Vertex>* border2,
29
                                  const BezierCurveT<Vertex>* border3)
30
      {
31
        init(patch);
32
      }
33

34
      __forceinline BSplinePatchT (const HalfEdge* edge, const char* vertices, size_t stride) {
35
        init(edge,vertices,stride);
36
      }
37

38
      __forceinline Vertex hard_corner(const                    Vertex& v01, const Vertex& v02, 
39
                                       const Vertex& v10, const Vertex& v11, const Vertex& v12, 
40
                                       const Vertex& v20, const Vertex& v21, const Vertex& v22)
41
      {
42
        return 4.0f*v11 - 2.0f*(v12+v21) + v22;
43
      }
44

45
      __forceinline Vertex soft_convex_corner( const                    Vertex& v01, const Vertex& v02, 
46
                                               const Vertex& v10, const Vertex& v11, const Vertex& v12, 
47
                                               const Vertex& v20, const Vertex& v21, const Vertex& v22)
48
      {
49
        return -8.0f*v11 + 4.0f*(v12+v21) + v22;
50
      }
51

52
      __forceinline Vertex convex_corner(const float vertex_crease_weight, 
53
                                         const                    Vertex& v01, const Vertex& v02, 
54
                                         const Vertex& v10, const Vertex& v11, const Vertex& v12, 
55
                                         const Vertex& v20, const Vertex& v21, const Vertex& v22)
56
      {
57
        if (std::isinf(vertex_crease_weight)) return hard_corner(v01,v02,v10,v11,v12,v20,v21,v22);
58
        else                                  return soft_convex_corner(v01,v02,v10,v11,v12,v20,v21,v22);
59
      }
60

61
      __forceinline Vertex load(const HalfEdge* edge, const char* vertices, size_t stride) {
62
        return Vertex_t::loadu(vertices+edge->getStartVertexIndex()*stride);
63
      }
64

65
      __forceinline void init_border(const CatmullClarkRing& edge0,
66
                                     Vertex& v01, Vertex& v02,
67
                                     const Vertex& v11, const Vertex& v12,
68
                                     const Vertex& v21, const Vertex& v22)
69
      {
70
        if (likely(edge0.has_opposite_back(0)))
71
        {
72
          v01 = edge0.back(2);
73
          v02 = edge0.back(1);
74
        } else {
75
          v01 = 2.0f*v11-v21;
76
          v02 = 2.0f*v12-v22;
77
        }
78
      }
79

80
      __forceinline void init_corner(const CatmullClarkRing& edge0,
81
                                     Vertex& v00,       const Vertex& v01, const Vertex& v02, 
82
                                     const Vertex& v10, const Vertex& v11, const Vertex& v12, 
83
                                     const Vertex& v20, const Vertex& v21, const Vertex& v22)
84
      {
85
        const bool MAYBE_UNUSED has_back1 = edge0.has_opposite_back(1);
86
        const bool has_back0 = edge0.has_opposite_back(0);
87
        const bool has_front1 = edge0.has_opposite_front(1);
88
        const bool MAYBE_UNUSED has_front2 = edge0.has_opposite_front(2);
89
        
90
        if (likely(has_back0)) {
91
          if (likely(has_front1)) { assert(has_back1 && has_front2); v00 = edge0.back(3); }
92
          else { assert(!has_back1); v00 = 2.0f*v01-v02; }
93
        }
94
        else {
95
          if (likely(has_front1)) { assert(!has_front2); v00 = 2.0f*v10-v20; }
96
          else v00 = convex_corner(edge0.vertex_crease_weight,v01,v02,v10,v11,v12,v20,v21,v22);
97
        }
98
      }
99
      
100
      void init(const CatmullClarkPatch& patch)
101
      {
102
        /* fill inner vertices */
103
        const Vertex v11 = v[1][1] = patch.ring[0].vtx;
104
        const Vertex v12 = v[1][2] = patch.ring[1].vtx;
105
        const Vertex v22 = v[2][2] = patch.ring[2].vtx; 
106
        const Vertex v21 = v[2][1] = patch.ring[3].vtx; 
107
        
108
        /* fill border vertices */
109
        init_border(patch.ring[0],v[0][1],v[0][2],v11,v12,v21,v22);
110
        init_border(patch.ring[1],v[1][3],v[2][3],v12,v22,v11,v21);
111
        init_border(patch.ring[2],v[3][2],v[3][1],v22,v21,v12,v11);
112
        init_border(patch.ring[3],v[2][0],v[1][0],v21,v11,v22,v12);
113
        
114
        /* fill corner vertices */
115
        init_corner(patch.ring[0],v[0][0],v[0][1],v[0][2],v[1][0],v11,v12,v[2][0],v21,v22);
116
        init_corner(patch.ring[1],v[0][3],v[1][3],v[2][3],v[0][2],v12,v22,v[0][1],v11,v21);
117
        init_corner(patch.ring[2],v[3][3],v[3][2],v[3][1],v[2][3],v22,v21,v[1][3],v12,v11);
118
        init_corner(patch.ring[3],v[3][0],v[2][0],v[1][0],v[3][1],v21,v11,v[3][2],v22,v12);
119
      }
120
      
121
      void init_border(const HalfEdge* edge0, const char* vertices, size_t stride,
122
                                     Vertex& v01, Vertex& v02,
123
                                     const Vertex& v11, const Vertex& v12,
124
                                     const Vertex& v21, const Vertex& v22)
125
      {
126
        if (likely(edge0->hasOpposite())) 
127
        {
128
          const HalfEdge* e = edge0->opposite()->next()->next(); 
129
          v01 = load(e,vertices,stride); 
130
          v02 = load(e->next(),vertices,stride);
131
        } else {
132
          v01 = 2.0f*v11-v21;
133
          v02 = 2.0f*v12-v22;
134
        }
135
      }
136
      
137
      void init_corner(const HalfEdge* edge0, const char* vertices, size_t stride,
138
                       Vertex& v00, const Vertex& v01, const Vertex& v02, 
139
                       const Vertex& v10, const Vertex& v11, const Vertex& v12, 
140
                       const Vertex& v20, const Vertex& v21, const Vertex& v22)
141
      {
142
        const bool has_back0 = edge0->hasOpposite();
143
        const bool has_front1 = edge0->prev()->hasOpposite();
144

145
        if (likely(has_back0))
146
        { 
147
          const HalfEdge* e = edge0->opposite()->next();
148
          if (likely(has_front1))
149
          {
150
            assert(e->hasOpposite());
151
            assert(edge0->prev()->opposite()->prev()->hasOpposite());
152
            v00 = load(e->opposite()->prev(),vertices,stride);
153
          } 
154
          else {
155
            assert(!e->hasOpposite());
156
            v00 = 2.0f*v01-v02;
157
          }
158
        }
159
        else
160
        {
161
          if (likely(has_front1)) {
162
            assert(!edge0->prev()->opposite()->prev()->hasOpposite());
163
            v00 = 2.0f*v10-v20;
164
          }
165
          else {
166
            assert(edge0->vertex_crease_weight == 0.0f || std::isinf(edge0->vertex_crease_weight));
167
            v00 = convex_corner(edge0->vertex_crease_weight,v01,v02,v10,v11,v12,v20,v21,v22);
168
          }
169
        }
170
      }
171
      
172
      void init(const HalfEdge* edge0, const char* vertices, size_t stride)
173
      {
174
        assert( edge0->isRegularFace() );
175
        
176
        /* fill inner vertices */
177
        const Vertex v11 = v[1][1] = load(edge0,vertices,stride); const HalfEdge* edge1 = edge0->next();
178
        const Vertex v12 = v[1][2] = load(edge1,vertices,stride); const HalfEdge* edge2 = edge1->next();
179
        const Vertex v22 = v[2][2] = load(edge2,vertices,stride); const HalfEdge* edge3 = edge2->next();
180
        const Vertex v21 = v[2][1] = load(edge3,vertices,stride); assert(edge0  == edge3->next());
181
        
182
        /* fill border vertices */
183
        init_border(edge0,vertices,stride,v[0][1],v[0][2],v11,v12,v21,v22);
184
        init_border(edge1,vertices,stride,v[1][3],v[2][3],v12,v22,v11,v21);
185
        init_border(edge2,vertices,stride,v[3][2],v[3][1],v22,v21,v12,v11);
186
        init_border(edge3,vertices,stride,v[2][0],v[1][0],v21,v11,v22,v12);
187
        
188
        /* fill corner vertices */
189
        init_corner(edge0,vertices,stride,v[0][0],v[0][1],v[0][2],v[1][0],v11,v12,v[2][0],v21,v22);
190
        init_corner(edge1,vertices,stride,v[0][3],v[1][3],v[2][3],v[0][2],v12,v22,v[0][1],v11,v21);
191
        init_corner(edge2,vertices,stride,v[3][3],v[3][2],v[3][1],v[2][3],v22,v21,v[1][3],v12,v11);
192
        init_corner(edge3,vertices,stride,v[3][0],v[2][0],v[1][0],v[3][1],v21,v11,v[3][2],v22,v12);
193
      }
194
      
195
      __forceinline BBox<Vertex> bounds() const
196
      {
197
        const Vertex* const cv = &v[0][0];
198
        BBox<Vertex> bounds (cv[0]);
199
        for (size_t i=1; i<16 ; i++)
200
          bounds.extend( cv[i] );
201
        return bounds;
202
      }
203
      
204
      __forceinline Vertex eval(const float uu, const float vv) const
205
      {
206
        const Vec4f v_n = BSplineBasis::eval(vv);
207
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
208
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
209
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
210
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
211
        
212
        const Vec4f u_n = BSplineBasis::eval(uu);
213
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
214
      }
215
      
216
      __forceinline Vertex eval_du(const float uu, const float vv) const
217
      {
218
        const Vec4f v_n = BSplineBasis::eval(vv);
219
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
220
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
221
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
222
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
223
        
224
        const Vec4f u_n = BSplineBasis::derivative(uu);
225
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
226
      }
227
      
228
      __forceinline Vertex eval_dv(const float uu, const float vv) const
229
      {
230
        const Vec4f v_n = BSplineBasis::derivative(vv);
231
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
232
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
233
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
234
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
235
        
236
        const Vec4f u_n = BSplineBasis::eval(uu);
237
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
238
      }
239

240
      __forceinline Vertex eval_dudu(const float uu, const float vv) const
241
      {
242
        const Vec4f v_n = BSplineBasis::eval(vv);
243
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
244
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
245
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
246
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
247
        
248
        const Vec4f u_n = BSplineBasis::derivative2(uu);
249
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
250
      }
251

252
      __forceinline Vertex eval_dvdv(const float uu, const float vv) const
253
      {
254
        const Vec4f v_n = BSplineBasis::derivative2(vv);
255
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
256
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
257
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
258
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
259
        
260
        const Vec4f u_n = BSplineBasis::eval(uu);
261
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
262
      }
263

264
      __forceinline Vertex eval_dudv(const float uu, const float vv) const
265
      {
266
        const Vec4f v_n = BSplineBasis::derivative(vv);
267
        const Vertex_t curve0 = madd(v_n[0],v[0][0],madd(v_n[1],v[1][0],madd(v_n[2],v[2][0],v_n[3] * v[3][0])));
268
        const Vertex_t curve1 = madd(v_n[0],v[0][1],madd(v_n[1],v[1][1],madd(v_n[2],v[2][1],v_n[3] * v[3][1])));
269
        const Vertex_t curve2 = madd(v_n[0],v[0][2],madd(v_n[1],v[1][2],madd(v_n[2],v[2][2],v_n[3] * v[3][2])));
270
        const Vertex_t curve3 = madd(v_n[0],v[0][3],madd(v_n[1],v[1][3],madd(v_n[2],v[2][3],v_n[3] * v[3][3])));
271
        
272
        const Vec4f u_n = BSplineBasis::derivative(uu);
273
        return madd(u_n[0],curve0,madd(u_n[1],curve1,madd(u_n[2],curve2,u_n[3] * curve3)));
274
      }
275
      
276
      __forceinline Vertex normal(const float uu, const float vv) const
277
      {
278
        const Vertex tu = eval_du(uu,vv);
279
        const Vertex tv = eval_dv(uu,vv);
280
        return cross(tu,tv);
281
      }   
282

283
      template<typename T>
284
      __forceinline Vec3<T> eval(const T& uu, const T& vv, const Vec4<T>& u_n, const Vec4<T>& v_n) const
285
      {
286
        const T curve0_x = madd(v_n[0],T(v[0][0].x),madd(v_n[1],T(v[1][0].x),madd(v_n[2],T(v[2][0].x),v_n[3] * T(v[3][0].x))));
287
        const T curve1_x = madd(v_n[0],T(v[0][1].x),madd(v_n[1],T(v[1][1].x),madd(v_n[2],T(v[2][1].x),v_n[3] * T(v[3][1].x))));
288
        const T curve2_x = madd(v_n[0],T(v[0][2].x),madd(v_n[1],T(v[1][2].x),madd(v_n[2],T(v[2][2].x),v_n[3] * T(v[3][2].x))));
289
        const T curve3_x = madd(v_n[0],T(v[0][3].x),madd(v_n[1],T(v[1][3].x),madd(v_n[2],T(v[2][3].x),v_n[3] * T(v[3][3].x))));
290
        const T x = madd(u_n[0],curve0_x,madd(u_n[1],curve1_x,madd(u_n[2],curve2_x,u_n[3] * curve3_x)));
291
                  
292
        const T curve0_y = madd(v_n[0],T(v[0][0].y),madd(v_n[1],T(v[1][0].y),madd(v_n[2],T(v[2][0].y),v_n[3] * T(v[3][0].y))));
293
        const T curve1_y = madd(v_n[0],T(v[0][1].y),madd(v_n[1],T(v[1][1].y),madd(v_n[2],T(v[2][1].y),v_n[3] * T(v[3][1].y))));
294
        const T curve2_y = madd(v_n[0],T(v[0][2].y),madd(v_n[1],T(v[1][2].y),madd(v_n[2],T(v[2][2].y),v_n[3] * T(v[3][2].y))));
295
        const T curve3_y = madd(v_n[0],T(v[0][3].y),madd(v_n[1],T(v[1][3].y),madd(v_n[2],T(v[2][3].y),v_n[3] * T(v[3][3].y))));
296
        const T y = madd(u_n[0],curve0_y,madd(u_n[1],curve1_y,madd(u_n[2],curve2_y,u_n[3] * curve3_y)));
297
          
298
        const T curve0_z = madd(v_n[0],T(v[0][0].z),madd(v_n[1],T(v[1][0].z),madd(v_n[2],T(v[2][0].z),v_n[3] * T(v[3][0].z))));
299
        const T curve1_z = madd(v_n[0],T(v[0][1].z),madd(v_n[1],T(v[1][1].z),madd(v_n[2],T(v[2][1].z),v_n[3] * T(v[3][1].z))));
300
        const T curve2_z = madd(v_n[0],T(v[0][2].z),madd(v_n[1],T(v[1][2].z),madd(v_n[2],T(v[2][2].z),v_n[3] * T(v[3][2].z))));
301
        const T curve3_z = madd(v_n[0],T(v[0][3].z),madd(v_n[1],T(v[1][3].z),madd(v_n[2],T(v[2][3].z),v_n[3] * T(v[3][3].z))));
302
        const T z = madd(u_n[0],curve0_z,madd(u_n[1],curve1_z,madd(u_n[2],curve2_z,u_n[3] * curve3_z)));
303
        
304
        return Vec3<T>(x,y,z);
305
      }
306
      
307
      template<typename T>
308
      __forceinline Vec3<T> eval(const T& uu, const T& vv) const
309
      {
310
        const Vec4<T> u_n = BSplineBasis::eval(uu);
311
        const Vec4<T> v_n = BSplineBasis::eval(vv);
312
        return eval(uu,vv,u_n,v_n);
313
      }
314

315
      template<typename T>
316
      __forceinline Vec3<T> eval_du(const T& uu, const T& vv) const
317
      {
318
        const Vec4<T> u_n = BSplineBasis::derivative(uu); 
319
        const Vec4<T> v_n = BSplineBasis::eval(vv); 
320
        return eval(uu,vv,u_n,v_n);      
321
      }
322
      
323
      template<typename T>
324
      __forceinline Vec3<T> eval_dv(const T& uu, const T& vv) const
325
      {
326
        const Vec4<T> u_n = BSplineBasis::eval(uu); 
327
        const Vec4<T> v_n = BSplineBasis::derivative(vv); 
328
        return eval(uu,vv,u_n,v_n);      
329
      }
330

331
      template<typename T>
332
      __forceinline Vec3<T> eval_dudu(const T& uu, const T& vv) const
333
      {
334
        const Vec4<T> u_n = BSplineBasis::derivative2(uu); 
335
        const Vec4<T> v_n = BSplineBasis::eval(vv); 
336
        return eval(uu,vv,u_n,v_n);      
337
      }
338

339
      template<typename T>
340
      __forceinline Vec3<T> eval_dvdv(const T& uu, const T& vv) const
341
      {
342
        const Vec4<T> u_n = BSplineBasis::eval(uu); 
343
        const Vec4<T> v_n = BSplineBasis::derivative2(vv); 
344
        return eval(uu,vv,u_n,v_n);      
345
      }
346

347
      template<typename T>
348
      __forceinline Vec3<T> eval_dudv(const T& uu, const T& vv) const
349
      {
350
        const Vec4<T> u_n = BSplineBasis::derivative(uu); 
351
        const Vec4<T> v_n = BSplineBasis::derivative(vv); 
352
        return eval(uu,vv,u_n,v_n);      
353
      }
354
      
355
      template<typename T>
356
      __forceinline Vec3<T> normal(const T& uu, const T& vv) const {
357
        return cross(eval_du(uu,vv),eval_dv(uu,vv));
358
      }
359

360
      void eval(const float u, const float v, 
361
                Vertex* P, Vertex* dPdu, Vertex* dPdv, Vertex* ddPdudu, Vertex* ddPdvdv, Vertex* ddPdudv, 
362
                const float dscale = 1.0f) const
363
      {
364
        if (P) {
365
          *P = eval(u,v); 
366
        }
367
        if (dPdu) {
368
          assert(dPdu); *dPdu = eval_du(u,v)*dscale; 
369
          assert(dPdv); *dPdv = eval_dv(u,v)*dscale; 
370
        }
371
        if (ddPdudu) {
372
          assert(ddPdudu); *ddPdudu = eval_dudu(u,v)*sqr(dscale); 
373
          assert(ddPdvdv); *ddPdvdv = eval_dvdv(u,v)*sqr(dscale); 
374
          assert(ddPdudv); *ddPdudv = eval_dudv(u,v)*sqr(dscale); 
375
        }
376
      }
377

378
      template<class vfloat>
379
      __forceinline vfloat eval(const size_t i, const vfloat& uu, const vfloat& vv, const Vec4<vfloat>& u_n, const Vec4<vfloat>& v_n) const
380
      {
381
        const vfloat curve0_x = madd(v_n[0],vfloat(v[0][0][i]),madd(v_n[1],vfloat(v[1][0][i]),madd(v_n[2],vfloat(v[2][0][i]),v_n[3] * vfloat(v[3][0][i]))));
382
        const vfloat curve1_x = madd(v_n[0],vfloat(v[0][1][i]),madd(v_n[1],vfloat(v[1][1][i]),madd(v_n[2],vfloat(v[2][1][i]),v_n[3] * vfloat(v[3][1][i]))));
383
        const vfloat curve2_x = madd(v_n[0],vfloat(v[0][2][i]),madd(v_n[1],vfloat(v[1][2][i]),madd(v_n[2],vfloat(v[2][2][i]),v_n[3] * vfloat(v[3][2][i]))));
384
        const vfloat curve3_x = madd(v_n[0],vfloat(v[0][3][i]),madd(v_n[1],vfloat(v[1][3][i]),madd(v_n[2],vfloat(v[2][3][i]),v_n[3] * vfloat(v[3][3][i]))));
385
        return madd(u_n[0],curve0_x,madd(u_n[1],curve1_x,madd(u_n[2],curve2_x,u_n[3] * curve3_x)));
386
      }
387
        
388
      template<typename vbool, typename vfloat>
389
      void eval(const vbool& valid, const vfloat& uu, const vfloat& vv, 
390
                float* P, float* dPdu, float* dPdv, float* ddPdudu, float* ddPdvdv, float* ddPdudv, 
391
                const float dscale, const size_t dstride, const size_t N) const
392
      {
393
        if (P) {
394
          const Vec4<vfloat> u_n = BSplineBasis::eval(uu); 
395
          const Vec4<vfloat> v_n = BSplineBasis::eval(vv); 
396
          for (size_t i=0; i<N; i++) vfloat::store(valid,P+i*dstride,eval(i,uu,vv,u_n,v_n));
397
        }
398
        if (dPdu) 
399
        {
400
          {
401
            assert(dPdu);
402
            const Vec4<vfloat> u_n = BSplineBasis::derivative(uu); 
403
            const Vec4<vfloat> v_n = BSplineBasis::eval(vv);
404
            for (size_t i=0; i<N; i++) vfloat::store(valid,dPdu+i*dstride,eval(i,uu,vv,u_n,v_n)*dscale);
405
          }
406
          {
407
            assert(dPdv);
408
            const Vec4<vfloat> u_n = BSplineBasis::eval(uu); 
409
            const Vec4<vfloat> v_n = BSplineBasis::derivative(vv);
410
            for (size_t i=0; i<N; i++) vfloat::store(valid,dPdv+i*dstride,eval(i,uu,vv,u_n,v_n)*dscale);
411
          }
412
        }
413
        if (ddPdudu) 
414
        {
415
          {
416
            assert(ddPdudu);
417
            const Vec4<vfloat> u_n = BSplineBasis::derivative2(uu); 
418
            const Vec4<vfloat> v_n = BSplineBasis::eval(vv);
419
            for (size_t i=0; i<N; i++) vfloat::store(valid,ddPdudu+i*dstride,eval(i,uu,vv,u_n,v_n)*sqr(dscale));
420
          }
421
          {
422
            assert(ddPdvdv);
423
            const Vec4<vfloat> u_n = BSplineBasis::eval(uu); 
424
            const Vec4<vfloat> v_n = BSplineBasis::derivative2(vv);
425
            for (size_t i=0; i<N; i++) vfloat::store(valid,ddPdvdv+i*dstride,eval(i,uu,vv,u_n,v_n)*sqr(dscale));
426
          }
427
          {
428
            assert(ddPdudv);
429
            const Vec4<vfloat> u_n = BSplineBasis::derivative(uu); 
430
            const Vec4<vfloat> v_n = BSplineBasis::derivative(vv);
431
            for (size_t i=0; i<N; i++) vfloat::store(valid,ddPdudv+i*dstride,eval(i,uu,vv,u_n,v_n)*sqr(dscale));
432
          }
433
        }
434
      }
435

436
      friend __forceinline embree_ostream operator<<(embree_ostream o, const BSplinePatchT& p)
437
      {
438
        for (size_t y=0; y<4; y++)
439
          for (size_t x=0; x<4; x++)
440
            o << "[" << y << "][" << x << "] " << p.v[y][x] << embree_endl;
441
        return o;
442
      } 
443

444
    public:
445
      Vertex v[4][4];
446
    };
447
  
448
  typedef BSplinePatchT<Vec3fa,Vec3fa_t> BSplinePatch3fa;
449
}
450

451
Product

Resources

Company
