FitToHeightMap

Description¶

This example uses vtkFitToHeightMapFilter to cover a height map with vtkPolyData generated by a vtkPlaneSource. The plane's origin, point1 and point2 are calculated from the DEM's bounding box. The z coordinate of each parameter is set to the height of the DEM. With this plane, the example uses vtkProbeFilter to apply the scalar elevation data from the DEM.

vtkWarpScalar is used to "elevate" the DEM mesh. Once the planes are draped over the warped DEM mesh, they are displayed side-by-side using the same vtkLookupTable.

The displayed meshes are original (left), point fit (middle) and cell fit (right).

Note

DEM files for the United States are available here.

Info

This example requires vtk version 8.2 or newer.

Other languages

See (PythonicAPI)

Question

If you have a question about this example, please use the VTK Discourse Forum

Code¶

FitToHeightMap.cxx

#include <vtkActor.h>
#include <vtkCamera.h>
#include <vtkDEMReader.h>
#include <vtkFitToHeightMapFilter.h>
#include <vtkImageData.h>
#include <vtkImageDataGeometryFilter.h>
#include <vtkLookupTable.h>
#include <vtkNamedColors.h>
#include <vtkNew.h>
#include <vtkPlaneSource.h>
#include <vtkPolyDataMapper.h>
#include <vtkProbeFilter.h>
#include <vtkProperty.h>
#include <vtkRenderWindow.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkRenderer.h>
#include <vtkWarpScalar.h>

int main(int argc, char* argv[])
{
  if (argc < 2)
  {
    std::cout << "Usage: " << argv[0] << " file.dem e.g. SainteHelens.dem"
              << std::endl;
    return EXIT_FAILURE;
  }
  // Create the RenderWindow, Renderer
  //
  vtkNew<vtkNamedColors> colors;

  vtkNew<vtkRenderer> ren0;
  ren0->SetViewport(0, 0, 1.0 / 3.0, 1);
  ren0->SetBackground(colors->GetColor3d("Wheat").GetData());

  vtkNew<vtkRenderer> ren1;
  ren1->SetViewport(1.0 / 3.0, 0, 2.0 / 3.0, 1);
  ren1->SetBackground(colors->GetColor3d("BurlyWood").GetData());

  vtkNew<vtkRenderer> ren2;
  ren2->SetViewport(2.0 / 3.0, 0, 1, 1);
  ren2->SetBackground(colors->GetColor3d("Tan").GetData());

  vtkNew<vtkRenderWindow> renderWindow;
  renderWindow->SetSize(1200, 400);
  renderWindow->AddRenderer(ren0);
  renderWindow->AddRenderer(ren1);
  renderWindow->AddRenderer(ren2);
  renderWindow->SetWindowName("FitToHeightMap");

  vtkNew<vtkRenderWindowInteractor> interactor;
  interactor->SetRenderWindow(renderWindow);

  // Create pipeline. Load terrain data.
  //
  vtkNew<vtkLookupTable> lut;
  lut->SetHueRange(0.6, 0);
  lut->SetSaturationRange(1.0, 0);
  lut->SetValueRange(0.5, 1.0);

  // Read the data: a height field results
  vtkNew<vtkDEMReader> demReader;
  demReader->SetFileName(argv[1]);
  demReader->Update(); // Update so we can get bounds

  vtkNew<vtkImageDataGeometryFilter> surface;
  surface->SetInputConnection(demReader->GetOutputPort());

  // Warp the surface in the vertical direction
  vtkNew<vtkWarpScalar> warp;
  warp->SetInputConnection(surface->GetOutputPort());
  warp->SetScaleFactor(1);
  warp->UseNormalOn();
  warp->SetNormal(0, 0, 1);
  warp->Update();

  // Show the terrain
  double lo = demReader->GetOutput()->GetScalarRange()[0];
  double hi = demReader->GetOutput()->GetScalarRange()[1];

  vtkNew<vtkPolyDataMapper> demMapper;
  demMapper->SetInputConnection(warp->GetOutputPort());
  demMapper->ScalarVisibilityOn();
  demMapper->SetScalarRange(lo, hi);
  demMapper->SetLookupTable(lut);

  vtkNew<vtkActor> demActor;
  demActor->SetMapper(demMapper);

  // Create polygon(s) to fit. z-values should be height of DEM
  // to accomodate ProbeFilter
  double zLevel = demReader->GetOutput()->GetBounds()[5];
  double* demBounds;
  demBounds = demReader->GetOutput()->GetBounds();

  vtkNew<vtkPlaneSource> plane;
  plane->SetOrigin(demBounds[0], demBounds[2], zLevel);
  plane->SetPoint1(demBounds[1], demBounds[2], zLevel);
  plane->SetPoint2(demBounds[0], demBounds[3], zLevel);
  plane->SetResolution(128, 128);
  plane->Update();

  // Get the scalrs from the DEM
  vtkNew<vtkProbeFilter> probeDem;
  probeDem->SetSourceData(demReader->GetOutput());
  probeDem->SetInputConnection(plane->GetOutputPort());
  probeDem->Update();

  // Fit polygons to surface using two strategies
  // Both will share the same lookup tables
  vtkNew<vtkFitToHeightMapFilter> pointFit;
  pointFit->SetInputConnection(probeDem->GetOutputPort());
  pointFit->SetHeightMapConnection(demReader->GetOutputPort());
  pointFit->SetFittingStrategyToPointProjection();
  pointFit->UseHeightMapOffsetOn();

  vtkNew<vtkPolyDataMapper> pointMapper;
  pointMapper->SetInputConnection(pointFit->GetOutputPort());
  pointMapper->ScalarVisibilityOn();
  pointMapper->SetScalarRange(lo, hi);
  pointMapper->SetLookupTable(lut);

  vtkNew<vtkActor> pointActor;
  pointActor->SetMapper(pointMapper);

  // Fit polygons to surface (cell strategy)
  vtkNew<vtkFitToHeightMapFilter> cellFit;
  cellFit->SetInputConnection(probeDem->GetOutputPort());
  cellFit->SetHeightMapConnection(demReader->GetOutputPort());
  cellFit->SetFittingStrategyToCellAverageHeight();
  cellFit->UseHeightMapOffsetOn();

  vtkNew<vtkPolyDataMapper> cellMapper;
  cellMapper->SetInputConnection(cellFit->GetOutputPort());
  cellMapper->ScalarVisibilityOn();
  cellMapper->SetScalarRange(lo, hi);
  cellMapper->SetLookupTable(lut);

  vtkNew<vtkActor> cellActor;
  cellActor->SetMapper(cellMapper);

  // Render the three representations
  ren0->AddActor(demActor);
  ren1->AddActor(pointActor);
  ren2->AddActor(cellActor);

  // Look down the x axis
  ren0->GetActiveCamera()->SetPosition(1, 0, 0);
  ren0->GetActiveCamera()->SetFocalPoint(0, 1, 0);
  ren0->GetActiveCamera()->SetViewUp(0, 0, 1);
  ren0->ResetCamera();

  // Rotate to an oblique view
  ren0->GetActiveCamera()->Azimuth(30.0);
  ren0->GetActiveCamera()->Elevation(60.0);

  // Share the cameras
  ren1->SetActiveCamera(ren0->GetActiveCamera());
  ren2->SetActiveCamera(ren0->GetActiveCamera());

  renderWindow->Render();
  interactor->Start();

  return EXIT_SUCCESS;
}

CMakeLists.txt¶

cmake_minimum_required(VERSION 3.12 FATAL_ERROR)

project(FitToHeightMap)

find_package(VTK COMPONENTS 
  CommonColor
  CommonCore
  CommonDataModel
  FiltersCore
  FiltersGeneral
  FiltersGeometry
  FiltersModeling
  FiltersSources
  IOImage
  InteractionStyle
  RenderingContextOpenGL2
  RenderingCore
  RenderingFreeType
  RenderingGL2PSOpenGL2
  RenderingOpenGL2
)

if (NOT VTK_FOUND)
  message(FATAL_ERROR "FitToHeightMap: Unable to find the VTK build folder.")
endif()

# Prevent a "command line is too long" failure in Windows.
set(CMAKE_NINJA_FORCE_RESPONSE_FILE "ON" CACHE BOOL "Force Ninja to use response files.")
add_executable(FitToHeightMap MACOSX_BUNDLE FitToHeightMap.cxx )
  target_link_libraries(FitToHeightMap PRIVATE ${VTK_LIBRARIES}
)
# vtk_module_autoinit is needed
vtk_module_autoinit(
  TARGETS FitToHeightMap
  MODULES ${VTK_LIBRARIES}
)

Download and Build FitToHeightMap¶

Click here to download FitToHeightMap and its CMakeLists.txt file. Once the tarball FitToHeightMap.tar has been downloaded and extracted,

cd FitToHeightMap/build

If VTK is installed:

cmake ..

If VTK is not installed but compiled on your system, you will need to specify the path to your VTK build:

cmake -DVTK_DIR:PATH=/home/me/vtk_build ..

Build the project:

make

and run it:

./FitToHeightMap

WINDOWS USERS

Be sure to add the VTK bin directory to your path. This will resolve the VTK dll's at run time.