source: vtkvis/trunk/Contour2D.cpp @ 4790

Last change on this file since 4790 was 3680, checked in by ldelgass, 11 years ago

Improved cloud support in vtkvis: handle ugrids with no cells as well as
polydata clouds, add cloudstyle options to some graphics objects.

  • Property svn:eol-style set to native
File size: 17.6 KB
Line 
1/* -*- mode: c++; c-basic-offset: 4; indent-tabs-mode: nil -*- */
2/*
3 * Copyright (C) 2004-2012  HUBzero Foundation, LLC
4 *
5 * Author: Leif Delgass <ldelgass@purdue.edu>
6 */
7
8#include <cfloat>
9
10#include <vtkVersion.h>
11#include <vtkDataSet.h>
12#include <vtkPointData.h>
13#include <vtkCellData.h>
14#include <vtkCellDataToPointData.h>
15#include <vtkContourFilter.h>
16#include <vtkStripper.h>
17#include <vtkPolyDataMapper.h>
18#include <vtkUnstructuredGrid.h>
19#include <vtkProperty.h>
20#include <vtkTransform.h>
21#include <vtkDelaunay2D.h>
22#include <vtkDelaunay3D.h>
23#include <vtkDataSetSurfaceFilter.h>
24
25#include "Contour2D.h"
26#include "Renderer.h"
27#include "Trace.h"
28
29using namespace VtkVis;
30
31Contour2D::Contour2D(int numContours) :
32    GraphicsObject(),
33    _numContours(numContours),
34    _pipelineInitialized(false),
35    _colorMap(NULL),
36    _colorMode(COLOR_BY_SCALAR),
37    _colorFieldType(DataSet::POINT_DATA),
38    _renderer(NULL)
39{
40    _edgeColor[0] = _color[0];
41    _edgeColor[1] = _color[0];
42    _edgeColor[2] = _color[0];
43    _colorFieldRange[0] = DBL_MAX;
44    _colorFieldRange[1] = -DBL_MAX;
45    _lighting = false;
46}
47
48Contour2D::Contour2D(const std::vector<double>& contours) :
49    GraphicsObject(),
50    _numContours(contours.size()),
51    _contours(contours),
52    _pipelineInitialized(false),
53    _colorMap(NULL),
54    _colorMode(COLOR_BY_SCALAR),
55    _colorFieldType(DataSet::POINT_DATA),
56    _renderer(NULL)
57{
58    _edgeColor[0] = _color[0];
59    _edgeColor[1] = _color[0];
60    _edgeColor[2] = _color[0];
61    _colorFieldRange[0] = DBL_MAX;
62    _colorFieldRange[1] = -DBL_MAX;
63    _lighting = false;
64}
65
66Contour2D::~Contour2D()
67{
68#ifdef WANT_TRACE
69    if (_dataSet != NULL)
70        TRACE("Deleting Contour2D for %s", _dataSet->getName().c_str());
71    else
72        TRACE("Deleting Contour2D with NULL DataSet");
73#endif
74}
75
76void Contour2D::setDataSet(DataSet *dataSet,
77                           Renderer *renderer)
78{
79    if (_dataSet != dataSet) {
80        _dataSet = dataSet;
81        _renderer = renderer;
82
83        if (renderer->getUseCumulativeRange()) {
84            renderer->getCumulativeDataRange(_dataRange,
85                                             _dataSet->getActiveScalarsName(),
86                                             1);
87            const char *activeVectors = _dataSet->getActiveVectorsName();
88            if (activeVectors != NULL) {
89                renderer->getCumulativeDataRange(_vectorMagnitudeRange,
90                                                 _dataSet->getActiveVectorsName(),
91                                                 3);
92                for (int i = 0; i < 3; i++) {
93                    renderer->getCumulativeDataRange(_vectorComponentRange[i],
94                                                     _dataSet->getActiveVectorsName(),
95                                                     3, i);
96                }
97            }
98        } else {
99            _dataSet->getScalarRange(_dataRange);
100            _dataSet->getVectorRange(_vectorMagnitudeRange);
101            for (int i = 0; i < 3; i++) {
102                _dataSet->getVectorRange(_vectorComponentRange[i], i);
103            }
104        }
105
106        update();
107    }
108}
109
110/**
111 * \brief Internal method to re-compute contours after a state change
112 */
113void Contour2D::update()
114{
115    if (_dataSet == NULL) {
116        return;
117    }
118    vtkDataSet *ds = _dataSet->getVtkDataSet();
119
120    // Contour filter to generate isolines
121    if (_contourFilter == NULL) {
122        _contourFilter = vtkSmartPointer<vtkContourFilter>::New();
123    }
124
125    if (!_pipelineInitialized) {
126        vtkSmartPointer<vtkCellDataToPointData> cellToPtData;
127
128        if (ds->GetPointData() == NULL ||
129            ds->GetPointData()->GetScalars() == NULL) {
130            TRACE("No scalar point data in dataset %s", _dataSet->getName().c_str());
131            if (ds->GetCellData() != NULL &&
132                ds->GetCellData()->GetScalars() != NULL) {
133                cellToPtData =
134                    vtkSmartPointer<vtkCellDataToPointData>::New();
135#ifdef USE_VTK6
136                cellToPtData->SetInputData(ds);
137#else
138                cellToPtData->SetInput(ds);
139#endif
140                //cellToPtData->PassCellDataOn();
141                cellToPtData->Update();
142                ds = cellToPtData->GetOutput();
143            } else {
144                USER_ERROR("No scalar field was found in the data set");
145            }
146        }
147
148        if (_dataSet->isCloud()) {
149            // DataSet is a point cloud
150            PrincipalPlane plane;
151            double offset;
152            if (_dataSet->is2D(&plane, &offset)) {
153                vtkSmartPointer<vtkDelaunay2D> mesher = vtkSmartPointer<vtkDelaunay2D>::New();
154                if (plane == PLANE_ZY) {
155                    vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
156                    trans->RotateWXYZ(90, 0, 1, 0);
157                    if (offset != 0.0) {
158                        trans->Translate(-offset, 0, 0);
159                    }
160                    mesher->SetTransform(trans);
161                } else if (plane == PLANE_XZ) {
162                    vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
163                    trans->RotateWXYZ(-90, 1, 0, 0);
164                    if (offset != 0.0) {
165                        trans->Translate(0, -offset, 0);
166                    }
167                    mesher->SetTransform(trans);
168                } else if (offset != 0.0) {
169                    // XY with Z offset
170                    vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
171                    trans->Translate(0, 0, -offset);
172                    mesher->SetTransform(trans);
173                }
174#ifdef USE_VTK6
175                mesher->SetInputData(ds);
176#else
177                mesher->SetInput(ds);
178#endif
179                _contourFilter->SetInputConnection(mesher->GetOutputPort());
180            } else {
181                vtkSmartPointer<vtkDelaunay3D> mesher = vtkSmartPointer<vtkDelaunay3D>::New();
182#ifdef USE_VTK6
183                mesher->SetInputData(ds);
184#else
185                mesher->SetInput(ds);
186#endif
187                vtkSmartPointer<vtkDataSetSurfaceFilter> gf = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
188                gf->SetInputConnection(mesher->GetOutputPort());
189                _contourFilter->SetInputConnection(gf->GetOutputPort());
190            }
191        } else if (vtkPolyData::SafeDownCast(ds) != NULL) {
192            // DataSet is a vtkPolyData with lines and/or polygons
193#ifdef USE_VTK6
194            _contourFilter->SetInputData(ds);
195#else
196            _contourFilter->SetInput(ds);
197#endif
198        } else {
199            TRACE("Generating surface for data set");
200            // DataSet is NOT a vtkPolyData
201            vtkSmartPointer<vtkDataSetSurfaceFilter> gf = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
202#ifdef USE_VTK6
203            gf->SetInputData(ds);
204#else
205            gf->SetInput(ds);
206#endif
207            _contourFilter->SetInputConnection(gf->GetOutputPort());
208        }
209    }
210
211    _pipelineInitialized = true;
212
213    _contourFilter->ComputeNormalsOff();
214    _contourFilter->ComputeGradientsOff();
215
216    // Speed up multiple contour computation at cost of extra memory use
217    if (_numContours > 1) {
218        _contourFilter->UseScalarTreeOn();
219    } else {
220        _contourFilter->UseScalarTreeOff();
221    }
222
223    _contourFilter->SetNumberOfContours(_numContours);
224
225    if (_contours.empty()) {
226        // Evenly spaced isovalues
227        _contourFilter->GenerateValues(_numContours, _dataRange[0], _dataRange[1]);
228    } else {
229        // User-supplied isovalues
230        for (int i = 0; i < _numContours; i++) {
231            _contourFilter->SetValue(i, _contours[i]);
232        }
233    }
234
235    initProp();
236
237    if (_mapper == NULL) {
238        _mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
239        _mapper->SetResolveCoincidentTopologyToPolygonOffset();
240        _mapper->ScalarVisibilityOff();
241        vtkSmartPointer<vtkStripper> stripper = vtkSmartPointer<vtkStripper>::New();
242        stripper->SetInputConnection(_contourFilter->GetOutputPort());
243        _mapper->SetInputConnection(stripper->GetOutputPort());
244        getActor()->SetMapper(_mapper);
245    }
246
247    if (_lut == NULL) {
248        setColorMap(ColorMap::getDefault());
249        setColorMode(_colorMode);
250    }
251
252    _mapper->Update();
253}
254
255void Contour2D::updateRanges(Renderer *renderer)
256{
257    if (_dataSet == NULL) {
258        ERROR("called before setDataSet");
259        return;
260    }
261
262    if (renderer->getUseCumulativeRange()) {
263        renderer->getCumulativeDataRange(_dataRange,
264                                         _dataSet->getActiveScalarsName(),
265                                         1);
266        const char *activeVectors = _dataSet->getActiveVectorsName();
267        if (activeVectors != NULL) {
268            renderer->getCumulativeDataRange(_vectorMagnitudeRange,
269                                             _dataSet->getActiveVectorsName(),
270                                             3);
271            for (int i = 0; i < 3; i++) {
272                renderer->getCumulativeDataRange(_vectorComponentRange[i],
273                                                 _dataSet->getActiveVectorsName(),
274                                                 3, i);
275            }
276        }
277    } else {
278        _dataSet->getScalarRange(_dataRange);
279        _dataSet->getVectorRange(_vectorMagnitudeRange);
280        for (int i = 0; i < 3; i++) {
281            _dataSet->getVectorRange(_vectorComponentRange[i], i);
282        }
283    }
284 
285    // Need to update color map ranges
286    double *rangePtr = _colorFieldRange;
287    if (_colorFieldRange[0] > _colorFieldRange[1]) {
288        rangePtr = NULL;
289    }
290    setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
291
292    if (_contours.empty() && _numContours > 0) {
293        // Contour isovalues need to be recomputed
294        update();
295    }
296}
297
298void Contour2D::setContourField(const char *name)
299{
300    if (_contourFilter != NULL) {
301        _contourFilter->SetInputArrayToProcess(0, 0, 0, vtkDataObject::FIELD_ASSOCIATION_POINTS, name);
302        update();
303    }
304}
305
306void Contour2D::setColorMode(ColorMode mode)
307{
308    _colorMode = mode;
309    if (_dataSet == NULL)
310        return;
311
312    switch (mode) {
313    case COLOR_BY_SCALAR:
314        setColorMode(mode,
315                     _dataSet->getActiveScalarsType(),
316                     _dataSet->getActiveScalarsName(),
317                     _dataRange);
318        break;
319    case COLOR_BY_VECTOR_MAGNITUDE:
320        setColorMode(mode,
321                     _dataSet->getActiveVectorsType(),
322                     _dataSet->getActiveVectorsName(),
323                     _vectorMagnitudeRange);
324        break;
325    case COLOR_BY_VECTOR_X:
326        setColorMode(mode,
327                     _dataSet->getActiveVectorsType(),
328                     _dataSet->getActiveVectorsName(),
329                     _vectorComponentRange[0]);
330        break;
331    case COLOR_BY_VECTOR_Y:
332        setColorMode(mode,
333                     _dataSet->getActiveVectorsType(),
334                     _dataSet->getActiveVectorsName(),
335                     _vectorComponentRange[1]);
336        break;
337    case COLOR_BY_VECTOR_Z:
338        setColorMode(mode,
339                     _dataSet->getActiveVectorsType(),
340                     _dataSet->getActiveVectorsName(),
341                     _vectorComponentRange[2]);
342        break;
343    case COLOR_CONSTANT:
344    default:
345        setColorMode(mode, DataSet::POINT_DATA, NULL, NULL);
346        break;
347    }
348}
349
350void Contour2D::setColorMode(ColorMode mode,
351                             const char *name, double range[2])
352{
353    if (_dataSet == NULL)
354        return;
355    DataSet::DataAttributeType type = DataSet::POINT_DATA;
356    int numComponents = 1;
357    if (name != NULL && strlen(name) > 0 &&
358        !_dataSet->getFieldInfo(name, &type, &numComponents)) {
359        ERROR("Field not found: %s", name);
360        return;
361    }
362    setColorMode(mode, type, name, range);
363}
364
365void Contour2D::setColorMode(ColorMode mode, DataSet::DataAttributeType type,
366                             const char *name, double range[2])
367{
368    _colorMode = mode;
369    _colorFieldType = type;
370    if (name == NULL)
371        _colorFieldName.clear();
372    else
373        _colorFieldName = name;
374    if (range == NULL) {
375        _colorFieldRange[0] = DBL_MAX;
376        _colorFieldRange[1] = -DBL_MAX;
377    } else {
378        memcpy(_colorFieldRange, range, sizeof(double)*2);
379    }
380
381    if (_dataSet == NULL || _mapper == NULL)
382        return;
383
384    switch (type) {
385    case DataSet::POINT_DATA:
386        _mapper->SetScalarModeToUsePointFieldData();
387        break;
388    case DataSet::CELL_DATA:
389        _mapper->SetScalarModeToUseCellFieldData();
390        break;
391    default:
392        ERROR("Unsupported DataAttributeType: %d", type);
393        return;
394    }
395
396    if (name != NULL && strlen(name) > 0) {
397        _mapper->SelectColorArray(name);
398    } else {
399        _mapper->SetScalarModeToDefault();
400    }
401
402    if (_lut != NULL) {
403        if (range != NULL) {
404            _lut->SetRange(range);
405        } else if (name != NULL && strlen(name) > 0) {
406            double r[2];
407            int comp = -1;
408            if (mode == COLOR_BY_VECTOR_X)
409                comp = 0;
410            else if (mode == COLOR_BY_VECTOR_Y)
411                comp = 1;
412            else if (mode == COLOR_BY_VECTOR_Z)
413                comp = 2;
414
415            if (_renderer->getUseCumulativeRange()) {
416                int numComponents;
417                if  (!_dataSet->getFieldInfo(name, type, &numComponents)) {
418                    ERROR("Field not found: %s, type: %d", name, type);
419                    return;
420                } else if (numComponents < comp+1) {
421                    ERROR("Request for component %d in field with %d components",
422                          comp, numComponents);
423                    return;
424                }
425                _renderer->getCumulativeDataRange(r, name, type, numComponents, comp);
426            } else {
427                _dataSet->getDataRange(r, name, type, comp);
428            }
429            _lut->SetRange(r);
430        }
431    }
432
433    switch (mode) {
434    case COLOR_BY_SCALAR:
435        _mapper->ScalarVisibilityOn();
436        break;
437    case COLOR_BY_VECTOR_MAGNITUDE:
438        _mapper->ScalarVisibilityOn();
439        if (_lut != NULL) {
440            _lut->SetVectorModeToMagnitude();
441        }
442        break;
443    case COLOR_BY_VECTOR_X:
444        _mapper->ScalarVisibilityOn();
445        if (_lut != NULL) {
446            _lut->SetVectorModeToComponent();
447            _lut->SetVectorComponent(0);
448        }
449        break;
450    case COLOR_BY_VECTOR_Y:
451        _mapper->ScalarVisibilityOn();
452        if (_lut != NULL) {
453            _lut->SetVectorModeToComponent();
454            _lut->SetVectorComponent(1);
455        }
456        break;
457    case COLOR_BY_VECTOR_Z:
458        _mapper->ScalarVisibilityOn();
459        if (_lut != NULL) {
460            _lut->SetVectorModeToComponent();
461            _lut->SetVectorComponent(2);
462        }
463        break;
464    case COLOR_CONSTANT:
465    default:
466        _mapper->ScalarVisibilityOff();
467        break;
468    }
469}
470
471/**
472 * \brief Called when the color map has been edited
473 */
474void Contour2D::updateColorMap()
475{
476    setColorMap(_colorMap);
477}
478
479/**
480 * \brief Associate a colormap lookup table with the DataSet
481 */
482void Contour2D::setColorMap(ColorMap *cmap)
483{
484    if (cmap == NULL)
485        return;
486
487    _colorMap = cmap;
488 
489    if (_lut == NULL) {
490        _lut = vtkSmartPointer<vtkLookupTable>::New();
491        if (_mapper != NULL) {
492            _mapper->UseLookupTableScalarRangeOn();
493            _mapper->SetLookupTable(_lut);
494        }
495        _lut->DeepCopy(cmap->getLookupTable());
496        switch (_colorMode) {
497        case COLOR_CONSTANT:
498        case COLOR_BY_SCALAR:
499            _lut->SetRange(_dataRange);
500            break;
501        case COLOR_BY_VECTOR_MAGNITUDE:
502            _lut->SetRange(_vectorMagnitudeRange);
503            break;
504        case COLOR_BY_VECTOR_X:
505            _lut->SetRange(_vectorComponentRange[0]);
506            break;
507        case COLOR_BY_VECTOR_Y:
508            _lut->SetRange(_vectorComponentRange[1]);
509            break;
510        case COLOR_BY_VECTOR_Z:
511            _lut->SetRange(_vectorComponentRange[2]);
512            break;
513        default:
514            break;
515        }
516    } else {
517        double range[2];
518        _lut->GetTableRange(range);
519        _lut->DeepCopy(cmap->getLookupTable());
520        _lut->SetRange(range);
521        _lut->Modified();
522    }
523
524    switch (_colorMode) {
525    case COLOR_BY_VECTOR_MAGNITUDE:
526        _lut->SetVectorModeToMagnitude();
527        break;
528    case COLOR_BY_VECTOR_X:
529        _lut->SetVectorModeToComponent();
530        _lut->SetVectorComponent(0);
531        break;
532    case COLOR_BY_VECTOR_Y:
533        _lut->SetVectorModeToComponent();
534        _lut->SetVectorComponent(1);
535        break;
536    case COLOR_BY_VECTOR_Z:
537        _lut->SetVectorModeToComponent();
538        _lut->SetVectorComponent(2);
539        break;
540    default:
541         break;
542    }
543}
544
545/**
546 * \brief Specify number of evenly spaced contour lines to render
547 *
548 * Will override any existing contours
549 */
550void Contour2D::setNumContours(int numContours)
551{
552    _contours.clear();
553    _numContours = numContours;
554
555    update();
556}
557
558/**
559 * \brief Specify an explicit list of contour isovalues
560 *
561 * Will override any existing contours
562 */
563void Contour2D::setContourList(const std::vector<double>& contours)
564{
565    _contours = contours;
566    _numContours = (int)_contours.size();
567
568    update();
569}
570
571/**
572 * \brief Get the number of contours
573 */
574int Contour2D::getNumContours() const
575{
576    return _numContours;
577}
578
579/**
580 * \brief Get the contour list (may be empty if number of contours
581 * was specified in place of a list)
582 */
583const std::vector<double>& Contour2D::getContourList() const
584{
585    return _contours;
586}
587
588/**
589 * \brief Set a group of world coordinate planes to clip rendering
590 *
591 * Passing NULL for planes will remove all cliping planes
592 */
593void Contour2D::setClippingPlanes(vtkPlaneCollection *planes)
594{
595    if (_mapper != NULL) {
596        _mapper->SetClippingPlanes(planes);
597    }
598}
Note: See TracBrowser for help on using the repository browser.