SAGA API  v9.6
tin_elements.cpp
Go to the documentation of this file.
1 
3 // //
4 // SAGA //
5 // //
6 // System for Automated Geoscientific Analyses //
7 // //
8 // Application Programming Interface //
9 // //
10 // Library: SAGA_API //
11 // //
12 //-------------------------------------------------------//
13 // //
14 // tin_elements.cpp //
15 // //
16 // Copyright (C) 2005 by Olaf Conrad //
17 // //
18 //-------------------------------------------------------//
19 // //
20 // This file is part of 'SAGA - System for Automated //
21 // Geoscientific Analyses'. //
22 // //
23 // This library is free software; you can redistribute //
24 // it and/or modify it under the terms of the GNU Lesser //
25 // General Public License as published by the Free //
26 // Software Foundation, either version 2.1 of the //
27 // License, or (at your option) any later version. //
28 // //
29 // This library is distributed in the hope that it will //
30 // be useful, but WITHOUT ANY WARRANTY; without even the //
31 // implied warranty of MERCHANTABILITY or FITNESS FOR A //
32 // PARTICULAR PURPOSE. See the GNU Lesser General Public //
33 // License for more details. //
34 // //
35 // You should have received a copy of the GNU Lesser //
36 // General Public License along with this program; if //
37 // not, see <http://www.gnu.org/licenses/>. //
38 // //
39 //-------------------------------------------------------//
40 // //
41 // contact: Olaf Conrad //
42 // Institute of Geography //
43 // University of Goettingen //
44 // Goldschmidtstr. 5 //
45 // 37077 Goettingen //
46 // Germany //
47 // //
48 // e-mail: oconrad@saga-gis.org //
49 // //
51 
52 //---------------------------------------------------------
53 #include "tin.h"
54 
55 
57 // //
58 // //
59 // //
61 
62 //---------------------------------------------------------
63 CSG_TIN_Node::CSG_TIN_Node(CSG_TIN *pOwner, sLong Index)
64  : CSG_Table_Record(pOwner, Index)
65 {
66  m_Point.x = m_Point.y = 0.;
67 
68  m_Neighbors = NULL;
69  m_nNeighbors = 0;
70 
71  m_Triangles = NULL;
72  m_nTriangles = 0;
73 }
74 
75 //---------------------------------------------------------
76 CSG_TIN_Node::~CSG_TIN_Node(void)
77 {
78  _Del_Relations();
79 }
80 
81 //---------------------------------------------------------
82 bool CSG_TIN_Node::_Add_Triangle(CSG_TIN_Triangle *pTriangle)
83 {
84  for(int i=0; i<m_nTriangles; i++)
85  {
86  if( m_Triangles[i] == pTriangle )
87  {
88  return( false );
89  }
90  }
91 
92  m_Triangles = (CSG_TIN_Triangle **)SG_Realloc(m_Triangles, ((uLong)m_nTriangles + 1) * sizeof(CSG_TIN_Triangle *));
93  m_Triangles[m_nTriangles++] = pTriangle;
94 
95 // _Add_Neighbor(pTriangle->Get_Point());
96 // _Add_Neighbor(pTriangle->Get_Point(1));
97 // _Add_Neighbor(pTriangle->Get_Point(2));
98 
99  return( true );
100 }
101 
102 //---------------------------------------------------------
103 bool CSG_TIN_Node::_Add_Neighbor(CSG_TIN_Node *pNeighbor)
104 {
105  if( pNeighbor == this )
106  {
107  return( false );
108  }
109 
110  for(int i=0; i<m_nNeighbors; i++)
111  {
112  if( m_Neighbors[i] == pNeighbor )
113  {
114  return( false );
115  }
116  }
117 
118  m_Neighbors = (CSG_TIN_Node **)SG_Realloc(m_Neighbors, ((uLong)m_nNeighbors + 1) * sizeof(CSG_TIN_Node *));
119  m_Neighbors[m_nNeighbors++] = pNeighbor;
120 
121  return( true );
122 }
123 
124 //---------------------------------------------------------
125 bool CSG_TIN_Node::_Del_Relations(void)
126 {
127  if( m_nTriangles > 0 )
128  {
129  SG_Free(m_Triangles);
130  m_Triangles = NULL;
131  m_nTriangles = 0;
132  }
133 
134  if( m_nNeighbors > 0 )
135  {
136  SG_Free(m_Neighbors);
137  m_Neighbors = NULL;
138  m_nNeighbors = 0;
139  }
140 
141  return( true );
142 }
143 
144 //---------------------------------------------------------
145 double CSG_TIN_Node::Get_Gradient(int iNeighbor, int iField)
146 {
147  double dx, dy, dz;
148  CSG_TIN_Node *pNeighbor;
149 
150  if( (pNeighbor = Get_Neighbor(iNeighbor)) != NULL )
151  {
152  dx = Get_X() - pNeighbor->Get_X();
153  dy = Get_Y() - pNeighbor->Get_Y();
154  dz = asDouble(iField) - pNeighbor->asDouble(iField);
155 
156  if( (dx = sqrt(dx*dx + dy*dy)) > 0.0 )
157  {
158  return( dz / dx );
159  }
160  }
161 
162  return( 0.0 );
163 }
164 
165 //---------------------------------------------------------
166 int SG_TIN_Compare_Triangle_Center(const void *pz1, const void *pz2)
167 {
168  double z1 = ((TSG_Point_3D *)pz1)->z;
169  double z2 = ((TSG_Point_3D *)pz2)->z;
170 
171  return( z1 < z2 ? -1 : z1 > z2 ? 1 : 0 );
172 }
173 
174 //---------------------------------------------------------
176 {
177  Points.Clear();
178 
179  if( m_nTriangles >= 3 )
180  {
181  int i; CSG_Points_3D p;
182 
183  for(i=0; i<m_nTriangles; i++)
184  {
185  TSG_Point c = m_Triangles[i]->Get_CircumCircle_Point();
186 
187  p.Add(c.x, c.y, atan2(m_Point.y - c.y, m_Point.x - c.x));
188  }
189 
190  qsort(&(p[0]), p.Get_Count(), sizeof(TSG_Point_3D), SG_TIN_Compare_Triangle_Center);
191 
192  for(i=0; i<m_nTriangles; i++)
193  {
194  Points.Add(p[i].x, p[i].y);
195  }
196 
197  return( true );
198  }
199 
200  return( false );
201 }
202 
203 //---------------------------------------------------------
205 {
206  CSG_Points Points;
207 
208  if( Get_Polygon(Points) )
209  {
210  return( SG_Get_Polygon_Area(Points) );
211  }
212 
213  return( 0.0 );
214 }
215 
216 
218 // //
219 // //
220 // //
222 
223 //---------------------------------------------------------
224 CSG_TIN_Edge::CSG_TIN_Edge(CSG_TIN_Node *a, CSG_TIN_Node *b)
225 {
226  m_Nodes[0] = a;
227  m_Nodes[1] = b;
228 }
229 
230 //---------------------------------------------------------
231 CSG_TIN_Edge::~CSG_TIN_Edge(void)
232 {
233 }
234 
235 
237 // //
238 // //
239 // //
241 
242 //---------------------------------------------------------
243 CSG_TIN_Triangle::CSG_TIN_Triangle(CSG_TIN_Node *a, CSG_TIN_Node *b, CSG_TIN_Node *c)
244 {
245  m_Nodes[0] = a;
246  m_Nodes[1] = b;
247  m_Nodes[2] = c;
248 
249  //-----------------------------------------------------
250  double xMin, yMin, xMax, yMax;
251 
252  xMin = xMax = a->Get_X();
253  yMin = yMax = a->Get_Y();
254 
255  if( xMin > b->Get_X() )
256  xMin = b->Get_X();
257  else if( xMax < b->Get_X() )
258  xMax = b->Get_X();
259 
260  if( yMin > b->Get_Y() )
261  yMin = b->Get_Y();
262  else if( yMax < b->Get_Y() )
263  yMax = b->Get_Y();
264 
265  if( xMin > c->Get_X() )
266  xMin = c->Get_X();
267  else if( xMax < c->Get_X() )
268  xMax = c->Get_X();
269 
270  if( yMin > c->Get_Y() )
271  yMin = c->Get_Y();
272  else if( yMax < c->Get_Y() )
273  yMax = c->Get_Y();
274 
275  m_Extent.Assign(xMin, yMin, xMax, yMax);
276 
277  //-----------------------------------------------------
278  m_Area = fabs( a->Get_X() * (b->Get_Y() - c->Get_Y())
279  + b->Get_X() * (c->Get_Y() - a->Get_Y())
280  + c->Get_X() * (a->Get_Y() - b->Get_Y()) ) / 2.0;
281 
282  //-----------------------------------------------------
283  TSG_Point Points[3];
284 
285  Points[0] = m_Nodes[0]->Get_Point();
286  Points[1] = m_Nodes[1]->Get_Point();
287  Points[2] = m_Nodes[2]->Get_Point();
288 
289  SG_Get_Triangle_CircumCircle(Points, m_Center, m_Radius);
290 }
291 
292 //---------------------------------------------------------
293 CSG_TIN_Triangle::~CSG_TIN_Triangle(void)
294 {}
295 
296 
297 //---------------------------------------------------------
299 {
300  return( is_Containing(Point.x, Point.y) );
301 }
302 
303 //---------------------------------------------------------
304 #define IS_ONLINE(A, B) (A.y == B.y && ((A.x <= x && x <= B.x) || (B.x <= x && x <= A.x)))
305 
306 //---------------------------------------------------------
307 bool CSG_TIN_Triangle::is_Containing(double x, double y)
308 {
309  if( m_Extent.Contains(x, y) )
310  {
311  int nCrossings;
312  TSG_Point A, B, C;
313 
314  if( (x == m_Nodes[0]->Get_Point().x && y == m_Nodes[0]->Get_Point().y)
315  || (x == m_Nodes[1]->Get_Point().x && y == m_Nodes[1]->Get_Point().y)
316  || (x == m_Nodes[2]->Get_Point().x && y == m_Nodes[2]->Get_Point().y) )
317  return( true );
318 
319  if( y == m_Extent.Get_YMin() || y == m_Extent.Get_YMax() )
320  {
321  if( IS_ONLINE(m_Nodes[0]->Get_Point(), m_Nodes[1]->Get_Point())
322  || IS_ONLINE(m_Nodes[1]->Get_Point(), m_Nodes[2]->Get_Point())
323  || IS_ONLINE(m_Nodes[2]->Get_Point(), m_Nodes[0]->Get_Point()) )
324  return( true );
325  }
326 
327  nCrossings = 0;
328 
329  if( (y == m_Nodes[0]->Get_Point().y && x > m_Nodes[0]->Get_Point().x)
330  || (y == m_Nodes[1]->Get_Point().y && x > m_Nodes[1]->Get_Point().x)
331  || (y == m_Nodes[2]->Get_Point().y && x > m_Nodes[2]->Get_Point().x) )
332  nCrossings = -1;
333 
334  A.x = m_Extent.xMin - 1.0;
335  B.x = x;
336  A.y = B.y = y;
337 
338  if( SG_Get_Crossing(C, m_Nodes[0]->Get_Point(), m_Nodes[1]->Get_Point(), A, B) )
339  nCrossings++;
340 
341  if( SG_Get_Crossing(C, m_Nodes[1]->Get_Point(), m_Nodes[2]->Get_Point(), A, B) )
342  nCrossings++;
343 
344  if( SG_Get_Crossing(C, m_Nodes[2]->Get_Point(), m_Nodes[0]->Get_Point(), A, B) )
345  nCrossings++;
346 
347  return( nCrossings == 1 );
348  }
349 
350  return( false );
351 }
352 
353 //---------------------------------------------------------
354 bool CSG_TIN_Triangle::Get_Value(int zField, const TSG_Point &p, double &z)
355 {
356  return( Get_Value(zField, p.x, p.y, z) );
357 }
358 
359 bool CSG_TIN_Triangle::Get_Value(int zField, double x, double y, double &z)
360 {
361  CSG_Vector B, Z(3);
362  CSG_Matrix M(3, 3), Mt;
363 
364  for(int i=0; i<3; i++)
365  {
366  M[i][0] = 1.0;
367  M[i][1] = m_Nodes[i]->Get_X();
368  M[i][2] = m_Nodes[i]->Get_Y();
369  Z[i] = m_Nodes[i]->asDouble(zField);
370  }
371 
372  Mt = M.Get_Transpose();
373  B = (Mt * M).Get_Inverse() * (Mt * Z);
374 
375  z = B[0] + B[1] * x + B[2] * y;
376 
377  return( true );
378 }
379 
380 //---------------------------------------------------------
381 bool CSG_TIN_Triangle::Get_Gradient(int zField, double &Decline, double &Azimuth)
382 {
383  int i;
384  double x[3], y[3], z[3], A, B, C;
385 
386  for(i=0; i<3; i++)
387  {
388  x[i] = m_Nodes[i]->Get_X();
389  y[i] = m_Nodes[i]->Get_Y();
390  z[i] = m_Nodes[i]->asDouble(zField);
391  }
392 
393  A = z[0] * (x[1] - x[2]) + z[1] * (x[2] - x[0]) + z[2] * (x[0] - x[1]);
394  B = y[0] * (z[1] - z[2]) + y[1] * (z[2] - z[0]) + y[2] * (z[0] - z[1]);
395  C = x[0] * (y[1] - y[2]) + x[1] * (y[2] - y[0]) + x[2] * (y[0] - y[1]);
396 
397  if( C != 0.0 )
398  {
399  A = - A / C;
400  B = - B / C;
401 
402  Decline = atan(sqrt(A*A + B*B));
403 
404  if( A != 0.0 )
405  Azimuth = M_PI_180 + atan2(B, A);
406  else
407  Azimuth = B > 0.0 ? M_PI_270 : (B < 0.0 ? M_PI_090 : -1.0);
408 
409  return( true );
410  }
411 
412  Decline = -1.0;
413  Azimuth = -1.0;
414 
415  return( false );
416 }
417 
418 
420 // //
421 // //
422 // //
424 
425 //---------------------------------------------------------
CSG_Table_Record::asDouble
double asDouble(int Field) const
Definition: table_record.cpp:527
CSG_Rect::Assign
CSG_Rect & Assign(double xMin, double yMin, double xMax, double yMax)
Definition: geo_classes.cpp:727
CSG_Points::Add
bool Add(double x, double y)
Definition: geo_classes.cpp:350
CSG_Points_3D
Definition: geo_tools.h:320
CSG_TIN::m_Triangles
CSG_TIN_Triangle ** m_Triangles
Definition: tin.h:299
SG_Get_Triangle_CircumCircle
bool SG_Get_Triangle_CircumCircle(TSG_Point Triangle[3], TSG_Point &Point, double &Radius)
Definition: geo_functions.cpp:465
CSG_Points_3D::Add
bool Add(double x, double y, double z)
Definition: geo_classes.cpp:582
CSG_TIN_Node::Get_X
double Get_X(void)
Definition: tin.h:101
A
#define A
CSG_Table_Record
Definition: table.h:130
SG_Get_Polygon_Area
double SG_Get_Polygon_Area(TSG_Point *Points, int nPoints)
Definition: geo_functions.cpp:514
C
#define C
CSG_TIN::m_nTriangles
sLong m_nTriangles
Definition: tin.h:293
SSG_Point_3D
Definition: geo_tools.h:264
CSG_Points_3D::Get_Count
sLong Get_Count(void) const
Definition: geo_tools.h:336
SG_Free
SAGA_API_DLL_EXPORT void SG_Free(void *memblock)
Definition: api_memory.cpp:83
IS_ONLINE
#define IS_ONLINE(A, B)
Definition: tin_elements.cpp:304
CSG_TIN_Node::Get_Polygon_Area
double Get_Polygon_Area(void)
Definition: tin_elements.cpp:204
tin.h
SSG_Point
Definition: geo_tools.h:128
SSG_Rect::xMin
double xMin
Definition: geo_tools.h:465
CSG_Rect::Contains
bool Contains(double x, double y) const
Definition: geo_classes.cpp:870
CSG_TIN_Triangle::Get_CircumCircle_Point
TSG_Point Get_CircumCircle_Point(void)
Definition: tin.h:194
CSG_TIN
Definition: tin.h:222
CSG_Vector
Definition: mat_tools.h:360
CSG_Rect::Get_YMin
double Get_YMin(void) const
Definition: geo_tools.h:504
SG_TIN_Compare_Triangle_Center
int SG_TIN_Compare_Triangle_Center(const void *pz1, const void *pz2)
Definition: tin_elements.cpp:166
sLong
signed long long sLong
Definition: api_core.h:158
CSG_TIN_Node::Get_Y
double Get_Y(void)
Definition: tin.h:102
M_PI_180
#define M_PI_180
Definition: mat_tools.h:102
CSG_Points
Definition: geo_tools.h:184
CSG_Matrix::Get_Transpose
CSG_Matrix Get_Transpose(void) const
Definition: mat_matrix.cpp:1829
M_PI_270
#define M_PI_270
Definition: mat_tools.h:104
CSG_TIN_Node::Get_Polygon
bool Get_Polygon(CSG_Points &Points)
Definition: tin_elements.cpp:175
CSG_TIN_Triangle::Get_Value
bool Get_Value(int zField, const TSG_Point &p, double &z)
Definition: tin_elements.cpp:354
CSG_TIN_Triangle
Definition: tin.h:173
CSG_TIN_Triangle::is_Containing
bool is_Containing(const TSG_Point &p)
Definition: tin_elements.cpp:298
B
#define B
M_PI_090
#define M_PI_090
Definition: mat_tools.h:100
CSG_TIN_Node::Get_Neighbor
CSG_TIN_Node * Get_Neighbor(int iNeighbor)
Definition: tin.h:105
SSG_Point::x
double x
Definition: geo_tools.h:129
CSG_TIN_Node
Definition: tin.h:93
CSG_Rect::Get_YMax
double Get_YMax(void) const
Definition: geo_tools.h:505
CSG_TIN_Node::Get_Gradient
double Get_Gradient(int iNeighbor, int iField)
Definition: tin_elements.cpp:145
SSG_Point::y
double y
Definition: geo_tools.h:129
CSG_Points::Clear
bool Clear(void)
Definition: geo_tools.h:190
SG_Realloc
SAGA_API_DLL_EXPORT void * SG_Realloc(void *memblock, size_t size)
Definition: api_memory.cpp:77
SG_Get_Crossing
bool SG_Get_Crossing(TSG_Point &Crossing, const TSG_Point &a1, const TSG_Point &a2, const TSG_Point &b1, const TSG_Point &b2, bool bExactMatch)
Definition: geo_functions.cpp:280
CSG_Matrix
Definition: mat_tools.h:478
CSG_TIN_Triangle::Get_Gradient
bool Get_Gradient(int zField, double &Decline, double &Azimuth)
Definition: tin_elements.cpp:381
uLong
unsigned long long uLong
Definition: api_core.h:159
CSG_TIN_Node::Get_Point
const TSG_Point & Get_Point(void)
Definition: tin.h:100