3D-printing/OrcaSlicer
1
0
Fork 0
mirror of https://github.com/SoftFever/OrcaSlicer.git synced 2026-02-19 11:02:15 -07:00
Code Issues Projects Releases Packages Wiki Activity Actions 5
OrcaSlicer/tests/libslic3r/test_marchingsquares.cpp
Donovan Baarda 2a3e761ab9
Some checks are pending
Build all / Build All (push) Waiting to run
Details
Build all / Flatpak (push) Waiting to run
Details
Optimize and simplify MarchingSquares.hpp. (#10747) 
* Optimize and simplify MarchingSquares.hpp, and fix it's test.

This changes the implementation to get the possible next directions for a cell
when building the tags and clearing them as the cells are visited during the
march, instead of adding the visited previous direction to the tags during the
march. The Dir enum has been turned into bit flags that for the possible next
directions with boolean operators for testing/setting/clearing them. This
simplifies and optimizes many operations during the march and building the
polygons.

The complicated/broken and unused partial support for cell overlap has been
removed, simplifying the overly confusing grid iteration logic.

The broken test has been fixed by removing the now gone `RasterBase` namespace
from `sla::RasterBase::Pixeldim` and `sla:RasterBase:Resolution`, and the
CMakeLists.txt entry uncommented.

make Dir into flags

* Further optimize MarchingSquares.hpp and improve comments.

* Switch from a single byte-vector containing tags and dirs for each cell to a
  m_tags vector of bit-packed tags for each grid corner and an m_dirs vector
  of packed 4bit dirs for each cell. Since each grid corner tag is shared by
  the 4 adjacent cells this significantly reduces storage space and avoids
  redundantly calculating each tag 4x. It also significantly improves memory
  locality with each phase of calculating tags, calculating dirs, calculating
  rings operating only on the tags or dirs data required without them being
  interleaved with the data they don't need.

* Change NEXT_CCW to be initialized with a static constexpr lambda instead of
  a manually entered table. This avoids typo errors manually building the
  table.

* Optimize search_start_cell() so it can efficiently skip over cleared blocks
  of 8 dirs in the packed m_dirs vector.

* Change the tags logical labeling to better suit the packed tags vector data.
  This makes it a tiny bit more efficient to extract from the m_tags bitmap.

* Remove the now unused SquareTag enum class.

* Add comments explaining the algorithm, including corner-cases in cell
  iteration.

* Remove unused Dir operators and get_dirs() argument, and clang-format.

* Fix some bugs and add stream output operators for debugging.

* Fix a bug building tags where `step(gcrd, Dir::right)` was not assigned to
  update the gcrd grid point. Perhaps this should be a mutating method, or
  even a += operator? Also when wrapping at the end of a row it was updating
  the gcrd grid point by mutating the p raster point instead of itself.
  Perhaps Grid and Raster points should be different types? Maybe even
  templated?

* Fix a bug in get_tags() when the second row tags are packed into any of the
  2 LSB's of the uint32_t blocks. In hind-sight obviously `>>(o - 2)` will not
  shift left when `o < 2`.

* Move interpolation of the edge-crossings into a `interpolate()` method, and
  make it shift bottom and right side points "out" by one to account for
  raster pixel width. This makes the results track the raster shapes much more
  accurately for very small windows.

* Make `interpolate_rings()` check for and remove duplicated points. It turns
  out it's pretty common that two edge-crossing-points at a corner interpolate
  to the same point. This can also happen for the first and last points.

* For Coord add `==` and `!=` operators, and use them wherever Coord's are
  compared.

* Add `<<` stream output operators for Coord, Ring, and Dir classes. Add
  `streamtags(<stream>)` and `streamdirs(<stream>)` methods for dumping the
  tags and dirs data in an easy to understand text format. These make
  print-debugging much easier.

* Add `assert(idx < m_gridlen)` in a bunch of places where grid-indexes are
  used.

* For test_clipper_utils.cpp fix three "ambiguous overloading" compiler errors.

This just adds three `Polygons` qualifications to fix compiler errors about
ambiguous overloaded methods.

Note this file was formated with a mixture of tabs and spaces and had lots of
trailing whitespace. My editor cleaned these up resulting in a large looking
diff, but if you use `git diff -w` to ignore the whitespace changes you will
see it is actually tiny.

errros

* Update SLA/RasterToPolygons.* for MarchingSquares.hpp improvements.

Change the minimum and default window size from 2x2 to 1x1. Also remove the
strange pixel size re-scaling by (resolution/resolution-1).

The old MarchingSquares implementation had complications around a default
minimum 1 pixel "overlap" between cells which messed with the scaling a tiny
bit and meant when you requested a 2x2 window size it actually used a 1x1
window. Both of these meant you had to specify a window 1 pixel larger than
you really wanted, and you needed to undo the strange scaling artifact for
accurate dimensions of your results.

This has been fixed/removed in the new implementation, so the window is the
window, there is no overlap, and no strange miss-scaling.

* Fix test_marchingsquares.cpp and add StreamUtils.hpp.

This fixes the MarchingSquares unittests to both pass and be more strict than
they were before.

It also adds libslic3r/StreamUtils.hpp which includes some handy streaming
operators for standard libslic3r classes used to show extracted polys in the
unittests.

* Change Format/SL1.cpp to support the min 1x1 window for MarchingSquares.

* Fix the ring-walk termination condition.

Terminate the ring-walk when we return to the starting cell instead of when we
reach a cell with no remaining directions. This ensures we don't merge two
polygons if we started on an ambiguous cell.

* Revert the removal of duplicate points in interpolate_rings().

It turns out that duplicate points are only relatively common when using a 1x1
window. These happen when the line passes through the corner pixel on a
top-left corner in the raster, and the probability of this rapidly declines as
the window increases, so in many cases this filtering is just overhead. It can
also be potentially useful to see the points for every edge crossing even if
they are duplicates. This kind of filtering is already done and done better in
the polygon post-processing.

* rename `interpolate()` to `interpolate_edge()`, make it update the point
  in-place, and add asserts to ensure the input point is a valid edge
  interpolation point.
* Remove the duplicate point filtering from `interpolate_rings()` and simplify
  it.

* Optimize directions building.

This optimizes `get_dirs_block8()` to rapidly skip over blocks where the tags
produce no directions (all tags are 1's or 0's), and also to build the
directions faster when it has to by fetching the whole blocks worth of tags at
once instead of cell-by-cell.

* Rename `get_tags()` to `get_tags9()` and make it fetch a row of nine tags
  instead of the tags for a single cell.

* Optimize `get_dirs_block8()` to use `get_tags9()` to get the next nine tags
  for the current and next rows and then shift through them to generate the
  tags and directions for each cell in the block. Also abort early and just
  return an empty block if the tags are all 0's or all 1's.

* Tiny optimization for `get_tags_block32()`.

This avoids using the `step()` method for a simple step-right that can be done
with a simple increment of the column. It also avoids re-calculating the
raster-coodinates for every corner, instead incrementing the column by
`m_window.c` until the end of a row.

* Fix svg output in test_marchingsquares.cpp for recreate_object_from_rasters.

These SVG's were not properly centered...

* Fix 2 static_casts for compiling on Windows.

Thanks to RF47 for pointing this out on the #10747 pull request.

* Make edge iteration use O(ln(N)) binary search instead of linear.

This should be much faster when the window size is large.

* Make `CellIt` into a `std::random_access_iterator_tag` so that
`std::lower_bound()` can use a binary search to find the point on the edge
instead of a linear search.

* Change `step()` to support an optional distance argument and make it modify
the `Coord` in-place instead of return a new one.

* Update tests for the `step()` change.

* Add Catch2 BENCHMARK tests for MarchingSquares.

This required enabling the benchmarks in the tests/CMakeLists.txt config.

* Add a _Loop<> specialization for parallel execution using ExecutionTBB.

This is something that could be added wherever you are going to use this, but
I intend on using this in multiple places so we might as add this once in one
place where it can be reused.

* Fix whitespace in messed up by tab-replacements.

My editor renders, and replaces, tabs as 8 spaces. This messed up the
indenting in tests/libslic3r/CMakeLists.txt and
tests/libslic3r/test_clipper_utils.cpp when I made tiny changes in them.
This fixes the indenting using 4 chars. Note it will still show as a diff
because it is replacing tabs with 4 spaces, and removing trailing whitespace.

But at least it's now indented correctly...

---------

Co-authored-by: Donovan Baarda <dbaarda@google.com>
Co-authored-by: SoftFever <softfeverever@gmail.com>
2025-10-19 19:27:55 +08:00

570 lines
20 KiB
C++
Raw Blame History

#define NOMINMAX
#include <catch2/catch.hpp>
#include <test_utils.hpp>
#include <fstream>
#include <libslic3r/MarchingSquares.hpp>
#include <libslic3r/SLA/RasterToPolygons.hpp>
#include <libslic3r/SLA/AGGRaster.hpp>
#include <libslic3r/MTUtils.hpp>
#include <libslic3r/SVG.hpp>
#include <libslic3r/ClipperUtils.hpp>
#include <libslic3r/TriangleMeshSlicer.hpp>
#include <libslic3r/TriangulateWall.hpp>
#include <libslic3r/Tesselate.hpp>
#include <libslic3r/SlicesToTriangleMesh.hpp>
#include <libslic3r/StreamUtils.hpp>
using namespace Slic3r;
using namespace Catch::Matchers;
// Note this tests SLA/RasterToPolygons.hpp, SLA/AGGRaster.hpp, and
// ClipperUtils.hpp at least as much as MarchingSquares.hpp.
// Get the Point corresponding to a raster column and row.
Point rstPoint(const sla::RasterGrayscaleAA& rst, const size_t c, const size_t r)
{
size_t rows = rst.resolution().height_px, cols = rst.resolution().width_px;
auto pxd = rst.pixel_dimensions();
auto tr = rst.trafo();
coord_t width = scaled(cols * pxd.h_mm), height = scaled(rows * pxd.w_mm);
Point p = Point::new_scale(c * pxd.w_mm, r * pxd.h_mm);
// reverse the raster transformations
if (tr.mirror_y)
p.y() = height - p.y();
if (tr.mirror_x)
p.x() = width - p.x();
p.x() -= tr.center_x;
p.y() -= tr.center_y;
if (tr.flipXY)
std::swap(p.x(), p.y());
return p;
}
// Get the size of a raster pixel in coord_t.
static Point rstPixel(const sla::RasterGrayscaleAA& rst)
{
auto pxd = rst.pixel_dimensions();
return Point::new_scale(pxd.w_mm, pxd.h_mm);
}
// Get the size of a raster in coord_t.
static Point rstSize(const sla::RasterGrayscaleAA& rst)
{
auto pxd = rst.pixel_dimensions();
auto res = rst.resolution();
return Point::new_scale(pxd.w_mm * res.width_px, pxd.h_mm * res.height_px);
}
// Get the bounding box of a raster.
static BoundingBox rstBBox(const sla::RasterGrayscaleAA& rst)
{
auto center = rst.trafo().get_center();
return BoundingBox(Point(0, 0) - center, rstSize(rst) - center);
}
// Get the ExPolygons directly corresponding to a raster.
static ExPolygons rstGetPolys(sla::RasterGrayscaleAA& rst)
{
size_t rows = rst.resolution().height_px, cols = rst.resolution().width_px;
Polygons polys;
for (auto r = 0; r < rows; r++) {
// use c0==cols as a sentinel marker for "no start column yet".
size_t c0 = cols;
for (auto c = 0; c <= cols; c++) {
if (c < cols && rst.read_pixel(c, r) > 128) {
// We have set pixels, set the c0 start column if it is not yet set.
if (c0 == cols)
c0 = c;
} else if (c0 < cols) {
// There is no pixel set, but we do have a c0 start column. Output a
// "row-rectangle" poly for this row between the start column c0 and
// the current column.
polys.push_back({rstPoint(rst, c0, r), rstPoint(rst, c0, r + 1), rstPoint(rst, c, r + 1), rstPoint(rst, c, r)});
// Make sure the poly is anti-clockwise, which it might not be
// depending on how rstPoint() reverses the raster transformations
// from (c,r) to (x,y) coordinates.
if (polys.back().is_clockwise())
polys.back().reverse();
// Clear the start column c0 for the next row-rectangle.
c0 = cols;
}
}
}
// Merge all the row-rectangle polys into contiguous raster ExPolygons.
return union_ex(polys);
}
// Get the length in mm of a "vector" Point.
static double len(const Point& v) { return unscaled(v.norm()); }
// Get the area in mm^2 of a box with corners at the origin and a Point.
static double area(const Point& v) { return unscaled(v.x()) * unscaled(v.y()); }
// Find the index of the nearest extracted ExPolygon for a reference ExPolygon.
static int find_closest_ext(const ExPolygons& exts, ExPolygon ref)
{
auto ref_center = ref.contour.bounding_box().center();
auto closest = std::min_element(exts.begin(), exts.end(), [&ref_center](auto a, auto b) {
auto a_center = a.contour.bounding_box().center();
auto b_center = b.contour.bounding_box().center();
return a_center.distance_to(ref_center) < b_center.distance_to(ref_center);
});
return std::distance(exts.begin(), closest);
}
static Slic3r::sla::RasterGrayscaleAA create_raster(const sla::Resolution& res, double disp_w = 100., double disp_h = 100.)
{
sla::PixelDim pixdim{disp_w / res.width_px, disp_h / res.height_px};
auto bb = BoundingBox({0, 0}, {scaled(disp_w), scaled(disp_h)});
sla::RasterBase::Trafo trafo;
trafo.center_x = bb.center().x();
trafo.center_y = bb.center().y();
return sla::RasterGrayscaleAA{res, pixdim, trafo, agg::gamma_threshold(.5)};
}
static ExPolygon square(double a, Point center = {0, 0})
{
ExPolygon poly;
coord_t V = scaled(a / 2.);
poly.contour.points = {{-V, -V}, {V, -V}, {V, V}, {-V, V}};
poly.translate(center.x(), center.y());
return poly;
}
static ExPolygon square_with_hole(double a, Point center = {0, 0})
{
ExPolygon poly = square(a);
poly.holes.emplace_back();
coord_t V = scaled(a / 4.);
poly.holes.front().points = {{-V, V}, {V, V}, {V, -V}, {-V, -V}};
poly.translate(center.x(), center.y());
return poly;
}
static ExPolygons circle_with_hole(double r, Point center = {0, 0})
{
ExPolygon poly;
std::vector<double> pis = linspace_vector(0., 2 * PI, 100);
coord_t rs = scaled(r);
for (double phi : pis) {
poly.contour.points.emplace_back(rs * std::cos(phi), rs * std::sin(phi));
}
poly.holes.emplace_back(poly.contour);
poly.holes.front().reverse();
for (auto& p : poly.holes.front().points)
p /= 2;
poly.translate(center.x(), center.y());
return {poly};
}
static const Vec2i32 W2x2 = {2, 2};
static const Vec2i32 W1x1 = {1, 1};
template<class Rst>
static void test_expolys(Rst&& rst, const ExPolygons& ref, Vec2i32 window, const std::string& name = "test", bool strict = true)
{
auto raster_bb = rstBBox(rst);
Point pixel_size = rstPixel(rst);
Point window_size{coord_t(pixel_size.x() * window.x()), coord_t(pixel_size.y() * window.y())};
double pixel_area = area(pixel_size);
double pixel_len = len(pixel_size);
double window_area = area(window_size);
double window_len = len(window_size);
for (const ExPolygon& expoly : ref)
rst.draw(expoly);
std::fstream out(name + ".png", std::ios::out);
out << rst.encode(sla::PNGRasterEncoder{});
out.close();
const ExPolygons bmp = rstGetPolys(rst);
const ExPolygons ext = sla::raster_to_polygons(rst, window);
SVG svg(name + ".svg", raster_bb);
svg.draw(bmp, "green");
if (pixel_size.x() >= scale_(0.5))
svg.draw_grid(raster_bb, "grey", scale_(0.05), pixel_size.x());
if (window_size.x() >= scale_(1.0))
svg.draw_grid(raster_bb, "grey", scale_(0.10), window_size.x());
svg.draw_outline(ref, "red", "red", scale_(0.3));
svg.draw_outline(ext, "blue", "blue");
svg.Close();
// Note all these areas are unscaled back to mm^2.
double raster_area = unscaled(unscaled(area(bmp)));
double reference_area = unscaled(unscaled(area(ref)));
double extracted_area = unscaled(unscaled(area(ext)));
// Note that errors accumulate with each step going from the reference
// polys to the extracted polys. The rendering of the reference polys to
// the raster does introduce pixelization errors too. This checks for
// acceptable errors going from reference to raster, and raster to
// reference.
for (size_t i = 0; i < ref.size(); ++i) {
if (ref[i].contour.size() < 20)
UNSCOPED_INFO("reference ref[" << i << "]: " << ref[i]);
}
CHECK_THAT(raster_area, WithinRel(reference_area, pixel_len * 0.05) || WithinAbs(reference_area, pixel_area));
for (size_t i = 0; i < ext.size(); ++i) {
if (ext[i].contour.size() < 20)
UNSCOPED_INFO("extracted ext[" << i << "]: " << ext[i]);
}
CHECK_THAT(extracted_area, WithinRel(raster_area, 0.05) || WithinAbs(raster_area, window_area));
for (auto i = 0; i < ext.size(); ++i) {
CHECK(ext[i].contour.is_counter_clockwise());
for (auto& h : ext[i].holes)
CHECK(h.is_clockwise());
}
BoundingBox ref_bb;
for (auto& expoly : ref)
ref_bb.merge(expoly.contour.bounding_box());
BoundingBox ext_bb;
for (auto& expoly : ext)
ext_bb.merge(expoly.contour.bounding_box());
CHECK(len(ext_bb.center() - ref_bb.center()) < pixel_len);
// In ambigous cases (when polygons just touch) there are multiple equally
// valid interpretations of the raster into polygons. Although
// MarchingSquares currently systematically selects the solution that
// breaks them into separate polygons, that might not always be true. Also,
// SLA/RasterToPolygons.hpp, and in particular union_ex() from
// ClipperUtils.hpp that it uses, can and does sometimes merge them back
// together. This means we cannot reliably make assertions about the
// extracted number of polygons and their shapes in these cases. So we skip
// the individual polygon checks for strict=false.
if (strict) {
CHECK(ext.size() == ref.size());
for (auto i = 0; i < ext.size(); ++i) {
auto j = find_closest_ext(ref, ext[i]);
INFO("Comparing ext[" << i << "] against closest ref[" << j << "]");
CHECK(ext[i].holes.size() == ref[j].holes.size());
double ext_i_area = unscaled(unscaled(ext[i].area()));
double ref_j_area = unscaled(unscaled(ref[j].area()));
CHECK_THAT(ext_i_area, WithinRel(ref_j_area, pixel_len * 0.05) || WithinAbs(ref_j_area, window_area));
auto ext_i_bb = ext[i].contour.bounding_box();
auto ref_j_bb = ref[j].contour.bounding_box();
CHECK(len(ext_i_bb.center() - ref_j_bb.center()) < pixel_len);
}
}
}
TEST_CASE("Empty raster should result in empty polygons", "[MarchingSquares]")
{
sla::RasterGrayscaleAAGammaPower rst{{}, {}, {}};
ExPolygons extracted = sla::raster_to_polygons(rst);
REQUIRE(extracted.size() == 0);
}
TEST_CASE("Marching squares directions", "[MarchingSquares]")
{
using namespace marchsq;
Coord crd{0, 0};
__impl::step(crd, __impl::Dir::left);
CHECK(crd == Coord(0, -1));
__impl::step(crd, __impl::Dir::down);
CHECK(crd == Coord(1, -1));
__impl::step(crd, __impl::Dir::right);
CHECK(crd == Coord(1, 0));
__impl::step(crd, __impl::Dir::up);
CHECK(crd == Coord(0, 0));
__impl::step(crd, __impl::Dir::left, 7);
CHECK(crd == Coord(0, -7));
__impl::step(crd, __impl::Dir::down, 7);
CHECK(crd == Coord(7, -7));
__impl::step(crd, __impl::Dir::right, 7);
CHECK(crd == Coord(7, 0));
__impl::step(crd, __impl::Dir::up, 7);
CHECK(crd == Coord(0, 0));
__impl::step(crd, __impl::Dir::left, -3);
CHECK(crd == Coord(0, 3));
__impl::step(crd, __impl::Dir::down, -3);
CHECK(crd == Coord(-3, 3));
__impl::step(crd, __impl::Dir::right, -3);
CHECK(crd == Coord(-3, 0));
__impl::step(crd, __impl::Dir::up, -3);
CHECK(crd == Coord(0, 0));
}
TEST_CASE("Fully covered raster should result in a rectangle", "[MarchingSquares]")
{
auto rst = create_raster({4, 4}, 4., 4.);
ExPolygon rect = square(4);
SECTION("Full accuracy") { test_expolys(rst, {rect}, W1x1, "fully_covered_full_acc"); }
SECTION("Half accuracy") { test_expolys(rst, {rect}, W2x2, "fully_covered_half_acc"); }
}
TEST_CASE("4x4 raster with one ring", "[MarchingSquares]")
{
sla::PixelDim pixdim{1, 1};
// We need one additional row and column to detect edges
sla::RasterGrayscaleAA rst{{4, 4}, pixdim, {}, agg::gamma_threshold(.5)};
ExPolygons one = {{{1, 1}, {3, 1}, {3, 3}, {2, 3}, {2, 2}, {1, 2}}};
for (ExPolygon& p : one)
p.scale(scaled(1.0));
test_expolys(rst, one, W1x1, "one_4x4");
}
TEST_CASE("10x10 raster with two rings", "[MarchingSquares]")
{
sla::PixelDim pixdim{1, 1};
// We need one additional row and column to detect edges
sla::RasterGrayscaleAA rst{{10, 10}, pixdim, {}, agg::gamma_threshold(.5)};
SECTION("Ambiguous case with 'bd' square")
{
ExPolygons ac = {{{1, 1}, {3, 1}, {3, 2}, {2, 2}, {2, 3}, {1, 3}}, {{4, 4}, {2, 4}, {2, 3}, {3, 3}, {3, 2}, {4, 2}}};
for (ExPolygon& p : ac)
p.scale(scaled(2.0));
test_expolys(rst, ac, W1x1, "bd_10x10", false);
}
SECTION("Ambiguous case with 'ac' square")
{
ExPolygons bd = {{{1, 4}, {1, 2}, {2, 2}, {2, 3}, {3, 3}, {3, 4}}, {{4, 1}, {4, 3}, {3, 3}, {3, 2}, {2, 2}, {2, 1}}};
for (ExPolygon& p : bd)
p.scale(scaled(2.0));
test_expolys(rst, bd, W1x1, "ac_10x10", false);
}
}
TEST_CASE("Square with hole in the middle", "[MarchingSquares]")
{
using namespace Slic3r;
ExPolygons inp = {square_with_hole(50.)};
SECTION("Proportional raster, 1x1 mm pixel size, full accuracy")
{
test_expolys(create_raster({100, 100}, 100., 100.), inp, W1x1, "square_with_hole_proportional_1x1_mm_px_full");
}
SECTION("Proportional raster, 1x1 mm pixel size, half accuracy")
{
test_expolys(create_raster({100, 100}, 100., 100.), inp, W2x2, "square_with_hole_proportional_1x1_mm_px_half");
}
SECTION("Landscape raster, 1x1 mm pixel size, full accuracy")
{
test_expolys(create_raster({150, 100}, 150., 100.), inp, W1x1, "square_with_hole_landsc_1x1_mm_px_full");
}
SECTION("Landscape raster, 1x1 mm pixel size, half accuracy")
{
test_expolys(create_raster({150, 100}, 150., 100.), inp, W2x2, "square_with_hole_landsc_1x1_mm_px_half");
}
SECTION("Portrait raster, 1x1 mm pixel size, full accuracy")
{
test_expolys(create_raster({100, 150}, 100., 150.), inp, W1x1, "square_with_hole_portrait_1x1_mm_px_full");
}
SECTION("Portrait raster, 1x1 mm pixel size, half accuracy")
{
test_expolys(create_raster({100, 150}, 100., 150.), inp, W2x2, "square_with_hole_portrait_1x1_mm_px_half");
}
SECTION("Proportional raster, 2x2 mm pixel size, full accuracy")
{
test_expolys(create_raster({50, 50}, 100., 100.), inp, W1x1, "square_with_hole_proportional_2x2_mm_px_full");
}
SECTION("Proportional raster, 2x2 mm pixel size, half accuracy")
{
test_expolys(create_raster({50, 50}, 100., 100.), inp, W2x2, "square_with_hole_proportional_2x2_mm_px_half");
}
SECTION("Proportional raster, 0.5x0.5 mm pixel size, full accuracy")
{
test_expolys(create_raster({200, 200}, 100., 100.), inp, W1x1, "square_with_hole_proportional_0.5x0.5_mm_px_full");
}
SECTION("Proportional raster, 0.5x0.5 mm pixel size, half accuracy")
{
test_expolys(create_raster({200, 200}, 100., 100.), inp, W2x2, "square_with_hole_proportional_0.5x0.5_mm_px_half");
}
}
TEST_CASE("Circle with hole in the middle", "[MarchingSquares]")
{
using namespace Slic3r;
test_expolys(create_raster({1000, 1000}), circle_with_hole(25.), W1x1, "circle_with_hole");
}
static void recreate_object_from_rasters(const std::string& objname, float lh)
{
TriangleMesh mesh = load_model(objname);
auto bb = mesh.bounding_box();
Vec3f tr = -bb.center().cast<float>();
mesh.translate(tr.x(), tr.y(), tr.z());
bb = mesh.bounding_box();
std::vector<ExPolygons> layers = slice_mesh_ex(mesh.its, grid(float(bb.min.z()) + lh, float(bb.max.z()), lh));
sla::Resolution res{2560, 1440};
double disp_w = 120.96;
double disp_h = 68.04;
#ifndef NDEBUG
size_t cntr = 0;
#endif
for (ExPolygons& layer : layers) {
auto rst = create_raster(res, disp_w, disp_h);
for (ExPolygon& island : layer) {
rst.draw(island);
}
#ifndef NDEBUG
std::fstream out(objname + std::to_string(cntr) + ".png", std::ios::out);
out << rst.encode(sla::PNGRasterEncoder{});
out.close();
#endif
ExPolygons layer_ = sla::raster_to_polygons(rst);
// float delta = scaled(std::min(rst.pixel_dimensions().h_mm,
// rst.pixel_dimensions().w_mm)) / 2;
// layer_ = expolygons_simplify(layer_, delta);
#ifndef NDEBUG
SVG svg(objname + std::to_string(cntr) + ".svg", rstBBox(rst));
svg.draw(layer_);
svg.draw(layer, "green");
svg.Close();
#endif
double layera = 0., layera_ = 0.;
for (auto& p : layer)
layera += p.area();
for (auto& p : layer_)
layera_ += p.area();
#ifndef NDEBUG
std::cout << cntr++ << std::endl;
#endif
double diff = std::abs(layera_ - layera);
REQUIRE((diff <= 0.1 * layera || diff < scaled<double>(1.) * scaled<double>(1.)));
layer = std::move(layer_);
}
indexed_triangle_set out = slices_to_mesh(layers, bb.min.z(), double(lh), double(lh));
its_write_obj(out, "out_from_rasters.obj");
}
TEST_CASE("Recreate object from rasters", "[SL1Import]") { recreate_object_from_rasters("frog_legs.obj", 0.05f); }
namespace marchsq {
static constexpr float layerf = 0.20; // layer height in mm (used for z values).
static constexpr float gsizef = 100.0; // grid size in mm (box volume side length).
static constexpr float wsizef = 0.50; // grid window size in mm (roughly line segment length).
static constexpr float psizef = 0.01; // raster pixel size in mm (roughly point accuracy).
static constexpr float isoval = 0.0; // iso value threshold to use.
static constexpr size_t wsize = std::round(wsizef / psizef);
static float period = 10.0; // gyroid "wavelength" in mm (2x line spacing).
static float freq = 2 * PI / period; // gyroid frequency in waves per mm.
void set_period(float len = 10.0)
{
period = len;
freq = 2 * PI / period;
}
static size_t layer_n;
static size_t ring_n;
static size_t point_n;
static size_t get_n;
void reset_stats()
{
layer_n = 0;
ring_n = 0;
point_n = 0;
get_n = 0;
}
using Rings = std::vector<Ring>;
template<> struct _RasterTraits<size_t>
{
// using Rst = Slic3r::sla::RasterGrayscaleAA;
// The type of pixel cell in the raster
using ValueType = float;
// Value at a given position
static float get(const size_t& layer, size_t row, size_t col)
{
get_n++;
const float x = col * psizef * freq;
const float y = row * psizef * freq;
const float z = layer * psizef * freq;
return sinf(x) * cosf(y) + sinf(y) * cosf(z) + sinf(z) * cosf(x);
}
// Number of rows and cols of the raster
static size_t rows(const size_t& layer) { return std::round(gsizef / psizef); }
static size_t cols(const size_t& layer) { return std::round(gsizef / psizef); }
};
Rings get_gyroids(size_t l)
{
size_t layer = l;
Rings rings = execute(layer, isoval, {wsize, wsize});
layer_n++;
ring_n += rings.size();
for (auto r : rings)
point_n += r.size();
return rings;
}
}; // namespace marchsq
void benchmark_gyroid(float period)
{
marchsq::reset_stats();
marchsq::set_period(period);
INFO("grid size: " << marchsq::gsizef << "mm\nlayer height: " << marchsq::layerf << "mm\n");
INFO("window size: " << marchsq::wsizef << "mm\npoint size: " << marchsq::psizef << "mm\n");
INFO("gyroid period: " << marchsq::period << "mm\n");
BENCHMARK("indexed", i) { return marchsq::get_gyroids(i); };
INFO("output avg rings/layer: " << float(marchsq::ring_n) / float(marchsq::layer_n) << "\n");
INFO("output avg points/layer: " << float(marchsq::point_n) / float(marchsq::layer_n) << "\n");
INFO("output avg gets/layer: " << float(marchsq::get_n) / float(marchsq::layer_n) << "\n");
REQUIRE(marchsq::layer_n > 0);
}
TEST_CASE("Benchmark gyroid cube period 10.0mm", "[MarchingSquares]") { benchmark_gyroid(10.0); }
TEST_CASE("Benchmark gyroid cube period 5.0mm", "[MarchingSquares]") { benchmark_gyroid(5.0); }
Reference in a new issue View git blame Copy permalink