source: branches/vtkvis_threaded/RpVtkDataSet.cpp @ 2549

Last change on this file since 2549 was 2478, checked in by ldelgass, 13 years ago

Add dataset bounding box outline color command. Legend - if no title/labels
are requested, try to make color bar fill the legend area (not currently pixel
exact because of coordinate system conversions and rounding in
vtkScalarBarActor)

  • Property svn:eol-style set to native
File size: 21.9 KB
Line 
1/* -*- mode: c++; c-basic-offset: 4; indent-tabs-mode: nil -*- */
2/*
3 * Copyright (C) 2011, Purdue Research Foundation
4 *
5 * Author: Leif Delgass <ldelgass@purdue.edu>
6 */
7
8#include <cassert>
9#include <cstring>
10#include <cfloat>
11#include <cmath>
12
13#include <vtkCharArray.h>
14#include <vtkDataSetReader.h>
15#include <vtkPolyData.h>
16#include <vtkStructuredPoints.h>
17#include <vtkStructuredGrid.h>
18#include <vtkRectilinearGrid.h>
19#include <vtkUnstructuredGrid.h>
20#include <vtkProperty.h>
21#include <vtkPointData.h>
22#include <vtkCellData.h>
23#include <vtkCell.h>
24#include <vtkLookupTable.h>
25
26#include "RpVtkDataSet.h"
27#include "Trace.h"
28
29using namespace Rappture::VtkVis;
30
31DataSet::DataSet(const std::string& name) :
32    _name(name),
33    _visible(true),
34    _opacity(1),
35    _cellSizeAverage(0)
36{
37    _cellSizeRange[0] = -1;
38    _cellSizeRange[1] = -1;
39}
40
41DataSet::~DataSet()
42{
43}
44
45/**
46 * \brief Create and initialize a VTK Prop to render the outline
47 */
48void DataSet::initProp()
49{
50    if (_prop == NULL) {
51        _prop = vtkSmartPointer<vtkActor>::New();
52        vtkProperty *property = _prop->GetProperty();
53        property->EdgeVisibilityOff();
54        property->SetOpacity(_opacity);
55        property->SetAmbient(.2);
56        property->LightingOff();
57        _prop->SetVisibility((_visible ? 1 : 0));
58    }
59}
60
61/**
62 * \brief Create and initialize a wireframe outline
63 */
64void DataSet::showOutline(bool state)
65{
66    if (state) {
67        if (_outlineFilter == NULL) {
68            _outlineFilter = vtkSmartPointer<vtkOutlineFilter>::New();
69            _outlineFilter->SetInput(_dataSet);
70        }
71        if (_outlineMapper == NULL) {
72            _outlineMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
73            _outlineMapper->SetInputConnection(_outlineFilter->GetOutputPort());
74        }
75        initProp();
76        _prop->SetMapper(_outlineMapper);
77    } else {
78        if (_prop != NULL) {
79            _prop->SetMapper(NULL);
80        }
81        if (_outlineMapper != NULL) {
82            _outlineMapper = NULL;
83        }
84        if (_outlineFilter != NULL) {
85            _outlineFilter = NULL;
86        }
87    }
88}
89
90/**
91 * \brief Set color of outline bounding box
92 */
93void DataSet::setOutlineColor(float color[3])
94{
95    if (_prop == NULL) {
96        initProp();
97    }
98    _prop->GetProperty()->SetColor(color[0], color[1], color[2]);
99}
100
101/**
102 * \brief Set opacity of DataSet outline
103 *
104 * This method is used for record-keeping and opacity of the
105 * DataSet bounds outline.  The renderer controls opacity
106 * of other related graphics objects.
107 */
108void DataSet::setOpacity(double opacity)
109{
110    _opacity = opacity;
111    if (_prop != NULL) {
112        _prop->GetProperty()->SetOpacity(opacity);
113    }
114}
115
116/**
117 * \brief Set visibility flag in DataSet
118 *
119 * This method is used for record-keeping and visibility of the
120 * DataSet bounds outline.  The renderer controls visibility
121 * of other related graphics objects.
122 */
123void DataSet::setVisibility(bool state)
124{
125    _visible = state;
126    if (_prop != NULL) {
127        _prop->SetVisibility((state ? 1 : 0));
128    }
129}
130
131/**
132 * \brief Get visibility flag in DataSet
133 *
134 * This method is used for record-keeping.  The renderer controls
135 * the visibility of related graphics objects.
136 */
137bool DataSet::getVisibility() const
138{
139    return _visible;
140}
141
142/**
143 * \brief Read a VTK data file
144 */
145bool DataSet::setDataFile(const char *filename)
146{
147    vtkSmartPointer<vtkDataSetReader> reader = vtkSmartPointer<vtkDataSetReader>::New();
148
149#if defined(WANT_TRACE) && defined(DEBUG)
150    reader->DebugOn();
151#endif
152
153    reader->SetFileName(filename);
154    reader->ReadAllNormalsOn();
155    //reader->ReadAllTCoordsOn();
156    reader->ReadAllScalarsOn();
157    //reader->ReadAllColorScalarsOn();
158    reader->ReadAllVectorsOn();
159    //reader->ReadAllTensorsOn();
160    reader->ReadAllFieldsOn();
161
162    return setData(reader);
163}
164
165/**
166 * \brief Read a VTK data file from a memory buffer
167 */
168bool DataSet::setData(char *data, int nbytes)
169{
170    vtkSmartPointer<vtkDataSetReader> reader = vtkSmartPointer<vtkDataSetReader>::New();
171    vtkSmartPointer<vtkCharArray> dataSetString = vtkSmartPointer<vtkCharArray>::New();
172
173#if defined(WANT_TRACE) && defined(DEBUG)
174    reader->DebugOn();
175    dataSetString->DebugOn();
176#endif
177
178    dataSetString->SetArray(data, nbytes, 1);
179    reader->SetInputArray(dataSetString);
180    reader->ReadFromInputStringOn();
181    reader->ReadAllNormalsOn();
182    //reader->ReadAllTCoordsOn();
183    reader->ReadAllScalarsOn();
184    //reader->ReadAllColorScalarsOn();
185    reader->ReadAllVectorsOn();
186    //reader->ReadAllTensorsOn();
187    reader->ReadAllFieldsOn();
188
189    return setData(reader);
190}
191
192void DataSet::print() const
193{
194    TRACE("DataSet class: %s", _dataSet->GetClassName());
195
196    TRACE("DataSet memory: %g MiB", _dataSet->GetActualMemorySize()/1024.);
197
198    double bounds[6];
199    getBounds(bounds);
200
201    // Topology
202    TRACE("DataSet bounds: %g %g %g %g %g %g",
203          bounds[0], bounds[1],
204          bounds[2], bounds[3],
205          bounds[4], bounds[5]);
206    TRACE("Points: %d Cells: %d", _dataSet->GetNumberOfPoints(), _dataSet->GetNumberOfCells());
207
208    double dataRange[2];
209    if (_dataSet->GetPointData() != NULL) {
210        TRACE("PointData arrays: %d", _dataSet->GetPointData()->GetNumberOfArrays());
211        for (int i = 0; i < _dataSet->GetPointData()->GetNumberOfArrays(); i++) {
212            _dataSet->GetPointData()->GetArray(i)->GetRange(dataRange, -1);
213            TRACE("PointData[%d]: '%s' comp:%d, (%g,%g)", i,
214                  _dataSet->GetPointData()->GetArrayName(i),
215                  _dataSet->GetPointData()->GetArray(i)->GetNumberOfComponents(),
216                  dataRange[0], dataRange[1]);
217        }
218    }
219    if (_dataSet->GetCellData() != NULL) {
220        TRACE("CellData arrays: %d", _dataSet->GetCellData()->GetNumberOfArrays());
221        for (int i = 0; i < _dataSet->GetCellData()->GetNumberOfArrays(); i++) {
222            _dataSet->GetCellData()->GetArray(i)->GetRange(dataRange, -1);
223            TRACE("CellData[%d]: '%s' comp:%d, (%g,%g)", i,
224                  _dataSet->GetCellData()->GetArrayName(i),
225                  _dataSet->GetCellData()->GetArray(i)->GetNumberOfComponents(),
226                  dataRange[0], dataRange[1]);
227        }
228    }
229    if (_dataSet->GetFieldData() != NULL) {
230        TRACE("FieldData arrays: %d", _dataSet->GetFieldData()->GetNumberOfArrays());
231        for (int i = 0; i < _dataSet->GetFieldData()->GetNumberOfArrays(); i++) {
232            _dataSet->GetFieldData()->GetArray(i)->GetRange(dataRange, -1);
233            TRACE("FieldData[%d]: '%s' comp:%d, tuples:%d (%g,%g)", i,
234                  _dataSet->GetFieldData()->GetArrayName(i),
235                  _dataSet->GetFieldData()->GetArray(i)->GetNumberOfComponents(),
236                  _dataSet->GetFieldData()->GetArray(i)->GetNumberOfTuples(),
237                  dataRange[0], dataRange[1]);
238        }
239    }
240}
241
242void DataSet::setDefaultArrays()
243{
244    if (_dataSet->GetPointData() != NULL &&
245        _dataSet->GetPointData()->GetScalars() == NULL &&
246        _dataSet->GetPointData()->GetNumberOfArrays() > 0) {
247        for (int i = 0; i < _dataSet->GetPointData()->GetNumberOfArrays(); i++) {
248            if (_dataSet->GetPointData()->GetArray(i)->GetNumberOfComponents() == 1) {
249                TRACE("Setting point scalars to '%s'", _dataSet->GetPointData()->GetArrayName(i));
250                _dataSet->GetPointData()->SetActiveScalars(_dataSet->GetPointData()->GetArrayName(i));
251                break;
252            }
253        }
254    }
255    if (_dataSet->GetPointData() != NULL &&
256        _dataSet->GetPointData()->GetVectors() == NULL &&
257        _dataSet->GetPointData()->GetNumberOfArrays() > 0) {
258        for (int i = 0; i < _dataSet->GetPointData()->GetNumberOfArrays(); i++) {
259            if (_dataSet->GetPointData()->GetArray(i)->GetNumberOfComponents() == 3) {
260                TRACE("Setting point vectors to '%s'", _dataSet->GetPointData()->GetArrayName(i));
261                _dataSet->GetPointData()->SetActiveVectors(_dataSet->GetPointData()->GetArrayName(i));
262                break;
263            }
264        }
265    }
266    if (_dataSet->GetCellData() != NULL &&
267        _dataSet->GetCellData()->GetScalars() == NULL &&
268        _dataSet->GetCellData()->GetNumberOfArrays() > 0) {
269        for (int i = 0; i < _dataSet->GetCellData()->GetNumberOfArrays(); i++) {
270            if (_dataSet->GetCellData()->GetArray(i)->GetNumberOfComponents() == 1) {
271                TRACE("Setting cell scalars to '%s'", _dataSet->GetCellData()->GetArrayName(i));
272                _dataSet->GetCellData()->SetActiveScalars(_dataSet->GetCellData()->GetArrayName(i));
273                break;
274            }
275        }
276    }
277    if (_dataSet->GetCellData() != NULL &&
278        _dataSet->GetCellData()->GetVectors() == NULL &&
279        _dataSet->GetCellData()->GetNumberOfArrays() > 0) {
280        for (int i = 0; i < _dataSet->GetCellData()->GetNumberOfArrays(); i++) {
281            if (_dataSet->GetCellData()->GetArray(i)->GetNumberOfComponents() == 3) {
282                TRACE("Setting cell vectors to '%s'", _dataSet->GetCellData()->GetArrayName(i));
283                _dataSet->GetCellData()->SetActiveVectors(_dataSet->GetCellData()->GetArrayName(i));
284                break;
285            }
286        }
287    }
288}
289
290/**
291 * \brief Read dataset using supplied reader
292 *
293 * Pipeline information is removed from the resulting
294 * vtkDataSet, so that the reader and its data can be
295 * released
296 */
297bool DataSet::setData(vtkDataSetReader *reader)
298{
299    TRACE("Enter");
300    // Force reading data set
301    reader->SetLookupTableName("");
302    reader->Update();
303
304    _dataSet = reader->GetOutput();
305    _dataSet->SetPipelineInformation(NULL);
306
307    if (_dataSet->GetPointData() != NULL &&
308        _dataSet->GetPointData()->GetScalars() != NULL &&
309        _dataSet->GetPointData()->GetScalars()->GetLookupTable() != NULL) {
310        ERROR("No lookup table should be specified in DataSets");
311    }
312
313    setDefaultArrays();
314
315#ifdef WANT_TRACE
316    print();
317#endif
318    TRACE("Leave");
319    return true;
320}
321
322/**
323 * \brief Set DataSet from existing vtkDataSet object
324 *
325 * Pipeline information is removed from the supplied vtkDataSet
326 */
327bool DataSet::setData(vtkDataSet *ds)
328{
329    _dataSet = ds;
330    _dataSet->SetPipelineInformation(NULL);
331
332    if (_dataSet->GetPointData() != NULL &&
333        _dataSet->GetPointData()->GetScalars() != NULL &&
334        _dataSet->GetPointData()->GetScalars()->GetLookupTable() != NULL) {
335        ERROR("No lookup table should be specified in DataSets");
336    }
337
338    setDefaultArrays();
339
340#ifdef WANT_TRACE
341    print();
342#endif
343    return true;
344}
345
346/**
347 * \brief Copy an existing vtkDataSet object
348 *
349 * Pipeline information is not copied from the supplied vtkDataSet
350 * into this DataSet, but pipeline information in the supplied
351 * vtkDataSet is not removed
352 */
353vtkDataSet *DataSet::copyData(vtkDataSet *ds)
354{
355    if (vtkPolyData::SafeDownCast(ds) != NULL) {
356        _dataSet = vtkSmartPointer<vtkPolyData>::New();
357        _dataSet->DeepCopy(vtkPolyData::SafeDownCast(ds));
358    } else if (vtkStructuredPoints::SafeDownCast(ds) != NULL) {
359        _dataSet = vtkSmartPointer<vtkStructuredPoints>::New();
360        _dataSet->DeepCopy(vtkStructuredPoints::SafeDownCast(ds));
361    } else if (vtkStructuredGrid::SafeDownCast(ds) != NULL) {
362        _dataSet = vtkSmartPointer<vtkStructuredGrid>::New();
363        _dataSet->DeepCopy(vtkStructuredGrid::SafeDownCast(ds));
364    } else if (vtkRectilinearGrid::SafeDownCast(ds) != NULL) {
365        _dataSet = vtkSmartPointer<vtkRectilinearGrid>::New();
366        _dataSet->DeepCopy(vtkRectilinearGrid::SafeDownCast(ds));
367    } else if (vtkUnstructuredGrid::SafeDownCast(ds) != NULL) {
368        _dataSet = vtkSmartPointer<vtkUnstructuredGrid>::New();
369        _dataSet->DeepCopy(vtkUnstructuredGrid::SafeDownCast(ds));
370    } else {
371        ERROR("Unknown data type");
372        return NULL;
373    }
374
375#ifdef WANT_TRACE
376    print();
377#endif   
378    return _dataSet;
379}
380
381/**
382 * \brief Does DataSet lie in the XY plane (z = 0)
383 */
384bool DataSet::isXY() const
385{
386    double bounds[6];
387    getBounds(bounds);
388    return (bounds[4] == 0. && bounds[4] == bounds[5]);
389}
390
391/**
392 * \brief Returns the dimensionality of the AABB
393 */
394int DataSet::numDimensions() const
395{
396    double bounds[6];
397    getBounds(bounds);
398    int numDims = 0;
399    if (bounds[0] != bounds[1])
400        numDims++;
401    if (bounds[2] != bounds[3])
402        numDims++;
403    if (bounds[4] != bounds[5])
404        numDims++;
405
406    return numDims;
407}
408
409/**
410 * \brief Determines if DataSet lies in a principal axis plane
411 * and if so, returns the plane normal and offset from origin
412 */
413bool DataSet::is2D(DataSet::PrincipalPlane *plane, double *offset) const
414{
415    double bounds[6];
416    getBounds(bounds);
417    if (bounds[4] == bounds[5]) {
418        // Z = 0, XY plane
419        if (plane != NULL) {
420            *plane = PLANE_XY;
421        }
422        if (offset != NULL)
423            *offset = bounds[4];
424        return true;
425    } else if (bounds[0] == bounds[1]) {
426        // X = 0, ZY plane
427        if (plane != NULL) {
428            *plane = PLANE_ZY;
429        }
430        if (offset != NULL)
431            *offset = bounds[0];
432        return true;
433    } else if (bounds[2] == bounds[3]) {
434        // Y = 0, XZ plane
435        if (plane != NULL) {
436            *plane = PLANE_XZ;
437         }
438        if (offset != NULL)
439            *offset = bounds[2];
440        return true;
441    }
442    return false;
443}
444
445/**
446 * \brief Determines a principal plane with the
447 * largest two dimensions of the AABB
448 */
449DataSet::PrincipalPlane DataSet::principalPlane() const
450{
451    double bounds[6];
452    getBounds(bounds);
453    double xlen = bounds[1] - bounds[0];
454    double ylen = bounds[3] - bounds[2];
455    double zlen = bounds[5] - bounds[4];
456    if (zlen <= xlen && zlen <= ylen) {
457        return PLANE_XY;
458    } else if (xlen <= ylen && xlen <= zlen) {
459        return PLANE_ZY;
460    } else {
461        return PLANE_XZ;
462    }
463}
464
465/**
466 * \brief Get the name/id of this dataset
467 */
468const std::string& DataSet::getName() const
469{
470    return _name;
471}
472
473/**
474 * \brief Get the underlying VTK DataSet object
475 */
476vtkDataSet *DataSet::getVtkDataSet()
477{
478    return _dataSet;
479}
480
481/**
482 * \brief Get the underlying VTK DataSet subclass class name
483 */
484const char *DataSet::getVtkType() const
485{
486    return _dataSet->GetClassName();
487}
488
489/**
490 * \brief Set the ative scalar array to the named field
491 */
492bool DataSet::setActiveScalars(const char *name)
493{
494    bool found = false;
495    if (_dataSet != NULL) {
496        if (_dataSet->GetPointData() != NULL) {
497            if (_dataSet->GetPointData()->SetActiveScalars(name) >= 0) {
498                TRACE("Set active point data scalars to %s", name);
499                found = true;
500            }
501        }
502        if (_dataSet->GetCellData() != NULL) {
503            if (_dataSet->GetCellData()->SetActiveScalars(name) >= 0) {
504                TRACE("Set active cell data scalars to %s", name);
505                found = true;
506            }
507        }
508    }
509#ifdef WANT_TRACE
510    if (_dataSet->GetPointData() != NULL) {
511        if (_dataSet->GetPointData()->GetScalars() != NULL) {
512            TRACE("Point data scalars: %s", _dataSet->GetPointData()->GetScalars()->GetName());
513        } else {
514            TRACE("NULL point data scalars");
515        }
516    }
517    if (_dataSet->GetCellData() != NULL) {
518        if (_dataSet->GetCellData()->GetScalars() != NULL) {
519            TRACE("Cell data scalars: %s", _dataSet->GetCellData()->GetScalars()->GetName());
520        } else {
521            TRACE("NULL cell data scalars");
522        }
523    }
524#endif
525    return found;
526}
527
528/**
529 * \brief Get the active scalar array field name
530 */
531const char *DataSet::getActiveScalarsName()
532{
533    if (_dataSet != NULL) {
534         if (_dataSet->GetPointData() != NULL &&
535             _dataSet->GetPointData()->GetScalars() != NULL) {
536            return _dataSet->GetPointData()->GetScalars()->GetName();
537        }
538        TRACE("No point scalars");
539        if (_dataSet->GetCellData() != NULL &&
540            _dataSet->GetCellData()->GetScalars() != NULL) {
541            return _dataSet->GetCellData()->GetScalars()->GetName();
542        }
543        TRACE("No cell scalars");
544    }
545    return NULL;
546}
547
548/**
549 * \brief Set the ative vector array to the named field
550 */
551bool DataSet::setActiveVectors(const char *name)
552{
553    bool found = false;
554    if (_dataSet != NULL) {
555        if (_dataSet->GetPointData() != NULL) {
556            if (_dataSet->GetPointData()->SetActiveVectors(name) >= 0) {
557                TRACE("Set active point data vectors to %s", name);
558                found = true;
559            }
560        }
561        if (_dataSet->GetCellData() != NULL) {
562            if (_dataSet->GetCellData()->SetActiveVectors(name) >= 0) {
563                TRACE("Set active cell data vectors to %s", name);
564                found = true;
565            }
566        }
567    }
568
569    return found;
570}
571
572/**
573 * \brief Get the active vector array field name
574 */
575const char *DataSet::getActiveVectorsName()
576{
577    if (_dataSet != NULL) {
578        if (_dataSet->GetPointData() != NULL &&
579            _dataSet->GetPointData()->GetVectors() != NULL) {
580            return _dataSet->GetPointData()->GetVectors()->GetName();
581        }
582        TRACE("No point vectors");
583        if (_dataSet->GetCellData() != NULL &&
584            _dataSet->GetCellData()->GetVectors() != NULL) {
585            return _dataSet->GetCellData()->GetVectors()->GetName();
586        }
587        TRACE("No cell vectors");
588    }
589    return NULL;
590}
591
592/**
593 * \brief Get the range of scalar values in the DataSet
594 */
595void DataSet::getScalarRange(double minmax[2]) const
596{
597    if (_dataSet == NULL)
598        return;
599    if (_dataSet->GetPointData() != NULL &&
600        _dataSet->GetPointData()->GetScalars() != NULL) {
601        _dataSet->GetPointData()->GetScalars()->GetRange(minmax, -1);
602    } else if (_dataSet->GetCellData() != NULL &&
603               _dataSet->GetCellData()->GetScalars() != NULL) {
604        _dataSet->GetCellData()->GetScalars()->GetRange(minmax, -1);
605    }
606}
607
608/**
609 * \brief Get the range of a vector component in the DataSet
610 *
611 * \param[out] minmax The data range
612 * \param[in] component The field component, -1 means magnitude
613 */
614void DataSet::getVectorRange(double minmax[2], int component) const
615{
616    if (_dataSet == NULL)
617        return;
618    if (_dataSet->GetPointData() != NULL &&
619        _dataSet->GetPointData()->GetVectors() != NULL) {
620        _dataSet->GetPointData()->GetVectors()->GetRange(minmax, component);
621    } else if (_dataSet->GetCellData() != NULL &&
622               _dataSet->GetCellData()->GetVectors() != NULL) {
623        _dataSet->GetCellData()->GetVectors()->GetRange(minmax, component);
624    }
625}
626
627/**
628 * \brief Get the range of values for the named field in the DataSet
629 *
630 * \param[out] minmax The data range
631 * \param[in] fieldName The array name
632 * \param[in] component The field component, -1 means magnitude
633 */
634void DataSet::getDataRange(double minmax[2], const char *fieldName, int component) const
635{
636    if (_dataSet == NULL)
637        return;
638    if (_dataSet->GetPointData() != NULL &&
639        _dataSet->GetPointData()->GetArray(fieldName) != NULL) {
640        _dataSet->GetPointData()->GetArray(fieldName)->GetRange(minmax, component);
641    } else if (_dataSet->GetCellData() != NULL &&
642        _dataSet->GetCellData()->GetArray(fieldName) != NULL) {
643        _dataSet->GetCellData()->GetArray(fieldName)->GetRange(minmax, component);
644    } else if (_dataSet->GetFieldData() != NULL &&
645        _dataSet->GetFieldData()->GetArray(fieldName) != NULL) {
646        _dataSet->GetFieldData()->GetArray(fieldName)->GetRange(minmax, component);
647    }
648}
649
650/**
651 * \brief Get the bounds the DataSet
652 */
653void DataSet::getBounds(double bounds[6]) const
654{
655    _dataSet->GetBounds(bounds);
656}
657
658/**
659 * \brief Get the range of cell AABB diagonal lengths in the DataSet
660 */
661void DataSet::getCellSizeRange(double minmax[2], double *average)
662{
663    if (_dataSet == NULL ||
664        _dataSet->GetNumberOfCells() < 1) {
665        minmax[0] = 1;
666        minmax[1] = 1;
667        *average = 1;
668        return;
669    }
670
671    if (_cellSizeRange[0] >= 0.0) {
672        minmax[0] = _cellSizeRange[0];
673        minmax[1] = _cellSizeRange[1];
674        *average = _cellSizeAverage;
675        return;
676    }
677
678    minmax[0] = DBL_MAX;
679    minmax[1] = -DBL_MAX;
680
681    *average = 0;
682    for (int i = 0; i < _dataSet->GetNumberOfCells(); i++) {
683        double length2 = _dataSet->GetCell(i)->GetLength2();
684        if (length2 < minmax[0])
685            minmax[0] = length2;
686        if (length2 > minmax[1])
687            minmax[1] = length2;
688        *average += length2;
689    }
690    if (minmax[0] == DBL_MAX)
691        minmax[0] = 1;
692    if (minmax[1] == -DBL_MAX)
693        minmax[1] = 1;
694
695    minmax[0] = sqrt(minmax[0]);
696    minmax[1] = sqrt(minmax[1]);
697    *average = sqrt(*average/((double)_dataSet->GetNumberOfCells()));
698    _cellSizeRange[0] = minmax[0];
699    _cellSizeRange[1] = minmax[1];
700    _cellSizeAverage = *average;
701}
702
703/**
704 * \brief Get nearest data value given world coordinates x,y,z
705 *
706 * Note: no interpolation is performed on data
707 *
708 * \return the value of the nearest point or 0 if no scalar data available
709 */
710bool DataSet::getScalarValue(double x, double y, double z, double *value) const
711{
712    if (_dataSet == NULL)
713        return false;
714    if (_dataSet->GetPointData() == NULL ||
715        _dataSet->GetPointData()->GetScalars() == NULL) {
716        return false;
717    }
718    vtkIdType pt = _dataSet->FindPoint(x, y, z);
719    if (pt < 0)
720        return false;
721    *value = _dataSet->GetPointData()->GetScalars()->GetComponent(pt, 0);
722    return true;
723}
724
725/**
726 * \brief Get nearest vector data value given world coordinates x,y,z
727 *
728 * Note: no interpolation is performed on data
729 *
730 * \param[in] x World x coordinate to probe
731 * \param[in] y World y coordinate to probe
732 * \param[in] z World z coordinate to probe
733 * \param[out] vector On success, contains the data values
734 * \return boolean indicating success or failure
735 */
736bool DataSet::getVectorValue(double x, double y, double z, double vector[3]) const
737{
738    if (_dataSet == NULL)
739        return false;
740    if (_dataSet->GetPointData() == NULL ||
741        _dataSet->GetPointData()->GetVectors() == NULL) {
742        return false;
743    }
744    vtkIdType pt = _dataSet->FindPoint(x, y, z);
745    if (pt < 0)
746        return false;
747    assert(_dataSet->GetPointData()->GetVectors()->GetNumberOfComponents() == 3);
748    _dataSet->GetPointData()->GetVectors()->GetTuple(pt, vector);
749    return true;
750}
Note: See TracBrowser for help on using the repository browser.