6#include <vtkInformation.h>
7#include <vtkInformationVector.h>
10#include <vtkCellData.h>
11#include <vtkDataSet.h>
12#include <vtkUnstructuredGrid.h>
14#include <vtkFloatArray.h>
15#include <vtkIdTypeArray.h>
24 vtkInformation *info) {
27 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(),
"vtkDataSet");
32 vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(),
"vtkUnstructuredGrid");
39 vtkInformation *info) {
43 info->Set(vtkDataObject::DATA_TYPE_NAME(),
"vtkUnstructuredGrid");
52 vtkInformationVector **iVec,
53 vtkInformationVector *oVec) {
55 _outFld = vtkUnstructuredGrid::GetData(oVec, 0);
56 _outPts = vtkUnstructuredGrid::GetData(oVec, 1);
61 vtkUnstructuredGrid::GetData(iVec[2])))
78 auto scalars = GetInputArrayToProcess(0, dataset);
82 const double radius = ui_spherical ? -1. : 0.;
87 printErr(
"super->setInputField failed with code "
88 + std::to_string(errorCode));
92 _triangTypeCode = triangulation->
getType();
93 _scalarTypeCode = scalars->GetDataType();
94 _scalarsName = scalars->GetName();
96 std::ostringstream stream;
97 stream <<
"Scalar type: " << scalars->GetDataTypeAsString() <<
" (code "
98 << _scalarTypeCode <<
")";
107 vtkUnstructuredGrid *arcs) {
110 auto points = nodes->GetPoints();
111 if(points->GetDataType() != VTK_FLOAT) {
112 printErr(
"The point coordinates must be of type float");
118 auto pData = nodes->GetPointData();
119 auto scalarBuf = getBuffer<float>(pData,
"Scalar", VTK_FLOAT,
"float");
120 auto codeBuf = getBuffer<int>(pData,
"CriticalType", VTK_INT,
"int");
121 if(!scalarBuf || !codeBuf)
127 auto cells = arcs->GetCells();
128 const auto maxNvPerC = cells->GetMaxCellSize();
131 "The points must come in pairs but there is at least one cell with "
132 + std::to_string(maxNvPerC) +
" points");
138 auto cData = arcs->GetCellData();
139 auto c2p = getBuffer<int>(cData,
"upNodeId", VTK_INT,
"int");
140 auto c2q = getBuffer<int>(cData,
"downNodeId", VTK_INT,
"int");
146 static constexpr int minCode = 0;
147 static constexpr int maxCode = 3;
148 const double sadFac = ui_extension * 0.01;
149 const double extFac = 1 - sadFac;
156 auto addPoint = [
this, coords, scalarBuf](
int p,
float isoval,
int code) {
157 const auto point = &coords[p * 3];
158 _coords.push_back(point[0]);
159 _coords.push_back(point[1]);
160 _coords.push_back(point[2]);
161 _scalars.push_back(scalarBuf[p]);
162 _isovals.push_back(isoval);
163 _flags.push_back(code == minCode ? 0 : 1);
166 const vtkIdType nc = arcs->GetNumberOfCells();
167 for(vtkIdType c = 0; c < nc; ++c) {
168 const auto p = c2p[c];
169 const auto q = c2q[c];
170 const auto pCode = codeBuf[p];
171 const auto qCode = codeBuf[q];
173 const bool pIsSad = pCode != minCode && pCode != maxCode;
174 const bool qIsSad = qCode != minCode && qCode != maxCode;
175 if(pIsSad || qIsSad) {
179 const auto ext = pIsSad ? q : p;
180 const auto sad = pIsSad ? p : q;
181 const float isoval = scalarBuf[ext] * extFac + scalarBuf[sad] * sadFac;
182 addPoint(ext, isoval, codeBuf[ext]);
184 printWrn(
"Arc " + std::to_string(c) +
" joins a minimum and a maximum");
185 const auto pVal = scalarBuf[p];
186 const auto qVal = scalarBuf[q];
187 const auto cVal = (pVal + qVal) / 2;
188 addPoint(p, pVal * extFac + cVal * sadFac, pCode);
189 addPoint(q, qVal * extFac + cVal * sadFac, qCode);
193 auto np = _scalars.size();
195 _coords.data(), _scalars.data(), _isovals.data(), _flags.data(), np);
197 printErr(
"setInputPoints failed with code " + std::to_string(errorCode));
207 switch(_scalarTypeCode) {
208 vtkTemplateMacro((errorCode = this->execute<VTK_TT>()));
215 printErr(
"super->execute failed with code " + std::to_string(errorCode));
230 std::vector<int> ctypes(nc);
232 vtkIdType cinfoCounter = 0;
235 assert(nvOfCell >= 2 && nvOfCell <= 3);
236 ctypes[c] = nvOfCell == 2 ? VTK_LINE : VTK_TRIANGLE;
237 cinfoCounter += nvOfCell + 1;
244 cells->SetCells(nc, cinfoArr);
245 _outFld->SetCells(ctypes.data(), cells);
250 printErr(
"The API has changed! We have expected the default "
251 "coordinate type to be float");
257 coordArr->SetNumberOfComponents(3);
260 points->SetData(coordArr);
261 _outFld->SetPoints(points);
263 auto scalarArr = vtkFloatArray::New();
266 scalarArr->SetName(_scalarsName);
267 _outFld->GetPointData()->AddArray(scalarArr);
269 auto flagArr = vtkIntArray::New();
272 flagArr->SetName(
"isMax");
273 _outFld->GetPointData()->AddArray(flagArr);
279 coordArr->SetNumberOfComponents(3);
282 points->SetData(coordArr);
283 _outPts->SetPoints(points);
285 scalarArr = vtkFloatArray::New();
288 scalarArr->SetName(_scalarsName);
289 _outPts->GetPointData()->AddArray(scalarArr);
291 flagArr = vtkIntArray::New();
294 flagArr->SetName(
"isMax");
295 _outPts->GetPointData()->AddArray(flagArr);
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
ttk::Triangulation * GetTriangulation(vtkDataSet *dataSet)
TTK VTK-filter that wraps the contourAroundPoint processing package.
bool postprocess()
Assemble the output object from the results of the TTK module.
int RequestData(vtkInformation *request, vtkInformationVector **iVec, vtkInformationVector *oVec) override
bool preprocessPts(vtkUnstructuredGrid *nodes, vtkUnstructuredGrid *arcs)
int FillOutputPortInformation(int port, vtkInformation *info) override
bool preprocessFld(vtkDataSet *dataset)
int FillInputPortInformation(int port, vtkInformation *info) override
static void * GetVoidPointer(vtkDataArray *array, vtkIdType start=0)
static void SetVoidArray(vtkDataArray *array, void *data, vtkIdType size, int save)
int setInputField(triangulationType *triangulation, void *scalars, double sizeFilter, double radius=0.)
std::vector< float > _outCentroidsCoords
std::vector< float > _outContoursScalars
std::vector< float > _outContoursCoords
std::vector< int > _outCentroidsFlags
std::vector< float > _outCentroidsScalars
std::vector< LongSimplexId > _outContoursCinfos
int setInputPoints(float *coords, float *scalars, float *isovals, int *flags, std::size_t np)
std::vector< int > _outContoursFlags
int printWrn(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
int printMsg(const std::string &msg, const debug::Priority &priority=debug::Priority::INFO, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cout) const
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Triangulation is a class that provides time and memory efficient traversal methods on triangulations ...
Triangulation::Type getType() const
int SimplexId
Identifier type for simplices of any dimension.
vtkStandardNewMacro(ttkContourAroundPoint)