www.mooseframework.org
ImageSampler.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 // MOOSE includes
11 #include "ImageSampler.h"
12 #include "MooseApp.h"
13 #include "ImageMesh.h"
14 
15 #include "libmesh/mesh_tools.h"
16 
19 {
20  // Define the general parameters
23 
24  params.addParam<Point>("origin", "Origin of the image (defaults to mesh origin)");
25  params.addParam<Point>("dimensions",
26  "x,y,z dimensions of the image (defaults to mesh dimensions)");
27  params.addParam<unsigned int>(
28  "component",
29  "The image RGB-component to return, leaving this blank will result in a greyscale value "
30  "for the image to be created. The component number is zero based, i.e. 0 returns the first "
31  "(RED) component of the image.");
32 
33  // Shift and Scale (application of these occurs prior to threshold)
34  params.addParam<double>("shift", 0, "Value to add to all pixels; occurs prior to scaling");
35  params.addParam<double>(
36  "scale", 1, "Multiplier to apply to all pixel values; occurs after shifting");
37  params.addParamNamesToGroup("shift scale", "Rescale");
38 
39  // Threshold parameters
40  params.addParam<double>("threshold", "The threshold value");
41  params.addParam<double>(
42  "upper_value", 1, "The value to set for data greater than the threshold value");
43  params.addParam<double>(
44  "lower_value", 0, "The value to set for data less than the threshold value");
45  params.addParamNamesToGroup("threshold upper_value lower_value", "Threshold");
46 
47  // Flip image
48  params.addParam<bool>("flip_x", false, "Flip the image along the x-axis");
49  params.addParam<bool>("flip_y", false, "Flip the image along the y-axis");
50  params.addParam<bool>("flip_z", false, "Flip the image along the z-axis");
51  params.addParamNamesToGroup("flip_x flip_y flip_z", "Flip");
52 
53  return params;
54 }
55 
57  : FileRangeBuilder(parameters),
58 #ifdef LIBMESH_HAVE_VTK
59  _data(nullptr),
60  _algorithm(nullptr),
61 #endif
62  _is_pars(parameters),
63  _is_console(
64  (parameters.getCheckedPointerParam<MooseApp *>("_moose_app"))->getOutputWarehouse()),
65  _flip({{_is_pars.get<bool>("flip_x"),
66  _is_pars.get<bool>("flip_y"),
67  _is_pars.get<bool>("flip_z")}})
68 
69 {
70 #ifndef LIBMESH_HAVE_VTK
71  // This should be impossible to reach, the registration of ImageSampler is also guarded with
72  // LIBMESH_HAVE_VTK
73  mooseError("libMesh must be configured with VTK enabled to utilize ImageSampler");
74 #endif
75 }
76 
77 void
79 {
80  // Don't warn that mesh or _is_pars are unused when VTK is not enabled.
83 
84 #ifdef LIBMESH_HAVE_VTK
85  // Get access to the Mesh object
86  BoundingBox bbox = MeshTools::create_bounding_box(mesh.getMesh());
87 
88  // Set the dimensions from the Mesh if not set by the User
89  if (_is_pars.isParamValid("dimensions"))
90  _physical_dims = _is_pars.get<Point>("dimensions");
91 
92  else
93  {
94  _physical_dims(0) = bbox.max()(0) - bbox.min()(0);
95 #if LIBMESH_DIM > 1
96  _physical_dims(1) = bbox.max()(1) - bbox.min()(1);
97 #endif
98 #if LIBMESH_DIM > 2
99  _physical_dims(2) = bbox.max()(2) - bbox.min()(2);
100 #endif
101  }
102 
103  // Set the origin from the Mesh if not set in the input file
104  if (_is_pars.isParamValid("origin"))
105  _origin = _is_pars.get<Point>("origin");
106  else
107  {
108  _origin(0) = bbox.min()(0);
109 #if LIBMESH_DIM > 1
110  _origin(1) = bbox.min()(1);
111 #endif
112 #if LIBMESH_DIM > 2
113  _origin(2) = bbox.min()(2);
114 #endif
115  }
116 
117  // An array of filenames, to be filled in
118  std::vector<std::string> filenames;
119 
120  // The file suffix, to be determined
121  std::string file_suffix;
122 
123  // Try to parse our own file range parameters. If that fails, then
124  // see if the associated Mesh is an ImageMesh and use its. If that
125  // also fails, then we have to throw an error...
126  //
127  // The parseFileRange method sets parameters, thus a writable reference to the InputParameters
128  // object must be obtained from the warehouse. Generally, this should be avoided, but
129  // this is a special case.
130  if (_status != 0)
131  {
132  // We don't have parameters, so see if we can get them from ImageMesh
133  ImageMesh * image_mesh = dynamic_cast<ImageMesh *>(&mesh);
134  if (!image_mesh)
135  mooseError("No file range parameters were provided and the Mesh is not an ImageMesh.");
136 
137  // Get the ImageMesh's parameters. This should work, otherwise
138  // errors would already have been thrown...
139  filenames = image_mesh->filenames();
140  file_suffix = image_mesh->fileSuffix();
141  }
142  else
143  {
144  // Use our own parameters (using 'this' b/c of conflicts with filenames the local variable)
145  filenames = this->filenames();
146  file_suffix = fileSuffix();
147  }
148 
149  // Storage for the file names
150  _files = vtkSmartPointer<vtkStringArray>::New();
151 
152  for (const auto & filename : filenames)
153  _files->InsertNextValue(filename);
154 
155  // Error if no files where located
156  if (_files->GetNumberOfValues() == 0)
157  mooseError("No image file(s) located");
158 
159  // Read the image stack. Hurray for VTK not using polymorphism in a
160  // smart way... we actually have to explicitly create the type of
161  // reader based on the file extension, using an if-statement...
162  if (file_suffix == "png")
163  _image = vtkSmartPointer<vtkPNGReader>::New();
164  else if (file_suffix == "tiff" || file_suffix == "tif")
165  _image = vtkSmartPointer<vtkTIFFReader>::New();
166  else
167  mooseError("Un-supported file type '", file_suffix, "'");
168 
169  // Now that _image is set up, actually read the images
170  // Indicate that data read has started
171  _is_console << "Reading image(s)..." << std::endl;
172 
173  // Extract the data
174  _image->SetFileNames(_files);
175  _image->Update();
176  _data = _image->GetOutput();
177  _algorithm = _image->GetOutputPort();
178 
179  // Set the image dimensions and voxel size member variable
180  int * dims = _data->GetDimensions();
181  for (unsigned int i = 0; i < 3; ++i)
182  {
183  _dims.push_back(dims[i]);
184  _voxel.push_back(_physical_dims(i) / _dims[i]);
185  }
186 
187  // Set the dimensions of the image and bounding box
188  _data->SetSpacing(_voxel[0], _voxel[1], _voxel[2]);
189  _data->SetOrigin(_origin(0), _origin(1), _origin(2));
190  _bounding_box.min() = _origin;
192 
193  // Indicate data read is completed
194  _is_console << " ...image read finished" << std::endl;
195 
196  // Set the component parameter
197  // If the parameter is not set then vtkMagnitude() will applied
198  if (_is_pars.isParamValid("component"))
199  {
200  unsigned int n = _data->GetNumberOfScalarComponents();
201  _component = _is_pars.get<unsigned int>("component");
202  if (_component >= n)
203  mooseError("'component' parameter must be empty or have a value of 0 to ", n - 1);
204  }
205  else
206  _component = 0;
207 
208  // Apply filters, the toggling on and off of each filter is handled internally
209  vtkMagnitude();
211  vtkThreshold();
212 #endif
213 }
214 
215 Real
216 ImageSampler::sample(const Point & p) const
217 {
218 #ifdef LIBMESH_HAVE_VTK
219 
220  // Do nothing if the point is outside of the image domain
221  if (!_bounding_box.contains_point(p))
222  return 0.0;
223 
224  // Determine pixel coordinates
225  std::vector<int> x(3, 0);
226  for (int i = 0; i < LIBMESH_DIM; ++i)
227  {
228  // Compute position, only if voxel size is greater than zero
229  if (_voxel[i] != 0)
230  {
231  x[i] = std::floor((p(i) - _origin(i)) / _voxel[i]);
232 
233  // If the point falls on the mesh extents the index needs to be decreased by one
234  if (x[i] == _dims[i])
235  x[i]--;
236 
237  // flip
238  if (_flip[i])
239  x[i] = _dims[i] - x[i] - 1;
240  }
241  }
242 
243  // Return the image data at the given point
244  return _data->GetScalarComponentAsDouble(x[0], x[1], x[2], _component);
245 
246 #else
247  libmesh_ignore(p); // avoid un-used parameter warnings
248  return 0.0;
249 #endif
250 }
251 
252 void
254 {
255 #ifdef LIBMESH_HAVE_VTK
256  // Do nothing if 'component' is set
257  if (_is_pars.isParamValid("component"))
258  return;
259 
260  // Apply the greyscale filtering
261  _magnitude_filter = vtkSmartPointer<vtkImageMagnitude>::New();
262  _magnitude_filter->SetInputConnection(_algorithm);
263  _magnitude_filter->Update();
264 
265  // Update the pointers
266  _data = _magnitude_filter->GetOutput();
267  _algorithm = _magnitude_filter->GetOutputPort();
268 #endif
269 }
270 
271 void
273 {
274 #ifdef LIBMESH_HAVE_VTK
275  // Capture the parameters
276  double shift = _is_pars.get<double>("shift");
277  double scale = _is_pars.get<double>("scale");
278 
279  // Do nothing if shift and scale are not set
280  if (shift == 0 && scale == 1)
281  return;
282 
283  // Perform the scaling and offset actions
284  _shift_scale_filter = vtkSmartPointer<vtkImageShiftScale>::New();
285  _shift_scale_filter->SetOutputScalarTypeToDouble();
286 
287  _shift_scale_filter->SetInputConnection(_algorithm);
288  _shift_scale_filter->SetShift(shift);
289  _shift_scale_filter->SetScale(scale);
290  _shift_scale_filter->Update();
291 
292  // Update the pointers
293  _data = _shift_scale_filter->GetOutput();
294  _algorithm = _shift_scale_filter->GetOutputPort();
295 #endif
296 }
297 
298 void
300 {
301 #ifdef LIBMESH_HAVE_VTK
302  // Do nothing if threshold not set
303  if (!_is_pars.isParamValid("threshold"))
304  return;
305 
306  // Error if both upper and lower are not set
307  if (!_is_pars.isParamValid("upper_value") || !_is_pars.isParamValid("lower_value"))
308  mooseError("When thresholding is applied, both the upper_value and lower_value parameters must "
309  "be set");
310 
311  // Create the thresholding object
312  _image_threshold = vtkSmartPointer<vtkImageThreshold>::New();
313 
314  // Set the data source
315  _image_threshold->SetInputConnection(_algorithm);
316 
317  // Setup the thresholding options
318  _image_threshold->ThresholdByUpper(_is_pars.get<Real>("threshold"));
319  _image_threshold->ReplaceInOn();
320  _image_threshold->SetInValue(_is_pars.get<Real>("upper_value"));
321  _image_threshold->ReplaceOutOn();
322  _image_threshold->SetOutValue(_is_pars.get<Real>("lower_value"));
323  _image_threshold->SetOutputScalarTypeToDouble();
324 
325  // Perform the thresholding
326  _image_threshold->Update();
327 
328  // Update the pointers
329  _data = _image_threshold->GetOutput();
330  _algorithm = _image_threshold->GetOutputPort();
331 #endif
332 }
unsigned int _component
Component to extract.
Definition: ImageSampler.h:130
vtkAlgorithmOutput * _algorithm
VTK-6 seems to work better in terms of "algorithm outputs" rather than vtkImageData pointers...
Definition: ImageSampler.h:101
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:284
std::vector< std::pair< R1, R2 > > get(const std::string &param1, const std::string &param2) const
Combine two vector parameters into a single vector of pairs.
void scale(MeshBase &mesh, const Real xs, const Real ys=0., const Real zs=0.)
virtual void setupImageSampler(MooseMesh &mesh)
Perform initialization of image data.
Definition: ImageSampler.C:78
const std::vector< std::string > & filenames()
MeshBase & mesh
Base class for MOOSE-based applications.
Definition: MooseApp.h:73
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
vtkSmartPointer< vtkImageReader2 > _image
Complete image data.
Definition: ImageSampler.h:104
ConsoleStream _is_console
Create a console stream object for this helper class.
Definition: ImageSampler.h:140
virtual Real sample(const Point &p) const
Return the pixel value for the given point.
Definition: ImageSampler.C:216
static InputParameters validParams()
BoundingBox _bounding_box
Bounding box for testing points.
Definition: ImageSampler.h:134
InputParameters emptyInputParameters()
vtkImageData * _data
Complete image data.
Definition: ImageSampler.h:98
std::vector< double > _voxel
Physical pixel size.
Definition: ImageSampler.h:126
std::string fileSuffix()
void libmesh_ignore(const Args &...)
vtkSmartPointer< vtkImageMagnitude > _magnitude_filter
Pointer to the magnitude filter.
Definition: ImageSampler.h:113
void vtkShiftAndScale()
Apply image re-scaling using the vtkImageShiftAndRescale object.
Definition: ImageSampler.C:272
Point _origin
Origin of image.
Definition: ImageSampler.h:117
vtkSmartPointer< vtkImageShiftScale > _shift_scale_filter
Pointer to the shift and scaling filter.
Definition: ImageSampler.h:110
const InputParameters & _is_pars
Parameters for interface.
Definition: ImageSampler.h:137
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
static InputParameters validParams()
Constructor.
Definition: ImageSampler.C:18
ImageSampler(const InputParameters &parameters)
Definition: ImageSampler.C:56
void vtkThreshold()
Perform thresholding.
Definition: ImageSampler.C:299
vtkSmartPointer< vtkStringArray > _files
List of file names to extract data.
Definition: ImageSampler.h:95
Point _physical_dims
Physical dimensions of image.
Definition: ImageSampler.h:123
To be called in the validParams functions of classes that need to operate on ranges of files...
std::vector< int > _dims
Pixel dimension of image.
Definition: ImageSampler.h:120
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void vtkMagnitude()
Convert the image to greyscale.
Definition: ImageSampler.C:253
vtkSmartPointer< vtkImageThreshold > _image_threshold
Pointer to thresholding filter.
Definition: ImageSampler.h:107
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
std::array< bool, 3 > _flip
image flip
Definition: ImageSampler.h:143
A 2D GeneratedMesh where xmin, xmax, etc.
Definition: ImageMesh.h:18
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.