TTK
Loading...
Searching...
No Matches
ttkTextureMapFromField.cpp
Go to the documentation of this file.
1#include <vtkDataArray.h>
2#include <vtkDataSet.h>
3#include <vtkFloatArray.h>
4#include <vtkInformation.h>
5#include <vtkObjectFactory.h>
6#include <vtkPointData.h>
7
8#include <array>
9
11
13
15 this->setDebugMsgPrefix("TextureMapFromField");
16 this->SetNumberOfInputPorts(1);
17 this->SetNumberOfOutputPorts(1);
18}
19
21 vtkInformation *info) {
22 if(port == 0) {
23 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
24 return 1;
25 }
26 return 0;
27};
29 vtkInformation *info) {
30 if(port == 0) {
32 return 1;
33 }
34 return 0;
35};
36
38 vtkInformationVector **inputVector,
39 vtkInformationVector *outputVector) {
40
41 const auto input = vtkDataSet::GetData(inputVector[0]);
42 auto output = vtkDataSet::GetData(outputVector);
43
44 ttk::Timer t;
45
46 output->ShallowCopy(input);
47
48 const auto inputScalarFieldU = this->GetInputArrayToProcess(0, input);
49 const auto inputScalarFieldV = this->GetInputArrayToProcess(1, input);
50
51 if(inputScalarFieldU == nullptr || inputScalarFieldV == nullptr) {
52 return -1;
53 }
54
55 vtkNew<vtkFloatArray> textureCoordinates{};
56 textureCoordinates->SetNumberOfComponents(2);
57 textureCoordinates->SetName("UV coordinates from field");
58
59 if(textureCoordinates->GetNumberOfTuples() != output->GetNumberOfPoints()) {
60 textureCoordinates->SetNumberOfTuples(output->GetNumberOfPoints());
61 }
62
63 double uRange[2], vRange[2];
64 inputScalarFieldU->GetRange(uRange);
65 inputScalarFieldV->GetRange(vRange);
66
67 std::vector<std::array<double, 2>> coordinates(threadNumber_);
68
69#ifdef TTK_ENABLE_OPENMP
70#pragma omp parallel for num_threads(threadNumber_)
71#endif
72 for(int i = 0; i < output->GetNumberOfPoints(); i++) {
73
74 ttk::ThreadId threadId = 0;
75
76#ifdef TTK_ENABLE_OPENMP
77 threadId = omp_get_thread_num();
78#endif
79
80 coordinates[threadId][0] = coordinates[threadId][1] = 0;
81
82 if(!OnlyVComponent) {
83 inputScalarFieldU->GetTuple(i, &(coordinates[threadId][0]));
84 if(!RepeatUTexture) {
85 coordinates[threadId][0]
86 = (coordinates[threadId][0] - uRange[0]) / (uRange[1] - uRange[0]);
87 }
88 }
89
90 if(!OnlyUComponent) {
91 inputScalarFieldV->GetTuple(i, &(coordinates[threadId][1]));
92 if(!RepeatVTexture) {
93 coordinates[threadId][1]
94 = (coordinates[threadId][1] - vRange[0]) / (vRange[1] - vRange[0]);
95 }
96 }
97
98 textureCoordinates->SetTuple(i, coordinates[threadId].data());
99 }
100
101 output->GetPointData()->SetTCoords(textureCoordinates);
102
103 this->printMsg(
104 "Computed texture map", 1.0, t.getElapsedTime(), this->threadNumber_);
105
106 return 1;
107}
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
Definition BaseClass.h:47
static vtkInformationIntegerKey * SAME_DATA_TYPE_AS_INPUT_PORT()
TTK VTK-filter which generates a texture map from one or two point data scalar fields.
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
int FillOutputPortInformation(int port, vtkInformation *info) override
int FillInputPortInformation(int port, vtkInformation *info) override
void setDebugMsgPrefix(const std::string &prefix)
Definition Debug.h:364
double getElapsedTime()
Definition Timer.h:15
int ThreadId
Identifier type for threads (i.e. with OpenMP).
Definition DataTypes.h:26
vtkStandardNewMacro(ttkTextureMapFromField)
printMsg(debug::output::BOLD+" | | | | | . \\ | | (__| | / __/| |_| / __/|__ _|"+debug::output::ENDCOLOR, debug::Priority::PERFORMANCE, debug::LineMode::NEW, stream)