15 #include "libmesh/mesh_tools.h" 24 params.
addParam<Point>(
"origin",
"Origin of the image (defaults to mesh origin)");
26 "x,y,z dimensions of the image (defaults to mesh dimensions)");
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.");
34 params.
addParam<
double>(
"shift", 0,
"Value to add to all pixels; occurs prior to scaling");
36 "scale", 1,
"Multiplier to apply to all pixel values; occurs after shifting");
40 params.
addParam<
double>(
"threshold",
"The threshold value");
42 "upper_value", 1,
"The value to set for data greater than the threshold value");
44 "lower_value", 0,
"The value to set for data less than the threshold value");
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");
58 #ifdef LIBMESH_HAVE_VTK
64 (parameters.getCheckedPointerParam<
MooseApp *>(
"_moose_app"))->getOutputWarehouse()),
70 #ifndef LIBMESH_HAVE_VTK 73 mooseError(
"libMesh must be configured with VTK enabled to utilize ImageSampler");
84 #ifdef LIBMESH_HAVE_VTK 86 BoundingBox bbox = MeshTools::create_bounding_box(
mesh.getMesh());
121 std::string file_suffix;
135 mooseError(
"No file range parameters were provided and the Mesh is not an ImageMesh.");
150 _files = vtkSmartPointer<vtkStringArray>::New();
153 _files->InsertNextValue(filename);
156 if (
_files->GetNumberOfValues() == 0)
162 if (file_suffix ==
"png")
163 _image = vtkSmartPointer<vtkPNGReader>::New();
164 else if (file_suffix ==
"tiff" || file_suffix ==
"tif")
165 _image = vtkSmartPointer<vtkTIFFReader>::New();
167 mooseError(
"Un-supported file type '", file_suffix,
"'");
180 int * dims =
_data->GetDimensions();
181 for (
unsigned int i = 0; i < 3; ++i)
183 _dims.push_back(dims[i]);
194 _is_console <<
" ...image read finished" << std::endl;
200 unsigned int n =
_data->GetNumberOfScalarComponents();
203 mooseError(
"'component' parameter must be empty or have a value of 0 to ", n - 1);
218 #ifdef LIBMESH_HAVE_VTK 225 std::vector<int> x(3, 0);
226 for (
int i = 0; i < LIBMESH_DIM; ++i)
234 if (x[i] ==
_dims[i])
239 x[i] =
_dims[i] - x[i] - 1;
244 return _data->GetScalarComponentAsDouble(x[0], x[1], x[2],
_component);
255 #ifdef LIBMESH_HAVE_VTK 274 #ifdef LIBMESH_HAVE_VTK 280 if (shift == 0 &&
scale == 1)
301 #ifdef LIBMESH_HAVE_VTK 308 mooseError(
"When thresholding is applied, both the upper_value and lower_value parameters must " unsigned int _component
Component to extract.
vtkAlgorithmOutput * _algorithm
VTK-6 seems to work better in terms of "algorithm outputs" rather than vtkImageData pointers...
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
virtual void setupImageSampler(MooseMesh &mesh)
Perform initialization of image data.
const std::vector< std::string > & filenames()
Base class for MOOSE-based applications.
vtkSmartPointer< vtkImageReader2 > _image
Complete image data.
ConsoleStream _is_console
Create a console stream object for this helper class.
virtual Real sample(const Point &p) const
Return the pixel value for the given point.
static InputParameters validParams()
BoundingBox _bounding_box
Bounding box for testing points.
vtkImageData * _data
Complete image data.
std::vector< double > _voxel
Physical pixel size.
void libmesh_ignore(const Args &...)
vtkSmartPointer< vtkImageMagnitude > _magnitude_filter
Pointer to the magnitude filter.
void vtkShiftAndScale()
Apply image re-scaling using the vtkImageShiftAndRescale object.
Point _origin
Origin of image.
vtkSmartPointer< vtkImageShiftScale > _shift_scale_filter
Pointer to the shift and scaling filter.
const InputParameters & _is_pars
Parameters for interface.
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
static InputParameters validParams()
Constructor.
ImageSampler(const InputParameters ¶meters)
void vtkThreshold()
Perform thresholding.
vtkSmartPointer< vtkStringArray > _files
List of file names to extract data.
Point _physical_dims
Physical dimensions of image.
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.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void vtkMagnitude()
Convert the image to greyscale.
vtkSmartPointer< vtkImageThreshold > _image_threshold
Pointer to thresholding filter.
std::array< bool, 3 > _flip
image flip
A 2D GeneratedMesh where xmin, xmax, etc.