Address #1101 comments

This commit is contained in:
Sterling Harper 2018-10-22 21:55:06 -04:00
parent fa419332c9
commit 321b255aea
35 changed files with 200 additions and 98 deletions

View file

@ -7,7 +7,7 @@
#include "openmc/hdf5_interface.h"
#include "openmc/particle.h"
#include "openmc/xml_interface.h"
#include "pugixml.hpp"
namespace openmc {
@ -43,10 +43,10 @@ public:
class Filter
{
public:
virtual std::string type() const = 0;
virtual ~Filter() = 0;
virtual std::string type() const = 0;
//! Uses an XML input to fill the filter's data fields.
virtual void from_xml(pugi::xml_node node) = 0;
@ -55,7 +55,7 @@ public:
//! \param[out] match will contain the matching bins and corresponding
//! weights; note that there may be zero matching bins
virtual void
get_all_bins(Particle* p, int estimator, FilterMatch& match) const = 0;
get_all_bins(const Particle* p, int estimator, FilterMatch& match) const = 0;
//! Writes data describing this filter to an HDF5 statepoint group.
virtual void

View file

@ -21,7 +21,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -24,14 +24,17 @@ public:
void initialize() override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
//! The indices of the cells binned by this filter.
std::vector<int32_t> cells_;
//! A map from cell indices to filter bin indices.
std::unordered_map<int32_t, int> map_;
};

View file

@ -16,7 +16,7 @@ class CellbornFilter : public CellFilter
public:
std::string type() const override {return "cellborn";}
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
std::string text_label(int bin) const override;

View file

@ -16,7 +16,7 @@ class CellFromFilter : public CellFilter
public:
std::string type() const override {return "cellfrom";}
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
std::string text_label(int bin) const override;

View file

@ -22,7 +22,7 @@ public:
void initialize() override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -27,14 +27,17 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
//! Incident neutron energy interpolation grid.
std::vector<double> energy_;
//! Interpolant values.
std::vector<double> y_;
};

View file

@ -20,7 +20,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -24,14 +24,17 @@ public:
void initialize() override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
//! The indices of the materials binned by this filter.
std::vector<int32_t> materials_;
//! A map from material indices to filter bin indices.
std::unordered_map<int32_t, int> map_;
};

View file

@ -22,7 +22,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -12,7 +12,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
std::string text_label(int bin) const override;

View file

@ -21,7 +21,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -21,7 +21,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;

View file

@ -25,7 +25,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
@ -33,6 +33,8 @@ public:
std::string text_label(int bin) const override;
int order_;
//! The type of angle that this filter measures when binning events.
SphericalHarmonicsCosine cosine_ {SphericalHarmonicsCosine::particle};
};

View file

@ -25,7 +25,7 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
@ -33,8 +33,15 @@ public:
std::string text_label(int bin) const override;
int order_;
//! The Cartesian coordinate axis that the Legendre expansion is applied to.
LegendreAxis axis_;
double min_, max_;
//! The minimum coordinate along the reference axis that the expansion covers.
double min_;
//! The maximum coordinate along the reference axis that the expansion covers.
double max_;
};
} // namespace openmc

View file

@ -24,14 +24,17 @@ public:
void initialize() override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
//! The indices of the surfaces binned by this filter.
std::vector<int32_t> surfaces_;
//! A map from surface indices to filter bin indices.
std::unordered_map<int32_t, int> map_;
};

View file

@ -24,14 +24,17 @@ public:
void initialize() override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
//! The indices of the universes binned by this filter.
std::vector<int32_t> universes_;
//! A map from universe indices to filter bin indices.
std::unordered_map<int32_t, int> map_;
};

View file

@ -20,17 +20,28 @@ public:
void from_xml(pugi::xml_node node) override;
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
void to_statepoint(hid_t filter_group) const override;
std::string text_label(int bin) const override;
virtual void calc_n_bins() {n_bins_ = ((order_+1) * (order_+2)) / 2;}
int order() const {return order_;}
virtual void set_order(int order);
//! Cartesian x coordinate for the origin of this expansion.
double x_;
//! Cartesian y coordinate for the origin of this expansion.
double y_;
//! Maximum radius from the origin covered by this expansion.
double r_;
protected:
int order_;
double x_, y_, r_;
};
//==============================================================================
@ -42,12 +53,12 @@ class ZernikeRadialFilter : public ZernikeFilter
public:
std::string type() const override {return "zernikeradial";}
void get_all_bins(Particle* p, int estimator, FilterMatch& match)
void get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const override;
std::string text_label(int bin) const override;
void calc_n_bins() override {n_bins_ = order_ / 2 + 1;}
void set_order(int order) override;
};
} // namespace openmc

View file

@ -6,6 +6,7 @@
#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/search.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -36,8 +37,8 @@ AzimuthalFilter::from_xml(pugi::xml_node node)
}
void
AzimuthalFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
AzimuthalFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
double phi;
if (estimator == ESTIMATOR_TRACKLENGTH) {
@ -46,10 +47,11 @@ const
phi = std::atan2(p->last_uvw[1], p->last_uvw[0]);
}
if (phi >= bins_[0] && phi <= bins_.back()) {
if (phi >= bins_.front() && phi <= bins_.back()) {
//TODO: off-by-one
auto bin = lower_bound_index(bins_.begin(), bins_.end(), phi) + 1;
match.bins_.push_back(bin);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}

View file

@ -5,6 +5,7 @@
#include "openmc/capi.h"
#include "openmc/cell.h"
#include "openmc/error.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -18,6 +19,7 @@ CellFilter::from_xml(pugi::xml_node node)
void
CellFilter::initialize()
{
// Convert cell IDs to indices of the global array.
for (auto& c : cells_) {
auto search = cell_map.find(c);
if (search != cell_map.end()) {
@ -29,19 +31,22 @@ CellFilter::initialize()
}
}
// Populate the index->bin map.
for (int i = 0; i < cells_.size(); i++) {
map_[cells_[i]] = i;
}
}
void
CellFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match) const
CellFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
for (int i = 0; i < p->n_coord; i++) {
auto search = map_.find(p->coord[i].cell);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}
}

View file

@ -5,13 +5,14 @@
namespace openmc {
void
CellbornFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
CellbornFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
auto search = map_.find(p->cell_born);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}

View file

@ -5,14 +5,15 @@
namespace openmc {
void
CellFromFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
CellFromFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
for (int i = 0; i < p->last_n_coord; i++) {
auto search = map_.find(p->last_cell[i]);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}
}

View file

@ -4,6 +4,7 @@
#include "openmc/error.h"
#include "openmc/geometry_aux.h" // For distribcell_path
#include "openmc/lattice.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -20,6 +21,7 @@ DistribcellFilter::from_xml(pugi::xml_node node)
void
DistribcellFilter::initialize()
{
// Convert the cell ID to an index of the global array.
auto search = cell_map.find(cell_);
if (search != cell_map.end()) {
cell_ = search->second;
@ -33,8 +35,8 @@ DistribcellFilter::initialize()
}
void
DistribcellFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
DistribcellFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
int offset = 0;
auto distribcell_index = cells[cell_]->distribcell_index_;
@ -52,8 +54,9 @@ const
}
}
if (cell_ == p->coord[i].cell) {
//TODO: off-by-one
match.bins_.push_back(offset + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
return;
}
}

View file

@ -7,6 +7,7 @@
#include "openmc/error.h"
#include "openmc/search.h"
#include "openmc/settings.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -29,7 +30,7 @@ EnergyFunctionFilter::from_xml(pugi::xml_node node)
}
void
EnergyFunctionFilter::get_all_bins(Particle* p, int estimator,
EnergyFunctionFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
if (p->last_E >= energy_.front() && p->last_E <= energy_.back()) {
@ -40,6 +41,7 @@ EnergyFunctionFilter::get_all_bins(Particle* p, int estimator,
double f = (p->last_E - energy_[i]) / (energy_[i+1] - energy_[i]);
// Interpolate on the lin-lin grid.
//TODO: off-by-one
match.bins_.push_back(1);
match.weights_.push_back((1-f) * y_[i] + f * y_[i+1]);
}

View file

@ -3,6 +3,7 @@
#include "openmc/capi.h"
#include "openmc/error.h"
#include "openmc/math_functions.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -14,12 +15,13 @@ LegendreFilter::from_xml(pugi::xml_node node)
}
void
LegendreFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
LegendreFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
double wgt[n_bins_];
calc_pn_c(order_, p->mu, wgt);
for (int i = 0; i < n_bins_; i++) {
//TODO: off-by-one
match.bins_.push_back(i + 1);
match.weights_.push_back(wgt[i]);
}

View file

@ -5,6 +5,7 @@
#include "openmc/capi.h"
#include "openmc/error.h"
#include "openmc/material.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -18,6 +19,7 @@ MaterialFilter::from_xml(pugi::xml_node node)
void
MaterialFilter::initialize()
{
// Convert material IDs to indices of the global array.
for (auto& m : materials_) {
auto search = material_map.find(m);
if (search != material_map.end()) {
@ -30,19 +32,21 @@ MaterialFilter::initialize()
}
}
// Populate the index->bin map.
for (int i = 0; i < materials_.size(); i++) {
map_[materials_[i]] = i;
}
}
void
MaterialFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
MaterialFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
auto search = map_.find(p->material - 1);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}

View file

@ -6,6 +6,7 @@
#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/mesh.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -33,12 +34,13 @@ MeshFilter::from_xml(pugi::xml_node node)
}
void
MeshFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
MeshFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const
{
if (estimator != ESTIMATOR_TRACKLENGTH) {
auto bin = meshes[mesh_]->get_bin(p->coord[0].xyz);
if (bin >= 0) {
//TODO: off-by-one
match.bins_.push_back(bin);
match.weights_.push_back(1.0);
}

View file

@ -16,8 +16,8 @@ MeshSurfaceFilter::from_xml(pugi::xml_node node)
}
void
MeshSurfaceFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
MeshSurfaceFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
meshes[mesh_]->surface_bins_crossed(p, match.bins_);
for (auto b : match.bins_) match.weights_.push_back(1.0);
@ -37,30 +37,43 @@ MeshSurfaceFilter::text_label(int bin) const
std::string out = MeshFilter::text_label(i_mesh);
// Get surface part of label.
if (i_surf == OUT_LEFT) {
out += " Outgoing, x-min";
} else if (i_surf == IN_LEFT) {
out += " Incoming, x-min";
} else if (i_surf == OUT_RIGHT) {
out += " Outgoing, x-max";
} else if (i_surf == IN_RIGHT) {
out += " Incoming, x-max";
} else if (i_surf == OUT_BACK) {
out += " Outgoing, y-min";
} else if (i_surf == IN_BACK) {
out += " Incoming, y-min";
} else if (i_surf == OUT_FRONT) {
out += " Outgoing, y-max";
} else if (i_surf == IN_FRONT) {
out += " Incoming, y-max";
} else if (i_surf == OUT_BOTTOM) {
out += " Outgoing, z-min";
} else if (i_surf == IN_BOTTOM) {
out += " Incoming, z-min";
} else if (i_surf == OUT_TOP) {
out += " Outgoing, z-max";
} else if (i_surf == IN_TOP) {
out += " Incoming, z-max";
switch (i_surf) {
case OUT_LEFT:
out += " Outgoing, x-min";
break;
case IN_LEFT:
out += " Incoming, x-min";
break;
case OUT_RIGHT:
out += " Outgoing, x-max";
break;
case IN_RIGHT:
out += " Incoming, x-max";
break;
case OUT_BACK:
out += " Outgoing, y-min";
break;
case IN_BACK:
out += " Incoming, y-min";
break;
case OUT_FRONT:
out += " Outgoing, y-max";
break;
case IN_FRONT:
out += " Incoming, y-max";
break;
case OUT_BOTTOM:
out += " Outgoing, z-min";
break;
case IN_BOTTOM:
out += " Incoming, z-min";
break;
case OUT_TOP:
out += " Outgoing, z-max";
break;
case IN_TOP:
out += " Incoming, z-max";
break;
}
return out;

View file

@ -4,6 +4,7 @@
#include "openmc/error.h"
#include "openmc/search.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -34,12 +35,14 @@ MuFilter::from_xml(pugi::xml_node node)
}
void
MuFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match) const
MuFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const
{
if (p->mu >= bins_[0] && p->mu <= bins_.back()) {
if (p->mu >= bins_.front() && p->mu <= bins_.back()) {
//TODO: off-by-one
auto bin = lower_bound_index(bins_.begin(), bins_.end(), p->mu) + 1;
match.bins_.push_back(bin);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}

View file

@ -5,6 +5,7 @@
#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/search.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -35,7 +36,7 @@ PolarFilter::from_xml(pugi::xml_node node)
}
void
PolarFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
PolarFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const
{
double theta;
@ -45,10 +46,11 @@ const
theta = std::acos(p->last_uvw[2]);
}
if (theta >= bins_[0] && theta <= bins_.back()) {
if (theta >= bins_.front() && theta <= bins_.back()) {
//TODO: off-by-one
auto bin = lower_bound_index(bins_.begin(), bins_.end(), theta) + 1;
match.bins_.push_back(bin);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}

View file

@ -3,6 +3,7 @@
#include "openmc/capi.h"
#include "openmc/error.h"
#include "openmc/math_functions.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -28,7 +29,7 @@ SphericalHarmonicsFilter::from_xml(pugi::xml_node node)
}
void
SphericalHarmonicsFilter::get_all_bins(Particle* p, int estimator,
SphericalHarmonicsFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
// Determine cosine term for scatter expansion if necessary
@ -51,6 +52,7 @@ SphericalHarmonicsFilter::get_all_bins(Particle* p, int estimator,
// Append the matching (bin,weight) for each moment
for (int i = 0; i < num_nm; i++) {
match.weights_.push_back(wgt[n] * rn[j]);
//TODO: off-by-one
match.bins_.push_back(++j);
}
}

View file

@ -3,6 +3,7 @@
#include "openmc/capi.h"
#include "openmc/error.h"
#include "openmc/math_functions.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -29,8 +30,8 @@ SpatialLegendreFilter::from_xml(pugi::xml_node node)
}
void
SpatialLegendreFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
SpatialLegendreFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
// Get the coordinate along the axis of interest.
double x;
@ -50,6 +51,7 @@ const
double wgt[order_ + 1];
calc_pn_c(order_, x_norm, wgt);
for (int i = 0; i < order_ + 1; i++) {
//TODO: off-by-one
match.bins_.push_back(i + 1);
match.weights_.push_back(wgt[i]);
}

View file

@ -4,6 +4,7 @@
#include "openmc/error.h"
#include "openmc/surface.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -17,6 +18,7 @@ SurfaceFilter::from_xml(pugi::xml_node node)
void
SurfaceFilter::initialize()
{
// Convert surface IDs to indices of the global array.
for (auto& s : surfaces_) {
auto search = surface_map.find(s);
if (search != surface_map.end()) {
@ -29,22 +31,24 @@ SurfaceFilter::initialize()
}
}
// Populate the index->bin map.
for (int i = 0; i < surfaces_.size(); i++) {
map_[surfaces_[i]] = i;
}
}
void
SurfaceFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
SurfaceFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
auto search = map_.find(std::abs(p->surface)-1);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
if (p->surface < 0) {
match.weights_.push_back(-1);
match.weights_.push_back(-1.0);
} else {
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}
}

View file

@ -4,6 +4,7 @@
#include "openmc/cell.h"
#include "openmc/error.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -17,6 +18,7 @@ UniverseFilter::from_xml(pugi::xml_node node)
void
UniverseFilter::initialize()
{
// Convert universe IDs to indices of the global array.
for (auto& u : universes_) {
auto search = universe_map.find(u);
if (search != universe_map.end()) {
@ -29,20 +31,22 @@ UniverseFilter::initialize()
}
}
// Populate the index->bin map.
for (int i = 0; i < universes_.size(); i++) {
map_[universes_[i]] = i;
}
}
void
UniverseFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
UniverseFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
for (int i = 0; i < p->n_coord; i++) {
auto search = map_.find(p->coord[i].universe);
if (search != map_.end()) {
//TODO: off-by-one
match.bins_.push_back(search->second + 1);
match.weights_.push_back(1);
match.weights_.push_back(1.0);
}
}
}

View file

@ -6,6 +6,7 @@
#include "openmc/capi.h"
#include "openmc/error.h"
#include "openmc/math_functions.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -16,16 +17,15 @@ namespace openmc {
void
ZernikeFilter::from_xml(pugi::xml_node node)
{
order_ = std::stoi(get_node_value(node, "order"));
set_order(std::stoi(get_node_value(node, "order")));
x_ = std::stod(get_node_value(node, "x"));
y_ = std::stod(get_node_value(node, "y"));
r_ = std::stod(get_node_value(node, "r"));
calc_n_bins();
}
void
ZernikeFilter::get_all_bins(Particle* p, int estimator, FilterMatch& match)
const
ZernikeFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
// Determine the normalized (r,theta) coordinates.
double x = p->coord[0].xyz[0] - x_;
@ -38,6 +38,7 @@ const
double zn[n_bins_];
calc_zn(order_, r, theta, zn);
for (int i = 0; i < n_bins_; i++) {
//TODO: off-by-one
match.bins_.push_back(i+1);
match.weights_.push_back(zn[i]);
}
@ -69,12 +70,19 @@ ZernikeFilter::text_label(int bin) const
}
}
void
ZernikeFilter::set_order(int order)
{
order_ = order;
n_bins_ = ((order+1) * (order+2)) / 2;
}
//==============================================================================
// ZernikeRadialFilter implementation
//==============================================================================
void
ZernikeRadialFilter::get_all_bins(Particle* p, int estimator,
ZernikeRadialFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
// Determine the normalized radius coordinate.
@ -87,6 +95,7 @@ ZernikeRadialFilter::get_all_bins(Particle* p, int estimator,
double zn[n_bins_];
calc_zn_rad(order_, r, zn);
for (int i = 0; i < n_bins_; i++) {
//TODO: off-by-one
match.bins_.push_back(i+1);
match.weights_.push_back(zn[i]);
}
@ -99,6 +108,13 @@ ZernikeRadialFilter::text_label(int bin) const
return "Zernike expansion, Z" + std::to_string(2*(bin-1)) + ",0";
}
void
ZernikeRadialFilter::set_order(int order)
{
order_ = order;
n_bins_ = order / 2 + 1;
}
//==============================================================================
// C-API functions
//==============================================================================
@ -116,7 +132,7 @@ openmc_zernike_filter_get_order(int32_t index, int* order)
}
auto z_filt = static_cast<ZernikeFilter*>(filt);
*order = z_filt->order_;
*order = z_filt->order();
return 0;
}
@ -153,8 +169,7 @@ openmc_zernike_filter_set_order(int32_t index, int order)
}
auto z_filt = static_cast<ZernikeFilter*>(filt);
z_filt->order_ = order;
z_filt->calc_n_bins();
z_filt->set_order(order);
filter_update_n_bins(index);
return 0;
}