1 /**********************************************************
2 * Version $Id$
3 *********************************************************/
4
5 ///////////////////////////////////////////////////////////
6 // //
7 // SAGA //
8 // //
9 // System for Automated Geoscientific Analyses //
10 // //
11 // Tool Library //
12 // shapes_polygons //
13 // //
14 //-------------------------------------------------------//
15 // //
16 // Polygon_Intersection.cpp //
17 // //
18 // Copyright (C) 2003 by //
19 // Olaf Conrad //
20 // //
21 //-------------------------------------------------------//
22 // //
23 // This file is part of 'SAGA - System for Automated //
24 // Geoscientific Analyses'. SAGA is free software; you //
25 // can redistribute it and/or modify it under the terms //
26 // of the GNU General Public License as published by the //
27 // Free Software Foundation, either version 2 of the //
28 // License, or (at your option) any later version. //
29 // //
30 // SAGA is distributed in the hope that it will be //
31 // useful, but WITHOUT ANY WARRANTY; without even the //
32 // implied warranty of MERCHANTABILITY or FITNESS FOR A //
33 // PARTICULAR PURPOSE. See the GNU General Public //
34 // License for more details. //
35 // //
36 // You should have received a copy of the GNU General //
37 // Public License along with this program; if not, see //
38 // <http://www.gnu.org/licenses/>. //
39 // //
40 //-------------------------------------------------------//
41 // //
42 // e-mail: oconrad@saga-gis.org //
43 // //
44 // contact: Olaf Conrad //
45 // Institute of Geography //
46 // University of Goettingen //
47 // Goldschmidtstr. 5 //
48 // 37077 Goettingen //
49 // Germany //
50 // //
51 ///////////////////////////////////////////////////////////
52
53 //---------------------------------------------------------
54
55
56 ///////////////////////////////////////////////////////////
57 // //
58 // //
59 // //
60 ///////////////////////////////////////////////////////////
61
62 //---------------------------------------------------------
63 #include "Polygon_Intersection.h"
64
65
66 ///////////////////////////////////////////////////////////
67 // //
68 // //
69 // //
70 ///////////////////////////////////////////////////////////
71
72 //---------------------------------------------------------
CPolygon_Overlay(const CSG_String & Name)73 CPolygon_Overlay::CPolygon_Overlay(const CSG_String &Name)
74 {
75 //-----------------------------------------------------
76 Set_Name (Name);
77
78 Set_Author ("O.Conrad (c) 2003-17");
79
80 Set_Description (_TW(
81 "Uses the free and open source software library <b>Clipper</b> created by Angus Johnson.\n"
82 "<a target=\"_blank\" href=\"http://www.angusj.com/delphi/clipper.php\">Clipper Homepage</a>\n"
83 "<a target=\"_blank\" href=\"http://sourceforge.net/projects/polyclipping/\">Clipper at SourceForge</a>\n"
84 ));
85
86 //-----------------------------------------------------
87 Parameters.Add_Shapes("",
88 "A" , _TL("Layer A"),
89 _TL(""),
90 PARAMETER_INPUT, SHAPE_TYPE_Polygon
91 );
92
93 Parameters.Add_Shapes("",
94 "B" , _TL("Layer B"),
95 _TL(""),
96 PARAMETER_INPUT, SHAPE_TYPE_Polygon
97 );
98
99 Parameters.Add_Shapes("",
100 "RESULT", Name,
101 _TL(""),
102 PARAMETER_OUTPUT, SHAPE_TYPE_Polygon
103 );
104
105 Parameters.Add_Bool("",
106 "SPLIT" , _TL("Split Parts"),
107 _TL("Set true if you want multipart polygons to become separate polygons."),
108 true
109 );
110 }
111
112
113 ///////////////////////////////////////////////////////////
114 // //
115 ///////////////////////////////////////////////////////////
116
117 //---------------------------------------------------------
Add_Description(const CSG_String & Text)118 bool CPolygon_Overlay::Add_Description(const CSG_String &Text)
119 {
120 Set_Description(Text + "\n" + Get_Description());
121
122 return( true );
123 }
124
125 //---------------------------------------------------------
Initialize(CSG_Shapes ** ppA,CSG_Shapes ** ppB,bool bBothAttributes)126 bool CPolygon_Overlay::Initialize(CSG_Shapes **ppA, CSG_Shapes **ppB, bool bBothAttributes)
127 {
128 *ppA = Parameters("A")->asShapes(); if( (*ppA)->Get_Type() != SHAPE_TYPE_Polygon || !(*ppA)->is_Valid() ) return( false );
129 *ppB = Parameters("B")->asShapes(); if( (*ppB)->Get_Type() != SHAPE_TYPE_Polygon || !(*ppB)->is_Valid() ) return( false );
130
131 m_bSplit = Parameters("SPLIT")->asBool();
132
133 m_pA = NULL;
134 m_pB = NULL;
135
136 m_pAB = Parameters("RESULT")->asShapes();
137
138 if( m_pAB == *ppA || m_pAB == *ppB )
139 {
140 Error_Set(_TL("Output layer must not be one of the input layers!"));
141
142 return( false );
143 }
144
145 m_pAB->Create(SHAPE_TYPE_Polygon, SG_T(""), *ppA);
146 m_pAB->Fmt_Name("%s [%s]-[%s]", Get_Name().c_str(), (*ppA)->Get_Name(), (*ppB)->Get_Name());
147
148 if( bBothAttributes )
149 {
150 for(int i=0; i<(*ppB)->Get_Field_Count(); i++)
151 {
152 m_pAB->Add_Field((*ppB)->Get_Field_Name(i), (*ppB)->Get_Field_Type(i));
153 }
154 }
155
156 return( true );
157 }
158
159
160 ///////////////////////////////////////////////////////////
161 // //
162 ///////////////////////////////////////////////////////////
163
164 //---------------------------------------------------------
Get_Intersection(CSG_Shapes * pA,CSG_Shapes * pB)165 bool CPolygon_Overlay::Get_Intersection(CSG_Shapes *pA, CSG_Shapes *pB)
166 {
167 m_bInvert = false;
168
169 m_pA = pA;
170 m_pB = pB;
171
172 CSG_Shapes Result(SHAPE_TYPE_Polygon);
173
174 CSG_Shape_Polygon *pResult = (CSG_Shape_Polygon *)Result.Add_Shape();
175
176 //-----------------------------------------------------
177 for(int id_A=0; id_A<m_pA->Get_Count() && Set_Progress(id_A, m_pA->Get_Count()); id_A++)
178 {
179 for(int id_B=0; id_B<m_pB->Get_Count(); id_B++)
180 {
181 if( SG_Polygon_Intersection(m_pA->Get_Shape(id_A), m_pB->Get_Shape(id_B), pResult) )
182 {
183 _Add_Polygon(pResult, id_A, id_B);
184 }
185 }
186 }
187
188 return( true );
189 }
190
191
192 ///////////////////////////////////////////////////////////
193 // //
194 ///////////////////////////////////////////////////////////
195
196 //---------------------------------------------------------
Get_Difference(CSG_Shapes * pA,CSG_Shapes * pB,bool bInvert)197 bool CPolygon_Overlay::Get_Difference(CSG_Shapes *pA, CSG_Shapes *pB, bool bInvert)
198 {
199 m_bInvert = bInvert;
200
201 m_pA = pA;
202 m_pB = pB;
203
204 CSG_Shapes Result(SHAPE_TYPE_Polygon);
205
206 CSG_Shape_Polygon *pResult = (CSG_Shape_Polygon *)Result.Add_Shape();
207
208 //-----------------------------------------------------
209 for(int id_A=0; id_A<m_pA->Get_Count() && Set_Progress(id_A, m_pA->Get_Count()); id_A++)
210 {
211 pResult->Assign(m_pA->Get_Shape(id_A), false);
212
213 for(int id_B=0; id_B<m_pB->Get_Count() && pResult->is_Valid(); id_B++)
214 {
215 switch( pResult->Intersects(m_pB->Get_Shape(id_B)) )
216 {
217 case INTERSECTION_None:
218 break;
219
220 case INTERSECTION_Identical:
221 case INTERSECTION_Contained:
222 pResult->Del_Parts();
223 break;
224
225 case INTERSECTION_Contains:
226 case INTERSECTION_Overlaps:
227 SG_Polygon_Difference(pResult, m_pB->Get_Shape(id_B));
228 break;
229 }
230 }
231
232 if( pResult->is_Valid() )
233 {
234 _Add_Polygon(pResult, id_A);
235 }
236 }
237
238 return( true );
239 }
240
241
242 ///////////////////////////////////////////////////////////
243 // //
244 ///////////////////////////////////////////////////////////
245
246 //---------------------------------------------------------
_Add_Polygon(CSG_Shape_Polygon * pPolygon,int id_A,int id_B)247 bool CPolygon_Overlay::_Add_Polygon(CSG_Shape_Polygon *pPolygon, int id_A, int id_B)
248 {
249 if( !_Fit_Polygon(pPolygon) )
250 {
251 return( false );
252 }
253
254 //-----------------------------------------------------
255 if( !m_bSplit || pPolygon->Get_Part_Count() <= 1 )
256 {
257 CSG_Shape_Polygon *pNew = _Add_Polygon(id_A, id_B);
258
259 if( pNew )
260 {
261 pNew->Assign(pPolygon, false);
262 }
263
264 return( true );
265 }
266
267 //-----------------------------------------------------
268 for(int iPart=0; iPart<pPolygon->Get_Part_Count(); iPart++)
269 {
270 CSG_Shape_Polygon *pNew = pPolygon->is_Lake(iPart) ? NULL : _Add_Polygon(id_A, id_B);
271
272 if( pNew )
273 {
274 pNew->Add_Part(pPolygon->Get_Part(iPart));
275
276 for(int jPart=0; jPart<pPolygon->Get_Part_Count(); jPart++)
277 {
278 if( pPolygon->is_Lake(jPart) && pNew->Contains(pPolygon->Get_Point(0, jPart)) )
279 {
280 pNew->Add_Part(pPolygon->Get_Part(jPart));
281 }
282 }
283 }
284 }
285
286 return( true );
287 }
288
289 //---------------------------------------------------------
_Add_Polygon(int id_A,int id_B)290 CSG_Shape_Polygon * CPolygon_Overlay::_Add_Polygon(int id_A, int id_B)
291 {
292 CSG_Shape *pOriginal, *pNew = m_pAB->Add_Shape();
293
294 if( pNew )
295 {
296 if( 1 )
297 {
298 for(int i=0; i<m_pAB->Get_Field_Count(); i++)
299 {
300 pNew->Set_NoData(i);
301 }
302 }
303
304 if( (pOriginal = m_pA->Get_Shape(id_A)) != NULL )
305 {
306 for(int i=0, j=m_bInvert ? m_pB->Get_Field_Count() : 0; i<m_pA->Get_Field_Count() && j<m_pAB->Get_Field_Count(); i++, j++)
307 {
308 if( pOriginal->is_NoData(i) ) pNew->Set_NoData(i); else *pNew->Get_Value(j) = *pOriginal->Get_Value(i);
309 }
310 }
311
312 if( (pOriginal = m_pB->Get_Shape(id_B)) != NULL )
313 {
314 for(int i=0, j=m_bInvert ? 0 : m_pA->Get_Field_Count(); i<m_pB->Get_Field_Count() && j<m_pAB->Get_Field_Count(); i++, j++)
315 {
316 if( pOriginal->is_NoData(i) ) pNew->Set_NoData(i); else *pNew->Get_Value(j) = *pOriginal->Get_Value(i);
317 }
318 }
319 }
320
321 return( (CSG_Shape_Polygon *)pNew );
322 }
323
324 //---------------------------------------------------------
_Fit_Polygon(CSG_Shape_Polygon * pPolygon)325 bool CPolygon_Overlay::_Fit_Polygon(CSG_Shape_Polygon *pPolygon)
326 {
327 for(int iPart=pPolygon->Get_Part_Count()-1; iPart>=0; iPart--)
328 {
329 if( pPolygon->Get_Area(iPart) <= 0.0 )
330 {
331 pPolygon->Del_Part(iPart);
332 }
333 else if( pPolygon->Get_Point_Count(iPart) <= 3 )
334 {
335 CSG_Point a(pPolygon->Get_Point(0, iPart)), b(pPolygon->Get_Point(1, iPart)), c(pPolygon->Get_Point(2, iPart));
336
337 if( a == b || b == c || c == a )
338 {
339 pPolygon->Del_Part(iPart);
340 }
341 }
342 }
343
344 return( pPolygon->is_Valid() );
345 }
346
347
348 ///////////////////////////////////////////////////////////
349 // //
350 // //
351 // //
352 ///////////////////////////////////////////////////////////
353
354 //---------------------------------------------------------
CPolygon_Intersection(void)355 CPolygon_Intersection::CPolygon_Intersection(void)
356 : CPolygon_Overlay(_TL("Intersect"))
357 {
358 Add_Description(_TW(
359 "Calculates the geometric intersection of the overlayed polygon layers, "
360 "i.e. layer A and layer B."
361 ));
362 }
363
364 //---------------------------------------------------------
On_Execute(void)365 bool CPolygon_Intersection::On_Execute(void)
366 {
367 CSG_Shapes *pA, *pB;
368
369 if( !CPolygon_Overlay::Initialize(&pA, &pB, true) )
370 {
371 return( false );
372 }
373
374 return( Get_Intersection(pA, pB) );
375 }
376
377
378 ///////////////////////////////////////////////////////////
379 // //
380 ///////////////////////////////////////////////////////////
381
382 //---------------------------------------------------------
CPolygon_Difference(void)383 CPolygon_Difference::CPolygon_Difference(void)
384 : CPolygon_Overlay(_TL("Difference"))
385 {
386 Add_Description(_TW(
387 "Calculates the geometric difference of the overlayed polygon layers, "
388 "i.e. layer A less layer B. Sometimes referred to as \'Erase\' command."
389 ));
390 }
391
392 //---------------------------------------------------------
On_Execute(void)393 bool CPolygon_Difference::On_Execute(void)
394 {
395 CSG_Shapes *pA, *pB;
396
397 if( !CPolygon_Overlay::Initialize(&pA, &pB, false) )
398 {
399 return( false );
400 }
401
402 return( Get_Difference(pA, pB) );
403 }
404
405
406 ///////////////////////////////////////////////////////////
407 // //
408 ///////////////////////////////////////////////////////////
409
410 //---------------------------------------------------------
CPolygon_SymDifference(void)411 CPolygon_SymDifference::CPolygon_SymDifference(void)
412 : CPolygon_Overlay(_TL("Symmetrical Difference"))
413 {
414 Add_Description(_TW(
415 "Calculates the symmetrical geometric difference of the overlayed polygon layers, "
416 "i.e. layer A less layer B plus layer B less layer A."
417 ));
418 }
419
420 //---------------------------------------------------------
On_Execute(void)421 bool CPolygon_SymDifference::On_Execute(void)
422 {
423 CSG_Shapes *pA, *pB;
424
425 if( !CPolygon_Overlay::Initialize(&pA, &pB, true) )
426 {
427 return( false );
428 }
429
430 return( Get_Difference(pA, pB, false)
431 && Get_Difference(pB, pA, true) );
432 }
433
434
435 ///////////////////////////////////////////////////////////
436 // //
437 ///////////////////////////////////////////////////////////
438
439 //---------------------------------------------------------
CPolygon_Union(void)440 CPolygon_Union::CPolygon_Union(void)
441 : CPolygon_Overlay(_TL("Union"))
442 {
443 Add_Description(_TW(
444 "Calculates the geometric union of the overlayed polygon layers, "
445 "i.e. the intersection plus the symmetrical difference of layers A and B."
446 ));
447 }
448
449 //---------------------------------------------------------
On_Execute(void)450 bool CPolygon_Union::On_Execute(void)
451 {
452 CSG_Shapes *pA, *pB;
453
454 if( !CPolygon_Overlay::Initialize(&pA, &pB, true) )
455 {
456 return( false );
457 }
458
459 return( Get_Intersection(pA, pB)
460 && Get_Difference (pA, pB, false)
461 && Get_Difference (pB, pA, true) );
462 }
463
464
465 ///////////////////////////////////////////////////////////
466 // //
467 ///////////////////////////////////////////////////////////
468
469 //---------------------------------------------------------
CPolygon_Identity(void)470 CPolygon_Identity::CPolygon_Identity(void)
471 : CPolygon_Overlay(_TL("Identity"))
472 {
473 Add_Description(_TW(
474 "Calculates the geometric intersection between both layers "
475 "and adds the difference of layer A less layer B."
476 ));
477 }
478
479 //---------------------------------------------------------
On_Execute(void)480 bool CPolygon_Identity::On_Execute(void)
481 {
482 CSG_Shapes *pA, *pB;
483
484 if( !CPolygon_Overlay::Initialize(&pA, &pB, true) )
485 {
486 return( false );
487 }
488
489 return( Get_Intersection(pA, pB)
490 && Get_Difference (pA, pB) );
491 }
492
493
494 ///////////////////////////////////////////////////////////
495 // //
496 ///////////////////////////////////////////////////////////
497
498 //---------------------------------------------------------
CPolygon_Update(void)499 CPolygon_Update::CPolygon_Update(void)
500 : CPolygon_Overlay(_TL("Update"))
501 {
502 Add_Description(_TW(
503 "Updates features of layer A with the features of layer B, "
504 "i.e. all features of layer B will be supplemented with the "
505 "difference of layer A less layer B plus. It is assumed, "
506 "that both input layers share the same attribute structure."
507 ));
508 }
509
510 //---------------------------------------------------------
On_Execute(void)511 bool CPolygon_Update::On_Execute(void)
512 {
513 CSG_Shapes *pA, *pB;
514
515 if( !CPolygon_Overlay::Initialize(&pA, &pB, false) )
516 {
517 return( false );
518 }
519
520 if( !Get_Difference(pA, pB) )
521 {
522 return( false );
523 }
524
525 //-----------------------------------------------------
526 CSG_Shapes *pAB = Parameters("RESULT")->asShapes();
527
528 for(int i=0; i<pB->Get_Count(); i++)
529 {
530 pAB->Add_Shape(pB->Get_Shape(i));
531 }
532
533 return( true );
534 }
535
536
537 ///////////////////////////////////////////////////////////
538 // //
539 // //
540 // //
541 ///////////////////////////////////////////////////////////
542
543 //---------------------------------------------------------
544