source: branches/blt4/packages/vizservers/vtkvis/RpVtkRenderer.cpp @ 2681

Last change on this file since 2681 was 2681, checked in by gah, 11 years ago
File size: 109.6 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 <cfloat>
9#include <cstring>
10#include <cassert>
11#include <cmath>
12
13#include <GL/gl.h>
14
15#ifdef WANT_TRACE
16#include <sys/time.h>
17#endif
18
19#include <vtkMath.h>
20#include <vtkCamera.h>
21#include <vtkLight.h>
22#include <vtkCoordinate.h>
23#include <vtkTransform.h>
24#include <vtkCharArray.h>
25#include <vtkAxisActor2D.h>
26#include <vtkCubeAxesActor.h>
27#ifdef USE_CUSTOM_AXES
28#include "vtkRpCubeAxesActor2D.h"
29#else
30#include <vtkCubeAxesActor2D.h>
31#endif
32#include <vtkDataSetReader.h>
33#include <vtkDataSetMapper.h>
34#include <vtkContourFilter.h>
35#include <vtkPolyDataMapper.h>
36#include <vtkProperty.h>
37#include <vtkProperty2D.h>
38#include <vtkPointData.h>
39#include <vtkLookupTable.h>
40#include <vtkTextProperty.h>
41#include <vtkOpenGLRenderWindow.h>
42#include <vtkVersion.h>
43
44#include "RpVtkRenderer.h"
45#include "RpVtkRendererGraphicsObjs.h"
46#include "RpMath.h"
47#include "ColorMap.h"
48#include "Trace.h"
49
50#define ELAPSED_TIME(t1, t2) \
51    ((t1).tv_sec == (t2).tv_sec ? (((t2).tv_usec - (t1).tv_usec)/1.0e+3f) : \
52     (((t2).tv_sec - (t1).tv_sec))*1.0e+3f + (float)((t2).tv_usec - (t1).tv_usec)/1.0e+3f)
53
54using namespace Rappture::VtkVis;
55
56Renderer::Renderer() :
57    _needsRedraw(true),
58    _windowWidth(500),
59    _windowHeight(500),
60    _imgCameraPlane(PLANE_XY),
61    _imgCameraOffset(0),
62    _cameraZoomRatio(1),
63    _useCumulativeRange(true),
64    _cumulativeRangeOnlyVisible(false),
65    _cameraMode(PERSPECTIVE)
66{
67    _bgColor[0] = 0;
68    _bgColor[1] = 0;
69    _bgColor[2] = 0;
70    _cameraPan[0] = 0;
71    _cameraPan[1] = 0;
72    _cameraOrientation[0] = 1.0;
73    _cameraOrientation[1] = 0.0;
74    _cameraOrientation[2] = 0.0;
75    _cameraOrientation[3] = 0.0;
76    // clipping planes to prevent overdrawing axes
77    _activeClipPlanes = vtkSmartPointer<vtkPlaneCollection>::New();
78    // bottom
79    _cameraClipPlanes[0] = vtkSmartPointer<vtkPlane>::New();
80    _cameraClipPlanes[0]->SetNormal(0, 1, 0);
81    _cameraClipPlanes[0]->SetOrigin(0, 0, 0);
82    if (_cameraMode == IMAGE)
83        _activeClipPlanes->AddItem(_cameraClipPlanes[0]);
84    // left
85    _cameraClipPlanes[1] = vtkSmartPointer<vtkPlane>::New();
86    _cameraClipPlanes[1]->SetNormal(1, 0, 0);
87    _cameraClipPlanes[1]->SetOrigin(0, 0, 0);
88    if (_cameraMode == IMAGE)
89        _activeClipPlanes->AddItem(_cameraClipPlanes[1]);
90    // top
91    _cameraClipPlanes[2] = vtkSmartPointer<vtkPlane>::New();
92    _cameraClipPlanes[2]->SetNormal(0, -1, 0);
93    _cameraClipPlanes[2]->SetOrigin(0, 1, 0);
94    if (_cameraMode == IMAGE)
95        _activeClipPlanes->AddItem(_cameraClipPlanes[2]);
96    // right
97    _cameraClipPlanes[3] = vtkSmartPointer<vtkPlane>::New();
98    _cameraClipPlanes[3]->SetNormal(-1, 0, 0);
99    _cameraClipPlanes[3]->SetOrigin(1, 0, 0);
100    if (_cameraMode == IMAGE)
101        _activeClipPlanes->AddItem(_cameraClipPlanes[3]);
102    _renderer = vtkSmartPointer<vtkRenderer>::New();
103    vtkSmartPointer<vtkLight> headlight = vtkSmartPointer<vtkLight>::New();
104    headlight->SetLightTypeToHeadlight();
105    headlight->PositionalOff();
106    //headlight->SetAmbientColor(1, 1, 1);
107    _renderer->SetAmbient(.2, .2, .2);
108    _renderer->AddLight(headlight);
109    vtkSmartPointer<vtkLight> skylight = vtkSmartPointer<vtkLight>::New();
110    skylight->SetLightTypeToCameraLight();
111    skylight->SetPosition(0, 1, 0);
112    skylight->SetFocalPoint(0, 0, 0);
113    skylight->PositionalOff();
114    //skylight->SetAmbientColor(1, 1, 1);
115    _renderer->AddLight(skylight);
116    _renderer->LightFollowCameraOn();
117    _renderWindow = vtkSmartPointer<vtkRenderWindow>::New();
118#ifdef USE_OFFSCREEN_RENDERING
119    _renderWindow->DoubleBufferOff();
120    _renderWindow->OffScreenRenderingOn();
121#else
122    _renderWindow->SwapBuffersOff();
123#endif
124    _renderWindow->SetSize(_windowWidth, _windowHeight);
125    // Next 2 options needed to support depth peeling
126    _renderWindow->SetAlphaBitPlanes(1);
127    _renderWindow->SetMultiSamples(0);
128    _renderer->SetMaximumNumberOfPeels(100);
129    _renderer->SetUseDepthPeeling(1);
130    _renderWindow->AddRenderer(_renderer);
131    setViewAngle(_windowHeight);
132    initAxes();
133    initCamera();
134    addColorMap("default", ColorMap::getDefault());
135    addColorMap("grayDefault", ColorMap::getGrayDefault());
136    addColorMap("volumeDefault", ColorMap::getVolumeDefault());
137    addColorMap("elementDefault", ColorMap::getElementDefault());
138}
139
140Renderer::~Renderer()
141{
142    TRACE("Enter");
143    TRACE("Deleting Contour2Ds");
144    for (Contour2DHashmap::iterator itr = _contour2Ds.begin();
145         itr != _contour2Ds.end(); ++itr) {
146        delete itr->second;
147    }
148    _contour2Ds.clear();
149    TRACE("Deleting Contour3Ds");
150    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
151         itr != _contour3Ds.end(); ++itr) {
152        delete itr->second;
153    }
154    _contour3Ds.clear();
155    TRACE("Deleting Cutplanes");
156    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
157         itr != _cutplanes.end(); ++itr) {
158        delete itr->second;
159    }
160    _cutplanes.clear();
161    TRACE("Deleting Glyphs");
162    for (GlyphsHashmap::iterator itr = _glyphs.begin();
163         itr != _glyphs.end(); ++itr) {
164        delete itr->second;
165    }
166    _glyphs.clear();
167    TRACE("Deleting HeightMaps");
168    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
169         itr != _heightMaps.end(); ++itr) {
170        delete itr->second;
171    }
172    _heightMaps.clear();
173    TRACE("Deleting LICs");
174    for (LICHashmap::iterator itr = _lics.begin();
175         itr != _lics.end(); ++itr) {
176        delete itr->second;
177    }
178    _lics.clear();
179    TRACE("Deleting Molecules");
180    for (MoleculeHashmap::iterator itr = _molecules.begin();
181         itr != _molecules.end(); ++itr) {
182        delete itr->second;
183    }
184    _molecules.clear();
185    TRACE("Deleting PolyDatas");
186    for (PolyDataHashmap::iterator itr = _polyDatas.begin();
187         itr != _polyDatas.end(); ++itr) {
188        delete itr->second;
189    }
190    _polyDatas.clear();
191    TRACE("Deleting PseudoColors");
192    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
193         itr != _pseudoColors.end(); ++itr) {
194        delete itr->second;
195    }
196    _pseudoColors.clear();
197    TRACE("Deleting Streamlines");
198    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
199         itr != _streamlines.end(); ++itr) {
200        delete itr->second;
201    }
202    _streamlines.clear();
203    TRACE("Deleting Volumes");
204    for (VolumeHashmap::iterator itr = _volumes.begin();
205         itr != _volumes.end(); ++itr) {
206        delete itr->second;
207    }
208    _volumes.clear();
209    TRACE("Deleting ColorMaps");
210    // Delete color maps and data sets last in case references still
211    // exist
212    for (ColorMapHashmap::iterator itr = _colorMaps.begin();
213         itr != _colorMaps.end(); ++itr) {
214        delete itr->second;
215    }
216    _colorMaps.clear();
217    TRACE("Deleting DataSets");
218    for (DataSetHashmap::iterator itr = _dataSets.begin();
219         itr != _dataSets.end(); ++itr) {
220        delete itr->second;
221    }
222    _dataSets.clear();
223
224    clearFieldRanges();
225
226    TRACE("Leave");
227}
228
229/**
230 * \brief Add a DataSet to this Renderer
231 *
232 * This just adds the DataSet to the Renderer's list of data sets.
233 * In order to render the data, a graphics object using the data
234 * set must be added to the Renderer.
235 */
236void Renderer::addDataSet(const DataSetId& id)
237{
238    if (getDataSet(id) != NULL) {
239        WARN("Replacing existing DataSet %s", id.c_str());
240        deleteDataSet(id);
241    }
242    _dataSets[id] = new DataSet(id);
243}
244
245/**
246 * \brief Remove the specified DataSet and associated rendering objects
247 *
248 * The underlying DataSet and any associated graphics
249 * objects are deleted, freeing the memory used.
250 */
251void Renderer::deleteDataSet(const DataSetId& id)
252{
253    DataSetHashmap::iterator itr;
254
255    bool doAll = false;
256
257    if (id.compare("all") == 0) {
258        itr = _dataSets.begin();
259        doAll = true;
260    } else {
261        itr = _dataSets.find(id);
262    }
263    if (itr == _dataSets.end()) {
264        ERROR("Unknown dataset %s", id.c_str());
265        return;
266    }
267
268    do {
269        TRACE("Deleting dataset %s", itr->second->getName().c_str());
270
271        deleteGraphicsObject<Contour2D>(itr->second->getName());
272        deleteGraphicsObject<Contour3D>(itr->second->getName());
273        deleteGraphicsObject<Cutplane>(itr->second->getName());
274        deleteGraphicsObject<Glyphs>(itr->second->getName());
275        deleteGraphicsObject<HeightMap>(itr->second->getName());
276        deleteGraphicsObject<LIC>(itr->second->getName());
277        deleteGraphicsObject<Molecule>(itr->second->getName());
278        deleteGraphicsObject<PolyData>(itr->second->getName());
279        deleteGraphicsObject<PseudoColor>(itr->second->getName());
280        deleteGraphicsObject<Streamlines>(itr->second->getName());
281        deleteGraphicsObject<Volume>(itr->second->getName());
282 
283        if (itr->second->getProp() != NULL) {
284            _renderer->RemoveViewProp(itr->second->getProp());
285        }
286
287        TRACE("After deleting graphics objects");
288
289        delete itr->second;
290        itr = _dataSets.erase(itr);
291    } while (doAll && itr != _dataSets.end());
292
293    // Update cumulative data range
294    initFieldRanges();
295    updateFieldRanges();
296
297    initCamera();
298    _needsRedraw = true;
299}
300
301/**
302 * \brief Get a list of DataSets this Renderer knows about
303 */
304void Renderer::getDataSetNames(std::vector<std::string>& names)
305{
306    names.clear();
307    for (DataSetHashmap::iterator itr = _dataSets.begin();
308         itr != _dataSets.end(); ++itr) {
309        names.push_back(itr->second->getName());
310    }
311}
312
313/**
314 * \brief Find the DataSet for the given DataSetId key
315 *
316 * \return A pointer to the DataSet, or NULL if not found
317 */
318DataSet *Renderer::getDataSet(const DataSetId& id)
319{
320    DataSetHashmap::iterator itr = _dataSets.find(id);
321    if (itr == _dataSets.end()) {
322#ifdef DEBUG
323        TRACE("DataSet not found: %s", id.c_str());
324#endif
325        return NULL;
326    } else
327        return itr->second;
328}
329
330/**
331 * \brief (Re-)load the data for the specified DataSet key from a file
332 */
333bool Renderer::setDataFile(const DataSetId& id, const char *filename)
334{
335    DataSet *ds = getDataSet(id);
336    if (ds) {
337        bool ret = ds->setDataFile(filename);
338        initFieldRanges();
339        updateFieldRanges();
340        _needsRedraw = true;
341        return ret;
342    } else
343        return false;
344}
345
346/**
347 * \brief (Re-)load the data for the specified DataSet key from a memory buffer
348 */
349bool Renderer::setData(const DataSetId& id, char *data, int nbytes)
350{
351    DataSet *ds = getDataSet(id);
352    if (ds) {
353        bool ret = ds->setData(data, nbytes);
354        initFieldRanges();
355        updateFieldRanges();
356        _needsRedraw = true;
357        return ret;
358    } else
359        return false;
360}
361
362/**
363 * \brief Set the active scalar field array by name for a DataSet
364 */
365bool Renderer::setDataSetActiveScalars(const DataSetId& id, const char *scalarName)
366{
367    DataSetHashmap::iterator itr;
368
369    bool doAll = false;
370
371    if (id.compare("all") == 0) {
372        itr = _dataSets.begin();
373        doAll = true;
374    } else {
375        itr = _dataSets.find(id);
376    }
377    if (itr == _dataSets.end()) {
378        ERROR("DataSet not found: %s", id.c_str());
379        return false;
380    }
381
382    bool ret = true;
383    bool needRangeUpdate = false;
384    do {
385        if (strcmp(itr->second->getActiveScalarsName(), scalarName) != 0) {
386            if (!itr->second->setActiveScalars(scalarName)) {
387                ret = false;
388            } else {
389                needRangeUpdate = true;
390            }
391        }
392    } while (doAll && ++itr != _dataSets.end());
393
394    if (needRangeUpdate) {
395         updateFieldRanges();
396        _needsRedraw = true;
397    }
398
399    return ret;
400}
401
402/**
403 * \brief Set the active vector field array by name for a DataSet
404 */
405bool Renderer::setDataSetActiveVectors(const DataSetId& id, const char *vectorName)
406{
407    DataSetHashmap::iterator itr;
408
409    bool doAll = false;
410
411    if (id.compare("all") == 0) {
412        itr = _dataSets.begin();
413        doAll = true;
414    } else {
415        itr = _dataSets.find(id);
416    }
417    if (itr == _dataSets.end()) {
418        ERROR("DataSet not found: %s", id.c_str());
419        return false;
420    }
421
422    bool ret = true;
423    bool needRangeUpdate = false;
424    do {
425        if (strcmp(itr->second->getActiveVectorsName(), vectorName) != 0) {
426            if (!itr->second->setActiveVectors(vectorName)) {
427                ret = false;
428            } else {
429                needRangeUpdate = true;
430            }
431        }
432    } while (doAll && ++itr != _dataSets.end());
433
434    if (needRangeUpdate) {
435        updateFieldRanges();
436        _needsRedraw = true;
437    }
438
439    return ret;
440}
441
442/**
443 * \brief Control whether the cumulative data range of datasets is used for
444 * colormapping and contours
445 */
446void Renderer::setUseCumulativeDataRange(bool state, bool onlyVisible)
447{
448    if (state != _useCumulativeRange) {
449        _useCumulativeRange = state;
450        _cumulativeRangeOnlyVisible = onlyVisible;
451        updateFieldRanges();
452        _needsRedraw = true;
453    }
454}
455
456void Renderer::resetAxes(double bounds[6])
457{
458    TRACE("Resetting axes");
459    if (_cubeAxesActor == NULL ||
460        _cubeAxesActor2D == NULL) {
461        initAxes();
462    }
463    if (_cameraMode == IMAGE) {
464        if (_renderer->HasViewProp(_cubeAxesActor)) {
465            TRACE("Removing 3D axes");
466            _renderer->RemoveViewProp(_cubeAxesActor);
467        }
468        if (!_renderer->HasViewProp(_cubeAxesActor2D)) {
469            TRACE("Adding 2D axes");
470            _renderer->AddViewProp(_cubeAxesActor2D);
471        }
472    } else {
473        if (_renderer->HasViewProp(_cubeAxesActor2D)) {
474            TRACE("Removing 2D axes");
475            _renderer->RemoveViewProp(_cubeAxesActor2D);
476        }
477        if (!_renderer->HasViewProp(_cubeAxesActor)) {
478            TRACE("Adding 3D axes");
479            _renderer->AddViewProp(_cubeAxesActor);
480        }
481        if (bounds == NULL) {
482            double newBounds[6];
483            collectBounds(newBounds, false);
484#if ((VTK_MAJOR_VERSION > 5) || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION >= 8))
485            _cubeAxesActor->SetXAxisRange(newBounds[0], newBounds[1]);
486            _cubeAxesActor->SetYAxisRange(newBounds[2], newBounds[3]);
487            _cubeAxesActor->SetZAxisRange(newBounds[4], newBounds[5]);
488#endif
489            _cubeAxesActor->SetBounds(newBounds);
490            TRACE("Bounds (computed): %g %g %g %g %g %g",
491                  newBounds[0],
492                  newBounds[1],
493                  newBounds[2],
494                  newBounds[3],
495                  newBounds[4],
496                  newBounds[5]);
497        } else {
498#if ((VTK_MAJOR_VERSION > 5) || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION >= 8))
499            _cubeAxesActor->SetXAxisRange(bounds[0], bounds[1]);
500            _cubeAxesActor->SetYAxisRange(bounds[2], bounds[3]);
501            _cubeAxesActor->SetZAxisRange(bounds[4], bounds[5]);
502#endif
503            _cubeAxesActor->SetBounds(bounds);
504            TRACE("Bounds (supplied): %g %g %g %g %g %g",
505                  bounds[0],
506                  bounds[1],
507                  bounds[2],
508                  bounds[3],
509                  bounds[4],
510                  bounds[5]);
511        }
512    }
513}
514
515/**
516 * \brief Set inital properties on the 2D Axes
517 */
518void Renderer::initAxes()
519{
520    TRACE("Initializing axes");
521    if (_cubeAxesActor == NULL)
522        _cubeAxesActor = vtkSmartPointer<vtkCubeAxesActor>::New();
523    _cubeAxesActor->SetCamera(_renderer->GetActiveCamera());
524    _cubeAxesActor->GetProperty()->LightingOff();
525    // Don't offset labels at origin
526    _cubeAxesActor->SetCornerOffset(0);
527    _cubeAxesActor->SetFlyModeToStaticTriad();
528
529#ifdef USE_CUSTOM_AXES
530    if (_cubeAxesActor2D == NULL)
531        _cubeAxesActor2D = vtkSmartPointer<vtkRpCubeAxesActor2D>::New();
532#else
533    if (_cubeAxesActor2D == NULL)
534        _cubeAxesActor2D = vtkSmartPointer<vtkCubeAxesActor2D>::New();
535#endif
536
537    _cubeAxesActor2D->SetCamera(_renderer->GetActiveCamera());
538    _cubeAxesActor2D->ZAxisVisibilityOff();
539    _cubeAxesActor2D->SetCornerOffset(0);
540    _cubeAxesActor2D->SetFlyModeToClosestTriad();
541
542    _cubeAxesActor2D->ScalingOff();
543    //_cubeAxesActor2D->SetShowActualBounds(0);
544    _cubeAxesActor2D->SetFontFactor(1.25);
545    // Use "nice" range and number of ticks/labels
546    _cubeAxesActor2D->GetXAxisActor2D()->AdjustLabelsOn();
547    _cubeAxesActor2D->GetYAxisActor2D()->AdjustLabelsOn();
548    _cubeAxesActor2D->GetZAxisActor2D()->AdjustLabelsOn();
549
550#ifdef USE_CUSTOM_AXES
551    _cubeAxesActor2D->SetAxisTitleTextProperty(NULL);
552    _cubeAxesActor2D->SetAxisLabelTextProperty(NULL);
553    //_cubeAxesActor2D->GetXAxisActor2D()->SizeFontRelativeToAxisOn();
554    _cubeAxesActor2D->GetXAxisActor2D()->GetTitleTextProperty()->BoldOn();
555    _cubeAxesActor2D->GetXAxisActor2D()->GetTitleTextProperty()->ItalicOff();
556    _cubeAxesActor2D->GetXAxisActor2D()->GetTitleTextProperty()->ShadowOn();
557    _cubeAxesActor2D->GetXAxisActor2D()->GetLabelTextProperty()->BoldOff();
558    _cubeAxesActor2D->GetXAxisActor2D()->GetLabelTextProperty()->ItalicOff();
559    _cubeAxesActor2D->GetXAxisActor2D()->GetLabelTextProperty()->ShadowOff();
560
561    //_cubeAxesActor2D->GetYAxisActor2D()->SizeFontRelativeToAxisOn();
562    _cubeAxesActor2D->GetYAxisActor2D()->GetTitleTextProperty()->BoldOn();
563    _cubeAxesActor2D->GetYAxisActor2D()->GetTitleTextProperty()->ItalicOff();
564    _cubeAxesActor2D->GetYAxisActor2D()->GetTitleTextProperty()->ShadowOn();
565    _cubeAxesActor2D->GetYAxisActor2D()->GetLabelTextProperty()->BoldOff();
566    _cubeAxesActor2D->GetYAxisActor2D()->GetLabelTextProperty()->ItalicOff();
567    _cubeAxesActor2D->GetYAxisActor2D()->GetLabelTextProperty()->ShadowOff();
568
569    //_cubeAxesActor2D->GetZAxisActor2D()->SizeFontRelativeToAxisOn();
570    _cubeAxesActor2D->GetZAxisActor2D()->GetTitleTextProperty()->BoldOn();
571    _cubeAxesActor2D->GetZAxisActor2D()->GetTitleTextProperty()->ItalicOff();
572    _cubeAxesActor2D->GetZAxisActor2D()->GetTitleTextProperty()->ShadowOn();
573    _cubeAxesActor2D->GetZAxisActor2D()->GetLabelTextProperty()->BoldOff();
574    _cubeAxesActor2D->GetZAxisActor2D()->GetLabelTextProperty()->ItalicOff();
575    _cubeAxesActor2D->GetZAxisActor2D()->GetLabelTextProperty()->ShadowOff();
576#else
577    _cubeAxesActor2D->GetAxisTitleTextProperty()->BoldOn();
578    _cubeAxesActor2D->GetAxisTitleTextProperty()->ItalicOff();
579    _cubeAxesActor2D->GetAxisTitleTextProperty()->ShadowOn();
580    _cubeAxesActor2D->GetAxisLabelTextProperty()->BoldOff();
581    _cubeAxesActor2D->GetAxisLabelTextProperty()->ItalicOff();
582    _cubeAxesActor2D->GetAxisLabelTextProperty()->ShadowOff();
583#endif
584
585    if (_cameraMode == IMAGE) {
586        if (!_renderer->HasViewProp(_cubeAxesActor2D))
587            _renderer->AddViewProp(_cubeAxesActor2D);
588    } else {
589        if (!_renderer->HasViewProp(_cubeAxesActor))
590            _renderer->AddViewProp(_cubeAxesActor);
591    }
592}
593
594/**
595 * \brief Set Fly mode of axes
596 */
597void Renderer::setAxesFlyMode(AxesFlyMode mode)
598{
599    if (_cubeAxesActor == NULL)
600        initAxes();
601    switch (mode) {
602    case FLY_STATIC_EDGES:
603        _cubeAxesActor->SetFlyModeToStaticEdges();
604        break;
605    case FLY_STATIC_TRIAD:
606        _cubeAxesActor->SetFlyModeToStaticTriad();
607        break;
608    case FLY_OUTER_EDGES:
609        _cubeAxesActor->SetFlyModeToOuterEdges();
610        break;
611    case FLY_FURTHEST_TRIAD:
612        _cubeAxesActor->SetFlyModeToFurthestTriad();
613        break;
614    case FLY_CLOSEST_TRIAD:
615    default:
616        _cubeAxesActor->SetFlyModeToClosestTriad();
617        break;
618    }
619    _needsRedraw = true;
620}
621
622/**
623 * \brief Set color of axes, ticks, labels, titles
624 */
625void Renderer::setAxesColor(double color[3])
626{
627    if (_cubeAxesActor != NULL) {
628        _cubeAxesActor->GetProperty()->SetColor(color);
629        _needsRedraw = true;
630    }
631    if (_cubeAxesActor2D != NULL) {
632        _cubeAxesActor2D->GetProperty()->SetColor(color);
633#ifdef USE_CUSTOM_AXES
634        _cubeAxesActor2D->GetXAxisActor2D()->GetTitleTextProperty()->SetColor(color);
635        _cubeAxesActor2D->GetXAxisActor2D()->GetLabelTextProperty()->SetColor(color);
636        _cubeAxesActor2D->GetYAxisActor2D()->GetTitleTextProperty()->SetColor(color);
637        _cubeAxesActor2D->GetYAxisActor2D()->GetLabelTextProperty()->SetColor(color);
638        _cubeAxesActor2D->GetZAxisActor2D()->GetTitleTextProperty()->SetColor(color);
639        _cubeAxesActor2D->GetZAxisActor2D()->GetLabelTextProperty()->SetColor(color);
640#else
641        _cubeAxesActor2D->GetAxisTitleTextProperty()->SetColor(color);
642        _cubeAxesActor2D->GetAxisLabelTextProperty()->SetColor(color);
643#endif
644        _needsRedraw = true;
645    }
646}
647
648/**
649 * \brief Turn on/off rendering of all enabled axes
650 */
651void Renderer::setAxesVisibility(bool state)
652{
653    if (_cubeAxesActor != NULL) {
654        _cubeAxesActor->SetVisibility((state ? 1 : 0));
655        _needsRedraw = true;
656    }
657    if (_cubeAxesActor2D != NULL) {
658        _cubeAxesActor2D->SetVisibility((state ? 1 : 0));
659        _needsRedraw = true;
660    }
661}
662
663/**
664 * \brief Turn on/off rendering of all axes gridlines
665 */
666void Renderer::setAxesGridVisibility(bool state)
667{
668    setAxisGridVisibility(X_AXIS, state);
669    setAxisGridVisibility(Y_AXIS, state);
670    setAxisGridVisibility(Z_AXIS, state);
671}
672
673/**
674 * \brief Turn on/off rendering of all axis labels
675 */
676void Renderer::setAxesLabelVisibility(bool state)
677{
678    setAxisLabelVisibility(X_AXIS, state);
679    setAxisLabelVisibility(Y_AXIS, state);
680    setAxisLabelVisibility(Z_AXIS, state);
681}
682
683/**
684 * \brief Turn on/off rendering of all axis ticks
685 */
686void Renderer::setAxesTickVisibility(bool state)
687{
688    setAxisTickVisibility(X_AXIS, state);
689    setAxisTickVisibility(Y_AXIS, state);
690    setAxisTickVisibility(Z_AXIS, state);
691}
692
693/**
694 * \brief Control position of ticks on 3D axes
695 */
696void Renderer::setAxesTickPosition(AxesTickPosition pos)
697{
698    if (_cubeAxesActor == NULL)
699        return;
700
701    switch (pos) {
702    case TICKS_BOTH:
703        _cubeAxesActor->SetTickLocationToBoth();
704        break;
705    case TICKS_OUTSIDE:
706        _cubeAxesActor->SetTickLocationToOutside();
707        break;
708    case TICKS_INSIDE:
709    default:
710        _cubeAxesActor->SetTickLocationToInside();
711        break;
712    }
713    _needsRedraw = true;
714}
715
716/**
717 * \brief Turn on/off rendering of the specified axis
718 */
719void Renderer::setAxisVisibility(Axis axis, bool state)
720{
721    if (_cubeAxesActor != NULL) {
722        if (axis == X_AXIS) {
723            _cubeAxesActor->SetXAxisVisibility((state ? 1 : 0));
724        } else if (axis == Y_AXIS) {
725            _cubeAxesActor->SetYAxisVisibility((state ? 1 : 0));
726        } else if (axis == Z_AXIS) {
727            _cubeAxesActor->SetZAxisVisibility((state ? 1 : 0));
728        }
729        _needsRedraw = true;
730    }
731    if (_cubeAxesActor2D != NULL) {
732        if (axis == X_AXIS) {
733            _cubeAxesActor2D->SetXAxisVisibility((state ? 1 : 0));
734        } else if (axis == Y_AXIS) {
735            _cubeAxesActor2D->SetYAxisVisibility((state ? 1 : 0));
736        } else if (axis == Z_AXIS) {
737            _cubeAxesActor2D->SetZAxisVisibility((state ? 1 : 0));
738        }
739        _needsRedraw = true;
740    }
741}
742
743/**
744 * \brief Turn on/off rendering of single axis gridlines
745 */
746void Renderer::setAxisGridVisibility(Axis axis, bool state)
747{
748    if (_cubeAxesActor != NULL) {
749        if (axis == X_AXIS) {
750            _cubeAxesActor->SetDrawXGridlines((state ? 1 : 0));
751        } else if (axis == Y_AXIS) {
752            _cubeAxesActor->SetDrawYGridlines((state ? 1 : 0));
753        } else if (axis == Z_AXIS) {
754            _cubeAxesActor->SetDrawZGridlines((state ? 1 : 0));
755        }
756        _needsRedraw = true;
757    }
758}
759
760/**
761 * \brief Toggle label visibility for the specified axis
762 */
763void Renderer::setAxisLabelVisibility(Axis axis, bool state)
764{
765    if (_cubeAxesActor != NULL) {
766        if (axis == X_AXIS) {
767            _cubeAxesActor->SetXAxisLabelVisibility((state ? 1 : 0));
768        } else if (axis == Y_AXIS) {
769            _cubeAxesActor->SetYAxisLabelVisibility((state ? 1 : 0));
770        } else if (axis == Z_AXIS) {
771            _cubeAxesActor->SetZAxisLabelVisibility((state ? 1 : 0));
772        }
773        _needsRedraw = true;
774    }
775    if (_cubeAxesActor2D != NULL) {
776        if (axis == X_AXIS) {
777            _cubeAxesActor2D->GetXAxisActor2D()->SetLabelVisibility((state ? 1 : 0));
778        } else if (axis == Y_AXIS) {
779            _cubeAxesActor2D->GetYAxisActor2D()->SetLabelVisibility((state ? 1 : 0));
780        } else if (axis == Z_AXIS) {
781            _cubeAxesActor2D->GetZAxisActor2D()->SetLabelVisibility((state ? 1 : 0));
782        }
783        _needsRedraw = true;
784    }
785}
786
787/**
788 * \brief Toggle tick visibility for the specified axis
789 */
790void Renderer::setAxisTickVisibility(Axis axis, bool state)
791{
792    if (_cubeAxesActor != NULL) {
793        if (axis == X_AXIS) {
794            _cubeAxesActor->SetXAxisTickVisibility((state ? 1 : 0));
795        } else if (axis == Y_AXIS) {
796            _cubeAxesActor->SetYAxisTickVisibility((state ? 1 : 0));
797        } else if (axis == Z_AXIS) {
798            _cubeAxesActor->SetZAxisTickVisibility((state ? 1 : 0));
799        }
800        _needsRedraw = true;
801    }
802    if (_cubeAxesActor2D != NULL) {
803        if (axis == X_AXIS) {
804            _cubeAxesActor2D->GetXAxisActor2D()->SetTickVisibility((state ? 1 : 0));
805        } else if (axis == Y_AXIS) {
806            _cubeAxesActor2D->GetYAxisActor2D()->SetTickVisibility((state ? 1 : 0));
807        } else if (axis == Z_AXIS) {
808            _cubeAxesActor2D->GetZAxisActor2D()->SetTickVisibility((state ? 1 : 0));
809        }
810        _needsRedraw = true;
811    }
812}
813
814/**
815 * \brief Set title of the specified axis
816 */
817void Renderer::setAxisTitle(Axis axis, const char *title)
818{
819    if (_cubeAxesActor != NULL) {
820        if (axis == X_AXIS) {
821            _cubeAxesActor->SetXTitle(title);
822        } else if (axis == Y_AXIS) {
823            _cubeAxesActor->SetYTitle(title);
824        } else if (axis == Z_AXIS) {
825            _cubeAxesActor->SetZTitle(title);
826        }
827        _needsRedraw = true;
828    }
829    if (_cubeAxesActor2D != NULL) {
830        if (axis == X_AXIS) {
831            _cubeAxesActor2D->SetXLabel(title);
832        } else if (axis == Y_AXIS) {
833            _cubeAxesActor2D->SetYLabel(title);
834        } else if (axis == Z_AXIS) {
835            _cubeAxesActor2D->SetZLabel(title);
836        }
837        _needsRedraw = true;
838    }
839}
840
841/**
842 * \brief Set units of the specified axis
843 */
844void Renderer::setAxisUnits(Axis axis, const char *units)
845{
846    if (_cubeAxesActor != NULL) {
847        if (axis == X_AXIS) {
848            _cubeAxesActor->SetXUnits(units);
849        } else if (axis == Y_AXIS) {
850            _cubeAxesActor->SetYUnits(units);
851        } else if (axis == Z_AXIS) {
852            _cubeAxesActor->SetZUnits(units);
853        }
854        _needsRedraw = true;
855    }
856#ifdef notdef
857    if (_cubeAxesActor2D != NULL) {
858        if (axis == X_AXIS) {
859            _cubeAxesActor2D->SetXUnits(units);
860        } else if (axis == Y_AXIS) {
861            _cubeAxesActor2D->SetYUnits(units);
862        } else if (axis == Z_AXIS) {
863            _cubeAxesActor2D->SetZUnits(units);
864        }
865        _needsRedraw = true;
866    }
867#endif
868}
869
870/**
871 * \brief Notify graphics objects that color map has changed
872 */
873void Renderer::updateColorMap(ColorMap *cmap)
874{
875    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
876         itr != _contour3Ds.end(); ++itr) {
877        if (itr->second->getColorMap() == cmap) {
878            itr->second->updateColorMap();
879            _needsRedraw = true;
880        }
881    }
882    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
883         itr != _cutplanes.end(); ++itr) {
884        if (itr->second->getColorMap() == cmap) {
885            itr->second->updateColorMap();
886            _needsRedraw = true;
887        }
888    }
889    for (GlyphsHashmap::iterator itr = _glyphs.begin();
890         itr != _glyphs.end(); ++itr) {
891        if (itr->second->getColorMap() == cmap) {
892            itr->second->updateColorMap();
893            _needsRedraw = true;
894        }
895    }
896    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
897         itr != _heightMaps.end(); ++itr) {
898        if (itr->second->getColorMap() == cmap) {
899            itr->second->updateColorMap();
900            _needsRedraw = true;
901        }
902    }
903    for (LICHashmap::iterator itr = _lics.begin();
904         itr != _lics.end(); ++itr) {
905        if (itr->second->getColorMap() == cmap) {
906            itr->second->updateColorMap();
907            _needsRedraw = true;
908        }
909    }
910    for (MoleculeHashmap::iterator itr = _molecules.begin();
911         itr != _molecules.end(); ++itr) {
912        if (itr->second->getColorMap() == cmap) {
913            itr->second->updateColorMap();
914            _needsRedraw = true;
915        }
916    }
917    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
918         itr != _pseudoColors.end(); ++itr) {
919        if (itr->second->getColorMap() == cmap) {
920            itr->second->updateColorMap();
921            _needsRedraw = true;
922        }
923    }
924    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
925         itr != _streamlines.end(); ++itr) {
926        if (itr->second->getColorMap() == cmap) {
927            itr->second->updateColorMap();
928            _needsRedraw = true;
929        }
930    }
931    for (VolumeHashmap::iterator itr = _volumes.begin();
932         itr != _volumes.end(); ++itr) {
933        if (itr->second->getColorMap() == cmap) {
934            itr->second->updateColorMap();
935            _needsRedraw = true;
936        }
937    }
938}
939
940/**
941 * \brief Check if a ColorMap is in use by graphics objects
942 */
943bool Renderer::colorMapUsed(ColorMap *cmap)
944{
945    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
946         itr != _contour3Ds.end(); ++itr) {
947        if (itr->second->getColorMap() == cmap)
948            return true;
949    }
950    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
951         itr != _cutplanes.end(); ++itr) {
952        if (itr->second->getColorMap() == cmap)
953            return true;
954    }
955    for (GlyphsHashmap::iterator itr = _glyphs.begin();
956         itr != _glyphs.end(); ++itr) {
957        if (itr->second->getColorMap() == cmap)
958            return true;
959    }
960    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
961         itr != _heightMaps.end(); ++itr) {
962        if (itr->second->getColorMap() == cmap)
963            return true;
964    }
965    for (LICHashmap::iterator itr = _lics.begin();
966         itr != _lics.end(); ++itr) {
967        if (itr->second->getColorMap() == cmap)
968            return true;
969    }
970    for (MoleculeHashmap::iterator itr = _molecules.begin();
971         itr != _molecules.end(); ++itr) {
972        if (itr->second->getColorMap() == cmap)
973            return true;
974    }
975    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
976         itr != _pseudoColors.end(); ++itr) {
977        if (itr->second->getColorMap() == cmap)
978            return true;
979    }
980    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
981         itr != _streamlines.end(); ++itr) {
982        if (itr->second->getColorMap() == cmap)
983            return true;
984    }
985    for (VolumeHashmap::iterator itr = _volumes.begin();
986         itr != _volumes.end(); ++itr) {
987        if (itr->second->getColorMap() == cmap)
988            return true;
989    }
990    return false;
991}
992
993/**
994 * \brief Add/replace a ColorMap for use in the Renderer
995 */
996void Renderer::addColorMap(const ColorMapId& id, ColorMap *colorMap)
997{
998    if (colorMap != NULL) {
999        colorMap->build();
1000        if (getColorMap(id) != NULL) {
1001            TRACE("Replacing existing ColorMap %s", id.c_str());
1002            // Copy to current colormap to avoid invalidating
1003            // pointers in graphics objects using the color map
1004            *_colorMaps[id] = *colorMap;
1005            delete colorMap;
1006            // Notify graphics objects of change
1007            updateColorMap(_colorMaps[id]);
1008        } else
1009            _colorMaps[id] = colorMap;
1010    } else {
1011        ERROR("NULL ColorMap");
1012    }
1013}
1014
1015/**
1016 * \brief Return the ColorMap associated with the colormap key given
1017 */
1018ColorMap *Renderer::getColorMap(const ColorMapId& id)
1019{
1020    ColorMapHashmap::iterator itr = _colorMaps.find(id);
1021
1022    if (itr == _colorMaps.end())
1023        return NULL;
1024    else
1025        return itr->second;
1026}
1027
1028/**
1029 * \brief Remove the colormap associated with the key given
1030 *
1031 * The underlying vtkLookupTable will be deleted if it is not referenced
1032 * by any other objects
1033 */
1034void Renderer::deleteColorMap(const ColorMapId& id)
1035{
1036    ColorMapHashmap::iterator itr;
1037
1038    bool doAll = false;
1039
1040    if (id.compare("all") == 0) {
1041        itr = _colorMaps.begin();
1042        doAll = true;
1043    } else {
1044        itr = _colorMaps.find(id);
1045    }
1046
1047    if (itr == _colorMaps.end()) {
1048        ERROR("Unknown ColorMap %s", id.c_str());
1049        return;
1050    }
1051
1052    do {
1053        if (itr->second->getName().compare("default") == 0 ||
1054            itr->second->getName().compare("grayDefault") == 0 ||
1055            itr->second->getName().compare("volumeDefault") == 0 ||
1056            itr->second->getName().compare("elementDefault") == 0) {
1057            if (id.compare("all") != 0) {
1058                WARN("Cannot delete a default color map");
1059            }
1060            continue;
1061        } else if (colorMapUsed(itr->second)) {
1062            WARN("Cannot delete color map '%s', it is in use", itr->second->getName().c_str());
1063            continue;
1064        }
1065
1066        TRACE("Deleting ColorMap %s", itr->second->getName().c_str());
1067
1068        delete itr->second;
1069        itr = _colorMaps.erase(itr);
1070    } while (doAll && itr != _colorMaps.end());
1071}
1072
1073/**
1074 * \brief Render a labelled legend image for the given colormap
1075 *
1076 * The field is assumed to be the active scalar or vector field
1077 * based on the legendType.
1078 *
1079 * \param[in] id ColorMap name
1080 * \param[in] dataSetID DataSet name
1081 * \param[in] legendType scalar or vector field legend
1082 * \param[in,out] title If supplied, draw title ("#auto" means to
1083 * fill in field name and draw).  If blank, do not draw title. 
1084 * If title was blank or "#auto", will be filled with field name on
1085 * return
1086 * \param[in,out] range Data range to use in legend.  Set min > max to have
1087 * range computed, will be filled with valid min and max values
1088 * \param[in] width Pixel width of legend (aspect controls orientation)
1089 * \param[in] height Pixel height of legend (aspect controls orientation)
1090 * \param[in] numLabels Number of labels to render (includes min/max)
1091 * \param[in,out] imgData Pointer to array to fill with image bytes. Array
1092 * will be resized if needed.
1093 * \return The image is rendered into the supplied array, false is
1094 * returned if the color map is not found
1095 */
1096bool Renderer::renderColorMap(const ColorMapId& id,
1097                              const DataSetId& dataSetID,
1098                              Renderer::LegendType legendType,
1099                              std::string& title,
1100                              double range[2],
1101                              int width, int height,
1102                              int numLabels,
1103                              vtkUnsignedCharArray *imgData)
1104{
1105    DataSet *dataSet = NULL;
1106    if (dataSetID.compare("all") == 0) {
1107        if (_dataSets.empty()) {
1108            WARN("No DataSets exist, can't fill title or range");
1109            return renderColorMap(id, dataSetID, legendType,
1110                                  NULL,
1111                                  DataSet::POINT_DATA,
1112                                  title, range, width, height, numLabels, imgData);
1113        } else {
1114            dataSet = _dataSets.begin()->second;
1115        }
1116    } else {
1117        dataSet = getDataSet(dataSetID);
1118        if (dataSet == NULL) {
1119            ERROR("DataSet '%s' not found", dataSetID.c_str());
1120            return false;
1121        }
1122    }
1123
1124    if (legendType == LEGEND_SCALAR) {
1125        return renderColorMap(id, dataSetID, legendType,
1126                              dataSet->getActiveScalarsName(),
1127                              dataSet->getActiveScalarsType(),
1128                              title, range, width, height, numLabels, imgData);
1129    } else {
1130        return renderColorMap(id, dataSetID, legendType,
1131                              dataSet->getActiveVectorsName(),
1132                              dataSet->getActiveVectorsType(),
1133                              title, range, width, height, numLabels, imgData);
1134    }
1135}
1136
1137/**
1138 * \brief Render a labelled legend image for the given colormap
1139 *
1140 * The field is assumed to be point data, if the field is not found
1141 * as point data, cell data is used.
1142 *
1143 * \param[in] id ColorMap name
1144 * \param[in] dataSetID DataSet name
1145 * \param[in] legendType scalar or vector field legend
1146 * \param[in] fieldName Name of the field array this legend is for
1147 * \param[in,out] title If supplied, draw title ("#auto" means to
1148 * fill in field name and draw).  If blank, do not draw title. 
1149 * If title was blank or "#auto", will be filled with field name on
1150 * return
1151 * \param[in,out] range Data range to use in legend.  Set min > max to have
1152 * range computed, will be filled with valid min and max values
1153 * \param[in] width Pixel width of legend (aspect controls orientation)
1154 * \param[in] height Pixel height of legend (aspect controls orientation)
1155 * \param[in] numLabels Number of labels to render (includes min/max)
1156 * \param[in,out] imgData Pointer to array to fill with image bytes. Array
1157 * will be resized if needed.
1158 * \return The image is rendered into the supplied array, false is
1159 * returned if the color map is not found
1160 */
1161bool Renderer::renderColorMap(const ColorMapId& id,
1162                              const DataSetId& dataSetID,
1163                              Renderer::LegendType legendType,
1164                              const char *fieldName,
1165                              std::string& title,
1166                              double range[2],
1167                              int width, int height,
1168                              int numLabels,
1169                              vtkUnsignedCharArray *imgData)
1170{
1171    DataSet *dataSet = NULL;
1172    if (dataSetID.compare("all") == 0) {
1173        if (_dataSets.empty()) {
1174            WARN("No DataSets exist, can't fill title or range");
1175            return renderColorMap(id, dataSetID, legendType,
1176                                  NULL,
1177                                  DataSet::POINT_DATA,
1178                                  title, range, width, height, numLabels, imgData);
1179        } else {
1180            dataSet = _dataSets.begin()->second;
1181        }
1182    } else {
1183        dataSet = getDataSet(dataSetID);
1184        if (dataSet == NULL) {
1185            ERROR("DataSet '%s' not found", dataSetID.c_str());
1186            return false;
1187        }
1188    }
1189
1190    DataSet::DataAttributeType attrType;
1191    int numComponents;
1192
1193    dataSet->getFieldInfo(fieldName, &attrType, &numComponents);
1194
1195    return renderColorMap(id, dataSetID, legendType,
1196                          fieldName,
1197                          attrType,
1198                          title, range, width, height, numLabels, imgData);
1199}
1200
1201/**
1202 * \brief Render a labelled legend image for the given colormap
1203 *
1204 * \param[in] id ColorMap name
1205 * \param[in] dataSetID DataSet name
1206 * \param[in] legendType scalar or vector field legend
1207 * \param[in] fieldName Name of the field array this legend is for
1208 * \param[in] type DataAttributeType of the field
1209 * \param[in,out] title If supplied, draw title ("#auto" means to
1210 * fill in field name and draw).  If blank, do not draw title. 
1211 * If title was blank or "#auto", will be filled with field name on
1212 * return
1213 * \param[in,out] range Data range to use in legend.  Set min > max to have
1214 * range computed, will be filled with valid min and max values
1215 * \param[in] width Pixel width of legend (aspect controls orientation)
1216 * \param[in] height Pixel height of legend (aspect controls orientation)
1217 * \param[in] numLabels Number of labels to render (includes min/max)
1218 * \param[in,out] imgData Pointer to array to fill with image bytes. Array
1219 * will be resized if needed.
1220 * \return The image is rendered into the supplied array, false is
1221 * returned if the color map is not found
1222 */
1223bool Renderer::renderColorMap(const ColorMapId& id,
1224                              const DataSetId& dataSetID,
1225                              Renderer::LegendType legendType,
1226                              const char *fieldName,
1227                              DataSet::DataAttributeType type,
1228                              std::string& title,
1229                              double range[2],
1230                              int width, int height,
1231                              int numLabels,
1232                              vtkUnsignedCharArray *imgData)
1233{
1234    TRACE("Enter");
1235    ColorMap *colorMap = getColorMap(id);
1236    if (colorMap == NULL)
1237        return false;
1238
1239    if (_legendRenderWindow == NULL) {
1240        _legendRenderWindow = vtkSmartPointer<vtkRenderWindow>::New();
1241        _legendRenderWindow->SetMultiSamples(0);
1242#ifdef USE_OFFSCREEN_RENDERING
1243        _legendRenderWindow->DoubleBufferOff();
1244        _legendRenderWindow->OffScreenRenderingOn();
1245#else
1246        _legendRenderWindow->DoubleBufferOn();
1247        _legendRenderWindow->SwapBuffersOff();
1248#endif
1249    }
1250
1251    _legendRenderWindow->SetSize(width, height);
1252
1253    if (_legendRenderer == NULL) {
1254        _legendRenderer = vtkSmartPointer<vtkRenderer>::New();
1255        _legendRenderWindow->AddRenderer(_legendRenderer);
1256    }
1257    _legendRenderer->SetBackground(_bgColor[0], _bgColor[1], _bgColor[2]);
1258
1259    if (_scalarBarActor == NULL) {
1260        _scalarBarActor = vtkSmartPointer<vtkScalarBarActor>::New();
1261        _scalarBarActor->UseOpacityOn();
1262        _legendRenderer->AddViewProp(_scalarBarActor);
1263    }
1264
1265    if (width > height) {
1266        _scalarBarActor->SetOrientationToHorizontal();
1267    } else {
1268        _scalarBarActor->SetOrientationToVertical();
1269    }
1270
1271    // Set viewport-relative width/height/pos
1272    if (title.empty() && numLabels == 0) {
1273        _scalarBarActor->SetPosition(0, 0);
1274        if (width > height) {
1275            // horizontal
1276            _scalarBarActor->SetHeight(2.5); // VTK: floor(actorHeight * .4)
1277            _scalarBarActor->SetWidth(1);
1278        } else {
1279            // vertical
1280            _scalarBarActor->SetHeight((double)height/(.86*height)); // VTK: floor(actorHeight * .86)
1281            _scalarBarActor->SetWidth(((double)(width+5))/((double)width)); // VTK: actorWidth - 4 pixels
1282        }
1283    } else {
1284        if (width > height) {
1285            // horizontal
1286            _scalarBarActor->SetPosition(.075, .1);
1287            _scalarBarActor->SetHeight(0.8);
1288            _scalarBarActor->SetWidth(0.85);
1289        } else {
1290            // vertical
1291            _scalarBarActor->SetPosition(.1, .05);
1292            _scalarBarActor->SetHeight(0.9);
1293            _scalarBarActor->SetWidth(0.8);
1294        }
1295    }
1296
1297    vtkSmartPointer<vtkLookupTable> lut = colorMap->getLookupTable();
1298    DataSet *dataSet = NULL;
1299    bool cumulative = _useCumulativeRange;
1300    if (dataSetID.compare("all") == 0) {
1301        if (_dataSets.empty()) {
1302            WARN("No DataSets exist, can't fill title or range");
1303        } else {
1304            dataSet = _dataSets.begin()->second;
1305        }
1306        cumulative = true;
1307    } else {
1308        dataSet = getDataSet(dataSetID);
1309        if (dataSet == NULL) {
1310            ERROR("DataSet '%s' not found", dataSetID.c_str());
1311            return false;
1312        }
1313    }
1314
1315    bool drawTitle = false;
1316    if (!title.empty()) {
1317        drawTitle = true;
1318        if (title.compare("#auto") == 0) {
1319            title.clear();
1320        }
1321    }
1322
1323    bool needRange = false;
1324    if (range[0] > range[1]) {
1325        range[0] = 0.0;
1326        range[1] = 1.0;
1327        needRange = true;
1328    }
1329
1330    switch (legendType) {
1331    case LEGEND_VECTOR_MAGNITUDE:
1332        if (needRange) {
1333            if (cumulative) {
1334                getCumulativeDataRange(range, fieldName, type, 3);
1335            } else if (dataSet != NULL) {
1336                dataSet->getDataRange(range, fieldName, type);
1337            }
1338        }
1339
1340        lut->SetRange(range);
1341
1342        if (title.empty() && dataSet != NULL) {
1343            if (fieldName != NULL) {
1344                title = fieldName;
1345                title.append("(mag)");
1346            }
1347        }
1348        break;
1349    case LEGEND_VECTOR_X:
1350        if (needRange) {
1351            if (cumulative) {
1352                getCumulativeDataRange(range, fieldName, type, 3, 0);
1353            } else if (dataSet != NULL) {
1354                dataSet->getDataRange(range, fieldName, type, 0);
1355            }
1356        }
1357
1358        lut->SetRange(range);
1359
1360        if (title.empty() && dataSet != NULL) {
1361            if (fieldName != NULL) {
1362                title = fieldName;
1363                title.append("(x)");
1364            }
1365        }
1366        break;
1367    case LEGEND_VECTOR_Y:
1368        if (needRange) {
1369            if (cumulative) {
1370                getCumulativeDataRange(range, fieldName, type, 3, 1);
1371            } else if (dataSet != NULL) {
1372                dataSet->getDataRange(range, fieldName, type, 1);
1373            }
1374        }
1375
1376        lut->SetRange(range);
1377
1378        if (title.empty() && dataSet != NULL) {
1379            if (fieldName != NULL) {
1380                title = fieldName;
1381                title.append("(y)");
1382            }
1383        }
1384        break;
1385    case LEGEND_VECTOR_Z:
1386        if (needRange) {
1387            if (cumulative) {
1388                getCumulativeDataRange(range, fieldName, type, 3, 2);
1389            } else if (dataSet != NULL) {
1390                dataSet->getDataRange(range, fieldName, type, 1);
1391            }
1392        }
1393
1394        lut->SetRange(range);
1395
1396        if (title.empty() && dataSet != NULL) {
1397            if (fieldName != NULL) {
1398                title = fieldName;
1399                title.append("(z)");
1400            }
1401        }
1402        break;
1403    case LEGEND_SCALAR:
1404    default:
1405        if (needRange) {
1406            if (cumulative) {
1407                getCumulativeDataRange(range, fieldName, type, 1);
1408            } else if (dataSet != NULL) {
1409                dataSet->getDataRange(range, fieldName, type);
1410            }
1411        }
1412
1413        lut->SetRange(range);
1414
1415        if (title.empty() && dataSet != NULL) {
1416            if (fieldName != NULL)
1417                title = fieldName;
1418        }
1419        break;
1420    }
1421
1422    _scalarBarActor->SetLookupTable(lut);
1423
1424    if (drawTitle) {
1425        _scalarBarActor->SetTitle(title.c_str());
1426    } else {
1427        _scalarBarActor->SetTitle("");
1428    }
1429    _scalarBarActor->GetTitleTextProperty()->ItalicOff();
1430    _scalarBarActor->SetNumberOfLabels(numLabels);
1431    _scalarBarActor->GetLabelTextProperty()->BoldOff();
1432    _scalarBarActor->GetLabelTextProperty()->ItalicOff();
1433    _scalarBarActor->GetLabelTextProperty()->ShadowOff();
1434
1435    _legendRenderWindow->Render();
1436
1437#ifdef RENDER_TARGA
1438    _legendRenderWindow->MakeCurrent();
1439    // Must clear previous errors first.
1440    while (glGetError() != GL_NO_ERROR){
1441        ;
1442    }
1443    int bytesPerPixel = TARGA_BYTES_PER_PIXEL;
1444    int size = bytesPerPixel * width * height;
1445
1446    if (imgData->GetMaxId() + 1 != size)
1447    {
1448        imgData->SetNumberOfComponents(bytesPerPixel);
1449        imgData->SetNumberOfValues(size);
1450    }
1451    glDisable(GL_TEXTURE_2D);
1452    if (_legendRenderWindow->GetDoubleBuffer()) {
1453        glReadBuffer(static_cast<GLenum>(vtkOpenGLRenderWindow::SafeDownCast(_legendRenderWindow)->GetBackLeftBuffer()));
1454    } else {
1455        glReadBuffer(static_cast<GLenum>(vtkOpenGLRenderWindow::SafeDownCast(_legendRenderWindow)->GetFrontLeftBuffer()));
1456    }
1457    glPixelStorei(GL_PACK_ALIGNMENT, 1);
1458    if (bytesPerPixel == 4) {
1459        glReadPixels(0, 0, width, height, GL_BGRA,
1460                     GL_UNSIGNED_BYTE, imgData->GetPointer(0));
1461    } else {
1462        glReadPixels(0, 0, width, height, GL_BGR,
1463                     GL_UNSIGNED_BYTE, imgData->GetPointer(0));
1464    }
1465    if (glGetError() != GL_NO_ERROR) {
1466        ERROR("glReadPixels");
1467    }
1468#else
1469    _legendRenderWindow->GetPixelData(0, 0, width-1, height-1,
1470                                      !_legendRenderWindow->GetDoubleBuffer(),
1471                                      imgData);
1472#endif
1473    TRACE("Leave");
1474    return true;
1475}
1476
1477/**
1478 * \brief Set camera FOV based on render window height
1479 *
1480 * Computes a field-of-view angle based on some assumptions about
1481 * viewer's distance to screen and pixel density
1482 */
1483void Renderer::setViewAngle(int height)
1484{
1485    // Distance of eyes from screen in inches
1486    double d = 20.0;
1487    // Assume 72 ppi screen
1488    double h = (double)height / 72.0;
1489
1490    double angle = vtkMath::DegreesFromRadians(2.0 * atan((h/2.0)/d));
1491    _renderer->GetActiveCamera()->SetViewAngle(angle);
1492
1493    TRACE("Setting view angle: %g", angle);
1494}
1495
1496/**
1497 * \brief Resize the render window (image size for renderings)
1498 */
1499void Renderer::setWindowSize(int width, int height)
1500{
1501    if (_windowWidth == width &&
1502        _windowHeight == height)
1503        return;
1504
1505    //setViewAngle(height);
1506
1507    // FIXME: Fix up panning on aspect change
1508#ifdef notdef
1509    if (_cameraPan[0] != 0.0) {
1510        _cameraPan[0] *= ((double)_windowWidth / width);
1511    }
1512    if (_cameraPan[1] != 0.0) {
1513        _cameraPan[1] *= ((double)_windowHeight / height);
1514    }
1515#endif
1516    _windowWidth = width;
1517    _windowHeight = height;
1518    _renderWindow->SetSize(_windowWidth, _windowHeight);
1519    if (_cameraMode == IMAGE) {
1520        setCameraZoomRegion(_imgWorldOrigin[0], _imgWorldOrigin[1],
1521                            _imgWorldDims[0], _imgWorldDims[1]);
1522    }
1523    _needsRedraw = true;
1524}
1525
1526/**
1527 * \brief Change the camera type: perspective, orthographic or image view
1528 *
1529 * Perspective mode is a normal 3D camera.
1530 *
1531 * Orthogrphic mode is parallel projection.
1532 *
1533 * Image mode is an orthographic camera with fixed axes and a clipping region
1534 * around the plot area, use setCameraZoomRegion to control the displayed area
1535 *
1536 * \param[in] mode Enum specifying camera type
1537 */
1538void Renderer::setCameraMode(CameraMode mode)
1539{
1540    if (_cameraMode == mode) return;
1541
1542    CameraMode origMode = _cameraMode;
1543    _cameraMode = mode;
1544    resetAxes();
1545
1546    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1547    switch (mode) {
1548    case ORTHO: {
1549        TRACE("Set camera to Ortho mode");
1550        camera->ParallelProjectionOn();
1551        if (origMode == IMAGE) {
1552            resetCamera(true);
1553        }
1554        break;
1555    }
1556    case PERSPECTIVE: {
1557        TRACE("Set camera to Perspective mode");
1558        camera->ParallelProjectionOff();
1559        if (origMode == IMAGE) {
1560            resetCamera(true);
1561        }
1562        break;
1563    }
1564    case IMAGE: {
1565        camera->ParallelProjectionOn();
1566        setCameraZoomRegion(_imgWorldOrigin[0], _imgWorldOrigin[1],
1567                            _imgWorldDims[0],_imgWorldDims[1]);
1568        TRACE("Set camera to Image mode");
1569        break;
1570    }
1571    default:
1572        ERROR("Unkown camera mode: %d", mode);
1573    }
1574    _needsRedraw = true;
1575}
1576
1577/**
1578 * \brief Get the current camera mode
1579 */
1580Renderer::CameraMode Renderer::getCameraMode() const
1581{
1582    return _cameraMode;
1583}
1584
1585/**
1586 * \brief Set the VTK camera parameters based on a 4x4 view matrix
1587 */
1588void Renderer::setCameraFromMatrix(vtkCamera *camera, vtkMatrix4x4& mat)
1589{
1590    double d = camera->GetDistance();
1591    double vu[3];
1592    vu[0] = mat[1][0];
1593    vu[1] = mat[1][1];
1594    vu[2] = mat[1][2];
1595    double trans[3];
1596    trans[0] = mat[0][3];
1597    trans[1] = mat[1][3];
1598    trans[2] = mat[2][3];
1599    mat[0][3] = 0;
1600    mat[1][3] = 0;
1601    mat[2][3] = 0;
1602    double vpn[3] = {mat[2][0], mat[2][1], mat[2][2]};
1603    double pos[3];
1604    // With translation removed, we have an orthogonal matrix,
1605    // so the inverse is the transpose
1606    mat.Transpose();
1607    mat.MultiplyPoint(trans, pos);
1608    pos[0] = -pos[0];
1609    pos[1] = -pos[1];
1610    pos[2] = -pos[2];
1611    double fp[3];
1612    fp[0] = pos[0] - vpn[0] * d;
1613    fp[1] = pos[1] - vpn[1] * d;
1614    fp[2] = pos[2] - vpn[2] * d;
1615    camera->SetPosition(pos);
1616    camera->SetFocalPoint(fp);
1617    camera->SetViewUp(vu);
1618}
1619
1620/**
1621 * \brief Set the orientation of the camera from a quaternion
1622 * rotation
1623 *
1624 * \param[in] quat A quaternion with scalar part first: w,x,y,z
1625 * \param[in] absolute Is rotation absolute or relative?
1626 */
1627void Renderer::setCameraOrientation(const double quat[4], bool absolute)
1628{
1629    if (_cameraMode == IMAGE)
1630        return;
1631    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1632    vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
1633    vtkSmartPointer<vtkMatrix4x4> rotMat = vtkSmartPointer<vtkMatrix4x4>::New();
1634
1635    double q[4];
1636    copyQuat(quat, q);
1637
1638    if (absolute) {
1639        double abs[4];
1640        // Save absolute rotation
1641        copyQuat(q, abs);
1642        // Compute relative rotation
1643        quatMult(quatReciprocal(_cameraOrientation), q, q);
1644        // Store absolute rotation
1645        copyQuat(abs, _cameraOrientation);
1646    } else {
1647        // Compute new absolute rotation
1648        quatMult(_cameraOrientation, q, _cameraOrientation);
1649    }
1650
1651    quaternionToTransposeMatrix4x4(q, *rotMat);
1652#ifdef DEBUG
1653    TRACE("Rotation matrix:\n %g %g %g\n %g %g %g\n %g %g %g",
1654          (*rotMat)[0][0], (*rotMat)[0][1], (*rotMat)[0][2],
1655          (*rotMat)[1][0], (*rotMat)[1][1], (*rotMat)[1][2],
1656          (*rotMat)[2][0], (*rotMat)[2][1], (*rotMat)[2][2]);
1657    vtkSmartPointer<vtkMatrix4x4> camMat = vtkSmartPointer<vtkMatrix4x4>::New();
1658    camMat->DeepCopy(camera->GetViewTransformMatrix());
1659    TRACE("Camera matrix:\n %g %g %g %g\n %g %g %g %g\n %g %g %g %g\n %g %g %g %g",
1660          (*camMat)[0][0], (*camMat)[0][1], (*camMat)[0][2], (*camMat)[0][3],
1661          (*camMat)[1][0], (*camMat)[1][1], (*camMat)[1][2], (*camMat)[1][3],
1662          (*camMat)[2][0], (*camMat)[2][1], (*camMat)[2][2], (*camMat)[2][3],
1663          (*camMat)[3][0], (*camMat)[3][1], (*camMat)[3][2], (*camMat)[3][3]);
1664    printCameraInfo(camera);
1665#endif
1666    trans->Translate(0, 0, -camera->GetDistance());
1667    trans->Concatenate(rotMat);
1668    trans->Translate(0, 0, camera->GetDistance());
1669    trans->Concatenate(camera->GetViewTransformMatrix());
1670    setCameraFromMatrix(camera, *trans->GetMatrix());
1671
1672    _renderer->ResetCameraClippingRange();
1673    printCameraInfo(camera);
1674    _needsRedraw = true;
1675}
1676
1677/**
1678 * \brief Set the position and orientation of the camera
1679 *
1680 * \param[in] position x,y,z position of camera in world coordinates
1681 * \param[in] focalPoint x,y,z look-at point in world coordinates
1682 * \param[in] viewUp x,y,z up vector of camera
1683 */
1684void Renderer::setCameraOrientationAndPosition(const double position[3],
1685                                               const double focalPoint[3],
1686                                               const double viewUp[3])
1687{
1688    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1689    camera->SetPosition(position);
1690    camera->SetFocalPoint(focalPoint);
1691    camera->SetViewUp(viewUp);
1692    _renderer->ResetCameraClippingRange();
1693    _needsRedraw = true;
1694}
1695
1696/**
1697 * \brief Get the position and orientation of the camera
1698 *
1699 * \param[out] position x,y,z position of camera in world coordinates
1700 * \param[out] focalPoint x,y,z look-at point in world coordinates
1701 * \param[out] viewUp x,y,z up vector of camera
1702 */
1703void Renderer::getCameraOrientationAndPosition(double position[3],
1704                                               double focalPoint[3],
1705                                               double viewUp[3])
1706{
1707    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1708    camera->GetPosition(position);
1709    camera->GetFocalPoint(focalPoint);
1710    camera->GetViewUp(viewUp);
1711}
1712
1713/**
1714 * \brief Reset pan, zoom, clipping planes and optionally rotation
1715 *
1716 * \param[in] resetOrientation Reset the camera rotation/orientation also
1717 */
1718void Renderer::resetCamera(bool resetOrientation)
1719{
1720    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1721    if (_cameraMode == IMAGE) {
1722        initCamera();
1723    } else {
1724        if (resetOrientation) {
1725            camera->SetPosition(0, 0, 1);
1726            camera->SetFocalPoint(0, 0, 0);
1727            camera->SetViewUp(0, 1, 0);
1728            _cameraOrientation[0] = 1.0;
1729            _cameraOrientation[1] = 0.0;
1730            _cameraOrientation[2] = 0.0;
1731            _cameraOrientation[3] = 0.0;
1732        }
1733        //setViewAngle(_windowHeight);
1734        double bounds[6];
1735        collectBounds(bounds, false);
1736        _renderer->ResetCamera(bounds);
1737        _renderer->ResetCameraClippingRange();
1738        //computeScreenWorldCoords();
1739    }
1740
1741    printCameraInfo(camera);
1742
1743    _cameraZoomRatio = 1;
1744    _cameraPan[0] = 0;
1745    _cameraPan[1] = 0;
1746
1747    _needsRedraw = true;
1748}
1749
1750void Renderer::resetCameraClippingRange()
1751{
1752    _renderer->ResetCameraClippingRange();
1753}
1754
1755/**
1756 * \brief Perform a relative rotation to current camera orientation
1757 *
1758 * Angles are in degrees, rotation is about focal point
1759 */
1760void Renderer::rotateCamera(double yaw, double pitch, double roll)
1761{
1762    if (_cameraMode == IMAGE)
1763        return;
1764    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1765    camera->Azimuth(yaw); // Rotate about object
1766    //camera->SetYaw(yaw); // Rotate about camera
1767    camera->Elevation(pitch); // Rotate about object
1768    //camera->SetPitch(pitch); // Rotate about camera
1769    camera->Roll(roll); // Roll about camera view axis
1770    _renderer->ResetCameraClippingRange();
1771    //computeScreenWorldCoords();
1772    _needsRedraw = true;
1773}
1774
1775/**
1776 * \brief Perform a 2D translation of the camera
1777 *
1778 * x,y pan amount are specified as signed absolute pan amount in viewport
1779 * units -- i.e. 0 is no pan, .5 is half the viewport, 2 is twice the viewport,
1780 * etc.
1781 *
1782 * \param[in] x Viewport coordinate horizontal panning (positive number pans
1783 * camera left, object right)
1784 * \param[in] y Viewport coordinate vertical panning (positive number pans
1785 * camera up, object down)
1786 * \param[in] absolute Control if pan amount is relative to current or absolute
1787 */
1788void Renderer::panCamera(double x, double y, bool absolute)
1789{
1790    TRACE("Enter panCamera: %g %g, current abs: %g %g",
1791          x, y, _cameraPan[0], _cameraPan[1]);
1792
1793    if (_cameraMode == IMAGE) {
1794        // Reverse x rather than y, since we are panning the camera, and client
1795        // expects to be panning/moving the object
1796        x = -x * _screenWorldCoords[2];
1797        y = y * _screenWorldCoords[3];
1798
1799        if (absolute) {
1800            double panAbs[2];
1801            panAbs[0] = x;
1802            panAbs[1] = y;
1803            x -= _cameraPan[0];
1804            y -= _cameraPan[1];
1805            _cameraPan[0] = panAbs[0];
1806            _cameraPan[1] = panAbs[1];
1807        } else {
1808            _cameraPan[0] += x;
1809            _cameraPan[1] += y;
1810        }
1811
1812        _imgWorldOrigin[0] += x;
1813        _imgWorldOrigin[1] += y;
1814        setCameraZoomRegion(_imgWorldOrigin[0], _imgWorldOrigin[1],
1815                            _imgWorldDims[0], _imgWorldDims[1]);
1816    } else {
1817        y = -y;
1818        if (absolute) {
1819            double panAbs[2];
1820            panAbs[0] = x;
1821            panAbs[1] = y;
1822            x -= _cameraPan[0];
1823            y -= _cameraPan[1];
1824            _cameraPan[0] = panAbs[0];
1825            _cameraPan[1] = panAbs[1];
1826        } else {
1827            _cameraPan[0] += x;
1828            _cameraPan[1] += y;
1829        }
1830
1831        if (x != 0.0 || y != 0.0) {
1832            vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1833            double viewFocus[4], focalDepth, viewPoint[3];
1834            double newPickPoint[4], oldPickPoint[4], motionVector[3];
1835
1836            camera->GetFocalPoint(viewFocus);
1837            computeWorldToDisplay(viewFocus[0], viewFocus[1], viewFocus[2],
1838                                  viewFocus);
1839            focalDepth = viewFocus[2];
1840
1841            computeDisplayToWorld(( x * 2. + 1.) * _windowWidth / 2.0,
1842                                  ( y * 2. + 1.) * _windowHeight / 2.0,
1843                                  focalDepth,
1844                                  newPickPoint);
1845
1846            computeDisplayToWorld(_windowWidth / 2.0,
1847                                  _windowHeight / 2.0,
1848                                  focalDepth,
1849                                  oldPickPoint);
1850
1851            // Camera motion is reversed
1852            motionVector[0] = oldPickPoint[0] - newPickPoint[0];
1853            motionVector[1] = oldPickPoint[1] - newPickPoint[1];
1854            motionVector[2] = oldPickPoint[2] - newPickPoint[2];
1855
1856            camera->GetFocalPoint(viewFocus);
1857            camera->GetPosition(viewPoint);
1858            camera->SetFocalPoint(motionVector[0] + viewFocus[0],
1859                                  motionVector[1] + viewFocus[1],
1860                                  motionVector[2] + viewFocus[2]);
1861
1862            camera->SetPosition(motionVector[0] + viewPoint[0],
1863                                motionVector[1] + viewPoint[1],
1864                                motionVector[2] + viewPoint[2]);
1865
1866            _renderer->ResetCameraClippingRange();
1867            //computeScreenWorldCoords();
1868        }
1869    }
1870
1871    TRACE("Leave panCamera: %g %g, current abs: %g %g",
1872          x, y, _cameraPan[0], _cameraPan[1]);
1873
1874    _needsRedraw = true;
1875}
1876
1877/**
1878 * \brief Dolly camera or set orthographic scaling based on camera type
1879 *
1880 * \param[in] z Ratio to change zoom (greater than 1 is zoom in, less than 1 is zoom out)
1881 * \param[in] absolute Control if zoom factor is relative to current setting or absolute
1882 */
1883void Renderer::zoomCamera(double z, bool absolute)
1884{
1885    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
1886    TRACE("Enter Zoom: current abs: %g, z: %g, view angle %g",
1887          _cameraZoomRatio, z, camera->GetViewAngle());
1888
1889    if (absolute) {
1890        assert(_cameraZoomRatio > 0.0);
1891        double zAbs = z;
1892        z *= 1.0/_cameraZoomRatio;
1893        _cameraZoomRatio = zAbs;
1894    } else {
1895        _cameraZoomRatio *= z;
1896    }
1897
1898    if (_cameraMode == IMAGE) {
1899        double dx = _imgWorldDims[0];
1900        double dy = _imgWorldDims[1];
1901        _imgWorldDims[0] /= z;
1902        _imgWorldDims[1] /= z;
1903        dx -= _imgWorldDims[0];
1904        dy -= _imgWorldDims[1];
1905        _imgWorldOrigin[0] += dx/2.0;
1906        _imgWorldOrigin[1] += dy/2.0;
1907        setCameraZoomRegion(_imgWorldOrigin[0], _imgWorldOrigin[1],
1908                            _imgWorldDims[0], _imgWorldDims[1]);
1909    } else {
1910        // Keep ortho and perspective modes in sync
1911        // Move camera forward/back for perspective camera
1912        camera->Dolly(z);
1913        // Change ortho parallel scale
1914        camera->SetParallelScale(camera->GetParallelScale()/z);
1915        _renderer->ResetCameraClippingRange();
1916        //computeScreenWorldCoords();
1917    }
1918
1919    TRACE("Leave Zoom: rel %g, new abs: %g, view angle %g",
1920          z, _cameraZoomRatio, camera->GetViewAngle());
1921
1922    _needsRedraw = true;
1923}
1924
1925/**
1926 * \brief Set the pan/zoom using a corner and dimensions in pixel coordinates
1927 *
1928 * \param[in] x left pixel coordinate
1929 * \param[in] y bottom  pixel coordinate (with y=0 at top of window)
1930 * \param[in] width Width of zoom region in pixel coordinates
1931 * \param[in] height Height of zoom region in pixel coordinates
1932 */
1933void Renderer::setCameraZoomRegionPixels(int x, int y, int width, int height)
1934{
1935    double wx, wy, ww, wh;
1936
1937    y = _windowHeight - y;
1938    double pxToWorldX = _screenWorldCoords[2] / (double)_windowWidth;
1939    double pxToWorldY = _screenWorldCoords[3] / (double)_windowHeight;
1940
1941    wx = _screenWorldCoords[0] + x * pxToWorldX;
1942    wy = _screenWorldCoords[1] + y * pxToWorldY;
1943    ww = abs(width) *  pxToWorldX;
1944    wh = abs(height) * pxToWorldY;
1945    setCameraZoomRegion(wx, wy, ww, wh);
1946
1947    TRACE("\npx: %d %d %d %d\nworld: %g %g %g %g",
1948          x, y, width, height,
1949          wx, wy, ww, wh);
1950}
1951
1952/**
1953 * \brief Set the pan/zoom using a corner and dimensions in world coordinates
1954 *
1955 * \param[in] x left world coordinate
1956 * \param[in] y bottom  world coordinate
1957 * \param[in] width Width of zoom region in world coordinates
1958 * \param[in] height Height of zoom region in world coordinates
1959 */
1960void Renderer::setCameraZoomRegion(double x, double y, double width, double height)
1961{
1962    double camPos[2];
1963
1964    int pxOffsetX = (int)(0.17 * (double)_windowWidth);
1965    pxOffsetX = (pxOffsetX > 100 ? 100 : pxOffsetX);
1966    int pxOffsetY = (int)(0.15 * (double)_windowHeight);
1967    pxOffsetY = (pxOffsetY > 75 ? 75 : pxOffsetY);
1968    int outerGutter = (int)(0.03 * (double)_windowWidth);
1969    outerGutter = (outerGutter > 15 ? 15 : outerGutter);
1970
1971    int imgHeightPx = _windowHeight - pxOffsetY - outerGutter;
1972    int imgWidthPx = _windowWidth - pxOffsetX - outerGutter;
1973
1974    double imgAspect = width / height;
1975    double winAspect = (double)_windowWidth / _windowHeight;
1976
1977    double pxToWorld;
1978
1979    if (imgAspect >= winAspect) {
1980        pxToWorld = width / imgWidthPx;
1981    } else {
1982        pxToWorld = height / imgHeightPx;
1983    }
1984
1985    double offsetX = pxOffsetX * pxToWorld;
1986    double offsetY = pxOffsetY * pxToWorld;
1987
1988    TRACE("Window: %d %d", _windowWidth, _windowHeight);
1989    TRACE("ZoomRegion: %g %g %g %g", x, y, width, height);
1990    TRACE("pxToWorld: %g", pxToWorld);
1991    TRACE("offset: %g %g", offsetX, offsetY);
1992
1993    setCameraMode(IMAGE);
1994
1995    _imgWorldOrigin[0] = x;
1996    _imgWorldOrigin[1] = y;
1997    _imgWorldDims[0] = width;
1998    _imgWorldDims[1] = height;
1999
2000    camPos[0] = _imgWorldOrigin[0] - offsetX + (_windowWidth * pxToWorld)/2.0;
2001    camPos[1] = _imgWorldOrigin[1] - offsetY + (_windowHeight * pxToWorld)/2.0;
2002
2003    vtkSmartPointer<vtkCamera> camera = _renderer->GetActiveCamera();
2004    camera->ParallelProjectionOn();
2005    camera->SetClippingRange(1, 2);
2006    // Half of world coordinate height of viewport (Documentation is wrong)
2007    camera->SetParallelScale(_windowHeight * pxToWorld / 2.0);
2008
2009    if (_imgCameraPlane == PLANE_XY) {
2010        // XY plane
2011        camera->SetPosition(camPos[0], camPos[1], _imgCameraOffset + 1.);
2012        camera->SetFocalPoint(camPos[0], camPos[1], _imgCameraOffset);
2013        camera->SetViewUp(0, 1, 0);
2014        // bottom
2015        _cameraClipPlanes[0]->SetOrigin(0, _imgWorldOrigin[1], 0);
2016        _cameraClipPlanes[0]->SetNormal(0, 1, 0);
2017        // left
2018        _cameraClipPlanes[1]->SetOrigin(_imgWorldOrigin[0], 0, 0);
2019        _cameraClipPlanes[1]->SetNormal(1, 0, 0);
2020        // top
2021        _cameraClipPlanes[2]->SetOrigin(0, _imgWorldOrigin[1] + _imgWorldDims[1], 0);
2022        _cameraClipPlanes[2]->SetNormal(0, -1, 0);
2023        // right
2024        _cameraClipPlanes[3]->SetOrigin(_imgWorldOrigin[0] + _imgWorldDims[0], 0, 0);
2025        _cameraClipPlanes[3]->SetNormal(-1, 0, 0);
2026        _cubeAxesActor2D->SetBounds(_imgWorldOrigin[0], _imgWorldOrigin[0] + _imgWorldDims[0],
2027                                    _imgWorldOrigin[1], _imgWorldOrigin[1] + _imgWorldDims[1],
2028                                    _imgCameraOffset, _imgCameraOffset);
2029        _cubeAxesActor2D->XAxisVisibilityOn();
2030        _cubeAxesActor2D->YAxisVisibilityOn();
2031        _cubeAxesActor2D->ZAxisVisibilityOff();
2032    } else if (_imgCameraPlane == PLANE_ZY) {
2033        // ZY plane
2034        camera->SetPosition(_imgCameraOffset - 1., camPos[1], camPos[0]);
2035        camera->SetFocalPoint(_imgCameraOffset, camPos[1], camPos[0]);
2036        camera->SetViewUp(0, 1, 0);
2037        // bottom
2038        _cameraClipPlanes[0]->SetOrigin(0, _imgWorldOrigin[1], 0);
2039        _cameraClipPlanes[0]->SetNormal(0, 1, 0);
2040        // left
2041        _cameraClipPlanes[1]->SetOrigin(0, 0, _imgWorldOrigin[0]);
2042        _cameraClipPlanes[1]->SetNormal(0, 0, 1);
2043        // top
2044        _cameraClipPlanes[2]->SetOrigin(0, _imgWorldOrigin[1] + _imgWorldDims[1], 0);
2045        _cameraClipPlanes[2]->SetNormal(0, -1, 0);
2046        // right
2047        _cameraClipPlanes[3]->SetOrigin(0, 0, _imgWorldOrigin[0] + _imgWorldDims[0]);
2048        _cameraClipPlanes[3]->SetNormal(0, 0, -1);
2049        _cubeAxesActor2D->SetBounds(_imgCameraOffset, _imgCameraOffset,
2050                                    _imgWorldOrigin[1], _imgWorldOrigin[1] + _imgWorldDims[1],
2051                                    _imgWorldOrigin[0], _imgWorldOrigin[0] + _imgWorldDims[0]);
2052        _cubeAxesActor2D->XAxisVisibilityOff();
2053        _cubeAxesActor2D->YAxisVisibilityOn();
2054        _cubeAxesActor2D->ZAxisVisibilityOn();
2055    } else {
2056        // XZ plane
2057        camera->SetPosition(camPos[0], _imgCameraOffset - 1., camPos[1]);
2058        camera->SetFocalPoint(camPos[0], _imgCameraOffset, camPos[1]);
2059        camera->SetViewUp(0, 0, 1);
2060        // bottom
2061        _cameraClipPlanes[0]->SetOrigin(0, 0, _imgWorldOrigin[1]);
2062        _cameraClipPlanes[0]->SetNormal(0, 0, 1);
2063        // left
2064        _cameraClipPlanes[1]->SetOrigin(_imgWorldOrigin[0], 0, 0);
2065        _cameraClipPlanes[1]->SetNormal(1, 0, 0);
2066        // top
2067        _cameraClipPlanes[2]->SetOrigin(0, 0, _imgWorldOrigin[1] + _imgWorldDims[1]);
2068        _cameraClipPlanes[2]->SetNormal(0, 0, -1);
2069        // right
2070        _cameraClipPlanes[3]->SetOrigin(_imgWorldOrigin[0] + _imgWorldDims[0], 0, 0);
2071        _cameraClipPlanes[3]->SetNormal(-1, 0, 0);
2072        _cubeAxesActor2D->SetBounds(_imgWorldOrigin[0], _imgWorldOrigin[0] + _imgWorldDims[0],
2073                                    _imgCameraOffset, _imgCameraOffset,
2074                                    _imgWorldOrigin[1], _imgWorldOrigin[1] + _imgWorldDims[1]);
2075        _cubeAxesActor2D->XAxisVisibilityOn();
2076        _cubeAxesActor2D->YAxisVisibilityOff();
2077        _cubeAxesActor2D->ZAxisVisibilityOn();
2078    }
2079
2080    // Compute screen world coordinates
2081    computeScreenWorldCoords();
2082
2083#ifdef DEBUG
2084    printCameraInfo(camera);
2085#endif
2086
2087    _needsRedraw = true;
2088}
2089
2090/**
2091 * \brief Convert pixel/display coordinates to world coordinates based on current camera
2092 */
2093void Renderer::computeDisplayToWorld(double x, double y, double z, double worldPt[4])
2094{
2095    _renderer->SetDisplayPoint(x, y, z);
2096    _renderer->DisplayToWorld();
2097    _renderer->GetWorldPoint(worldPt);
2098    if (worldPt[3]) {
2099        worldPt[0] /= worldPt[3];
2100        worldPt[1] /= worldPt[3];
2101        worldPt[2] /= worldPt[3];
2102        worldPt[3] = 1.0;
2103    }
2104}
2105
2106/**
2107 * \brief Convert world coordinates to pixel/display coordinates based on current camera
2108 */
2109void Renderer::computeWorldToDisplay(double x, double y, double z, double displayPt[3])
2110{
2111    _renderer->SetWorldPoint(x, y, z, 1.0);
2112    _renderer->WorldToDisplay();
2113    _renderer->GetDisplayPoint(displayPt);
2114}
2115
2116/**
2117 * \brief Compute the world coordinate bounds of the display rectangle
2118 */
2119void Renderer::computeScreenWorldCoords()
2120{
2121    // Start with viewport coords [-1,1]
2122    double x0 = -1;
2123    double y0 = -1;
2124    double z0 = -1;
2125    double x1 = 1;
2126    double y1 = 1;
2127    double z1 = -1;
2128
2129    vtkMatrix4x4 *mat = vtkMatrix4x4::New();
2130    double result[4];
2131
2132    // get the perspective transformation from the active camera
2133    mat->DeepCopy(_renderer->GetActiveCamera()->
2134                  GetCompositeProjectionTransformMatrix(_renderer->GetTiledAspectRatio(),0,1));
2135
2136    // use the inverse matrix
2137    mat->Invert();
2138
2139    // Transform point to world coordinates
2140    result[0] = x0;
2141    result[1] = y0;
2142    result[2] = z0;
2143    result[3] = 1.0;
2144
2145    mat->MultiplyPoint(result, result);
2146
2147    // Get the transformed vector & set WorldPoint
2148    // while we are at it try to keep w at one
2149    if (result[3]) {
2150        x0 = result[0] / result[3];
2151        y0 = result[1] / result[3];
2152        z0 = result[2] / result[3];
2153    }
2154
2155    result[0] = x1;
2156    result[1] = y1;
2157    result[2] = z1;
2158    result[3] = 1.0;
2159
2160    mat->MultiplyPoint(result, result);
2161
2162    if (result[3]) {
2163        x1 = result[0] / result[3];
2164        y1 = result[1] / result[3];
2165        z1 = result[2] / result[3];
2166    }
2167
2168    mat->Delete();
2169
2170    if (_imgCameraPlane == PLANE_XZ) {
2171        _screenWorldCoords[0] = x0;
2172        _screenWorldCoords[1] = z0;
2173        _screenWorldCoords[2] = x1 - x0;
2174        _screenWorldCoords[3] = z1 - z0;
2175    } else if (_imgCameraPlane == PLANE_ZY) {
2176        _screenWorldCoords[0] = z0;
2177        _screenWorldCoords[1] = y0;
2178        _screenWorldCoords[2] = z1 - z0;
2179        _screenWorldCoords[3] = y1 - y0;
2180    } else {
2181        // XY
2182        _screenWorldCoords[0] = x0;
2183        _screenWorldCoords[1] = y0;
2184        _screenWorldCoords[2] = x1 - x0;
2185        _screenWorldCoords[3] = y1 - y0;
2186    }
2187}
2188
2189/**
2190 * \brief Get the world coordinates of the image camera plot area
2191 *
2192 * \param[out] xywh Array to hold x,y,width,height world coordinates
2193 */
2194void Renderer::getCameraZoomRegion(double xywh[4]) const
2195{
2196    xywh[0] = _imgWorldOrigin[0];
2197    xywh[1] = _imgWorldOrigin[1];
2198    xywh[2] = _imgWorldDims[0];
2199    xywh[3] = _imgWorldDims[1];
2200}
2201
2202/**
2203 * \brief Get the world origin and dimensions of the screen
2204 *
2205 * \param[out] xywh Array to hold x,y,width,height world coordinates
2206 */
2207void Renderer::getScreenWorldCoords(double xywh[4]) const
2208{
2209    memcpy(xywh, _screenWorldCoords, sizeof(double)*4);
2210}
2211
2212/**
2213 * \brief Compute bounding box containing the two input bounding boxes
2214 *
2215 * \param[out] boundsDest Union of the two input bounding boxes
2216 * \param[in] bounds1 Input bounding box
2217 * \param[in] bounds2 Input bounding box
2218 */
2219void Renderer::mergeBounds(double *boundsDest,
2220                           const double *bounds1, const double *bounds2)
2221{
2222    assert(boundsDest != NULL);
2223    assert(bounds1 != NULL);
2224    if (bounds2 == NULL) {
2225        WARN("NULL bounds2 array");
2226        return;
2227    }
2228    for (int i = 0; i < 6; i++) {
2229        if (i % 2 == 0)
2230            boundsDest[i] = min2(bounds1[i], bounds2[i]);
2231        else
2232            boundsDest[i] = max2(bounds1[i], bounds2[i]);
2233    }
2234}
2235
2236/**
2237 * \brief Collect bounds of all graphics objects
2238 *
2239 * \param[out] bounds Bounds of all scene objects
2240 * \param[in] onlyVisible Only collect bounds of visible objects
2241 */
2242void Renderer::collectBounds(double *bounds, bool onlyVisible)
2243{
2244    bounds[0] = DBL_MAX;
2245    bounds[1] = -DBL_MAX;
2246    bounds[2] = DBL_MAX;
2247    bounds[3] = -DBL_MAX;
2248    bounds[4] = DBL_MAX;
2249    bounds[5] = -DBL_MAX;
2250
2251    for (DataSetHashmap::iterator itr = _dataSets.begin();
2252             itr != _dataSets.end(); ++itr) {
2253        if ((!onlyVisible || itr->second->getVisibility()) &&
2254            itr->second->getProp() != NULL)
2255            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2256    }
2257    for (Contour2DHashmap::iterator itr = _contour2Ds.begin();
2258             itr != _contour2Ds.end(); ++itr) {
2259        if ((!onlyVisible || itr->second->getVisibility()) &&
2260            itr->second->getProp() != NULL)
2261            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2262    }
2263    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
2264             itr != _contour3Ds.end(); ++itr) {
2265        if ((!onlyVisible || itr->second->getVisibility()) &&
2266            itr->second->getProp() != NULL)
2267            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2268    }
2269    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
2270             itr != _cutplanes.end(); ++itr) {
2271        if ((!onlyVisible || itr->second->getVisibility()) &&
2272            itr->second->getProp() != NULL)
2273            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2274    }
2275    for (GlyphsHashmap::iterator itr = _glyphs.begin();
2276             itr != _glyphs.end(); ++itr) {
2277        if ((!onlyVisible || itr->second->getVisibility()) &&
2278            itr->second->getProp() != NULL)
2279            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2280    }
2281    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
2282             itr != _heightMaps.end(); ++itr) {
2283        if ((!onlyVisible || itr->second->getVisibility()) &&
2284            itr->second->getProp() != NULL)
2285            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2286    }
2287    for (LICHashmap::iterator itr = _lics.begin();
2288             itr != _lics.end(); ++itr) {
2289        if ((!onlyVisible || itr->second->getVisibility()) &&
2290            itr->second->getProp() != NULL)
2291            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2292    }
2293    for (MoleculeHashmap::iterator itr = _molecules.begin();
2294             itr != _molecules.end(); ++itr) {
2295        if ((!onlyVisible || itr->second->getVisibility()) &&
2296            itr->second->getProp() != NULL)
2297            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2298    }
2299    for (PolyDataHashmap::iterator itr = _polyDatas.begin();
2300             itr != _polyDatas.end(); ++itr) {
2301        if ((!onlyVisible || itr->second->getVisibility()) &&
2302            itr->second->getProp() != NULL)
2303            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2304    }
2305    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
2306             itr != _pseudoColors.end(); ++itr) {
2307        if ((!onlyVisible || itr->second->getVisibility()) &&
2308            itr->second->getProp() != NULL)
2309            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2310    }
2311    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
2312             itr != _streamlines.end(); ++itr) {
2313        if ((!onlyVisible || itr->second->getVisibility()) &&
2314            itr->second->getProp() != NULL)
2315            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2316    }
2317    for (VolumeHashmap::iterator itr = _volumes.begin();
2318             itr != _volumes.end(); ++itr) {
2319        if ((!onlyVisible || itr->second->getVisibility()) &&
2320            itr->second->getProp() != NULL)
2321            mergeBounds(bounds, bounds, itr->second->getProp()->GetBounds());
2322    }
2323
2324    for (int i = 0; i < 6; i += 2) {
2325        if (bounds[i+1] < bounds[i]) {
2326            bounds[i] = -0.5;
2327            bounds[i+1] = 0.5;
2328        }
2329    }
2330
2331    int numDims = 0;
2332    if (bounds[0] != bounds[1])
2333        numDims++;
2334    if (bounds[2] != bounds[3])
2335        numDims++;
2336    if (bounds[4] != bounds[5])
2337        numDims++;
2338
2339    if (numDims == 0) {
2340        bounds[0] -= .5;
2341        bounds[1] += .5;
2342        bounds[2] -= .5;
2343        bounds[3] += .5;
2344    }
2345
2346    TRACE("Bounds: %g %g %g %g %g %g",
2347          bounds[0],
2348          bounds[1],
2349          bounds[2],
2350          bounds[3],
2351          bounds[4],
2352          bounds[5]);
2353}
2354
2355/**
2356 * \brief Update data ranges for color-mapping and contours
2357 */
2358void Renderer::updateFieldRanges()
2359{
2360    collectDataRanges();
2361
2362    for (Contour2DHashmap::iterator itr = _contour2Ds.begin();
2363         itr != _contour2Ds.end(); ++itr) {
2364        itr->second->updateRanges(this);
2365    }
2366    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
2367         itr != _contour3Ds.end(); ++itr) {
2368        itr->second->updateRanges(this);
2369    }
2370    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
2371         itr != _cutplanes.end(); ++itr) {
2372        itr->second->updateRanges(this);
2373    }
2374    for (GlyphsHashmap::iterator itr = _glyphs.begin();
2375         itr != _glyphs.end(); ++itr) {
2376        itr->second->updateRanges(this);
2377    }
2378    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
2379         itr != _heightMaps.end(); ++itr) {
2380        itr->second->updateRanges(this);
2381    }
2382    for (LICHashmap::iterator itr = _lics.begin();
2383         itr != _lics.end(); ++itr) {
2384        itr->second->updateRanges(this);
2385    }
2386    for (MoleculeHashmap::iterator itr = _molecules.begin();
2387         itr != _molecules.end(); ++itr) {
2388        itr->second->updateRanges(this);
2389    }
2390    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
2391         itr != _pseudoColors.end(); ++itr) {
2392        itr->second->updateRanges(this);
2393    }
2394    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
2395         itr != _streamlines.end(); ++itr) {
2396        itr->second->updateRanges(this);
2397    }
2398    for (VolumeHashmap::iterator itr = _volumes.begin();
2399         itr != _volumes.end(); ++itr) {
2400        itr->second->updateRanges(this);
2401    }
2402}
2403
2404/**
2405 * \brief Collect cumulative data range of all DataSets
2406 *
2407 * \param[out] range Data range of all DataSets
2408 * \param[in] name Field name
2409 * \param[in] type Attribute type: e.g. POINT_DATA, CELL_DATA
2410 * \param[in] component Array component or -1 for magnitude
2411 * \param[in] onlyVisible Only collect range of visible DataSets
2412 */
2413void Renderer::collectDataRanges(double *range, const char *name,
2414                                 DataSet::DataAttributeType type,
2415                                 int component, bool onlyVisible)
2416{
2417    range[0] = DBL_MAX;
2418    range[1] = -DBL_MAX;
2419
2420    for (DataSetHashmap::iterator itr = _dataSets.begin();
2421         itr != _dataSets.end(); ++itr) {
2422        if (!onlyVisible || itr->second->getVisibility()) {
2423            double r[2];
2424            itr->second->getDataRange(r, name, type, component);
2425            range[0] = min2(range[0], r[0]);
2426            range[1] = max2(range[1], r[1]);
2427        }
2428    }
2429    if (range[0] == DBL_MAX)
2430        range[0] = 0;
2431    if (range[1] == -DBL_MAX)
2432        range[1] = 1;
2433   
2434    TRACE("n:'%s' t:%d c:%d [%g,%g]", name, type, component,
2435          range[0], range[1]);
2436}
2437
2438/**
2439 * \brief Clear field range hashtables and free memory
2440 */
2441void Renderer::clearFieldRanges()
2442{
2443    TRACE("Deleting Field Ranges");
2444    for (FieldRangeHashmap::iterator itr = _scalarPointDataRange.begin();
2445         itr != _scalarPointDataRange.end(); ++itr) {
2446        delete [] itr->second;
2447    }
2448    _scalarPointDataRange.clear();
2449    for (FieldRangeHashmap::iterator itr = _scalarCellDataRange.begin();
2450         itr != _scalarCellDataRange.end(); ++itr) {
2451        delete [] itr->second;
2452    }
2453    _scalarCellDataRange.clear();
2454    for (FieldRangeHashmap::iterator itr = _vectorPointDataRange.begin();
2455         itr != _vectorPointDataRange.end(); ++itr) {
2456        delete [] itr->second;
2457    }
2458    _vectorPointDataRange.clear();
2459    for (int i = 0; i < 3; i++) {
2460        for (FieldRangeHashmap::iterator itr = _vectorCompPointDataRange[i].begin();
2461             itr != _vectorCompPointDataRange[i].end(); ++itr) {
2462            delete [] itr->second;
2463        }
2464        _vectorCompPointDataRange[i].clear();
2465    }
2466    for (FieldRangeHashmap::iterator itr = _vectorCellDataRange.begin();
2467         itr != _vectorCellDataRange.end(); ++itr) {
2468        delete [] itr->second;
2469    }
2470    _vectorCellDataRange.clear();
2471    for (int i = 0; i < 3; i++) {
2472        for (FieldRangeHashmap::iterator itr = _vectorCompCellDataRange[i].begin();
2473             itr != _vectorCompCellDataRange[i].end(); ++itr) {
2474            delete [] itr->second;
2475        }
2476        _vectorCompCellDataRange[i].clear();
2477    }
2478}
2479
2480void Renderer::initFieldRanges()
2481{
2482    clearFieldRanges();
2483
2484    for (DataSetHashmap::iterator itr = _dataSets.begin();
2485         itr != _dataSets.end(); ++itr) {
2486        DataSet *ds = itr->second;
2487        std::vector<std::string> names;
2488        ds->getFieldNames(names, DataSet::POINT_DATA, 1);
2489        for (std::vector<std::string>::iterator itr = names.begin();
2490             itr != names.end(); ++itr) {
2491            FieldRangeHashmap::iterator fritr =
2492                _scalarPointDataRange.find(*itr);
2493            if (fritr == _scalarPointDataRange.end()) {
2494                _scalarPointDataRange[*itr] = new double[2];
2495                _scalarPointDataRange[*itr][0] = 0;
2496                _scalarPointDataRange[*itr][1] = 1;
2497            }
2498        }
2499        names.clear();
2500        ds->getFieldNames(names, DataSet::CELL_DATA, 1);
2501        for (std::vector<std::string>::iterator itr = names.begin();
2502             itr != names.end(); ++itr) {
2503            FieldRangeHashmap::iterator fritr =
2504                _scalarCellDataRange.find(*itr);
2505            if (fritr == _scalarCellDataRange.end()) {
2506                _scalarCellDataRange[*itr] = new double[2];
2507                _scalarCellDataRange[*itr][0] = 0;
2508                _scalarCellDataRange[*itr][1] = 1;
2509            }
2510        }
2511        names.clear();
2512        ds->getFieldNames(names, DataSet::POINT_DATA, 3);
2513        for (std::vector<std::string>::iterator itr = names.begin();
2514             itr != names.end(); ++itr) {
2515            FieldRangeHashmap::iterator fritr =
2516                _vectorPointDataRange.find(*itr);
2517            if (fritr == _vectorPointDataRange.end()) {
2518                _vectorPointDataRange[*itr] = new double[2];
2519                _vectorPointDataRange[*itr][0] = 0;
2520                _vectorPointDataRange[*itr][1] = 1;
2521            }
2522            for (int i = 0; i < 3; i++) {
2523                fritr = _vectorCompPointDataRange[i].find(*itr);
2524                if (fritr == _vectorCompPointDataRange[i].end()) {
2525                    _vectorCompPointDataRange[i][*itr] = new double[2];
2526                    _vectorCompPointDataRange[i][*itr][0] = 0;
2527                    _vectorCompPointDataRange[i][*itr][1] = 1;
2528                }
2529            }
2530        }
2531        names.clear();
2532        ds->getFieldNames(names, DataSet::CELL_DATA, 3);
2533        for (std::vector<std::string>::iterator itr = names.begin();
2534             itr != names.end(); ++itr) {
2535            FieldRangeHashmap::iterator fritr =
2536                _vectorCellDataRange.find(*itr);
2537            if (fritr == _vectorCellDataRange.end()) {
2538                _vectorCellDataRange[*itr] = new double[2];
2539                _vectorCellDataRange[*itr][0] = 0;
2540                _vectorCellDataRange[*itr][1] = 1;
2541            }
2542            for (int i = 0; i < 3; i++) {
2543                fritr = _vectorCompCellDataRange[i].find(*itr);
2544                if (fritr == _vectorCompCellDataRange[i].end()) {
2545                    _vectorCompCellDataRange[i][*itr] = new double[2];
2546                    _vectorCompCellDataRange[i][*itr][0] = 0;
2547                    _vectorCompCellDataRange[i][*itr][1] = 1;
2548                }
2549            }
2550        }
2551    }
2552}
2553
2554/**
2555 * \brief Returns boolean flag indicating if cumulative data ranges
2556 * should be used
2557 */
2558bool Renderer::getUseCumulativeRange()
2559{
2560    return _useCumulativeRange;
2561}
2562
2563/**
2564 * \brief Get the cumulative range across all DataSets for a point
2565 * data field if it exists, otherwise a cell data field if it exists
2566 *
2567 * \param[out] range Pointer to an array of 2 doubles
2568 * \param[in] name Field name
2569 * \param[in] numComponents Number of components in field
2570 * \param[in] component Index of component or -1 for magnitude/scalar
2571 * \return boolean indicating if field was found
2572 */
2573bool Renderer::getCumulativeDataRange(double *range, const char *name,
2574                                      int numComponents,
2575                                      int component)
2576{
2577    bool ret;
2578    if (!(ret = getCumulativeDataRange(range, name, DataSet::POINT_DATA,
2579                                       numComponents, component)))
2580        ret = getCumulativeDataRange(range, name, DataSet::CELL_DATA,
2581                                     numComponents, component);
2582    return ret;
2583}
2584
2585/**
2586 * \brief Get the cumulative range across all DataSets for a field
2587 *
2588 * \param[out] range Pointer to an array of 2 doubles
2589 * \param[in] name Field name
2590 * \param[in] type DataAttributeType of field
2591 * \param[in] numComponents Number of components in field
2592 * \param[in] component Index of component or -1 for magnitude/scalar
2593 * \return boolean indicating if field was found
2594 */
2595bool Renderer::getCumulativeDataRange(double *range, const char *name,
2596                                      DataSet::DataAttributeType type,
2597                                      int numComponents,
2598                                      int component)
2599{
2600    if (name == NULL)
2601        return false;
2602    switch (type) {
2603    case DataSet::POINT_DATA: {
2604        FieldRangeHashmap::iterator itr;
2605        if (numComponents == 1) {
2606            itr = _scalarPointDataRange.find(name);
2607            if (itr == _scalarPointDataRange.end()) {
2608                return false;
2609            } else {
2610                memcpy(range, itr->second, sizeof(double)*2);
2611                return true;
2612            }
2613        } else if (numComponents == 3) {
2614            if (component == -1) {
2615                itr = _vectorPointDataRange.find(name);
2616                if (itr == _vectorPointDataRange.end()) {
2617                    return false;
2618                } else {
2619                    memcpy(range, itr->second, sizeof(double)*2);
2620                    return true;
2621                }
2622            } else if (component >= 0 && component <= 3) {
2623                itr = _vectorCompPointDataRange[component].find(name);
2624                if (itr == _vectorCompPointDataRange[component].end()) {
2625                    return false;
2626                } else {
2627                    memcpy(range, itr->second, sizeof(double)*2);
2628                    return true;
2629                }
2630            }
2631        }
2632    }
2633        break;
2634    case DataSet::CELL_DATA: {
2635        FieldRangeHashmap::iterator itr;
2636        if (numComponents == 1) {
2637            itr = _scalarCellDataRange.find(name);
2638            if (itr == _scalarCellDataRange.end()) {
2639                return false;
2640            } else {
2641                memcpy(range, itr->second, sizeof(double)*2);
2642                return true;
2643            }
2644        } else if (numComponents == 3) {
2645            if (component == -1) {
2646                itr = _vectorCellDataRange.find(name);
2647                if (itr == _vectorCellDataRange.end()) {
2648                    return false;
2649                } else {
2650                    memcpy(range, itr->second, sizeof(double)*2);
2651                    return true;
2652                }
2653            } else if (component >= 0 && component <= 3) {
2654                itr = _vectorCompCellDataRange[component].find(name);
2655                if (itr == _vectorCompCellDataRange[component].end()) {
2656                    return false;
2657                } else {
2658                    memcpy(range, itr->second, sizeof(double)*2);
2659                    return true;
2660                }
2661            }
2662        }
2663    }
2664        break;
2665    default:
2666        break;
2667    }
2668    return false;
2669}
2670
2671void Renderer::collectDataRanges()
2672{
2673    for (FieldRangeHashmap::iterator itr = _scalarPointDataRange.begin();
2674         itr != _scalarPointDataRange.end(); ++itr) {
2675        collectDataRanges(itr->second, itr->first.c_str(),
2676                          DataSet::POINT_DATA, -1,
2677                          _cumulativeRangeOnlyVisible);
2678    }
2679    for (FieldRangeHashmap::iterator itr = _scalarCellDataRange.begin();
2680         itr != _scalarCellDataRange.end(); ++itr) {
2681        collectDataRanges(itr->second, itr->first.c_str(),
2682                          DataSet::CELL_DATA, -1,
2683                          _cumulativeRangeOnlyVisible);
2684    }
2685    for (FieldRangeHashmap::iterator itr = _vectorPointDataRange.begin();
2686         itr != _vectorPointDataRange.end(); ++itr) {
2687        collectDataRanges(itr->second, itr->first.c_str(),
2688                          DataSet::POINT_DATA, -1,
2689                          _cumulativeRangeOnlyVisible);
2690    }
2691    for (int i = 0; i < 3; i++) {
2692        for (FieldRangeHashmap::iterator itr = _vectorCompPointDataRange[i].begin();
2693             itr != _vectorCompPointDataRange[i].end(); ++itr) {
2694            collectDataRanges(itr->second, itr->first.c_str(),
2695                              DataSet::POINT_DATA, i,
2696                              _cumulativeRangeOnlyVisible);
2697        }
2698    }
2699    for (FieldRangeHashmap::iterator itr = _vectorCellDataRange.begin();
2700         itr != _vectorCellDataRange.end(); ++itr) {
2701        collectDataRanges(itr->second, itr->first.c_str(),
2702                          DataSet::CELL_DATA, -1,
2703                          _cumulativeRangeOnlyVisible);
2704    }
2705    for (int i = 0; i < 3; i++) {
2706        for (FieldRangeHashmap::iterator itr = _vectorCompCellDataRange[i].begin();
2707             itr != _vectorCompCellDataRange[i].end(); ++itr) {
2708            collectDataRanges(itr->second, itr->first.c_str(),
2709                              DataSet::CELL_DATA, i,
2710                              _cumulativeRangeOnlyVisible);
2711        }
2712    }
2713}
2714
2715/**
2716 * \brief Determines if AABB lies in a principal axis plane
2717 * and if so, returns the plane normal
2718 */
2719bool Renderer::is2D(const double bounds[6],
2720                    PrincipalPlane *plane,
2721                    double *offset) const
2722{
2723    if (bounds[4] == bounds[5]) {
2724        // Z = 0, XY plane
2725        if (plane)
2726            *plane = PLANE_XY;
2727        if (offset)
2728            *offset = bounds[4];
2729        return true;
2730    } else if (bounds[0] == bounds[1]) {
2731        // X = 0, ZY plane
2732        if (plane)
2733            *plane = PLANE_ZY;
2734        if (offset)
2735            *offset = bounds[0];
2736        return true;
2737    } else if (bounds[2] == bounds[3]) {
2738        // Y = 0, XZ plane
2739        if (plane)
2740            *plane = PLANE_XZ;
2741        if (offset)
2742            *offset = bounds[2];
2743        return true;
2744    }
2745    *plane = PLANE_XY;
2746    *offset = 0;
2747    return false;
2748}
2749
2750/**
2751 * \brief Initialize the camera zoom region to include the bounding volume given
2752 */
2753void Renderer::initCamera()
2754{
2755#ifdef WANT_TRACE
2756    switch (_cameraMode) {
2757    case IMAGE:
2758        TRACE("Image camera");
2759        break;
2760    case ORTHO:
2761        TRACE("Ortho camera");
2762        break;
2763    case PERSPECTIVE:
2764        TRACE("Perspective camera");
2765        break;
2766    default:
2767        TRACE("Unknown camera mode");
2768    }
2769#endif
2770    double bounds[6];
2771    collectBounds(bounds, false);
2772    bool twod = is2D(bounds, &_imgCameraPlane, &_imgCameraOffset);
2773    if (twod) {
2774        _cameraMode = IMAGE;
2775        if (_imgCameraPlane == PLANE_ZY) {
2776            _imgWorldOrigin[0] = bounds[4];
2777            _imgWorldOrigin[1] = bounds[2];
2778            _imgWorldDims[0] = bounds[5] - bounds[4];
2779            _imgWorldDims[1] = bounds[3] - bounds[2];
2780        } else if (_imgCameraPlane == PLANE_XZ) {
2781            _imgWorldOrigin[0] = bounds[0];
2782            _imgWorldOrigin[1] = bounds[4];
2783            _imgWorldDims[0] = bounds[1] - bounds[0];
2784            _imgWorldDims[1] = bounds[5] - bounds[4];
2785        } else {
2786            _imgWorldOrigin[0] = bounds[0];
2787            _imgWorldOrigin[1] = bounds[2];
2788            _imgWorldDims[0] = bounds[1] - bounds[0];
2789            _imgWorldDims[1] = bounds[3] - bounds[2];
2790        }
2791    } else {
2792        _imgWorldOrigin[0] = bounds[0];
2793        _imgWorldOrigin[1] = bounds[2];
2794        _imgWorldDims[0] = bounds[1] - bounds[0];
2795        _imgWorldDims[1] = bounds[3] - bounds[2];
2796    }
2797
2798    _cameraPan[0] = 0;
2799    _cameraPan[1] = 0;
2800    _cameraZoomRatio = 1;
2801
2802    switch (_cameraMode) {
2803    case IMAGE:
2804        _renderer->ResetCamera(bounds);
2805        setCameraZoomRegion(_imgWorldOrigin[0], _imgWorldOrigin[1],
2806                            _imgWorldDims[0], _imgWorldDims[1]);
2807        resetAxes();
2808        break;
2809    case ORTHO:
2810        _renderer->GetActiveCamera()->ParallelProjectionOn();
2811        resetAxes(bounds);
2812        _renderer->ResetCamera(bounds);
2813        //computeScreenWorldCoords();
2814        break;
2815    case PERSPECTIVE:
2816        _renderer->GetActiveCamera()->ParallelProjectionOff();
2817        resetAxes(bounds);
2818        _renderer->ResetCamera(bounds);
2819        //computeScreenWorldCoords();
2820        break;
2821    default:
2822        ERROR("Unknown camera mode");
2823    }
2824
2825#ifdef WANT_TRACE
2826    printCameraInfo(_renderer->GetActiveCamera());
2827#endif
2828}
2829
2830/**
2831 * \brief Print debugging info about a vtkCamera
2832 */
2833void Renderer::printCameraInfo(vtkCamera *camera)
2834{
2835    TRACE("pscale: %g, angle: %g, d: %g pos: %g %g %g, fpt: %g %g %g, vup: %g %g %g, clip: %g %g",
2836          camera->GetParallelScale(),
2837          camera->GetViewAngle(),
2838          camera->GetDistance(),
2839          camera->GetPosition()[0],
2840          camera->GetPosition()[1],
2841          camera->GetPosition()[2],
2842          camera->GetFocalPoint()[0],
2843          camera->GetFocalPoint()[1],
2844          camera->GetFocalPoint()[2],
2845          camera->GetViewUp()[0],
2846          camera->GetViewUp()[1],
2847          camera->GetViewUp()[2],
2848          camera->GetClippingRange()[0],
2849          camera->GetClippingRange()[1]);
2850}
2851
2852/**
2853 * \brief Set the RGB background color to render into the image
2854 */
2855void Renderer::setBackgroundColor(float color[3])
2856{
2857    _bgColor[0] = color[0];
2858    _bgColor[1] = color[1];
2859    _bgColor[2] = color[2];
2860    _renderer->SetBackground(_bgColor[0], _bgColor[1], _bgColor[2]);
2861    _needsRedraw = true;
2862}
2863
2864/**
2865 * \brief Set the opacity of the specified DataSet's associated graphics objects
2866 */
2867void Renderer::setDataSetOpacity(const DataSetId& id, double opacity)
2868{
2869    DataSetHashmap::iterator itr;
2870
2871    bool doAll = false;
2872
2873    if (id.compare("all") == 0) {
2874        itr = _dataSets.begin();
2875        doAll = true;
2876    } else {
2877        itr = _dataSets.find(id);
2878    }
2879    if (itr == _dataSets.end()) {
2880        ERROR("Unknown dataset %s", id.c_str());
2881        return;
2882    }
2883
2884    do {
2885        itr->second->setOpacity(opacity);
2886    } while (doAll && ++itr != _dataSets.end());
2887
2888    if (id.compare("all") == 0 || getGraphicsObject<Contour2D>(id) != NULL)
2889        setGraphicsObjectOpacity<Contour2D>(id, opacity);
2890    if (id.compare("all") == 0 || getGraphicsObject<Contour3D>(id) != NULL)
2891        setGraphicsObjectOpacity<Contour3D>(id, opacity);
2892    if (id.compare("all") == 0 || getGraphicsObject<Cutplane>(id) != NULL)
2893        setGraphicsObjectOpacity<Cutplane>(id, opacity);
2894    if (id.compare("all") == 0 || getGraphicsObject<Glyphs>(id) != NULL)
2895        setGraphicsObjectOpacity<Glyphs>(id, opacity);
2896    if (id.compare("all") == 0 || getGraphicsObject<HeightMap>(id) != NULL)
2897        setGraphicsObjectOpacity<HeightMap>(id, opacity);
2898    if (id.compare("all") == 0 || getGraphicsObject<LIC>(id) != NULL)
2899        setGraphicsObjectOpacity<LIC>(id, opacity);
2900    if (id.compare("all") == 0 || getGraphicsObject<Molecule>(id) != NULL)
2901        setGraphicsObjectOpacity<Molecule>(id, opacity);
2902    if (id.compare("all") == 0 || getGraphicsObject<PolyData>(id) != NULL)
2903        setGraphicsObjectOpacity<PolyData>(id, opacity);
2904    if (id.compare("all") == 0 || getGraphicsObject<PseudoColor>(id) != NULL)
2905        setGraphicsObjectOpacity<PseudoColor>(id, opacity);
2906    if (id.compare("all") == 0 || getGraphicsObject<Streamlines>(id) != NULL)
2907        setGraphicsObjectOpacity<Streamlines>(id, opacity);
2908    if (id.compare("all") == 0 || getGraphicsObject<Volume>(id) != NULL)
2909        setGraphicsObjectOpacity<Volume>(id, opacity);
2910
2911    _needsRedraw = true;
2912}
2913
2914/**
2915 * \brief Turn on/off rendering of the specified DataSet's associated graphics objects
2916 */
2917void Renderer::setDataSetVisibility(const DataSetId& id, bool state)
2918{
2919    DataSetHashmap::iterator itr;
2920
2921    bool doAll = false;
2922
2923    if (id.compare("all") == 0) {
2924        itr = _dataSets.begin();
2925        doAll = true;
2926    } else {
2927        itr = _dataSets.find(id);
2928    }
2929    if (itr == _dataSets.end()) {
2930        ERROR("Unknown dataset %s", id.c_str());
2931        return;
2932    }
2933
2934    do {
2935        itr->second->setVisibility(state);
2936    } while (doAll && ++itr != _dataSets.end());
2937
2938    if (id.compare("all") == 0 || getGraphicsObject<Contour2D>(id) != NULL)
2939        setGraphicsObjectVisibility<Contour2D>(id, state);
2940    if (id.compare("all") == 0 || getGraphicsObject<Contour3D>(id) != NULL)
2941        setGraphicsObjectVisibility<Contour3D>(id, state);
2942    if (id.compare("all") == 0 || getGraphicsObject<Cutplane>(id) != NULL)
2943        setGraphicsObjectVisibility<Cutplane>(id, state);
2944    if (id.compare("all") == 0 || getGraphicsObject<Glyphs>(id) != NULL)
2945        setGraphicsObjectVisibility<Glyphs>(id, state);
2946    if (id.compare("all") == 0 || getGraphicsObject<HeightMap>(id) != NULL)
2947        setGraphicsObjectVisibility<HeightMap>(id, state);
2948    if (id.compare("all") == 0 || getGraphicsObject<LIC>(id) != NULL)
2949        setGraphicsObjectVisibility<LIC>(id, state);
2950    if (id.compare("all") == 0 || getGraphicsObject<Molecule>(id) != NULL)
2951        setGraphicsObjectVisibility<Molecule>(id, state);
2952    if (id.compare("all") == 0 || getGraphicsObject<PolyData>(id) != NULL)
2953        setGraphicsObjectVisibility<PolyData>(id, state);
2954    if (id.compare("all") == 0 || getGraphicsObject<PseudoColor>(id) != NULL)
2955        setGraphicsObjectVisibility<PseudoColor>(id, state);
2956    if (id.compare("all") == 0 || getGraphicsObject<Streamlines>(id) != NULL)
2957        setGraphicsObjectVisibility<Streamlines>(id, state);
2958    if (id.compare("all") == 0 || getGraphicsObject<Volume>(id) != NULL)
2959        setGraphicsObjectVisibility<Volume>(id, state);
2960
2961    _needsRedraw = true;
2962}
2963
2964/**
2965 * \brief Toggle rendering of actors' bounding box
2966 */
2967void Renderer::setDataSetShowBounds(const DataSetId& id, bool state)
2968{
2969    DataSetHashmap::iterator itr;
2970
2971    bool doAll = false;
2972
2973    if (id.compare("all") == 0) {
2974        itr = _dataSets.begin();
2975        doAll = true;
2976    } else {
2977        itr = _dataSets.find(id);
2978    }
2979    if (itr == _dataSets.end()) {
2980        ERROR("Unknown dataset %s", id.c_str());
2981        return;
2982    }
2983
2984    do {
2985        if (!state && itr->second->getProp()) {
2986            _renderer->RemoveViewProp(itr->second->getProp());
2987        }
2988
2989        itr->second->showOutline(state);
2990
2991        if (state && !_renderer->HasViewProp(itr->second->getProp())) {
2992            _renderer->AddViewProp(itr->second->getProp());
2993        }
2994    } while (doAll && ++itr != _dataSets.end());
2995
2996    initCamera();
2997    _needsRedraw = true;
2998}
2999
3000/**
3001 * \brief Set color of outline bounding box
3002 */
3003void Renderer::setDataSetOutlineColor(const DataSetId& id, float color[3])
3004{
3005    DataSetHashmap::iterator itr;
3006
3007    bool doAll = false;
3008
3009    if (id.compare("all") == 0) {
3010        itr = _dataSets.begin();
3011        doAll = true;
3012    } else {
3013        itr = _dataSets.find(id);
3014    }
3015    if (itr == _dataSets.end()) {
3016        ERROR("Unknown dataset %s", id.c_str());
3017        return;
3018    }
3019
3020    do {
3021        itr->second->setOutlineColor(color);
3022    } while (doAll && ++itr != _dataSets.end());
3023
3024    _needsRedraw = true;
3025}
3026
3027/**
3028 * \brief Set a user clipping plane
3029 *
3030 * TODO: Fix clip plane positions after a change in actor bounds
3031 */
3032void Renderer::setClipPlane(Axis axis, double ratio, int direction)
3033{
3034    double bounds[6];
3035    collectBounds(bounds, false);
3036
3037    switch (axis) {
3038    case X_AXIS:
3039        if (direction > 0) {
3040            if (ratio > 0.0) {
3041                if (_userClipPlanes[0] == NULL) {
3042                    _userClipPlanes[0] = vtkSmartPointer<vtkPlane>::New();
3043                    _userClipPlanes[0]->SetNormal(1, 0, 0);
3044                }
3045                _userClipPlanes[0]->SetOrigin(bounds[0] + (bounds[1]-bounds[0])*ratio, 0, 0);
3046            } else {
3047                _userClipPlanes[0] = NULL;
3048            }
3049        } else {
3050            if (ratio < 1.0) {
3051                if (_userClipPlanes[1] == NULL) {
3052                    _userClipPlanes[1] = vtkSmartPointer<vtkPlane>::New();
3053                    _userClipPlanes[1]->SetNormal(-1, 0, 0);
3054                }
3055                _userClipPlanes[1]->SetOrigin(bounds[0] + (bounds[1]-bounds[0])*ratio, 0, 0);
3056            } else {
3057                _userClipPlanes[1] = NULL;
3058            }
3059        }
3060        break;
3061    case Y_AXIS:
3062        if (direction > 0) {
3063            if (ratio > 0.0) {
3064                if (_userClipPlanes[2] == NULL) {
3065                    _userClipPlanes[2] = vtkSmartPointer<vtkPlane>::New();
3066                    _userClipPlanes[2]->SetNormal(0, 1, 0);
3067                }
3068                _userClipPlanes[2]->SetOrigin(0, bounds[2] + (bounds[3]-bounds[2])*ratio, 0);
3069            } else {
3070                _userClipPlanes[2] = NULL;
3071            }
3072        } else {
3073            if (ratio < 1.0) {
3074                if (_userClipPlanes[3] == NULL) {
3075                    _userClipPlanes[3] = vtkSmartPointer<vtkPlane>::New();
3076                    _userClipPlanes[3]->SetNormal(0, -1, 0);
3077                }
3078                _userClipPlanes[3]->SetOrigin(0, bounds[2] + (bounds[3]-bounds[2])*ratio, 0);
3079            } else {
3080                _userClipPlanes[3] = NULL;
3081            }
3082        }
3083        break;
3084    case Z_AXIS:
3085        if (direction > 0) {
3086            if (ratio > 0.0) {
3087                if (_userClipPlanes[4] == NULL) {
3088                    _userClipPlanes[4] = vtkSmartPointer<vtkPlane>::New();
3089                    _userClipPlanes[4]->SetNormal(0, 0, 1);
3090                }
3091                _userClipPlanes[4]->SetOrigin(0, 0, bounds[4] + (bounds[5]-bounds[4])*ratio);
3092            } else {
3093                _userClipPlanes[4] = NULL;
3094            }
3095        } else {
3096            if (ratio < 1.0) {
3097                if (_userClipPlanes[5] == NULL) {
3098                    _userClipPlanes[5] = vtkSmartPointer<vtkPlane>::New();
3099                    _userClipPlanes[5]->SetNormal(0, 0, -1);
3100                }
3101                _userClipPlanes[5]->SetOrigin(0, 0, bounds[4] + (bounds[5]-bounds[4])*ratio);
3102            } else {
3103                _userClipPlanes[5] = NULL;
3104            }
3105        }
3106        break;
3107    default:
3108        ;
3109    }
3110
3111    _needsRedraw = true;
3112}
3113
3114/**
3115 * \brief Set up clipping planes for image camera mode if needed
3116 */
3117void Renderer::setCameraClippingPlanes()
3118{
3119    /* XXX: Note that there appears to be a bug with setting the
3120     * clipping plane collection to NULL in the VTK Mappers --
3121     * the old clip planes are still applied.  The workaround here
3122     * is to keep the PlaneCollection and add/remove the planes
3123     * to/from the PlaneCollection as needed.
3124     */
3125    if (_cameraMode == IMAGE) {
3126        if (_activeClipPlanes->GetNumberOfItems() == 0) {
3127            for (int i = 0; i < 4; i++)
3128                _activeClipPlanes->AddItem(_cameraClipPlanes[i]);
3129        }
3130    } else {
3131        if (_activeClipPlanes->GetNumberOfItems() > 0)
3132            _activeClipPlanes->RemoveAllItems();
3133        for (int i = 0; i < 6; i++) {
3134            if (_userClipPlanes[i] != NULL) {
3135                _activeClipPlanes->AddItem(_userClipPlanes[i]);
3136            }
3137        }
3138    }
3139
3140    /* Ensure all Mappers are using the PlaneCollection
3141     * This will not change the state or timestamp of
3142     * Mappers already using the PlaneCollection
3143     */
3144    for (Contour2DHashmap::iterator itr = _contour2Ds.begin();
3145         itr != _contour2Ds.end(); ++itr) {
3146        itr->second->setClippingPlanes(_activeClipPlanes);
3147    }
3148    for (Contour3DHashmap::iterator itr = _contour3Ds.begin();
3149         itr != _contour3Ds.end(); ++itr) {
3150        itr->second->setClippingPlanes(_activeClipPlanes);
3151    }
3152    for (CutplaneHashmap::iterator itr = _cutplanes.begin();
3153         itr != _cutplanes.end(); ++itr) {
3154        itr->second->setClippingPlanes(_activeClipPlanes);
3155    }
3156    for (GlyphsHashmap::iterator itr = _glyphs.begin();
3157         itr != _glyphs.end(); ++itr) {
3158        itr->second->setClippingPlanes(_activeClipPlanes);
3159    }
3160    for (HeightMapHashmap::iterator itr = _heightMaps.begin();
3161         itr != _heightMaps.end(); ++itr) {
3162        itr->second->setClippingPlanes(_activeClipPlanes);
3163    }
3164    for (LICHashmap::iterator itr = _lics.begin();
3165         itr != _lics.end(); ++itr) {
3166        itr->second->setClippingPlanes(_activeClipPlanes);
3167    }
3168    for (MoleculeHashmap::iterator itr = _molecules.begin();
3169         itr != _molecules.end(); ++itr) {
3170        itr->second->setClippingPlanes(_activeClipPlanes);
3171    }
3172    for (PolyDataHashmap::iterator itr = _polyDatas.begin();
3173         itr != _polyDatas.end(); ++itr) {
3174        itr->second->setClippingPlanes(_activeClipPlanes);
3175    }
3176    for (PseudoColorHashmap::iterator itr = _pseudoColors.begin();
3177         itr != _pseudoColors.end(); ++itr) {
3178        itr->second->setClippingPlanes(_activeClipPlanes);
3179    }
3180    for (StreamlinesHashmap::iterator itr = _streamlines.begin();
3181         itr != _streamlines.end(); ++itr) {
3182        itr->second->setClippingPlanes(_activeClipPlanes);
3183    }
3184    for (VolumeHashmap::iterator itr = _volumes.begin();
3185         itr != _volumes.end(); ++itr) {
3186        itr->second->setClippingPlanes(_activeClipPlanes);
3187    }
3188}
3189
3190/**
3191 * \brief Control the use of two sided lighting
3192 */
3193void Renderer::setUseTwoSidedLighting(bool state)
3194{
3195    _renderer->SetTwoSidedLighting(state ? 1 : 0);
3196    _needsRedraw = true;
3197}
3198
3199/**
3200 * \brief Control the use of the depth peeling algorithm for transparency
3201 */
3202void Renderer::setUseDepthPeeling(bool state)
3203{
3204    _renderer->SetUseDepthPeeling(state ? 1 : 0);
3205    _needsRedraw = true;
3206}
3207
3208/**
3209 * \brief Sets flag to trigger rendering next time render() is called
3210 */
3211void Renderer::eventuallyRender()
3212{
3213    _needsRedraw = true;
3214}
3215
3216/**
3217 * \brief Cause the rendering to render a new image if needed
3218 *
3219 * The _needsRedraw flag indicates if a state change has occured since
3220 * the last rendered frame
3221 */
3222bool Renderer::render()
3223{
3224    if (_needsRedraw) {
3225        setCameraClippingPlanes();
3226        _renderWindow->Render();
3227        _needsRedraw = false;
3228        return true;
3229    } else
3230        return false;
3231}
3232
3233/// Get the pixel width of the render window/image
3234int Renderer::getWindowWidth() const
3235{
3236    return _windowWidth;
3237}
3238
3239/// Get the pixel height of the render window/image
3240int Renderer::getWindowHeight() const
3241{
3242    return _windowHeight;
3243}
3244
3245/**
3246 * \brief Read back the rendered framebuffer image
3247 */
3248void Renderer::getRenderedFrame(vtkUnsignedCharArray *imgData)
3249{
3250#ifdef RENDER_TARGA
3251    _renderWindow->MakeCurrent();
3252    // Must clear previous errors first.
3253    while (glGetError() != GL_NO_ERROR){
3254        ;
3255    }
3256    int bytesPerPixel = TARGA_BYTES_PER_PIXEL;
3257    int size = bytesPerPixel * _windowWidth * _windowHeight;
3258
3259    if (imgData->GetMaxId() + 1 != size)
3260    {
3261        imgData->SetNumberOfComponents(bytesPerPixel);
3262        imgData->SetNumberOfValues(size);
3263    }
3264    glDisable(GL_TEXTURE_2D);
3265    if (_renderWindow->GetDoubleBuffer()) {
3266        glReadBuffer(static_cast<GLenum>(vtkOpenGLRenderWindow::SafeDownCast(_renderWindow)->GetBackLeftBuffer()));
3267    } else {
3268        glReadBuffer(static_cast<GLenum>(vtkOpenGLRenderWindow::SafeDownCast(_renderWindow)->GetFrontLeftBuffer()));
3269    }
3270    glPixelStorei(GL_PACK_ALIGNMENT, 1);
3271#ifdef WANT_TRACE
3272    struct timeval t1, t2;
3273    glFinish();
3274    gettimeofday(&t1, 0);
3275#endif
3276    if (bytesPerPixel == 4) {
3277        glReadPixels(0, 0, _windowWidth, _windowHeight, GL_BGRA,
3278                     GL_UNSIGNED_BYTE, imgData->GetPointer(0));
3279    } else {
3280        glReadPixels(0, 0, _windowWidth, _windowHeight, GL_BGR,
3281                     GL_UNSIGNED_BYTE, imgData->GetPointer(0));
3282    }
3283#ifdef WANT_TRACE
3284    gettimeofday(&t2, 0);
3285    static unsigned int numFrames = 0;
3286    static double accum = 0;
3287    numFrames++;
3288    accum += ELAPSED_TIME(t1, t2);
3289#endif
3290    TRACE("Readback time: %g ms", ELAPSED_TIME(t1, t2));
3291    TRACE("Readback avg: %g ms", accum/numFrames);
3292    if (glGetError() != GL_NO_ERROR) {
3293        ERROR("glReadPixels");
3294    }
3295#else
3296    _renderWindow->GetPixelData(0, 0, _windowWidth-1, _windowHeight-1,
3297                                !_renderWindow->GetDoubleBuffer(), imgData);
3298#endif
3299    TRACE("Image data size: %d", imgData->GetSize());
3300}
3301
3302/**
3303 * \brief Get nearest data value given display coordinates x,y
3304 *
3305 * Note: no interpolation is performed on data
3306 */
3307bool Renderer::getScalarValueAtPixel(const DataSetId& id, int x, int y, double *value)
3308{
3309    vtkSmartPointer<vtkCoordinate> coord = vtkSmartPointer<vtkCoordinate>::New();
3310    coord->SetCoordinateSystemToDisplay();
3311    coord->SetValue(x, _windowHeight - y, 0);
3312    double *worldCoords = coord->GetComputedWorldValue(_renderer);
3313
3314    TRACE("Pixel coords: %d, %d\nWorld coords: %g, %g, %g", x, y,
3315          worldCoords[0],
3316          worldCoords[1],
3317          worldCoords[2]);
3318
3319    return getScalarValue(id, worldCoords[0], worldCoords[1], worldCoords[2], value);
3320}
3321
3322/**
3323 * \brief Get nearest data value given world coordinates x,y,z
3324 *
3325 * Note: no interpolation is performed on data
3326 */
3327bool Renderer::getScalarValue(const DataSetId& id, double x, double y, double z, double *value)
3328{
3329    DataSet *ds = getDataSet(id);
3330    if (ds == NULL)
3331        return false;
3332
3333    return ds->getScalarValue(x, y, z, value);
3334}
3335
3336/**
3337 * \brief Get nearest data value given display coordinates x,y
3338 *
3339 * Note: no interpolation is performed on data
3340 */
3341bool Renderer::getVectorValueAtPixel(const DataSetId& id, int x, int y, double vector[3])
3342{
3343    vtkSmartPointer<vtkCoordinate> coord = vtkSmartPointer<vtkCoordinate>::New();
3344    coord->SetCoordinateSystemToDisplay();
3345    coord->SetValue(x, _windowHeight - y, 0);
3346    double *worldCoords = coord->GetComputedWorldValue(_renderer);
3347
3348    TRACE("Pixel coords: %d, %d\nWorld coords: %g, %g, %g", x, y,
3349          worldCoords[0],
3350          worldCoords[1],
3351          worldCoords[2]);
3352
3353    return getVectorValue(id, worldCoords[0], worldCoords[1], worldCoords[2], vector);
3354}
3355
3356/**
3357 * \brief Get nearest data value given world coordinates x,y,z
3358 *
3359 * Note: no interpolation is performed on data
3360 */
3361bool Renderer::getVectorValue(const DataSetId& id, double x, double y, double z, double vector[3])
3362{
3363    DataSet *ds = getDataSet(id);
3364    if (ds == NULL)
3365        return false;
3366
3367    return ds->getVectorValue(x, y, z, vector);
3368}
Note: See TracBrowser for help on using the repository browser.