1 // Copyright 2009-2021 Intel Corporation
2 // SPDX-License-Identifier: Apache-2.0
3 
4 #pragma once
5 
6 #include "catmullclark_ring.h"
7 #include "bezier_curve.h"
8 
9 namespace embree
10 {
11   template<typename Vertex, typename Vertex_t = Vertex>
12     class __aligned(64) CatmullClarkPatchT
13     {
14     public:
15     typedef CatmullClark1RingT<Vertex,Vertex_t> CatmullClark1Ring;
16     typedef typename CatmullClark1Ring::Type Type;
17 
18     array_t<CatmullClark1RingT<Vertex,Vertex_t>,4> ring;
19 
20     public:
CatmullClarkPatchT()21     __forceinline CatmullClarkPatchT () {}
22 
CatmullClarkPatchT(const HalfEdge * first_half_edge,const char * vertices,size_t stride)23     __forceinline CatmullClarkPatchT (const HalfEdge* first_half_edge, const char* vertices, size_t stride) {
24       init(first_half_edge,vertices,stride);
25     }
26 
CatmullClarkPatchT(const HalfEdge * first_half_edge,const BufferView<Vec3fa> & vertices)27     __forceinline CatmullClarkPatchT (const HalfEdge* first_half_edge, const BufferView<Vec3fa>& vertices) {
28       init(first_half_edge,vertices.getPtr(),vertices.getStride());
29     }
30 
init(const HalfEdge * first_half_edge,const char * vertices,size_t stride)31     __forceinline void init (const HalfEdge* first_half_edge, const char* vertices, size_t stride)
32     {
33       for (unsigned i=0; i<4; i++)
34         ring[i].init(first_half_edge+i,vertices,stride);
35 
36       assert(verify());
37     }
38 
bytes()39     __forceinline size_t bytes() const {
40       return ring[0].bytes()+ring[1].bytes()+ring[2].bytes()+ring[3].bytes();
41     }
42 
serialize(void * ptr,size_t & ofs)43     __forceinline void serialize(void* ptr, size_t& ofs) const
44     {
45       for (size_t i=0; i<4; i++)
46         ring[i].serialize((char*)ptr,ofs);
47     }
48 
deserialize(void * ptr)49     __forceinline void deserialize(void* ptr)
50     {
51       size_t ofs = 0;
52       for (size_t i=0; i<4; i++)
53         ring[i].deserialize((char*)ptr,ofs);
54     }
55 
bounds()56     __forceinline BBox3fa bounds() const
57     {
58       BBox3fa bounds (ring[0].bounds());
59       for (size_t i=1; i<4; i++)
60 	bounds.extend(ring[i].bounds());
61       return bounds;
62     }
63 
type()64     __forceinline Type type() const
65     {
66       const int ty0 = ring[0].type() ^ CatmullClark1Ring::TYPE_CREASES;
67       const int ty1 = ring[1].type() ^ CatmullClark1Ring::TYPE_CREASES;
68       const int ty2 = ring[2].type() ^ CatmullClark1Ring::TYPE_CREASES;
69       const int ty3 = ring[3].type() ^ CatmullClark1Ring::TYPE_CREASES;
70       return (Type) ((ty0 & ty1 & ty2 & ty3) ^ CatmullClark1Ring::TYPE_CREASES);
71     }
72 
isFinalResolution(float res)73     __forceinline bool isFinalResolution(float res) const {
74       return ring[0].isFinalResolution(res) && ring[1].isFinalResolution(res) && ring[2].isFinalResolution(res) && ring[3].isFinalResolution(res);
75     }
76 
init_regular(const CatmullClark1RingT<Vertex,Vertex_t> & p0,const CatmullClark1RingT<Vertex,Vertex_t> & p1,CatmullClark1RingT<Vertex,Vertex_t> & dest0,CatmullClark1RingT<Vertex,Vertex_t> & dest1)77     static __forceinline void init_regular(const CatmullClark1RingT<Vertex,Vertex_t>& p0,
78 					   const CatmullClark1RingT<Vertex,Vertex_t>& p1,
79 					   CatmullClark1RingT<Vertex,Vertex_t>& dest0,
80 					   CatmullClark1RingT<Vertex,Vertex_t>& dest1)
81     {
82       assert(p1.face_valence > 2);
83       dest1.vertex_level = dest0.vertex_level = p0.edge_level;
84       dest1.face_valence = dest0.face_valence = 4;
85       dest1.edge_valence = dest0.edge_valence = 8;
86       dest1.border_index = dest0.border_index = -1;
87       dest1.vtx = dest0.vtx = (Vertex_t)p0.ring[0];
88       dest1.vertex_crease_weight = dest0.vertex_crease_weight = 0.0f;
89 
90       dest1.ring[2] = dest0.ring[0] = (Vertex_t)p0.ring[1];
91       dest1.ring[1] = dest0.ring[7] = (Vertex_t)p1.ring[0];
92       dest1.ring[0] = dest0.ring[6] = (Vertex_t)p1.vtx;
93       dest1.ring[7] = dest0.ring[5] = (Vertex_t)p1.ring[4];
94       dest1.ring[6] = dest0.ring[4] = (Vertex_t)p0.ring[p0.edge_valence-1];
95       dest1.ring[5] = dest0.ring[3] = (Vertex_t)p0.ring[p0.edge_valence-2];
96       dest1.ring[4] = dest0.ring[2] = (Vertex_t)p0.vtx;
97       dest1.ring[3] = dest0.ring[1] = (Vertex_t)p0.ring[2];
98 
99       dest1.crease_weight[1] = dest0.crease_weight[0] = 0.0f;
100       dest1.crease_weight[0] = dest0.crease_weight[3] = p1.crease_weight[1];
101       dest1.crease_weight[3] = dest0.crease_weight[2] = 0.0f;
102       dest1.crease_weight[2] = dest0.crease_weight[1] = p0.crease_weight[0];
103 
104       if (p0.eval_unique_identifier <= p1.eval_unique_identifier)
105       {
106         dest0.eval_start_index = 3;
107         dest1.eval_start_index = 0;
108         dest0.eval_unique_identifier = p0.eval_unique_identifier;
109         dest1.eval_unique_identifier = p0.eval_unique_identifier;
110       }
111       else
112       {
113         dest0.eval_start_index = 1;
114         dest1.eval_start_index = 2;
115         dest0.eval_unique_identifier = p1.eval_unique_identifier;
116         dest1.eval_unique_identifier = p1.eval_unique_identifier;
117       }
118     }
119 
init_border(const CatmullClark1RingT<Vertex,Vertex_t> & p0,const CatmullClark1RingT<Vertex,Vertex_t> & p1,CatmullClark1RingT<Vertex,Vertex_t> & dest0,CatmullClark1RingT<Vertex,Vertex_t> & dest1)120     static __forceinline void init_border(const CatmullClark1RingT<Vertex,Vertex_t> &p0,
121                                           const CatmullClark1RingT<Vertex,Vertex_t> &p1,
122                                           CatmullClark1RingT<Vertex,Vertex_t> &dest0,
123                                           CatmullClark1RingT<Vertex,Vertex_t> &dest1)
124     {
125       dest1.vertex_level = dest0.vertex_level = p0.edge_level;
126       dest1.face_valence = dest0.face_valence = 3;
127       dest1.edge_valence = dest0.edge_valence = 6;
128       dest0.border_index = 2;
129       dest1.border_index = 4;
130       dest1.vtx  = dest0.vtx = (Vertex_t)p0.ring[0];
131       dest1.vertex_crease_weight = dest0.vertex_crease_weight = 0.0f;
132 
133       dest1.ring[2] = dest0.ring[0] = (Vertex_t)p0.ring[1];
134       dest1.ring[1] = dest0.ring[5] = (Vertex_t)p1.ring[0];
135       dest1.ring[0] = dest0.ring[4] = (Vertex_t)p1.vtx;
136       dest1.ring[5] = dest0.ring[3] = (Vertex_t)p0.ring[p0.border_index+1]; // dummy
137       dest1.ring[4] = dest0.ring[2] = (Vertex_t)p0.vtx;
138       dest1.ring[3] = dest0.ring[1] = (Vertex_t)p0.ring[2];
139 
140       dest1.crease_weight[1] = dest0.crease_weight[0] = 0.0f;
141       dest1.crease_weight[0] = dest0.crease_weight[2] = p1.crease_weight[1];
142       dest1.crease_weight[2] = dest0.crease_weight[1] = p0.crease_weight[0];
143 
144       if (p0.eval_unique_identifier <= p1.eval_unique_identifier)
145       {
146         dest0.eval_start_index = 1;
147         dest1.eval_start_index = 2;
148         dest0.eval_unique_identifier = p0.eval_unique_identifier;
149         dest1.eval_unique_identifier = p0.eval_unique_identifier;
150       }
151       else
152       {
153         dest0.eval_start_index = 2;
154         dest1.eval_start_index = 0;
155         dest0.eval_unique_identifier = p1.eval_unique_identifier;
156         dest1.eval_unique_identifier = p1.eval_unique_identifier;
157       }
158     }
159 
init_regular(const Vertex_t & center,const Vertex_t center_ring[8],const unsigned int offset,CatmullClark1RingT<Vertex,Vertex_t> & dest)160     static __forceinline void init_regular(const Vertex_t &center, const Vertex_t center_ring[8], const unsigned int offset, CatmullClark1RingT<Vertex,Vertex_t> &dest)
161     {
162       dest.vertex_level = 0.0f;
163       dest.face_valence = 4;
164       dest.edge_valence = 8;
165       dest.border_index = -1;
166       dest.vtx     = (Vertex_t)center;
167       dest.vertex_crease_weight = 0.0f;
168       for (size_t i=0; i<8; i++)
169 	dest.ring[i] = (Vertex_t)center_ring[(offset+i)%8];
170       for (size_t i=0; i<4; i++)
171         dest.crease_weight[i] = 0.0f;
172 
173       dest.eval_start_index = (8-offset)>>1;
174       if (dest.eval_start_index >= dest.face_valence) dest.eval_start_index -= dest.face_valence;
175       assert( dest.eval_start_index < dest.face_valence );
176       dest.eval_unique_identifier = 0;
177     }
178 
subdivide(array_t<CatmullClarkPatchT,4> & patch)179     __noinline void subdivide(array_t<CatmullClarkPatchT,4>& patch) const
180     {
181       ring[0].subdivide(patch[0].ring[0]);
182       ring[1].subdivide(patch[1].ring[1]);
183       ring[2].subdivide(patch[2].ring[2]);
184       ring[3].subdivide(patch[3].ring[3]);
185 
186       patch[0].ring[0].edge_level = 0.5f*ring[0].edge_level;
187       patch[0].ring[1].edge_level = 0.25f*(ring[1].edge_level+ring[3].edge_level);
188       patch[0].ring[2].edge_level = 0.25f*(ring[0].edge_level+ring[2].edge_level);
189       patch[0].ring[3].edge_level = 0.5f*ring[3].edge_level;
190 
191       patch[1].ring[0].edge_level = 0.5f*ring[0].edge_level;
192       patch[1].ring[1].edge_level = 0.5f*ring[1].edge_level;
193       patch[1].ring[2].edge_level = 0.25f*(ring[0].edge_level+ring[2].edge_level);
194       patch[1].ring[3].edge_level = 0.25f*(ring[1].edge_level+ring[3].edge_level);
195 
196       patch[2].ring[0].edge_level = 0.25f*(ring[0].edge_level+ring[2].edge_level);
197       patch[2].ring[1].edge_level = 0.5f*ring[1].edge_level;
198       patch[2].ring[2].edge_level = 0.5f*ring[2].edge_level;
199       patch[2].ring[3].edge_level = 0.25f*(ring[1].edge_level+ring[3].edge_level);
200 
201       patch[3].ring[0].edge_level = 0.25f*(ring[0].edge_level+ring[2].edge_level);
202       patch[3].ring[1].edge_level = 0.25f*(ring[1].edge_level+ring[3].edge_level);
203       patch[3].ring[2].edge_level = 0.5f*ring[2].edge_level;
204       patch[3].ring[3].edge_level = 0.5f*ring[3].edge_level;
205 
206       const bool regular0 = ring[0].has_last_face() && ring[1].face_valence > 2;
207       if (likely(regular0))
208         init_regular(patch[0].ring[0],patch[1].ring[1],patch[0].ring[1],patch[1].ring[0]);
209       else
210         init_border(patch[0].ring[0],patch[1].ring[1],patch[0].ring[1],patch[1].ring[0]);
211 
212       const bool regular1 = ring[1].has_last_face() && ring[2].face_valence > 2;
213       if (likely(regular1))
214         init_regular(patch[1].ring[1],patch[2].ring[2],patch[1].ring[2],patch[2].ring[1]);
215       else
216         init_border(patch[1].ring[1],patch[2].ring[2],patch[1].ring[2],patch[2].ring[1]);
217 
218       const bool regular2 = ring[2].has_last_face() && ring[3].face_valence > 2;
219       if (likely(regular2))
220         init_regular(patch[2].ring[2],patch[3].ring[3],patch[2].ring[3],patch[3].ring[2]);
221       else
222         init_border(patch[2].ring[2],patch[3].ring[3],patch[2].ring[3],patch[3].ring[2]);
223 
224       const bool regular3 = ring[3].has_last_face() && ring[0].face_valence > 2;
225       if (likely(regular3))
226         init_regular(patch[3].ring[3],patch[0].ring[0],patch[3].ring[0],patch[0].ring[3]);
227       else
228         init_border(patch[3].ring[3],patch[0].ring[0],patch[3].ring[0],patch[0].ring[3]);
229 
230       Vertex_t center = (ring[0].vtx + ring[1].vtx + ring[2].vtx + ring[3].vtx) * 0.25f;
231 
232       Vertex_t center_ring[8];
233       center_ring[0] = (Vertex_t)patch[3].ring[3].ring[0];
234       center_ring[7] = (Vertex_t)patch[3].ring[3].vtx;
235       center_ring[6] = (Vertex_t)patch[2].ring[2].ring[0];
236       center_ring[5] = (Vertex_t)patch[2].ring[2].vtx;
237       center_ring[4] = (Vertex_t)patch[1].ring[1].ring[0];
238       center_ring[3] = (Vertex_t)patch[1].ring[1].vtx;
239       center_ring[2] = (Vertex_t)patch[0].ring[0].ring[0];
240       center_ring[1] = (Vertex_t)patch[0].ring[0].vtx;
241 
242       init_regular(center,center_ring,0,patch[0].ring[2]);
243       init_regular(center,center_ring,2,patch[1].ring[3]);
244       init_regular(center,center_ring,4,patch[2].ring[0]);
245       init_regular(center,center_ring,6,patch[3].ring[1]);
246 
247       assert(patch[0].verify());
248       assert(patch[1].verify());
249       assert(patch[2].verify());
250       assert(patch[3].verify());
251     }
252 
verify()253     bool verify() const {
254       return ring[0].hasValidPositions() && ring[1].hasValidPositions() && ring[2].hasValidPositions() && ring[3].hasValidPositions();
255     }
256 
init(FinalQuad & quad)257     __forceinline void init( FinalQuad& quad ) const
258     {
259       quad.vtx[0] = (Vertex_t)ring[0].vtx;
260       quad.vtx[1] = (Vertex_t)ring[1].vtx;
261       quad.vtx[2] = (Vertex_t)ring[2].vtx;
262       quad.vtx[3] = (Vertex_t)ring[3].vtx;
263     };
264 
265     friend __forceinline embree_ostream operator<<(embree_ostream o, const CatmullClarkPatchT &p)
266     {
267       o << "CatmullClarkPatch { " << embree_endl;
268       for (size_t i=0; i<4; i++)
269 	o << "ring" << i << ": " << p.ring[i] << embree_endl;
270       o << "}" << embree_endl;
271       return o;
272     }
273     };
274 
275   typedef CatmullClarkPatchT<Vec3fa,Vec3fa_t> CatmullClarkPatch3fa;
276 
277   template<typename Vertex, typename Vertex_t = Vertex>
278     class __aligned(64) GeneralCatmullClarkPatchT
279     {
280     public:
281     typedef CatmullClarkPatchT<Vertex,Vertex_t> CatmullClarkPatch;
282     typedef CatmullClark1RingT<Vertex,Vertex_t> CatmullClark1Ring;
283     typedef BezierCurveT<Vertex> BezierCurve;
284 
285     static const unsigned SIZE = MAX_PATCH_VALENCE;
286     DynamicStackArray<GeneralCatmullClark1RingT<Vertex,Vertex_t>,8,SIZE> ring;
287     unsigned N;
288 
GeneralCatmullClarkPatchT()289     __forceinline GeneralCatmullClarkPatchT ()
290     : N(0) {}
291 
GeneralCatmullClarkPatchT(const HalfEdge * h,const char * vertices,size_t stride)292     GeneralCatmullClarkPatchT (const HalfEdge* h, const char* vertices, size_t stride) {
293       init(h,vertices,stride);
294     }
295 
GeneralCatmullClarkPatchT(const HalfEdge * first_half_edge,const BufferView<Vec3fa> & vertices)296     __forceinline GeneralCatmullClarkPatchT (const HalfEdge* first_half_edge, const BufferView<Vec3fa>& vertices) {
297       init(first_half_edge,vertices.getPtr(),vertices.getStride());
298     }
299 
init(const HalfEdge * h,const char * vertices,size_t stride)300     __forceinline void init (const HalfEdge* h, const char* vertices, size_t stride)
301     {
302       unsigned int i = 0;
303       const HalfEdge* edge = h;
304       do {
305         ring[i].init(edge,vertices,stride);
306         edge = edge->next();
307         i++;
308       } while ((edge != h) && (i < SIZE));
309       N = i;
310     }
311 
size()312     __forceinline unsigned size() const {
313       return N;
314     }
315 
isQuadPatch()316     __forceinline bool isQuadPatch() const {
317       return (N == 4) && ring[0].only_quads && ring[1].only_quads && ring[2].only_quads && ring[3].only_quads;
318     }
319 
init_regular(const CatmullClark1RingT<Vertex,Vertex_t> & p0,const CatmullClark1RingT<Vertex,Vertex_t> & p1,CatmullClark1RingT<Vertex,Vertex_t> & dest0,CatmullClark1RingT<Vertex,Vertex_t> & dest1)320     static __forceinline void init_regular(const CatmullClark1RingT<Vertex,Vertex_t>& p0,
321 					   const CatmullClark1RingT<Vertex,Vertex_t>& p1,
322 					   CatmullClark1RingT<Vertex,Vertex_t>& dest0,
323 					   CatmullClark1RingT<Vertex,Vertex_t>& dest1)
324     {
325       assert(p1.face_valence > 2);
326       dest1.vertex_level = dest0.vertex_level = p0.edge_level;
327       dest1.face_valence = dest0.face_valence = 4;
328       dest1.edge_valence = dest0.edge_valence = 8;
329       dest1.border_index = dest0.border_index = -1;
330       dest1.vtx = dest0.vtx = (Vertex_t)p0.ring[0];
331       dest1.vertex_crease_weight = dest0.vertex_crease_weight = 0.0f;
332 
333       dest1.ring[2] = dest0.ring[0] = (Vertex_t)p0.ring[1];
334       dest1.ring[1] = dest0.ring[7] = (Vertex_t)p1.ring[0];
335       dest1.ring[0] = dest0.ring[6] = (Vertex_t)p1.vtx;
336       dest1.ring[7] = dest0.ring[5] = (Vertex_t)p1.ring[4];
337       dest1.ring[6] = dest0.ring[4] = (Vertex_t)p0.ring[p0.edge_valence-1];
338       dest1.ring[5] = dest0.ring[3] = (Vertex_t)p0.ring[p0.edge_valence-2];
339       dest1.ring[4] = dest0.ring[2] = (Vertex_t)p0.vtx;
340       dest1.ring[3] = dest0.ring[1] = (Vertex_t)p0.ring[2];
341 
342       dest1.crease_weight[1] = dest0.crease_weight[0] = 0.0f;
343       dest1.crease_weight[0] = dest0.crease_weight[3] = p1.crease_weight[1];
344       dest1.crease_weight[3] = dest0.crease_weight[2] = 0.0f;
345       dest1.crease_weight[2] = dest0.crease_weight[1] = p0.crease_weight[0];
346 
347       if (p0.eval_unique_identifier <= p1.eval_unique_identifier)
348       {
349         dest0.eval_start_index = 3;
350         dest1.eval_start_index = 0;
351         dest0.eval_unique_identifier = p0.eval_unique_identifier;
352         dest1.eval_unique_identifier = p0.eval_unique_identifier;
353       }
354       else
355       {
356         dest0.eval_start_index = 1;
357         dest1.eval_start_index = 2;
358         dest0.eval_unique_identifier = p1.eval_unique_identifier;
359         dest1.eval_unique_identifier = p1.eval_unique_identifier;
360       }
361     }
362 
363 
init_border(const CatmullClark1RingT<Vertex,Vertex_t> & p0,const CatmullClark1RingT<Vertex,Vertex_t> & p1,CatmullClark1RingT<Vertex,Vertex_t> & dest0,CatmullClark1RingT<Vertex,Vertex_t> & dest1)364     static __forceinline void init_border(const CatmullClark1RingT<Vertex,Vertex_t> &p0,
365                                           const CatmullClark1RingT<Vertex,Vertex_t> &p1,
366                                           CatmullClark1RingT<Vertex,Vertex_t> &dest0,
367                                           CatmullClark1RingT<Vertex,Vertex_t> &dest1)
368     {
369       dest1.vertex_level = dest0.vertex_level = p0.edge_level;
370       dest1.face_valence = dest0.face_valence = 3;
371       dest1.edge_valence = dest0.edge_valence = 6;
372       dest0.border_index = 2;
373       dest1.border_index = 4;
374       dest1.vtx  = dest0.vtx = (Vertex_t)p0.ring[0];
375       dest1.vertex_crease_weight = dest0.vertex_crease_weight = 0.0f;
376 
377       dest1.ring[2] = dest0.ring[0] = (Vertex_t)p0.ring[1];
378       dest1.ring[1] = dest0.ring[5] = (Vertex_t)p1.ring[0];
379       dest1.ring[0] = dest0.ring[4] = (Vertex_t)p1.vtx;
380       dest1.ring[5] = dest0.ring[3] = (Vertex_t)p0.ring[p0.border_index+1]; // dummy
381       dest1.ring[4] = dest0.ring[2] = (Vertex_t)p0.vtx;
382       dest1.ring[3] = dest0.ring[1] = (Vertex_t)p0.ring[2];
383 
384       dest1.crease_weight[1] = dest0.crease_weight[0] = 0.0f;
385       dest1.crease_weight[0] = dest0.crease_weight[2] = p1.crease_weight[1];
386       dest1.crease_weight[2] = dest0.crease_weight[1] = p0.crease_weight[0];
387 
388       if (p0.eval_unique_identifier <= p1.eval_unique_identifier)
389       {
390         dest0.eval_start_index = 1;
391         dest1.eval_start_index = 2;
392         dest0.eval_unique_identifier = p0.eval_unique_identifier;
393         dest1.eval_unique_identifier = p0.eval_unique_identifier;
394       }
395       else
396       {
397         dest0.eval_start_index = 2;
398         dest1.eval_start_index = 0;
399         dest0.eval_unique_identifier = p1.eval_unique_identifier;
400         dest1.eval_unique_identifier = p1.eval_unique_identifier;
401       }
402     }
403 
init_regular(const Vertex_t & center,const array_t<Vertex_t,2* SIZE> & center_ring,const float vertex_level,const unsigned int N,const unsigned int offset,CatmullClark1RingT<Vertex,Vertex_t> & dest)404     static __forceinline void init_regular(const Vertex_t &center, const array_t<Vertex_t,2*SIZE>& center_ring, const float vertex_level, const unsigned int N, const unsigned int offset, CatmullClark1RingT<Vertex,Vertex_t> &dest)
405     {
406       assert(N<(MAX_RING_FACE_VALENCE));
407       assert(2*N<(MAX_RING_EDGE_VALENCE));
408       dest.vertex_level = vertex_level;
409       dest.face_valence = N;
410       dest.edge_valence = 2*N;
411       dest.border_index = -1;
412       dest.vtx     = (Vertex_t)center;
413       dest.vertex_crease_weight = 0.0f;
414       for (unsigned i=0; i<2*N; i++) {
415         dest.ring[i] = (Vertex_t)center_ring[(2*N+offset+i-1)%(2*N)];
416         assert(isvalid(dest.ring[i]));
417       }
418       for (unsigned i=0; i<N; i++)
419         dest.crease_weight[i] = 0.0f;
420 
421       assert(offset <= 2*N);
422       dest.eval_start_index = (2*N-offset)>>1;
423       if (dest.eval_start_index >= dest.face_valence) dest.eval_start_index -= dest.face_valence;
424 
425       assert( dest.eval_start_index < dest.face_valence );
426       dest.eval_unique_identifier = 0;
427     }
428 
subdivide(array_t<CatmullClarkPatch,SIZE> & patch,unsigned & N_o)429     __noinline void subdivide(array_t<CatmullClarkPatch,SIZE>& patch, unsigned& N_o) const
430     {
431       N_o = N;
432       assert( N );
433       for (unsigned i=0; i<N; i++) {
434         unsigned ip1 = (i+1)%N; // FIXME: %
435         ring[i].subdivide(patch[i].ring[0]);
436         patch[i]  .ring[0].edge_level = 0.5f*ring[i].edge_level;
437         patch[ip1].ring[3].edge_level = 0.5f*ring[i].edge_level;
438 
439 	assert( patch[i].ring[0].hasValidPositions() );
440 
441       }
442       assert(N < 2*SIZE);
443       Vertex_t center = Vertex_t(0.0f);
444       array_t<Vertex_t,2*SIZE> center_ring;
445       float center_vertex_level = 2.0f; // guarantees that irregular vertices get always isolated also for non-quads
446 
447       for (unsigned i=0; i<N; i++)
448       {
449         unsigned ip1 = (i+1)%N; // FIXME: %
450         unsigned im1 = (i+N-1)%N; // FIXME: %
451         bool regular = ring[i].has_last_face() && ring[ip1].face_valence > 2;
452         if (likely(regular)) init_regular(patch[i].ring[0],patch[ip1].ring[0],patch[i].ring[1],patch[ip1].ring[3]);
453         else                 init_border (patch[i].ring[0],patch[ip1].ring[0],patch[i].ring[1],patch[ip1].ring[3]);
454 
455 	assert( patch[i].ring[1].hasValidPositions() );
456 	assert( patch[ip1].ring[3].hasValidPositions() );
457 
458 	float level = 0.25f*(ring[im1].edge_level+ring[ip1].edge_level);
459         patch[i].ring[1].edge_level = patch[ip1].ring[2].edge_level = level;
460 	center_vertex_level = max(center_vertex_level,level);
461 
462         center += ring[i].vtx;
463         center_ring[2*i+0] = (Vertex_t)patch[i].ring[0].vtx;
464         center_ring[2*i+1] = (Vertex_t)patch[i].ring[0].ring[0];
465       }
466       center /= float(N);
467 
468       for (unsigned int i=0; i<N; i++) {
469         init_regular(center,center_ring,center_vertex_level,N,2*i,patch[i].ring[2]);
470 
471 	assert( patch[i].ring[2].hasValidPositions() );
472       }
473     }
474 
init(CatmullClarkPatch & patch)475     void init(CatmullClarkPatch& patch) const
476     {
477       assert(size() == 4);
478       ring[0].convert(patch.ring[0]);
479       ring[1].convert(patch.ring[1]);
480       ring[2].convert(patch.ring[2]);
481       ring[3].convert(patch.ring[3]);
482     }
483 
fix_quad_ring_order(array_t<CatmullClarkPatch,GeneralCatmullClarkPatchT::SIZE> & patches)484     static void fix_quad_ring_order (array_t<CatmullClarkPatch,GeneralCatmullClarkPatchT::SIZE>& patches)
485     {
486       CatmullClark1Ring patches1ring1 = patches[1].ring[1];
487       patches[1].ring[1] = patches[1].ring[0]; // FIXME: optimize these assignments
488       patches[1].ring[0] = patches[1].ring[3];
489       patches[1].ring[3] = patches[1].ring[2];
490       patches[1].ring[2] = patches1ring1;
491 
492       CatmullClark1Ring patches2ring2 = patches[2].ring[2];
493       patches[2].ring[2] = patches[2].ring[0];
494       patches[2].ring[0] = patches2ring2;
495       CatmullClark1Ring patches2ring3 = patches[2].ring[3];
496       patches[2].ring[3] = patches[2].ring[1];
497       patches[2].ring[1] = patches2ring3;
498 
499       CatmullClark1Ring patches3ring3 = patches[3].ring[3];
500       patches[3].ring[3] = patches[3].ring[0];
501       patches[3].ring[0] = patches[3].ring[1];
502       patches[3].ring[1] = patches[3].ring[2];
503       patches[3].ring[2] = patches3ring3;
504     }
505 
getLimitBorder(BezierCurve curves[GeneralCatmullClarkPatchT::SIZE])506     __forceinline void getLimitBorder(BezierCurve curves[GeneralCatmullClarkPatchT::SIZE]) const
507     {
508       Vertex P0 = ring[0].getLimitVertex();
509       for (unsigned i=0; i<N; i++)
510       {
511         const unsigned i0 = i, i1 = i+1==N ? 0 : i+1;
512         const Vertex P1 = madd(1.0f/3.0f,ring[i0].getLimitTangent(),P0);
513         const Vertex P3 = ring[i1].getLimitVertex();
514         const Vertex P2 = madd(1.0f/3.0f,ring[i1].getSecondLimitTangent(),P3);
515         new (&curves[i]) BezierCurve(P0,P1,P2,P3);
516         P0 = P3;
517       }
518     }
519 
getLimitBorder(BezierCurve curves[2],const unsigned subPatch)520     __forceinline void getLimitBorder(BezierCurve curves[2], const unsigned subPatch) const
521     {
522       const unsigned i0 = subPatch;
523       const Vertex t0_p = ring[i0].getLimitTangent();
524       const Vertex t0_m = ring[i0].getSecondLimitTangent();
525 
526       const unsigned i1 = subPatch+1 == N ? 0 : subPatch+1;
527       const Vertex t1_p = ring[i1].getLimitTangent();
528       const Vertex t1_m = ring[i1].getSecondLimitTangent();
529 
530       const unsigned i2 = subPatch == 0 ? N-1 : subPatch-1;
531       const Vertex t2_p = ring[i2].getLimitTangent();
532       const Vertex t2_m = ring[i2].getSecondLimitTangent();
533 
534       const Vertex b00 = ring[i0].getLimitVertex();
535       const Vertex b03 = ring[i1].getLimitVertex();
536       const Vertex b33 = ring[i2].getLimitVertex();
537 
538       const Vertex b01 = madd(1.0/3.0f,t0_p,b00);
539       const Vertex b11 = madd(1.0/3.0f,t0_m,b00);
540 
541       //const Vertex b13 = madd(1.0/3.0f,t1_p,b03);
542       const Vertex b02 = madd(1.0/3.0f,t1_m,b03);
543 
544       const Vertex b22 = madd(1.0/3.0f,t2_p,b33);
545       const Vertex b23 = madd(1.0/3.0f,t2_m,b33);
546 
547       new (&curves[0]) BezierCurve(b00,b01,b02,b03);
548       new (&curves[1]) BezierCurve(b33,b22,b11,b00);
549     }
550 
551     friend __forceinline embree_ostream operator<<(embree_ostream o, const GeneralCatmullClarkPatchT &p)
552     {
553       o << "GeneralCatmullClarkPatch { " << embree_endl;
554       for (unsigned i=0; i<p.N; i++)
555 	o << "ring" << i << ": " << p.ring[i] << embree_endl;
556       o << "}" << embree_endl;
557       return o;
558     }
559     };
560 
561   typedef GeneralCatmullClarkPatchT<Vec3fa,Vec3fa_t> GeneralCatmullClarkPatch3fa;
562 }
563