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. |
Regression linéaire multiple
Répondre avec citation
Partager