From 1db2ac83a553c3a29d1de1ea70ef3eea9c341111 Mon Sep 17 00:00:00 2001
From: Karsten Rink <karsten.rink@ufz.de>
Date: Fri, 24 Nov 2017 10:15:42 +0100
Subject: [PATCH] Raster to mesh conversion (#1991)
* reimplemented raster to mesh conversion
* fixed index issues and added 3d index calculation
* added necessary changes to gui dialogue
* added tests
* clang format
* added ifdef for data explorer/vtk includes
* modified dereferencing of unique pointers, formatting
* changed pointer to vector in vtk prep
* using const& for nodes vec
* fixed remaining occurences of dereferencing raster object
---
.../DataView/MeshElementRemovalDialog.cpp | 2 +-
.../NetCdfDialog/NetCdfConfigureDialog.cpp | 3 +-
.../DataExplorer/VtkVis/CMakeLists.txt | 3 +-
.../DataExplorer/VtkVis/MeshFromRaster.ui | 10 +
.../VtkVis/MeshFromRasterDialog.cpp | 18 +-
.../VtkVis/MeshFromRasterDialog.h | 2 +
.../VtkCompositeTextureOnSurfaceFilter.cpp | 4 -
.../DataExplorer/VtkVis/VtkRaster.cpp | 8 +-
Applications/FileIO/AsciiRasterInterface.cpp | 3 +
GeoLib/Raster.cpp | 3 +-
GeoLib/Raster.h | 4 +-
MeshLib/MeshGenerators/RasterToMesh.cpp | 311 +++++++++---------
MeshLib/MeshGenerators/RasterToMesh.h | 59 +---
MeshLib/MeshSearch/ElementSearch.h | 11 +-
Tests/CMakeLists.txt | 2 +-
Tests/MeshLib/TestMeshValidation.cpp | 4 +-
Tests/MeshLib/TestNodeSearch.cpp | 4 +-
Tests/MeshLib/TestRasterToMesh.cpp | 304 +++++++++++++++++
Tests/lfs-data/MeshLib/testraster_selke.asc | 39 +++
19 files changed, 566 insertions(+), 228 deletions(-)
create mode 100644 Tests/MeshLib/TestRasterToMesh.cpp
create mode 100644 Tests/lfs-data/MeshLib/testraster_selke.asc
diff --git a/Applications/DataExplorer/DataView/MeshElementRemovalDialog.cpp b/Applications/DataExplorer/DataView/MeshElementRemovalDialog.cpp
index dfcd451863a..7e52a1a6366 100644
--- a/Applications/DataExplorer/DataView/MeshElementRemovalDialog.cpp
+++ b/Applications/DataExplorer/DataView/MeshElementRemovalDialog.cpp
@@ -74,7 +74,7 @@ void MeshElementRemovalDialog::accept()
{
QList<QListWidgetItem*> items = this->materialListWidget->selectedItems();
for (auto& item : items)
- ex.searchByPropertyValue(item->text().toInt());
+ ex.searchByPropertyValue(item->text().toInt(), "MaterialIDs");
anything_checked = true;
}
if (this->boundingBoxCheckBox->isChecked())
diff --git a/Applications/DataExplorer/NetCdfDialog/NetCdfConfigureDialog.cpp b/Applications/DataExplorer/NetCdfDialog/NetCdfConfigureDialog.cpp
index 88b9095500b..f4897339419 100644
--- a/Applications/DataExplorer/NetCdfDialog/NetCdfConfigureDialog.cpp
+++ b/Applications/DataExplorer/NetCdfDialog/NetCdfConfigureDialog.cpp
@@ -353,7 +353,8 @@ void NetCdfConfigureDialog::createDataObject()
if (originLat > lastLat) // reverse lines in vertical direction if the original file has its origin in the northwest corner
this->reverseNorthSouth(data_array, sizeLon, sizeLat);
- GeoLib::RasterHeader header = {sizeLon,sizeLat,origin,resolution,-9999};
+ GeoLib::RasterHeader const header = {sizeLon, sizeLat, 1,
+ origin, resolution, -9999};
if (this->radioMesh->isChecked())
{
MeshLib::MeshElemType meshElemType = MeshLib::MeshElemType::QUAD;
diff --git a/Applications/DataExplorer/VtkVis/CMakeLists.txt b/Applications/DataExplorer/VtkVis/CMakeLists.txt
index 1ff4e6f778c..b746c07548b 100644
--- a/Applications/DataExplorer/VtkVis/CMakeLists.txt
+++ b/Applications/DataExplorer/VtkVis/CMakeLists.txt
@@ -127,9 +127,8 @@ if(GEOTIFF_FOUND)
endif() # GEOTIFF_FOUND
target_link_libraries(VtkVis
- PUBLIC BaseLib GeoLib MeshLib DataHolderLib QtBase
+ PUBLIC BaseLib GeoLib MeshLib DataHolderLib QtBase ${VTK_LIBRARIES} NetCdfDialogLib
PRIVATE MathLib ApplicationsFileIO Qt5::Gui logog
- INTERFACE ${VTK_LIBRARIES}
)
set_property(TARGET VtkVis PROPERTY FOLDER "DataExplorer")
diff --git a/Applications/DataExplorer/VtkVis/MeshFromRaster.ui b/Applications/DataExplorer/VtkVis/MeshFromRaster.ui
index 61a16073870..746fd2e790e 100644
--- a/Applications/DataExplorer/VtkVis/MeshFromRaster.ui
+++ b/Applications/DataExplorer/VtkVis/MeshFromRaster.ui
@@ -125,6 +125,16 @@
</property>
</widget>
</item>
+ <item>
+ <widget class="QRadioButton" name="prismButton">
+ <property name="enabled">
+ <bool>false</bool>
+ </property>
+ <property name="text">
+ <string>set of two prisms</string>
+ </property>
+ </widget>
+ </item>
<item>
<widget class="QRadioButton" name="hexButton">
<property name="enabled">
diff --git a/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.cpp b/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.cpp
index d06af5490d0..9823b8cd289 100644
--- a/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.cpp
+++ b/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.cpp
@@ -17,7 +17,7 @@
#include "MeshGenerators/VtkMeshConverter.h"
MeshFromRasterDialog::MeshFromRasterDialog(QDialog* parent)
-: QDialog(parent), _mesh_name("mesh"), _array_name("Colour")
+: QDialog(parent), _mesh_name("mesh"), _array_name("MaterialIDs")
{
setupUi(this);
@@ -28,6 +28,20 @@ MeshFromRasterDialog::MeshFromRasterDialog(QDialog* parent)
MeshFromRasterDialog::~MeshFromRasterDialog() = default;
+void MeshFromRasterDialog::on_elevationButton_toggled(bool isChecked)
+{
+ if (isChecked)
+ {
+ if (this->prismButton->isChecked())
+ this->triButton->setChecked(true);
+ if (this->hexButton->isChecked())
+ this->quadButton->setChecked(true);
+ }
+
+ this->prismButton->setEnabled(!isChecked);
+ this->hexButton->setEnabled(!isChecked);
+}
+
void MeshFromRasterDialog::accept()
{
if (this->mshNameEdit->text().isEmpty())
@@ -54,6 +68,8 @@ void MeshFromRasterDialog::accept()
}
_element_selection = MeshLib::MeshElemType::TRIANGLE;
if (this->quadButton->isChecked()) _element_selection = MeshLib::MeshElemType::QUAD;
+ else if (this->prismButton->isChecked())
+ _element_selection = MeshLib::MeshElemType::PRISM;
else if (this->hexButton->isChecked()) _element_selection = MeshLib::MeshElemType::HEXAHEDRON;
this->done(QDialog::Accepted);
diff --git a/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.h b/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.h
index 2933ef04c09..68bff68f873 100644
--- a/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.h
+++ b/Applications/DataExplorer/VtkVis/MeshFromRasterDialog.h
@@ -43,6 +43,8 @@ public:
MeshLib::UseIntensityAs getIntensitySelection() const { return _intensity_selection; }
private slots:
+ void on_elevationButton_toggled(bool isChecked);
+
/// Instructions if the OK-Button has been pressed.
void accept() override;
diff --git a/Applications/DataExplorer/VtkVis/VtkCompositeTextureOnSurfaceFilter.cpp b/Applications/DataExplorer/VtkVis/VtkCompositeTextureOnSurfaceFilter.cpp
index c6d5aaf410e..c253f0c6cc4 100644
--- a/Applications/DataExplorer/VtkVis/VtkCompositeTextureOnSurfaceFilter.cpp
+++ b/Applications/DataExplorer/VtkVis/VtkCompositeTextureOnSurfaceFilter.cpp
@@ -90,10 +90,6 @@ void VtkCompositeTextureOnSurfaceFilter::init()
image->GetOutput()->GetOrigin(origin);
double spacing[3];
image->GetOutput()->GetSpacing(spacing);
-/*
- VtkCompositeColormapToImageFilter* cm = new VtkCompositeColormapToImageFilter(image);
- vtkImageAlgorithm* img = dynamic_cast<vtkImageAlgorithm*>(cm->GetOutputAlgorithm());
-*/
surface->SetRaster(image, origin[0], origin[1], spacing[0]);
surface->Update();
}
diff --git a/Applications/DataExplorer/VtkVis/VtkRaster.cpp b/Applications/DataExplorer/VtkVis/VtkRaster.cpp
index 8a65b927c4e..1ee2b08499e 100644
--- a/Applications/DataExplorer/VtkVis/VtkRaster.cpp
+++ b/Applications/DataExplorer/VtkVis/VtkRaster.cpp
@@ -70,7 +70,7 @@ vtkImageAlgorithm* VtkRaster::loadImage(const std::string &fileName,
}
vtkImageImport* VtkRaster::loadImageFromArray(double const*const data_array, GeoLib::RasterHeader header)
{
- const unsigned length = header.n_rows*header.n_cols;
+ const unsigned length = header.n_rows * header.n_cols * header.n_depth;
auto* data = new float[length * 2];
float max_val = *std::max_element(data_array, data_array+length);
for (unsigned j=0; j<length; ++j)
@@ -88,8 +88,10 @@ vtkImageImport* VtkRaster::loadImageFromArray(double const*const data_array, Geo
vtkImageImport* image = vtkImageImport::New();
image->SetDataSpacing(header.cell_size, header.cell_size, header.cell_size);
image->SetDataOrigin(header.origin[0]+(header.cell_size/2.0), header.origin[1]+(header.cell_size/2.0), 0); // translate whole mesh by half a pixel in x and y
- image->SetWholeExtent(0, header.n_cols-1, 0, header.n_rows-1, 0, 0);
- image->SetDataExtent(0, header.n_cols-1, 0, header.n_rows-1, 0, 0);
+ image->SetWholeExtent(0, header.n_cols - 1, 0, header.n_rows - 1, 0,
+ header.n_depth - 1);
+ image->SetDataExtent(0, header.n_cols - 1, 0, header.n_rows - 1, 0,
+ header.n_depth - 1);
image->SetDataExtentToWholeExtent();
image->SetDataScalarTypeToFloat();
image->SetNumberOfScalarComponents(2);
diff --git a/Applications/FileIO/AsciiRasterInterface.cpp b/Applications/FileIO/AsciiRasterInterface.cpp
index fd7d1468776..560192ab3ee 100644
--- a/Applications/FileIO/AsciiRasterInterface.cpp
+++ b/Applications/FileIO/AsciiRasterInterface.cpp
@@ -86,6 +86,8 @@ bool AsciiRasterInterface::readASCHeader(std::ifstream &in, GeoLib::RasterHeader
header.n_rows = atoi(value.c_str());
} else return false;
+ header.n_depth = 1;
+
in >> tag;
if (tag == "xllcorner" || tag == "xllcenter")
{
@@ -192,6 +194,7 @@ bool AsciiRasterInterface::readSurferHeader(
ERR("Error in readSurferHeader() - Anisotropic cellsize detected.");
return false;
}
+ header.n_depth = 1;
header.no_data = min - 1;
in >> min >> max;
diff --git a/GeoLib/Raster.cpp b/GeoLib/Raster.cpp
index 625d37dd808..2a2fdff78dd 100644
--- a/GeoLib/Raster.cpp
+++ b/GeoLib/Raster.cpp
@@ -85,7 +85,8 @@ Raster* Raster::getRasterFromSurface(Surface const& sfc, double cell_size, doubl
}
}
- RasterHeader header = {std::size_t(n_cols), std::size_t(n_rows), MathLib::Point3d(ll), cell_size, static_cast<double>(-9999)};
+ RasterHeader header = { std::size_t(n_cols), std::size_t(n_rows), 1,
+ MathLib::Point3d(ll), cell_size, static_cast<double>(-9999) };
return new Raster(header, z_vals, z_vals+n_cols*n_rows);
}
diff --git a/GeoLib/Raster.h b/GeoLib/Raster.h
index 07a5914c199..32e6229b320 100644
--- a/GeoLib/Raster.h
+++ b/GeoLib/Raster.h
@@ -13,17 +13,19 @@
#pragma once
+#include <array>
#include <utility>
#include "Surface.h"
namespace GeoLib {
-/// Contains the relevant information when handling with geoscientific raster data
+/// Contains the relevant information when storing a geoscientific raster data
struct RasterHeader
{
std::size_t n_cols; // width
std::size_t n_rows; // height
+ std::size_t n_depth; // depth (for 3d image)
MathLib::Point3d origin; // lower left corner
double cell_size; // edge length of each pixel
double no_data; // no data value
diff --git a/MeshLib/MeshGenerators/RasterToMesh.cpp b/MeshLib/MeshGenerators/RasterToMesh.cpp
index dba1407c6c2..5b4e8089959 100644
--- a/MeshLib/MeshGenerators/RasterToMesh.cpp
+++ b/MeshLib/MeshGenerators/RasterToMesh.cpp
@@ -12,6 +12,9 @@
#include "MeshLib/Elements/Elements.h"
#include "MeshLib/Mesh.h"
#include "MeshLib/Node.h"
+#include "MeshLib/MeshEditing/RemoveMeshComponents.h"
+#include "MeshLib/MeshGenerators/MeshGenerator.h"
+#include "MeshLib/MeshSearch/ElementSearch.h"
#include <vtkImageData.h>
#include <vtkDataArray.h>
@@ -40,14 +43,25 @@ MeshLib::Mesh* RasterToMesh::convert(
UseIntensityAs intensity_type,
std::string const& array_name)
{
- if ((elem_type != MeshElemType::TRIANGLE) && (elem_type != MeshElemType::QUAD))
+ if ((elem_type != MeshElemType::TRIANGLE) &&
+ (elem_type != MeshElemType::QUAD) &&
+ (elem_type != MeshElemType::HEXAHEDRON) &&
+ (elem_type != MeshElemType::PRISM))
{
- ERR("VtkMeshConverter::convertImgToMesh(): Invalid Mesh Element Type.");
+ ERR("Invalid Mesh Element Type.");
return nullptr;
}
- vtkSmartPointer<vtkDataArray> pixelData = vtkSmartPointer<vtkDataArray>(img->GetPointData()->GetScalars());
int* dims = img->GetDimensions();
+ if (((elem_type == MeshElemType::TRIANGLE) || (elem_type == MeshElemType::QUAD)) &&
+ dims[2] != 1)
+ {
+ ERR("Triangle or Quad elements cannot be used to construct meshes from 3D rasters.");
+ return nullptr;
+ }
+
+ vtkSmartPointer<vtkDataArray> pixelData =
+ vtkSmartPointer<vtkDataArray>(img->GetPointData()->GetScalars());
int nTuple = pixelData->GetNumberOfComponents();
if (nTuple < 1 || nTuple > 4)
{
@@ -55,193 +69,168 @@ MeshLib::Mesh* RasterToMesh::convert(
return nullptr;
}
- MathLib::Point3d const orig (std::array<double,3>{{origin[0], origin[1], origin[2]}});
- GeoLib::RasterHeader const header =
- {static_cast<std::size_t>(dims[0]), static_cast<std::size_t>(dims[1]), orig, scalingFactor, -9999};
- const std::size_t incHeight = header.n_rows+1;
- const std::size_t incWidth = header.n_cols+1;
- std::vector<double> pix_val (incHeight * incWidth, std::numeric_limits<double>::max());
- std::vector<bool> pix_vis (header.n_rows * header.n_cols, false);
+ MathLib::Point3d const orig(std::array<double, 3>{ {origin[0], origin[1], origin[2]}});
+ GeoLib::RasterHeader const header = {static_cast<std::size_t>(dims[0]),
+ static_cast<std::size_t>(dims[1]),
+ static_cast<std::size_t>(dims[2]),
+ orig,
+ scalingFactor,
+ -9999};
- for (std::size_t i = 0; i < header.n_rows; i++)
- for (std::size_t j = 0; j < header.n_cols; j++)
+ std::vector<double> pix(header.n_cols * header.n_rows * header.n_depth, 0);
+ for (std::size_t k = 0; k < header.n_depth; k++)
+ {
+ std::size_t const layer_idx = (k*header.n_rows*header.n_cols);
+ for (std::size_t i = 0; i < header.n_rows; i++)
{
- std::size_t const img_idx = i*header.n_cols + j;
- std::size_t const fld_idx = i*incWidth + j;
-
- // colour of current pixel
- double* colour = pixelData->GetTuple(img_idx);
- // is current pixel visible?
- bool const visible = (nTuple == 2 || nTuple == 4) ? (colour[nTuple-1] != 0) : true;
- if (!visible)
- continue;
-
- double const value = (nTuple < 3) ?
- colour[0] : // grey (+ alpha)
- (0.3 * colour[0] + 0.6 * colour[1] + 0.1 * colour[2]); // rgb(a)
- pix_vis[img_idx] = true;
- pix_val[fld_idx] = value;
- pix_val[fld_idx+1] = value;
- pix_val[fld_idx+incWidth] = value;
- pix_val[fld_idx+incWidth+1] = value;
+ std::size_t const idx = i * header.n_cols + layer_idx;
+ for (std::size_t j = 0; j < header.n_cols; j++)
+ {
+ double* colour = pixelData->GetTuple(idx + j);
+ bool const visible = (nTuple == 2 || nTuple == 4) ? (colour[nTuple - 1] != 0) : true;
+ if (!visible)
+ pix[idx + j] = header.no_data;
+ else
+ pix[idx + j] = (nTuple < 3) ?
+ colour[0] : // grey (+ alpha)
+ (0.3 * colour[0] + 0.6 * colour[1] + 0.1 * colour[2]); // rgb(a)
+ }
}
+ }
- return constructMesh(pix_val, array_name, pix_vis, header, elem_type, intensity_type);
+ MeshLib::Mesh* mesh = convert(pix.data(), header, elem_type, intensity_type, array_name);
+ return mesh;
}
MeshLib::Mesh* RasterToMesh::convert(
- double const* img,
+ double const*const img,
GeoLib::RasterHeader const& header,
MeshElemType elem_type,
UseIntensityAs intensity_type,
std::string const& array_name)
{
- if ((elem_type != MeshElemType::TRIANGLE) && (elem_type != MeshElemType::QUAD))
+ if ((elem_type != MeshElemType::TRIANGLE) &&
+ (elem_type != MeshElemType::QUAD) &&
+ (elem_type != MeshElemType::HEXAHEDRON) &&
+ (elem_type != MeshElemType::PRISM))
{
ERR("Invalid Mesh Element Type.");
return nullptr;
}
- std::size_t const incHeight (header.n_rows+1);
- std::size_t const incWidth (header.n_cols+1);
- std::vector<double> pix_val (incHeight * incWidth, std::numeric_limits<double>::max());
- std::vector<bool> pix_vis (incHeight * incWidth, false);
-
- for (std::size_t i = 0; i < header.n_rows; i++)
- for (std::size_t j = 0; j < header.n_cols; j++)
- {
- std::size_t const img_idx = i*header.n_cols + j;
- std::size_t const fld_idx = i*incWidth + j;
- if (img[img_idx] == -9999)
- continue;
-
- pix_vis[img_idx] = true;
- pix_val[fld_idx] = img[img_idx];
- pix_val[fld_idx+1] = img[img_idx];
- pix_val[fld_idx+incWidth] = img[img_idx];
- pix_val[fld_idx+incWidth+1] = img[img_idx];
- }
-
- return constructMesh(pix_val, array_name, pix_vis, header, elem_type, intensity_type);
-}
-
-MeshLib::Mesh* RasterToMesh::constructMesh(
- std::vector<double> const& pix_val,
- std::string const& array_name,
- std::vector<bool> const& pix_vis,
- GeoLib::RasterHeader const& header,
- MeshLib::MeshElemType elem_type,
- MeshLib::UseIntensityAs intensity_type)
-{
- std::vector<int> node_idx_map ((header.n_rows+1) * (header.n_cols+1), -1);
- bool const use_elevation (intensity_type == MeshLib::UseIntensityAs::ELEVATION);
- std::vector<MeshLib::Node*> nodes (createNodeVector(pix_val, node_idx_map, header, use_elevation));
- if (nodes.empty())
+ if (((elem_type == MeshElemType::TRIANGLE) || (elem_type == MeshElemType::QUAD)) &&
+ header.n_depth != 1)
+ {
+ ERR("Triangle or Quad elements cannot be used to construct meshes from 3D rasters.");
return nullptr;
+ }
- std::vector<MeshLib::Element*> elements(createElementVector(
- pix_vis, nodes, node_idx_map, header.n_rows, header.n_cols, elem_type));
- if (elements.empty())
+ if (intensity_type == UseIntensityAs::ELEVATION &&
+ ((elem_type == MeshElemType::PRISM) || (elem_type == MeshElemType::HEXAHEDRON)))
+ {
+ ERR("Elevation mapping can only be performed for 2D meshes.");
return nullptr;
+ }
- MeshLib::Properties properties;
- if (intensity_type == MeshLib::UseIntensityAs::MATERIALS)
+ MathLib::Point3d mesh_origin({header.origin[0] - (header.cell_size / 2.0),
+ header.origin[1] - (header.cell_size / 2.0),
+ header.origin[2]});
+ std::unique_ptr<MeshLib::Mesh> mesh (nullptr);
+ if (elem_type == MeshElemType::TRIANGLE)
+ mesh.reset(
+ MeshLib::MeshGenerator::generateRegularTriMesh(header.n_cols,
+ header.n_rows,
+ header.cell_size,
+ mesh_origin,
+ "RasterDataMesh"));
+ else if (elem_type == MeshElemType::QUAD)
+ mesh.reset(
+ MeshLib::MeshGenerator::generateRegularQuadMesh(header.n_cols,
+ header.n_rows,
+ header.cell_size,
+ mesh_origin,
+ "RasterDataMesh"));
+ else if (elem_type == MeshElemType::PRISM)
+ mesh.reset(
+ MeshLib::MeshGenerator::generateRegularPrismMesh(header.n_cols,
+ header.n_rows,
+ header.n_depth,
+ header.cell_size,
+ mesh_origin,
+ "RasterDataMesh"));
+ else if (elem_type == MeshElemType::HEXAHEDRON)
+ mesh.reset(
+ MeshLib::MeshGenerator::generateRegularHexMesh(header.n_cols,
+ header.n_rows,
+ header.n_depth,
+ header.cell_size,
+ mesh_origin,
+ "RasterDataMesh"));
+
+ MeshLib::Mesh* new_mesh(nullptr);
+ std::vector<std::size_t> elements_to_remove;
+ if (intensity_type == UseIntensityAs::ELEVATION)
{
- auto* const prop_vec = properties.createNewPropertyVector<int>(
- "MaterialIDs", MeshLib::MeshItemType::Cell, 1);
- fillPropertyVector<int>(*prop_vec, pix_val, pix_vis, header.n_rows, header.n_cols, elem_type);
+ std::vector<MeshLib::Node*> const& nodes(mesh->getNodes());
+ std::vector<MeshLib::Element*> const& elems(mesh->getElements());
+ std::size_t const n_nodes(elems[0]->getNumberOfNodes());
+ bool const double_idx = (elem_type == MeshElemType::TRIANGLE) || (elem_type == MeshElemType::PRISM);
+ std::size_t const m = (double_idx) ? 2 : 1;
+ for (std::size_t k = 0; k < header.n_depth; k++)
+ {
+ std::size_t const layer_idx = (k*header.n_rows*header.n_cols);
+ for (std::size_t i = 0; i < header.n_cols; i++)
+ {
+ std::size_t const idx(i * header.n_rows + layer_idx);
+ for (std::size_t j = 0; j < header.n_rows; j++)
+ {
+ double const val(img[idx + j]);
+ if (val == header.no_data)
+ {
+ elements_to_remove.push_back(m * (idx + j));
+ if (double_idx)
+ elements_to_remove.push_back(m * (idx + j) + 1);
+ continue;
+ }
+ for (std::size_t n = 0; n < n_nodes; ++n)
+ {
+ (*(nodes[elems[m * (idx + j)]->getNodeIndex(n)]))[2] = val;
+ if (double_idx)
+ (*(nodes[elems[m * (idx + j) + 1]->getNodeIndex(n)]))[2] = val;
+ }
+ }
+ }
+ }
}
- else if (intensity_type == MeshLib::UseIntensityAs::DATAVECTOR)
+ else
{
- auto* const prop_vec = properties.createNewPropertyVector<double>(
- array_name, MeshLib::MeshItemType::Cell, 1);
- fillPropertyVector<double>(*prop_vec, pix_val, pix_vis, header.n_rows, header.n_cols, elem_type);
- }
-
- return new MeshLib::Mesh("RasterDataMesh", nodes, elements, properties);
-}
-
-std::vector<MeshLib::Node*> RasterToMesh::createNodeVector(
- std::vector<double> const& elevation,
- std::vector<int> & node_idx_map,
- GeoLib::RasterHeader const& header,
- bool use_elevation)
-{
- std::size_t node_idx_count(0);
- double const x_offset(header.origin[0] - header.cell_size/2.0);
- double const y_offset(header.origin[1] - header.cell_size/2.0);
- std::vector<MeshLib::Node*> nodes;
- for (std::size_t i = 0; i < (header.n_rows+1); i++)
- for (std::size_t j = 0; j < (header.n_cols+1); j++)
+ MeshLib::Properties &properties = mesh->getProperties();
+ MeshLib::ElementSearch ex(*mesh);
+ if (array_name == "MaterialIDs")
{
- std::size_t const index = i * (header.n_cols+1) + j;
- if (elevation[index] == std::numeric_limits<double>::max())
- continue;
-
- double const zValue = (use_elevation) ? elevation[index] : 0;
- auto* node(new MeshLib::Node(x_offset + (header.cell_size * j),
- y_offset + (header.cell_size * i),
- zValue));
- nodes.push_back(node);
- node_idx_map[index] = node_idx_count;
- node_idx_count++;
+ auto* const prop_vec = properties.createNewPropertyVector<int>(
+ array_name, MeshLib::MeshItemType::Cell, 1);
+ fillPropertyVector<int>(*prop_vec, img, header, elem_type);
+ ex.searchByPropertyValue<int>(header.no_data, array_name);
}
- return nodes;
-}
-
-std::vector<MeshLib::Element*> RasterToMesh::createElementVector(
- std::vector<bool> const& pix_vis,
- std::vector<MeshLib::Node*> const& nodes,
- std::vector<int> const&node_idx_map,
- std::size_t const imgHeight,
- std::size_t const imgWidth,
- MeshElemType elem_type)
-{
- std::vector<MeshLib::Element*> elements;
- std::size_t const incWidth (imgWidth+1);
- for (std::size_t i = 0; i < imgHeight; i++)
- for (std::size_t j = 0; j < imgWidth; j++)
+ else
{
- if (!pix_vis[i*imgWidth+j])
- continue;
-
- int const idx = i * incWidth + j;
- if (elem_type == MeshElemType::TRIANGLE)
- {
- auto** tri1_nodes = new MeshLib::Node*[3];
- tri1_nodes[0] = nodes[node_idx_map[idx]];
- tri1_nodes[1] = nodes[node_idx_map[idx+1]];
- tri1_nodes[2] = nodes[node_idx_map[idx+incWidth]];
-
- auto** tri2_nodes = new MeshLib::Node*[3];
- tri2_nodes[0] = nodes[node_idx_map[idx+1]];
- tri2_nodes[1] = nodes[node_idx_map[idx+incWidth+1]];
- tri2_nodes[2] = nodes[node_idx_map[idx+incWidth]];
-
- elements.push_back(new MeshLib::Tri(tri1_nodes)); // upper left triangle
- elements.push_back(new MeshLib::Tri(tri2_nodes)); // lower right triangle
- }
- else if (elem_type == MeshElemType::QUAD)
- {
- auto** quad_nodes = new MeshLib::Node*[4];
- quad_nodes[0] = nodes[node_idx_map[idx]];
- quad_nodes[1] = nodes[node_idx_map[idx + 1]];
- quad_nodes[2] = nodes[node_idx_map[idx + incWidth + 1]];
- quad_nodes[3] = nodes[node_idx_map[idx + incWidth]];
- elements.push_back(new MeshLib::Quad(quad_nodes));
- }
+ auto* const prop_vec = properties.createNewPropertyVector<double>(
+ array_name, MeshLib::MeshItemType::Cell, 1);
+ fillPropertyVector<double>(*prop_vec, img, header, elem_type);
+ ex.searchByPropertyValue<double>(header.no_data, array_name);
}
- return elements;
-}
-
-double RasterToMesh::getExistingValue(const double* img, std::size_t length)
-{
- for (std::size_t i=0; i<length; i++)
- {
- if (img[i] != -9999)
- return img[i];
+ elements_to_remove = ex.getSearchedElementIDs();
}
- return -9999;
+ if (!elements_to_remove.empty())
+ new_mesh = MeshLib::removeElements(*mesh, elements_to_remove, mesh->getName());
+ else
+ new_mesh = mesh.release();
+
+ if (intensity_type == UseIntensityAs::NONE)
+ new_mesh->getProperties().removePropertyVector(array_name);
+
+ return new_mesh;
}
} // end namespace MeshLib
diff --git a/MeshLib/MeshGenerators/RasterToMesh.h b/MeshLib/MeshGenerators/RasterToMesh.h
index 8ea537a7493..85487f11e00 100644
--- a/MeshLib/MeshGenerators/RasterToMesh.h
+++ b/MeshLib/MeshGenerators/RasterToMesh.h
@@ -75,61 +75,36 @@ public:
* \param array_name mesh property name, defaults to "Colour" if not
* given.
*/
- static MeshLib::Mesh* convert(const double* img,
- GeoLib::RasterHeader const& header,
- MeshElemType elem_type,
- UseIntensityAs intensity_type,
- std::string const& array_name = "Colour");
-
-private:
- static MeshLib::Mesh* constructMesh(
- std::vector<double> const& pix_val,
- std::string const& array_name,
- std::vector<bool> const& pix_vis,
- GeoLib::RasterHeader const& header,
- MeshLib::MeshElemType elem_type,
- MeshLib::UseIntensityAs intensity_type);
-
- static std::vector<MeshLib::Node*> createNodeVector(
- std::vector<double> const& elevation,
- std::vector<int> &node_idx_map,
+ static MeshLib::Mesh* convert(const double*const img,
GeoLib::RasterHeader const& header,
- bool use_elevation);
-
- static std::vector<MeshLib::Element*> createElementVector(
- std::vector<bool> const& pix_vis,
- std::vector<MeshLib::Node*> const& nodes,
- std::vector<int> const& node_idx_map,
- std::size_t const imgHeight,
- std::size_t const imgWidth,
- MeshElemType elem_type);
+ MeshElemType elem_type,
+ UseIntensityAs intensity_type,
+ std::string const& array_name = "Colour");
+private:
template<typename T>
static void fillPropertyVector(
MeshLib::PropertyVector<T> &prop_vec,
- std::vector<double> const& pix_val,
- std::vector<bool> const& pix_vis,
- const std::size_t &imgHeight,
- const std::size_t &imgWidth,
+ double const*const img,
+ GeoLib::RasterHeader const& header,
MeshElemType elem_type)
{
- for (std::size_t i = 0; i < imgHeight; i++)
- for (std::size_t j = 0; j < imgWidth; j++)
+ for (std::size_t k = 0; k < header.n_depth; k++)
+ {
+ std::size_t const layer_idx = (k*header.n_rows*header.n_cols);
+ for (std::size_t i = 0; i < header.n_cols; i++)
{
- if (!pix_vis[i*imgWidth+j])
- continue;
- auto val(static_cast<T>(pix_val[i * (imgWidth + 1) + j]));
- if (elem_type == MeshElemType::TRIANGLE)
+ std::size_t idx(i * header.n_rows + layer_idx);
+ for (std::size_t j = 0; j < header.n_rows; j++)
{
+ auto val(static_cast<T>(img[idx + j]));
prop_vec.push_back(val);
- prop_vec.push_back(val);
+ if (elem_type == MeshElemType::TRIANGLE || elem_type == MeshElemType::PRISM)
+ prop_vec.push_back(val); // because each pixel is represented by two cells
}
- else if (elem_type == MeshElemType::QUAD)
- prop_vec.push_back(val);
}
+ }
}
-
- static double getExistingValue(const double* img, std::size_t length);
};
} // end namespace MeshLib
diff --git a/MeshLib/MeshSearch/ElementSearch.h b/MeshLib/MeshSearch/ElementSearch.h
index d18d9f10d5a..d0aed4843c0 100644
--- a/MeshLib/MeshSearch/ElementSearch.h
+++ b/MeshLib/MeshSearch/ElementSearch.h
@@ -45,17 +45,15 @@ public:
template <typename PROPERTY_TYPE>
std::size_t searchByPropertyValue(
PROPERTY_TYPE const property_value,
- std::string const& property_name = "MaterialIDs")
+ std::string const& property_name)
{
- if (!_mesh.getProperties().existsPropertyVector<PROPERTY_TYPE>(
- property_name))
+ if (!_mesh.getProperties().existsPropertyVector<PROPERTY_TYPE>(property_name))
{
WARN("Property \"%s\" not found in mesh.", property_name.c_str());
return 0;
}
auto const* const pv =
- _mesh.getProperties().getPropertyVector<PROPERTY_TYPE>(
- property_name);
+ _mesh.getProperties().getPropertyVector<PROPERTY_TYPE>(property_name);
if (pv->getMeshItemType() != MeshLib::MeshItemType::Cell)
{
@@ -65,7 +63,8 @@ public:
}
std::vector<std::size_t> matchedIDs;
- for (std::size_t i(0); i < pv->getNumberOfTuples(); ++i) {
+ for (std::size_t i(0); i < pv->getNumberOfTuples(); ++i)
+ {
if ((*pv)[i] == property_value)
matchedIDs.push_back(i);
}
diff --git a/Tests/CMakeLists.txt b/Tests/CMakeLists.txt
index fb6b261ecdf..a623a41f84a 100644
--- a/Tests/CMakeLists.txt
+++ b/Tests/CMakeLists.txt
@@ -63,7 +63,7 @@ if(OGS_INSITU)
endif()
if(Qt5XmlPatterns_FOUND)
- target_link_libraries(testrunner Qt5::Core Qt5::Gui Qt5::Xml Qt5::Network)
+ target_link_libraries(testrunner Qt5::Core Qt5::Gui Qt5::Xml Qt5::Network VtkVis vtkGUISupportQt QtDataView)
if(CMAKE_CROSSCOMPILING)
target_link_libraries(testrunner
${QT_XML_DEPS_LIBRARIES}
diff --git a/Tests/MeshLib/TestMeshValidation.cpp b/Tests/MeshLib/TestMeshValidation.cpp
index acbf18911eb..4be198a0ffa 100644
--- a/Tests/MeshLib/TestMeshValidation.cpp
+++ b/Tests/MeshLib/TestMeshValidation.cpp
@@ -46,8 +46,8 @@ detectHoles(MeshLib::Mesh const& mesh,
TEST(MeshValidation, DetectHolesTri)
{
std::array<double, 12> pix = {{0,0.1,0.2,0.1,0,0,0.1,0,0,0,-0.1,0}};
- GeoLib::RasterHeader const header =
- {4,3,MathLib::Point3d(std::array<double,3>{{0,0,0}}),1,-9999};
+ GeoLib::RasterHeader const header = {
+ 4, 3, 1, MathLib::Point3d(std::array<double, 3>{{0, 0, 0}}), 1, -9999};
GeoLib::Raster const raster(header ,pix.begin(), pix.end());
std::unique_ptr<MeshLib::Mesh> mesh (MeshLib::RasterToMesh::convert(
raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION));
diff --git a/Tests/MeshLib/TestNodeSearch.cpp b/Tests/MeshLib/TestNodeSearch.cpp
index 923e6344a92..7cb79a3557a 100644
--- a/Tests/MeshLib/TestNodeSearch.cpp
+++ b/Tests/MeshLib/TestNodeSearch.cpp
@@ -23,8 +23,8 @@
TEST(NodeSearch, UnusedNodes)
{
std::array<double, 12> pix = {{0,0.1,0.2,0.1,0,0,0.1,0,0,0,-0.1,0}};
- GeoLib::RasterHeader const header =
- {4,3,MathLib::Point3d(std::array<double,3>{{0,0,0}}),1,-9999};
+ GeoLib::RasterHeader const header = {
+ 4, 3, 1, MathLib::Point3d(std::array<double, 3>{{0, 0, 0}}), 1, -9999};
GeoLib::Raster const raster(header,pix.begin(), pix.end());
MeshLib::Mesh* mesh (MeshLib::RasterToMesh::convert(raster, MeshLib::MeshElemType::TRIANGLE, MeshLib::UseIntensityAs::ELEVATION));
MeshLib::NodeSearch ns(*mesh);
diff --git a/Tests/MeshLib/TestRasterToMesh.cpp b/Tests/MeshLib/TestRasterToMesh.cpp
new file mode 100644
index 00000000000..2fac62f2b8a
--- /dev/null
+++ b/Tests/MeshLib/TestRasterToMesh.cpp
@@ -0,0 +1,304 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2017, 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 <cstdio>
+#include <memory>
+
+#include "gtest/gtest.h"
+
+#include "BaseLib/BuildInfo.h"
+
+#include "Applications/FileIO/AsciiRasterInterface.h"
+#include "GeoLib/Raster.h"
+#include "MeshLib/Mesh.h"
+#include "MeshLib/MeshGenerators/RasterToMesh.h"
+#include "MeshLib/MeshInformation.h"
+#include "MeshLib/Node.h"
+#include "MeshLib/IO/VtkIO/VtuInterface.h"
+
+#ifdef OGS_BUILD_GUI
+#include "Applications/DataExplorer/VtkVis/VtkGeoImageSource.h"
+#include "Applications/DataExplorer/VtkVis/VtkRaster.h"
+
+#include <vtkImageData.h>
+#endif
+
+class RasterToMeshTest : public ::testing::Test
+{
+public:
+ RasterToMeshTest()
+ : _file_name(BaseLib::BuildInfo::data_path + "/MeshLib/testraster_selke.asc"),
+ _mesh_name(BaseLib::BuildInfo::data_path + "/MeshLib/testraster_selke.vtu")
+ {
+ _raster.reset(FileIO::AsciiRasterInterface::readRaster(_file_name));
+ }
+
+ ~RasterToMeshTest() { std::remove(_mesh_name.c_str()); }
+
+protected:
+ std::size_t const _n_pix = 542;
+ std::size_t const _n_nodes = 626;
+ double _spacing = 1000;
+ std::string const _file_name;
+ std::string const _mesh_name;
+ std::unique_ptr<GeoLib::Raster> _raster;
+};
+
+TEST_F(RasterToMeshTest, convertRasterToTriMeshElevation)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::TRIANGLE,
+ MeshLib::UseIntensityAs::ELEVATION, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ MeshLib::IO::VtuInterface vtkio(mesh.get(), 0, false);
+ std::string name(BaseLib::BuildInfo::data_path +
+ "/MeshLib/testraster_selke.vtu");
+ vtkio.writeToFile(name);
+
+ std::cout << mesh->getNodes().size() << ", " << mesh->getNumberOfNodes()
+ << std::endl;
+ ASSERT_EQ(_n_nodes, mesh->getNodes().size());
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_TRUE(names.empty());
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(2 * _n_pix, n_types[1]);
+
+ GeoLib::AABB aabb = MeshLib::MeshInformation::getBoundingBox(*mesh);
+ ASSERT_NEAR(aabb.getMinPoint()[2], 0,
+ std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(aabb.getMaxPoint()[2], 0.07,
+ std::numeric_limits<double>::epsilon());
+}
+
+TEST_F(RasterToMeshTest, convertRasterToQuadMeshElevation)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::QUAD, MeshLib::UseIntensityAs::ELEVATION, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names = mesh->getProperties().getPropertyVectorNames();
+ ASSERT_TRUE(names.empty());
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(_n_pix, n_types[2]);
+
+ GeoLib::AABB aabb = MeshLib::MeshInformation::getBoundingBox(*mesh);
+ ASSERT_NEAR(aabb.getMinPoint()[2], 0,
+ std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(aabb.getMaxPoint()[2], 0.07,
+ std::numeric_limits<double>::epsilon());
+}
+
+TEST_F(RasterToMeshTest, convertRasterTo3DMeshElevation)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::PRISM,
+ MeshLib::UseIntensityAs::ELEVATION, "test"));
+ ASSERT_TRUE(mesh == nullptr);
+
+ std::unique_ptr<MeshLib::Mesh> mesh2(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::HEXAHEDRON,
+ MeshLib::UseIntensityAs::ELEVATION, "test"));
+ ASSERT_TRUE(mesh2 == nullptr);
+}
+
+TEST_F(RasterToMeshTest, convertRasterToTriMeshValue)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::TRIANGLE,
+ MeshLib::UseIntensityAs::DATAVECTOR, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_EQ(1, names.size());
+
+ MeshLib::PropertyVector<double>* prop =
+ mesh->getProperties().getPropertyVector<double>("test");
+ ASSERT_EQ(2 * _n_pix, prop->size());
+
+ std::pair<double, double> const& bounds =
+ MeshLib::MeshInformation::getValueBounds<double>(*mesh, "test");
+ ASSERT_NEAR(0, bounds.first, std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(0.07, bounds.second, std::numeric_limits<double>::epsilon());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_NEAR(0, (*n)[2], std::numeric_limits<double>::epsilon());
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(2 * _n_pix, n_types[1]);
+}
+
+TEST_F(RasterToMeshTest, convertRasterToQuadMeshValue)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::QUAD,
+ MeshLib::UseIntensityAs::DATAVECTOR, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_EQ(1, names.size());
+
+ MeshLib::PropertyVector<double>* prop =
+ mesh->getProperties().getPropertyVector<double>("test");
+ ASSERT_EQ(_n_pix, prop->size());
+
+ std::pair<double, double> const& bounds =
+ MeshLib::MeshInformation::getValueBounds<double>(*mesh, "test");
+ ASSERT_NEAR(0, bounds.first, std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(0.07, bounds.second, std::numeric_limits<double>::epsilon());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_TRUE((*n)[2] == 0);
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(_n_pix, n_types[2]);
+}
+
+TEST_F(RasterToMeshTest, convertRasterToPrismMeshValue)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::PRISM,
+ MeshLib::UseIntensityAs::DATAVECTOR, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(2 * _n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_EQ(1, names.size());
+
+ MeshLib::PropertyVector<double>* prop =
+ mesh->getProperties().getPropertyVector<double>("test");
+ ASSERT_EQ(2 * _n_pix, prop->size());
+
+ std::pair<double, double> const& bounds =
+ MeshLib::MeshInformation::getValueBounds<double>(*mesh, "test");
+ ASSERT_NEAR(0, bounds.first, std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(0.07, bounds.second, std::numeric_limits<double>::epsilon());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_TRUE(((*n)[2] == 0) || ((*n)[2] == _spacing));
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(2 * _n_pix, n_types[6]);
+}
+
+TEST_F(RasterToMeshTest, convertRasterToHexMeshValue)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::HEXAHEDRON,
+ MeshLib::UseIntensityAs::MATERIALS, "MaterialIDs"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(2 * _n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_EQ(1, names.size());
+
+ MeshLib::PropertyVector<int>* prop =
+ mesh->getProperties().getPropertyVector<int>("MaterialIDs");
+ ASSERT_EQ(_n_pix, prop->size());
+
+ std::pair<int, int> const& bounds =
+ MeshLib::MeshInformation::getValueBounds<int>(*mesh, "MaterialIDs");
+ ASSERT_NEAR(0, bounds.first, std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(0, bounds.second, std::numeric_limits<double>::epsilon());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_TRUE(((*n)[2] == 0) || ((*n)[2] == _spacing));
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(_n_pix, n_types[4]);
+}
+
+TEST_F(RasterToMeshTest, convertRasterToQuadMeshNone)
+{
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ *_raster, MeshLib::MeshElemType::QUAD,
+ MeshLib::UseIntensityAs::NONE, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_TRUE(names.empty());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_TRUE((*n)[2] == 0);
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(_n_pix, n_types[2]);
+}
+
+#ifdef OGS_BUILD_GUI
+TEST_F(RasterToMeshTest, vtkImage)
+{
+ double x0, y0, spacing;
+ vtkImageAlgorithm* raster =
+ VtkRaster::loadImage(_file_name, x0, y0, spacing);
+ double origin[3];
+ raster->GetOutput()->GetOrigin(origin);
+
+ std::unique_ptr<MeshLib::Mesh> mesh(MeshLib::RasterToMesh::convert(
+ raster->GetOutput(), origin, spacing, MeshLib::MeshElemType::TRIANGLE,
+ MeshLib::UseIntensityAs::DATAVECTOR, "test"));
+ ASSERT_TRUE(mesh != nullptr);
+
+ ASSERT_EQ(_n_nodes, mesh->getNumberOfBaseNodes());
+
+ std::vector<std::string> names =
+ mesh->getProperties().getPropertyVectorNames();
+ ASSERT_EQ(1, names.size());
+
+ MeshLib::PropertyVector<double>* prop =
+ mesh->getProperties().getPropertyVector<double>("test");
+ ASSERT_EQ(2 * _n_pix, prop->size());
+
+ std::pair<double, double> const& bounds =
+ MeshLib::MeshInformation::getValueBounds<double>(*mesh, "test");
+ ASSERT_NEAR(0, bounds.first, std::numeric_limits<double>::epsilon());
+ ASSERT_NEAR(0.07, bounds.second, std::numeric_limits<float>::epsilon());
+
+ std::vector<MeshLib::Node*> const& nodes = mesh->getNodes();
+ for (MeshLib::Node* n : nodes)
+ ASSERT_TRUE((*n)[2] == 0);
+
+ std::array<unsigned, 7> n_types =
+ MeshLib::MeshInformation::getNumberOfElementTypes(*mesh);
+ ASSERT_EQ(2 * _n_pix, n_types[1]);
+}
+#endif
diff --git a/Tests/lfs-data/MeshLib/testraster_selke.asc b/Tests/lfs-data/MeshLib/testraster_selke.asc
new file mode 100644
index 00000000000..52713da098b
--- /dev/null
+++ b/Tests/lfs-data/MeshLib/testraster_selke.asc
@@ -0,0 +1,39 @@
+ncols 43
+nrows 33
+xllcorner 4423656.0
+yllcorner 5716944.0
+cellsize 1000
+NODATA_value -9999
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.02 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.02 0.01 0.01 0.01 0.01 0.01 0.02 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.02 0.02 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.00 0.01 0.00 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.00 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.00 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.01 0.01 0.01 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.00 0.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.00 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.00 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.01 0.03 -9999.00 -9999.00 0.02 -9999.00 -9999.00 -9999.00 0.02 0.02 0.02 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.03 0.03 0.04 0.04 0.02 0.03 0.02 0.02 0.02 0.02 0.02 0.01 0.01 0.02 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.01 0.02 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.04 0.04 0.04 0.04 0.04 0.03 0.03 0.02 0.02 0.02 0.02 0.02 0.02 0.02 0.02 0.01 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.05 0.05 0.04 0.04 0.04 0.04 0.03 0.03 0.03 0.02 0.02 0.02 0.02 0.02 0.02 0.01 0.02 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 0.06 0.05 0.05 0.05 0.05 0.06 0.05 0.04 0.05 0.04 0.03 0.03 0.03 0.03 0.03 0.02 0.02 0.02 0.02 0.01 0.01 0.01 0.01 0.01 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ 0.06 0.06 0.05 0.05 0.05 0.06 0.05 0.05 0.04 0.05 0.04 0.04 0.03 0.03 0.02 0.01 0.02 0.02 0.01 0.02 0.02 0.02 0.02 0.02 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ 0.06 0.06 0.05 0.05 0.05 0.03 0.04 0.04 0.04 0.04 0.02 0.04 0.03 0.03 0.02 0.02 0.01 0.00 0.02 0.02 0.02 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ 0.06 0.06 0.06 0.06 0.06 0.05 0.05 0.04 0.04 0.04 0.04 0.03 0.03 0.02 0.02 0.03 0.03 0.03 0.03 0.02 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ 0.06 0.06 0.06 0.05 0.06 0.05 0.05 0.04 0.04 0.05 0.05 0.04 0.03 0.02 0.03 0.03 0.03 0.03 0.03 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ 0.07 0.06 0.07 0.06 0.06 0.05 0.05 0.04 0.04 0.04 0.02 0.04 0.03 0.03 0.03 0.04 0.04 0.02 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 0.06 0.06 0.05 0.05 0.05 0.04 0.04 0.05 0.04 0.04 0.03 0.03 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.05 0.05 0.04 0.04 0.05 0.04 0.04 0.04 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
+ -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 0.04 0.05 0.05 0.05 0.04 0.04 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00 -9999.00
--
GitLab