diff --git a/Applications/DataExplorer/VtkVis/CMakeLists.txt b/Applications/DataExplorer/VtkVis/CMakeLists.txt index b746c07548ba55ff12501db005deba09ab9cab41..ef1894a69343f8762e53030e4fe1be8aa770c7eb 100644 --- a/Applications/DataExplorer/VtkVis/CMakeLists.txt +++ b/Applications/DataExplorer/VtkVis/CMakeLists.txt @@ -19,15 +19,21 @@ set(SOURCES VtkCompositeElementSelectionFilter.cpp VtkCompositeGeoObjectFilter.cpp VtkCompositeImageToCylindersFilter.cpp + VtkCompositeImageToPointCloudFilter.cpp + VtkCompositeImageToSurfacePointsFilter.cpp VtkCompositeLineToTubeFilter.cpp VtkCompositeNodeSelectionFilter.cpp VtkCompositePointToGlyphFilter.cpp VtkCompositeTextureOnSurfaceFilter.cpp VtkCompositeThresholdFilter.cpp VtkConsoleOutputWindow.cpp + VtkCustomInteractorStyle.cpp VtkFilterFactory.cpp VtkGeoImageSource.cpp VtkImageDataToLinePolyDataFilter.cpp + VtkImageDataToPointCloudFilter.cpp + VtkImageDataToSurfacePointsFilter.cpp + VtkPickCallback.cpp VtkPolylinesSource.cpp VtkPointsSource.cpp VtkRaster.cpp @@ -41,11 +47,18 @@ set(SOURCES VtkVisPipelineView.cpp VtkVisPointSetItem.cpp VtkVisTabWidget.cpp - VtkPickCallback.cpp - VtkCustomInteractorStyle.cpp ) set(HEADERS + MeshFromRasterDialog.h + QVtkDataSetMapper.h + VisPrefsDialog.h + VisualizationWidget.h + VtkAddFilterDialog.h + VtkAlgorithmProperties.h + VtkAlgorithmPropertyLineEdit.h + VtkAlgorithmPropertyCheckbox.h + VtkAlgorithmPropertyVectorEdit.h VtkAppendArrayFilter.h VtkBGImageSource.h VtkColorByHeightFilter.h @@ -57,14 +70,21 @@ set(HEADERS VtkCompositeElementSelectionFilter.h VtkCompositeGeoObjectFilter.h VtkCompositeImageToCylindersFilter.h + VtkCompositeImageToPointCloudFilter.h + VtkCompositeImageToSurfacePointsFilter.h VtkCompositeLineToTubeFilter.h VtkCompositeNodeSelectionFilter.h VtkCompositePointToGlyphFilter.h VtkCompositeTextureOnSurfaceFilter.h VtkCompositeThresholdFilter.h + VtkConsoleOutputWindow.h + VtkCustomInteractorStyle.h VtkFilterFactory.h VtkGeoImageSource.h VtkImageDataToLinePolyDataFilter.h + VtkImageDataToPointCloudFilter.h + VtkImageDataToSurfacePointsFilter.h + VtkPickCallback.h VtkPolylinesSource.h VtkPointsSource.h VtkRaster.h @@ -73,23 +93,11 @@ set(HEADERS VtkTextureOnSurfaceFilter.h VtkVisHelper.h VtkVisImageItem.h - VtkVisPipelineItem.h - VtkVisPointSetItem.h - MeshFromRasterDialog.h - QVtkDataSetMapper.h - VisPrefsDialog.h - VisualizationWidget.h - VtkAddFilterDialog.h - VtkAlgorithmProperties.h - VtkAlgorithmPropertyLineEdit.h - VtkAlgorithmPropertyCheckbox.h - VtkAlgorithmPropertyVectorEdit.h VtkVisPipeline.h + VtkVisPipelineItem.h VtkVisPipelineView.h + VtkVisPointSetItem.h VtkVisTabWidget.h - VtkConsoleOutputWindow.h - VtkPickCallback.h - VtkCustomInteractorStyle.h ) # Visual Studio folder diff --git a/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.cpp b/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.cpp new file mode 100644 index 0000000000000000000000000000000000000000..2c7bb9c8706b668f37ef4252306d249a910edd64 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.cpp @@ -0,0 +1,101 @@ +/** + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#include "VtkCompositeImageToPointCloudFilter.h" + +#include <vtkPointData.h> +#include <vtkSmartPointer.h> + +#include <QMap> +#include <QString> +#include <QVariant> + +#include "VtkImageDataToPointCloudFilter.h" + +VtkCompositeImageToPointCloudFilter::VtkCompositeImageToPointCloudFilter(vtkAlgorithm* inputAlgorithm) + : VtkCompositeFilter(inputAlgorithm) +{ + init(); +} + +void VtkCompositeImageToPointCloudFilter::init() +{ + _inputDataObjectType = VTK_IMAGE_DATA; + _outputDataObjectType = VTK_POLY_DATA; + + VtkImageDataToPointCloudFilter* point_cloud_filter = + VtkImageDataToPointCloudFilter::New(); + point_cloud_filter->SetInputConnection(_inputAlgorithm->GetOutputPort()); + _inputAlgorithm->Update(); + + QList<QVariant> n_points_range_list; + n_points_range_list.push_back(point_cloud_filter->GetMinNumberOfPointsPerCell()); + n_points_range_list.push_back(point_cloud_filter->GetMaxNumberOfPointsPerCell()); + (*_algorithmUserVectorProperties)["Number of points range"] = n_points_range_list; + QList<QVariant> vertical_extent_list; + vertical_extent_list.push_back(point_cloud_filter->GetMinHeight()); + vertical_extent_list.push_back(point_cloud_filter->GetMaxHeight()); + (*_algorithmUserVectorProperties)["Vertical extent"] = vertical_extent_list; + (*_algorithmUserProperties)["Logarithmic interpolation"] = !point_cloud_filter->GetIsLinear(); + (*_algorithmUserProperties)["Gamma value"] = point_cloud_filter->GetGamma(); + + point_cloud_filter->Update(); + _outputAlgorithm = point_cloud_filter; +} + +void VtkCompositeImageToPointCloudFilter::SetUserProperty(QString name, QVariant value) +{ + VtkAlgorithmProperties::SetUserProperty(name, value); + + if ((name == "Gamma value") && (value.toDouble() > 0)) + { + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm)->SetGamma(value.toDouble()); + } + if (name == "Logarithmic interpolation") + { + if (value.toBool() == true) + { + double const gamma = + VtkAlgorithmProperties::GetUserProperty("Gamma value").toDouble(); + if (gamma > 0) + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->useLogarithmicInterpolation(gamma); + } + else + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->useLinearInterpolation(); + } +} + +void VtkCompositeImageToPointCloudFilter::SetUserVectorProperty(QString name, QList<QVariant> values) +{ + VtkAlgorithmProperties::SetUserVectorProperty(name, values); + + if (name == "Number of points range") + { + if (values[0].toInt() >= 0 && values[1].toInt() >= 0 && + values[0].toInt() <= values[1].toInt()) + { + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->SetMinNumberOfPointsPerCell(values[0].toInt()); + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->SetMaxNumberOfPointsPerCell(values[1].toInt()); + } + } + else if (name == "Vertical extent") + { + if (values[0].toDouble() <= values[1].toDouble()) + { + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->SetMinHeight(values[0].toDouble()); + static_cast<VtkImageDataToPointCloudFilter*>(_outputAlgorithm) + ->SetMaxHeight(values[1].toDouble()); + } + } +} diff --git a/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.h b/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.h new file mode 100644 index 0000000000000000000000000000000000000000..d0a78a63bb1b26506fb42633906f847f2ee3373e --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkCompositeImageToPointCloudFilter.h @@ -0,0 +1,27 @@ +/** + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#pragma once + +#include "VtkCompositeFilter.h" + +class VtkImageDataToPointCloudFilter; + +class VtkCompositeImageToPointCloudFilter : public VtkCompositeFilter +{ +public: + VtkCompositeImageToPointCloudFilter(vtkAlgorithm* inputAlgorithm); + ~VtkCompositeImageToPointCloudFilter() = default; + + void init() override; + + void SetUserProperty(QString name, QVariant value) override; + + void SetUserVectorProperty(QString name, QList<QVariant> values) override; +}; diff --git a/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.cpp b/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.cpp new file mode 100644 index 0000000000000000000000000000000000000000..1a5c0cf15f0f93364a05d0ddcd8cc24f22653725 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.cpp @@ -0,0 +1,57 @@ +/** + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#include "VtkCompositeImageToSurfacePointsFilter.h" + +#include <vtkPointData.h> +#include <vtkSmartPointer.h> + +#include <QMap> +#include <QString> +#include <QVariant> + +#include "VtkImageDataToSurfacePointsFilter.h" + +VtkCompositeImageToSurfacePointsFilter::VtkCompositeImageToSurfacePointsFilter( + vtkAlgorithm* inputAlgorithm) + : VtkCompositeFilter(inputAlgorithm) +{ + init(); +} + +void VtkCompositeImageToSurfacePointsFilter::init() +{ + _inputDataObjectType = VTK_IMAGE_DATA; + _outputDataObjectType = VTK_POLY_DATA; + + VtkImageDataToSurfacePointsFilter* point_cloud_filter = + VtkImageDataToSurfacePointsFilter::New(); + point_cloud_filter->SetInputConnection(_inputAlgorithm->GetOutputPort()); + _inputAlgorithm->Update(); + + (*_algorithmUserProperties)["Points per pixel"] = + static_cast<int>(point_cloud_filter->GetPointsPerPixel()); + point_cloud_filter->Update(); + _outputAlgorithm = point_cloud_filter; +} + +void VtkCompositeImageToSurfacePointsFilter::SetUserProperty(QString name, QVariant value) +{ + VtkAlgorithmProperties::SetUserProperty(name, value); + if ((name == "Points per pixel") && (value.toInt() > 0)) + static_cast<VtkImageDataToSurfacePointsFilter*>(_outputAlgorithm) + ->SetPointsPerPixel(static_cast<vtkIdType>(value.toInt())); +} + +void VtkCompositeImageToSurfacePointsFilter::SetUserVectorProperty( + QString name, QList<QVariant> values) +{ + Q_UNUSED(name); + Q_UNUSED(values); +} diff --git a/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.h b/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.h new file mode 100644 index 0000000000000000000000000000000000000000..feed48bf5948dffeb9482e8345939b44f9927a74 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkCompositeImageToSurfacePointsFilter.h @@ -0,0 +1,27 @@ +/** + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#pragma once + +#include "VtkCompositeFilter.h" + +class VtkImageDataToSurfacePointsFilter; + +class VtkCompositeImageToSurfacePointsFilter : public VtkCompositeFilter +{ +public: + VtkCompositeImageToSurfacePointsFilter(vtkAlgorithm* inputAlgorithm); + ~VtkCompositeImageToSurfacePointsFilter() = default; + + void init() override; + + void SetUserProperty(QString name, QVariant value) override; + + void SetUserVectorProperty(QString name, QList<QVariant> values) override; +}; diff --git a/Applications/DataExplorer/VtkVis/VtkFilterFactory.cpp b/Applications/DataExplorer/VtkVis/VtkFilterFactory.cpp index f915ea178e7f262f9a9948b7e8262022302185f8..580f5b9e7c587de89142a92bafa481e09986eea9 100644 --- a/Applications/DataExplorer/VtkVis/VtkFilterFactory.cpp +++ b/Applications/DataExplorer/VtkVis/VtkFilterFactory.cpp @@ -21,6 +21,8 @@ #include "VtkCompositeElementSelectionFilter.h" #include "VtkCompositeGeoObjectFilter.h" #include "VtkCompositeImageToCylindersFilter.h" +#include "VtkCompositeImageToPointCloudFilter.h" +#include "VtkCompositeImageToSurfacePointsFilter.h" #include "VtkCompositeLineToTubeFilter.h" #include "VtkCompositeNodeSelectionFilter.h" #include "VtkCompositePointToGlyphFilter.h" @@ -83,6 +85,18 @@ const QVector<VtkFilterInfo> VtkFilterFactory::GetFilterList() "Visualisation of contour-lines/-planes within dense scalar fields.", VTK_UNSTRUCTURED_GRID, VTK_UNSTRUCTURED_GRID)); + filterList.push_back(VtkFilterInfo( + "VtkCompositeImageToPointCloudFilter", "Image to point cloud", + "This filter creates a point cloud with a density based on the first " + "component of pixel values.", + VTK_IMAGE_DATA, VTK_POLY_DATA)); + + filterList.push_back(VtkFilterInfo( + "VtkCompositeImageToSurfacePointsFilter", "Image to surface points", + "This filter creates a point cloud representing a surface based on the " + "first component of pixel values.", + VTK_IMAGE_DATA, VTK_POLY_DATA)); + // Simple filters filterList.push_back(VtkFilterInfo( "VtkImageDataToLinePolyDataFilter", @@ -131,6 +145,10 @@ VtkCompositeFilter* VtkFilterFactory::CreateCompositeFilter( QString type, return new VtkCompositeContourFilter(inputAlgorithm); if (type.compare(QString("VtkCompositeGeoObjectFilter")) == 0) return new VtkCompositeGeoObjectFilter(inputAlgorithm); + if (type.compare(QString("VtkCompositeImageToPointCloudFilter")) == 0) + return new VtkCompositeImageToPointCloudFilter(inputAlgorithm); + if (type.compare(QString("VtkCompositeImageToSurfacePointsFilter")) == 0) + return new VtkCompositeImageToSurfacePointsFilter(inputAlgorithm); return nullptr; } diff --git a/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.cpp b/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.cpp new file mode 100644 index 0000000000000000000000000000000000000000..5b9d0a643f9a668126eff5040498c7c253228d73 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.cpp @@ -0,0 +1,173 @@ +/** + * \file + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +// ** INCLUDES ** +#include "VtkImageDataToPointCloudFilter.h" + +#include <algorithm> +#include <vector> + +#include <vtkIdList.h> +#include <vtkImageData.h> +#include <vtkInformation.h> +#include <vtkInformationVector.h> +#include <vtkLine.h> +#include <vtkObjectFactory.h> +#include <vtkPointData.h> +#include <vtkPolyData.h> +#include <vtkSmartPointer.h> + +vtkStandardNewMacro(VtkImageDataToPointCloudFilter); + +VtkImageDataToPointCloudFilter::VtkImageDataToPointCloudFilter() + : Gamma(1.0), + PointScaleFactor(1.0), + MinHeight(0), + MaxHeight(1000), + MinNumberOfPointsPerCell(1), + MaxNumberOfPointsPerCell(20), + IsLinear(true) +{ +} + +void VtkImageDataToPointCloudFilter::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); +} + +int VtkImageDataToPointCloudFilter::FillInputPortInformation(int, vtkInformation* info) +{ + info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkImageData"); + return 1; +} + +int VtkImageDataToPointCloudFilter::RequestData( + vtkInformation*, + vtkInformationVector** inputVector, + vtkInformationVector* outputVector) +{ + vtkDebugMacro(<< "Executing VtkImageDataToPointCloudFilter"); + + vtkInformation* input_info = inputVector[0]->GetInformationObject(0); + vtkInformation* output_info = outputVector->GetInformationObject(0); + vtkImageData* input = vtkImageData::SafeDownCast( + input_info->Get(vtkDataObject::DATA_OBJECT())); + vtkPolyData* output = vtkPolyData::SafeDownCast( + output_info->Get(vtkDataObject::DATA_OBJECT())); + + void* pixvals = input->GetScalarPointer(); + int const n_comp = input->GetNumberOfScalarComponents(); + if (n_comp < 1) + vtkDebugMacro("Error reading number of components."); + if (n_comp > 2) + vtkDebugMacro("RGB colours detected. Using only first channel for computation."); + + vtkIdType const n_points = input->GetNumberOfPoints(); + if (n_points == 0) + { + vtkDebugMacro("No data found!"); + return -1; + } + + double spacing[3]; + input->GetSpacing(spacing); + double range[2]; + vtkPointData* input_data = input->GetPointData(); + input_data->GetArray(0)->GetRange(range); + + std::vector<std::size_t> density; + density.reserve(n_points); + for (vtkIdType i = 0; i < n_points; ++i) + { + // Skip transparent pixels + if (n_comp == 2 || n_comp == 4) + { + if (((float*)pixvals)[(i + 1) * n_comp - 1] < 0.00000001f) + { + density.push_back(0); + continue; + } + } + float const val(((float*)pixvals)[i * n_comp]); + double const calc_gamma = (IsLinear) ? 1 : Gamma; + double const pnts_per_cell = interpolate(range[0], range[1], val, calc_gamma); + density.push_back(static_cast<std::size_t>( + std::floor(pnts_per_cell * GetPointScaleFactor()))); + } + + // Allocate the space needed for output objects + std::size_t const sum = std::accumulate(density.begin(), density.end(), 0); + vtkSmartPointer<vtkPoints> new_points = vtkSmartPointer<vtkPoints>::New(); + new_points->SetNumberOfPoints(sum); + vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New(); + cells->Allocate(sum); + + vtkPointData* output_data = output->GetPointData(); + double const half_cellsize(spacing[0] / 2.0); + std::size_t pnt_idx(0); + for (std::size_t i = 0; i < static_cast<std::size_t>(n_points); ++i) + { + if (density[i] == 0) + continue; + double p[3]; + input->GetPoint(i, p); + + MathLib::Point3d min_pnt{std::array<double, 3>{ + {p[0] - half_cellsize, p[1] - half_cellsize, MinHeight}}}; + MathLib::Point3d max_pnt{std::array<double, 3>{ + {p[0] + half_cellsize, p[1] + half_cellsize, MaxHeight}}}; + createPoints(new_points, cells, pnt_idx, density[i], min_pnt, max_pnt); + pnt_idx += density[i]; + } + + output->SetPoints(new_points); + output->SetVerts(cells); + output->Squeeze(); + + vtkDebugMacro(<< "Created " << new_points->GetNumberOfPoints() << " points."); + return 1; +} + +void VtkImageDataToPointCloudFilter::createPoints( + vtkSmartPointer<vtkPoints>& points, + vtkSmartPointer<vtkCellArray>& cells, + std::size_t pnt_idx, + std::size_t n_points, + MathLib::Point3d const& min_pnt, + MathLib::Point3d const& max_pnt) const +{ + for (std::size_t i = 0; i < n_points; ++i) + { + double p[3]; + p[0] = getRandomNumber(min_pnt[0], max_pnt[0]); + p[1] = getRandomNumber(min_pnt[1], max_pnt[1]); + p[2] = getRandomNumber(min_pnt[2], max_pnt[2]); + points->SetPoint(pnt_idx + i, p); + cells->InsertNextCell(1); + cells->InsertCellPoint(pnt_idx + i); + } +} + +double VtkImageDataToPointCloudFilter::getRandomNumber(double const& min, double const& max) const +{ + return ((double)std::rand() / RAND_MAX) * (max - min) + min; +} + +std::size_t VtkImageDataToPointCloudFilter::interpolate(double low, + double high, + double p, + double gamma) const +{ + assert(p >= low && p <= high); + double const r = (p - low) / (high - low); + return static_cast<std::size_t>( + (MaxNumberOfPointsPerCell - MinNumberOfPointsPerCell) * std::pow(r, gamma) + + MinNumberOfPointsPerCell); +} diff --git a/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.h b/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.h new file mode 100644 index 0000000000000000000000000000000000000000..f008fef21e6cbc8d568dff1a051a2d0db22f8e34 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkImageDataToPointCloudFilter.h @@ -0,0 +1,97 @@ +/** + * \file + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#pragma once + +#include <vtkPolyDataAlgorithm.h> + +#include "MathLib/Point3d.h" +#include "VtkAlgorithmProperties.h" + +class vtkInformation; +class vtkInformationVector; + +/// A VTK Filter that will create a point cloud with local densities based on +/// pixel values. +class VtkImageDataToPointCloudFilter : public vtkPolyDataAlgorithm +{ +public: + /// Create a new objects (required because of VTKs reference counting) + static VtkImageDataToPointCloudFilter* New(); + + vtkTypeMacro(VtkImageDataToPointCloudFilter, vtkPolyDataAlgorithm); + + /// Prints information about the filter + void PrintSelf(ostream& os, vtkIndent indent) override; + + void useLinearInterpolation() { SetIsLinear(true); } + void useLogarithmicInterpolation(double gamma) + { + SetIsLinear(false); + SetGamma(gamma); + } + + vtkGetMacro(IsLinear, bool); + vtkSetMacro(IsLinear, bool); + + vtkGetMacro(MinNumberOfPointsPerCell, vtkIdType); + vtkSetMacro(MinNumberOfPointsPerCell, vtkIdType); + + vtkGetMacro(MaxNumberOfPointsPerCell, vtkIdType); + vtkSetMacro(MaxNumberOfPointsPerCell, vtkIdType); + + vtkGetMacro(Gamma, double); + vtkSetMacro(Gamma, double); + + vtkGetMacro(PointScaleFactor, double); + vtkSetMacro(PointScaleFactor, double); + + vtkGetMacro(MinHeight, double); + vtkSetMacro(MinHeight, double); + + vtkGetMacro(MaxHeight, double); + vtkSetMacro(MaxHeight, double); + +protected: + VtkImageDataToPointCloudFilter(); + ~VtkImageDataToPointCloudFilter() = default; + + /// Sets input port to vtkImageData. + int FillInputPortInformation(int port, vtkInformation* info) override; + + /// Updates the graphical object + int RequestData(vtkInformation* request, vtkInformationVector** inputVector, + vtkInformationVector* outputVector) override; + + double Gamma; + double PointScaleFactor; + double MinHeight; + double MaxHeight; + vtkIdType MinNumberOfPointsPerCell; + vtkIdType MaxNumberOfPointsPerCell; + bool IsLinear; + + +private: + /// Creates the point objects based on the parameters set by the user + void createPoints(vtkSmartPointer<vtkPoints>& points, + vtkSmartPointer<vtkCellArray>& cells, std::size_t pnt_idx, + std::size_t n_points, MathLib::Point3d const& min_pnt, + MathLib::Point3d const& max_pnt) const; + + /// Returns a random number in [min, max] + double getRandomNumber(double const& min, double const& max) const; + + /// Interpolates the required number of points. Using gamma = 1 this is a linear + /// interpolation, for values smaller/greater than 1 the curve becomes + /// logarithmic/exponential, resulting in a smaller/larger difference of point + /// densities given the same input values. + std::size_t interpolate(double low, double high, double p, double gamma) const; +}; diff --git a/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.cpp b/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.cpp new file mode 100644 index 0000000000000000000000000000000000000000..e5324ab3b0614f11976040db0cfd52af29470773 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.cpp @@ -0,0 +1,156 @@ +/** + * \file + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +// ** INCLUDES ** +#include "VtkImageDataToSurfacePointsFilter.h" + +#include <algorithm> +#include <vector> + +#include <vtkIdList.h> +#include <vtkImageData.h> +#include <vtkInformation.h> +#include <vtkInformationVector.h> +#include <vtkLine.h> +#include <vtkObjectFactory.h> +#include <vtkPointData.h> +#include <vtkPolyData.h> +#include <vtkSmartPointer.h> + +vtkStandardNewMacro(VtkImageDataToSurfacePointsFilter); + +VtkImageDataToSurfacePointsFilter::VtkImageDataToSurfacePointsFilter() + : PointsPerPixel(20) +{ +} + +void VtkImageDataToSurfacePointsFilter::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); +} + +int VtkImageDataToSurfacePointsFilter::FillInputPortInformation(int, vtkInformation* info) +{ + info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkImageData"); + return 1; +} + +int VtkImageDataToSurfacePointsFilter::RequestData( + vtkInformation*, + vtkInformationVector** inputVector, + vtkInformationVector* outputVector) +{ + vtkDebugMacro(<< "Executing VtkImageDataToSurfacePointsFilter"); + + vtkInformation* input_info = inputVector[0]->GetInformationObject(0); + vtkInformation* output_info = outputVector->GetInformationObject(0); + vtkImageData* input = vtkImageData::SafeDownCast( + input_info->Get(vtkDataObject::DATA_OBJECT())); + vtkPolyData* output = vtkPolyData::SafeDownCast( + output_info->Get(vtkDataObject::DATA_OBJECT())); + + void* pixvals = input->GetScalarPointer(); + int n_comp = input->GetNumberOfScalarComponents(); + if (n_comp < 1) + vtkDebugMacro("Error reading number of components."); + if (n_comp > 2) + vtkDebugMacro("RGB colours detected. Using only first channel for computation."); + + std::size_t const n_points = static_cast<std::size_t>(input->GetNumberOfPoints()); + if (n_points == 0) + { + vtkDebugMacro("No data found!"); + return -1; + } + + double spacing[3]; + input->GetSpacing(spacing); + int dimensions[3]; + input->GetDimensions(dimensions); + double origin[3]; + input->GetOrigin(origin); + MathLib::Point3d const ll(std::array<double, 3>{{origin[0], origin[1], origin[2]}}); + GeoLib::RasterHeader const header = { + static_cast<std::size_t>(dimensions[0]), + static_cast<std::size_t>(dimensions[1]), + 1, ll, spacing[0], -9999}; + std::vector<double> pixels; + pixels.reserve(n_points); + for (std::size_t i = 0; i < n_points; ++i) + { + if ((n_comp == 2 || n_comp == 4) && + (((float*)pixvals)[(i + 1) * n_comp - 1] < 0.00000001f)) + pixels.push_back(-9999); + else + pixels.push_back(((float*)pixvals)[i * n_comp]); + } + GeoLib::Raster const* const raster(new GeoLib::Raster(header, pixels.begin(), pixels.end())); + + vtkSmartPointer<vtkPoints> new_points = vtkSmartPointer<vtkPoints>::New(); + new_points->SetNumberOfPoints(PointsPerPixel * n_points); + vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New(); + cells->Allocate(PointsPerPixel * n_points); + + vtkPointData* output_data = output->GetPointData(); + double const half_cellsize(spacing[0] / 2.0); + std::size_t pnt_idx(0); + for (std::size_t i = 0; i < static_cast<std::size_t>(n_points); ++i) + { + // Skip transparent pixels + if (n_comp == 2 || n_comp == 4) + { + if (((float*)pixvals)[(i + 1) * n_comp - 1] < 0.00000001f) + continue; + } + + double p[3]; + input->GetPoint(i, p); + MathLib::Point3d min_pnt{std::array<double, 3>{ + {p[0] - half_cellsize, p[1] - half_cellsize, 0}}}; + MathLib::Point3d max_pnt{std::array<double, 3>{ + {p[0] + half_cellsize, p[1] + half_cellsize, 0}}}; + createPointSurface( + new_points, cells, pnt_idx, min_pnt, max_pnt, *raster); + pnt_idx += PointsPerPixel; + } + + output->SetPoints(new_points); + output->SetVerts(cells); + output->Squeeze(); + + vtkDebugMacro(<< "Created " << new_points->GetNumberOfPoints() << " points."); + return 1; +} + +void VtkImageDataToSurfacePointsFilter::createPointSurface( + vtkSmartPointer<vtkPoints>& points, + vtkSmartPointer<vtkCellArray>& cells, + std::size_t pnt_idx, + MathLib::Point3d const& min_pnt, + MathLib::Point3d const& max_pnt, + GeoLib::Raster const& raster) +{ + std::size_t const n_points(static_cast<std::size_t>(this->GetPointsPerPixel())); + for (std::size_t i = 0; i < n_points; ++i) + { + double p[3]; + p[0] = getRandomNumber(min_pnt[0], max_pnt[0]); + p[1] = getRandomNumber(min_pnt[1], max_pnt[1]); + p[2] = raster.interpolateValueAtPoint(MathLib::Point3d({p[0], p[1], 0})); + points->SetPoint(pnt_idx + i, p); + cells->InsertNextCell(1); + cells->InsertCellPoint(pnt_idx + i); + } +} + +double VtkImageDataToSurfacePointsFilter::getRandomNumber(double const& min, double const& max) const +{ + return ((double)std::rand() / RAND_MAX) * (max - min) + min; +} diff --git a/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.h b/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.h new file mode 100644 index 0000000000000000000000000000000000000000..6b038770f8fef7b27bd48c631394499d29088e73 --- /dev/null +++ b/Applications/DataExplorer/VtkVis/VtkImageDataToSurfacePointsFilter.h @@ -0,0 +1,61 @@ +/** + * \file + * \copyright + * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org) + * Distributed under a Modified BSD License. + * See accompanying file LICENSE.txt or + * http://www.opengeosys.org/project/license + * + */ + +#pragma once + +#include <vtkPolyDataAlgorithm.h> + +#include "GeoLib/Raster.h" +#include "VtkAlgorithmProperties.h" + +class vtkInformation; +class vtkInformationVector; + +/// A VTK filter that creates a point cloud representing a surface defined by +/// pixel values +class VtkImageDataToSurfacePointsFilter : public vtkPolyDataAlgorithm +{ +public: + /// Create a new objects (required because of VTKs reference counting) + static VtkImageDataToSurfacePointsFilter* New(); + + vtkTypeMacro(VtkImageDataToSurfacePointsFilter, vtkPolyDataAlgorithm); + + /// @brief Prints information about itself. + void PrintSelf(ostream& os, vtkIndent indent) override; + + vtkGetMacro(PointsPerPixel, vtkIdType); + vtkSetMacro(PointsPerPixel, vtkIdType); + +protected: + VtkImageDataToSurfacePointsFilter(); + ~VtkImageDataToSurfacePointsFilter() = default; + + /// Sets input port to vtkImageData. + int FillInputPortInformation(int port, vtkInformation* info) override; + + /// Updates the graphical object + int RequestData(vtkInformation* request, vtkInformationVector** inputVector, + vtkInformationVector* outputVector) override; + +private: + /// Creates the point objects based on the parameters set by the user + void createPointSurface(vtkSmartPointer<vtkPoints>& points, + vtkSmartPointer<vtkCellArray>& cells, + std::size_t pnt_idx, + MathLib::Point3d const& min_pnt, + MathLib::Point3d const& max_pnt, + GeoLib::Raster const& raster); + + /// Returns a random number in [min, max] + double getRandomNumber(double const& min, double const& max) const; + + vtkIdType PointsPerPixel; +};