CoCalc -- bspline

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

4
#pragma once
5

6
#include "../common/default.h"
7
#include "bezier_curve.h"
8

9
namespace embree
10
{
11
  class BSplineBasis
12
  {
13
  public:
14

15
    template<typename T>
16
      static __forceinline Vec4<T> eval(const T& u) 
17
    {
18
      const T t  = u;
19
      const T s  = T(1.0f) - u;
20
      const T n0 = s*s*s;
21
      const T n1 = (4.0f*(s*s*s)+(t*t*t)) + (12.0f*((s*t)*s) + 6.0f*((t*s)*t));
22
      const T n2 = (4.0f*(t*t*t)+(s*s*s)) + (12.0f*((t*s)*t) + 6.0f*((s*t)*s));
23
      const T n3 = t*t*t;
24
      return T(1.0f/6.0f)*Vec4<T>(n0,n1,n2,n3);
25
    }
26
    
27
    template<typename T>
28
      static __forceinline Vec4<T>  derivative(const T& u)
29
    {
30
      const T t  =  u;
31
      const T s  =  1.0f - u;
32
      const T n0 = -s*s;
33
      const T n1 = -t*t - 4.0f*(t*s);
34
      const T n2 =  s*s + 4.0f*(s*t);
35
      const T n3 =  t*t;
36
      return T(0.5f)*Vec4<T>(n0,n1,n2,n3);
37
    }
38

39
    template<typename T>
40
      static __forceinline Vec4<T>  derivative2(const T& u)
41
    {
42
      const T t  =  u;
43
      const T s  =  1.0f - u;
44
      const T n0 = s;
45
      const T n1 = t - 2.0f*s;
46
      const T n2 = s - 2.0f*t;
47
      const T n3 = t;
48
      return Vec4<T>(n0,n1,n2,n3);
49
    }
50
  };
51
  
52
  struct PrecomputedBSplineBasis
53
  {
54
    enum { N = 16 };
55
  public:
56
    PrecomputedBSplineBasis() {}
57
    PrecomputedBSplineBasis(int shift);
58

59
    /* basis for bspline evaluation */
60
  public:
61
    float c0[N+1][N+1];
62
    float c1[N+1][N+1];
63
    float c2[N+1][N+1];
64
    float c3[N+1][N+1];
65
    
66
    /* basis for bspline derivative evaluation */
67
  public:
68
    float d0[N+1][N+1];
69
    float d1[N+1][N+1];
70
    float d2[N+1][N+1];
71
    float d3[N+1][N+1];
72
  };
73
  extern PrecomputedBSplineBasis bspline_basis0;
74
  extern PrecomputedBSplineBasis bspline_basis1;
75

76
  template<typename Vertex>
77
    struct BSplineCurveT
78
    {
79
      Vertex v0,v1,v2,v3;
80
      
81
      __forceinline BSplineCurveT() {}
82
      
83
      __forceinline BSplineCurveT(const Vertex& v0, const Vertex& v1, const Vertex& v2, const Vertex& v3)
84
        : v0(v0), v1(v1), v2(v2), v3(v3) {}
85

86
      __forceinline Vertex begin() const {
87
        return madd(1.0f/6.0f,v0,madd(2.0f/3.0f,v1,1.0f/6.0f*v2));
88
      }
89

90
      __forceinline Vertex end() const {
91
        return madd(1.0f/6.0f,v1,madd(2.0f/3.0f,v2,1.0f/6.0f*v3));
92
      }
93

94
      __forceinline Vertex center() const {
95
        return 0.25f*(v0+v1+v2+v3);
96
      }
97

98
      __forceinline BBox<Vertex> bounds() const {
99
        return merge(BBox<Vertex>(v0),BBox<Vertex>(v1),BBox<Vertex>(v2),BBox<Vertex>(v3));
100
      }
101
      
102
      __forceinline friend BSplineCurveT operator -( const BSplineCurveT& a, const Vertex& b ) {
103
        return BSplineCurveT(a.v0-b,a.v1-b,a.v2-b,a.v3-b);
104
      }
105

106
      __forceinline BSplineCurveT<Vec3ff> xfm_pr(const LinearSpace3fa& space, const Vec3fa& p) const
107
      {
108
        const Vec3ff q0(xfmVector(space,(Vec3fa)v0-p), v0.w);
109
        const Vec3ff q1(xfmVector(space,(Vec3fa)v1-p), v1.w);
110
        const Vec3ff q2(xfmVector(space,(Vec3fa)v2-p), v2.w);
111
        const Vec3ff q3(xfmVector(space,(Vec3fa)v3-p), v3.w);
112
        return BSplineCurveT<Vec3ff>(q0,q1,q2,q3);
113
      }
114
      
115
      __forceinline Vertex eval(const float t) const 
116
      {
117
        const Vec4<float> b = BSplineBasis::eval(t);
118
        return madd(b.x,v0,madd(b.y,v1,madd(b.z,v2,b.w*v3)));
119
      }
120
      
121
      __forceinline Vertex eval_du(const float t) const
122
      {
123
        const Vec4<float> b = BSplineBasis::derivative(t);
124
        return madd(b.x,v0,madd(b.y,v1,madd(b.z,v2,b.w*v3)));
125
      }
126
      
127
      __forceinline Vertex eval_dudu(const float t) const 
128
      {
129
        const Vec4<float> b = BSplineBasis::derivative2(t);
130
        return madd(b.x,v0,madd(b.y,v1,madd(b.z,v2,b.w*v3)));
131
      }
132

133
      __forceinline void eval(const float t, Vertex& p, Vertex& dp) const
134
      {
135
        p = eval(t);
136
        dp = eval_du(t);
137
      }
138
      
139
      __forceinline void eval(const float t, Vertex& p, Vertex& dp, Vertex& ddp) const
140
      {
141
        p = eval(t);
142
        dp = eval_du(t);
143
        ddp = eval_dudu(t);
144
      }
145

146
      template<int M>
147
      __forceinline Vec4vf<M> veval(const vfloat<M>& t) const 
148
      {
149
        const Vec4vf<M> b = BSplineBasis::eval(t);
150
        return madd(b.x, Vec4vf<M>(v0), madd(b.y, Vec4vf<M>(v1), madd(b.z, Vec4vf<M>(v2), b.w * Vec4vf<M>(v3))));
151
      }
152

153
      template<int M>
154
      __forceinline Vec4vf<M> veval_du(const vfloat<M>& t) const 
155
      {
156
        const Vec4vf<M> b = BSplineBasis::derivative(t);
157
        return madd(b.x, Vec4vf<M>(v0), madd(b.y, Vec4vf<M>(v1), madd(b.z, Vec4vf<M>(v2), b.w * Vec4vf<M>(v3))));
158
      }
159

160
      template<int M>
161
      __forceinline Vec4vf<M> veval_dudu(const vfloat<M>& t) const 
162
      {
163
        const Vec4vf<M> b = BSplineBasis::derivative2(t);
164
        return madd(b.x, Vec4vf<M>(v0), madd(b.y, Vec4vf<M>(v1), madd(b.z, Vec4vf<M>(v2), b.w * Vec4vf<M>(v3))));
165
      }
166

167
      template<int M>
168
      __forceinline void veval(const vfloat<M>& t, Vec4vf<M>& p, Vec4vf<M>& dp) const
169
      {
170
        p = veval<M>(t);
171
        dp = veval_du<M>(t);
172
      }
173
      
174
      template<int M>
175
      __forceinline Vec4vf<M> eval0(const int ofs, const int size) const
176
      {
177
        assert(size <= PrecomputedBSplineBasis::N);
178
        assert(ofs <= size);
179
        return madd(vfloat<M>::loadu(&bspline_basis0.c0[size][ofs]), Vec4vf<M>(v0),
180
                    madd(vfloat<M>::loadu(&bspline_basis0.c1[size][ofs]), Vec4vf<M>(v1),
181
                         madd(vfloat<M>::loadu(&bspline_basis0.c2[size][ofs]), Vec4vf<M>(v2),
182
                              vfloat<M>::loadu(&bspline_basis0.c3[size][ofs]) * Vec4vf<M>(v3))));
183
      }
184
      
185
      template<int M>
186
      __forceinline Vec4vf<M> eval1(const int ofs, const int size) const
187
      {
188
        assert(size <= PrecomputedBSplineBasis::N);
189
        assert(ofs <= size);
190
        return madd(vfloat<M>::loadu(&bspline_basis1.c0[size][ofs]), Vec4vf<M>(v0), 
191
                    madd(vfloat<M>::loadu(&bspline_basis1.c1[size][ofs]), Vec4vf<M>(v1),
192
                         madd(vfloat<M>::loadu(&bspline_basis1.c2[size][ofs]), Vec4vf<M>(v2),
193
                              vfloat<M>::loadu(&bspline_basis1.c3[size][ofs]) * Vec4vf<M>(v3))));
194
      }
195
      
196
      template<int M>
197
      __forceinline Vec4vf<M> derivative0(const int ofs, const int size) const
198
      {
199
        assert(size <= PrecomputedBSplineBasis::N);
200
        assert(ofs <= size);
201
        return madd(vfloat<M>::loadu(&bspline_basis0.d0[size][ofs]), Vec4vf<M>(v0),
202
                    madd(vfloat<M>::loadu(&bspline_basis0.d1[size][ofs]), Vec4vf<M>(v1),
203
                         madd(vfloat<M>::loadu(&bspline_basis0.d2[size][ofs]), Vec4vf<M>(v2),
204
                              vfloat<M>::loadu(&bspline_basis0.d3[size][ofs]) * Vec4vf<M>(v3))));
205
      }
206
      
207
      template<int M>
208
      __forceinline Vec4vf<M> derivative1(const int ofs, const int size) const
209
      {
210
        assert(size <= PrecomputedBSplineBasis::N);
211
        assert(ofs <= size);
212
        return madd(vfloat<M>::loadu(&bspline_basis1.d0[size][ofs]), Vec4vf<M>(v0),
213
                    madd(vfloat<M>::loadu(&bspline_basis1.d1[size][ofs]), Vec4vf<M>(v1),
214
                         madd(vfloat<M>::loadu(&bspline_basis1.d2[size][ofs]), Vec4vf<M>(v2),
215
                              vfloat<M>::loadu(&bspline_basis1.d3[size][ofs]) * Vec4vf<M>(v3))));
216
      }
217
      
218
      /* calculates bounds of bspline curve geometry */
219
      __forceinline BBox3fa accurateRoundBounds() const
220
      {
221
        const int N = 7;
222
        const float scale = 1.0f/(3.0f*(N-1));
223
        Vec4vfx pl(pos_inf), pu(neg_inf);
224
        for (int i=0; i<=N; i+=VSIZEX)
225
        {
226
          vintx vi = vintx(i)+vintx(step);
227
          vboolx valid = vi <= vintx(N);
228
          const Vec4vfx p  = eval0<VSIZEX>(i,N);
229
          const Vec4vfx dp = derivative0<VSIZEX>(i,N);
230
          const Vec4vfx pm = p-Vec4vfx(scale)*select(vi!=vintx(0),dp,Vec4vfx(zero));
231
          const Vec4vfx pp = p+Vec4vfx(scale)*select(vi!=vintx(N),dp,Vec4vfx(zero));
232
          pl = select(valid,min(pl,p,pm,pp),pl); // FIXME: use masked min
233
          pu = select(valid,max(pu,p,pm,pp),pu); // FIXME: use masked min
234
        }
235
        const Vec3fa lower(reduce_min(pl.x),reduce_min(pl.y),reduce_min(pl.z));
236
        const Vec3fa upper(reduce_max(pu.x),reduce_max(pu.y),reduce_max(pu.z));
237
        const float r_min = reduce_min(pl.w);
238
        const float r_max = reduce_max(pu.w);
239
        const Vec3fa upper_r = Vec3fa(max(abs(r_min),abs(r_max)));
240
        return enlarge(BBox3fa(lower,upper),upper_r);
241
      }
242
      
243
      /* calculates bounds when tessellated into N line segments */
244
      __forceinline BBox3fa accurateFlatBounds(int N) const
245
      {
246
        if (likely(N == 4))
247
        {
248
          const Vec4vf4 pi = eval0<4>(0,4);
249
          const Vec3fa lower(reduce_min(pi.x),reduce_min(pi.y),reduce_min(pi.z));
250
          const Vec3fa upper(reduce_max(pi.x),reduce_max(pi.y),reduce_max(pi.z));
251
          const Vec3fa upper_r = Vec3fa(reduce_max(abs(pi.w)));
252
          const Vec3ff pe = end();
253
          return enlarge(BBox3fa(min(lower,pe),max(upper,pe)),max(upper_r,Vec3fa(abs(pe.w))));
254
        } 
255
        else
256
        {
257
          Vec3vfx pl(pos_inf), pu(neg_inf); vfloatx ru(0.0f);
258
          for (int i=0; i<=N; i+=VSIZEX)
259
          {
260
            vboolx valid = vintx(i)+vintx(step) <= vintx(N);
261
            const Vec4vfx pi = eval0<VSIZEX>(i,N);
262
            
263
            pl.x = select(valid,min(pl.x,pi.x),pl.x); // FIXME: use masked min
264
            pl.y = select(valid,min(pl.y,pi.y),pl.y); 
265
            pl.z = select(valid,min(pl.z,pi.z),pl.z); 
266
            
267
            pu.x = select(valid,max(pu.x,pi.x),pu.x); // FIXME: use masked min
268
            pu.y = select(valid,max(pu.y,pi.y),pu.y); 
269
            pu.z = select(valid,max(pu.z,pi.z),pu.z); 
270
            
271
            ru = select(valid,max(ru,abs(pi.w)),ru); 
272
          }
273
          const Vec3fa lower(reduce_min(pl.x),reduce_min(pl.y),reduce_min(pl.z));
274
          const Vec3fa upper(reduce_max(pu.x),reduce_max(pu.y),reduce_max(pu.z));
275
          const Vec3fa upper_r(reduce_max(ru));
276
          return enlarge(BBox3fa(lower,upper),upper_r);
277
        }
278
      }
279
      
280
      friend __forceinline embree_ostream operator<<(embree_ostream cout, const BSplineCurveT& curve) {
281
        return cout << "BSplineCurve { v0 = " << curve.v0 << ", v1 = " << curve.v1 << ", v2 = " << curve.v2 << ", v3 = " << curve.v3 << " }";
282
      }
283
    };
284
  
285
  template<typename Vertex>
286
    __forceinline void convert(const BezierCurveT<Vertex>& icurve, BezierCurveT<Vertex>& ocurve) {
287
    ocurve = icurve;
288
  }
289
  
290
  template<typename Vertex>
291
    __forceinline void convert(const BSplineCurveT<Vertex>& icurve, BSplineCurveT<Vertex>& ocurve) {
292
    ocurve = icurve;
293
  }
294
  
295
  template<typename Vertex>
296
    __forceinline void convert(const BezierCurveT<Vertex>& icurve, BSplineCurveT<Vertex>& ocurve)
297
  {
298
    const Vertex v0 = madd(6.0f,icurve.v0,madd(-7.0f,icurve.v1,2.0f*icurve.v2));
299
    const Vertex v1 = msub(2.0f,icurve.v1,icurve.v2);
300
    const Vertex v2 = msub(2.0f,icurve.v2,icurve.v1);
301
    const Vertex v3 = madd(2.0f,icurve.v1,madd(-7.0f,icurve.v2,6.0f*icurve.v3));
302
    ocurve = BSplineCurveT<Vertex>(v0,v1,v2,v3);
303
  }
304
  
305
  template<typename Vertex>
306
    __forceinline void convert(const BSplineCurveT<Vertex>& icurve, BezierCurveT<Vertex>& ocurve)
307
  {
308
    const Vertex v0 = madd(1.0f/6.0f,icurve.v0,madd(2.0f/3.0f,icurve.v1,1.0f/6.0f*icurve.v2));
309
    const Vertex v1 = madd(2.0f/3.0f,icurve.v1,1.0f/3.0f*icurve.v2);
310
    const Vertex v2 = madd(1.0f/3.0f,icurve.v1,2.0f/3.0f*icurve.v2);
311
    const Vertex v3 = madd(1.0f/6.0f,icurve.v1,madd(2.0f/3.0f,icurve.v2,1.0f/6.0f*icurve.v3));
312
    ocurve = BezierCurveT<Vertex>(v0,v1,v2,v3);
313
  }
314

315
  template<typename CurveGeometry>
316
  __forceinline BSplineCurveT<Vec3ff> enlargeRadiusToMinWidth(const RayQueryContext* context, const CurveGeometry* geom, const Vec3fa& ray_org, const BSplineCurveT<Vec3ff>& curve)
317
  {
318
    return BSplineCurveT<Vec3ff>(enlargeRadiusToMinWidth(context,geom,ray_org,curve.v0),
319
                                 enlargeRadiusToMinWidth(context,geom,ray_org,curve.v1),
320
                                 enlargeRadiusToMinWidth(context,geom,ray_org,curve.v2),
321
                                 enlargeRadiusToMinWidth(context,geom,ray_org,curve.v3));
322
  }
323
  
324
  typedef BSplineCurveT<Vec3fa> BSplineCurve3fa;
325
}
326

327

328
Product

Resources

Company
