source: vtkvis/trunk/HeightMap.cpp @ 5835

Last change on this file since 5835 was 5835, checked in by ldelgass, 8 years ago

Require VTK >= 6.0.0

  • Property svn:eol-style set to native
File size: 41.8 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 <cassert>
9
10#include <vtkDataSet.h>
11#include <vtkPointData.h>
12#include <vtkCellData.h>
13#include <vtkTrivialProducer.h>
14#include <vtkCellDataToPointData.h>
15#include <vtkPolyDataMapper.h>
16#include <vtkUnstructuredGrid.h>
17#include <vtkProperty.h>
18#include <vtkImageData.h>
19#include <vtkLookupTable.h>
20#include <vtkTransform.h>
21#include <vtkDelaunay2D.h>
22#include <vtkDelaunay3D.h>
23#include <vtkGaussianSplatter.h>
24#include <vtkExtractVOI.h>
25#include <vtkDataSetSurfaceFilter.h>
26#include <vtkContourFilter.h>
27#include <vtkStripper.h>
28#include <vtkWarpScalar.h>
29#include <vtkPropAssembly.h>
30#include <vtkCutter.h>
31#include <vtkPlane.h>
32
33#include "HeightMap.h"
34#include "Renderer.h"
35#include "Trace.h"
36
37using namespace VtkVis;
38
39#define TRANSLATE_TO_ZORIGIN
40
41HeightMap::HeightMap(int numContours, double heightScale) :
42    GraphicsObject(),
43    _numContours(numContours),
44    _contourColorMap(false),
45    _contourEdgeWidth(1.0),
46    _warpScale(heightScale),
47    _sliceAxis(Z_AXIS),
48    _pipelineInitialized(false),
49    _cloudStyle(CLOUD_MESH),
50    _colorMap(NULL),
51    _colorMode(COLOR_BY_SCALAR),
52    _colorFieldType(DataSet::POINT_DATA),
53    _renderer(NULL)
54{
55    _contourEdgeColor[0] = 0.0f;
56    _contourEdgeColor[1] = 0.0f;
57    _contourEdgeColor[2] = 0.0f;
58    _colorFieldRange[0] = DBL_MAX;
59    _colorFieldRange[1] = -DBL_MAX;
60}
61
62HeightMap::HeightMap(const std::vector<double>& contours, double heightScale) :
63    GraphicsObject(),
64    _numContours(contours.size()),
65    _contours(contours),
66    _contourColorMap(false),
67    _contourEdgeWidth(1.0),
68    _warpScale(heightScale),
69    _sliceAxis(Z_AXIS),
70    _pipelineInitialized(false),
71    _cloudStyle(CLOUD_MESH),
72    _colorMap(NULL),
73    _colorMode(COLOR_BY_SCALAR),
74    _colorFieldType(DataSet::POINT_DATA),
75    _renderer(NULL)
76{
77    _contourEdgeColor[0] = 0.0f;
78    _contourEdgeColor[1] = 0.0f;
79    _contourEdgeColor[2] = 0.0f;
80    _colorFieldRange[0] = DBL_MAX;
81    _colorFieldRange[1] = -DBL_MAX;
82}
83
84HeightMap::~HeightMap()
85{
86#ifdef WANT_TRACE
87    if (_dataSet != NULL)
88        TRACE("Deleting HeightMap for %s", _dataSet->getName().c_str());
89    else
90        TRACE("Deleting HeightMap with NULL DataSet");
91#endif
92}
93
94void HeightMap::setDataSet(DataSet *dataSet,
95                           Renderer *renderer)
96{
97    if (_dataSet != dataSet) {
98        _dataSet = dataSet;
99        _renderer = renderer;
100
101        if (_dataSet != NULL) {
102            TRACE("DataSet name: '%s' type: %s",
103                  _dataSet->getName().c_str(),
104                  _dataSet->getVtkType());
105
106            if (renderer->getUseCumulativeRange()) {
107                renderer->getCumulativeDataRange(_dataRange,
108                                                 _dataSet->getActiveScalarsName(),
109                                                 1);
110                const char *activeVectors = _dataSet->getActiveVectorsName();
111                if (activeVectors != NULL) {
112                    renderer->getCumulativeDataRange(_vectorMagnitudeRange,
113                                                     activeVectors,
114                                                     3);
115                    for (int i = 0; i < 3; i++) {
116                        renderer->getCumulativeDataRange(_vectorComponentRange[i],
117                                                         activeVectors,
118                                                         3, i);
119                    }
120                }
121            } else {
122                _dataSet->getScalarRange(_dataRange);
123                _dataSet->getVectorRange(_vectorMagnitudeRange);
124                for (int i = 0; i < 3; i++) {
125                    _dataSet->getVectorRange(_vectorComponentRange[i], i);
126                }
127            }
128        }
129
130        update();
131    }
132}
133
134/**
135 * Compute a data scaling factor to make maximum height (at default _warpScale)
136 * equivalent to largest dimension of data set
137 */
138void HeightMap::computeDataScale()
139{
140    if (_dataSet == NULL)
141        return;
142
143    double bounds[6];
144    double boundsRange = 0.0;
145    _dataSet->getBounds(bounds);
146    for (int i = 0; i < 6; i+=2) {
147        double r = bounds[i+1] - bounds[i];
148        if (r > boundsRange)
149            boundsRange = r;
150    }
151    double datRange = _dataRange[1] - _dataRange[0];
152    if (datRange < 1.0e-17) {
153        _dataScale = 1.0;
154    } else {
155        _dataScale = boundsRange / datRange;
156    }
157    TRACE("Bounds range: %g data range: %g scaling: %g",
158          boundsRange, datRange, _dataScale);
159}
160
161/**
162 * \brief Create and initialize VTK Props to render the colormapped dataset
163 */
164void HeightMap::initProp()
165{
166    if (_dsActor == NULL) {
167        _dsActor = vtkSmartPointer<vtkActor>::New();
168        _dsActor->GetProperty()->SetOpacity(_opacity);
169        _dsActor->GetProperty()->SetColor(_color[0],
170                                          _color[1],
171                                          _color[2]);
172        _dsActor->GetProperty()->SetEdgeColor(_edgeColor[0],
173                                              _edgeColor[1],
174                                              _edgeColor[2]);
175        _dsActor->GetProperty()->SetLineWidth(_edgeWidth);
176        _dsActor->GetProperty()->EdgeVisibilityOff();
177        _dsActor->GetProperty()->SetAmbient(0.2);
178        _dsActor->GetProperty()->LightingOn();
179    }
180    if (_contourActor == NULL) {
181        _contourActor = vtkSmartPointer<vtkActor>::New();
182        _contourActor->GetProperty()->SetOpacity(_opacity);
183        _contourActor->GetProperty()->SetColor(_contourEdgeColor[0],
184                                               _contourEdgeColor[1],
185                                               _contourEdgeColor[2]);
186        _contourActor->GetProperty()->SetEdgeColor(_contourEdgeColor[0],
187                                                   _contourEdgeColor[1],
188                                                   _contourEdgeColor[2]);
189        _contourActor->GetProperty()->SetLineWidth(_contourEdgeWidth);
190        _contourActor->GetProperty()->EdgeVisibilityOff();
191        _contourActor->GetProperty()->SetAmbient(1.0);
192        _contourActor->GetProperty()->LightingOff();
193    }
194    // Contour actor is added in update() if contours are produced
195    if (_prop == NULL) {
196        _prop = vtkSmartPointer<vtkAssembly>::New();
197        getAssembly()->AddPart(_dsActor);
198    } else {
199        getAssembly()->RemovePart(_contourActor);
200    }
201}
202
203/**
204 * \brief Internal method to set up pipeline after a state change
205 */
206void HeightMap::update()
207{
208    if (_dataSet == NULL)
209        return;
210
211    TRACE("DataSet: %s", _dataSet->getName().c_str());
212
213    vtkDataSet *ds = _dataSet->getVtkDataSet();
214
215    computeDataScale();
216
217    // Contour filter to generate isolines
218    if (_contourFilter == NULL) {
219        _contourFilter = vtkSmartPointer<vtkContourFilter>::New();
220    }
221
222    // Mapper, actor to render color-mapped data set
223    if (_mapper == NULL) {
224        _mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
225        // Map scalars through lookup table regardless of type
226        _mapper->SetColorModeToMapScalars();
227    }
228
229    if (!_pipelineInitialized) {
230        vtkAlgorithmOutput *dsOutput = NULL;
231        vtkSmartPointer<vtkCellDataToPointData> cellToPtData;
232        vtkSmartPointer<vtkTrivialProducer> prod;
233
234        if (ds->GetPointData() == NULL ||
235            ds->GetPointData()->GetScalars() == NULL) {
236            TRACE("No scalar point data in dataset %s", _dataSet->getName().c_str());
237            if (ds->GetCellData() != NULL &&
238                ds->GetCellData()->GetScalars() != NULL) {
239                cellToPtData =
240                    vtkSmartPointer<vtkCellDataToPointData>::New();
241                cellToPtData->SetInputData(ds);
242                cellToPtData->PassCellDataOn();
243                dsOutput = cellToPtData->GetOutputPort();
244            } else {
245                USER_ERROR("No scalar field was found in the data set.");
246                return;
247            }
248        } else {
249            prod = vtkSmartPointer<vtkTrivialProducer>::New();
250            prod->SetOutput(ds);
251            dsOutput = prod->GetOutputPort();
252        }
253
254        _splatter = NULL;
255
256        if (_dataSet->isCloud()) {
257            // DataSet is a point cloud
258            PrincipalPlane plane;
259            double offset;
260            if (_dataSet->is2D(&plane, &offset)) {
261                if (_cloudStyle == CLOUD_MESH) {
262                    // Result of Delaunay2D is a PolyData
263                    vtkSmartPointer<vtkDelaunay2D> mesher = vtkSmartPointer<vtkDelaunay2D>::New();
264                    if (plane == PLANE_ZY) {
265                        vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
266                        trans->RotateWXYZ(90, 0, 1, 0);
267                        if (offset != 0.0) {
268                            trans->Translate(-offset, 0, 0);
269                        }
270                        mesher->SetTransform(trans);
271                        _sliceAxis = X_AXIS;
272                    } else if (plane == PLANE_XZ) {
273                        vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
274                        trans->RotateWXYZ(-90, 1, 0, 0);
275                        if (offset != 0.0) {
276                            trans->Translate(0, -offset, 0);
277                        }
278                        mesher->SetTransform(trans);
279                        _sliceAxis = Y_AXIS;
280                    } else if (offset != 0.0) {
281                        // XY with Z offset
282                        vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
283                        trans->Translate(0, 0, -offset);
284                        mesher->SetTransform(trans);
285                    }
286                    mesher->SetInputConnection(dsOutput);
287                    vtkAlgorithmOutput *warpOutput = initWarp(mesher->GetOutputPort());
288                    _mapper->SetInputConnection(warpOutput);
289                    _contourFilter->SetInputConnection(warpOutput);
290                } else {
291                    // _cloudStyle == CLOUD_SPLAT
292                    if (_splatter == NULL)
293                        _splatter = vtkSmartPointer<vtkGaussianSplatter>::New();
294                    if (_volumeSlicer == NULL)
295                        _volumeSlicer = vtkSmartPointer<vtkExtractVOI>::New();
296                    _splatter->SetInputConnection(dsOutput);
297                    int dims[3];
298                    _splatter->GetSampleDimensions(dims);
299                    TRACE("Sample dims: %d %d %d", dims[0], dims[1], dims[2]);
300                    if (plane == PLANE_ZY) {
301                        dims[0] = 3;
302                        _volumeSlicer->SetVOI(1, 1, 0, dims[1]-1, 0, dims[1]-1);
303                        _sliceAxis = X_AXIS;
304                    } else if (plane == PLANE_XZ) {
305                        dims[1] = 3;
306                        _volumeSlicer->SetVOI(0, dims[0]-1, 1, 1, 0, dims[2]-1);
307                        _sliceAxis = Y_AXIS;
308                    } else {
309                        dims[2] = 3;
310                        _volumeSlicer->SetVOI(0, dims[0]-1, 0, dims[1]-1, 1, 1);
311                    }
312                    _splatter->SetSampleDimensions(dims);
313                    double bounds[6];
314                    _splatter->Update();
315                    _splatter->GetModelBounds(bounds);
316                    TRACE("Model bounds: %g %g %g %g %g %g",
317                          bounds[0], bounds[1],
318                          bounds[2], bounds[3],
319                          bounds[4], bounds[5]);
320                    _volumeSlicer->SetInputConnection(_splatter->GetOutputPort());
321                    _volumeSlicer->SetSampleRate(1, 1, 1);
322                    vtkSmartPointer<vtkDataSetSurfaceFilter> gf = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
323                    gf->UseStripsOn();
324                    gf->SetInputConnection(_volumeSlicer->GetOutputPort());
325                    vtkAlgorithmOutput *warpOutput = initWarp(gf->GetOutputPort());
326                    _mapper->SetInputConnection(warpOutput);
327                    _contourFilter->SetInputConnection(warpOutput);
328                }
329            } else {
330                // 3D point cloud
331                vtkSmartPointer<vtkDataSetSurfaceFilter> gf = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
332                gf->UseStripsOn();
333                if (_cloudStyle == CLOUD_MESH) {
334                     // Result of Delaunay3D mesher is unstructured grid
335                    vtkSmartPointer<vtkDelaunay3D> mesher = vtkSmartPointer<vtkDelaunay3D>::New();
336                    mesher->SetInputConnection(dsOutput);
337                    // Run the mesher
338                    mesher->Update();
339                    // Get bounds of resulting grid
340                    double bounds[6];
341                    mesher->GetOutput()->GetBounds(bounds);
342                    // Sample a plane within the grid bounding box
343                    vtkSmartPointer<vtkCutter> cutter = vtkSmartPointer<vtkCutter>::New();
344                    cutter->SetInputConnection(mesher->GetOutputPort());
345                    if (_cutPlane == NULL) {
346                        _cutPlane = vtkSmartPointer<vtkPlane>::New();
347                    }
348                    _cutPlane->SetNormal(0, 0, 1);
349                    _cutPlane->SetOrigin(0,
350                                         0,
351                                         bounds[4] + (bounds[5]-bounds[4])/2.);
352                    cutter->SetCutFunction(_cutPlane);
353                    gf->SetInputConnection(cutter->GetOutputPort());
354                } else {
355                    // _cloudStyle == CLOUD_SPLAT
356                    if (_splatter == NULL)
357                        _splatter = vtkSmartPointer<vtkGaussianSplatter>::New();
358                    _splatter->SetInputConnection(dsOutput);
359                    int dims[3];
360                    _splatter->GetSampleDimensions(dims);
361                    TRACE("Sample dims: %d %d %d", dims[0], dims[1], dims[2]);
362                    dims[2] = 3;
363                    _splatter->SetSampleDimensions(dims);
364                    double bounds[6];
365                    _splatter->Update();
366                    _splatter->GetModelBounds(bounds);
367                    TRACE("Model bounds: %g %g %g %g %g %g",
368                          bounds[0], bounds[1],
369                          bounds[2], bounds[3],
370                          bounds[4], bounds[5]);
371                    if (_volumeSlicer == NULL)
372                        _volumeSlicer = vtkSmartPointer<vtkExtractVOI>::New();
373                    _volumeSlicer->SetInputConnection(_splatter->GetOutputPort());
374                    _volumeSlicer->SetVOI(0, dims[0]-1, 0, dims[1]-1, 1, 1);
375                    _volumeSlicer->SetSampleRate(1, 1, 1);
376                    gf->SetInputConnection(_volumeSlicer->GetOutputPort());
377                }
378                vtkAlgorithmOutput *warpOutput = initWarp(gf->GetOutputPort());
379                _mapper->SetInputConnection(warpOutput);
380                _contourFilter->SetInputConnection(warpOutput);
381            }
382        } else if (vtkPolyData::SafeDownCast(ds) != NULL) {
383            // DataSet is a vtkPolyData with lines and/or polygons
384            vtkAlgorithmOutput *warpOutput = initWarp(dsOutput);
385            _mapper->SetInputConnection(warpOutput);
386            _contourFilter->SetInputConnection(warpOutput);
387        } else {
388            // DataSet is NOT a vtkPolyData
389            // Can be: image/volume/uniform grid, structured grid, unstructured grid, rectilinear grid
390            vtkSmartPointer<vtkDataSetSurfaceFilter> gf = vtkSmartPointer<vtkDataSetSurfaceFilter>::New();
391            gf->UseStripsOn();
392            vtkImageData *imageData = vtkImageData::SafeDownCast(ds);
393            if (!_dataSet->is2D() && imageData != NULL) {
394                // 3D image/volume/uniform grid
395                if (_volumeSlicer == NULL)
396                    _volumeSlicer = vtkSmartPointer<vtkExtractVOI>::New();
397                int dims[3];
398                imageData->GetDimensions(dims);
399                TRACE("Image data dimensions: %d %d %d", dims[0], dims[1], dims[2]);
400                _volumeSlicer->SetInputConnection(dsOutput);
401                _volumeSlicer->SetVOI(0, dims[0]-1, 0, dims[1]-1, (dims[2]-1)/2, (dims[2]-1)/2);
402                _volumeSlicer->SetSampleRate(1, 1, 1);
403                gf->SetInputConnection(_volumeSlicer->GetOutputPort());
404            } else if (!_dataSet->is2D() && imageData == NULL) {
405                // 3D structured grid, unstructured grid, or rectilinear grid
406                double bounds[6];
407                ds->GetBounds(bounds);
408                // Sample a plane within the grid bounding box
409                vtkSmartPointer<vtkCutter> cutter = vtkSmartPointer<vtkCutter>::New();
410                cutter->SetInputConnection(dsOutput);
411                if (_cutPlane == NULL) {
412                    _cutPlane = vtkSmartPointer<vtkPlane>::New();
413                }
414                _cutPlane->SetNormal(0, 0, 1);
415                _cutPlane->SetOrigin(0,
416                                     0,
417                                     bounds[4] + (bounds[5]-bounds[4])/2.);
418                cutter->SetCutFunction(_cutPlane);
419                gf->SetInputConnection(cutter->GetOutputPort());
420            } else {
421                // 2D data
422                gf->SetInputConnection(dsOutput);
423            }
424            vtkAlgorithmOutput *warpOutput = initWarp(gf->GetOutputPort());
425            _mapper->SetInputConnection(warpOutput);
426            _contourFilter->SetInputConnection(warpOutput);
427        }
428    }
429
430    _contourFilter->ComputeNormalsOff();
431    _contourFilter->ComputeGradientsOff();
432
433    // Speed up multiple contour computation at cost of extra memory use
434    if (_numContours > 1) {
435        _contourFilter->UseScalarTreeOn();
436    } else {
437        _contourFilter->UseScalarTreeOff();
438    }
439
440    _contourFilter->SetNumberOfContours(_numContours);
441
442    if (_contours.empty()) {
443        // Evenly spaced isovalues
444        if (_colorFieldRange[0] <= _colorFieldRange[1]) {
445            _contourFilter->GenerateValues(_numContours, _colorFieldRange[0], _colorFieldRange[1]);
446        } else {
447            _contourFilter->GenerateValues(_numContours, _dataRange[0], _dataRange[1]);
448        }
449    } else {
450        // User-supplied isovalues
451        for (int i = 0; i < _numContours; i++) {
452            _contourFilter->SetValue(i, _contours[i]);
453        }
454    }
455
456    initProp();
457
458    if (_contourMapper == NULL) {
459        _contourMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
460        _contourMapper->ScalarVisibilityOff();
461        _contourMapper->SetResolveCoincidentTopologyToPolygonOffset();
462        vtkSmartPointer<vtkStripper> stripper = vtkSmartPointer<vtkStripper>::New();
463        stripper->SetInputConnection(_contourFilter->GetOutputPort());
464        _contourMapper->SetInputConnection(stripper->GetOutputPort());
465        _contourActor->SetMapper(_contourMapper);
466        //_contourActor->InterpolateScalarsBeforeMappingOn();
467    }
468
469    setInterpolateBeforeMapping(true);
470
471    if (_lut == NULL) {
472        setColorMap(ColorMap::getDefault());
473        setColorMode(_colorMode);
474    } else if (!_pipelineInitialized) {
475        double *rangePtr = _colorFieldRange;
476        if (_colorFieldRange[0] > _colorFieldRange[1]) {
477            rangePtr = NULL;
478        }
479        setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
480    }
481
482    _pipelineInitialized = true;
483
484    TRACE("DataSet active scalars: %s", _dataSet->getActiveScalarsName());
485    // Ensure updated dataScale is applied
486    if (_warp != NULL) {
487        _warp->SetScaleFactor(_warpScale * _dataScale);
488    }
489#ifdef TRANSLATE_TO_ZORIGIN
490    double pos[3];
491    pos[0] = 0;
492    pos[1] = 0;
493    pos[2] = - _warpScale * _dataScale * _dataRange[0];
494    setPosition(pos);
495    TRACE("Z_POS: %g", pos[2]);
496#endif
497    _dsActor->SetMapper(_mapper);
498
499    _mapper->Update();
500    _contourMapper->Update();
501
502    // Only add contour actor to assembly if contours were
503    // produced, in order to prevent messing up assembly bounds
504    double bounds[6];
505    _contourActor->GetBounds(bounds);
506    if (bounds[0] <= bounds[1]) {
507        getAssembly()->AddPart(_contourActor);
508    }
509}
510
511void HeightMap::setCloudStyle(CloudStyle style)
512{
513    if (style != _cloudStyle) {
514        _cloudStyle = style;
515        if (_dataSet != NULL) {
516            _pipelineInitialized = false;
517            update();
518        }
519    }
520}
521
522void HeightMap::setInterpolateBeforeMapping(bool state)
523{
524    if (_mapper != NULL) {
525        _mapper->SetInterpolateScalarsBeforeMapping((state ? 1 : 0));
526    }
527}
528
529void HeightMap::setAspect(double aspect)
530{
531    double bounds[6];
532    vtkDataSet *ds = _dataSet->getVtkDataSet();
533    ds->GetBounds(bounds);
534    double size[3];
535    size[0] = bounds[1] - bounds[0];
536    size[1] = bounds[3] - bounds[2];
537    size[2] = bounds[5] - bounds[4];
538    double scale[3];
539    scale[0] = scale[1] = scale[2] = 1.;
540
541    if (aspect == 1.0) {
542        // Square
543        switch (_sliceAxis) {
544        case X_AXIS: {
545            if (size[1] > size[2] && size[2] > 0.0) {
546                scale[2] = size[1] / size[2];
547            } else if (size[2] > size[1] && size[1] > 0.0) {
548                scale[1] = size[2] / size[1];
549            }
550        }
551            break;
552        case Y_AXIS: {
553            if (size[0] > size[2] && size[2] > 0.0) {
554                scale[2] = size[0] / size[2];
555            } else if (size[2] > size[0] && size[0] > 0.0) {
556                scale[0] = size[2] / size[0];
557            }
558        }
559            break;
560        case Z_AXIS: {
561            if (size[0] > size[1] && size[1] > 0.0) {
562                scale[1] = size[0] / size[1];
563            } else if (size[1] > size[0] && size[0] > 0.0) {
564                scale[0] = size[1] / size[0];
565            }
566        }
567            break;
568        }
569    } else if (aspect != 0.0) {
570        switch (_sliceAxis) {
571        case X_AXIS: {
572            if (aspect > 1.0) {
573                if (size[2] > size[1] && size[1] > 0.0) {
574                    scale[1] = (size[2] / aspect) / size[1];
575                } else if (size[2] > 0.0) {
576                    scale[2] = (size[1] * aspect) / size[2];
577                }
578            } else {
579                if (size[1] > size[2] && size[2] > 0.0) {
580                    scale[2] = (size[1] * aspect) / size[2];
581                } else if (size[1] > 0.0) {
582                    scale[1] = (size[2] / aspect) / size[1];
583                }
584            }
585        }
586            break;
587        case Y_AXIS: {
588            if (aspect > 1.0) {
589                if (size[0] > size[2] && size[2] > 0.0) {
590                    scale[2] = (size[0] / aspect) / size[2];
591                } else if (size[0] > 0.0) {
592                    scale[0] = (size[2] * aspect) / size[0];
593                }
594            } else {
595                if (size[2] > size[0] && size[0] > 0.0) {
596                    scale[0] = (size[2] * aspect) / size[0];
597                } else if (size[2] > 0.0) {
598                    scale[2] = (size[0] / aspect) / size[2];
599                }
600            }
601        }
602            break;
603        case Z_AXIS: {
604            if (aspect > 1.0) {
605                if (size[0] > size[1] && size[1] > 0.0) {
606                    scale[1] = (size[0] / aspect) / size[1];
607                } else if (size[0] > 0.0) {
608                    scale[0] = (size[1] * aspect) / size[0];
609                }
610            } else {
611                if (size[1] > size[0] && size[0] > 0.0) {
612                    scale[0] = (size[1] * aspect) / size[0];
613                } else if (size[1] > 0.0) {
614                    scale[1] = (size[0] / aspect) / size[1];
615                }
616            }
617        }
618            break;
619        }
620    }
621
622    TRACE("%s dims %g,%g,%g", _dataSet->getName().c_str(),
623          size[0], size[1], size[2]);
624    TRACE("Setting scale to %g,%g,%g", scale[0], scale[1], scale[2]);
625    setScale(scale);
626}
627
628vtkAlgorithmOutput *HeightMap::initWarp(vtkAlgorithmOutput *input)
629{
630    TRACE("Warp scale: %g", _warpScale);
631    if (_warpScale == 0.0) {
632        _warp = NULL;
633        _normalsGenerator = NULL;
634        return input;
635    } else if (input == NULL) {
636        ERROR("NULL input");
637        return input;
638    } else {
639        if (_warp == NULL)
640            _warp = vtkSmartPointer<vtkWarpScalar>::New();
641        if (_normalsGenerator == NULL)
642            _normalsGenerator = vtkSmartPointer<vtkPolyDataNormals>::New();
643        switch (_sliceAxis) {
644        case X_AXIS:
645            _warp->SetNormal(1, 0, 0);
646            break;
647        case Y_AXIS:
648            _warp->SetNormal(0, 1, 0);
649            break;
650        default:
651            _warp->SetNormal(0, 0, 1);
652        }
653        _warp->UseNormalOn();
654        _warp->SetScaleFactor(_warpScale * _dataScale);
655        _warp->SetInputConnection(input);
656        _normalsGenerator->SetInputConnection(_warp->GetOutputPort());
657        _normalsGenerator->SetFeatureAngle(90.);
658        _normalsGenerator->AutoOrientNormalsOff();
659        _normalsGenerator->ComputePointNormalsOn();
660        return _normalsGenerator->GetOutputPort();
661    }
662}
663
664/**
665 * \brief Controls relative scaling of height of mountain plot
666 */
667void HeightMap::setHeightScale(double scale)
668{
669    if (_warpScale == scale)
670        return;
671
672    _warpScale = scale;
673    if (_warp == NULL && scale != 0.0) {
674        vtkAlgorithmOutput *warpOutput = initWarp(_mapper->GetInputConnection(0, 0));
675        _mapper->SetInputConnection(warpOutput);
676        _contourFilter->SetInputConnection(warpOutput);
677    } else if (scale == 0.0) {
678        vtkAlgorithmOutput *warpInput = _warp->GetInputConnection(0, 0);
679        _mapper->SetInputConnection(warpInput);
680        _contourFilter->SetInputConnection(warpInput);
681        _warp = NULL;
682    } else {
683        _warp->SetScaleFactor(_warpScale * _dataScale);
684    }
685#ifdef TRANSLATE_TO_ZORIGIN
686    double pos[3];
687    pos[0] = 0;
688    pos[1] = 0;
689    pos[2] = - _warpScale * _dataScale * _dataRange[0];
690    setPosition(pos);
691    TRACE("Z_POS: %g", pos[2]);
692#endif
693    if (_mapper != NULL)
694        _mapper->Update();
695    if (_contourMapper != NULL)
696        _contourMapper->Update();
697}
698
699/**
700 * \brief Select a 2D slice plane from a 3D DataSet
701 *
702 * \param[in] axis Axis of slice plane
703 * \param[in] ratio Position [0,1] of slice plane along axis
704 */
705void HeightMap::selectVolumeSlice(Axis axis, double ratio)
706{
707    if (_dataSet->is2D()) {
708        WARN("DataSet not 3D, returning");
709        return;
710    }
711
712    if (_volumeSlicer == NULL &&
713        _cutPlane == NULL) {
714        WARN("Called before update() or DataSet is not a volume");
715        return;
716    }
717
718    _sliceAxis = axis;
719    if (_warp != NULL) {
720        switch (axis) {
721        case X_AXIS:
722            _warp->SetNormal(1, 0, 0);
723            break;
724        case Y_AXIS:
725            _warp->SetNormal(0, 1, 0);
726            break;
727        case Z_AXIS:
728            _warp->SetNormal(0, 0, 1);
729            break;
730        default:
731            ERROR("Invalid Axis");
732            return;
733        }
734    }
735
736    if (_cutPlane != NULL) {
737        double bounds[6];
738        _dataSet->getBounds(bounds);
739        switch (axis) {
740        case X_AXIS:
741            _cutPlane->SetNormal(1, 0, 0);
742            _cutPlane->SetOrigin(bounds[0] + (bounds[1]-bounds[0]) * ratio,
743                                 0,
744                                 0);
745            break;
746        case Y_AXIS:
747            _cutPlane->SetNormal(0, 1, 0);
748            _cutPlane->SetOrigin(0,
749                                 bounds[2] + (bounds[3]-bounds[2]) * ratio,
750                                 0);
751            break;
752        case Z_AXIS:
753            _cutPlane->SetNormal(0, 0, 1);
754            _cutPlane->SetOrigin(0,
755                                 0,
756                                 bounds[4] + (bounds[5]-bounds[4]) * ratio);
757            break;
758        default:
759            ERROR("Invalid Axis");
760            return;
761        }
762    } else {
763        int dims[3];
764        if (_splatter != NULL) {
765            _splatter->GetSampleDimensions(dims);
766        } else {
767            vtkImageData *imageData = vtkImageData::SafeDownCast(_dataSet->getVtkDataSet());
768            if (imageData == NULL) {
769                ERROR("Not a volume data set");
770                return;
771            }
772            imageData->GetDimensions(dims);
773        }
774        int voi[6];
775
776        switch (axis) {
777        case X_AXIS:
778            voi[0] = voi[1] = (int)((dims[0]-1) * ratio);
779            voi[2] = 0;
780            voi[3] = dims[1]-1;
781            voi[4] = 0;
782            voi[5] = dims[2]-1;
783            break;
784        case Y_AXIS:
785            voi[0] = 0;
786            voi[1] = dims[0]-1;
787            voi[2] = voi[3] = (int)((dims[1]-1) * ratio);
788            voi[4] = 0;
789            voi[5] = dims[2]-1;
790            break;
791        case Z_AXIS:
792            voi[0] = 0;
793            voi[1] = dims[0]-1;
794            voi[2] = 0;
795            voi[3] = dims[1]-1;
796            voi[4] = voi[5] = (int)((dims[2]-1) * ratio);
797            break;
798        default:
799            ERROR("Invalid Axis");
800            return;
801        }
802
803        _volumeSlicer->SetVOI(voi);
804    }
805
806    if (_mapper != NULL)
807        _mapper->Update();
808    if (_contourMapper != NULL)
809        _contourMapper->Update();
810}
811
812void HeightMap::updateRanges(Renderer *renderer)
813{
814    TRACE("Enter");
815
816    if (_dataSet == NULL) {
817        ERROR("called before setDataSet");
818        return;
819    }
820
821    if (renderer->getUseCumulativeRange()) {
822        renderer->getCumulativeDataRange(_dataRange,
823                                         _dataSet->getActiveScalarsName(),
824                                         1);
825        const char *activeVectors = _dataSet->getActiveVectorsName();
826        if (activeVectors != NULL) {
827            renderer->getCumulativeDataRange(_vectorMagnitudeRange,
828                                             activeVectors,
829                                             3);
830            for (int i = 0; i < 3; i++) {
831                renderer->getCumulativeDataRange(_vectorComponentRange[i],
832                                                 activeVectors,
833                                                 3, i);
834            }
835        }
836    } else {
837        _dataSet->getScalarRange(_dataRange);
838        _dataSet->getVectorRange(_vectorMagnitudeRange);
839        for (int i = 0; i < 3; i++) {
840            _dataSet->getVectorRange(_vectorComponentRange[i], i);
841        }
842    }
843 
844    // Need to update color map ranges
845    double *rangePtr = _colorFieldRange;
846    if (_colorFieldRange[0] > _colorFieldRange[1]) {
847        rangePtr = NULL;
848    }
849    setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
850
851    if ((_contours.empty() && _numContours > 0) ||
852        _warpScale > 0.0) {
853        // Contour isovalues and/or warp need to be recomputed
854        update();
855    } else {
856        computeDataScale();
857    }
858
859    TRACE("Leave");
860}
861
862void HeightMap::setContourLineColorMapEnabled(bool mode)
863{
864    _contourColorMap = mode;
865
866    double *rangePtr = _colorFieldRange;
867    if (_colorFieldRange[0] > _colorFieldRange[1]) {
868        rangePtr = NULL;
869    }
870    setColorMode(_colorMode, _colorFieldType, _colorFieldName.c_str(), rangePtr);
871}
872
873void HeightMap::setColorMode(ColorMode mode)
874{
875    _colorMode = mode;
876    if (_dataSet == NULL)
877        return;
878
879    switch (mode) {
880    case COLOR_BY_SCALAR:
881        setColorMode(mode,
882                     _dataSet->getActiveScalarsType(),
883                     _dataSet->getActiveScalarsName());
884        break;
885    case COLOR_BY_VECTOR_MAGNITUDE:
886        setColorMode(mode,
887                     _dataSet->getActiveVectorsType(),
888                     _dataSet->getActiveVectorsName());
889        break;
890    case COLOR_BY_VECTOR_X:
891        setColorMode(mode,
892                     _dataSet->getActiveVectorsType(),
893                     _dataSet->getActiveVectorsName());
894        break;
895    case COLOR_BY_VECTOR_Y:
896        setColorMode(mode,
897                     _dataSet->getActiveVectorsType(),
898                     _dataSet->getActiveVectorsName());
899        break;
900    case COLOR_BY_VECTOR_Z:
901        setColorMode(mode,
902                     _dataSet->getActiveVectorsType(),
903                     _dataSet->getActiveVectorsName());
904        break;
905    case COLOR_CONSTANT:
906    default:
907        setColorMode(mode, DataSet::POINT_DATA, NULL, NULL);
908        break;
909    }
910}
911
912void HeightMap::setColorMode(ColorMode mode,
913                             const char *name, double range[2])
914{
915    if (_dataSet == NULL)
916        return;
917    DataSet::DataAttributeType type = DataSet::POINT_DATA;
918    int numComponents = 1;
919    if (name != NULL && strlen(name) > 0 &&
920        !_dataSet->getFieldInfo(name, &type, &numComponents)) {
921        ERROR("Field not found: %s", name);
922        return;
923    }
924    setColorMode(mode, type, name, range);
925}
926
927void HeightMap::setColorMode(ColorMode mode, DataSet::DataAttributeType type,
928                             const char *name, double range[2])
929{
930    _colorMode = mode;
931    _colorFieldType = type;
932    if (name == NULL)
933        _colorFieldName.clear();
934    else
935        _colorFieldName = name;
936    if (range == NULL) {
937        _colorFieldRange[0] = DBL_MAX;
938        _colorFieldRange[1] = -DBL_MAX;
939    } else {
940        memcpy(_colorFieldRange, range, sizeof(double)*2);
941    }
942
943    if (_dataSet == NULL || _mapper == NULL)
944        return;
945
946    switch (type) {
947    case DataSet::POINT_DATA:
948        _mapper->SetScalarModeToUsePointFieldData();
949        _contourMapper->SetScalarModeToUsePointFieldData();
950        break;
951    case DataSet::CELL_DATA:
952        _mapper->SetScalarModeToUseCellFieldData();
953        _contourMapper->SetScalarModeToUseCellFieldData();
954        break;
955    default:
956        ERROR("Unsupported DataAttributeType: %d", type);
957        return;
958    }
959
960    if (_splatter != NULL) {
961        if (name != NULL && strlen(name) > 0) {
962            _splatter->SetInputArrayToProcess(0, 0, 0, vtkDataObject::FIELD_ASSOCIATION_POINTS, name);
963        }
964        _mapper->SelectColorArray("SplatterValues");
965        _contourMapper->SelectColorArray("SplatterValues");
966    } else if (name != NULL && strlen(name) > 0) {
967        _mapper->SelectColorArray(name);
968        _contourMapper->SelectColorArray(name);
969    } else {
970        _mapper->SetScalarModeToDefault();
971        _contourMapper->SetScalarModeToDefault();
972    }
973
974    if (_lut != NULL) {
975        if (range != NULL) {
976            _lut->SetRange(range);
977            TRACE("range %g, %g", range[0], range[1]);
978        } else if (name != NULL && strlen(name) > 0) {
979            double r[2];
980            int comp = -1;
981            if (mode == COLOR_BY_VECTOR_X)
982                comp = 0;
983            else if (mode == COLOR_BY_VECTOR_Y)
984                comp = 1;
985            else if (mode == COLOR_BY_VECTOR_Z)
986                comp = 2;
987
988            if (_renderer->getUseCumulativeRange()) {
989                int numComponents;
990                if  (!_dataSet->getFieldInfo(name, type, &numComponents)) {
991                    ERROR("Field not found: %s, type: %d", name, type);
992                    return;
993                } else if (numComponents < comp+1) {
994                    ERROR("Request for component %d in field with %d components",
995                          comp, numComponents);
996                    return;
997                }
998                _renderer->getCumulativeDataRange(r, name, type, numComponents, comp);
999            } else {
1000                _dataSet->getDataRange(r, name, type, comp);
1001            }
1002            _lut->SetRange(r);
1003            TRACE("%s range %g, %g", name, r[0], r[1]);
1004        }
1005        _lut->Modified();
1006    }
1007
1008    if (_contourColorMap) {
1009        _contourMapper->ScalarVisibilityOn();
1010    } else {
1011        _contourMapper->ScalarVisibilityOff();
1012    }
1013
1014    switch (mode) {
1015    case COLOR_BY_SCALAR:
1016        _mapper->ScalarVisibilityOn();
1017        break;
1018    case COLOR_BY_VECTOR_MAGNITUDE:
1019        _mapper->ScalarVisibilityOn();
1020        if (_lut != NULL) {
1021            _lut->SetVectorModeToMagnitude();
1022        }
1023        break;
1024    case COLOR_BY_VECTOR_X:
1025        _mapper->ScalarVisibilityOn();
1026        if (_lut != NULL) {
1027            _lut->SetVectorModeToComponent();
1028            _lut->SetVectorComponent(0);
1029        }
1030        break;
1031    case COLOR_BY_VECTOR_Y:
1032        _mapper->ScalarVisibilityOn();
1033        if (_lut != NULL) {
1034            _lut->SetVectorModeToComponent();
1035            _lut->SetVectorComponent(1);
1036        }
1037        break;
1038    case COLOR_BY_VECTOR_Z:
1039        _mapper->ScalarVisibilityOn();
1040        if (_lut != NULL) {
1041            _lut->SetVectorModeToComponent();
1042            _lut->SetVectorComponent(2);
1043        }
1044        break;
1045    case COLOR_CONSTANT:
1046    default:
1047        _mapper->ScalarVisibilityOff();
1048        break;
1049    }
1050}
1051
1052/**
1053 * \brief Called when the color map has been edited
1054 */
1055void HeightMap::updateColorMap()
1056{
1057    setColorMap(_colorMap);
1058}
1059
1060/**
1061 * \brief Associate a colormap lookup table with the DataSet
1062 */
1063void HeightMap::setColorMap(ColorMap *cmap)
1064{
1065    if (cmap == NULL)
1066        return;
1067
1068    _colorMap = cmap;
1069
1070    if (_lut == NULL) {
1071        _lut = vtkSmartPointer<vtkLookupTable>::New();
1072        if (_mapper != NULL) {
1073            _mapper->UseLookupTableScalarRangeOn();
1074            _mapper->SetLookupTable(_lut);
1075        }
1076        if (_contourMapper != NULL) {
1077            _contourMapper->UseLookupTableScalarRangeOn();
1078            _contourMapper->SetLookupTable(_lut);
1079        }
1080        _lut->DeepCopy(cmap->getLookupTable());
1081        switch (_colorMode) {
1082        case COLOR_CONSTANT:
1083        case COLOR_BY_SCALAR:
1084            _lut->SetRange(_dataRange);
1085            break;
1086        case COLOR_BY_VECTOR_MAGNITUDE:
1087            _lut->SetRange(_vectorMagnitudeRange);
1088            break;
1089        case COLOR_BY_VECTOR_X:
1090            _lut->SetRange(_vectorComponentRange[0]);
1091            break;
1092        case COLOR_BY_VECTOR_Y:
1093            _lut->SetRange(_vectorComponentRange[1]);
1094            break;
1095        case COLOR_BY_VECTOR_Z:
1096            _lut->SetRange(_vectorComponentRange[2]);
1097            break;
1098        default:
1099            break;
1100        }
1101    } else {
1102        double range[2];
1103        _lut->GetTableRange(range);
1104        _lut->DeepCopy(cmap->getLookupTable());
1105        _lut->SetRange(range);
1106        _lut->Modified();
1107    }
1108
1109    switch (_colorMode) {
1110    case COLOR_BY_VECTOR_MAGNITUDE:
1111        _lut->SetVectorModeToMagnitude();
1112        break;
1113    case COLOR_BY_VECTOR_X:
1114        _lut->SetVectorModeToComponent();
1115        _lut->SetVectorComponent(0);
1116        break;
1117    case COLOR_BY_VECTOR_Y:
1118        _lut->SetVectorModeToComponent();
1119        _lut->SetVectorComponent(1);
1120        break;
1121    case COLOR_BY_VECTOR_Z:
1122        _lut->SetVectorModeToComponent();
1123        _lut->SetVectorComponent(2);
1124        break;
1125    default:
1126         break;
1127    }
1128}
1129
1130/**
1131 * \brief Specify number of evenly spaced contour lines to render
1132 *
1133 * Will override any existing contours
1134 */
1135void HeightMap::setNumContours(int numContours)
1136{
1137    _contours.clear();
1138    _numContours = numContours;
1139
1140    update();
1141}
1142
1143/**
1144 * \brief Specify an explicit list of contour isovalues
1145 *
1146 * Will override any existing contours
1147 */
1148void HeightMap::setContourList(const std::vector<double>& contours)
1149{
1150    _contours = contours;
1151    _numContours = (int)_contours.size();
1152
1153    update();
1154}
1155
1156/**
1157 * \brief Get the number of contours
1158 */
1159int HeightMap::getNumContours() const
1160{
1161    return _numContours;
1162}
1163
1164/**
1165 * \brief Get the contour list (may be empty if number of contours
1166 * was specified in place of a list)
1167 */
1168const std::vector<double>& HeightMap::getContourList() const
1169{
1170    return _contours;
1171}
1172
1173/**
1174 * \brief Turn on/off lighting of this object
1175 */
1176void HeightMap::setLighting(bool state)
1177{
1178    _lighting = state;
1179    if (_dsActor != NULL)
1180        _dsActor->GetProperty()->SetLighting((state ? 1 : 0));
1181}
1182
1183/**
1184 * \brief Turn on/off rendering of mesh edges
1185 */
1186void HeightMap::setEdgeVisibility(bool state)
1187{
1188    if (_dsActor != NULL) {
1189        _dsActor->GetProperty()->SetEdgeVisibility((state ? 1 : 0));
1190    }
1191}
1192
1193/**
1194 * \brief Turn on/off rendering of colormaped surface
1195 */
1196void HeightMap::setContourSurfaceVisibility(bool state)
1197{
1198    if (_dsActor != NULL) {
1199        _dsActor->SetVisibility((state ? 1 : 0));
1200    }
1201}
1202
1203/**
1204 * \brief Turn on/off rendering of contour isolines
1205 */
1206void HeightMap::setContourLineVisibility(bool state)
1207{
1208    if (_contourActor != NULL) {
1209        _contourActor->SetVisibility((state ? 1 : 0));
1210    }
1211}
1212
1213/**
1214 * \brief Set RGB color of mesh
1215 */
1216void HeightMap::setColor(float color[3])
1217{
1218    _color[0] = color[0];
1219    _color[1] = color[1];
1220    _color[2] = color[2];
1221    if (_dsActor != NULL)
1222        _dsActor->GetProperty()->SetColor(_color[0], _color[1], _color[2]);
1223}
1224
1225/**
1226 * \brief Set RGB color of polygon edges
1227 */
1228void HeightMap::setEdgeColor(float color[3])
1229{
1230    _edgeColor[0] = color[0];
1231    _edgeColor[1] = color[1];
1232    _edgeColor[2] = color[2];
1233    if (_dsActor != NULL)
1234        _dsActor->GetProperty()->SetEdgeColor(_edgeColor[0], _edgeColor[1], _edgeColor[2]);
1235}
1236
1237/**
1238 * \brief Set pixel width of polygon edges (may be a no-op)
1239 */
1240void HeightMap::setEdgeWidth(float edgeWidth)
1241{
1242    _edgeWidth = edgeWidth;
1243    if (_dsActor != NULL)
1244        _dsActor->GetProperty()->SetLineWidth(_edgeWidth);
1245}
1246
1247/**
1248 * \brief Set RGB color of contour isolines
1249 */
1250void HeightMap::setContourEdgeColor(float color[3])
1251{
1252    _contourEdgeColor[0] = color[0];
1253    _contourEdgeColor[1] = color[1];
1254    _contourEdgeColor[2] = color[2];
1255    if (_contourActor != NULL) {
1256        _contourActor->GetProperty()->SetColor(_contourEdgeColor[0],
1257                                               _contourEdgeColor[1],
1258                                               _contourEdgeColor[2]);
1259        _contourActor->GetProperty()->SetEdgeColor(_contourEdgeColor[0],
1260                                                   _contourEdgeColor[1],
1261                                                   _contourEdgeColor[2]);
1262    }
1263}
1264
1265/**
1266 * \brief Set pixel width of contour isolines (may be a no-op)
1267 */
1268void HeightMap::setContourEdgeWidth(float edgeWidth)
1269{
1270    _contourEdgeWidth = edgeWidth;
1271    if (_contourActor != NULL)
1272        _contourActor->GetProperty()->SetLineWidth(_contourEdgeWidth);
1273}
1274
1275/**
1276 * \brief Set a group of world coordinate planes to clip rendering
1277 *
1278 * Passing NULL for planes will remove all cliping planes
1279 */
1280void HeightMap::setClippingPlanes(vtkPlaneCollection *planes)
1281{
1282    if (_mapper != NULL) {
1283        _mapper->SetClippingPlanes(planes);
1284    }
1285    if (_contourMapper != NULL) {
1286        _contourMapper->SetClippingPlanes(planes);
1287    }
1288}
1289
Note: See TracBrowser for help on using the repository browser.