Logo ROOT   6.16/01
Reference Guide
TGLMarchingCubes.cxx
Go to the documentation of this file.
1// @(#)root/gl:$Id$
2// Author: Timur Pocheptsov 06/01/2009
3
4/*************************************************************************
5 * Copyright (C) 1995-2009, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include <algorithm>
13#include <cmath>
14
15#include "TError.h"
16#include "TF3.h"
17
18#include "TGLMarchingCubes.h"
19
20/*
21Implementation of "marching cubes" algortihm for GL module. Used by
22TF3GLPainter, TGLIsoPainter, TGL5DPainter.
23Good and clear algorithm explanation can be found here:
24http://local.wasp.uwa.edu.au/~pbourke/geometry/polygonise/
25*/
26
27namespace Rgl {
28namespace Mc {
29
30/*
31Some routines, values and tables for marching cube method.
32*/
33extern const UInt_t eInt[256];
34extern const Float_t vOff[8][3];
35extern const UChar_t eConn[12][2];
36extern const Float_t eDir[12][3];
37extern const Int_t conTbl[256][16];
38
39namespace {
40
41enum ECubeBitMasks {
42 k0 = 0x1,
43 k1 = 0x2,
44 k2 = 0x4,
45 k3 = 0x8,
46 k4 = 0x10,
47 k5 = 0x20,
48 k6 = 0x40,
49 k7 = 0x80,
50 k8 = 0x100,
51 k9 = 0x200,
52 k10 = 0x400,
53 k11 = 0x800,
54 //
55 k1_5 = k1 | k5,
56 k2_6 = k2 | k6,
57 k3_7 = k3 | k7,
58 k4_5_6_7 = k4 | k5 | k6 | k7,
59 k5_6 = k5 | k6,
60 k0_1_2_3_7_8_11 = k0 | k1 | k2 | k3 | k7 | k8 | k11,
61 k6_7 = k6 | k7
62};
63
64////////////////////////////////////////////////////////////////////////////////
65
66template<class E, class V>
67void ConnectTriangles(TCell<E> &cell, TIsoMesh<V> *mesh, V eps)
68{
69 UInt_t t[3];
70 for (UInt_t i = 0; i < 5; ++i) {
71 if (conTbl[cell.fType][3 * i] < 0)
72 break;
73 for (Int_t j = 2; j >= 0; --j)
74 t[j] = cell.fIds[conTbl[cell.fType][3 * i + j]];
75
76 const V *v0 = &mesh->fVerts[t[0] * 3];
77 const V *v1 = &mesh->fVerts[t[1] * 3];
78 const V *v2 = &mesh->fVerts[t[2] * 3];
79
80 if (std::abs(v0[0] - v1[0]) < eps &&
81 std::abs(v0[1] - v1[1]) < eps &&
82 std::abs(v0[2] - v1[2]) < eps)
83 continue;
84
85 if (std::abs(v2[0] - v1[0]) < eps &&
86 std::abs(v2[1] - v1[1]) < eps &&
87 std::abs(v2[2] - v1[2]) < eps)
88 continue;
89
90 if (std::abs(v0[0] - v2[0]) < eps &&
91 std::abs(v0[1] - v2[1]) < eps &&
92 std::abs(v0[2] - v2[2]) < eps)
93 continue;
94
95 mesh->AddTriangle(t);
96 }
97}
98
99}//unnamed namespace.
100
101/*
102TF3Adapter.
103*/
104////////////////////////////////////////////////////////////////////////////////
105
107{
108 fTF3 = f3;
109 fW = f3->GetXaxis()->GetNbins();//f3->GetNpx();
110 fH = f3->GetYaxis()->GetNbins();//f3->GetNpy();
111 fD = f3->GetZaxis()->GetNbins();//f3->GetNpz();
112}
113
114////////////////////////////////////////////////////////////////////////////////
115
117{
121}
122
123/*
124TF3 split edge implementation.
125*/
126////////////////////////////////////////////////////////////////////////////////
127///Split the edge and find normal in a new vertex.
128
130 Double_t x, Double_t y, Double_t z, Double_t iso)const
131{
132 Double_t v[3] = {};
133 const Double_t ofst = GetOffset(cell.fVals[eConn[i][0]], cell.fVals[eConn[i][1]], iso);
134 v[0] = x + (vOff[eConn[i][0]][0] + ofst * eDir[i][0]) * fStepX;
135 v[1] = y + (vOff[eConn[i][0]][1] + ofst * eDir[i][1]) * fStepY;
136 v[2] = z + (vOff[eConn[i][0]][2] + ofst * eDir[i][2]) * fStepZ;
137 cell.fIds[i] = mesh->AddVertex(v);
138
139 const Double_t stepXU = fStepX * fXScaleInverted;
140 const Double_t xU = x * fXScaleInverted;
141 const Double_t stepYU = fStepY * fYScaleInverted;
142 const Double_t yU = y * fYScaleInverted;
143 const Double_t stepZU = fStepZ * fZScaleInverted;
144 const Double_t zU = z * fZScaleInverted;
145
146 Double_t vU[3] = {};//U - unscaled.
147 vU[0] = xU + (vOff[eConn[i][0]][0] + ofst * eDir[i][0]) * stepXU;
148 vU[1] = yU + (vOff[eConn[i][0]][1] + ofst * eDir[i][1]) * stepYU;
149 vU[2] = zU + (vOff[eConn[i][0]][2] + ofst * eDir[i][2]) * stepZU;
150 //Find normals.
151 Double_t n[3];
152 n[0] = fTF3->Eval(vU[0] - 0.1 * stepXU, vU[1], vU[2]) -
153 fTF3->Eval(vU[0] + 0.1 * stepXU, vU[1], vU[2]);
154 n[1] = fTF3->Eval(vU[0], vU[1] - 0.1 * stepYU, vU[2]) -
155 fTF3->Eval(vU[0], vU[1] + 0.1 * stepYU, vU[2]);
156 n[2] = fTF3->Eval(vU[0], vU[1], vU[2] - 0.1 * stepZU) -
157 fTF3->Eval(vU[0], vU[1], vU[2] + 0.1 * stepZU);
158
159 const Double_t len = std::sqrt(n[0] * n[0] + n[1] * n[1] + n[2] * n[2]);
160 if (len > 1e-7) {
161 n[0] /= len;
162 n[1] /= len;
163 n[2] /= len;
164 }
165
166 mesh->AddNormal(n);
167}
168/*
169TMeshBuilder's implementation.
170"this->" is used with type-dependant names
171in templates.
172*/
173////////////////////////////////////////////////////////////////////////////////
174///Build iso-mesh using marching cubes.
175
176template<class D, class V>
178 MeshType_t *m, V iso)
179{
180 static_cast<TGridGeometry<V> &>(*this) = g;
181
182 this->SetDataSource(s);
183
184 if (GetW() < 2 || GetH() < 2 || GetD() < 2) {
185 Error("TMeshBuilder::BuildMesh",
186 "Bad grid size, one of dimensions is less than 2");
187 return;
188 }
189
190 fSlices[0].ResizeSlice(GetW() - 1, GetH() - 1);
191 fSlices[1].ResizeSlice(GetW() - 1, GetH() - 1);
192
193 this->SetNormalEvaluator(s);
194
195 fMesh = m;
196 fIso = iso;
197
198 SliceType_t *slice1 = fSlices;
199 SliceType_t *slice2 = fSlices + 1;
200
201 this->FetchDensities();
202 NextStep(0, 0, slice1);
203
204 for (UInt_t i = 1, e = GetD(); i < e - 1; ++i) {
205 NextStep(i, slice1, slice2);
206 std::swap(slice1, slice2);
207 }
208
209 if(fAvgNormals)
210 BuildNormals();
211}
212
213////////////////////////////////////////////////////////////////////////////////
214///Fill slice with vertices and triangles.
215
216template<class D, class V>
218 SliceType_t *curr)const
219{
220 if (!prevSlice) {
221 //The first slice in mc grid.
222 BuildFirstCube(curr);
223 BuildRow(curr);
224 BuildCol(curr);
225 BuildSlice(curr);
226 } else {
227 BuildFirstCube(depth, prevSlice, curr);
228 BuildRow(depth, prevSlice, curr);
229 BuildCol(depth, prevSlice, curr);
230 BuildSlice(depth, prevSlice, curr);
231 }
232}
233
234////////////////////////////////////////////////////////////////////////////////
235///The first cube in a grid. nx == 0, ny == 0, nz ==0.
236
237template<class D, class V>
239{
240 CellType_t & cell = s->fCells[0];
241 cell.fVals[0] = GetData(0, 0, 0);
242 cell.fVals[1] = GetData(1, 0, 0);
243 cell.fVals[2] = GetData(1, 1, 0);
244 cell.fVals[3] = GetData(0, 1, 0);
245 cell.fVals[4] = GetData(0, 0, 1);
246 cell.fVals[5] = GetData(1, 0, 1);
247 cell.fVals[6] = GetData(1, 1, 1);
248 cell.fVals[7] = GetData(0, 1, 1);
249
250 cell.fType = 0;
251 for (UInt_t i = 0; i < 8; ++i) {
252 if (cell.fVals[i] <= fIso)
253 cell.fType |= 1 << i;
254 }
255
256 for (UInt_t i = 0, edges = eInt[cell.fType]; i < 12; ++i) {
257 if (edges & (1 << i))
258 SplitEdge(cell, fMesh, i, this->fMinX, this->fMinY, this->fMinZ, fIso);
259 }
260
261 ConnectTriangles(cell, fMesh, fEpsilon);
262}
263
264////////////////////////////////////////////////////////////////////////////////
265///The first row (along x) in the first slice:
266///ny == 0, nz == 0, nx : [1, W - 1].
267///Each cube has previous cube.
268///Values 0, 3, 4, 7 are taken from the previous cube.
269///Edges 3, 7, 8, 11 are taken from the previous cube.
270
271template<class D, class V>
273{
274 for (UInt_t i = 1, e = GetW() - 1; i < e; ++i) {
275 const CellType_t &prev = s->fCells[i - 1];
276 CellType_t &cell = s->fCells[i];
277 cell.fType = 0;
278
279 cell.fVals[0] = prev.fVals[1], cell.fVals[4] = prev.fVals[5];
280 cell.fVals[7] = prev.fVals[6], cell.fVals[3] = prev.fVals[2];
281 cell.fType |= (prev.fType & k1_5) >> 1;
282 cell.fType |= (prev.fType & k2_6) << 1;
283
284 if ((cell.fVals[1] = GetData(i + 1, 0, 0)) <= fIso)
285 cell.fType |= k1;
286 if ((cell.fVals[2] = GetData(i + 1, 1, 0)) <= fIso)
287 cell.fType |= k2;
288 if ((cell.fVals[5] = GetData(i + 1, 0, 1)) <= fIso)
289 cell.fType |= k5;
290 if ((cell.fVals[6] = GetData(i + 1, 1, 1)) <= fIso)
291 cell.fType |= k6;
292
293 const UInt_t edges = eInt[cell.fType];
294 if (!edges)
295 continue;
296 //1. Take edges 3, 7, 8, 11 from the previous cube.
297 if (edges & k3)
298 cell.fIds[3] = prev.fIds[1];
299 if (edges & k7)
300 cell.fIds[7] = prev.fIds[5];
301 if (edges & k8)
302 cell.fIds[8] = prev.fIds[9];
303 if (edges & k11)
304 cell.fIds[11] = prev.fIds[10];
305 //2. Intersect edges 0, 1, 2, 4, 5, 6, 9, 10.
306 const V x = this->fMinX + i * this->fStepX;
307 if (edges & k0)
308 SplitEdge(cell, fMesh, 0, x, this->fMinY, this->fMinZ, fIso);
309 if (edges & k1)
310 SplitEdge(cell, fMesh, 1, x, this->fMinY, this->fMinZ, fIso);
311 if (edges & k2)
312 SplitEdge(cell, fMesh, 2, x, this->fMinY, this->fMinZ, fIso);
313 if (edges & k4)
314 SplitEdge(cell, fMesh, 4, x, this->fMinY, this->fMinZ, fIso);
315 if (edges & k5)
316 SplitEdge(cell, fMesh, 5, x, this->fMinY, this->fMinZ, fIso);
317 if (edges & k6)
318 SplitEdge(cell, fMesh, 6, x, this->fMinY, this->fMinZ, fIso);
319 if (edges & k9)
320 SplitEdge(cell, fMesh, 9, x, this->fMinY, this->fMinZ, fIso);
321 if (edges & k10)
322 SplitEdge(cell, fMesh, 10, x, this->fMinY, this->fMinZ, fIso);
323 //3. Connect new triangles.
324 ConnectTriangles(cell, fMesh, fEpsilon);
325 }
326}
327
328////////////////////////////////////////////////////////////////////////////////
329///"Col" (column) consists of cubes along y axis
330///on the first slice (nx == 0, nz == 0).
331///Each cube has a previous cube and shares values:
332///0, 1, 4, 5 (in prev.: 3, 2, 7, 6); and edges:
333///0, 4, 8, 9 (in prev.: 2, 6, 10, 11).
334
335template<class D, class V>
337{
338 const UInt_t w = GetW();
339 const UInt_t h = GetH();
340
341 for (UInt_t i = 1; i < h - 1; ++i) {
342 const CellType_t &prev = s->fCells[(i - 1) * (w - 1)];
343 CellType_t &cell = s->fCells[i * (w - 1)];
344 cell.fType = 0;
345 //Take values 0, 1, 4, 5 from the prev. cube.
346 cell.fVals[0] = prev.fVals[3], cell.fVals[1] = prev.fVals[2];
347 cell.fVals[4] = prev.fVals[7], cell.fVals[5] = prev.fVals[6];
348 cell.fType |= (prev.fType & k2_6) >> 1;
349 cell.fType |= (prev.fType & k3_7) >> 3;
350 //Calculate values 2, 3, 6, 7.
351 if((cell.fVals[2] = GetData(1, i + 1, 0)) <= fIso)
352 cell.fType |= k2;
353 if((cell.fVals[3] = GetData(0, i + 1, 0)) <= fIso)
354 cell.fType |= k3;
355 if((cell.fVals[6] = GetData(1, i + 1, 1)) <= fIso)
356 cell.fType |= k6;
357 if((cell.fVals[7] = GetData(0, i + 1, 1)) <= fIso)
358 cell.fType |= k7;
359
360 const UInt_t edges = eInt[cell.fType];
361 if(!edges)
362 continue;
363 //Take edges from the previous cube.
364 if (edges & k0)
365 cell.fIds[0] = prev.fIds[2];
366 if (edges & k4)
367 cell.fIds[4] = prev.fIds[6];
368 if (edges & k9)
369 cell.fIds[9] = prev.fIds[10];
370 if (edges & k8)
371 cell.fIds[8] = prev.fIds[11];
372 //Find the remaining edges.
373 const V y = this->fMinY + i * this->fStepY;
374
375 if (edges & k1)
376 SplitEdge(cell, fMesh, 1, this->fMinX, y, this->fMinZ, fIso);
377 if (edges & k2)
378 SplitEdge(cell, fMesh, 2, this->fMinX, y, this->fMinZ, fIso);
379 if (edges & k3)
380 SplitEdge(cell, fMesh, 3, this->fMinX, y, this->fMinZ, fIso);
381 if (edges & k5)
382 SplitEdge(cell, fMesh, 5, this->fMinX, y, this->fMinZ, fIso);
383 if (edges & k6)
384 SplitEdge(cell, fMesh, 6, this->fMinX, y, this->fMinZ, fIso);
385 if (edges & k7)
386 SplitEdge(cell, fMesh, 7, this->fMinX, y, this->fMinZ, fIso);
387 if (edges & k10)
388 SplitEdge(cell, fMesh, 10, this->fMinX, y, this->fMinZ, fIso);
389 if (edges & k11)
390 SplitEdge(cell, fMesh, 11, this->fMinX, y, this->fMinZ, fIso);
391
392 ConnectTriangles(cell, fMesh, fEpsilon);
393 }
394}
395
396////////////////////////////////////////////////////////////////////////////////
397///Slice with nz == 0.
398///nx : [1, W - 1], ny : [1, H - 1].
399///nx increased inside inner loop, ny - enclosing loop.
400///Each cube has two neighbours: ny - 1 => "left",
401///nx - 1 => "right".
402
403template<class D, class V>
405{
406 const UInt_t w = GetW();
407 const UInt_t h = GetH();
408
409 for (UInt_t i = 1; i < h - 1; ++i) {
410 const V y = this->fMinY + i * this->fStepY;
411
412 for (UInt_t j = 1; j < w - 1; ++j) {
413 const CellType_t &left = s->fCells[(i - 1) * (w - 1) + j];
414 const CellType_t &right = s->fCells[i * (w - 1) + j - 1];
415 CellType_t &cell = s->fCells[i * (w - 1) + j];
416 cell.fType = 0;
417 //Take values 0, 1, 4, 5 from left cube.
418 cell.fVals[1] = left.fVals[2];
419 cell.fVals[0] = left.fVals[3];
420 cell.fVals[5] = left.fVals[6];
421 cell.fVals[4] = left.fVals[7];
422 cell.fType |= (left.fType & k2_6) >> 1;
423 cell.fType |= (left.fType & k3_7) >> 3;
424 //3, 7 from right cube.
425 cell.fVals[3] = right.fVals[2];
426 cell.fVals[7] = right.fVals[6];
427 cell.fType |= (right.fType & k2_6) << 1;
428 //Calculate values 2, 6.
429 if((cell.fVals[2] = GetData(j + 1, i + 1, 0)) <= fIso)
430 cell.fType |= k2;
431 if((cell.fVals[6] = GetData(j + 1, i + 1, 1)) <= fIso)
432 cell.fType |= k6;
433
434 const UInt_t edges = eInt[cell.fType];
435 if(!edges)
436 continue;
437 //Take edges 0, 4, 8, 9 from the "left" cube.
438 //In left cube their indices are 2, 6, 11, 10.
439 if(edges & k0)
440 cell.fIds[0] = left.fIds[2];
441 if(edges & k4)
442 cell.fIds[4] = left.fIds[6];
443 if(edges & k8)
444 cell.fIds[8] = left.fIds[11];
445 if(edges & k9)
446 cell.fIds[9] = left.fIds[10];
447 //Take edges 3, 7, 11 from the "right" cube.
448 //Their "right" indices are 1, 5, 10.
449 if(edges & k3)
450 cell.fIds[3] = right.fIds[1];
451 if(edges & k7)
452 cell.fIds[7] = right.fIds[5];
453 if(edges & k11)
454 cell.fIds[11] = right.fIds[10];
455 //Calculate the remaining intersections: edges
456 //1, 2, 5, 6, 10.
457 const V x = this->fMinX + j * this->fStepX;
458 if (edges & k1)
459 SplitEdge(cell, fMesh, 1, x, y, this->fMinZ, fIso);
460 if (edges & k2)
461 SplitEdge(cell, fMesh, 2, x, y, this->fMinZ, fIso);
462 if (edges & k5)
463 SplitEdge(cell, fMesh, 5, x, y, this->fMinZ, fIso);
464 if (edges & k6)
465 SplitEdge(cell, fMesh, 6, x, y, this->fMinZ, fIso);
466 if (edges & k10)
467 SplitEdge(cell, fMesh, 10, x, y, this->fMinZ, fIso);
468
469 ConnectTriangles(cell, fMesh, fEpsilon);
470 }
471 }
472}
473
474////////////////////////////////////////////////////////////////////////////////
475///The first cube in a slice with nz == depth.
476///Neighbour is the first cube in the previous slice.
477///Four values and four edges come from the previous cube.
478
479template<class D, class V>
481 SliceType_t *slice)const
482{
483 const CellType_t &prevCell = prevSlice->fCells[0];
484 CellType_t &cell = slice->fCells[0];
485 cell.fType = 0;
486 //Values 0, 1, 2, 3 are 4, 5, 6, 7
487 //in the previous cube.
488 cell.fVals[0] = prevCell.fVals[4];
489 cell.fVals[1] = prevCell.fVals[5];
490 cell.fVals[2] = prevCell.fVals[6];
491 cell.fVals[3] = prevCell.fVals[7];
492 cell.fType |= (prevCell.fType & k4_5_6_7) >> 4;
493 //Calculate 4, 5, 6, 7.
494 if((cell.fVals[4] = GetData(0, 0, depth + 1)) <= fIso)
495 cell.fType |= k4;
496 if((cell.fVals[5] = GetData(1, 0, depth + 1)) <= fIso)
497 cell.fType |= k5;
498 if((cell.fVals[6] = GetData(1, 1, depth + 1)) <= fIso)
499 cell.fType |= k6;
500 if((cell.fVals[7] = GetData(0, 1, depth + 1)) <= fIso)
501 cell.fType |= k7;
502
503 const UInt_t edges = eInt[cell.fType];
504 if(!edges)
505 return;
506
507 //Edges 0, 1, 2, 3 taken from the prev. cube -
508 //they have indices 4, 5, 6, 7 there.
509 if(edges & k0)
510 cell.fIds[0] = prevCell.fIds[4];
511 if(edges & k1)
512 cell.fIds[1] = prevCell.fIds[5];
513 if(edges & k2)
514 cell.fIds[2] = prevCell.fIds[6];
515 if(edges & k3)
516 cell.fIds[3] = prevCell.fIds[7];
517
518 const V z = this->fMinZ + depth * this->fStepZ;
519
520 if(edges & k4)
521 SplitEdge(cell, fMesh, 4, this->fMinX, this->fMinY, z, fIso);
522 if(edges & k5)
523 SplitEdge(cell, fMesh, 5, this->fMinX, this->fMinY, z, fIso);
524 if(edges & k6)
525 SplitEdge(cell, fMesh, 6, this->fMinX, this->fMinY, z, fIso);
526 if(edges & k7)
527 SplitEdge(cell, fMesh, 7, this->fMinX, this->fMinY, z, fIso);
528 if(edges & k8)
529 SplitEdge(cell, fMesh, 8, this->fMinX, this->fMinY, z, fIso);
530 if(edges & k9)
531 SplitEdge(cell, fMesh, 9, this->fMinX, this->fMinY, z, fIso);
532 if(edges & k10)
533 SplitEdge(cell, fMesh, 10, this->fMinX, this->fMinY, z, fIso);
534 if(edges & k11)
535 SplitEdge(cell, fMesh, 11, this->fMinX, this->fMinY, z, fIso);
536
537 ConnectTriangles(cell, fMesh, fEpsilon);
538}
539
540////////////////////////////////////////////////////////////////////////////////
541///Row with ny == 0 and nz == depth, nx : [1, W - 1].
542///Two neighbours: one from previous slice (called bottom cube here),
543///the second is the previous cube in a row.
544
545template<class D, class V>
547 SliceType_t *slice)const
548{
549 const V z = this->fMinZ + depth * this->fStepZ;
550 const UInt_t w = GetW();
551
552 for (UInt_t i = 1; i < w - 1; ++i) {
553 const CellType_t &prevCell = slice->fCells[i - 1];
554 const CellType_t &bottCell = prevSlice->fCells[i];
555 CellType_t &cell = slice->fCells[i];
556 cell.fType = 0;
557 //Value 0 is not required,
558 //only bit number 0 in fType is interesting.
559 //3, 4, 7 come from the previous box (2, 5, 6)
560 cell.fVals[3] = prevCell.fVals[2];
561 cell.fVals[4] = prevCell.fVals[5];
562 cell.fVals[7] = prevCell.fVals[6];
563 cell.fType |= (prevCell.fType & k1_5) >> 1;
564 cell.fType |= (prevCell.fType & k2_6) << 1;
565 //1, 2 can be taken from the bottom cube (5, 6).
566 cell.fVals[1] = bottCell.fVals[5];
567 cell.fVals[2] = bottCell.fVals[6];
568 cell.fType |= (bottCell.fType & k5_6) >> 4;
569 //5, 6 must be calculated.
570 if((cell.fVals[5] = GetData(i + 1, 0, depth + 1)) <= fIso)
571 cell.fType |= k5;
572 if((cell.fVals[6] = GetData(i + 1, 1, depth + 1)) <= fIso)
573 cell.fType |= k6;
574
575 UInt_t edges = eInt[cell.fType];
576
577 if(!edges)
578 continue;
579 //Take edges 3, 7, 8, 11 from the previous cube (1, 5, 9, 10).
580 if(edges & k3)
581 cell.fIds[3] = prevCell.fIds[1];
582 if(edges & k7)
583 cell.fIds[7] = prevCell.fIds[5];
584 if(edges & k8)
585 cell.fIds[8] = prevCell.fIds[9];
586 if(edges & k11)
587 cell.fIds[11] = prevCell.fIds[10];
588 //Take edges 0, 1, 2 from the bottom cube (4, 5, 6).
589 if(edges & k0)
590 cell.fIds[0] = bottCell.fIds[4];
591 if(edges & k1)
592 cell.fIds[1] = bottCell.fIds[5];
593 if(edges & k2)
594 cell.fIds[2] = bottCell.fIds[6];
595
596 edges &= ~k0_1_2_3_7_8_11;
597
598 if (edges) {
599 const V x = this->fMinX + i * this->fStepX;
600
601 if(edges & k4)
602 SplitEdge(cell, fMesh, 4, x, this->fMinY, z, fIso);
603 if(edges & k5)
604 SplitEdge(cell, fMesh, 5, x, this->fMinY, z, fIso);
605 if(edges & k6)
606 SplitEdge(cell, fMesh, 6, x, this->fMinY, z, fIso);
607 if(edges & k9)
608 SplitEdge(cell, fMesh, 9, x, this->fMinY, z, fIso);
609 if(edges & k10)
610 SplitEdge(cell, fMesh, 10, x, this->fMinY, z, fIso);
611 }
612
613 ConnectTriangles(cell, fMesh, fEpsilon);
614 }
615}
616
617////////////////////////////////////////////////////////////////////////////////
618///nz == depth, nx == 0, ny : [1, H - 1].
619///Two neighbours - from previous slice ("bottom" cube)
620///and previous cube in a column.
621
622template<class D, class V>
624 SliceType_t *slice)const
625{
626 const V z = this->fMinZ + depth * this->fStepZ;
627 const UInt_t w = GetW();
628 const UInt_t h = GetH();
629
630 for (UInt_t i = 1; i < h - 1; ++i) {
631 const CellType_t &left = slice->fCells[(i - 1) * (w - 1)];
632 const CellType_t &bott = prevSlice->fCells[i * (w - 1)];
633 CellType_t &cell = slice->fCells[i * (w - 1)];
634 cell.fType = 0;
635 //Value 0 is not required, only bit.
636 //Take 1, 4, 5 from left cube.
637 cell.fVals[1] = left.fVals[2];
638 cell.fVals[4] = left.fVals[7];
639 cell.fVals[5] = left.fVals[6];
640 cell.fType |= (left.fType & k2_6) >> 1;
641 cell.fType |= (left.fType & k3_7) >> 3;
642 //2, 3 from bottom.
643 cell.fVals[2] = bott.fVals[6];
644 cell.fVals[3] = bott.fVals[7];
645 cell.fType |= (bott.fType & k6_7) >> 4;
646 //Calculate 6, 7.
647 if((cell.fVals[6] = GetData(1, i + 1, depth + 1)) <= fIso)
648 cell.fType |= k6;
649 if((cell.fVals[7] = GetData(0, i + 1, depth + 1)) <= fIso)
650 cell.fType |= k7;
651
652 const UInt_t edges = eInt[cell.fType];
653 if(!edges)
654 continue;
655
656 if(edges & k0)
657 cell.fIds[0] = left.fIds[2];
658 if(edges & k4)
659 cell.fIds[4] = left.fIds[6];
660 if(edges & k8)
661 cell.fIds[8] = left.fIds[11];
662 if(edges & k9)
663 cell.fIds[9] = left.fIds[10];
664
665 if(edges & k1)
666 cell.fIds[1] = bott.fIds[5];
667 if(edges & k2)
668 cell.fIds[2] = bott.fIds[6];
669 if(edges & k3)
670 cell.fIds[3] = bott.fIds[7];
671
672 const V y = this->fMinY + i * this->fStepY;
673
674 if(edges & k5)
675 SplitEdge(cell, fMesh, 5, this->fMinX, y, z, fIso);
676 if(edges & k6)
677 SplitEdge(cell, fMesh, 6, this->fMinX, y, z, fIso);
678 if(edges & k7)
679 SplitEdge(cell, fMesh, 7, this->fMinX, y, z, fIso);
680 if(edges & k10)
681 SplitEdge(cell, fMesh, 10, this->fMinX, y, z, fIso);
682 if(edges & k11)
683 SplitEdge(cell, fMesh, 11, this->fMinX, y, z, fIso);
684
685 ConnectTriangles(cell, fMesh, fEpsilon);
686 }
687}
688
689
690////////////////////////////////////////////////////////////////////////////////
691///nz == depth, nx : [1, W - 1], ny : [1, H - 1].
692///Each cube has 3 neighbours, "bottom" cube from
693///the previous slice, "left" and "right" from the
694///current slice.
695
696template<class D, class V>
698 SliceType_t *slice)const
699{
700 const V z = this->fMinZ + depth * this->fStepZ;
701 const UInt_t h = GetH();
702 const UInt_t w = GetW();
703
704 for (UInt_t i = 1; i < h - 1; ++i) {
705 const V y = this->fMinY + i * this->fStepY;
706 for (UInt_t j = 1; j < w - 1; ++j) {
707 const CellType_t &left = slice->fCells[(i - 1) * (w - 1) + j];
708 const CellType_t &right = slice->fCells[i * (w - 1) + j - 1];
709 const CellType_t &bott = prevSlice->fCells[i * (w - 1) + j];
710 CellType_t &cell = slice->fCells[i * (w - 1) + j];
711 cell.fType = 0;
712
713 cell.fVals[1] = left.fVals[2];
714 cell.fVals[4] = left.fVals[7];
715 cell.fVals[5] = left.fVals[6];
716 cell.fType |= (left.fType & k2_6) >> 1;
717 cell.fType |= (left.fType & k3_7) >> 3;
718
719 cell.fVals[2] = bott.fVals[6];
720 cell.fVals[3] = bott.fVals[7];
721 cell.fType |= (bott.fType & k6_7) >> 4;
722
723 cell.fVals[7] = right.fVals[6];
724 cell.fType |= (right.fType & k6) << 1;
725
726 if ((cell.fVals[6] = GetData(j + 1, i + 1, depth + 1)) <= fIso)
727 cell.fType |= k6;
728
729 const UInt_t edges = eInt[cell.fType];
730 if (!edges)
731 continue;
732
733 if(edges & k0)
734 cell.fIds[0] = left.fIds[2];
735 if(edges & k4)
736 cell.fIds[4] = left.fIds[6];
737 if(edges & k8)
738 cell.fIds[8] = left.fIds[11];
739 if(edges & k9)
740 cell.fIds[9] = left.fIds[10];
741
742 if(edges & k3)
743 cell.fIds[3] = right.fIds[1];
744 if(edges & k7)
745 cell.fIds[7] = right.fIds[5];
746 if(edges & k11)
747 cell.fIds[11] = right.fIds[10];
748
749 if(edges & k1)
750 cell.fIds[1] = bott.fIds[5];
751 if(edges & k2)
752 cell.fIds[2] = bott.fIds[6];
753
754 const V x = this->fMinX + j * this->fStepX;
755 if(edges & k5)
756 SplitEdge(cell, fMesh, 5, x, y, z, fIso);
757 if(edges & k6)
758 SplitEdge(cell, fMesh, 6, x, y, z, fIso);
759 if(edges & k10)
760 SplitEdge(cell, fMesh, 10, x, y, z, fIso);
761
762 ConnectTriangles(cell, fMesh, fEpsilon);
763 }
764 }
765}
766
767////////////////////////////////////////////////////////////////////////////////
768///Build averaged normals using vertices and
769///trinagles.
770
771template<class D, class V>
773{
774 typedef std::vector<UInt_t>::size_type size_type;
775 const UInt_t *t;
776 V *p1, *p2, *p3;
777 V v1[3], v2[3], n[3];
778
779 fMesh->fNorms.assign(fMesh->fVerts.size(), V());
780
781 for (size_type i = 0, e = fMesh->fTris.size() / 3; i < e; ++i) {
782 t = &fMesh->fTris[i * 3];
783 p1 = &fMesh->fVerts[t[0] * 3];
784 p2 = &fMesh->fVerts[t[1] * 3];
785 p3 = &fMesh->fVerts[t[2] * 3];
786 v1[0] = p2[0] - p1[0];
787 v1[1] = p2[1] - p1[1];
788 v1[2] = p2[2] - p1[2];
789 v2[0] = p3[0] - p1[0];
790 v2[1] = p3[1] - p1[1];
791 v2[2] = p3[2] - p1[2];
792 n[0] = v1[1] * v2[2] - v1[2] * v2[1];
793 n[1] = v1[2] * v2[0] - v1[0] * v2[2];
794 n[2] = v1[0] * v2[1] - v1[1] * v2[0];
795
796 const V len = std::sqrt(n[0] * n[0] + n[1] * n[1] + n[2] * n[2]);
797
798 if (len < fEpsilon)//degenerated triangle
799 continue;
800
801 n[0] /= len;
802 n[1] /= len;
803 n[2] /= len;
804 UInt_t ind = t[0] * 3;
805 fMesh->fNorms[ind] += n[0];
806 fMesh->fNorms[ind + 1] += n[1];
807 fMesh->fNorms[ind + 2] += n[2];
808 ind = t[1] * 3;
809 fMesh->fNorms[ind] += n[0];
810 fMesh->fNorms[ind + 1] += n[1];
811 fMesh->fNorms[ind + 2] += n[2];
812 ind = t[2] * 3;
813 fMesh->fNorms[ind] += n[0];
814 fMesh->fNorms[ind + 1] += n[1];
815 fMesh->fNorms[ind + 2] += n[2];
816 }
817
818 for (size_type i = 0, e = fMesh->fNorms.size() / 3; i < e; ++i) {
819 V * nn = &fMesh->fNorms[i * 3];
820 const V len = std::sqrt(nn[0] * nn[0] + nn[1] * nn[1] + nn[2] * nn[2]);
821 if (len < fEpsilon)
822 continue;
823 fMesh->fNorms[i * 3] /= len;
824 fMesh->fNorms[i * 3 + 1] /= len;
825 fMesh->fNorms[i * 3 + 2] /= len;
826 }
827}
828
829/////////////////////////////////////////////////////////////////////////
830//****************************TABLES***********************************//
831/////////////////////////////////////////////////////////////////////////
832
833const UInt_t eInt[256] =
834{
835 0x000, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c,
836 0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
837 0x190, 0x099, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c,
838 0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
839 0x230, 0x339, 0x033, 0x13a, 0x636, 0x73f, 0x435, 0x53c,
840 0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
841 0x3a0, 0x2a9, 0x1a3, 0x0aa, 0x7a6, 0x6af, 0x5a5, 0x4ac,
842 0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
843 0x460, 0x569, 0x663, 0x76a, 0x066, 0x16f, 0x265, 0x36c,
844 0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
845 0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0x0ff, 0x3f5, 0x2fc,
846 0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
847 0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x055, 0x15c,
848 0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
849 0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0x0cc,
850 0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
851 0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc,
852 0x0cc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
853 0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c,
854 0x15c, 0x055, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
855 0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc,
856 0x2fc, 0x3f5, 0x0ff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
857 0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c,
858 0x36c, 0x265, 0x16f, 0x066, 0x76a, 0x663, 0x569, 0x460,
859 0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac,
860 0x4ac, 0x5a5, 0x6af, 0x7a6, 0x0aa, 0x1a3, 0x2a9, 0x3a0,
861 0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c,
862 0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x033, 0x339, 0x230,
863 0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c,
864 0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x099, 0x190,
865 0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c,
866 0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x000
867};
868
869const Float_t vOff[8][3] =
870{
871 {0.f, 0.f, 0.f}, {1.f, 0.f, 0.f}, {1.f, 1.f, 0.f},
872 {0.f, 1.f, 0.f}, {0.f, 0.f, 1.f}, {1.f, 0.f, 1.f},
873 {1.f, 1.f, 1.f}, {0.f, 1.f, 1.f}
874};
875
876const UChar_t eConn[12][2] =
877{
878 {0, 1}, {1, 2}, {2, 3}, {3, 0},
879 {4, 5}, {5, 6}, {6, 7}, {7, 4},
880 {0, 4}, {1, 5}, {2, 6}, {3, 7}
881};
882
883const Float_t eDir[12][3] =
884{
885 { 1.f, 0.f, 0.f}, {0.f, 1.f, 0.f}, {-1.f, 0.f, 0.f},
886 { 0.f, -1.f, 0.f}, {1.f, 0.f, 0.f}, { 0.f, 1.f, 0.f},
887 {-1.f, 0.f, 0.f}, {0.f, -1.f, 0.f}, { 0.f, 0.f, 1.f},
888 { 0.f, 0.f, 1.f}, {0.f, 0.f, 1.f}, { 0.f, 0.f, 1.f}
889};
890
891
892const Int_t conTbl[256][16] =
893{
894 {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
895 {0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
896 {0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
897 {1, 8, 3, 9, 8, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
898 {1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
899 {0, 8, 3, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
900 {9, 2, 10, 0, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
901 {2, 8, 3, 2, 10, 8, 10, 9, 8, -1, -1, -1, -1, -1, -1, -1},
902 {3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
903 {0, 11, 2, 8, 11, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
904 {1, 9, 0, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
905 {1, 11, 2, 1, 9, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1},
906 {3, 10, 1, 11, 10, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
907 {0, 10, 1, 0, 8, 10, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1},
908 {3, 9, 0, 3, 11, 9, 11, 10, 9, -1, -1, -1, -1, -1, -1, -1},
909 {9, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
910 {4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
911 {4, 3, 0, 7, 3, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
912 {0, 1, 9, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
913 {4, 1, 9, 4, 7, 1, 7, 3, 1, -1, -1, -1, -1, -1, -1, -1},
914 {1, 2, 10, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
915 {3, 4, 7, 3, 0, 4, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1},
916 {9, 2, 10, 9, 0, 2, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
917 {2, 10, 9, 2, 9, 7, 2, 7, 3, 7, 9, 4, -1, -1, -1, -1},
918 {8, 4, 7, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
919 {11, 4, 7, 11, 2, 4, 2, 0, 4, -1, -1, -1, -1, -1, -1, -1},
920 {9, 0, 1, 8, 4, 7, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
921 {4, 7, 11, 9, 4, 11, 9, 11, 2, 9, 2, 1, -1, -1, -1, -1},
922 {3, 10, 1, 3, 11, 10, 7, 8, 4, -1, -1, -1, -1, -1, -1, -1},
923 {1, 11, 10, 1, 4, 11, 1, 0, 4, 7, 11, 4, -1, -1, -1, -1},
924 {4, 7, 8, 9, 0, 11, 9, 11, 10, 11, 0, 3, -1, -1, -1, -1},
925 {4, 7, 11, 4, 11, 9, 9, 11, 10, -1, -1, -1, -1, -1, -1, -1},
926 {9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
927 {9, 5, 4, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
928 {0, 5, 4, 1, 5, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
929 {8, 5, 4, 8, 3, 5, 3, 1, 5, -1, -1, -1, -1, -1, -1, -1},
930 {1, 2, 10, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
931 {3, 0, 8, 1, 2, 10, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
932 {5, 2, 10, 5, 4, 2, 4, 0, 2, -1, -1, -1, -1, -1, -1, -1},
933 {2, 10, 5, 3, 2, 5, 3, 5, 4, 3, 4, 8, -1, -1, -1, -1},
934 {9, 5, 4, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
935 {0, 11, 2, 0, 8, 11, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
936 {0, 5, 4, 0, 1, 5, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
937 {2, 1, 5, 2, 5, 8, 2, 8, 11, 4, 8, 5, -1, -1, -1, -1},
938 {10, 3, 11, 10, 1, 3, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1},
939 {4, 9, 5, 0, 8, 1, 8, 10, 1, 8, 11, 10, -1, -1, -1, -1},
940 {5, 4, 0, 5, 0, 11, 5, 11, 10, 11, 0, 3, -1, -1, -1, -1},
941 {5, 4, 8, 5, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1},
942 {9, 7, 8, 5, 7, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
943 {9, 3, 0, 9, 5, 3, 5, 7, 3, -1, -1, -1, -1, -1, -1, -1},
944 {0, 7, 8, 0, 1, 7, 1, 5, 7, -1, -1, -1, -1, -1, -1, -1},
945 {1, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
946 {9, 7, 8, 9, 5, 7, 10, 1, 2, -1, -1, -1, -1, -1, -1, -1},
947 {10, 1, 2, 9, 5, 0, 5, 3, 0, 5, 7, 3, -1, -1, -1, -1},
948 {8, 0, 2, 8, 2, 5, 8, 5, 7, 10, 5, 2, -1, -1, -1, -1},
949 {2, 10, 5, 2, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1},
950 {7, 9, 5, 7, 8, 9, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1},
951 {9, 5, 7, 9, 7, 2, 9, 2, 0, 2, 7, 11, -1, -1, -1, -1},
952 {2, 3, 11, 0, 1, 8, 1, 7, 8, 1, 5, 7, -1, -1, -1, -1},
953 {11, 2, 1, 11, 1, 7, 7, 1, 5, -1, -1, -1, -1, -1, -1, -1},
954 {9, 5, 8, 8, 5, 7, 10, 1, 3, 10, 3, 11, -1, -1, -1, -1},
955 {5, 7, 0, 5, 0, 9, 7, 11, 0, 1, 0, 10, 11, 10, 0, -1},
956 {11, 10, 0, 11, 0, 3, 10, 5, 0, 8, 0, 7, 5, 7, 0, -1},
957 {11, 10, 5, 7, 11, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
958 {10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
959 {0, 8, 3, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
960 {9, 0, 1, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
961 {1, 8, 3, 1, 9, 8, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
962 {1, 6, 5, 2, 6, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
963 {1, 6, 5, 1, 2, 6, 3, 0, 8, -1, -1, -1, -1, -1, -1, -1},
964 {9, 6, 5, 9, 0, 6, 0, 2, 6, -1, -1, -1, -1, -1, -1, -1},
965 {5, 9, 8, 5, 8, 2, 5, 2, 6, 3, 2, 8, -1, -1, -1, -1},
966 {2, 3, 11, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
967 {11, 0, 8, 11, 2, 0, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
968 {0, 1, 9, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
969 {5, 10, 6, 1, 9, 2, 9, 11, 2, 9, 8, 11, -1, -1, -1, -1},
970 {6, 3, 11, 6, 5, 3, 5, 1, 3, -1, -1, -1, -1, -1, -1, -1},
971 {0, 8, 11, 0, 11, 5, 0, 5, 1, 5, 11, 6, -1, -1, -1, -1},
972 {3, 11, 6, 0, 3, 6, 0, 6, 5, 0, 5, 9, -1, -1, -1, -1},
973 {6, 5, 9, 6, 9, 11, 11, 9, 8, -1, -1, -1, -1, -1, -1, -1},
974 {5, 10, 6, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
975 {4, 3, 0, 4, 7, 3, 6, 5, 10, -1, -1, -1, -1, -1, -1, -1},
976 {1, 9, 0, 5, 10, 6, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
977 {10, 6, 5, 1, 9, 7, 1, 7, 3, 7, 9, 4, -1, -1, -1, -1},
978 {6, 1, 2, 6, 5, 1, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1},
979 {1, 2, 5, 5, 2, 6, 3, 0, 4, 3, 4, 7, -1, -1, -1, -1},
980 {8, 4, 7, 9, 0, 5, 0, 6, 5, 0, 2, 6, -1, -1, -1, -1},
981 {7, 3, 9, 7, 9, 4, 3, 2, 9, 5, 9, 6, 2, 6, 9, -1},
982 {3, 11, 2, 7, 8, 4, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
983 {5, 10, 6, 4, 7, 2, 4, 2, 0, 2, 7, 11, -1, -1, -1, -1},
984 {0, 1, 9, 4, 7, 8, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1},
985 {9, 2, 1, 9, 11, 2, 9, 4, 11, 7, 11, 4, 5, 10, 6, -1},
986 {8, 4, 7, 3, 11, 5, 3, 5, 1, 5, 11, 6, -1, -1, -1, -1},
987 {5, 1, 11, 5, 11, 6, 1, 0, 11, 7, 11, 4, 0, 4, 11, -1},
988 {0, 5, 9, 0, 6, 5, 0, 3, 6, 11, 6, 3, 8, 4, 7, -1},
989 {6, 5, 9, 6, 9, 11, 4, 7, 9, 7, 11, 9, -1, -1, -1, -1},
990 {10, 4, 9, 6, 4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
991 {4, 10, 6, 4, 9, 10, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1},
992 {10, 0, 1, 10, 6, 0, 6, 4, 0, -1, -1, -1, -1, -1, -1, -1},
993 {8, 3, 1, 8, 1, 6, 8, 6, 4, 6, 1, 10, -1, -1, -1, -1},
994 {1, 4, 9, 1, 2, 4, 2, 6, 4, -1, -1, -1, -1, -1, -1, -1},
995 {3, 0, 8, 1, 2, 9, 2, 4, 9, 2, 6, 4, -1, -1, -1, -1},
996 {0, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
997 {8, 3, 2, 8, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1},
998 {10, 4, 9, 10, 6, 4, 11, 2, 3, -1, -1, -1, -1, -1, -1, -1},
999 {0, 8, 2, 2, 8, 11, 4, 9, 10, 4, 10, 6, -1, -1, -1, -1},
1000 {3, 11, 2, 0, 1, 6, 0, 6, 4, 6, 1, 10, -1, -1, -1, -1},
1001 {6, 4, 1, 6, 1, 10, 4, 8, 1, 2, 1, 11, 8, 11, 1, -1},
1002 {9, 6, 4, 9, 3, 6, 9, 1, 3, 11, 6, 3, -1, -1, -1, -1},
1003 {8, 11, 1, 8, 1, 0, 11, 6, 1, 9, 1, 4, 6, 4, 1, -1},
1004 {3, 11, 6, 3, 6, 0, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1},
1005 {6, 4, 8, 11, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1006 {7, 10, 6, 7, 8, 10, 8, 9, 10, -1, -1, -1, -1, -1, -1, -1},
1007 {0, 7, 3, 0, 10, 7, 0, 9, 10, 6, 7, 10, -1, -1, -1, -1},
1008 {10, 6, 7, 1, 10, 7, 1, 7, 8, 1, 8, 0, -1, -1, -1, -1},
1009 {10, 6, 7, 10, 7, 1, 1, 7, 3, -1, -1, -1, -1, -1, -1, -1},
1010 {1, 2, 6, 1, 6, 8, 1, 8, 9, 8, 6, 7, -1, -1, -1, -1},
1011 {2, 6, 9, 2, 9, 1, 6, 7, 9, 0, 9, 3, 7, 3, 9, -1},
1012 {7, 8, 0, 7, 0, 6, 6, 0, 2, -1, -1, -1, -1, -1, -1, -1},
1013 {7, 3, 2, 6, 7, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1014 {2, 3, 11, 10, 6, 8, 10, 8, 9, 8, 6, 7, -1, -1, -1, -1},
1015 {2, 0, 7, 2, 7, 11, 0, 9, 7, 6, 7, 10, 9, 10, 7, -1},
1016 {1, 8, 0, 1, 7, 8, 1, 10, 7, 6, 7, 10, 2, 3, 11, -1},
1017 {11, 2, 1, 11, 1, 7, 10, 6, 1, 6, 7, 1, -1, -1, -1, -1},
1018 {8, 9, 6, 8, 6, 7, 9, 1, 6, 11, 6, 3, 1, 3, 6, -1},
1019 {0, 9, 1, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1020 {7, 8, 0, 7, 0, 6, 3, 11, 0, 11, 6, 0, -1, -1, -1, -1},
1021 {7, 11, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1022 {7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1023 {3, 0, 8, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1024 {0, 1, 9, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1025 {8, 1, 9, 8, 3, 1, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
1026 {10, 1, 2, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1027 {1, 2, 10, 3, 0, 8, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
1028 {2, 9, 0, 2, 10, 9, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
1029 {6, 11, 7, 2, 10, 3, 10, 8, 3, 10, 9, 8, -1, -1, -1, -1},
1030 {7, 2, 3, 6, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1031 {7, 0, 8, 7, 6, 0, 6, 2, 0, -1, -1, -1, -1, -1, -1, -1},
1032 {2, 7, 6, 2, 3, 7, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1},
1033 {1, 6, 2, 1, 8, 6, 1, 9, 8, 8, 7, 6, -1, -1, -1, -1},
1034 {10, 7, 6, 10, 1, 7, 1, 3, 7, -1, -1, -1, -1, -1, -1, -1},
1035 {10, 7, 6, 1, 7, 10, 1, 8, 7, 1, 0, 8, -1, -1, -1, -1},
1036 {0, 3, 7, 0, 7, 10, 0, 10, 9, 6, 10, 7, -1, -1, -1, -1},
1037 {7, 6, 10, 7, 10, 8, 8, 10, 9, -1, -1, -1, -1, -1, -1, -1},
1038 {6, 8, 4, 11, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1039 {3, 6, 11, 3, 0, 6, 0, 4, 6, -1, -1, -1, -1, -1, -1, -1},
1040 {8, 6, 11, 8, 4, 6, 9, 0, 1, -1, -1, -1, -1, -1, -1, -1},
1041 {9, 4, 6, 9, 6, 3, 9, 3, 1, 11, 3, 6, -1, -1, -1, -1},
1042 {6, 8, 4, 6, 11, 8, 2, 10, 1, -1, -1, -1, -1, -1, -1, -1},
1043 {1, 2, 10, 3, 0, 11, 0, 6, 11, 0, 4, 6, -1, -1, -1, -1},
1044 {4, 11, 8, 4, 6, 11, 0, 2, 9, 2, 10, 9, -1, -1, -1, -1},
1045 {10, 9, 3, 10, 3, 2, 9, 4, 3, 11, 3, 6, 4, 6, 3, -1},
1046 {8, 2, 3, 8, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1},
1047 {0, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1048 {1, 9, 0, 2, 3, 4, 2, 4, 6, 4, 3, 8, -1, -1, -1, -1},
1049 {1, 9, 4, 1, 4, 2, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1},
1050 {8, 1, 3, 8, 6, 1, 8, 4, 6, 6, 10, 1, -1, -1, -1, -1},
1051 {10, 1, 0, 10, 0, 6, 6, 0, 4, -1, -1, -1, -1, -1, -1, -1},
1052 {4, 6, 3, 4, 3, 8, 6, 10, 3, 0, 3, 9, 10, 9, 3, -1},
1053 {10, 9, 4, 6, 10, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1054 {4, 9, 5, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1055 {0, 8, 3, 4, 9, 5, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
1056 {5, 0, 1, 5, 4, 0, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
1057 {11, 7, 6, 8, 3, 4, 3, 5, 4, 3, 1, 5, -1, -1, -1, -1},
1058 {9, 5, 4, 10, 1, 2, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
1059 {6, 11, 7, 1, 2, 10, 0, 8, 3, 4, 9, 5, -1, -1, -1, -1},
1060 {7, 6, 11, 5, 4, 10, 4, 2, 10, 4, 0, 2, -1, -1, -1, -1},
1061 {3, 4, 8, 3, 5, 4, 3, 2, 5, 10, 5, 2, 11, 7, 6, -1},
1062 {7, 2, 3, 7, 6, 2, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1},
1063 {9, 5, 4, 0, 8, 6, 0, 6, 2, 6, 8, 7, -1, -1, -1, -1},
1064 {3, 6, 2, 3, 7, 6, 1, 5, 0, 5, 4, 0, -1, -1, -1, -1},
1065 {6, 2, 8, 6, 8, 7, 2, 1, 8, 4, 8, 5, 1, 5, 8, -1},
1066 {9, 5, 4, 10, 1, 6, 1, 7, 6, 1, 3, 7, -1, -1, -1, -1},
1067 {1, 6, 10, 1, 7, 6, 1, 0, 7, 8, 7, 0, 9, 5, 4, -1},
1068 {4, 0, 10, 4, 10, 5, 0, 3, 10, 6, 10, 7, 3, 7, 10, -1},
1069 {7, 6, 10, 7, 10, 8, 5, 4, 10, 4, 8, 10, -1, -1, -1, -1},
1070 {6, 9, 5, 6, 11, 9, 11, 8, 9, -1, -1, -1, -1, -1, -1, -1},
1071 {3, 6, 11, 0, 6, 3, 0, 5, 6, 0, 9, 5, -1, -1, -1, -1},
1072 {0, 11, 8, 0, 5, 11, 0, 1, 5, 5, 6, 11, -1, -1, -1, -1},
1073 {6, 11, 3, 6, 3, 5, 5, 3, 1, -1, -1, -1, -1, -1, -1, -1},
1074 {1, 2, 10, 9, 5, 11, 9, 11, 8, 11, 5, 6, -1, -1, -1, -1},
1075 {0, 11, 3, 0, 6, 11, 0, 9, 6, 5, 6, 9, 1, 2, 10, -1},
1076 {11, 8, 5, 11, 5, 6, 8, 0, 5, 10, 5, 2, 0, 2, 5, -1},
1077 {6, 11, 3, 6, 3, 5, 2, 10, 3, 10, 5, 3, -1, -1, -1, -1},
1078 {5, 8, 9, 5, 2, 8, 5, 6, 2, 3, 8, 2, -1, -1, -1, -1},
1079 {9, 5, 6, 9, 6, 0, 0, 6, 2, -1, -1, -1, -1, -1, -1, -1},
1080 {1, 5, 8, 1, 8, 0, 5, 6, 8, 3, 8, 2, 6, 2, 8, -1},
1081 {1, 5, 6, 2, 1, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1082 {1, 3, 6, 1, 6, 10, 3, 8, 6, 5, 6, 9, 8, 9, 6, -1},
1083 {10, 1, 0, 10, 0, 6, 9, 5, 0, 5, 6, 0, -1, -1, -1, -1},
1084 {0, 3, 8, 5, 6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1085 {10, 5, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1086 {11, 5, 10, 7, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1087 {11, 5, 10, 11, 7, 5, 8, 3, 0, -1, -1, -1, -1, -1, -1, -1},
1088 {5, 11, 7, 5, 10, 11, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1},
1089 {10, 7, 5, 10, 11, 7, 9, 8, 1, 8, 3, 1, -1, -1, -1, -1},
1090 {11, 1, 2, 11, 7, 1, 7, 5, 1, -1, -1, -1, -1, -1, -1, -1},
1091 {0, 8, 3, 1, 2, 7, 1, 7, 5, 7, 2, 11, -1, -1, -1, -1},
1092 {9, 7, 5, 9, 2, 7, 9, 0, 2, 2, 11, 7, -1, -1, -1, -1},
1093 {7, 5, 2, 7, 2, 11, 5, 9, 2, 3, 2, 8, 9, 8, 2, -1},
1094 {2, 5, 10, 2, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1},
1095 {8, 2, 0, 8, 5, 2, 8, 7, 5, 10, 2, 5, -1, -1, -1, -1},
1096 {9, 0, 1, 5, 10, 3, 5, 3, 7, 3, 10, 2, -1, -1, -1, -1},
1097 {9, 8, 2, 9, 2, 1, 8, 7, 2, 10, 2, 5, 7, 5, 2, -1},
1098 {1, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1099 {0, 8, 7, 0, 7, 1, 1, 7, 5, -1, -1, -1, -1, -1, -1, -1},
1100 {9, 0, 3, 9, 3, 5, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1},
1101 {9, 8, 7, 5, 9, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1102 {5, 8, 4, 5, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1},
1103 {5, 0, 4, 5, 11, 0, 5, 10, 11, 11, 3, 0, -1, -1, -1, -1},
1104 {0, 1, 9, 8, 4, 10, 8, 10, 11, 10, 4, 5, -1, -1, -1, -1},
1105 {10, 11, 4, 10, 4, 5, 11, 3, 4, 9, 4, 1, 3, 1, 4, -1},
1106 {2, 5, 1, 2, 8, 5, 2, 11, 8, 4, 5, 8, -1, -1, -1, -1},
1107 {0, 4, 11, 0, 11, 3, 4, 5, 11, 2, 11, 1, 5, 1, 11, -1},
1108 {0, 2, 5, 0, 5, 9, 2, 11, 5, 4, 5, 8, 11, 8, 5, -1},
1109 {9, 4, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1110 {2, 5, 10, 3, 5, 2, 3, 4, 5, 3, 8, 4, -1, -1, -1, -1},
1111 {5, 10, 2, 5, 2, 4, 4, 2, 0, -1, -1, -1, -1, -1, -1, -1},
1112 {3, 10, 2, 3, 5, 10, 3, 8, 5, 4, 5, 8, 0, 1, 9, -1},
1113 {5, 10, 2, 5, 2, 4, 1, 9, 2, 9, 4, 2, -1, -1, -1, -1},
1114 {8, 4, 5, 8, 5, 3, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1},
1115 {0, 4, 5, 1, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1116 {8, 4, 5, 8, 5, 3, 9, 0, 5, 0, 3, 5, -1, -1, -1, -1},
1117 {9, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1118 {4, 11, 7, 4, 9, 11, 9, 10, 11, -1, -1, -1, -1, -1, -1, -1},
1119 {0, 8, 3, 4, 9, 7, 9, 11, 7, 9, 10, 11, -1, -1, -1, -1},
1120 {1, 10, 11, 1, 11, 4, 1, 4, 0, 7, 4, 11, -1, -1, -1, -1},
1121 {3, 1, 4, 3, 4, 8, 1, 10, 4, 7, 4, 11, 10, 11, 4, -1},
1122 {4, 11, 7, 9, 11, 4, 9, 2, 11, 9, 1, 2, -1, -1, -1, -1},
1123 {9, 7, 4, 9, 11, 7, 9, 1, 11, 2, 11, 1, 0, 8, 3, -1},
1124 {11, 7, 4, 11, 4, 2, 2, 4, 0, -1, -1, -1, -1, -1, -1, -1},
1125 {11, 7, 4, 11, 4, 2, 8, 3, 4, 3, 2, 4, -1, -1, -1, -1},
1126 {2, 9, 10, 2, 7, 9, 2, 3, 7, 7, 4, 9, -1, -1, -1, -1},
1127 {9, 10, 7, 9, 7, 4, 10, 2, 7, 8, 7, 0, 2, 0, 7, -1},
1128 {3, 7, 10, 3, 10, 2, 7, 4, 10, 1, 10, 0, 4, 0, 10, -1},
1129 {1, 10, 2, 8, 7, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1130 {4, 9, 1, 4, 1, 7, 7, 1, 3, -1, -1, -1, -1, -1, -1, -1},
1131 {4, 9, 1, 4, 1, 7, 0, 8, 1, 8, 7, 1, -1, -1, -1, -1},
1132 {4, 0, 3, 7, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1133 {4, 8, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1134 {9, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1135 {3, 0, 9, 3, 9, 11, 11, 9, 10, -1, -1, -1, -1, -1, -1, -1},
1136 {0, 1, 10, 0, 10, 8, 8, 10, 11, -1, -1, -1, -1, -1, -1, -1},
1137 {3, 1, 10, 11, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1138 {1, 2, 11, 1, 11, 9, 9, 11, 8, -1, -1, -1, -1, -1, -1, -1},
1139 {3, 0, 9, 3, 9, 11, 1, 2, 9, 2, 11, 9, -1, -1, -1, -1},
1140 {0, 2, 11, 8, 0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1141 {3, 2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1142 {2, 3, 8, 2, 8, 10, 10, 8, 9, -1, -1, -1, -1, -1, -1, -1},
1143 {9, 10, 2, 0, 9, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1144 {2, 3, 8, 2, 8, 10, 0, 1, 8, 1, 10, 8, -1, -1, -1, -1},
1145 {1, 10, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1146 {1, 3, 8, 9, 1, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1147 {0, 9, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1148 {0, 3, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
1149 {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
1150};
1151
1152template class TMeshBuilder<TH3C, Float_t>;
1153template class TMeshBuilder<TH3S, Float_t>;
1154template class TMeshBuilder<TH3I, Float_t>;
1155template class TMeshBuilder<TH3F, Float_t>;
1156template class TMeshBuilder<TH3D, Float_t>;
1157template class TMeshBuilder<TF3, Double_t>;
1158//TMeshBuilder does not need any detail from TKDEFGT.
1159//TKDEFGT only helps to select correct implementation.
1160//Forward class declaration is enough for TKDEFGT.
1161template class TMeshBuilder<TKDEFGT, Float_t>;
1162
1163}//namespace Mc
1164}//namespace Rgl
SVector< double, 2 > v
Definition: Dict.h:5
#define g(i)
Definition: RSha256.hxx:105
#define h(i)
Definition: RSha256.hxx:106
#define e(i)
Definition: RSha256.hxx:103
static double p3(double t, double a, double b, double c, double d)
static double p1(double t, double a, double b)
static double p2(double t, double a, double b, double c)
int Int_t
Definition: RtypesCore.h:41
unsigned char UChar_t
Definition: RtypesCore.h:34
unsigned int UInt_t
Definition: RtypesCore.h:42
double Double_t
Definition: RtypesCore.h:55
float Float_t
Definition: RtypesCore.h:53
void Error(const char *location, const char *msgfmt,...)
double sqrt(double)
Double_t GetData(UInt_t i, UInt_t j, UInt_t k) const
void SetDataSource(const TF3 *f)
void SplitEdge(TCell< Double_t > &cell, TIsoMesh< Double_t > *mesh, UInt_t i, Double_t x, Double_t y, Double_t z, Double_t iso) const
Split the edge and find normal in a new vertex.
void AddNormal(const V *n)
Definition: TGLIsoMesh.h:50
UInt_t AddVertex(const V *v)
Definition: TGLIsoMesh.h:40
UInt_t AddTriangle(const UInt_t *t)
Definition: TGLIsoMesh.h:57
std::vector< V > fVerts
Definition: TGLIsoMesh.h:81
void BuildMesh(const DataSource *src, const TGridGeometry< ValueType > &geom, MeshType_t *mesh, ValueType iso)
Build iso-mesh using marching cubes.
void NextStep(UInt_t depth, const SliceType_t *prevSlice, SliceType_t *curr) const
Fill slice with vertices and triangles.
void BuildRow(SliceType_t *slice) const
The first row (along x) in the first slice: ny == 0, nz == 0, nx : [1, W - 1].
void BuildCol(SliceType_t *slice) const
"Col" (column) consists of cubes along y axis on the first slice (nx == 0, nz == 0).
void BuildSlice(SliceType_t *slice) const
Slice with nz == 0.
void BuildNormals() const
Build averaged normals using vertices and trinagles.
void BuildFirstCube(SliceType_t *slice) const
The first cube in a grid. nx == 0, ny == 0, nz ==0.
std::vector< TCell< V > > fCells
Int_t GetNbins() const
Definition: TAxis.h:121
TAxis * GetYaxis() const
Get y axis of the function.
Definition: TF1.cxx:2387
TAxis * GetZaxis() const
Get z axis of the function. (In case this object is a TF2 or TF3)
Definition: TF1.cxx:2398
virtual Double_t Eval(Double_t x, Double_t y=0, Double_t z=0, Double_t t=0) const
Evaluate this function.
Definition: TF1.cxx:1424
TAxis * GetXaxis() const
Get x axis of the function.
Definition: TF1.cxx:2376
A 3-Dim function with parameters.
Definition: TF3.h:28
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
const UChar_t eConn[12][2]
const Int_t conTbl[256][16]
const Float_t vOff[8][3]
V GetOffset(E val1, E val2, V iso)
const Float_t eDir[12][3]
const UInt_t eInt[256]
static constexpr double s
void swap(nlohmann::json &j1, nlohmann::json &j2) noexcept(is_nothrow_move_constructible< nlohmann::json >::value and is_nothrow_move_assignable< nlohmann::json >::value)
exchanges the values of two JSON objects
Definition: json.hpp:12929
auto * m
Definition: textangle.C:8