1 /*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkMeshQuality.cxx
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14 Copyright 2003-2008 Sandia Corporation.
15 Under the terms of Contract DE-AC04-94AL85000, there is a non-exclusive
16 license for use of this work by or on behalf of the
17 U.S. Government. Redistribution and use in source and binary forms, with
18 or without modification, are permitted provided that this Notice and any
19 statement of authorship are reproduced on all copies.
20
21 Contact: dcthomp@sandia.gov,pppebay@sandia.gov
22
23 =========================================================================*/
24 #include "vtkMeshQuality.h"
25 #include "vtkInformation.h"
26 #include "vtkInformationVector.h"
27 #include "vtkObjectFactory.h"
28 #include "vtkDataSet.h"
29 #include "vtkCellData.h"
30 #include "vtkDataArray.h"
31 #include "vtkDoubleArray.h"
32 #include "vtkCell.h"
33 #include "vtkCellTypes.h"
34 #include "vtkPoints.h"
35 #include "vtkMath.h"
36 #include "vtkTetra.h"
37 #include "vtkTriangle.h"
38
39 #include "vtk_verdict.h"
40
41 vtkStandardNewMacro(vtkMeshQuality);
42
43 typedef double (*CellQualityType)( vtkCell* );
44
45 double TetVolume( vtkCell* cell );
46
47 static const char* QualityMeasureNames[] =
48 {
49 "EdgeRatio",
50 "AspectRatio",
51 "RadiusRatio",
52 "AspectFrobenius",
53 "MedAspectFrobenius",
54 "MaxAspectFrobenius",
55 "MinAngle",
56 "CollapseRatio",
57 "MaxAngle",
58 "Condition",
59 "ScaledJacobian",
60 "Shear",
61 "RelativeSizeSquared",
62 "Shape",
63 "ShapeAndSize",
64 "Distortion",
65 "MaxEdgeRatio",
66 "Skew",
67 "Taper",
68 "Volume",
69 "Stretch",
70 "Diagonal",
71 "Dimension",
72 "Oddy",
73 "ShearAndSize",
74 "Jacobian",
75 "Warpage",
76 "AspectGamma",
77 "Area",
78 "AspectBeta"
79 };
80
81 double vtkMeshQuality::CurrentTriNormal[3];
82
PrintSelf(ostream & os,vtkIndent indent)83 void vtkMeshQuality::PrintSelf(ostream& os, vtkIndent indent )
84 {
85 const char onStr[] = "On";
86 const char offStr[] = "Off";
87
88 this->Superclass::PrintSelf( os, indent );
89
90 os << indent << "SaveCellQuality: "
91 << (this->SaveCellQuality ? onStr : offStr) << endl;
92 os << indent << "TriangleQualityMeasure: "
93 << QualityMeasureNames[this->TriangleQualityMeasure] << endl;
94 os << indent << "QuadQualityMeasure: "
95 << QualityMeasureNames[this->QuadQualityMeasure] << endl;
96 os << indent << "TetQualityMeasure: "
97 << QualityMeasureNames[this->TetQualityMeasure] << endl;
98 os << indent << "HexQualityMeasure: "
99 << QualityMeasureNames[this->HexQualityMeasure] << endl;
100 os << indent << "Volume: "
101 << (this->Volume ? onStr : offStr) << endl;
102 os << indent << "CompatibilityMode: "
103 << (this->CompatibilityMode ? onStr : offStr) << endl;
104 }
105
vtkMeshQuality()106 vtkMeshQuality::vtkMeshQuality()
107 {
108 this->SaveCellQuality = 1; // Default is On
109 this->TriangleQualityMeasure = VTK_QUALITY_ASPECT_RATIO;
110 this->QuadQualityMeasure = VTK_QUALITY_EDGE_RATIO;
111 this->TetQualityMeasure = VTK_QUALITY_ASPECT_RATIO;
112 this->HexQualityMeasure = VTK_QUALITY_MAX_ASPECT_FROBENIUS;
113 this->Volume = 0;
114 this->CompatibilityMode = 0;
115 }
116
~vtkMeshQuality()117 vtkMeshQuality::~vtkMeshQuality()
118 {
119 // Nothing yet.
120 }
121
RequestData(vtkInformation * vtkNotUsed (request),vtkInformationVector ** inputVector,vtkInformationVector * outputVector)122 int vtkMeshQuality::RequestData(
123 vtkInformation *vtkNotUsed(request),
124 vtkInformationVector **inputVector,
125 vtkInformationVector *outputVector)
126 {
127 // get the info objects
128 vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
129 vtkInformation *outInfo = outputVector->GetInformationObject(0);
130
131 // get the input and output
132 vtkDataSet *in = vtkDataSet::SafeDownCast(
133 inInfo->Get(vtkDataObject::DATA_OBJECT()));
134 vtkDataSet *out = vtkDataSet::SafeDownCast(
135 outInfo->Get(vtkDataObject::DATA_OBJECT()));
136
137 CellQualityType TriangleQuality,QuadQuality,TetQuality,HexQuality;
138 vtkDoubleArray* quality = nullptr;
139 vtkDoubleArray* volume = nullptr;
140 vtkIdType N = in->GetNumberOfCells();
141 double qtrim,qtriM,Eqtri,Eqtri2;
142 double qquam,qquaM,Eqqua,Eqqua2;
143 double qtetm,qtetM,Eqtet,Eqtet2;
144 double qhexm,qhexM,Eqhex,Eqhex2;
145 double q;
146 double V = 0.;
147 vtkIdType ntri = 0;
148 vtkIdType nqua = 0;
149 vtkIdType ntet = 0;
150 vtkIdType nhex = 0;
151 vtkCell* cell;
152 int progressNumer = 0;
153 double progressDenom = 20.;
154
155 this->CellNormals = in->GetCellData()->GetNormals();
156
157 if ( this->CellNormals )
158 v_set_tri_normal_func(reinterpret_cast<ComputeNormal>(vtkMeshQuality::GetCurrentTriangleNormal));
159 else
160 v_set_tri_normal_func( nullptr );
161
162 // Initialize the min and max values, std deviations, etc.
163 qtriM = qquaM = qtetM = qhexM = VTK_DOUBLE_MIN;
164 qtrim = qquam = qtetm = qhexm = VTK_DOUBLE_MAX;
165 Eqtri = Eqtri2 = Eqqua = Eqqua2 = Eqtet = Eqtet2 = Eqhex = Eqhex2 = 0.;
166
167 switch ( this->GetTriangleQualityMeasure() )
168 {
169 case VTK_QUALITY_AREA:
170 TriangleQuality = TriangleArea;
171 break;
172 case VTK_QUALITY_EDGE_RATIO:
173 TriangleQuality = TriangleEdgeRatio;
174 break;
175 case VTK_QUALITY_ASPECT_RATIO:
176 TriangleQuality = TriangleAspectRatio;
177 break;
178 case VTK_QUALITY_RADIUS_RATIO:
179 TriangleQuality = TriangleRadiusRatio;
180 break;
181 case VTK_QUALITY_ASPECT_FROBENIUS:
182 TriangleQuality = TriangleAspectFrobenius;
183 break;
184 case VTK_QUALITY_MIN_ANGLE:
185 TriangleQuality = TriangleMinAngle;
186 break;
187 case VTK_QUALITY_MAX_ANGLE:
188 TriangleQuality = TriangleMaxAngle;
189 break;
190 case VTK_QUALITY_CONDITION:
191 TriangleQuality = TriangleCondition;
192 break;
193 case VTK_QUALITY_SCALED_JACOBIAN:
194 TriangleQuality = TriangleScaledJacobian;
195 break;
196 case VTK_QUALITY_RELATIVE_SIZE_SQUARED:
197 TriangleQuality = TriangleRelativeSizeSquared;
198 break;
199 case VTK_QUALITY_SHAPE:
200 TriangleQuality = TriangleShape;
201 break;
202 case VTK_QUALITY_SHAPE_AND_SIZE:
203 TriangleQuality = TriangleShapeAndSize;
204 break;
205 case VTK_QUALITY_DISTORTION:
206 TriangleQuality = TriangleDistortion;
207 break;
208 default:
209 vtkWarningMacro( "Bad TriangleQualityMeasure ("
210 << this->GetTriangleQualityMeasure() << "), using RadiusRatio instead");
211 TriangleQuality = TriangleRadiusRatio;
212 break;
213 }
214
215 switch ( this->GetQuadQualityMeasure() )
216 {
217 case VTK_QUALITY_EDGE_RATIO:
218 QuadQuality = QuadEdgeRatio;
219 break;
220 case VTK_QUALITY_ASPECT_RATIO:
221 QuadQuality = QuadAspectRatio;
222 break;
223 case VTK_QUALITY_RADIUS_RATIO:
224 QuadQuality = QuadRadiusRatio;
225 break;
226 case VTK_QUALITY_MED_ASPECT_FROBENIUS:
227 QuadQuality = QuadMedAspectFrobenius;
228 break;
229 case VTK_QUALITY_MAX_ASPECT_FROBENIUS:
230 QuadQuality = QuadMaxAspectFrobenius;
231 break;
232 case VTK_QUALITY_MIN_ANGLE:
233 QuadQuality = QuadMinAngle;
234 break;
235 case VTK_QUALITY_MAX_EDGE_RATIO:
236 QuadQuality = QuadMaxEdgeRatios;
237 break;
238 case VTK_QUALITY_SKEW:
239 QuadQuality = QuadSkew;
240 break;
241 case VTK_QUALITY_TAPER:
242 QuadQuality = QuadTaper;
243 break;
244 case VTK_QUALITY_WARPAGE:
245 QuadQuality = QuadWarpage;
246 break;
247 case VTK_QUALITY_AREA:
248 QuadQuality = QuadArea;
249 break;
250 case VTK_QUALITY_STRETCH:
251 QuadQuality = QuadStretch;
252 break;
253 //case VTK_QUALITY_MIN_ANGLE:
254 case VTK_QUALITY_MAX_ANGLE:
255 QuadQuality = QuadMaxAngle;
256 break;
257 case VTK_QUALITY_ODDY:
258 QuadQuality = QuadOddy;
259 break;
260 case VTK_QUALITY_CONDITION:
261 QuadQuality = QuadCondition;
262 break;
263 case VTK_QUALITY_JACOBIAN:
264 QuadQuality = QuadJacobian;
265 break;
266 case VTK_QUALITY_SCALED_JACOBIAN:
267 QuadQuality = QuadScaledJacobian;
268 break;
269 case VTK_QUALITY_SHEAR:
270 QuadQuality = QuadShear;
271 break;
272 case VTK_QUALITY_SHAPE:
273 QuadQuality = QuadShape;
274 break;
275 case VTK_QUALITY_RELATIVE_SIZE_SQUARED:
276 QuadQuality = QuadRelativeSizeSquared;
277 break;
278 case VTK_QUALITY_SHAPE_AND_SIZE:
279 QuadQuality = QuadShapeAndSize;
280 break;
281 case VTK_QUALITY_SHEAR_AND_SIZE:
282 QuadQuality = QuadShearAndSize;
283 break;
284 case VTK_QUALITY_DISTORTION:
285 QuadQuality = QuadDistortion;
286 break;
287 default:
288 vtkWarningMacro( "Bad QuadQualityMeasure ("
289 << this->GetQuadQualityMeasure() << "), using EdgeRatio instead");
290 QuadQuality = QuadEdgeRatio;
291 break;
292 }
293
294 switch ( this->GetTetQualityMeasure() )
295 {
296 case VTK_QUALITY_EDGE_RATIO:
297 TetQuality = TetEdgeRatio;
298 break;
299 case VTK_QUALITY_ASPECT_RATIO:
300 TetQuality = TetAspectRatio;
301 break;
302 case VTK_QUALITY_RADIUS_RATIO:
303 TetQuality = TetRadiusRatio;
304 break;
305 case VTK_QUALITY_ASPECT_FROBENIUS:
306 TetQuality = TetAspectFrobenius;
307 break;
308 case VTK_QUALITY_MIN_ANGLE:
309 TetQuality = TetMinAngle;
310 break;
311 case VTK_QUALITY_COLLAPSE_RATIO:
312 TetQuality = TetCollapseRatio;
313 break;
314 case VTK_QUALITY_ASPECT_BETA:
315 TetQuality = TetAspectBeta;
316 break;
317 case VTK_QUALITY_ASPECT_GAMMA:
318 TetQuality = TetAspectGamma;
319 break;
320 case VTK_QUALITY_VOLUME:
321 TetQuality = TetVolume;
322 break;
323 case VTK_QUALITY_CONDITION:
324 TetQuality = TetCondition;
325 break;
326 case VTK_QUALITY_JACOBIAN:
327 TetQuality = TetJacobian;
328 break;
329 case VTK_QUALITY_SCALED_JACOBIAN:
330 TetQuality = TetScaledJacobian;
331 break;
332 case VTK_QUALITY_SHAPE:
333 TetQuality = TetShape;
334 break;
335 case VTK_QUALITY_RELATIVE_SIZE_SQUARED:
336 TetQuality = TetRelativeSizeSquared;
337 break;
338 case VTK_QUALITY_SHAPE_AND_SIZE:
339 TetQuality = TetShapeandSize;
340 break;
341 case VTK_QUALITY_DISTORTION:
342 TetQuality = TetDistortion;
343 break;
344 default:
345 vtkWarningMacro( "Bad TetQualityMeasure ("
346 << this->GetTetQualityMeasure() << "), using RadiusRatio instead");
347 TetQuality = TetRadiusRatio;
348 break;
349 }
350
351 switch ( this->GetHexQualityMeasure() )
352 {
353 case VTK_QUALITY_EDGE_RATIO:
354 HexQuality = HexEdgeRatio;
355 break;
356 case VTK_QUALITY_MED_ASPECT_FROBENIUS:
357 HexQuality = HexMedAspectFrobenius;
358 break;
359 case VTK_QUALITY_MAX_ASPECT_FROBENIUS:
360 HexQuality = HexMaxAspectFrobenius;
361 break;
362 case VTK_QUALITY_MAX_EDGE_RATIO:
363 HexQuality = HexMaxEdgeRatio;
364 break;
365 case VTK_QUALITY_SKEW:
366 HexQuality = HexSkew;
367 break;
368 case VTK_QUALITY_TAPER:
369 HexQuality = HexTaper;
370 break;
371 case VTK_QUALITY_VOLUME:
372 HexQuality = HexVolume;
373 break;
374 case VTK_QUALITY_STRETCH:
375 HexQuality = HexStretch;
376 break;
377 case VTK_QUALITY_DIAGONAL:
378 HexQuality = HexDiagonal;
379 break;
380 case VTK_QUALITY_DIMENSION:
381 HexQuality = HexDimension;
382 break;
383 case VTK_QUALITY_ODDY:
384 HexQuality = HexOddy;
385 break;
386 case VTK_QUALITY_CONDITION:
387 HexQuality = HexCondition;
388 break;
389 case VTK_QUALITY_JACOBIAN:
390 HexQuality = HexJacobian;
391 break;
392 case VTK_QUALITY_SCALED_JACOBIAN:
393 HexQuality = HexScaledJacobian;
394 break;
395 case VTK_QUALITY_SHEAR:
396 HexQuality = HexShear;
397 break;
398 case VTK_QUALITY_SHAPE:
399 HexQuality = HexShape;
400 break;
401 case VTK_QUALITY_RELATIVE_SIZE_SQUARED:
402 HexQuality = HexRelativeSizeSquared;
403 break;
404 case VTK_QUALITY_SHAPE_AND_SIZE:
405 HexQuality = HexShapeAndSize;
406 break;
407 case VTK_QUALITY_SHEAR_AND_SIZE:
408 HexQuality = HexShearAndSize;
409 break;
410 case VTK_QUALITY_DISTORTION:
411 HexQuality = HexDistortion;
412 break;
413 default:
414 vtkWarningMacro( "Bad HexQualityMeasure ("
415 << this->GetTetQualityMeasure() << "), using MaxAspectFrobenius instead");
416 HexQuality = HexMaxAspectFrobenius;
417 break;
418 }
419
420 out->ShallowCopy( in );
421
422 if ( this->SaveCellQuality )
423 {
424 quality = vtkDoubleArray::New();
425 if ( this->CompatibilityMode )
426 {
427 if ( this->Volume )
428 {
429 quality->SetNumberOfComponents(2);
430 }
431 else
432 {
433 quality->SetNumberOfComponents(1);
434 }
435 }
436 else
437 {
438 quality->SetNumberOfComponents(1);
439 }
440 quality->SetNumberOfTuples( N );
441 quality->SetName( "Quality" );
442 out->GetCellData()->AddArray( quality );
443 out->GetCellData()->SetActiveAttribute( "Quality", vtkDataSetAttributes::SCALARS );
444 quality->Delete();
445
446 if ( ! this->CompatibilityMode )
447 {
448 if ( this->Volume )
449 {
450 volume = vtkDoubleArray::New();
451 volume->SetNumberOfComponents(1);
452 volume->SetNumberOfTuples( N );
453 volume->SetName( "Volume" );
454 out->GetCellData()->AddArray( volume );
455 volume->Delete();
456 }
457 }
458 }
459
460 // These measures require the average area/volume for all cells of the same type in the mesh.
461 // Either use the hinted value (computed by a previous vtkMeshQuality filter) or compute it.
462 if ( this->GetTriangleQualityMeasure() == VTK_QUALITY_RELATIVE_SIZE_SQUARED ||
463 this->GetTriangleQualityMeasure() == VTK_QUALITY_SHAPE_AND_SIZE ||
464 this->GetQuadQualityMeasure() == VTK_QUALITY_RELATIVE_SIZE_SQUARED ||
465 this->GetQuadQualityMeasure() == VTK_QUALITY_SHAPE_AND_SIZE ||
466 this->GetQuadQualityMeasure() == VTK_QUALITY_SHEAR_AND_SIZE ||
467 this->GetTetQualityMeasure() == VTK_QUALITY_RELATIVE_SIZE_SQUARED ||
468 this->GetTetQualityMeasure() == VTK_QUALITY_SHAPE_AND_SIZE ||
469 this->GetHexQualityMeasure() == VTK_QUALITY_RELATIVE_SIZE_SQUARED ||
470 this->GetHexQualityMeasure() == VTK_QUALITY_SHAPE_AND_SIZE ||
471 this->GetHexQualityMeasure() == VTK_QUALITY_SHEAR_AND_SIZE )
472 {
473 vtkDataArray* triAreaHint = in->GetFieldData()->GetArray( "TriArea" );
474 vtkDataArray* quadAreaHint = in->GetFieldData()->GetArray( "QuadArea" );
475 vtkDataArray* tetVolHint = in->GetFieldData()->GetArray( "TetVolume" );
476 vtkDataArray* hexVolHint = in->GetFieldData()->GetArray( "HexVolume" );
477
478 double triAreaTuple[5];
479 double quadAreaTuple[5];
480 double tetVolTuple[5];
481 double hexVolTuple[5];
482
483 if ( triAreaHint && triAreaHint->GetNumberOfTuples() > 0 && triAreaHint->GetNumberOfComponents() == 5 &&
484 quadAreaHint && quadAreaHint->GetNumberOfTuples() > 0 && quadAreaHint->GetNumberOfComponents() == 5 &&
485 tetVolHint && tetVolHint->GetNumberOfTuples() > 0 && tetVolHint->GetNumberOfComponents() == 5 &&
486 hexVolHint && hexVolHint->GetNumberOfTuples() > 0 && hexVolHint->GetNumberOfComponents() == 5 )
487 {
488 triAreaHint->GetTuple( 0, triAreaTuple );
489 quadAreaHint->GetTuple( 0, quadAreaTuple );
490 tetVolHint->GetTuple( 0, tetVolTuple );
491 hexVolHint->GetTuple( 0, hexVolTuple );
492 v_set_tri_size( triAreaTuple[1] / triAreaTuple[4] );
493 v_set_quad_size( quadAreaTuple[1] / quadAreaTuple[4] );
494 v_set_tet_size( tetVolTuple[1] / tetVolTuple[4] );
495 v_set_hex_size( hexVolTuple[1] / hexVolTuple[4] );
496 }
497 else
498 {
499 for ( int i = 0; i < 5; ++i ) {
500 triAreaTuple[i] = 0;
501 quadAreaTuple[i] = 0;
502 tetVolTuple[i] = 0;
503 hexVolTuple[i] = 0;
504 }
505 for ( vtkIdType c = 0; c < N; ++c )
506 {
507 double a, v; // area and volume
508 cell = out->GetCell( c );
509 switch ( cell->GetCellType() )
510 {
511 case VTK_TRIANGLE:
512 a = TriangleArea( cell );
513 if ( a > triAreaTuple[2] )
514 {
515 if ( triAreaTuple[0] == triAreaTuple[2] )
516 { // min == max => min has not been set
517 triAreaTuple[0] = a;
518 }
519 triAreaTuple[2] = a;
520 }
521 else if ( a < triAreaTuple[0] )
522 {
523 triAreaTuple[0] = a;
524 }
525 triAreaTuple[1] += a;
526 triAreaTuple[3] += a * a;
527 ntri++;
528 break;
529 case VTK_QUAD:
530 a = QuadArea( cell );
531 if ( a > quadAreaTuple[2] )
532 {
533 if ( quadAreaTuple[0] == quadAreaTuple[2] )
534 { // min == max => min has not been set
535 quadAreaTuple[0] = a;
536 }
537 quadAreaTuple[2] = a;
538 }
539 else if ( a < quadAreaTuple[0] )
540 {
541 quadAreaTuple[0] = a;
542 }
543 quadAreaTuple[1] += a;
544 quadAreaTuple[3] += a * a;
545 nqua++;
546 break;
547 case VTK_TETRA:
548 v = TetVolume( cell );
549 if ( v > tetVolTuple[2] )
550 {
551 if ( tetVolTuple[0] == tetVolTuple[2] )
552 { // min == max => min has not been set
553 tetVolTuple[0] = v;
554 }
555 tetVolTuple[2] = v;
556 }
557 else if ( v < tetVolTuple[0] )
558 {
559 tetVolTuple[0] = v;
560 }
561 tetVolTuple[1] += v;
562 tetVolTuple[3] += v * v;
563 ntet++;
564 break;
565 case VTK_HEXAHEDRON:
566 v = HexVolume( cell );
567 if ( v > hexVolTuple[2] )
568 {
569 if ( hexVolTuple[0] == hexVolTuple[2] )
570 { // min == max => min has not been set
571 hexVolTuple[0] = v;
572 }
573 hexVolTuple[2] = v;
574 }
575 else if ( v < hexVolTuple[0] )
576 {
577 hexVolTuple[0] = v;
578 }
579 hexVolTuple[1] += v;
580 hexVolTuple[3] += v * v;
581 nhex++;
582 break;
583 }
584 }
585 triAreaTuple[4] = ntri;
586 quadAreaTuple[4] = nqua;
587 tetVolTuple[4] = ntet;
588 hexVolTuple[4] = nhex;
589 v_set_tri_size( triAreaTuple[1] / triAreaTuple[4] );
590 v_set_quad_size( quadAreaTuple[1] / quadAreaTuple[4] );
591 v_set_tet_size( tetVolTuple[1] / tetVolTuple[4] );
592 v_set_hex_size( hexVolTuple[1] / hexVolTuple[4] );
593 progressNumer = 20;
594 progressDenom = 40.;
595 ntri = 0;
596 nqua = 0;
597 ntet = 0;
598 nhex = 0;
599
600 // Save info as field data for downstream filters
601 triAreaHint = vtkDoubleArray::New();
602 triAreaHint->SetName( "TriArea" );
603 triAreaHint->SetNumberOfComponents( 5 );
604 triAreaHint->InsertNextTuple( triAreaTuple );
605 out->GetFieldData()->AddArray( triAreaHint );
606 triAreaHint->Delete();
607
608 quadAreaHint = vtkDoubleArray::New();
609 quadAreaHint->SetName( "QuadArea" );
610 quadAreaHint->SetNumberOfComponents( 5 );
611 quadAreaHint->InsertNextTuple( quadAreaTuple );
612 out->GetFieldData()->AddArray( quadAreaHint );
613 quadAreaHint->Delete();
614
615 tetVolHint = vtkDoubleArray::New();
616 tetVolHint->SetName( "TetVolume" );
617 tetVolHint->SetNumberOfComponents( 5 );
618 tetVolHint->InsertNextTuple( tetVolTuple );
619 out->GetFieldData()->AddArray( tetVolHint );
620 tetVolHint->Delete();
621
622 hexVolHint = vtkDoubleArray::New();
623 hexVolHint->SetName( "HexVolume" );
624 hexVolHint->SetNumberOfComponents( 5 );
625 hexVolHint->InsertNextTuple( hexVolTuple );
626 out->GetFieldData()->AddArray( hexVolHint );
627 hexVolHint->Delete();
628 }
629 }
630
631 int p;
632 vtkIdType c = 0;
633 vtkIdType sz = N / 20 + 1;
634 vtkIdType inner;
635 this->UpdateProgress( progressNumer/progressDenom + 0.01 );
636 for ( p = 0; p < 20; ++p )
637 {
638 for ( inner = 0; (inner < sz && c < N); ++c, ++inner )
639 {
640 cell = out->GetCell( c );
641 V = 0.;
642 switch ( cell->GetCellType() )
643 {
644 case VTK_TRIANGLE:
645 if ( this->CellNormals )
646 this->CellNormals->GetTuple( c, vtkMeshQuality::CurrentTriNormal );
647 q = TriangleQuality( cell );
648 if ( q > qtriM )
649 {
650 if ( qtrim > qtriM )
651 {
652 qtrim = q;
653 }
654 qtriM = q;
655 }
656 else if ( q < qtrim )
657 {
658 qtrim = q;
659 }
660 Eqtri += q;
661 Eqtri2 += q * q;
662 ++ ntri;
663 break;
664 case VTK_QUAD:
665 q = QuadQuality( cell );
666 if ( q > qquaM )
667 {
668 if ( qquam > qquaM )
669 {
670 qquam = q;
671 }
672 qquaM = q;
673 }
674 else if ( q < qquam )
675 {
676 qquam = q;
677 }
678 Eqqua += q;
679 Eqqua2 += q * q;
680 ++ nqua;
681 break;
682 case VTK_TETRA:
683 q = TetQuality( cell );
684 if ( q > qtetM )
685 {
686 if ( qtetm > qtetM )
687 {
688 qtetm = q;
689 }
690 qtetM = q;
691 }
692 else if ( q < qtetm )
693 {
694 qtetm = q;
695 }
696 Eqtet += q;
697 Eqtet2 += q * q;
698 ++ ntet;
699 if ( this->Volume )
700 {
701 V = TetVolume( cell );
702 if ( ! this->CompatibilityMode )
703 {
704 volume->SetTuple1( 0, V );
705 }
706 }
707 break;
708 case VTK_HEXAHEDRON:
709 q = HexQuality( cell );
710 if ( q > qhexM )
711 {
712 if ( qhexm > qhexM )
713 {
714 qhexm = q;
715 }
716 qhexM = q;
717 }
718 else if ( q < qhexm )
719 {
720 qhexm = q;
721 }
722 Eqhex += q;
723 Eqhex2 += q * q;
724 ++ nhex;
725 break;
726 default:
727 q = 0.;
728 }
729
730 if ( this->SaveCellQuality )
731 {
732 if ( this->CompatibilityMode && this->Volume )
733 {
734 quality->SetTuple2( c, V, q );
735 }
736 else
737 {
738 quality->SetTuple1( c, q );
739 }
740 }
741 }
742 this->UpdateProgress( double(p+1+progressNumer)/progressDenom );
743 }
744
745 if ( ntri )
746 {
747 Eqtri /= static_cast<double>( ntri );
748 double multFactor = 1. / static_cast<double>( ntri > 1 ? ntri - 1 : ntri );
749 Eqtri2 = multFactor * ( Eqtri2 - static_cast<double>( ntri ) * Eqtri * Eqtri );
750 }
751 else
752 {
753 qtrim = Eqtri = qtriM = Eqtri2 = 0.;
754 }
755
756 if ( nqua )
757 {
758 Eqqua /= static_cast<double>( nqua );
759 double multFactor = 1. / static_cast<double>( nqua > 1 ? nqua - 1 : nqua );
760 Eqqua2 = multFactor * ( Eqqua2 - static_cast<double>( nqua ) * Eqqua * Eqqua );
761 }
762 else
763 {
764 qquam = Eqqua = qquaM = Eqqua2 = 0.;
765 }
766
767 if ( ntet )
768 {
769 Eqtet /= static_cast<double>( ntet );
770 double multFactor = 1. / static_cast<double>( ntet > 1 ? ntet - 1 : ntet );
771 Eqtet2 = multFactor * ( Eqtet2 - static_cast<double>( ntet ) * Eqtet * Eqtet );
772 }
773 else
774 {
775 qtetm = Eqtet = qtetM = Eqtet2 = 0.;
776 }
777
778 if ( nhex )
779 {
780 Eqhex /= static_cast<double>( nhex );
781 double multFactor = 1. / static_cast<double>( nhex > 1 ? nhex - 1 : nhex );
782 Eqhex2 = multFactor * ( Eqhex2 - static_cast<double>( nhex ) * Eqhex * Eqhex );
783 }
784 else
785 {
786 qhexm = Eqhex = qhexM = Eqhex2 = 0.;
787 }
788
789 double tuple[5];
790 quality = vtkDoubleArray::New();
791 quality->SetName( "Mesh Triangle Quality" );
792 quality->SetNumberOfComponents(5);
793 tuple[0] = qtrim;
794 tuple[1] = Eqtri;
795 tuple[2] = qtriM;
796 tuple[3] = Eqtri2;
797 tuple[4] = ntri;
798 quality->InsertNextTuple( tuple );
799 out->GetFieldData()->AddArray( quality );
800 quality->Delete();
801
802 quality = vtkDoubleArray::New();
803 quality->SetName( "Mesh Quadrilateral Quality" );
804 quality->SetNumberOfComponents(5);
805 tuple[0] = qquam;
806 tuple[1] = Eqqua;
807 tuple[2] = qquaM;
808 tuple[3] = Eqqua2;
809 tuple[4] = nqua;
810 quality->InsertNextTuple( tuple );
811 out->GetFieldData()->AddArray( quality );
812 quality->Delete();
813
814 quality = vtkDoubleArray::New();
815 quality->SetName( "Mesh Tetrahedron Quality" );
816 quality->SetNumberOfComponents(5);
817 tuple[0] = qtetm;
818 tuple[1] = Eqtet;
819 tuple[2] = qtetM;
820 tuple[3] = Eqtet2;
821 tuple[4] = ntet;
822 quality->InsertNextTuple( tuple );
823 out->GetFieldData()->AddArray( quality );
824 quality->Delete();
825
826 quality = vtkDoubleArray::New();
827 quality->SetName( "Mesh Hexahedron Quality" );
828 quality->SetNumberOfComponents(5);
829 tuple[0] = qhexm;
830 tuple[1] = Eqhex;
831 tuple[2] = qhexM;
832 tuple[3] = Eqhex2;
833 tuple[4] = nhex;
834 quality->InsertNextTuple( tuple );
835 out->GetFieldData()->AddArray( quality );
836 quality->Delete();
837
838 return 1;
839 }
840
GetCurrentTriangleNormal(double point[3],double normal[3])841 int vtkMeshQuality::GetCurrentTriangleNormal( double point[3], double normal[3] )
842 {
843 // ignore the location where the normal should be evaluated.
844 (void) point;
845
846 // copy the cell normal
847 for ( int i = 0; i < 3; ++i )
848 normal[i] = vtkMeshQuality::CurrentTriNormal[i];
849 return 1;
850 }
851
852 // Triangle quality metrics
853
TriangleArea(vtkCell * cell)854 double vtkMeshQuality::TriangleArea( vtkCell* cell )
855 {
856 double pc[3][3];
857
858 vtkPoints *p = cell->GetPoints();
859 p->GetPoint(0, pc[0]);
860 p->GetPoint(1, pc[1]);
861 p->GetPoint(2, pc[2]);
862
863 return v_tri_area( 3, pc );
864 }
865
TriangleEdgeRatio(vtkCell * cell)866 double vtkMeshQuality::TriangleEdgeRatio( vtkCell* cell )
867 {
868 double pc[3][3];
869
870 vtkPoints *p = cell->GetPoints();
871 p->GetPoint(0, pc[0]);
872 p->GetPoint(1, pc[1]);
873 p->GetPoint(2, pc[2]);
874
875 return v_tri_edge_ratio( 3, pc );
876 }
877
TriangleAspectRatio(vtkCell * cell)878 double vtkMeshQuality::TriangleAspectRatio( vtkCell* cell )
879 {
880 double pc[3][3];
881
882 vtkPoints *p = cell->GetPoints();
883 p->GetPoint(0, pc[0]);
884 p->GetPoint(1, pc[1]);
885 p->GetPoint(2, pc[2]);
886
887 return v_tri_aspect_ratio( 3, pc );
888 }
889
TriangleRadiusRatio(vtkCell * cell)890 double vtkMeshQuality::TriangleRadiusRatio( vtkCell* cell )
891 {
892 double pc[3][3];
893
894 vtkPoints *p = cell->GetPoints();
895 p->GetPoint(0, pc[0]);
896 p->GetPoint(1, pc[1]);
897 p->GetPoint(2, pc[2]);
898
899 return v_tri_radius_ratio( 3, pc );
900 }
901
TriangleAspectFrobenius(vtkCell * cell)902 double vtkMeshQuality::TriangleAspectFrobenius( vtkCell* cell )
903 {
904 double pc[3][3];
905
906 vtkPoints *p = cell->GetPoints();
907 p->GetPoint(0, pc[0]);
908 p->GetPoint(1, pc[1]);
909 p->GetPoint(2, pc[2]);
910
911 return v_tri_aspect_frobenius( 3, pc );
912 }
913
TriangleMinAngle(vtkCell * cell)914 double vtkMeshQuality::TriangleMinAngle( vtkCell* cell )
915 {
916 double pc[3][3];
917
918 vtkPoints *p = cell->GetPoints();
919 p->GetPoint(0, pc[0]);
920 p->GetPoint(1, pc[1]);
921 p->GetPoint(2, pc[2]);
922
923 return v_tri_minimum_angle( 3, pc );
924 }
925
TriangleMaxAngle(vtkCell * cell)926 double vtkMeshQuality::TriangleMaxAngle( vtkCell* cell )
927 {
928 double pc[3][3];
929
930 vtkPoints *p = cell->GetPoints();
931 p->GetPoint(0, pc[0]);
932 p->GetPoint(1, pc[1]);
933 p->GetPoint(2, pc[2]);
934
935 return v_tri_maximum_angle( 3, pc );
936 }
937
TriangleCondition(vtkCell * cell)938 double vtkMeshQuality::TriangleCondition( vtkCell* cell )
939 {
940 double pc[3][3];
941
942 vtkPoints *p = cell->GetPoints();
943 p->GetPoint(0, pc[0]);
944 p->GetPoint(1, pc[1]);
945 p->GetPoint(2, pc[2]);
946
947 return v_tri_condition( 3, pc );
948 }
949
TriangleScaledJacobian(vtkCell * cell)950 double vtkMeshQuality::TriangleScaledJacobian( vtkCell* cell )
951 {
952 double pc[3][3];
953
954 vtkPoints *p = cell->GetPoints();
955 p->GetPoint(0, pc[0]);
956 p->GetPoint(1, pc[1]);
957 p->GetPoint(2, pc[2]);
958
959 return v_tri_scaled_jacobian( 3, pc );
960 }
961
TriangleRelativeSizeSquared(vtkCell * cell)962 double vtkMeshQuality::TriangleRelativeSizeSquared( vtkCell* cell )
963 {
964 double pc[3][3];
965
966 vtkPoints *p = cell->GetPoints();
967 p->GetPoint(0, pc[0]);
968 p->GetPoint(1, pc[1]);
969 p->GetPoint(2, pc[2]);
970
971 return v_tri_relative_size_squared( 3, pc );
972 }
973
TriangleShape(vtkCell * cell)974 double vtkMeshQuality::TriangleShape( vtkCell* cell )
975 {
976 double pc[3][3];
977
978 vtkPoints *p = cell->GetPoints();
979 p->GetPoint(0, pc[0]);
980 p->GetPoint(1, pc[1]);
981 p->GetPoint(2, pc[2]);
982
983 return v_tri_shape( 3, pc );
984 }
985
TriangleShapeAndSize(vtkCell * cell)986 double vtkMeshQuality::TriangleShapeAndSize( vtkCell* cell )
987 {
988 double pc[3][3];
989
990 vtkPoints *p = cell->GetPoints();
991 p->GetPoint(0, pc[0]);
992 p->GetPoint(1, pc[1]);
993 p->GetPoint(2, pc[2]);
994
995 return v_tri_shape_and_size( 3, pc );
996 }
997
TriangleDistortion(vtkCell * cell)998 double vtkMeshQuality::TriangleDistortion( vtkCell* cell )
999 {
1000 double pc[3][3];
1001
1002 vtkPoints *p = cell->GetPoints();
1003 p->GetPoint(0, pc[0]);
1004 p->GetPoint(1, pc[1]);
1005 p->GetPoint(2, pc[2]);
1006
1007 return v_tri_distortion( 3, pc );
1008 }
1009
1010 // Quadrangle quality metrics
1011
QuadEdgeRatio(vtkCell * cell)1012 double vtkMeshQuality::QuadEdgeRatio( vtkCell* cell )
1013 {
1014 double pc[4][3];
1015
1016 vtkPoints *p = cell->GetPoints();
1017 for ( int i = 0; i < 4; ++i )
1018 p->GetPoint( i, pc[i] );
1019
1020 return v_quad_edge_ratio( 4, pc );
1021 }
1022
QuadAspectRatio(vtkCell * cell)1023 double vtkMeshQuality::QuadAspectRatio( vtkCell* cell )
1024 {
1025 double pc[4][3];
1026
1027 vtkPoints *p = cell->GetPoints();
1028 for ( int i = 0; i < 4; ++i )
1029 p->GetPoint( i, pc[i] );
1030
1031 return v_quad_aspect_ratio( 4, pc );
1032 }
1033
QuadRadiusRatio(vtkCell * cell)1034 double vtkMeshQuality::QuadRadiusRatio( vtkCell* cell )
1035 {
1036 double pc[4][3];
1037
1038 vtkPoints *p = cell->GetPoints();
1039 for ( int i = 0; i < 4; ++i )
1040 p->GetPoint( i, pc[i] );
1041
1042 return v_quad_radius_ratio( 4, pc );
1043 }
1044
QuadMedAspectFrobenius(vtkCell * cell)1045 double vtkMeshQuality::QuadMedAspectFrobenius( vtkCell* cell )
1046 {
1047 double pc[4][3];
1048
1049 vtkPoints *p = cell->GetPoints();
1050 for ( int i = 0; i < 4; ++i )
1051 p->GetPoint( i, pc[i] );
1052
1053 return v_quad_med_aspect_frobenius( 4, pc );
1054 }
1055
QuadMaxAspectFrobenius(vtkCell * cell)1056 double vtkMeshQuality::QuadMaxAspectFrobenius( vtkCell* cell )
1057 {
1058 double pc[4][3];
1059
1060 vtkPoints *p = cell->GetPoints();
1061 for ( int i = 0; i < 4; ++i )
1062 p->GetPoint( i, pc[i] );
1063
1064 return v_quad_max_aspect_frobenius( 4, pc );
1065 }
1066
QuadMinAngle(vtkCell * cell)1067 double vtkMeshQuality::QuadMinAngle( vtkCell* cell )
1068 {
1069 double pc[4][3];
1070
1071 vtkPoints *p = cell->GetPoints();
1072 for ( int i = 0; i < 4; ++i )
1073 p->GetPoint( i, pc[i] );
1074
1075 return v_quad_minimum_angle( 4, pc );
1076 }
1077
QuadMaxEdgeRatios(vtkCell * cell)1078 double vtkMeshQuality::QuadMaxEdgeRatios( vtkCell* cell )
1079 {
1080 double pc[4][3];
1081
1082 vtkPoints *p = cell->GetPoints();
1083 for ( int i = 0; i < 4; ++i )
1084 p->GetPoint( i, pc[i] );
1085
1086 return v_quad_max_edge_ratio( 4, pc );
1087 }
1088
QuadSkew(vtkCell * cell)1089 double vtkMeshQuality::QuadSkew( vtkCell* cell )
1090 {
1091 double pc[4][3];
1092
1093 vtkPoints *p = cell->GetPoints();
1094 for ( int i = 0; i < 4; ++i )
1095 p->GetPoint( i, pc[i] );
1096
1097 return v_quad_skew( 4, pc );
1098 }
1099
QuadTaper(vtkCell * cell)1100 double vtkMeshQuality::QuadTaper( vtkCell* cell )
1101 {
1102 double pc[4][3];
1103
1104 vtkPoints *p = cell->GetPoints();
1105 for ( int i = 0; i < 4; ++i )
1106 p->GetPoint( i, pc[i] );
1107
1108 return v_quad_taper( 4, pc );
1109 }
1110
QuadWarpage(vtkCell * cell)1111 double vtkMeshQuality::QuadWarpage( vtkCell* cell )
1112 {
1113 double pc[4][3];
1114
1115 vtkPoints *p = cell->GetPoints();
1116 for ( int i = 0; i < 4; ++i )
1117 p->GetPoint( i, pc[i] );
1118
1119 return v_quad_warpage( 4, pc );
1120 }
1121
QuadArea(vtkCell * cell)1122 double vtkMeshQuality::QuadArea( vtkCell* cell )
1123 {
1124 double pc[4][3];
1125
1126 vtkPoints *p = cell->GetPoints();
1127 for ( int i = 0; i < 4; ++i )
1128 p->GetPoint( i, pc[i] );
1129
1130 return v_quad_area( 4, pc );
1131 }
1132
QuadStretch(vtkCell * cell)1133 double vtkMeshQuality::QuadStretch( vtkCell* cell )
1134 {
1135 double pc[4][3];
1136
1137 vtkPoints *p = cell->GetPoints();
1138 for ( int i = 0; i < 4; ++i )
1139 p->GetPoint( i, pc[i] );
1140
1141 return v_quad_stretch( 4, pc );
1142 }
1143
1144 #if 0
1145 // FIXME
1146 double vtkMeshQuality::QuadMinAngle( vtkCell* cell )
1147 {
1148 double pc[4][3];
1149
1150 vtkPoints *p = cell->GetPoints();
1151 for ( int i = 0; i < 4; ++i )
1152 p->GetPoint( i, pc[i] );
1153
1154 return v_quad_minimum_angle( 4, pc );
1155 }
1156 #endif // 0
1157
QuadMaxAngle(vtkCell * cell)1158 double vtkMeshQuality::QuadMaxAngle( vtkCell* cell )
1159 {
1160 double pc[4][3];
1161
1162 vtkPoints *p = cell->GetPoints();
1163 for ( int i = 0; i < 4; ++i )
1164 p->GetPoint( i, pc[i] );
1165
1166 return v_quad_maximum_angle( 4, pc );
1167 }
1168
QuadOddy(vtkCell * cell)1169 double vtkMeshQuality::QuadOddy( vtkCell* cell )
1170 {
1171 double pc[4][3];
1172
1173 vtkPoints *p = cell->GetPoints();
1174 for ( int i = 0; i < 4; ++i )
1175 p->GetPoint( i, pc[i] );
1176
1177 return v_quad_oddy( 4, pc );
1178 }
1179
QuadCondition(vtkCell * cell)1180 double vtkMeshQuality::QuadCondition( vtkCell* cell )
1181 {
1182 double pc[4][3];
1183
1184 vtkPoints *p = cell->GetPoints();
1185 for ( int i = 0; i < 4; ++i )
1186 p->GetPoint( i, pc[i] );
1187
1188 return v_quad_condition( 4, pc );
1189 }
1190
QuadJacobian(vtkCell * cell)1191 double vtkMeshQuality::QuadJacobian( vtkCell* cell )
1192 {
1193 double pc[4][3];
1194
1195 vtkPoints *p = cell->GetPoints();
1196 for ( int i = 0; i < 4; ++i )
1197 p->GetPoint( i, pc[i] );
1198
1199 return v_quad_jacobian( 4, pc );
1200 }
1201
QuadScaledJacobian(vtkCell * cell)1202 double vtkMeshQuality::QuadScaledJacobian( vtkCell* cell )
1203 {
1204 double pc[4][3];
1205
1206 vtkPoints *p = cell->GetPoints();
1207 for ( int i = 0; i < 4; ++i )
1208 p->GetPoint( i, pc[i] );
1209
1210 return v_quad_scaled_jacobian( 4, pc );
1211 }
1212
QuadShear(vtkCell * cell)1213 double vtkMeshQuality::QuadShear( vtkCell* cell )
1214 {
1215 double pc[4][3];
1216
1217 vtkPoints *p = cell->GetPoints();
1218 for ( int i = 0; i < 4; ++i )
1219 p->GetPoint( i, pc[i] );
1220
1221 return v_quad_shear( 4, pc );
1222 }
1223
QuadShape(vtkCell * cell)1224 double vtkMeshQuality::QuadShape( vtkCell* cell )
1225 {
1226 double pc[4][3];
1227
1228 vtkPoints *p = cell->GetPoints();
1229 for ( int i = 0; i < 4; ++i )
1230 p->GetPoint( i, pc[i] );
1231
1232 return v_quad_shape( 4, pc );
1233 }
1234
QuadRelativeSizeSquared(vtkCell * cell)1235 double vtkMeshQuality::QuadRelativeSizeSquared( vtkCell* cell )
1236 {
1237 double pc[4][3];
1238
1239 vtkPoints *p = cell->GetPoints();
1240 for ( int i = 0; i < 4; ++i )
1241 p->GetPoint( i, pc[i] );
1242
1243 return v_quad_relative_size_squared( 4, pc );
1244 }
1245
QuadShapeAndSize(vtkCell * cell)1246 double vtkMeshQuality::QuadShapeAndSize( vtkCell* cell )
1247 {
1248 double pc[4][3];
1249
1250 vtkPoints *p = cell->GetPoints();
1251 for ( int i = 0; i < 4; ++i )
1252 p->GetPoint( i, pc[i] );
1253
1254 return v_quad_shape_and_size( 4, pc );
1255 }
1256
QuadShearAndSize(vtkCell * cell)1257 double vtkMeshQuality::QuadShearAndSize( vtkCell* cell )
1258 {
1259 double pc[4][3];
1260
1261 vtkPoints *p = cell->GetPoints();
1262 for ( int i = 0; i < 4; ++i )
1263 p->GetPoint( i, pc[i] );
1264
1265 return v_quad_shear_and_size( 4, pc );
1266 }
1267
QuadDistortion(vtkCell * cell)1268 double vtkMeshQuality::QuadDistortion( vtkCell* cell )
1269 {
1270 double pc[4][3];
1271
1272 vtkPoints *p = cell->GetPoints();
1273 for ( int i = 0; i < 4; ++i )
1274 p->GetPoint( i, pc[i] );
1275
1276 return v_quad_distortion( 4, pc );
1277 }
1278
1279
1280 // Volume of a tetrahedron, for compatibility with the original vtkMeshQuality
1281
TetVolume(vtkCell * cell)1282 double TetVolume( vtkCell* cell )
1283 {
1284 double x0[3];
1285 double x1[3];
1286 double x2[3];
1287 double x3[3];
1288 cell->Points->GetPoint( 0, x0 );
1289 cell->Points->GetPoint( 1, x1 );
1290 cell->Points->GetPoint( 2, x2 );
1291 cell->Points->GetPoint( 3, x3 );
1292 return vtkTetra::ComputeVolume( x0, x1, x2, x3 );
1293 }
1294
1295 // Tetrahedral quality metrics
1296
TetEdgeRatio(vtkCell * cell)1297 double vtkMeshQuality::TetEdgeRatio( vtkCell* cell )
1298 {
1299 double pc[4][3];
1300
1301 vtkPoints *p = cell->GetPoints();
1302 for ( int i = 0; i < 4; ++i )
1303 p->GetPoint( i, pc[i] );
1304
1305 return v_tet_edge_ratio( 4, pc );
1306 }
1307
TetAspectRatio(vtkCell * cell)1308 double vtkMeshQuality::TetAspectRatio( vtkCell* cell )
1309 {
1310 double pc[4][3];
1311
1312 vtkPoints *p = cell->GetPoints();
1313 for ( int i = 0; i < 4; ++i )
1314 p->GetPoint( i, pc[i] );
1315
1316 return v_tet_aspect_ratio( 4, pc );
1317 }
1318
TetRadiusRatio(vtkCell * cell)1319 double vtkMeshQuality::TetRadiusRatio( vtkCell* cell )
1320 {
1321 double pc[4][3];
1322
1323 vtkPoints *p = cell->GetPoints();
1324 for ( int i = 0; i < 4; ++i )
1325 p->GetPoint( i, pc[i] );
1326
1327 return v_tet_radius_ratio( 4, pc );
1328 }
1329
TetAspectBeta(vtkCell * cell)1330 double vtkMeshQuality::TetAspectBeta( vtkCell* cell )
1331 {
1332 double pc[4][3];
1333
1334 vtkPoints *p = cell->GetPoints();
1335 for ( int i = 0; i < 4; ++i )
1336 p->GetPoint( i, pc[i] );
1337
1338 return v_tet_aspect_beta( 4, pc );
1339 }
1340
TetAspectFrobenius(vtkCell * cell)1341 double vtkMeshQuality::TetAspectFrobenius( vtkCell* cell )
1342 {
1343 double pc[4][3];
1344
1345 vtkPoints *p = cell->GetPoints();
1346 for ( int i = 0; i < 4; ++i )
1347 p->GetPoint( i, pc[i] );
1348
1349 return v_tet_aspect_frobenius( 4, pc );
1350 }
1351
TetMinAngle(vtkCell * cell)1352 double vtkMeshQuality::TetMinAngle( vtkCell* cell )
1353 {
1354 double pc[4][3];
1355
1356 vtkPoints *p = cell->GetPoints();
1357 for ( int i = 0; i < 4; ++i )
1358 p->GetPoint( i, pc[i] );
1359
1360 return v_tet_minimum_angle( 4, pc );
1361 }
1362
TetCollapseRatio(vtkCell * cell)1363 double vtkMeshQuality::TetCollapseRatio( vtkCell* cell )
1364 {
1365 double pc[4][3];
1366
1367 vtkPoints *p = cell->GetPoints();
1368 for ( int i = 0; i < 4; ++i )
1369 p->GetPoint( i, pc[i] );
1370
1371 return v_tet_collapse_ratio( 4, pc );
1372 }
1373
TetAspectGamma(vtkCell * cell)1374 double vtkMeshQuality::TetAspectGamma( vtkCell* cell )
1375 {
1376 double pc[4][3];
1377
1378 vtkPoints *p = cell->GetPoints();
1379 for ( int i = 0; i < 4; ++i )
1380 p->GetPoint( i, pc[i] );
1381
1382 return v_tet_aspect_gamma( 4, pc );
1383 }
1384
TetVolume(vtkCell * cell)1385 double vtkMeshQuality::TetVolume( vtkCell* cell )
1386 {
1387 double pc[4][3];
1388
1389 vtkPoints* p = cell->GetPoints();
1390 for ( int i = 0; i < 4; ++i )
1391 p->GetPoint( i, pc[i] );
1392
1393 return v_tet_volume( 4, pc );
1394 }
1395
TetCondition(vtkCell * cell)1396 double vtkMeshQuality::TetCondition( vtkCell* cell )
1397 {
1398 double pc[4][3];
1399
1400 vtkPoints* p = cell->GetPoints();
1401 for ( int i = 0; i < 4; ++i )
1402 p->GetPoint( i, pc[i] );
1403
1404 return v_tet_condition( 4, pc );
1405 }
1406
TetJacobian(vtkCell * cell)1407 double vtkMeshQuality::TetJacobian( vtkCell* cell )
1408 {
1409 double pc[4][3];
1410
1411 vtkPoints* p = cell->GetPoints();
1412 for ( int i = 0; i < 4; ++i )
1413 p->GetPoint( i, pc[i] );
1414
1415 return v_tet_jacobian( 4, pc );
1416 }
1417
TetScaledJacobian(vtkCell * cell)1418 double vtkMeshQuality::TetScaledJacobian( vtkCell* cell )
1419 {
1420 double pc[4][3];
1421
1422 vtkPoints* p = cell->GetPoints();
1423 for ( int i = 0; i < 4; ++i )
1424 p->GetPoint( i, pc[i] );
1425
1426 return v_tet_scaled_jacobian( 4, pc );
1427 }
1428
TetShape(vtkCell * cell)1429 double vtkMeshQuality::TetShape( vtkCell* cell )
1430 {
1431 double pc[4][3];
1432
1433 vtkPoints* p = cell->GetPoints();
1434 for ( int i = 0; i < 4; ++i )
1435 p->GetPoint( i, pc[i] );
1436
1437 return v_tet_shape( 4, pc );
1438 }
1439
TetRelativeSizeSquared(vtkCell * cell)1440 double vtkMeshQuality::TetRelativeSizeSquared( vtkCell* cell )
1441 {
1442 double pc[4][3];
1443
1444 vtkPoints* p = cell->GetPoints();
1445 for ( int i = 0; i < 4; ++i )
1446 p->GetPoint( i, pc[i] );
1447
1448 return v_tet_relative_size_squared( 4, pc );
1449 }
1450
TetShapeandSize(vtkCell * cell)1451 double vtkMeshQuality::TetShapeandSize( vtkCell* cell )
1452 {
1453 double pc[4][3];
1454
1455 vtkPoints* p = cell->GetPoints();
1456 for ( int i = 0; i < 4; ++i )
1457 p->GetPoint( i, pc[i] );
1458
1459 return v_tet_shape_and_size( 4, pc );
1460 }
1461
TetDistortion(vtkCell * cell)1462 double vtkMeshQuality::TetDistortion( vtkCell* cell )
1463 {
1464 double pc[4][3];
1465
1466 vtkPoints* p = cell->GetPoints();
1467 for ( int i = 0; i < 4; ++i )
1468 p->GetPoint( i, pc[i] );
1469
1470 return v_tet_distortion( 4, pc );
1471 }
1472
1473 // Hexahedral quality metrics
1474
HexEdgeRatio(vtkCell * cell)1475 double vtkMeshQuality::HexEdgeRatio( vtkCell* cell)
1476 {
1477 double pc[8][3];
1478 vtkPoints* p = cell->GetPoints();
1479
1480 for ( int i = 0; i < 8; ++i )
1481 p->GetPoint( i, pc[i] );
1482
1483 return v_hex_edge_ratio( 8, pc );
1484 }
1485
HexMedAspectFrobenius(vtkCell * cell)1486 double vtkMeshQuality::HexMedAspectFrobenius( vtkCell* cell )
1487 {
1488 double pc[8][3];
1489
1490 vtkPoints* p = cell->GetPoints();
1491 for ( int i = 0; i < 8; ++i )
1492 p->GetPoint( i, pc[i] );
1493
1494 return v_hex_med_aspect_frobenius( 8, pc );
1495 }
1496
HexMaxAspectFrobenius(vtkCell * cell)1497 double vtkMeshQuality::HexMaxAspectFrobenius( vtkCell* cell )
1498 {
1499 double pc[8][3];
1500
1501 vtkPoints* p = cell->GetPoints();
1502 for ( int i = 0; i < 8; ++i )
1503 p->GetPoint( i, pc[i] );
1504
1505 return v_hex_max_aspect_frobenius( 8, pc );
1506 }
1507
HexMaxEdgeRatio(vtkCell * cell)1508 double vtkMeshQuality::HexMaxEdgeRatio( vtkCell* cell )
1509 {
1510 double pc[8][3];
1511
1512 vtkPoints* p = cell->GetPoints();
1513 for ( int i = 0; i < 8; ++i )
1514 p->GetPoint( i, pc[i] );
1515
1516 return v_hex_max_edge_ratio( 8, pc );
1517 }
1518
HexSkew(vtkCell * cell)1519 double vtkMeshQuality::HexSkew( vtkCell* cell )
1520 {
1521 double pc[8][3];
1522 vtkPoints* p = cell->GetPoints();
1523
1524 for ( int i = 0; i < 8; ++i )
1525 p->GetPoint( i, pc[i] );
1526
1527 return v_hex_skew( 8, pc );
1528 }
1529
HexTaper(vtkCell * cell)1530 double vtkMeshQuality::HexTaper( vtkCell* cell )
1531 {
1532 double pc[8][3];
1533 vtkPoints* p = cell->GetPoints();
1534
1535 for ( int i = 0; i < 8; ++i )
1536 p->GetPoint( i, pc[i] );
1537
1538 return v_hex_taper( 8, pc );
1539 }
1540
HexVolume(vtkCell * cell)1541 double vtkMeshQuality::HexVolume( vtkCell* cell )
1542 {
1543 double pc[8][3];
1544 vtkPoints* p = cell->GetPoints();
1545
1546 for ( int i = 0; i < 8; ++i )
1547 p->GetPoint( i, pc[i] );
1548
1549 return v_hex_volume( 8, pc );
1550 }
1551
HexStretch(vtkCell * cell)1552 double vtkMeshQuality::HexStretch( vtkCell* cell )
1553 {
1554 double pc[8][3];
1555 vtkPoints* p = cell->GetPoints();
1556
1557 for ( int i = 0; i < 8; ++i )
1558 p->GetPoint( i, pc[i] );
1559
1560 return v_hex_stretch( 8, pc );
1561 }
1562
HexDiagonal(vtkCell * cell)1563 double vtkMeshQuality::HexDiagonal( vtkCell* cell )
1564 {
1565 double pc[8][3];
1566 vtkPoints* p = cell->GetPoints();
1567
1568 for ( int i = 0; i < 8; ++i )
1569 p->GetPoint( i, pc[i] );
1570
1571 return v_hex_diagonal( 8, pc );
1572 }
1573
HexDimension(vtkCell * cell)1574 double vtkMeshQuality::HexDimension( vtkCell* cell )
1575 {
1576 double pc[8][3];
1577 vtkPoints* p = cell->GetPoints();
1578
1579 for ( int i = 0; i < 8; ++i )
1580 p->GetPoint( i, pc[i] );
1581
1582 return v_hex_dimension( 8, pc );
1583 }
1584
HexOddy(vtkCell * cell)1585 double vtkMeshQuality::HexOddy( vtkCell* cell )
1586 {
1587 double pc[8][3];
1588 vtkPoints* p = cell->GetPoints();
1589
1590 for ( int i = 0; i < 8; ++i )
1591 p->GetPoint( i, pc[i] );
1592
1593 return v_hex_oddy( 8, pc );
1594 }
1595
HexCondition(vtkCell * cell)1596 double vtkMeshQuality::HexCondition( vtkCell* cell )
1597 {
1598 double pc[8][3];
1599 vtkPoints* p = cell->GetPoints();
1600
1601 for ( int i = 0; i < 8; ++i )
1602 p->GetPoint( i, pc[i] );
1603
1604 return v_hex_condition( 8, pc );
1605 }
1606
HexJacobian(vtkCell * cell)1607 double vtkMeshQuality::HexJacobian( vtkCell* cell )
1608 {
1609 double pc[8][3];
1610 vtkPoints* p = cell->GetPoints();
1611
1612 for ( int i = 0; i < 8; ++i )
1613 p->GetPoint( i, pc[i] );
1614
1615 return v_hex_jacobian( 8, pc );
1616 }
1617
HexScaledJacobian(vtkCell * cell)1618 double vtkMeshQuality::HexScaledJacobian( vtkCell* cell )
1619 {
1620 double pc[8][3];
1621 vtkPoints* p = cell->GetPoints();
1622
1623 for ( int i = 0; i < 8; ++i )
1624 p->GetPoint( i, pc[i] );
1625
1626 return v_hex_scaled_jacobian( 8, pc );
1627 }
1628
HexShear(vtkCell * cell)1629 double vtkMeshQuality::HexShear( vtkCell* cell )
1630 {
1631 double pc[8][3];
1632 vtkPoints* p = cell->GetPoints();
1633
1634 for ( int i = 0; i < 8; ++i )
1635 p->GetPoint( i, pc[i] );
1636
1637 return v_hex_shear( 8, pc );
1638 }
1639
HexShape(vtkCell * cell)1640 double vtkMeshQuality::HexShape( vtkCell* cell )
1641 {
1642 double pc[8][3];
1643 vtkPoints* p = cell->GetPoints();
1644
1645 for ( int i = 0; i < 8; ++i )
1646 p->GetPoint( i, pc[i] );
1647
1648 return v_hex_shape( 8, pc );
1649 }
1650
HexRelativeSizeSquared(vtkCell * cell)1651 double vtkMeshQuality::HexRelativeSizeSquared( vtkCell* cell )
1652 {
1653 double pc[8][3];
1654 vtkPoints* p = cell->GetPoints();
1655
1656 for ( int i = 0; i < 8; ++i )
1657 p->GetPoint( i, pc[i] );
1658
1659 return v_hex_relative_size_squared( 8, pc );
1660 }
1661
HexShapeAndSize(vtkCell * cell)1662 double vtkMeshQuality::HexShapeAndSize( vtkCell* cell )
1663 {
1664 double pc[8][3];
1665 vtkPoints* p = cell->GetPoints();
1666
1667 for ( int i = 0; i < 8; ++i )
1668 p->GetPoint( i, pc[i] );
1669
1670 return v_hex_shape_and_size( 8, pc );
1671 }
1672
HexShearAndSize(vtkCell * cell)1673 double vtkMeshQuality::HexShearAndSize( vtkCell* cell )
1674 {
1675 double pc[8][3];
1676 vtkPoints* p = cell->GetPoints();
1677
1678 for ( int i = 0; i < 8; ++i )
1679 p->GetPoint( i, pc[i] );
1680
1681 return v_hex_shear_and_size( 8, pc );
1682 }
1683
HexDistortion(vtkCell * cell)1684 double vtkMeshQuality::HexDistortion( vtkCell* cell )
1685 {
1686 double pc[8][3];
1687 vtkPoints* p = cell->GetPoints();
1688
1689 for ( int i = 0; i < 8; ++i )
1690 p->GetPoint( i, pc[i] );
1691
1692 return v_hex_distortion( 8, pc );
1693 }
1694