3

scipy (0.9.0) と matplotlib (1.0.1) の Delaunay 三角形分割ルーチンを比較すると、説明のつかない動作に気付きました。私のポイントは、に保存されている UTM 座標numpy.array([[easting, northing], [easting, northing], [easting, northing]])です。Scipy のエッジには私のポイントの一部がありませんが、matplotlib のエッジはすべてそこにあります。修正はありますか、それとも何か間違っていますか?

import scipy
import numpy
from scipy.spatial import Delaunay
import matplotlib.delaunay


def delaunay_edges(points):
    d = scipy.spatial.Delaunay(points)
    s = d.vertices
    return numpy.vstack((s[:,:2], s[:,1:], s[:,::-2]))


def delaunay_edges_matplotlib(points):
        cens, edges, tri, neig = matplotlib.delaunay.delaunay(points[:,0], points[:,1])
        return edges


points = numpy.array([[500000.25, 6220000.25],[500000.5, 6220000.5],[500001.0, 6220001.0],[500002.0, 6220003.0],[500003.0, 6220005.0]])

edges1 = delaunay_edges(points)
edges2 = delaunay_edges_matplotlib(points)

numpy.unique(edges1).shape # Some points missing, presumably nearby ones
numpy.unique(edges2).shape # Includes all points
4

1 に答える 1

6

のこの動作はscipy.spatial.Delaunay、浮動小数点演算の不正確さに関連している可能性があります。

ご存じのとおり、scipy.spatial.DelaunayCqhullライブラリを使用して Delaunay 三角形分割を計算します。は、アルゴリズムQhullの実装であり、この論文 (1)で著者によって詳細に説明されています。コンピュータで使用される浮動小数点演算が IEEE 754 標準を使用して実行されることもご存じかもしれません (たとえば、ウィキペディアで読むことができます)。標準によると、各有限数は 3 つの整数で最も簡単に記述されます。Quickhullscq= 指数。これらの整数を表すために使用されるビット数は、データ型によって異なります。したがって、数値軸上の浮動小数点数の分布の密度が一定ではないことは明らかです。数値が大きいほど、分散が緩くなります。Google 電卓でも見ることができます - 3333333333333333 を 3333333333333334 から引いて 0 を得ることができます。これは、3333333333333333 と 3333333333333334 の両方が同じ浮動小数点数に丸められるために発生します。

ここで、丸め誤差について知っているので、Copying with impresicionというタイトルの論文 (1) の第 4 章を読みたいと思うかもしれません。この章では、丸め誤差を処理す​​るアルゴリズムについて説明します。

Quickhull partitions a point and determines its horizon facets by computing 
whether the point is above or below a hyperplane. We have assumed that 
computations return consistent results ... With floating-point arithmetic, we 
cannot prevent errors from occurring, but we can repair the damage after 
processing a point. We use brute force: if adjacent facets are nonconvex, one of 
the facets is merged into a neighbor. Quickhull merges the facet that minimizes 
the maximum distance of a vertex to the neighbor.

これが起こる可能性があります -Quickhullは近くにある 2 つのポイントを区別できないため、2 つのファセットをマージして、そのうちの 1 つを正常に削除します。これらのエラーを解消するには、たとえば、座標の原点を移動してみてください。

co = points[0]
points = points - co

edges1 = delaunay_edges(points)
edges2 = delaunay_edges_matplotlib(points)

print numpy.unique(edges1) 
>>> [0 1 2 3 4]
print numpy.unique(edges2)
>>> [0 1 2 3 4]

これにより、浮動小数点数が適度に密集している領域に計算が移動します。

于 2011-11-21T10:35:39.437 に答える