1 /* Polygon.cpp
2 *
3 * Copyright (C) 1992-2012,2014-2021 Paul Boersma
4 *
5 * This code is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 2 of the License, or (at
8 * your option) any later version.
9 *
10 * This code is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 * See the GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this work. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19 #include "Polygon.h"
20
21 #include "oo_DESTROY.h"
22 #include "Polygon_def.h"
23 #include "oo_COPY.h"
24 #include "Polygon_def.h"
25 #include "oo_EQUAL.h"
26 #include "Polygon_def.h"
27 #include "oo_CAN_WRITE_AS_ENCODING.h"
28 #include "Polygon_def.h"
29 #include "oo_WRITE_BINARY.h"
30 #include "Polygon_def.h"
31 #include "oo_READ_BINARY.h"
32 #include "Polygon_def.h"
33 #include "oo_DESCRIPTION.h"
34 #include "Polygon_def.h"
35
36 Thing_implement (Polygon, Daata, 1);
37
v_info()38 void structPolygon :: v_info () {
39 our structDaata :: v_info ();
40 MelderInfo_writeLine (U"Number of points: ", our numberOfPoints);
41 MelderInfo_writeLine (U"Perimeter: ", Melder_single (Polygon_perimeter (this)));
42 }
43
v_writeText(MelderFile file)44 void structPolygon :: v_writeText (MelderFile file) {
45 texputi32 (file, our numberOfPoints, U"numberOfPoints");
46 for (integer i = 1; i <= our numberOfPoints; i ++) {
47 texputr64 (file, our x [i], U"x [", Melder_integer (i), U"]");
48 texputr64 (file, our y [i], U"y [", Melder_integer (i), U"]");
49 }
50 }
51
v_readText(MelderReadText text,int)52 void structPolygon :: v_readText (MelderReadText text, int /*formatVersion*/) {
53 our numberOfPoints = texgeti32 (text);
54 if (our numberOfPoints < 1)
55 Melder_throw (U"Cannot read a Polygon with only ", our numberOfPoints, U" points.");
56 our x = raw_VEC (our numberOfPoints);
57 our y = raw_VEC (our numberOfPoints);
58 for (integer i = 1; i <= our numberOfPoints; i ++) {
59 our x [i] = texgetr64 (text);
60 our y [i] = texgetr64 (text);
61 }
62 }
63
Polygon_create(integer numberOfPoints)64 autoPolygon Polygon_create (integer numberOfPoints) {
65 try {
66 autoPolygon me = Thing_new (Polygon);
67 my numberOfPoints = numberOfPoints;
68 my x = zero_VEC (numberOfPoints);
69 my y = zero_VEC (numberOfPoints);
70 return me;
71 } catch (MelderError) {
72 Melder_throw (U"Polygon not created.");
73 }
74 }
75
Polygon_randomize(Polygon me)76 void Polygon_randomize (Polygon me) {
77 for (integer i = 1; i <= my numberOfPoints; i ++) {
78 integer j = NUMrandomInteger (i, my numberOfPoints);
79 double xdum = my x [i];
80 double ydum = my y [i];
81 my x [i] = my x [j];
82 my y [i] = my y [j];
83 my x [j] = xdum;
84 my y [j] = ydum;
85 }
86 }
87
Polygon_perimeter(Polygon me)88 double Polygon_perimeter (Polygon me) {
89 if (my numberOfPoints < 1) return 0.0;
90 double dx = my x [1] - my x [my numberOfPoints], dy = my y [1] - my y [my numberOfPoints];
91 double result = hypot (dx, dy);
92 for (integer i = 1; i <= my numberOfPoints - 1; i ++) {
93 dx = my x [i] - my x [i + 1];
94 dy = my y [i] - my y [i + 1];
95 result += hypot (dx, dy);
96 }
97 return result;
98 }
99
computeDistanceTable(Polygon me,INTMAT const & table)100 static void computeDistanceTable (Polygon me, INTMAT const& table) {
101 for (integer i = 1; i <= my numberOfPoints - 1; i ++)
102 for (integer j = i + 1; j <= my numberOfPoints; j ++) {
103 double dx = my x [i] - my x [j], dy = my y [i] - my y [j];
104 table [i] [j] = table [j] [i] =
105 Melder_ifloor (hypot (dx, dy)); // round to zero
106 }
107 }
108
computeTotalDistance(constINTMATVU const & distance,constINTVECVU const & path)109 static integer computeTotalDistance (constINTMATVU const& distance, constINTVECVU const& path) {
110 integer result = 0;
111 for (integer i = 2; i <= path.size; i ++)
112 result += distance [path [i - 1]] [path [i]];
113 return result;
114 }
115
shuffle(INTVECVU const & path)116 static void shuffle (INTVECVU const& path) {
117 for (integer i = 1; i < path.size; i ++) {
118 integer j = NUMrandomInteger (i, path.size - 1);
119 std::swap (path [i], path [j]);
120 }
121 path [path.size] = path [1];
122 }
123
tryExchange(constINTMATVU const & distance,INTVECVU const & path,integer * totalDistance)124 static bool tryExchange (constINTMATVU const& distance, INTVECVU const& path, integer *totalDistance) {
125 bool result = false;
126 integer b1 = path [1];
127 integer b2nr = 1;
128 while (b2nr < path.size - 2) {
129 integer b2 = path [b2nr + 1];
130 integer distance_b1_b2 = distance [b1] [b2];
131 integer d2nr = b2nr + 2;
132 integer d1 = path [d2nr];
133 bool cont = true;
134 while (d2nr < path.size && cont) {
135 integer d2 = path [d2nr + 1];
136 integer gain = distance_b1_b2 + distance [d1] [d2] - distance [b1] [d1] - distance [b2] [d2];
137 if (gain > 0) {
138 integer below = b2nr + 1, above = d2nr;
139 cont = false;
140 do {
141 integer help = path [below];
142 path [below ++] = path [above];
143 path [above --] = help;
144 } while (below < above);
145 *totalDistance -= gain;
146 }
147 d1 = d2;
148 d2nr ++;
149 }
150 if (cont) {
151 b1 = b2;
152 b2nr ++;
153 } else
154 result = true;
155 }
156 return result;
157 }
158
tryAdoption(constINTMATVU const & distance,INTVECVU const & path,integer * totalDistance)159 static bool tryAdoption (constINTMATVU const& distance, INTVECVU const& path, integer *totalDistance)
160 {
161 autoINTVEC help = zero_INTVEC (path.size);
162 bool result = false;
163
164 /*
165 Compute the maximum distance between two successive cities.
166 */
167 integer city1 = path [1], city2 = path [2];
168 integer maximumDistance = distance [city1] [city2];
169 for (integer i = 2; i < path.size; i ++) {
170 city1 = city2;
171 city2 = path [i + 1];
172 if (distance [city1] [city2] > maximumDistance)
173 maximumDistance = distance [city1] [city2];
174 }
175
176 integer maximumGainLeft = maximumDistance;
177 for (integer i = 1; i < path.size; i ++) {
178 bool cont = true;
179 integer b1, b2, distance_b1_b2, d1nr = 3, cc, e1nrMax = 6;
180 integer numberOfCitiesMinus1 = path.size - 2;
181 for (integer j = 1; j < path.size; j ++)
182 path [j] = path [j + 1];
183 path [path.size] = path [1];
184 b1 = path [1];
185 b2 = path [2];
186 distance_b1_b2 = distance [b1] [b2];
187 cc = path [3];
188 while (d1nr < numberOfCitiesMinus1 && cont) {
189 integer d1 = path [d1nr + 1];
190 integer gain1 = distance_b1_b2 + distance [d1] [cc] - distance [d1] [b2];
191 if (gain1 + maximumGainLeft > 0) {
192 integer e1nr = d1nr + 1;
193 integer dn = path [d1nr + 1];
194 Melder_clipRight (& e1nrMax, numberOfCitiesMinus1);
195 while (e1nr < e1nrMax && cont) {
196 integer e1 = path [e1nr + 1];
197 integer gain = gain1 + distance [dn] [e1] - distance [dn] [b1] - distance [cc] [e1];
198 if (gain > 0) {
199 integer nAdoption = e1nr - d1nr;
200 integer dnnr = e1nr - 1;
201 cont = false;
202 *totalDistance -= gain;
203 for (integer j = 1; j <= dnnr; j ++)
204 help [j] = path [j + 1];
205 for (integer j = 1; j <= nAdoption; j ++)
206 path [j + 1] = help [dnnr - j + 1];
207 for (integer j = 1; j <= d1nr - 1; j ++)
208 path [nAdoption + j + 1] = help [j];
209 }
210 dn = e1;
211 e1nr ++;
212 }
213 }
214 e1nrMax ++;
215 cc = d1;
216 d1nr ++;
217 }
218 result |= ! cont;
219 }
220 return result;
221 }
222
Polygon_salesperson(Polygon me,integer numberOfIterations)223 void Polygon_salesperson (Polygon me, integer numberOfIterations) {
224 try {
225 integer numberOfShortest = 1, totalDistance, shortestDistance = 0;
226
227 integer numberOfCities = my numberOfPoints;
228 if (numberOfCities < 1)
229 Melder_throw (U"No points.");
230 autoINTMAT distance = zero_INTMAT (numberOfCities, numberOfCities);
231 computeDistanceTable (me, distance.get());
232 autoINTVEC path = zero_INTVEC (numberOfCities + 1);
233 for (integer i = 1; i <= numberOfCities; i ++)
234 path [i] = i;
235 path [numberOfCities + 1] = 1; // close path
236 autoINTVEC shortestPath = copy_INTVEC (path.all());
237 for (integer iteration = 1; iteration <= numberOfIterations; iteration ++) {
238 if (iteration > 1)
239 shuffle (path.all());
240 totalDistance = computeTotalDistance (distance.get(), path.all());
241 if (iteration == 1)
242 shortestDistance = totalDistance;
243 do {
244 do {
245 } while (tryExchange (distance.get(), path.all(), & totalDistance));
246 } while (tryAdoption (distance.get(), path.all(), & totalDistance));
247 if (totalDistance < shortestDistance) { // new shortest path
248 numberOfShortest = 1;
249 shortestPath.all() <<= path.all();
250 shortestDistance = totalDistance;
251 }
252 else if (totalDistance == shortestDistance) // shortest path confirmed
253 numberOfShortest ++;
254 }
255 if (numberOfIterations > 1)
256 Melder_casual (U"Polygon_salesperson:"
257 U" found ", numberOfShortest,
258 U" times the same shortest path.");
259
260 /* Change me: I will follow the shortest path found. */
261
262 autoPolygon help = Data_copy (me);
263 for (integer i = 1; i <= numberOfCities; i ++) {
264 my x [i] = help -> x [shortestPath [i]];
265 my y [i] = help -> y [shortestPath [i]];
266 }
267 } catch (MelderError) {
268 Melder_throw (me, U": shortest path not found.");
269 }
270 }
271
setWindow(Polygon me,Graphics graphics,double xmin,double xmax,double ymin,double ymax)272 static void setWindow (Polygon me, Graphics graphics,
273 double xmin, double xmax, double ymin, double ymax)
274 {
275 Melder_assert (me);
276 if (xmax == xmin) { // autoscaling along x axis
277 xmax = xmin = my x [1];
278 for (integer i = 2; i <= my numberOfPoints; i ++) {
279 if (my x [i] < xmin)
280 xmin = my x [i];
281 if (my x [i] > xmax)
282 xmax = my x [i];
283 }
284 if (xmin == xmax) {
285 xmin -= 1.0;
286 xmax += 1.0;
287 }
288 }
289 if (ymax == ymin) { // autoscaling along y axis
290 ymax = ymin = my y [1];
291 for (integer i = 2; i <= my numberOfPoints; i ++) {
292 if (my y [i] < ymin)
293 ymin = my y [i];
294 if (my y [i] > ymax)
295 ymax = my y [i];
296 }
297 if (ymin == ymax) {
298 ymin -= 1.0;
299 ymax += 1.0;
300 }
301 }
302 Graphics_setWindow (graphics, xmin, xmax, ymin, ymax);
303 }
304
Polygon_draw(Polygon me,Graphics g,double xmin,double xmax,double ymin,double ymax)305 void Polygon_draw (Polygon me, Graphics g, double xmin, double xmax, double ymin, double ymax) {
306 Graphics_setInner (g);
307 setWindow (me, g, xmin, xmax, ymin, ymax);
308 Graphics_polyline (g, my numberOfPoints, & my x [1], & my y [1]);
309 Graphics_unsetInner (g);
310 }
311
Polygon_drawClosed(Polygon me,Graphics g,double xmin,double xmax,double ymin,double ymax)312 void Polygon_drawClosed (Polygon me, Graphics g, double xmin, double xmax, double ymin, double ymax) {
313 Graphics_setInner (g);
314 setWindow (me, g, xmin, xmax, ymin, ymax);
315 Graphics_polyline_closed (g, my numberOfPoints, & my x [1], & my y [1]);
316 Graphics_unsetInner (g);
317 }
318
Polygon_paint(Polygon me,Graphics g,MelderColour colour,double xmin,double xmax,double ymin,double ymax)319 void Polygon_paint (Polygon me, Graphics g, MelderColour colour, double xmin, double xmax, double ymin, double ymax) {
320 Graphics_setInner (g);
321 setWindow (me, g, xmin, xmax, ymin, ymax);
322 Graphics_setColour (g, colour);
323 Graphics_fillArea (g, my numberOfPoints, & my x [1], & my y [1]);
324 Graphics_unsetInner (g);
325 }
326
Polygon_drawCircles(Polygon me,Graphics g,double xmin,double xmax,double ymin,double ymax,double diameter_mm)327 void Polygon_drawCircles (Polygon me, Graphics g,
328 double xmin, double xmax, double ymin, double ymax, double diameter_mm)
329 {
330 Graphics_setInner (g);
331 setWindow (me, g, xmin, xmax, ymin, ymax);
332 for (integer i = 1; i <= my numberOfPoints; i ++)
333 Graphics_circle_mm (g, my x [i], my y [i], diameter_mm);
334 Graphics_unsetInner (g);
335 }
336
Polygon_paintCircles(Polygon me,Graphics g,double xmin,double xmax,double ymin,double ymax,double diameter)337 void Polygon_paintCircles (Polygon me, Graphics g,
338 double xmin, double xmax, double ymin, double ymax, double diameter)
339 {
340 Graphics_setInner (g);
341 setWindow (me, g, xmin, xmax, ymin, ymax);
342 for (integer i = 1; i <= my numberOfPoints; i ++)
343 Graphics_fillCircle_mm (g, my x [i], my y [i], diameter);
344 Graphics_unsetInner (g);
345 }
346
Polygons_drawConnection(Polygon me,Polygon thee,Graphics g,double xmin,double xmax,double ymin,double ymax,int hasArrow,double relativeLength)347 void Polygons_drawConnection (Polygon me, Polygon thee, Graphics g,
348 double xmin, double xmax, double ymin, double ymax, int hasArrow, double relativeLength)
349 {
350 double w2 = 0.5 * (1 - relativeLength), w1 = 1 - w2;
351 integer n = my numberOfPoints;
352 if (thy numberOfPoints < n) n = thy numberOfPoints;
353 Graphics_setInner (g);
354 setWindow (me, g, xmin, xmax, ymin, ymax);
355 for (integer i = 1; i <= n; i ++) {
356 double x1 = my x [i], x2 = thy x [i], y1 = my y [i], y2 = thy y [i];
357 double dummy = w1 * x1 + w2 * x2;
358 x2 = w1 * x2 + w2 * x1;
359 x1 = dummy;
360 dummy = w1 * y1 + w2 * y2;
361 y2 = w1 * y2 + w2 * y1;
362 y1 = dummy;
363 if (hasArrow)
364 Graphics_arrow (g, x1, y1, x2, y2);
365 else
366 Graphics_line (g, x1, y1, x2, y2);
367 }
368 Graphics_unsetInner (g);
369 }
370
371 /* End of file Polygon.cpp */
372