OGS
VtkVisImageItem.cpp
Go to the documentation of this file.
1
15// ** INCLUDES **
16#include "VtkVisImageItem.h"
17
18#include <vtkActor.h>
19#include <vtkDataArray.h>
20#include <vtkDataSetMapper.h>
21#include <vtkImageAlgorithm.h>
22#include <vtkImageChangeInformation.h>
23#include <vtkImageData.h>
24#include <vtkImageShiftScale.h>
25#include <vtkPointData.h>
26#include <vtkRenderer.h>
27#include <vtkSmartPointer.h>
28
29#include <QFileDialog>
30
32#include "Base/OGSError.h"
33#include "BaseLib/FileTools.h"
34#include "BaseLib/Logging.h"
36#include "GeoLib/Raster.h"
38#include "VtkGeoImageSource.h"
39
40// export
41#include <vtkImageActor.h>
42#include <vtkXMLImageDataWriter.h>
43
45 vtkAlgorithm* algorithm, TreeItem* parentItem,
46 const QList<QVariant> data /*= QList<QVariant>()*/)
47 : VtkVisPipelineItem(algorithm, parentItem, data), _transformFilter(nullptr)
48{
49}
50
52 VtkCompositeFilter* compositeFilter, TreeItem* parentItem,
53 const QList<QVariant> data /*= QList<QVariant>()*/)
54 : VtkVisPipelineItem(compositeFilter, parentItem, data),
55 _transformFilter(nullptr)
56{
57}
58
63
64void VtkVisImageItem::Initialize(vtkRenderer* renderer)
65{
66 auto* image = dynamic_cast<vtkImageAlgorithm*>(_algorithm);
67 if (!image)
68 {
69 OGSError::box("Cast to image failed.");
70 return;
71 }
72 image->Update();
73
74 double origin[3];
75 double spacing[3];
76 double range[2];
77 vtkImageData* image_data = image->GetOutput();
78 if (!image_data)
79 {
80 OGSError::box("GetOutput() - Image could not be initialized.");
81 return;
82 }
83 image_data->GetOrigin(origin);
84 image_data->GetSpacing(spacing);
85 vtkPointData* point_data = image_data->GetPointData();
86 if (point_data)
87 {
88 auto* scalars = point_data->GetScalars();
89 if (scalars)
90 {
91 scalars->GetRange(range);
92 }
93 }
94 vtkImageShiftScale* scale = vtkImageShiftScale::New();
95 scale->SetOutputScalarTypeToUnsignedChar();
96 scale->SetInputConnection(image->GetOutputPort());
97 scale->SetShift(-range[0]);
98 scale->SetScale(255.0 / (range[1] - range[0]));
99
100 _transformFilter = vtkImageChangeInformation::New();
101 _transformFilter->SetInputConnection(scale->GetOutputPort());
102 // double origin[3];
103 // img->getOrigin(origin);
104 // double spacing = img->getSpacing();
105 //_transformFilter->SetOutputOrigin(origin2);
106 //_transformFilter->SetOutputSpacing(spacing2);
107 _transformFilter->Update();
108
109 _renderer = renderer;
110
111 // Use a special vtkImageActor instead of vtkActor
112 vtkImageActor* imageActor = vtkImageActor::New();
113 imageActor->SetInputData(_transformFilter->GetOutput());
114 _actor = imageActor;
115 _renderer->AddActor(_actor);
116
117 // Set pre-set properties
118 auto* vtkProps = dynamic_cast<VtkAlgorithmProperties*>(_algorithm);
119 if (vtkProps)
120 {
121 setVtkProperties(vtkProps);
122 }
123
124 auto* parentItem = dynamic_cast<VtkVisPipelineItem*>(this->parentItem());
125 while (parentItem)
126 {
127 auto* parentProps =
128 dynamic_cast<VtkAlgorithmProperties*>(parentItem->algorithm());
129 if (parentProps)
130 {
131 auto* newProps = new VtkAlgorithmProperties();
132 newProps->SetScalarVisibility(parentProps->GetScalarVisibility());
133 newProps->SetTexture(parentProps->GetTexture());
134 setVtkProperties(newProps);
135 vtkProps = newProps;
136 parentItem = nullptr;
137 }
138 else
139 {
140 parentItem =
141 dynamic_cast<VtkVisPipelineItem*>(parentItem->parentItem());
142 }
143 }
144
145 // Set active scalar to the desired one from VtkAlgorithmProperties
146 // or to match those of the parent.
147 if (vtkProps)
148 {
149 if (vtkProps->GetActiveAttribute().length() > 0)
150 {
151 this->SetActiveAttribute(vtkProps->GetActiveAttribute());
152 }
153 else
154 {
155 auto* visParentItem =
156 dynamic_cast<VtkVisPipelineItem*>(this->parentItem());
157 if (visParentItem)
158 {
159 this->SetActiveAttribute(visParentItem->GetActiveAttribute());
160 }
161 if (vtkProps->GetTexture() != nullptr)
162 {
163 this->SetActiveAttribute("Solid Color");
164 }
165 }
166 }
167}
168
170{
171 // todo implementation
172 (void)vtkProps;
173}
174
175int VtkVisImageItem::callVTKWriter(vtkAlgorithm* algorithm,
176 const std::string& filename) const
177{
178 std::string file_name_cpy(filename);
179 auto* algID = dynamic_cast<vtkImageAlgorithm*>(algorithm);
180 if (algID)
181 {
182 vtkSmartPointer<vtkXMLImageDataWriter> iWriter =
183 vtkSmartPointer<vtkXMLImageDataWriter>::New();
184 iWriter->SetInputData(algID->GetOutputDataObject(0));
185 if (BaseLib::getFileExtension(filename) != ".vti")
186 {
187 file_name_cpy.append(".vti");
188 }
189 iWriter->SetFileName(file_name_cpy.c_str());
190 return iWriter->Write();
191 }
192 ERR("VtkVisPipelineItem::writeToFile() - Unknown data type.");
193 return 0;
194}
195
196void VtkVisImageItem::setTranslation(double x, double y, double z) const
197{
198 _transformFilter->SetOriginTranslation(x, y, z);
199}
200
202{
203 return _transformFilter;
204}
205
207{
208 vtkSmartPointer<VtkGeoImageSource> imageSource =
209 VtkGeoImageSource::SafeDownCast(_algorithm);
210
211 auto const raster = VtkGeoImageSource::convertToRaster(imageSource);
212
213 if (!raster)
214 {
215 OGSError::box("Image could not be converted to a raster.");
216 return false;
217 }
218
219 QFileInfo const info(this->data(0).toString());
220 QString const rastername = info.completeBaseName() + ".asc";
221 QString const filetype("ESRI ASCII raster file (*.asc)");
222 QString const filename = QFileDialog::getSaveFileName(
223 nullptr, "Save raster as",
224 LastSavedFileDirectory::getDir() + rastername, filetype);
227 filename.toStdString());
228 return true;
229}
Definition of the AsciiRasterInterface class.
Filename manipulation routines.
Manages the last directory used for saving a file.
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:45
Definition of the OGSError class.
Definition of the GeoLib::Raster class.
Definition of the VtkAlgorithmProperties class.
Definition of the VtkGeoImageSource class.
Definition of the VtkVisImageItem class.
static void writeRasterAsASC(GeoLib::Raster const &raster, std::string const &file_name)
Writes an Esri asc-file.
static void setDir(const QString &path)
Sets the directory last used for saving a file.
static const QString getDir()
Returns the directory last used for saving a file.
static void box(const QString &e)
Definition OGSError.cpp:23
Objects nodes for the TreeModel.
Definition TreeItem.h:28
TreeItem * parentItem() const
Definition TreeItem.cpp:115
Contains properties for the visualization of objects as VtkVisPipelineItems.
Is used to combine several filter in one VtkVisPipelineItem. You can use vtk filter and custom filter...
static std::optional< GeoLib::Raster > convertToRaster(VtkGeoImageSource *const source)
vtkImageChangeInformation * _transformFilter
VtkVisImageItem(vtkAlgorithm *algorithm, TreeItem *parentItem, const QList< QVariant > data=QList< QVariant >())
Constructor for a source/filter object.
bool writeAsRaster()
Allows writing this item's source object as an ASCII raster file.
void Initialize(vtkRenderer *renderer) override
Initializes vtkMapper and vtkActor necessary for visualization of the item and sets the item's proper...
vtkAlgorithm * transformFilter() const override
int callVTKWriter(vtkAlgorithm *algorithm, const std::string &filename) const override
Selects the appropriate VTK-Writer object and writes the object to a file with the given name.
~VtkVisImageItem() override
void setVtkProperties(VtkAlgorithmProperties *vtkProps)
void setTranslation(double x, double y, double z) const override
Translates the item in visualisation-space. This function is empty and needs to be implemented by der...
An item in the VtkVisPipeline containing a graphic object to be visualized.
vtkAlgorithm * algorithm() const
Returns the algorithm object.
vtkAlgorithm * _algorithm
QVariant data(int column) const override
virtual void SetActiveAttribute(const QString &str)
std::string getFileExtension(const std::string &path)