[GGL][Polygon][Union] Union returns hole (inner ring) of the expected result as the output polygon (outer ring)
Hi, Union of two polygon_2d returns a hole (inner ring) as result. Where - Input polygons - Polygon 1: A big rectangle with a hexagonal hole inside Polygon 2: A triangle that intersects the hole from inside leaving portion of it in the hole Output polygon - The hexagonal hole without the triangle-portion. Expected output polygon - The big rectangle with current output polygon as hole (inner)as described in this picture. The Union of a rectangle with a hexagonal hole and a triangle intersecting the first polygon ***************************** ***************************** ***************************** ******** ********* ******* **X***** ****** XXX***** ***** XXX****** **** XX ***** ***** ****** ****** ******* ******* ******** ******** ********* ***************************** ***************************** ***************************** returns this polygon, instead of the correct result. ++++++++++++ ++++++++++++++ +++++++++++++++ +++++++++++++++ +++++++++++++++ +++ ++++++++++++++++++ ++++++++++++++++ ++++++++++++++ ++++++++++++ You can use the code bellow to generate the case and some PPM images. You may use a WSQ Viewer to see the polygons. Please let me know any workaround or anything thats I am missing. Thank and regards, Helal // the code void uniontest() // Your main should call it { double L = 200.0; const double outloop[][2] = { {L, L}, {L, -L}, {-L, -L}, {-L, L}, {L, L} }; double R = 100.0; double H = R * sqrt(3.0) / 2; const double inloop[][2] = { {R, 0}, {R/2, -H}, {-R/2, -H}, {-R, 0}, {-R/2, H}, {R/2, H}, {R, 0} }; polygon_2d BIGpoly; //A rectangle with a hexagonal hole in it assign(BIGpoly, outloop); //the rectangle BIGpoly.inners().resize(1); linear_ring<point_2d>& inner = BIGpoly.inners().back(); assign(inner, inloop); //the inner hexagon correct(BIGpoly); //correct the order iff not double a = 20.0; double h = a * sqrt(3.0) / 2; const double t1[][2] = { { 30+0, 10+0}, { 30+a/2, 10+h}, { 30+a, 10+0} }; const double t2[][2] = { { 0, 0}, { a/2, h}, { a, 0} }; const double t3[][2] = { {-10+0, -10+0}, {-10+a/2, -10+h}, {-10+a, -10+0} }; const double t4[][2] = { { 0, -15+0}, { a/2, -15+h}, { a, -15+0} }; const double t5[][2] = { {-50, -40}, {-110, 0}, {-100, 30}, {-50, -40} }; polygon_2d T1, T2, T3, T4, T5, T6; //the triangles assign(T1, t1); correct(T1); assign(T2, t2); correct(T2); assign(T3, t3); correct(T3); assign(T4, t4); correct(T4); assign(T5, t5); correct(T5); std::vector<polygon_2d > P_list, P_union, P_tmp; //the outer shape P_union.push_back(BIGpoly); DumbPPMPrint("c:\\BIGpoly.ppm", P_union.begin(), P_union.size()); //Case 1: floating triangle in hole P_union.clear(); ggl::union_inserter<ggl::polygon_2d >(BIGpoly, T1, std::back_inserter(P_union)); DumbPPMPrint("c:\\case1.ppm", P_union.begin(), P_union.size()); //Case 2: three solid triangle creating a small hole P_union.clear(); P_tmp.clear(); ggl::union_inserter<ggl::polygon_2d >(T2, T3, std::back_inserter(P_tmp)); DumbPPMPrint("c:\\case2_i.ppm", P_tmp.begin(), P_tmp.size()); ggl::union_inserter<ggl::polygon_2d >( (*(P_tmp.begin())), T4, std::back_inserter(P_union)); DumbPPMPrint("c:\\case2.ppm", P_union.begin(), P_union.size()); //Case 3: floating holed-polygon in hole //P_union.clear(); //contains result of case 2 P_tmp.clear(); ggl::union_inserter<ggl::polygon_2d >(BIGpoly, (*(P_union.begin())), std::back_inserter(P_tmp)); DumbPPMPrint("c:\\case3.ppm", P_tmp.begin(), P_tmp.size()); //Case 4: a triangle intersecting the hole from inside P_union.clear(); P_union.push_back(BIGpoly); P_union.push_back(T5); DumbPPMPrint("c:\\case4_pre.ppm", P_union.begin(), P_union.size()); P_union.clear(); ggl::union_inserter<ggl::polygon_2d >( BIGpoly, T5, std::back_inserter(P_union)); DumbPPMPrint("c:\\case4.ppm", P_union.begin(), P_union.size()); P_union.clear(); ggl::union_inserter<ggl::polygon_2d >( T5, BIGpoly, std::back_inserter(P_union)); DumbPPMPrint("c:\\case4_rev_order.ppm", P_union.begin(), P_union.size()); } void DumbPPMPrint(char *name, std::vector<polygon_2d >::iterator it, unsigned int n) { //WSQ Viewer can be of help //http://www.cognaxon.com/index.php?page=wsqview //the holes are black. //the more number of polygons are present at a pixel the brighter the pixel is. unsigned int imgwidth = 512; //pixel width of image unsigned int imgheight = 512; //pixel height of image float width = 400.0; //width of image: 10 units (e.g. meters) float height = 400.0; //height of image: 10 units (e.g. meters) //float width = 6000.0; //width of image: 10 units (e.g. meters) //float height = 6000.0; //height of image: 10 units (e.g. meters) //header FILE *signal = fopen(name, "wt"); fprintf(signal, "P3\n"); fprintf(signal, "# %s\n", name); fprintf(signal, "%d %d\n", imgwidth, imgheight); fprintf(signal, "255\n"); unsigned int x, y, i; float xpos, ypos; // PPM images carry red, green, blue color channels. int val0, val1, val2; std::vector<polygon_2d >::iterator poly; for (y = 0; y < imgheight; y++) { ypos = height/2 - height * y / imgheight; for (x = 0; x < imgwidth; x++) { xpos = width/2 - width * x / imgwidth; val0 = val1 = val2 = 0; //black ggl::point_xy<float> point(xpos, ypos); poly = it; for (i = 0; i < n; i++, poly++) { if( ggl::within( point, (*poly) ) ) { val0 += 255/n; val1 += 255/n; val2 += 255/n; } } //fi_poly // clamp to allowed range if (val0 < 0) val0 = 0; else if (val0 > 255) val0 = 255; if (val1 < 0) val1 = 0; else if (val1 > 255) val1 = 255; if (val2 < 0) val2 = 0; else if (val2 > 255) val2 = 255; fprintf(signal, "%3d %3d %3d ", val0, val1, val2); if (x%16 == 15) fprintf(signal, "\n"); } //fi_x } //fi_y fclose(signal); } //@End_of_Code
Hi Helal , This was not a known problem, I will investigate this, and answer you tomorrow or Friday (I cannot start immediately...).
The Union of a rectangle with a hexagonal hole and a triangle intersecting the first polygon
***************************** ***************************** ***************************** ******** ********* ******* **X***** ****** XXX***** ***** XXX****** **** XX ***** ***** ****** ****** ******* ******* ******** ******** ********* ***************************** ***************************** *****************************
returns this polygon, instead of the correct result.
++++++++++++ ++++++++++++++ +++++++++++++++ +++++++++++++++ +++++++++++++++ +++ ++++++++++++++++++ ++++++++++++++++ ++++++++++++++ ++++++++++++
You can use the code bellow to generate the case and some PPM images. You may use a WSQ Viewer to see the polygons.
Please let me know any workaround or anything thats I am missing.
The problem is completely clear, thanks for the useful pictures and source code. Regards, Barend
Hi,
Union of two polygon_2d returns a hole (inner ring) as result. Where -
Input polygons - Polygon 1: A big rectangle with a hexagonal hole inside Polygon 2: A triangle that intersects the hole from inside leaving portion of it in the hole
Output polygon - The hexagonal hole without the triangle-portion.
Expected output polygon - The big rectangle with current output polygon as hole (inner)as described in this picture.
The exterior ring was omitted from the output. This happened earlier with Luke's test in November, it was solved provisionally then, but apparently this is the same problem. In the meantime the problem has been solved and I hope to deploy it this weekend to Boost sandbox, together with the multi-polygon union (w.r.t. Kef's message). Regards, Barend
Hi Kef, Helal, others, Kef Schecter wrote:
When I try to use equals() on multi_polygons, the compiler blows up, giving a template error.
Mohammad Morshed Alam Helal wrote:
Expected output polygon - The big rectangle with current output polygon as hole (inner)as described in this picture.
These problems are solved now and more is updated. Boost.Geometry (headers & examples) is now publicly available in the boost.sandbox, https://svn.boost.org/svn/boost/sandbox/geometry Documentation can (still) be read here: http://geometrylibrary.geodan.nl/index.html , will be revised later (to Quickbook). Regards, Barend
Hello Barend, Union of multi_polygon has some bugs. I must admit I could not pinpoint the problem. Please let me know if I am missing something. These are the test cases I tried and most of them failed: Case 5: both list has equal number of polygon (PASS) Case 6: one of the list has a single polygon (FAIL) Case 7: both list has one polygon (intersecting) (PASS) Case 8: Different number of polygons in the lists (FAIL) Case 9: 2nd list has a single polygon (FAIL) Case 10: both list has one polygon (non-intersecting) (PASS) Case 11: one list has same polygon twice (FAIL) Case 12A: using same list twice. Its a mess. (FAIL) Case 12B: keeping one of the list empty. Returns empty. (FAIL) Images from Case 9: Before: Multi_polygon 1: **************************** **************************** **************************** ******** ******** ******* X ******* ****** XXX ****** ***** ***** **** Y **** ***** YYYS ***** ****** S SS ****** ******* XXXX ******* ******** ******** **************************** **************************** **************************** Multi_polygon 1: (only the X marked area) ............................ ............................ ............................ ........ ........ ....... ..X.... ...... XXX.... ..... XXX..... .... XX .... ..... ..... ...... ...... ....... ....... ........ ........ ............................ ............................ ............................ After: Multi_polygon: **************************** **************************** **************************** ******** ******** ******* ******* ****** ******* ***** ******** **** ** **** ***** ***** ****** ****** ******* ******* ******** ******** **************************** **************************** **************************** As you see, the floating polygons are missing. You can use the code below to regenerate the cases. I have the PPM images. Please let me know if you need them. Thank you for the fixes for the previous problem with union of polygons. It worked. Best regards, Helal //@_Begin_of_Code void uniontest() { bool case_1 = false, case_2 = false, case_2_3 = false, case_4 = false, case_5 = true, case_6 = true, case_7 = true, case_8 = true, case_9 = true, case_10 = true, case_11 = true, case_12 = true, isDump = false; double L = 200.0; const double outloop[][2] = { {L, L}, {L, -L}, {-L, -L}, {-L, L}, {L, L} }; double R = 100.0; double H = R * sqrt(3.0) / 2; const double inloop[][2] = { {R, 0}, {R/2, -H}, {-R/2, -H}, {-R, 0}, {-R/2, H}, {R/2, H}, {R, 0} }; polygon_2d BIGpoly; //A rectangle with a hexagonal hole in it assign(BIGpoly, outloop); //the rectangle BIGpoly.inners().resize(1); linear_ring<point_2d>& inner = BIGpoly.inners().back(); assign(inner, inloop); //the inner hexagon correct(BIGpoly); //correct the order iff not double a = 20.0; double h = a * sqrt(3.0) / 2; const double t1[][2] = { { 30+0, 10+0}, { 30+a/2, 10+h}, { 30+a, 10+0} }; const double t2[][2] = { { 0, 0}, { a/2, h}, { a, 0} }; const double t3[][2] = { {-10+0, -10+0}, {-10+a/2, -10+h}, {-10+a, -10+0} }; const double t4[][2] = { { 0, -15+0}, { a/2, -15+h}, { a, -15+0} }; const double t5[][2] = { {-50, -40}, {-110, 0}, {-100, 30}, {-50, -40} }; polygon_2d T1, T2, T3, T4, T5, T6; //the triangles assign(T1, t1); correct(T1); assign(T2, t2); correct(T2); assign(T3, t3); correct(T3); assign(T4, t4); correct(T4); assign(T5, t5); correct(T5); std::vector<polygon_2d > P_list, P_union, P_tmp; //the outer shape P_union.push_back(BIGpoly); if (isDump) DumbPPMPrint("c:\\BIGpoly.ppm", P_union.begin(), P_union.size()); //Case 1: floating triangle in hole if (case_1) { P_union.clear(); union_inserter<polygon_2d >(BIGpoly, T1, std::back_inserter(P_union)); if (isDump) DumbPPMPrint("c:\\case1.ppm", P_union.begin(), P_union.size()); } //Case 2: three solid triangle creating a small hole if (case_2) { P_union.clear(); P_tmp.clear(); union_inserter<polygon_2d >(T2, T3, std::back_inserter(P_tmp)); if (isDump) DumbPPMPrint("c:\\case2_i.ppm", P_tmp.begin(), P_tmp.size()); union_inserter<polygon_2d >( (*(P_tmp.begin())), T4, std::back_inserter(P_union)); if (isDump) DumbPPMPrint("c:\\case2.ppm", P_union.begin(), P_union.size()); } //Case 3: floating holed-polygon in hole if (case_2_3) { //P_union.clear(); //contains result of case 2 P_tmp.clear(); union_inserter<polygon_2d >(BIGpoly, (*(P_union.begin())), std::back_inserter(P_tmp)); if (isDump) DumbPPMPrint("c:\\case3.ppm", P_tmp.begin(), P_tmp.size()); } //Case 4: a triangle intersecting the hole from inside if (case_4) { P_union.clear(); P_union.push_back(BIGpoly); P_union.push_back(T5); if (isDump) DumbPPMPrint("c:\\case4_pre.ppm", P_union.begin(), P_union.size()); P_union.clear(); union_inserter<polygon_2d >( BIGpoly, T5, std::back_inserter(P_union)); if (isDump) DumbPPMPrint("c:\\case4.ppm", P_union.begin(), P_union.size()); P_union.clear(); union_inserter<polygon_2d >( T5, BIGpoly, std::back_inserter(P_union)); if (isDump) DumbPPMPrint("c:\\case4_rev_order.ppm", P_union.begin(), P_union.size()); } //Case 5: using multi_polygon (PASS) // both list has equal number of polygon if (case_5) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps1.push_back(T1); ps1.push_back(T2); ps2.push_back(T3); ps2.push_back(T4); ps2.push_back(T5); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case5.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case5_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case5_B.ppm", ps2.begin(), ps2.size()); } //Case 6: using multi_polygon (FAIL) // one of the list has a single polygon if (case_6) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps2.push_back(T1); ps2.push_back(T2); ps2.push_back(T3); ps2.push_back(T4); ps2.push_back(T5); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case6.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case6_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case6_B.ppm", ps2.begin(), ps2.size()); } //Case 7: using multi_polygon (PASS) // both list has one polygon (intersecting) if (case_7) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps2.push_back(T5); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case7.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case7_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case7_B.ppm", ps2.begin(), ps2.size()); } //Case 8: using multi_polygon (FAIL) // Different number of polygons in the lists if (case_8) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps1.push_back(T1); ps2.push_back(T2); ps2.push_back(T3); ps2.push_back(T4); ps2.push_back(T5); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case8.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case8_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case8_B.ppm", ps2.begin(), ps2.size()); } //Case 9: using multi_polygon (FAIL) // 2nd list has a single polygon if (case_9) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps1.push_back(T1); ps1.push_back(T2); ps1.push_back(T3); ps1.push_back(T4); ps2.push_back(T5); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case9.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case9_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case9_B.ppm", ps2.begin(), ps2.size()); } //Case 10: using multi_polygon (PASS) // both list has one polygon (non-intersecting) if (case_10) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps2.push_back(T4); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case10.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case10_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case10_B.ppm", ps2.begin(), ps2.size()); } //Case 11: using multi_polygon (FAIL) // one list has same polygon twice if (case_11) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps2.push_back(T4); ps2.push_back(T4); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); if (isDump) DumbPPMPrint("c:\\case11.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case11_A.ppm", ps1.begin(), ps1.size()); if (isDump) DumbPPMPrint("c:\\case11_B.ppm", ps2.begin(), ps2.size()); } //Case 12: using multi_polygon (FAIL) //A: using same list twice. Its a mess //B: keeping one of the list empty. Returns empty. if (case_12) { typedef multi_polygon<polygon_2d> polygon_set; polygon_set ps1, ps2, ps3; ps1.push_back(BIGpoly); ps1.push_back(T1); ps1.push_back(T2); ps1.push_back(T3); ps1.push_back(T4); ps1.push_back(T5); union_inserter<polygon_2d>(ps1, ps1, std::back_inserter(ps3)); //A //union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); //B if (isDump) DumbPPMPrint("c:\\case12.ppm", ps3.begin(), ps3.size()); if (isDump) DumbPPMPrint("c:\\case12_A.ppm", ps1.begin(), ps1.size()); //if (isDump) */DumbPPMPrint("c:\\case12_B.ppm", ps2.begin(), ps2.size()); } } void DumbPPMPrint(char *name, std::vector<polygon_2d >::iterator it, unsigned int n) { //WSQ Viewer can be of help //http://www.cognaxon.com/index.php?page=wsqview //the holes are black. //the more number of polygons are present at a pixel the brighter the pixel is. unsigned int imgwidth = 256; //pixel width of image unsigned int imgheight = 256; //pixel height of image float width = 400.0; //width of image: 10 units (e.g. meters) float height = 400.0; //height of image: 10 units (e.g. meters) //header FILE *signal = fopen(name, "wt"); fprintf(signal, "P3\n"); fprintf(signal, "# %s\n", name); fprintf(signal, "%d %d\n", imgwidth, imgheight); fprintf(signal, "255\n"); unsigned int x, y, i; float xpos, ypos; // PPM images carry red, green, blue color channels. int val0, val1, val2, weight; std::vector<polygon_2d >::iterator poly; weight = 255/n; if (weight < 101) weight = 101; for (y = 0; y < imgheight; y++) { ypos = height/2 - height * y / imgheight; for (x = 0; x < imgwidth; x++) { xpos = width/2 - width * x / imgwidth; val0 = val1 = val2 = 0; //black point_xy<float> point(xpos, ypos); poly = it; for (i = 0; i < n; i++, poly++) { if( within( point, (*poly) ) ) { val0 += weight; val1 += weight; val2 += weight; } } //fi_poly // clamp to allowed range if (val0 < 0) val0 = 0; else if (val0 > 255) val0 = 255; if (val1 < 0) val1 = 0; else if (val1 > 255) val1 = 255; if (val2 < 0) val2 = 0; else if (val2 > 255) val2 = 255; fprintf(signal, "%3d %3d %3d ", val0, val1, val2); if (x%16 == 15) fprintf(signal, "\n"); } //fi_x } //fi_y fclose(signal); } //@_End_of_Code
Hi Helal,
Union of multi_polygon has some bugs. I must admit I could not pinpoint the problem. Please let me know if I am missing something. These are the test cases I tried and most of them failed:
Thanks for your report, I will have a look soon. Regards, Barend
Hi Helal,
Union of multi_polygon has some bugs. I must admit I could not pinpoint the problem. Please let me know if I am missing something.
Yes, there was indeed a problem. Thanks again for your report. However, I could not reproduce some of your tests marked as failed.
These are the test cases I tried and most of them failed: Case 6: one of the list has a single polygon (FAIL) Case 8: Different number of polygons in the lists (FAIL) Case 9: 2nd list has a single polygon (FAIL) They indeed fail, all for the same reason
Case 11: one list has same polygon twice (FAIL) IMO it passes (see below)
Case 12A: using same list twice. Its a mess. (FAIL) IMO passes perfectly. Case 12B: keeping one of the list empty. Returns empty. (FAIL) IMO it passes and it does not return empty. It is in comments in your code, so I cannot see what was going wrong when it was not commented.
We have to define what is failing and what is passing. The overlay algorithms using multi-polygons currently do not remove internal boundaries. So, in case 12B, you have one multi-polygon of several overlapping polygons, and one empty. You get back the input (in my test). This is my adapted version of your code: ps3.clear(); union_inserter<polygon_2d>(ps1, ps2, std::back_inserter(ps3)); //B dump_wkt("case12b.wkt", ps1, ps2, ps3); (will explain dump_wkt below). (Note that I made a change here last week, it might indeed be that it failed earlier, didn't check). Case 11 passes the same polygon twice, and it is not overlapping with the other polygon, so you get it back like you entered it, combining the two geometries and listing that same polygon twice. I had a look at case 6,8,9 because they should not fail. It is solved and I will update boost.sandbox soon, I will notify the list if that is done. I will come back to the union of two multi-polygons which do have self intersections in another email, it will either be implemented or there will be a possibility to solve it differently. I added a "dump_wkt" which produces well-known text (WKT). Those are strings which I can easily use in my unit tests. In case you find more errors, maybe you can include the WKT's and I will probably not need the source code. #include <boost/geometry/extensions/gis/io/wkt/wkt.hpp> template <typename G1, typename G2, typename G3> inline void dump_wkt(std::string const& s, G1 const& a, G2 const& b, G3 const& c) { std::ofstream out(s.c_str()); out << s << std::endl; out << area(a) << " " << area(b) << " -> " << area(c); out << std::endl; out << wkt(a) << std::endl; out << wkt(b) << std::endl; out << wkt(c) << std::endl; } Regards, Barend
Hi Helal,
Union of multi_polygon has some bugs. I must admit I could not pinpoint the problem. Please let me know if I am missing something. These are the test cases I tried and most of them failed:
Case 6: one of the list has a single polygon (FAIL) Case 8: Different number of polygons in the lists (FAIL) Case 9: 2nd list has a single polygon (FAIL) These are now passing, I just committed the adapted sources. There was a failure in the assemble process which is again revised (and a bit simplified) now. Thanks again for the bug report.
Case 11: one list has same polygon twice (FAIL) Case 12A: using same list twice. Its a mess. (FAIL) Case 12B: keeping one of the list empty. Returns empty. (FAIL)
As said, these input-multi-polygons do have self-intersections which are not solved (not even detected) during the overlay process. According to OGC (Open Geospatial Consortium), these kind of input geometries are not valid. So spatial databases like PostGreSQL or SQL Server 2008 reject this input. There is a solution of course, we are implementing a "dissolve" algorithm, which removes self intersections (internal borders), it can be used for a polygon (e.g. to convert a keyholed interior ring to a normal interior ring) or for a multi-polygon (to solve things like this). This one is not yet ready, but if yes and you call it before union, it will give you the expected results. It probably avoids calls to union in many cases, dissolve does the process of unioning a multi-polygon already. An option is to call dissolve internally, from union. This however would halve the speed, in many cases the input is known to be valid, so this step is not necessary. Besides that, if you enter an input geometry A (multi-polygon with self-overlapping triangles) and B (empty), you will get back A immediatly for a union, you probably would be surprised that A is touched by the process. However, yet another option is to give the union a boolean flag to denote to dissolve first. This is not yet decided. Note that the dissolve is not yet in Boost.Sandbox because it is not yet ready, it is expected to be there one of the coming weeks. Regards, Barend
participants (2)
-
Barend Gehrels -
Mohammad Morshed Alam Helal