boost::polygon はまった

三角形分割するのに、boost::polygon を使ってたんですが、残念ながら実数版は実装されているもののサポート対象外みたいで、ロバスト性が保証されていないようです。 Android をはじめとする GLES では、Tesselation（三角形分割）がサポートされていないので、Trapezoid (台形）化してポリゴンを塗り潰してたんですけども、残念です。 以下、削ったコード。

#include <boost/cstdint.hpp>
#include <boost/polygon/polygon.hpp>
#include <vector>

using boost::int32_t;
using boost::uint8_t;

namespace bpl = boost::polygon;
typedef float point_value_type;
typedef bpl::polygon_with_holes_data<point_value_type> Polygon;
typedef bpl::polygon_data<point_value_type> Ring;
typedef bpl::polygon_traits<Polygon>::point_type Point;
typedef std::vector<Point> PointSet;
typedef std::vector<Polygon> PolygonSet;

void set_points( Polygon& polygon, const point_value_type* pts, int32_t count ) {
  polygon.self_.coords_.resize( count );
  for( int32_t i = 0; i < count; ++i ) {
    polygon.self_.coords_[i].x( *pts++ );
    polygon.self_.coords_[i].y( *pts++ );
  }
}

void construct( PolygonSet& polygons, const point_value_type* pts, int32_t count ) {
  using namespace bpl::operators;
  Polygon polygon;
  set_points( polygon, pts, count );
  polygons |= polygon;  //// SIGSEGV
}


int main() {
  float pts[] = {
   -30162.277344, -21684.941406,
   -30154.300781, -21684.707031,
   -30147.291016, -21768.253906,
   -30148.601562, -21752.630859,
   -30156.675781, -21753.292969,
   -30155.425781, -21768.546875,
  };
  
  PolygonSet polygons;
  construct( polygons, pts, 6 );

  
  return 0;
}

追記：どんな図形か調べてみた。ま、往々にしてありがちな図形やね。平行な直線が折り返しているとダメなのであろうか？普通、平行な直線は交点の計算ができないので、交点を求めたりはしない。これをやろうとすると泥沼に陥る事が多い。

2014/03/06 追記: 整数の同じ構成にしてテストを行うと、ちゃんとパスする。

#ifndef SEGF
  typedef int32_t point_value_type;
#else
  typedef float point_value_type;
#endif

にして

  point_value_type pts[] = {
#ifndef SEGF
   -22277344,  -4941406,
   -14300781,  -4707031,
    -7291016, -88253906,
    -8601562, -72630859,
   -16675781, -73292969,
   -15425781, -88546875,
#else
   -22.277344,  -4.941406,
   -14.300781,  -4.707031,
    -7.291016, -88.253906,
    -8.601562, -72.630859,
   -16.675781, -73.292969,
   -15.425781, -88.546875,
#endif
  };

ここは、難しいようで、 polygon_arbitrary_format.hpp にもデバッグ用のコードがコメントアウトされる形で埋め込まれていた。 両者で実行して、違いを見比べてみたところ セーフバージョン

processEvent_
loop
current Y -72630881
scanline size 2
( -16675781, -73292969), ( -8601822, -72630881); ( -14300781, -4707031), ( -8601822, -72630881); 
found element in scanline 1
finding elements in tree
first iter y is -7.26309e+007
loop2
loop2
counts_from_scanline size 2
aggregating
loop3
loop3

落ちるバージョン

processEvent_
loop
current Y -72.6309
scanline size 2
( -16.6758, -73.293), ( -8.60156, -72.6309); ( -14.3008, -4.70703), ( -7.29102, -88.2539); 
found element in scanline 1
finding elements in tree
first iter y is -72.6309
loop2
counts_from_scanline size 1
aggregating
loop3
loop3

更に、loop2 の条件判定文のコメントを polygon_arbitrary_format.hpp に追加して追試してみたところ、

        if(iter != scanData_.end())
          std::cout << "first iter y is " << iter->first.evalAtX(x_) << "\n";
        while(iter != scanData_.end() &&
              ((iter->first.pt.x() == x_ && iter->first.pt.y() == currentY) ||
               (iter->first.other_pt.x() == x_ && iter->first.other_pt.y() == currentY))) {
                //iter->first.evalAtX(x_) == (high_precision)currentY) {
          std::cout << "loop2\n";
          elementIters.push_back(iter);
          counts_from_scanline.push_back(std::pair<std::pair<std::pair<Point, Point>, int>, active_tail_arbitrary*>
                                         (std::pair<std::pair<Point, Point>, int>(std::pair<Point, Point>(iter->first.pt,
                                                                                                          iter->first.other_pt),
                                                                                  iter->first.count),
                                          iter->second));
          ++iter;
          if( iter != scanData_.end() ) {
            std::cout << "CHK loop2 " << iter->first.pt.x() << " == " << x_ << " && " << iter->first.pt.y() << " == " << currentY << "\n";
            std::cout << "CHK loop2 " << iter->first.other_pt.x() << " == " << x_ << " && " << iter->first.other_pt.y() << " == " << currentY << "\n";
          }
        }

セーフバージョンと落ちるバージョンでは、そもそもscanline の構成が異なっている事がわかった

processEvent_
loop
current Y -72630881
scanline size 2
( -16675781, -73292969), ( -8601822, -72630881); ( -14300781, -4707031), ( -8601822, -72630881); 
found element in scanline 1
finding elements in tree
first iter y is -7.26309e+007
loop2
CHK loop2 -14300781 == -8601822 && -4707031 == -72630881
CHK loop2 -8601822 == -8601822 && -72630881 == -72630881
loop2
counts_from_scanline size 2
aggregating

processEvent_
loop
current Y -72.6309
scanline size 2
( -16.6758, -73.293), ( -8.60156, -72.6309); ( -14.3008, -4.70703), ( -7.29102, -88.2539); 
found element in scanline 1
finding elements in tree
first iter y is -72.6309
loop2
CHK loop2 -14.3008 == -8.60156 && -4.70703 == -72.6309
CHK loop2 -7.29102 == -8.60156 && -88.2539 == -72.6309
counts_from_scanline size 1
aggregating
loop3
loop3

浮動小数点型の演算は、難しいものだ。やはり、ロバスト性を確保するのは並大抵のことではない感じだ。