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