• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

openmc-dev / openmc / 28819286649

06 Jul 2026 07:57PM UTC coverage: 81.26% (-0.03%) from 81.289%
28819286649

Pull #3998

github

web-flow
Merge 9239c236e into 0c6b3fb83
Pull Request #3998: Tracks to vtk Function for Lost Particles

18189 of 26398 branches covered (68.9%)

Branch coverage included in aggregate %.

3 of 40 new or added lines in 1 file covered. (7.5%)

138 existing lines in 6 files now uncovered.

59396 of 69079 relevant lines covered (85.98%)

48606606.3 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

61.61
/include/openmc/plot.h
1
#ifndef OPENMC_PLOT_H
2
#define OPENMC_PLOT_H
3

4
#include <cmath>
5
#include <sstream>
6
#include <unordered_map>
7
#include <unordered_set>
8
#include <vector>
9

10
#include "openmc/tensor.h"
11
#include "pugixml.hpp"
12

13
#include "hdf5.h"
14
#include "openmc/cell.h"
15
#include "openmc/constants.h"
16
#include "openmc/error.h"
17
#include "openmc/geometry.h"
18
#include "openmc/particle.h"
19
#include "openmc/position.h"
20
#include "openmc/random_lcg.h"
21
#include "openmc/ray.h"
22
#include "openmc/tallies/filter.h"
23
#include "openmc/tallies/filter_match.h"
24
#include "openmc/xml_interface.h"
25

26
namespace openmc {
27

28
//===============================================================================
29
// Global variables
30
//===============================================================================
31

32
class PlottableInterface;
33

34
namespace model {
35

36
extern std::unordered_map<int, int> plot_map; //!< map of plot ids to index
37
extern vector<std::unique_ptr<PlottableInterface>>
38
  plots; //!< Plot instance container
39

40
extern uint64_t plotter_seed; // Stream index used by the plotter
41

42
} // namespace model
43

44
//===============================================================================
45
// RGBColor holds color information for plotted objects
46
//===============================================================================
47

48
struct RGBColor {
49
  // Constructors
50
  RGBColor() : red(0), green(0), blue(0) {};
387,710!
51
  RGBColor(const int v[3]) : red(v[0]), green(v[1]), blue(v[2]) {};
52
  RGBColor(int r, int g, int b) : red(r), green(g), blue(b) {};
22✔
53

54
  RGBColor(const vector<int>& v)
308✔
55
  {
308✔
56
    if (v.size() != 3) {
308!
UNCOV
57
      throw std::out_of_range("Incorrect vector size for RGBColor.");
×
58
    }
59
    red = v[0];
308✔
60
    green = v[1];
308✔
61
    blue = v[2];
308✔
62
  }
308✔
63

64
  bool operator==(const RGBColor& other)
6,882✔
65
  {
66
    return red == other.red && green == other.green && blue == other.blue;
3,441!
67
  }
68

69
  RGBColor& operator*=(const double x)
741,631✔
70
  {
71
    red *= x;
741,631✔
72
    green *= x;
741,631✔
73
    blue *= x;
741,631✔
74
    return *this;
741,631✔
75
  }
76

77
  // Members
78
  uint8_t red, green, blue;
79
};
80

81
// some default colors
82
const RGBColor WHITE {255, 255, 255};
83
const RGBColor RED {255, 0, 0};
84
const RGBColor BLACK {0, 0, 0};
85

86
/**
87
 * \class PlottableInterface
88
 * \brief Interface for plottable objects.
89
 *
90
 * PlottableInterface classes must have unique IDs. If no ID (or -1) is
91
 * provided, the next available ID is assigned automatically. They guarantee
92
 * the ability to create output in some form. This interface is designed to be
93
 * implemented by classes that produce plot-relevant data which can be
94
 * visualized.
95
 */
96

97
typedef tensor::Tensor<RGBColor> ImageData;
98
class PlottableInterface {
99
public:
100
  PlottableInterface() = default;
11✔
101

102
  void set_default_colors();
103

104
private:
105
  void set_id(pugi::xml_node plot_node);
106
  int id_ {C_NONE}; // unique plot ID
107

108
  void set_bg_color(pugi::xml_node plot_node);
109
  void set_universe(pugi::xml_node plot_node);
110
  void set_color_by(pugi::xml_node plot_node);
111
  void set_user_colors(pugi::xml_node plot_node);
112
  void set_overlap_color(pugi::xml_node plot_node);
113
  void set_mask(pugi::xml_node plot_node);
114

115
protected:
116
  // Plot output filename, derived classes have logic to set it
117
  std::string path_plot_;
118

119
public:
120
  enum class PlotColorBy { cells = 0, mats = 1 };
121

122
  // Generates image data based on plot parameters and returns it
123
  virtual ImageData create_image() const = 0;
124

125
  // Creates the output image named path_plot_
126
  virtual void create_output() const = 0;
127

128
  // Write populated image data to file
129
  void write_image(const ImageData& data) const;
130

131
  // Print useful info to the terminal
132
  virtual void print_info() const = 0;
133

134
  const std::string& path_plot() const { return path_plot_; }
231✔
135
  std::string& path_plot() { return path_plot_; }
253✔
136
  int id() const { return id_; }
2,235!
137
  void set_id(int id = C_NONE);
138
  int level() const { return level_; }
242✔
139
  PlotColorBy color_by() const { return color_by_; }
140

141
  // Public color-related data
142
  PlottableInterface(pugi::xml_node plot_node);
143
  virtual ~PlottableInterface() = default;
981✔
144
  int level_ {-1};                           // Universe level to plot
145
  bool color_overlaps_ {false};              // Show overlapping cells?
146
  PlotColorBy color_by_ {PlotColorBy::mats}; // Plot coloring (cell/material)
147
  RGBColor not_found_ {WHITE};               // Plot background color
148
  RGBColor overlap_color_ {RED};             // Plot overlap color
149
  vector<RGBColor> colors_;                  // Plot colors
150
};
151

152
struct IdData {
9,746✔
153
  // Constructor
154
  IdData(size_t h_res, size_t v_res, bool include_filter = false);
155

156
  // Methods
157
  void set_value(size_t y, size_t x, const Particle& p, int level,
158
    Filter* filter = nullptr, FilterMatch* match = nullptr);
159
  void set_overlap(size_t y, size_t x, int overlap_idx);
160

161
  // Members
162
  tensor::Tensor<int32_t> data_; //!< 2D array of cell & material ids
163
};
164

UNCOV
165
struct PropertyData {
×
166
  // Constructor
167
  PropertyData(size_t h_res, size_t v_res, bool include_filter = false);
168

169
  // Methods
170
  void set_value(size_t y, size_t x, const Particle& p, int level,
171
    Filter* filter = nullptr, FilterMatch* match = nullptr);
172
  void set_overlap(size_t y, size_t x, int overlap_idx);
173

174
  // Members
175
  tensor::Tensor<double> data_; //!< 2D array of temperature & density data
176
};
177

178
struct RasterData {
179
  // Constructor
180
  RasterData(size_t h_res, size_t v_res, bool include_filter = false);
181

182
  // Methods
183
  void set_value(size_t y, size_t x, const Particle& p, int level,
184
    Filter* filter = nullptr, FilterMatch* match = nullptr);
185
  void set_overlap(size_t y, size_t x, int overlap_idx);
186

187
  // Members
188
  tensor::Tensor<int32_t>
189
    id_data_; //!< [v_res, h_res, 3 or 4]: cell, instance, mat, [filter_bin]
190
  tensor::Tensor<double>
191
    property_data_;     //!< [v_res, h_res, 2]: temperature, density
192
  bool include_filter_; //!< Whether filter bin index is included
193
};
194

195
//===============================================================================
196
// Plot class
197
//===============================================================================
198

199
class SlicePlotBase {
1,367✔
200
public:
201
  template<class T>
202
  T get_map(int32_t filter_index = -1) const;
203

204
  enum class PlotBasis { xy = 1, xz = 2, yz = 3 };
205

206
  // Accessors
207

208
  const std::array<size_t, 3>& pixels() const { return pixels_; }
116,347!
209
  std::array<size_t, 3>& pixels() { return pixels_; }
882✔
210

211
  // Members
212
public:
213
  Position origin_;         //!< Plot origin in geometry
214
  Direction u_span_;        //!< Full-width span vector in geometry
215
  Direction v_span_;        //!< Full-height span vector in geometry
216
  array<size_t, 3> pixels_; //!< Plot size in pixels
217
  bool show_overlaps_;      //!< Show overlapping cells?
218
  int slice_level_ {-1};    //!< Plot universe level
219
private:
220
};
221

222
template<class T>
223
T SlicePlotBase::get_map(int32_t filter_index) const
5,294✔
224
{
225

226
  size_t width = pixels_[0];
5,294✔
227
  size_t height = pixels_[1];
5,294✔
228

229
  // Determine if filter is being used
230
  bool include_filter = (filter_index >= 0);
5,294✔
231
  Filter* filter = nullptr;
5,294✔
232
  if (include_filter) {
5,294✔
233
    filter = model::tally_filters[filter_index].get();
22✔
234
  }
235

236
  // size data array
237
  T data(width, height, include_filter);
5,294✔
238

239
  // compute pixel steps and top-left pixel center
240
  Direction u_step = u_span_ / static_cast<double>(width);
5,294✔
241
  Direction v_step = v_span_ / static_cast<double>(height);
5,294✔
242

243
  Position start =
244
    origin_ - 0.5 * u_span_ + 0.5 * v_span_ + 0.5 * u_step - 0.5 * v_step;
5,294✔
245

246
  // Validate that span vectors define a valid plane
247
  Position cross = u_span_.cross(v_span_);
5,294✔
248
  if (cross.norm() == 0.0) {
5,294!
UNCOV
249
    fatal_error("Slice span vectors are invalid (zero area).");
×
250
  }
251

252
  // Use an arbitrary direction that is not aligned with any coordinate axis.
253
  // The direction has no physical meaning for plotting but is used by
254
  // Surface::sense() to break ties when a pixel is coincident with a surface.
255
  Direction dir = {1.0 / std::sqrt(2.0), 1.0 / std::sqrt(2.0), 0.0};
5,294✔
256

257
#pragma omp parallel
2,889✔
258
  {
259
    Particle p;
2,405✔
260
    p.r() = start;
2,405✔
261
    p.u() = dir;
2,405✔
262
    p.coord(0).universe() = model::root_universe;
2,405✔
263
    int level = slice_level_;
2,405✔
264
    int j {};
2,405✔
265
    FilterMatch match;
2,405✔
266

267
#pragma omp for
268
    for (int y = 0; y < height; y++) {
444,245✔
269
      Position row = start - v_step * static_cast<double>(y);
441,840✔
270
      for (int x = 0; x < width; x++) {
87,381,240✔
271
        p.r() = row + u_step * static_cast<double>(x);
86,939,400✔
272
        p.n_coord() = 1;
86,939,400✔
273
        // local variables
274
        bool found_cell = exhaustive_find_cell(p);
86,939,400✔
275
        j = p.n_coord() - 1;
86,939,400✔
276
        if (level >= 0) {
86,939,400✔
277
          j = level;
12,500✔
278
        }
279
        if (found_cell) {
86,939,400✔
280
          data.set_value(y, x, p, j, filter, &match);
17,745,190✔
281
        }
282
        if (show_overlaps_) {
86,939,400✔
283
          int overlap_idx = check_cell_overlap(p, false);
655,500✔
284
          if (overlap_idx >= 0) {
655,500✔
285
            data.set_overlap(y, x, overlap_idx);
179,100✔
286
          }
287
        }
288
      } // inner for
289
    }
290
  }
2,405✔
291

292
  return data;
5,294✔
UNCOV
293
}
×
294

295
// Represents either a voxel or pixel plot
296
class Plot : public PlottableInterface, public SlicePlotBase {
297

298
public:
299
  enum class PlotType { slice = 1, voxel = 2 };
300

301
  Plot(pugi::xml_node plot, PlotType type);
302

303
private:
304
  void set_output_path(pugi::xml_node plot_node);
305
  void set_basis(pugi::xml_node plot_node);
306
  void set_origin(pugi::xml_node plot_node);
307
  void set_width(pugi::xml_node plot_node);
308
  void set_meshlines(pugi::xml_node plot_node);
309

310
public:
311
  // Add mesh lines to ImageData
312
  void draw_mesh_lines(ImageData& data) const;
313
  ImageData create_image() const override;
314
  void create_voxel() const;
315

316
  void create_output() const override;
317
  void print_info() const override;
318

319
  PlotType type_;                 //!< Plot type (Slice/Voxel)
320
  Position width_;                //!< Axis-aligned width from plot.xml
321
  PlotBasis basis_;               //!< Basis from plot.xml for slice plots
322
  int meshlines_width_;           //!< Width of lines added to the plot
323
  int index_meshlines_mesh_ {-1}; //!< Index of the mesh to draw on the plot
324
  RGBColor meshlines_color_;      //!< Color of meshlines on the plot
325
};
326

327
/**
328
 * \class RaytracePlot
329
 * \brief Base class for plots that generate images through ray tracing.
330
 *
331
 * This class serves as a base for plots that create their visuals by tracing
332
 * rays from a camera through the problem geometry. It inherits from
333
 * PlottableInterface, ensuring that it provides an implementation for
334
 * generating output specific to ray-traced visualization. WireframeRayTracePlot
335
 * and SolidRayTracePlot provide concrete implementations of this class.
336
 */
UNCOV
337
class RayTracePlot : public PlottableInterface {
×
338
public:
339
  RayTracePlot() = default;
11✔
340
  RayTracePlot(pugi::xml_node plot);
341

342
  // Standard getters. No setting since it's done from XML.
343
  const Position& camera_position() const { return camera_position_; }
344
  Position& camera_position() { return camera_position_; }
11✔
345
  const Position& look_at() const { return look_at_; }
346
  Position& look_at() { return look_at_; }
11✔
347

348
  const double& horizontal_field_of_view() const
349
  {
350
    return horizontal_field_of_view_;
351
  }
352
  double& horizontal_field_of_view() { return horizontal_field_of_view_; }
353

354
  void print_info() const override;
355

356
  const std::array<int, 2>& pixels() const { return pixels_; }
9,190,381!
357
  std::array<int, 2>& pixels() { return pixels_; }
154!
358

359
  const Direction& up() const { return up_; }
360
  Direction& up() { return up_; }
11✔
361

362
  //! brief Updates the cached camera-to-model matrix after changes to
363
  //! camera parameters.
364
  void update_view();
365

366
protected:
367
  Direction camera_x_axis() const
440,000✔
368
  {
369
    return {camera_to_model_[0], camera_to_model_[3], camera_to_model_[6]};
440,000✔
370
  }
371

372
  Direction camera_y_axis() const
440,000✔
373
  {
374
    return {camera_to_model_[1], camera_to_model_[4], camera_to_model_[7]};
440,000✔
375
  }
376

377
  Direction camera_z_axis() const
440,000✔
378
  {
379
    return {camera_to_model_[2], camera_to_model_[5], camera_to_model_[8]};
440,000✔
380
  }
381

382
  void set_output_path(pugi::xml_node plot_node);
383

384
  /*
385
   * Gets the starting position and direction for the pixel corresponding
386
   * to this horizontal and vertical position.
387
   */
388
  std::pair<Position, Direction> get_pixel_ray(int horiz, int vert) const;
389

390
private:
391
  void set_look_at(pugi::xml_node node);
392
  void set_camera_position(pugi::xml_node node);
393
  void set_field_of_view(pugi::xml_node node);
394
  void set_pixels(pugi::xml_node node);
395
  void set_orthographic_width(pugi::xml_node node);
396

397
  double horizontal_field_of_view_ {70.0}; // horiz. f.o.v. in degrees
398
  Position camera_position_;               // where camera is
399
  Position look_at_;                     // point camera is centered looking at
400
  std::array<int, 2> pixels_ {100, 100}; // pixel dimension of resulting image
401
  Direction up_ {0.0, 0.0, 1.0};         // which way is up
402

403
  /* The horizontal thickness, if using an orthographic projection.
404
   * If set to zero, we assume using a perspective projection.
405
   */
406
  double orthographic_width_ {C_NONE};
407

408
  /*
409
   * Cached camera-to-model matrix with column vectors of axes. The x-axis is
410
   * the vector between the camera_position_ and look_at_; the y-axis is the
411
   * cross product of the x-axis with the up_ vector, and the z-axis is the
412
   * cross product of the x and y axes.
413
   */
414
  std::array<double, 9> camera_to_model_;
415
};
416

417
class ProjectionRay;
418

419
/**
420
 * \class WireframeRayTracePlot
421
 * \brief Creates plots that are like colorful x-ray imaging
422
 *
423
 * WireframeRayTracePlot is a specialized form of RayTracePlot designed for
424
 * creating projection plots. This involves tracing rays from a camera through
425
 * the problem geometry and rendering the results based on depth of penetration
426
 * through materials or cells and their colors.
427
 */
428
class WireframeRayTracePlot : public RayTracePlot {
429

430
  friend class ProjectionRay;
431

432
public:
433
  WireframeRayTracePlot(pugi::xml_node plot);
434

435
  ImageData create_image() const override;
436
  void create_output() const override;
437
  void print_info() const override;
438

439
private:
440
  void set_opacities(pugi::xml_node node);
441
  void set_wireframe_thickness(pugi::xml_node node);
442
  void set_wireframe_ids(pugi::xml_node node);
443
  void set_wireframe_color(pugi::xml_node node);
444

445
  /* Checks if a vector of two TrackSegments is equivalent. We define this
446
   * to mean not having matching intersection lengths, but rather having
447
   * a matching sequence of surface/cell/material intersections.
448
   */
449
  struct TrackSegment;
450
  bool trackstack_equivalent(const vector<TrackSegment>& track1,
451
    const vector<TrackSegment>& track2) const;
452

453
  /* Used for drawing wireframe and colors. We record the list of
454
   * surface/cell/material intersections and the corresponding lengths as a ray
455
   * traverses the geometry, then color by iterating in reverse.
456
   */
457
  struct TrackSegment {
458
    int id;        // material or cell ID (which is being colored)
459
    double length; // length of this track intersection
460

461
    /* Recording this allows us to draw edges on the wireframe. For instance
462
     * if two surfaces bound a single cell, it allows drawing that sharp edge
463
     * where the surfaces intersect.
464
     */
465
    int surface_index {-1}; // last surface index intersected in this segment
466
    TrackSegment(int id_a, double length_a, int surface_a)
2,359,148✔
467
      : id(id_a), length(length_a), surface_index(surface_a)
2,359,148✔
468
    {}
469
  };
470

471
  // which color IDs should be wireframed. If empty, all cells are wireframed.
472
  vector<int> wireframe_ids_;
473

474
  // Thickness of the wireframe lines. Can set to zero for no wireframe.
475
  int wireframe_thickness_ {1};
476

477
  RGBColor wireframe_color_ {BLACK}; // wireframe color
478
  vector<double> xs_; // macro cross section values for cell volume rendering
479
};
480

481
/**
482
 * \class SolidRayTracePlot
483
 * \brief Plots 3D objects as the eye might see them.
484
 *
485
 * Plots a geometry with single-scattered Phong lighting plus a diffuse lighting
486
 * contribution. The result is a physically reasonable, aesthetic 3D view of a
487
 * geometry.
488
 */
489
class SolidRayTracePlot : public RayTracePlot {
490
  friend class PhongRay;
491

492
public:
493
  SolidRayTracePlot() = default;
11✔
494

495
  SolidRayTracePlot(pugi::xml_node plot);
496

497
  ImageData create_image() const override;
498
  void create_output() const override;
499
  void print_info() const override;
500

501
  const std::unordered_set<int>& opaque_ids() const { return opaque_ids_; }
502
  std::unordered_set<int>& opaque_ids() { return opaque_ids_; }
22✔
503

504
  const Position& light_location() const { return light_location_; }
505
  Position& light_location() { return light_location_; }
11✔
506

507
  const double& diffuse_fraction() const { return diffuse_fraction_; }
508
  double& diffuse_fraction() { return diffuse_fraction_; }
509

510
private:
511
  void set_opaque_ids(pugi::xml_node node);
512
  void set_light_position(pugi::xml_node node);
513
  void set_diffuse_fraction(pugi::xml_node node);
514

515
  std::unordered_set<int> opaque_ids_;
516

517
  double diffuse_fraction_ {0.1};
518

519
  // By default, the light is at the camera unless otherwise specified.
520
  Position light_location_;
521
};
522

523
class ProjectionRay : public Ray {
1,000,000✔
524
public:
525
  ProjectionRay(Position r, Direction u, const WireframeRayTracePlot& plot,
1,000,000✔
526
    vector<WireframeRayTracePlot::TrackSegment>& line_segments)
527
    : Ray(r, u), plot_(plot), line_segments_(line_segments)
1,000,000✔
528
  {}
529

530
  void on_intersection() override;
531

532
private:
533
  /* Store a reference to the plot object which is running this ray, in order
534
   * to access some of the plot settings which influence the behavior where
535
   * intersections are.
536
   */
537
  const WireframeRayTracePlot& plot_;
538

539
  /* The ray runs through the geometry, and records the lengths of ray segments
540
   * and cells they lie in along the way.
541
   */
542
  vector<WireframeRayTracePlot::TrackSegment>& line_segments_;
543
};
544

545
class PhongRay : public Ray {
600,480✔
546
public:
547
  PhongRay(Position r, Direction u, const SolidRayTracePlot& plot)
1,321,056✔
548
    : Ray(r, u), plot_(plot)
1,321,056✔
549
  {
550
    result_color_ = plot_.not_found_;
1,321,056✔
551
  }
1,321,056✔
552

553
  void on_intersection() override;
554

555
  const RGBColor& result_color() { return result_color_; }
556

557
private:
558
  const SolidRayTracePlot& plot_;
559

560
  /* After the ray is reflected, it is moving towards the
561
   * camera. It does that in order to see if the exposed surface
562
   * is shadowed by something else.
563
   */
564
  bool reflected_ {false};
565

566
  // Have to record the first hit ID, so that if the region
567
  // does get shadowed, we recall what its color should be
568
  // when tracing from the surface to the light.
569
  int orig_hit_id_ {-1};
570

571
  RGBColor result_color_;
572
};
573

574
//===============================================================================
575
// Non-member functions
576
//===============================================================================
577

578
/* Write a PPM image
579
 * filename - name of output file
580
 * data - image data to write
581
 */
582
void output_ppm(const std::string& filename, const ImageData& data);
583

584
#ifdef USE_LIBPNG
585
/* Write a PNG image
586
 * filename - name of output file
587
 * data - image data to write
588
 */
589
void output_png(const std::string& filename, const ImageData& data);
590
#endif
591

592
//! Initialize a voxel file
593
//! \param[in] id of an open hdf5 file
594
//! \param[in] dimensions of the voxel file (dx, dy, dz)
595
//! \param[out] dataspace pointer to voxel data
596
//! \param[out] dataset pointer to voxesl data
597
//! \param[out] pointer to memory space of voxel data
598
void voxel_init(hid_t file_id, const hsize_t* dims, hid_t* dspace, hid_t* dset,
599
  hid_t* memspace);
600

601
//! Write a section of the voxel data to hdf5
602
//! \param[in] voxel slice
603
//! \param[out] dataspace pointer to voxel data
604
//! \param[out] dataset pointer to voxesl data
605
//! \param[out] pointer to data to write
606
void voxel_write_slice(
607
  int x, hid_t dspace, hid_t dset, hid_t memspace, void* buf);
608

609
//! Close voxel file entities
610
//! \param[in] data space to close
611
//! \param[in] dataset to close
612
//! \param[in] memory space to close
613
void voxel_finalize(hid_t dspace, hid_t dset, hid_t memspace);
614

615
//===============================================================================
616
// External functions
617
//===============================================================================
618

619
//! Read plot specifications from a plots.xml file
620
void read_plots_xml();
621

622
//! Read plot specifications from an XML Node
623
//! \param[in] XML node containing plot info
624
void read_plots_xml(pugi::xml_node root);
625

626
//! Clear memory
627
void free_memory_plot();
628

629
//! Create a randomly generated RGB color
630
//! \return RGBColor with random value
631
RGBColor random_color();
632

633
} // namespace openmc
634
#endif // OPENMC_PLOT_H
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2026 Coveralls, Inc