1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198
| unit DouglasPeuckers ;
{ Implementation de l'algorithme de Douglas-Peucker qui sert à simplifier
une polyligne ou un polygone (élimination des points presque alignés à
une tolérance prêt).
References:
David Douglas & Thomas Peucker, "Algorithms for the reduction of the number of
points required to represent a digitized line or its caricature", The Canadian
Cartographer 10(2), 112-122 (1973)
Code original de Nils Haeck
Utilisation :
PolySimplify(Tol: TFloat; const Orig: array of TPointFloat; var Simple: TPointFloatAr): integer;
Orig[] : Polyligne de départ
Simple[] : Polyligne résultante
Tol : tolérance d'alignement (en distance) pour les points à supprimer
Renvoie le nombre de points de Simple[]
}
interface
uses
Windows ;
type
// Type Float
TFloat = double ;
// Point de coordonnées de type double
TPointFloat = packed record
X: TFloat ;
Y: TFloat ;
end;
TPointFloatAr = array of TPointFloat ;
{ PolySimplify:
La polyligne stockée dans Orig est simplifiée et retournée dans Simple.
La distance maximum des points à supprimer est donnée dans Tol
Entrée: Tol = approximation tolerance
Orig[] = Tableau de points de la polyligne initiale
Sortie: Simple[] = Tableau de points de la polyligne simplifiée
(de longueur initiale à celle de Orig[]
Retour: Nombre de points de Simple[]
}
function PolySimplify(Tol: TFloat; const Orig: array of TPointFloat; var Simple: TPointFloatAr): integer;
implementation
uses Math ;
{ ========================================================================= }
function Vecteur(const A, B: TPointFloat): TPointFloat;
// Resultat Vecteur AB
begin
Result.X := B.X - A.X ;
Result.Y := B.Y - A.Y ;
end;
{ ========================================================================= }
function ProdScal(const A, B: TPointFloat): TFloat;
// Produit scalaire = A . B
begin
Result := A.X * B.X + A.Y * B.Y ;
end;
{ ========================================================================= }
function Norme2(const A: TPointFloat): TFloat;
// Carré de la norme |A|
begin
Result := A.X * A.X + A.Y * A.Y ;
end;
{ ========================================================================= }
function Distance2(const A, B: TPointFloat): TFloat;
// Carré de la distance entre A et B
begin
Result := Norme2(Vecteur(A, B)) ;
end;
{ ========================================================================= }
procedure Simplify(var Tol2: TFloat; const Orig: array of TPointFloat; var Marqueur: array of boolean; j, k: integer);
// Simplifie la polyligne Orig entre j et k en marquant les points a true pour les conserver
// par Marquer[i]
// Tol2 = carré de la tolérance pour les points à supprimer
//
var
i, MaxI: integer ; // Index maximum
MaxD2: TFloat ; // Maxi de la distance au carré du point au segment
CU, CW, B: TFloat ;
DV2: TFloat ; // Distance au carré
P0, P1, PB, U, W: TPointFloat ;
begin
// Polyligne composée de + de 2 points ==> on sort
if k <= j + 1 then
exit ;
// Traitement du segment j k
P0 := Orig[j] ;
P1 := Orig[k] ;
U := Vecteur(P0, P1) ; // Segment j k ou PO P1
CU := ProdScal(U, U) ; // Carré de la longueur du segment P0 P1
MaxD2 := 0 ;
MaxI := 0 ;
// Recherche du point le plus éloigné du segment
for i := j + 1 to k - 1 do
begin
W := Vecteur(P0, Orig[i]) ;
CW := ProdScal(W, U) ;
// Distance de ce point au segment
if CW <= 0 then
begin
// Projection avant le segment
DV2 := Distance2(Orig[i], P0)
end
else
begin
if CW > CU then
begin
// Projection après le segment
DV2 := Distance2(Orig[i], P1) ;
end
else
begin
// Projection sur le segment
if CU=0 then
B := 0
else
B := CW / CU ;
PB.X := P0.X + B * U.X;
PB.Y := P0.Y + B * U.Y ;
DV2 := Distance2(Orig[i], PB);
end ;
end ;
// test par rapport à la distance max
if DV2 > MaxD2 then
begin
// Orig[i] est un sommet à conserver
MaxI := i ;
MaxD2 := DV2 ;
end;
end;
// si le point le plus loin est hors tolérance, on découpe
if MaxD2 > Tol2 then
begin
// découpage de la polyligne
Marqueur[MaxI] := True ; // marquer Orig[maxi] pour la polyligne simplifiée
// simplification récursive des 2 polylignes découpées
Simplify(Tol2, Orig, Marqueur, j, MaxI) ; // polyligne Orig[j] à Orig[maxi]
Simplify(Tol2, Orig, Marqueur, MaxI, k) ; // polyligne Orig[maxi] à Orig[k]
end;
end;
{ ========================================================================= }
function PolySimplify(Tol: TFloat; const Orig: array of TPointFloat; var Simple: TPointFloatAr): integer;
// Simplifie la polyligne Orig[] entre 0 et N en supprimant les points presque alignés à Tol près
// Résultats dans Simple[]
// Tol = tolérance pour la distance des points à supprimer
// Renvoie le nombre de points de Simple[]
var
i : integer;
N: integer; // Taille de la polyligne initiale
Marqueur1: array of boolean; // à True si on doit conserver le point
Tol2: TFloat; // Tolérance pour éliminer les points
begin
Result := 0;
if length(Orig) < 2 then
exit;
Tol2 := sqr(Tol); // Car Simplify utilise le carré de la distance
// Création et initialisation du marqueur booléen des points à conserver
N := length(Orig);
SetLength(Marqueur1, N);
// les points début et fin font partie de la polyligne simplifiée
Marqueur1[0] := True;
Marqueur1[N - 1] := True;
// Exclusion des autres points de la polyligne
for i := 1 to N - 2 do Marqueur1[i] := False;
// Simplification de la polyligne
Simplify(Tol2, Orig, Marqueur1, 0, N - 1);
// remplissage de la liste des points résultats
SetLength(Simple, N);
for i := 0 to N - 1 do
begin
if Marqueur1[i] then
begin
Simple[Result] := Orig[i];
inc(Result);
end;
end;
// redimension de la liste des points résultats
SetLength(Simple, Result);
end;
{ ========================================================================= }
end. |
Partager