Обрезать часть многоугольника (кольца) из многоугольника (кольца) по геометрии буста - PullRequest
0 голосов
/ 25 мая 2018

Попытка обрезать кольцо от кольца снизу вверх (50%), но результаты не работают, как я ожидал.

Мое решение

using bg_point_type = boost::geometry::model::d2::point_xy<double>;
using bg_polygon_type = boost::geometry::model::ring<bg_point_type>;
//Find the envelope of the region of the object
bg_polygon_type obj_region = //... it is a rectangle
boost::geometry::model::box<bg_point_type> envelope_box;
boost::geometry::envelope(obj_region, envelope_box);

auto const top_y = envelope_box.min_corner().y();
auto const bottom_y = envelope_box.max_corner().y();
auto const left = envelope_box.min_corner().x();
auto const right = envelope_box.max_corner().x();
//with min_corner and max_corner, we can know the top of the ring
auto const line_y = bottom_y - (bottom_y - top_y) * ratio;

//find the intersection of the line and the polygon
boost::geometry::model::linestring<bg_point_type> const line{{left, line_y}, {right, line_y}};
std::vector<bg_point_type> intersect;
boost::geometry::intersection(obj_region, line, intersect);

//remove those points higher than top
bg_polygon_type result = obj_region;
auto it = std::remove_if(std::begin(result ), std::end(result ), [=](auto const &pt)
{
    return pt.y() < line_y;
});
result .erase(it, std::end(result ));  
//insert the intersect points into the results
std::copy(std::begin(intersect), std::end(intersect), std::back_inserter(result));    
boost::geometry::correct(result);

Это решение дает мне нужную точку, но оно не дает мне нужную форму.

С точками, равными "0,1, 0,1,0.1,0.5,0.5,0.5,0.5,0.1,0.1,0.1 ", у меня есть такое кольцо (для простоты я выбираю прямоугольник в качестве примера)

enter image description here

После того, как я обрежу кольцо, я надеюсь, что это может быть (с точками, 0,1,0.5,0.1,0.3,0.5,0.3,0.5,0.5,0.1,0.5)

enter image description here

Но это дало мне 0,1,0,5,0,5,0,5,0,1,0.3,0.5,0.3,0.1,0.5

enter image description here

Проблема в том, что «правильный» многоугольник не является уникальным, как я могу получить желаемые результаты, если он слишком длинный / слишком сложный, чтобы сказать здесь, какие ключевые слова мне следует искать?Я пробовал "обрезать многоугольник", "обрезать часть многоугольника", "обрезать многоугольник из многоугольника" и т. Д. Спасибо

Редактировать 1: obj_region - это простой многоугольник

Редактировать 2: ПроблемаРешение sehe

Я изменяю функцию crop_box на кадрирование снизу вверх и ожидаю, что api пересечений даст мне тот же результат, но вместо этого он выдаст странные результаты.

template <typename G>
bg_box_type crop_box(G const& geom, double ratio) {
    bg_box_type env;
    bg::envelope(geom, env);

    auto const miny = env.min_corner().y();
    auto const height = env.max_corner().y() - miny;

    env.min_corner().set<1>(env.max_corner().y() - height*ratio);
    return env;
}

Input:POLYGON ((100,0 100,0 100,0 200,0,200,0 200,0,200,0 100,0 100,0 100,0)) Выходные данные: POLYGON ((100 200 100 150 200 200 200 150 150 100 200)) - недействительный Выход после правильного: POLYGON ((100 200 100 150 200 200 200 150 150 100 200))--invalid

полные исходные коды

#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/multi_polygon.hpp>
#include <boost/geometry/geometries/linestring.hpp>
#include <boost/geometry/geometries/multi_polygon.hpp>
#include <boost/geometry/core/cs.hpp>
#include <boost/geometry/arithmetic/arithmetic.hpp>
#include <boost/geometry/algorithms/for_each.hpp>
#include <boost/geometry/algorithms/envelope.hpp>
#include <boost/geometry/algorithms/intersection.hpp>
#include <boost/geometry/io/io.hpp>
#include <iostream>
#include <fstream>

namespace bg = boost::geometry;

using bg_point_type = bg::model::d2::point_xy<double>;
using bg_polygon_type = bg::model::ring<bg_point_type>;
using bg_box_type = bg::model::box<bg_point_type>;

template <typename G>
bg_box_type crop_box(G const& geom, double ratio) {
    bg_box_type env;
    bg::envelope(geom, env);

    auto const miny = env.min_corner().y();
    auto const height = env.max_corner().y() - miny;

    env.min_corner().set<1>(env.max_corner().y() - height*ratio);
    return env;
}

template <typename G>
void diags(std::string name, G& geom) {
    std::cout << name << ": " << bg::wkt(geom) << "\n";

    std::string reason;
    if (!bg::is_valid(geom, reason)) {
        std::cout << name << ": " << reason << "\n";

        bg::correct(geom);

        std::cout << bg::wkt(geom) << "\n";
        if (!bg::is_valid(geom, reason)) {
            std::cout << name << " corrected: " << reason << "\n";
        }
    }
}

int main() {
    bool const visualize = true;

    //Find the envelope of the region of the object
    bg_polygon_type obj_region;
    bg::read_wkt("POLYGON((100.0 100.0,100.0 200.0,200.0 200.0,200.0 100.0,100.0 100.0))", obj_region);

    diags("Input", obj_region);

    bg_polygon_type out;
    bg::intersection(crop_box(obj_region, 0.5), obj_region, out);

    diags("output", out);

    std::cout << "Output: " << bg::wkt(out) << "\n";

    if (visualize) {
        std::ofstream svg("svg.svg");
        boost::geometry::svg_mapper<bg_point_type> mapper(svg, 600, 600);
        mapper.add(obj_region);
        mapper.add(out);

        mapper.map(obj_region, "fill-opacity:0.5;fill:rgb(204,153,0);stroke:rgb(204,153,0);stroke-width:2");
        mapper.map(out, "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
    }
}

enter image description here

1 Ответ

0 голосов
/ 25 мая 2018

Вы просто хотите пересечь оригинал region_obj с нижней половиной его envelope_box.

Вместо этого, находя только верхнюю линию пересечения и принудительно удаляяуказывает на случайно внутри полигона, вы создавали недопустимый полигон:

Live On Coliru

#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/linestring.hpp>
#include <boost/geometry/geometries/multi_polygon.hpp>
#include <boost/geometry/algorithms/envelope.hpp>
#include <boost/geometry/algorithms/intersection.hpp>
#include <boost/geometry/io/io.hpp>
#include <iostream>

namespace bg = boost::geometry;

using bg_point_type = bg::model::d2::point_xy<double>;
using bg_polygon_type = bg::model::ring<bg_point_type>;

template <typename G>
void diags(std::string name, G& geom) {
    std::cout << name << ": " << bg::wkt(geom) << "\n";

    std::string reason;
    if (!bg::is_valid(geom, reason)) {
        std::cout << name << " invalid: " << reason << "\n";

        bg::correct(geom);

        std::cout << name << " corrected: " << bg::wkt(geom) << "\n";
        if (!bg::is_valid(geom, reason)) {
            std::cout << name << " invalid: " << reason << "\n";
        }
    }
}

int main() {
    double const ratio = 0.5;

    //Find the envelope of the region of the object
    bg_polygon_type obj_region;
    bg::read_wkt("POLYGON((0.1 0.1,0.1 0.5,0.5 0.5,0.5 0.1,0.1 0.1))", obj_region);

    diags("Input", obj_region);

    bg::model::box<bg_point_type> envelope_box;
    bg::envelope(obj_region, envelope_box);

    auto const top_y    = envelope_box.min_corner().y();
    auto const bottom_y = envelope_box.max_corner().y();
    auto const left     = envelope_box.min_corner().x();
    auto const right    = envelope_box.max_corner().x();

    //with min_corner and max_corner, we can know the top of the ring
    auto const line_y   = bottom_y - (bottom_y - top_y) * ratio;

    //find the intersection of the line and the polygon
    bg::model::linestring<bg_point_type> const line{{left, line_y}, {right, line_y}};
    bg::model::multi_point<bg_point_type> intersect;
    bg::intersection(obj_region, line, intersect);

    std::cout << bg::wkt(intersect) << "\n";

    //remove those points higher than top
    bg_polygon_type result = obj_region;
    auto it = std::remove_if(std::begin(result ), std::end(result ), [=](auto const &pt) { return pt.y() < line_y; });
    result.erase(it, std::end(result ));  

    //insert the intersect points into the results
    std::copy(std::begin(intersect), std::end(intersect), std::back_inserter(result));    
    diags("Result", result);
}

Печать

Input: POLYGON((0.1 0.1,0.1 0.5,0.5 0.5,0.5 0.1,0.1 0.1))
MULTIPOINT((0.1 0.3),(0.5 0.3))
Result: POLYGON((0.1 0.5,0.5 0.5,0.1 0.3,0.5 0.3))
Result invalid: Geometry is defined as closed but is open
Result corrected: POLYGON((0.1 0.5,0.5 0.5,0.1 0.3,0.5 0.3,0.1 0.5))
Result invalid: Geometry has invalid self-intersections. A self-intersection point was found at (0.3, 0.4); method: i; operations: u/i; segment IDs {source, multi, ring, segment}: {0, -1, -1, 1}/{0, -1, -1, 3}

Исправление

Просто спросив, что вы имеете в виду:

bg_polygon_type obj_region;
bg::read_wkt("POLYGON((0.1 0.1,0.1 0.5,0.5 0.5,0.5 0.1,0.1 0.1))", obj_region);

diags("Input", obj_region); // check validity/attempt correction

bg_polygon_type out;
bg::intersection(crop_box(obj_region, 0.5), obj_region, out);

std::cout << "Output: " << bg::wkt(out) << "\n";

Конечно, определите crop_box аналогично тому, как вы это делали (только, получите конверт обрезки, а не многоточечный):

template <typename G>
bg_box_type crop_box(G const& geom, double ratio) {
    bg_box_type env;
    bg::envelope(geom, env);

    auto miny = env.min_corner().y();
    auto height = env.max_corner().y() - miny;

    env.max_corner().set<1>(miny + height*ratio);
    return env;
}

Демо

Добавление визуализации тоже:

Для визуализации мы масштабируем весь ввод в 100x

Live On Coliru

#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/multi_polygon.hpp>
#include <boost/geometry/geometries/linestring.hpp>
#include <boost/geometry/geometries/multi_polygon.hpp>
#include <boost/geometry/core/cs.hpp>
#include <boost/geometry/arithmetic/arithmetic.hpp>
#include <boost/geometry/algorithms/for_each.hpp>
#include <boost/geometry/algorithms/envelope.hpp>
#include <boost/geometry/algorithms/intersection.hpp>
#include <boost/geometry/io/io.hpp>
#include <iostream>
#include <fstream>

namespace bg = boost::geometry;

using bg_point_type = bg::model::d2::point_xy<double>;
using bg_ring_type = bg::model::ring<bg_point_type>;
using bg_polygon_type = bg::model::polygon<bg_point_type>;
using bg_multipolygon_type = bg::model::multi_polygon<bg_polygon_type>;
using bg_box_type = bg::model::box<bg_point_type>;

template <typename G>
bg_box_type crop_box(G const& geom, double ratio) {
    bg_box_type env;
    bg::envelope(geom, env);

    auto const maxy   = env.max_corner().y();
    auto const height = maxy - env.min_corner().y() ;

    env.min_corner().set<1>(maxy - height*ratio);
    return env;
}

template <typename G>
void diags(std::string name, G& geom) {
    std::cout << name << ": " << bg::wkt(geom) << "\n";

    std::string reason;
    if (!bg::is_valid(geom, reason)) {
        std::cout << name << ": " << reason << "\n";

        bg::correct(geom);

        std::cout << bg::wkt(geom) << "\n";
        if (!bg::is_valid(geom, reason)) {
            std::cout << name << " corrected: " << reason << "\n";
        }
    }
}

int main(int argc, char** argv) {
    bool const visualize = argc>1 && argv[1]==std::string("-v");

    //Find the envelope of the region of the object
    bg_polygon_type obj_region;
    bg::read_wkt("POLYGON((0.1 0.1,0.1 0.5,0.5 0.5,0.5 0.1,0.1 0.1))", obj_region);

    if (visualize) {
        bg::for_each_point(obj_region, [](bg_point_type& p) { bg::multiply_value(p, 100.0); });
    }

    diags("Input", obj_region);

    bg_multipolygon_type out;
    bg::intersection(crop_box(obj_region, 0.5), obj_region, out);
    diags("Output", out);

    if (visualize) {
        std::ofstream svg("svg.svg");
        boost::geometry::svg_mapper<bg_point_type> mapper(svg, 600, 600);
        mapper.add(obj_region);
        mapper.add(out);

        mapper.map(obj_region, "fill-opacity:0.5;fill:rgb(204,153,0);stroke:rgb(204,153,0);stroke-width:2");
        mapper.map(out, "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
    }
}

Отпечатки

Input: POLYGON((0.1 0.1,0.1 0.5,0.5 0.5,0.5 0.1,0.1 0.1))
Output: MULTIPOLYGON(((0.1 0.5,0.5 0.5,0.5 0.3,0.1 0.3,0.1 0.5)))

Или когда Включение визуализации :

Input: POLYGON((10 10,10 50,50 50,50 10,10 10))
Output: MULTIPOLYGON(((10 50,50 50,50 30,10 30,10 50)))

И следующий вывод SVG:

enter image description here

...