GeographicLib  1.37
PolygonArea.cpp
Go to the documentation of this file.
1 /**
2  * \file PolygonArea.cpp
3  * \brief Implementation for GeographicLib::PolygonArea class
4  *
5  * Copyright (c) Charles Karney (2010-2014) <charles@karney.com> and licensed
6  * under the MIT/X11 License. For more information, see
7  * http://geographiclib.sourceforge.net/
8  **********************************************************************/
9 
11 
12 namespace GeographicLib {
13 
14  using namespace std;
15 
16  template <class GeodType>
17  void PolygonAreaT<GeodType>::AddPoint(real lat, real lon) {
18  lon = Math::AngNormalize(lon);
19  if (_num == 0) {
20  _lat0 = _lat1 = lat;
21  _lon0 = _lon1 = lon;
22  } else {
23  real s12, S12, t;
24  _earth.GenInverse(_lat1, _lon1, lat, lon, _mask, s12, t, t, t, t, t, S12);
25  _perimetersum += s12;
26  if (!_polyline) {
27  _areasum += S12;
28  _crossings += transit(_lon1, lon);
29  }
30  _lat1 = lat; _lon1 = lon;
31  }
32  ++_num;
33  }
34 
35  template <class GeodType>
36  void PolygonAreaT<GeodType>::AddEdge(real azi, real s) {
37  if (_num) { // Do nothing if _num is zero
38  real lat, lon, S12, t;
39  _earth.GenDirect(_lat1, _lon1, azi, false, s, _mask,
40  lat, lon, t, t, t, t, t, S12);
41  _perimetersum += s;
42  if (!_polyline) {
43  _areasum += S12;
44  _crossings += transit(_lon1, lon);
45  }
46  _lat1 = lat; _lon1 = lon;
47  ++_num;
48  }
49  }
50 
51  template <class GeodType>
52  unsigned PolygonAreaT<GeodType>::Compute(bool reverse, bool sign,
53  real& perimeter, real& area) const {
54  real s12, S12, t;
55  if (_num < 2) {
56  perimeter = 0;
57  if (!_polyline)
58  area = 0;
59  return _num;
60  }
61  if (_polyline) {
62  perimeter = _perimetersum();
63  return _num;
64  }
65  _earth.GenInverse(_lat1, _lon1, _lat0, _lon0, _mask,
66  s12, t, t, t, t, t, S12);
67  perimeter = _perimetersum(s12);
68  Accumulator<> tempsum(_areasum);
69  tempsum += S12;
70  int crossings = _crossings + transit(_lon1, _lon0);
71  if (crossings & 1)
72  tempsum += (tempsum < 0 ? 1 : -1) * _area0/2;
73  // area is with the clockwise sense. If !reverse convert to
74  // counter-clockwise convention.
75  if (!reverse)
76  tempsum *= -1;
77  // If sign put area in (-area0/2, area0/2], else put area in [0, area0)
78  if (sign) {
79  if (tempsum > _area0/2)
80  tempsum -= _area0;
81  else if (tempsum <= -_area0/2)
82  tempsum += _area0;
83  } else {
84  if (tempsum >= _area0)
85  tempsum -= _area0;
86  else if (tempsum < 0)
87  tempsum += _area0;
88  }
89  area = 0 + tempsum();
90  return _num;
91  }
92 
93  template <class GeodType>
94  unsigned PolygonAreaT<GeodType>::TestPoint(real lat, real lon,
95  bool reverse, bool sign,
96  real& perimeter, real& area) const
97  {
98  if (_num == 0) {
99  perimeter = 0;
100  if (!_polyline)
101  area = 0;
102  return 1;
103  }
104  perimeter = _perimetersum();
105  real tempsum = _polyline ? 0 : _areasum();
106  int crossings = _crossings;
107  unsigned num = _num + 1;
108  for (int i = 0; i < (_polyline ? 1 : 2); ++i) {
109  real s12, S12, t;
110  _earth.GenInverse(i == 0 ? _lat1 : lat, i == 0 ? _lon1 : lon,
111  i != 0 ? _lat0 : lat, i != 0 ? _lon0 : lon,
112  _mask, s12, t, t, t, t, t, S12);
113  perimeter += s12;
114  if (!_polyline) {
115  tempsum += S12;
116  crossings += transit(i == 0 ? _lon1 : lon,
117  i != 0 ? _lon0 : lon);
118  }
119  }
120 
121  if (_polyline)
122  return num;
123 
124  if (crossings & 1)
125  tempsum += (tempsum < 0 ? 1 : -1) * _area0/2;
126  // area is with the clockwise sense. If !reverse convert to
127  // counter-clockwise convention.
128  if (!reverse)
129  tempsum *= -1;
130  // If sign put area in (-area0/2, area0/2], else put area in [0, area0)
131  if (sign) {
132  if (tempsum > _area0/2)
133  tempsum -= _area0;
134  else if (tempsum <= -_area0/2)
135  tempsum += _area0;
136  } else {
137  if (tempsum >= _area0)
138  tempsum -= _area0;
139  else if (tempsum < 0)
140  tempsum += _area0;
141  }
142  area = 0 + tempsum;
143  return num;
144  }
145 
146  template <class GeodType>
147  unsigned PolygonAreaT<GeodType>::TestEdge(real azi, real s,
148  bool reverse, bool sign,
149  real& perimeter, real& area) const {
150  if (_num == 0) { // we don't have a starting point!
151  perimeter = Math::NaN();
152  if (!_polyline)
153  area = Math::NaN();
154  return 0;
155  }
156  unsigned num = _num + 1;
157  perimeter = _perimetersum() + s;
158  if (_polyline)
159  return num;
160 
161  real tempsum = _areasum();
162  int crossings = _crossings;
163  {
164  real lat, lon, s12, S12, t;
165  _earth.GenDirect(_lat1, _lon1, azi, false, s, _mask,
166  lat, lon, t, t, t, t, t, S12);
167  tempsum += S12;
168  crossings += transit(_lon1, lon);
169  _earth.GenInverse(lat, lon, _lat0, _lon0, _mask, s12, t, t, t, t, t, S12);
170  perimeter += s12;
171  tempsum += S12;
172  crossings += transit(lon, _lon0);
173  }
174 
175  if (crossings & 1)
176  tempsum += (tempsum < 0 ? 1 : -1) * _area0/2;
177  // area is with the clockwise sense. If !reverse convert to
178  // counter-clockwise convention.
179  if (!reverse)
180  tempsum *= -1;
181  // If sign put area in (-area0/2, area0/2], else put area in [0, area0)
182  if (sign) {
183  if (tempsum > _area0/2)
184  tempsum -= _area0;
185  else if (tempsum <= -_area0/2)
186  tempsum += _area0;
187  } else {
188  if (tempsum >= _area0)
189  tempsum -= _area0;
190  else if (tempsum < 0)
191  tempsum += _area0;
192  }
193  area = 0 + tempsum;
194  return num;
195  }
196 
199 
200 } // namespace GeographicLib
static T AngNormalize(T x)
Definition: Math.hpp:399
unsigned TestEdge(real azi, real s, bool reverse, bool sign, real &perimeter, real &area) const
static T NaN()
Definition: Math.hpp:460
unsigned TestPoint(real lat, real lon, bool reverse, bool sign, real &perimeter, real &area) const
Definition: PolygonArea.cpp:94
#define GEOGRAPHICLIB_EXPORT
Definition: Constants.hpp:70
void AddEdge(real azi, real s)
Definition: PolygonArea.cpp:36
An accumulator for sums.
Definition: Accumulator.hpp:40
unsigned Compute(bool reverse, bool sign, real &perimeter, real &area) const
Definition: PolygonArea.cpp:52
Namespace for GeographicLib.
Definition: Accumulator.cpp:12
void AddPoint(real lat, real lon)
Definition: PolygonArea.cpp:17
Header for GeographicLib::PolygonArea class.