Line data Source code
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 : #ifdef MOOSE_HAVE_LIBPNG
11 :
12 : #include <fstream>
13 : #include "PNGOutput.h"
14 : #include "FEProblemBase.h"
15 : #include "NonlinearSystem.h"
16 : #include "AuxiliarySystem.h"
17 : #include "libmesh/mesh_tools.h"
18 :
19 : registerMooseObject("MooseApp", PNGOutput);
20 :
21 : InputParameters
22 3123 : PNGOutput::validParams()
23 : {
24 3123 : InputParameters params = FileOutput::validParams();
25 9369 : params.addParam<bool>("transparent_background",
26 6246 : false,
27 : "Determination of whether the background will be transparent.");
28 12492 : params.addRequiredParam<VariableName>("variable",
29 : "The name of the variable to use when creating the image");
30 12492 : params.addParam<Real>("max",
31 : "The maximum for the variable we want to use. If unspecified, the maximum "
32 : "of the variable's DOF values is used");
33 12492 : params.addParam<Real>("min",
34 : "The minimum for the variable we want to use. If unspecified, the minimum "
35 : "of the variable's DOF values is used");
36 :
37 9369 : params.addParam<Point>("frame_center",
38 6246 : Point(0, 0, 0),
39 : "Center of the frame when defining the rectangular plotting region "
40 : "dimensions. Applied after the frame rotation");
41 9369 : params.addParam<Point>(
42 : "first_axis",
43 6246 : Point(1, 0, 0),
44 : "First axis to generate the rectangular plotting region. The axis vector will be normalized");
45 9369 : params.addParam<Point>("second_axis",
46 6246 : Point(0, 1, 0),
47 : "Second axis to generate the rectangular plotting region. The axis vector "
48 : "will be normalized");
49 12492 : params.addParam<Real>(
50 : "min_x",
51 : "Minimum coordinate along the first axis. Defaults to mesh bounding box min X value");
52 12492 : params.addParam<Real>(
53 : "max_x",
54 : "Maximum coordinate along the first axis. Defaults to mesh bounding box max X value");
55 12492 : params.addParam<Real>(
56 : "min_y",
57 : "Minimum coordinate along the second axis. Default to mesh bounding box min Y value");
58 12492 : params.addParam<Real>(
59 : "max_y",
60 : "Maximum coordinate along the second axis. Default to mesh bounding box max Y value");
61 12492 : params.addParamNamesToGroup("frame_center first_axis second_axis min_x max_x min_y max_y",
62 : "Rectangular box extent");
63 :
64 12492 : MooseEnum color("GRAY BRYW BWR RWB BR", "BR");
65 12492 : params.addRequiredParam<MooseEnum>("color", color, "Choose the color scheme to use.");
66 15615 : params.addRangeCheckedParam<unsigned int>(
67 6246 : "resolution", 25, "resolution>0", "The length of the longest side of the image in pixels.");
68 15615 : params.addRangeCheckedParam<Real>("out_bounds_shade",
69 6246 : .5,
70 : "out_bounds_shade>=0 & out_bounds_shade<=1",
71 : "Color for the parts of the image that are out of bounds."
72 : "Value is between 1 and 0.");
73 18738 : params.addRangeCheckedParam<Real>("transparency",
74 : 1,
75 : "transparency>=0 & transparency<=1",
76 : "Value is between 0 and 1"
77 : "where 0 is completely transparent and 1 is completely opaque. "
78 : "Default transparency of the image is no transparency.");
79 12492 : params.addParamNamesToGroup("color resolution out_bounds_shade transparency", "Plotting");
80 3123 : params.addClassDescription("Output data in the PNG format");
81 6246 : return params;
82 3123 : }
83 :
84 31 : PNGOutput::PNGOutput(const InputParameters & parameters)
85 : : FileOutput(parameters),
86 31 : _resolution(getParam<unsigned int>("resolution")),
87 31 : _color(parameters.get<MooseEnum>("color")),
88 93 : _transparent_background(getParam<bool>("transparent_background")),
89 62 : _transparency(getParam<Real>("transparency")),
90 31 : _nl_sys_num(libMesh::invalid_uint),
91 62 : _variable(getParam<VariableName>("variable")),
92 62 : _max(isParamValid("max") ? getParam<Real>("max") : std::numeric_limits<Real>::max()),
93 62 : _min(isParamValid("min") ? getParam<Real>("min") : -std::numeric_limits<Real>::max()),
94 155 : _out_bounds_shade(getParam<Real>("out_bounds_shade"))
95 : {
96 31 : if (!MooseUtils::absoluteFuzzyEqual(
97 155 : getParam<Point>("first_axis") * getParam<Point>("second_axis"), 0))
98 0 : paramWarning("first_axis", "The two axis are not orthogonal");
99 31 : }
100 :
101 : // Function for making the _mesh_function object.
102 : void
103 320 : PNGOutput::makeMeshFunc()
104 : {
105 : // The number assigned to the variable. Used to build the correct mesh. Default is 0.
106 320 : unsigned int variable_number = 0;
107 :
108 : // PNGOutput does not currently scale for running in parallel.
109 320 : if (processor_id() != 0)
110 84 : mooseInfo("PNGOutput is not currently scalable.");
111 :
112 320 : bool var_found = false;
113 320 : if (_problem_ptr->getAuxiliarySystem().hasVariable(_variable))
114 : {
115 77 : variable_number = _problem_ptr->getAuxiliarySystem().getVariable(0, _variable).number();
116 77 : var_found = true;
117 : }
118 :
119 : else
120 486 : for (const auto nl_sys_num : make_range(_problem_ptr->numNonlinearSystems()))
121 243 : if (_problem_ptr->getNonlinearSystem(nl_sys_num).hasVariable(_variable))
122 : {
123 : variable_number =
124 243 : _problem_ptr->getNonlinearSystem(nl_sys_num).getVariable(0, _variable).number();
125 243 : _nl_sys_num = nl_sys_num;
126 243 : var_found = true;
127 : }
128 :
129 320 : if (!var_found)
130 0 : paramError("variable", "This doesn't exist.");
131 :
132 640 : const std::vector<unsigned int> var_nums = {variable_number};
133 :
134 : // If we want the background to be transparent, we need a number over 1.
135 320 : if (_transparent_background)
136 0 : _out_bounds_shade = 2;
137 :
138 : // Find the values that will be used for rescaling purposes.
139 320 : calculateRescalingValues();
140 :
141 : // Set up the mesh_function
142 320 : if (_nl_sys_num == libMesh::invalid_uint)
143 77 : _mesh_function = std::make_unique<libMesh::MeshFunction>(
144 77 : *_es_ptr,
145 77 : _problem_ptr->getAuxiliarySystem().serializedSolution(),
146 77 : _problem_ptr->getAuxiliarySystem().dofMap(),
147 77 : var_nums);
148 : else
149 243 : _mesh_function = std::make_unique<libMesh::MeshFunction>(
150 243 : *_es_ptr,
151 243 : _problem_ptr->getNonlinearSystem(_nl_sys_num).serializedSolution(),
152 243 : _problem_ptr->getNonlinearSystem(_nl_sys_num).dofMap(),
153 243 : var_nums);
154 320 : _mesh_function->init();
155 :
156 : // Need to enable out of mesh with the given control color scaled in reverse
157 : // so when scaling is done, this value retains it's original value.
158 320 : _mesh_function->enable_out_of_mesh_mode(reverseScale(_out_bounds_shade));
159 320 : }
160 :
161 : // Function to find the min and max values so that all the values can be scaled between the two.
162 : void
163 320 : PNGOutput::calculateRescalingValues()
164 : {
165 : // The max and min.
166 : // If the max value wasn't specified in the input file, find it from the system.
167 960 : if (!_pars.isParamValid("max"))
168 : {
169 320 : if (_nl_sys_num == libMesh::invalid_uint)
170 77 : _scaling_max = _problem_ptr->getAuxiliarySystem().serializedSolution().max();
171 : else
172 243 : _scaling_max = _problem_ptr->getNonlinearSystem(_nl_sys_num).serializedSolution().max();
173 : }
174 : else
175 0 : _scaling_max = _max;
176 :
177 : // If the min value wasn't specified in the input file, find it from the system.
178 960 : if (!_pars.isParamValid("min"))
179 : {
180 320 : if (_nl_sys_num == libMesh::invalid_uint)
181 77 : _scaling_min = _problem_ptr->getAuxiliarySystem().serializedSolution().min();
182 : else
183 243 : _scaling_min = _problem_ptr->getNonlinearSystem(_nl_sys_num).serializedSolution().min();
184 : }
185 : else
186 0 : _scaling_min = _min;
187 :
188 : // The amount the values will need to be shifted.
189 320 : _shift_value = 0;
190 :
191 : // Get the shift value.
192 320 : if (_scaling_min != 0)
193 : {
194 : // Shiftvalue will be the same magnitude, but
195 : // going in the opposite direction of the scalingMin
196 71 : _shift_value -= _scaling_min;
197 : }
198 :
199 : // Shift the max.
200 320 : _scaling_max += _shift_value;
201 320 : }
202 :
203 : // Function to apply the scale to the data points.
204 : // Needed to be able to see accurate images that cover the appropriate color spectrum.
205 : inline Real
206 134888 : PNGOutput::applyScale(Real value_to_scale)
207 : {
208 134888 : return ((value_to_scale + _shift_value) / _scaling_max);
209 : }
210 :
211 : // Function to reverse the scaling that happens to a value.
212 : // Needed to be able to accurately control the _out_bounds_shade.
213 : inline Real
214 320 : PNGOutput::reverseScale(Real value_to_unscale)
215 : {
216 320 : return ((value_to_unscale * _scaling_max) - _shift_value);
217 : }
218 :
219 : // Function that controls the colorization of the png image for non-grayscale images.
220 : void
221 134888 : PNGOutput::setRGB(png_byte * rgb, Real selection)
222 : {
223 : // With this system we have a color we start with when the value is 0 and another it approaches as
224 : // the value increases all the way to 255. If we want it to approach another color from that new
225 : // color, it will do so for the next 255, so the transition is from 256 - 511. For each
226 : // additional color we want to transition to, we need another 255. Transitioning from no color, or
227 : // black to Red then Green then Blue then the values of from black as it becomes Red would be 0 -
228 : // 255, Red to Green as 256 - 511 and then Green to Blue as 512 - 767 which gives us our total
229 : // colorSpectrum of 0 - 767, which includes those colors and each of their states in the
230 : // transistion.
231 134888 : unsigned int number_of_destination_colors = 1;
232 134888 : switch (_color)
233 : {
234 : // BRYW. Three destination colors (R,Y,W).
235 0 : case 1:
236 0 : number_of_destination_colors = 3;
237 0 : break;
238 :
239 : // BWR. Two destination colors (W,R).
240 0 : case 2:
241 0 : number_of_destination_colors = 2;
242 0 : break;
243 :
244 : // RWB. Two destination colors (W,B).
245 134888 : case 3:
246 134888 : number_of_destination_colors = 2;
247 134888 : break;
248 :
249 : // BR. One destination color (R).
250 0 : case 4:
251 0 : number_of_destination_colors = 1;
252 0 : break;
253 : }
254 :
255 : // We need to convert the number of colors into the spectrum max, then convert the value from the
256 : // mesh to a point somewhere in the range of 0 to color_spectrum_max.
257 134888 : auto color_spectrum_max = (256 * number_of_destination_colors) - 1;
258 134888 : auto color = (unsigned int)(selection * color_spectrum_max);
259 :
260 : // Unless we specifically say some part is transparent, we want the whole image to be opaque.
261 134888 : auto tran = (unsigned int)(_transparency * 255);
262 :
263 : // Make sure everything is within our colorSpectrum. If it's bigger, then we want a
264 : // transparent background.
265 134888 : if (color > color_spectrum_max)
266 : {
267 0 : color = color_spectrum_max;
268 0 : tran = 0;
269 : }
270 :
271 134888 : auto magnitude = color % 256;
272 :
273 134888 : switch (_color)
274 : {
275 : // Current color scheme: Blue->Red->Yellow->White
276 0 : case 1:
277 : // Blue->Red
278 0 : if (color < 256)
279 : {
280 0 : rgb[0] = magnitude;
281 0 : rgb[1] = 0;
282 0 : rgb[2] = 50; // 255 - magnitude;
283 : }
284 : // Red->Yellow
285 0 : else if (color < 512)
286 : {
287 0 : rgb[0] = 255;
288 0 : rgb[1] = magnitude;
289 0 : rgb[2] = 0;
290 : }
291 : // Yellow->White
292 : else
293 : {
294 0 : rgb[0] = 255;
295 0 : rgb[1] = 255;
296 0 : rgb[2] = magnitude;
297 : }
298 0 : break;
299 :
300 : // Color Scheme: Blue->White->Red
301 : // Using the RGB values found in Paraview
302 0 : case 2:
303 : // Blue->White
304 0 : if (color < 256)
305 : {
306 0 : rgb[0] = (int)(255.0 * (0.231373 + (0.002485 * (float)magnitude)));
307 0 : rgb[1] = (int)(255.0 * (0.298039 + (0.002223 * (float)magnitude)));
308 0 : rgb[2] = (int)(255.0 * (0.752941 + (0.000439 * (float)magnitude)));
309 : }
310 : // White->Red
311 : else
312 : {
313 0 : rgb[0] = (int)(255.0 * (0.865003 - (0.000624 * (float)magnitude)));
314 0 : rgb[1] = (int)(255.0 * (0.865003 - (0.003331 * (float)magnitude)));
315 0 : rgb[2] = (int)(255.0 * (0.865003 - (0.002808 * (float)magnitude)));
316 : }
317 0 : break;
318 :
319 : // Red->White->Blue
320 134888 : case 3:
321 : // Red->White
322 134888 : if (color < 256)
323 : {
324 102893 : rgb[0] = 255;
325 102893 : rgb[1] = magnitude;
326 102893 : rgb[2] = magnitude;
327 : }
328 : // White->Blue
329 : else
330 : {
331 31995 : rgb[0] = 255 - magnitude;
332 31995 : rgb[1] = 255 - magnitude;
333 31995 : rgb[2] = 255;
334 : }
335 134888 : break;
336 :
337 : // Blue->Red
338 0 : case 4:
339 : // Blue->Red
340 0 : rgb[0] = magnitude;
341 0 : rgb[1] = 0;
342 0 : rgb[2] = 255 - magnitude;
343 0 : break;
344 : }
345 : // Add any transparency.
346 134888 : rgb[3] = tran;
347 134888 : }
348 :
349 : void
350 320 : PNGOutput::output()
351 : {
352 320 : makeMeshFunc();
353 320 : _box = MeshTools::create_bounding_box(*_mesh_ptr);
354 :
355 : // Make sure this happens on processor 0
356 320 : if (processor_id() == 0)
357 236 : makePNG();
358 320 : }
359 :
360 : // Function the writes the PNG out to the appropriate filename.
361 : void
362 236 : PNGOutput::makePNG()
363 : {
364 : // Get the max and min of the BoundingBox
365 236 : Point max_point = _box.max();
366 236 : Point min_point = _box.min();
367 :
368 : // Modify the bounding box if the user passed different bounds
369 708 : if (isParamValid("min_x"))
370 0 : min_point(0) = getParam<Real>("min_x");
371 708 : if (isParamValid("min_y"))
372 0 : min_point(1) = getParam<Real>("min_y");
373 708 : if (isParamValid("max_x"))
374 0 : max_point(0) = getParam<Real>("max_x");
375 708 : if (isParamValid("max_y"))
376 0 : max_point(1) = getParam<Real>("max_y");
377 :
378 : // The total distance on the two plotting axes.
379 236 : Real dist_x = max_point(0) - min_point(0);
380 236 : Real dist_y = max_point(1) - min_point(1);
381 :
382 : // Width and height for the PNG image.
383 : Real width;
384 : Real height;
385 :
386 : // The longer dimension becomes the value to which we scale the other.
387 236 : if (dist_x > dist_y)
388 : {
389 12 : width = _resolution;
390 12 : height = (_resolution / dist_x) * dist_y;
391 : }
392 : else
393 : {
394 224 : height = _resolution;
395 224 : width = (_resolution / dist_y) * dist_x;
396 : }
397 :
398 : // Create the filename based on base and the test step number.
399 236 : std::ostringstream png_file;
400 236 : png_file << _file_base << "_" << std::setfill('0') << std::setw(3) << _t_step << ".png";
401 :
402 : // libpng is built on C, so by default it takes FILE*.
403 236 : FILE * fp = nullptr;
404 236 : png_structp pngp = nullptr;
405 236 : png_infop infop = nullptr;
406 : // Required depth for proper image clarity.
407 236 : Real depth = 8;
408 : // Allocate resources.
409 236 : std::vector<png_byte> row((width + 1) * 4);
410 :
411 : // Check if we can open and write to the file.
412 236 : MooseUtils::checkFileWriteable(png_file.str());
413 :
414 : // Open the file with write and bit modes.
415 236 : fp = fopen(png_file.str().c_str(), "wb");
416 :
417 : // Verify the libpng we are linked against at runtime matches the headers we
418 : // were compiled against. PNG_LIBPNG_VER is major*10000 + minor*100 + release;
419 : // libpng only guarantees ABI within the same major.minor, so a mismatch there
420 : // can cause png_create_write_struct to segfault rather than fail cleanly.
421 236 : const png_uint_32 runtime_libpng_ver = png_access_version_number();
422 236 : if (runtime_libpng_ver / 100 != PNG_LIBPNG_VER / 100)
423 : {
424 0 : mooseWarning("PNG output skipped for ",
425 0 : png_file.str(),
426 : ": runtime libpng version ",
427 : runtime_libpng_ver,
428 : " is ABI-incompatible with the version ",
429 0 : static_cast<png_uint_32>(PNG_LIBPNG_VER),
430 : " used at compile time.");
431 0 : if (fp != nullptr)
432 0 : fclose(fp);
433 0 : return;
434 : }
435 :
436 236 : pngp = png_create_write_struct(PNG_LIBPNG_VER_STRING, nullptr, nullptr, nullptr);
437 236 : if (!pngp)
438 0 : mooseError("Failed to make the pointer string for the png.");
439 :
440 236 : infop = png_create_info_struct(pngp);
441 236 : if (!infop)
442 0 : mooseError("Failed to make an info pointer for the png.");
443 :
444 : // Initializes the IO for the png. Needs FILE* to compile.
445 236 : png_init_io(pngp, fp);
446 :
447 : // Set up the PNG header.
448 236 : png_set_IHDR(pngp,
449 : infop,
450 : width,
451 : height,
452 : depth,
453 236 : (_color ? PNG_COLOR_TYPE_RGBA : PNG_COLOR_TYPE_GRAY),
454 : PNG_INTERLACE_NONE,
455 : PNG_COMPRESSION_TYPE_DEFAULT,
456 : PNG_FILTER_TYPE_DEFAULT);
457 :
458 236 : png_write_info(pngp, infop);
459 :
460 : // Initializing the point that will be used for populating the mesh values.
461 236 : Point pt(0, 0, 0);
462 :
463 : // Pre-compute the transformations
464 236 : Point center(0, 0, 0);
465 708 : if (isParamValid("frame_center"))
466 708 : center = getParam<Point>("frame_center");
467 : else
468 0 : center = Point((min_point(0) + max_point(0)) / 2, (min_point(1) + max_point(1)) / 2, 0);
469 :
470 472 : const auto first_axis = getParam<Point>("first_axis").unit();
471 472 : const auto second_axis = getParam<Point>("second_axis").unit();
472 :
473 : // Dense vector that we can pass into the _mesh_function to fill with a value for a given point.
474 236 : DenseVector<Number> dv(0);
475 :
476 : // Loop through to create the image.
477 6124 : for (const auto iy : make_range((unsigned int)height))
478 : {
479 140776 : for (const auto ix : make_range((unsigned int)width + 1))
480 : {
481 : // Move along the axis from the center
482 : // Center the point in the pixel
483 134888 : pt = center +
484 134888 : (min_point(0) + Real(ix) / _resolution * (max_point(0) - min_point(0))) * first_axis +
485 269776 : (min_point(1) + Real(iy) / _resolution * (max_point(1) - min_point(1))) * second_axis;
486 :
487 134888 : (*_mesh_function)(pt, _time, dv, nullptr);
488 :
489 : // Determine whether to create the PNG in color or grayscale
490 134888 : if (_color)
491 134888 : setRGB(&row.data()[ix * 4], applyScale(dv(0)));
492 : else
493 0 : row.data()[ix] = applyScale(dv(0)) * 255;
494 : }
495 5888 : png_write_row(pngp, row.data());
496 : }
497 :
498 : // Close the file and take care of some other png end stuff.
499 236 : png_write_end(pngp, nullptr);
500 236 : if (fp != nullptr)
501 236 : fclose(fp);
502 236 : if (infop != nullptr)
503 236 : png_free_data(pngp, infop, PNG_FREE_ALL, -1);
504 236 : if (pngp != nullptr)
505 236 : png_destroy_write_struct(&pngp, &infop);
506 236 : }
507 :
508 : #endif
|
