7#include <hrleSparseIterator.hpp>
11using namespace viennacore;
17 using hrleIndexType = viennahrle::IndexType;
18 SmartPointer<Domain<T, 2>> inputLevelSet =
nullptr;
19 SmartPointer<Domain<T, 3>> outputLevelSet =
nullptr;
20 Vec2D<T> extent = {0., 0.};
22 std::array<BoundaryConditionEnum, 3> boundaryConds = {};
28 SmartPointer<
Domain<T, 3>> passedOutputLS, Vec2D<T> passedExtent,
29 const int passedExtrudeDim = 2,
31 BoundaryConditionEnum::INFINITE_BOUNDARY)
32 : inputLevelSet(passedInputLS), outputLevelSet(passedOutputLS),
33 extent(passedExtent), extrudeDim(passedExtrudeDim),
34 boundaryConds{{passedInputLS->getGrid().getBoundaryConditions(0),
35 passedInputLS->getGrid().getBoundaryConditions(1),
36 passedBoundaryCond}} {}
39 SmartPointer<
Domain<T, 3>> passedOutputLS, Vec2D<T> passedExtent,
40 const int passedExtrudeDim,
41 std::array<BoundaryConditionEnum, 3> passedBoundaryConds)
42 : inputLevelSet(passedInputLS), outputLevelSet(passedOutputLS),
43 extent(passedExtent), extrudeDim(passedExtrudeDim),
44 boundaryConds(passedBoundaryConds) {}
47 inputLevelSet = passedInputLS;
52 outputLevelSet = passedOutputLS;
56 void setExtent(Vec2D<T> passedExtent) { extent = passedExtent; }
60 extrudeDim = passedExtrudeDim;
64 std::array<BoundaryConditionEnum, 3> passedBoundaryConds) {
65 boundaryConds = passedBoundaryConds;
69 for (
int i = 0; i < 3; i++)
70 boundaryConds[i] = passedBoundaryConds[i];
74 if (inputLevelSet ==
nullptr) {
76 .addWarning(
"No input Level Set supplied to Extrude! Not converting.")
79 if (outputLevelSet ==
nullptr) {
82 "No output Level Set supplied to Extrude! Not converting.")
88 const auto extrudeDims = getExtrudeDims();
90 std::vector<std::pair<viennahrle::Index<3>,
T>> points3D;
92 auto &domain2D = inputLevelSet->getDomain();
93 auto &grid2D = inputLevelSet->getGrid();
94 const T gridDelta = grid2D.getGridDelta();
95 auto minBounds = grid2D.getMinBounds();
96 auto maxBounds = grid2D.getMaxBounds();
98 double domainBounds[2 * 3];
99 domainBounds[2 * extrudeDim] = extent[0];
100 domainBounds[2 * extrudeDim + 1] = extent[1];
101 domainBounds[2 * extrudeDims[0]] = gridDelta * minBounds[0];
102 domainBounds[2 * extrudeDims[0] + 1] = gridDelta * maxBounds[0];
103 domainBounds[2 * extrudeDims[1]] = gridDelta * minBounds[1];
104 domainBounds[2 * extrudeDims[1] + 1] = gridDelta * maxBounds[1];
106 auto tmpLevelSet = SmartPointer<Domain<T, 3>>::New(
107 domainBounds, boundaryConds.data(), gridDelta);
108 outputLevelSet->deepCopy(tmpLevelSet);
112 !it.isFinished(); ++it) {
116 const auto index2D = it.getStartIndices();
117 const T value = it.getValue();
119 const hrleIndexType zStart = std::floor(extent[0] / gridDelta) - 1;
120 const hrleIndexType zEnd = std::ceil(extent[1] / gridDelta) + 1;
121 for (hrleIndexType z = zStart; z <= zEnd; ++z) {
122 viennahrle::Index<3> index3D;
123 index3D[0] = index2D[0];
124 index3D[1] = index2D[1];
127 points3D.emplace_back(index3D, value);
131 outputLevelSet->insertPoints(points3D);
132 outputLevelSet->finalize(2);
136 inline std::array<unsigned, 2> getExtrudeDims()
const {
137 assert(extrudeDim < 3);
138 if (extrudeDim == 0) {
140 }
else if (extrudeDim == 1) {
Class containing all information about the level set, including the dimensions of the domain,...
Definition lsDomain.hpp:27
viennahrle::Domain< T, D > DomainType
Definition lsDomain.hpp:32
void setBoundaryConditions(BoundaryConditionEnum passedBoundaryConds[3])
Definition lsExtrude.hpp:68
void setExtent(Vec2D< T > passedExtent)
Definition lsExtrude.hpp:56
void setBoundaryConditions(std::array< BoundaryConditionEnum, 3 > passedBoundaryConds)
Definition lsExtrude.hpp:63
Extrude(SmartPointer< Domain< T, 2 > > passedInputLS, SmartPointer< Domain< T, 3 > > passedOutputLS, Vec2D< T > passedExtent, const int passedExtrudeDim, std::array< BoundaryConditionEnum, 3 > passedBoundaryConds)
Definition lsExtrude.hpp:38
void apply()
Definition lsExtrude.hpp:73
Extrude(SmartPointer< Domain< T, 2 > > passedInputLS, SmartPointer< Domain< T, 3 > > passedOutputLS, Vec2D< T > passedExtent, const int passedExtrudeDim=2, BoundaryConditionEnum passedBoundaryCond=BoundaryConditionEnum::INFINITE_BOUNDARY)
Definition lsExtrude.hpp:27
void setExtrudeDimension(const int passedExtrudeDim)
Definition lsExtrude.hpp:59
void setOutputLevelSet(SmartPointer< Domain< T, 3 > > &passedOutputLS)
Definition lsExtrude.hpp:51
void setInputLevelSet(SmartPointer< Domain< T, 2 > > passedInputLS)
Definition lsExtrude.hpp:46
Definition lsAdvect.hpp:36
viennahrle::BoundaryType BoundaryConditionEnum
Definition lsDomain.hpp:23
double T
Definition pyWrap.cpp:68