GeographicLib  2.1.2
Planimeter.cpp
Go to the documentation of this file.
1 /**
2  * \file Planimeter.cpp
3  * \brief Command line utility for measuring the area of geodesic polygons
4  *
5  * Copyright (c) Charles Karney (2010-2022) <charles@karney.com> and licensed
6  * under the MIT/X11 License. For more information, see
7  * https://geographiclib.sourceforge.io/
8  *
9  * See the <a href="Planimeter.1.html">man page</a> for usage information.
10  **********************************************************************/
11 
12 #include <iostream>
13 #include <string>
14 #include <sstream>
15 #include <fstream>
17 #include <GeographicLib/DMS.hpp>
21 
22 #if defined(_MSC_VER)
23 // Squelch warnings about constant conditional expressions
24 # pragma warning (disable: 4127)
25 #endif
26 
27 #include "Planimeter.usage"
28 
29 int main(int argc, const char* const argv[]) {
30  try {
31  using namespace GeographicLib;
32  typedef Math::real real;
34  enum { GEODESIC, EXACT, AUTHALIC, RHUMB };
35  real
36  a = Constants::WGS84_a(),
37  f = Constants::WGS84_f();
38  bool reverse = false, sign = true, polyline = false, longfirst = false,
39  geoconvert_compat = false;
40  int linetype = GEODESIC;
41  int prec = 6;
42  std::string istring, ifile, ofile, cdelim;
43  char lsep = ';';
44 
45  for (int m = 1; m < argc; ++m) {
46  std::string arg(argv[m]);
47  if (arg == "-r")
48  reverse = !reverse;
49  else if (arg == "-s")
50  sign = !sign;
51  else if (arg == "-l")
52  polyline = !polyline;
53  else if (arg == "-e") {
54  if (m + 2 >= argc) return usage(1, true);
55  try {
56  a = Utility::val<real>(std::string(argv[m + 1]));
57  f = Utility::fract<real>(std::string(argv[m + 2]));
58  }
59  catch (const std::exception& e) {
60  std::cerr << "Error decoding arguments of -e: " << e.what() << "\n";
61  return 1;
62  }
63  m += 2;
64  } else if (arg == "-w")
65  longfirst = !longfirst;
66  else if (arg == "-p") {
67  if (++m == argc) return usage(1, true);
68  try {
69  prec = Utility::val<int>(std::string(argv[m]));
70  }
71  catch (const std::exception&) {
72  std::cerr << "Precision " << argv[m] << " is not a number\n";
73  return 1;
74  }
75  } else if (arg == "-G")
76  linetype = GEODESIC;
77  else if (arg == "-E")
78  linetype = EXACT;
79  else if (arg == "-Q")
80  linetype = AUTHALIC;
81  else if (arg == "-R")
82  linetype = RHUMB;
83  else if (arg == "--geoconvert-input")
84  geoconvert_compat = true;
85  else if (arg == "--input-string") {
86  if (++m == argc) return usage(1, true);
87  istring = argv[m];
88  } else if (arg == "--input-file") {
89  if (++m == argc) return usage(1, true);
90  ifile = argv[m];
91  } else if (arg == "--output-file") {
92  if (++m == argc) return usage(1, true);
93  ofile = argv[m];
94  } else if (arg == "--line-separator") {
95  if (++m == argc) return usage(1, true);
96  if (std::string(argv[m]).size() != 1) {
97  std::cerr << "Line separator must be a single character\n";
98  return 1;
99  }
100  lsep = argv[m][0];
101  } else if (arg == "--comment-delimiter") {
102  if (++m == argc) return usage(1, true);
103  cdelim = argv[m];
104  } else if (arg == "--version") {
105  std::cout << argv[0] << ": GeographicLib version "
106  << GEOGRAPHICLIB_VERSION_STRING << "\n";
107  return 0;
108  } else
109  return usage(!(arg == "-h" || arg == "--help"), arg != "--help");
110  }
111 
112  if (!ifile.empty() && !istring.empty()) {
113  std::cerr << "Cannot specify --input-string and --input-file together\n";
114  return 1;
115  }
116  if (ifile == "-") ifile.clear();
117  std::ifstream infile;
118  std::istringstream instring;
119  if (!ifile.empty()) {
120  infile.open(ifile.c_str());
121  if (!infile.is_open()) {
122  std::cerr << "Cannot open " << ifile << " for reading\n";
123  return 1;
124  }
125  } else if (!istring.empty()) {
126  std::string::size_type m = 0;
127  while (true) {
128  m = istring.find(lsep, m);
129  if (m == std::string::npos)
130  break;
131  istring[m] = '\n';
132  }
133  instring.str(istring);
134  }
135  std::istream* input = !ifile.empty() ? &infile :
136  (!istring.empty() ? &instring : &std::cin);
137 
138  std::ofstream outfile;
139  if (ofile == "-") ofile.clear();
140  if (!ofile.empty()) {
141  outfile.open(ofile.c_str());
142  if (!outfile.is_open()) {
143  std::cerr << "Cannot open " << ofile << " for writing\n";
144  return 1;
145  }
146  }
147  std::ostream* output = !ofile.empty() ? &outfile : &std::cout;
148 
149  const Ellipsoid ellip(a, f);
150  if (linetype == AUTHALIC) {
151  using std::sqrt;
152  a = sqrt(ellip.Area() / (4 * Math::pi()));
153  f = 0;
154  }
155  const Geodesic geod(a, f);
156  const GeodesicExact geode(a, f);
157  const Rhumb rhumb(a, f);
158  PolygonArea poly(geod, polyline);
159  PolygonAreaExact polye(geode, polyline);
160  PolygonAreaRhumb polyr(rhumb, polyline);
161  GeoCoords p;
162 
163  // Max precision = 10: 0.1 nm in distance, 10^-15 deg (= 0.11 nm),
164  // 10^-11 sec (= 0.3 nm).
165  prec = std::min(10 + Math::extra_digits(), std::max(0, prec));
166  std::string s, eol("\n");
167  real perimeter, area;
168  unsigned num;
169  std::istringstream str;
170  std::string slat, slon, junk;
171  real lat = 0, lon = 0;
172  while (std::getline(*input, s)) {
173  if (!cdelim.empty()) {
174  std::string::size_type m = s.find(cdelim);
175  if (m != std::string::npos) {
176  eol = " " + s.substr(m) + "\n";
177  s = s.substr(0, m);
178  }
179  }
180  bool endpoly = s.empty();
181  if (!endpoly) {
182  try {
183  using std::isnan;
184  if (geoconvert_compat) {
185  p.Reset(s, true, longfirst);
186  lat = p.Latitude(); lon = p.Longitude();
187  } else {
188  str.clear(); str.str(s);
189  if (!(str >> slat >> slon))
190  throw GeographicErr("incomplete input");
191  if (str >> junk)
192  throw GeographicErr("extra input");
193  DMS::DecodeLatLon(slat, slon, lat, lon, longfirst);
194  }
195  if (isnan(lat) || isnan(lon))
196  endpoly = true;
197  }
198  catch (const GeographicErr&) {
199  endpoly = true;
200  }
201  }
202  if (endpoly) {
203  num =
204  linetype == EXACT ? polye.Compute(reverse, sign, perimeter, area) :
205  linetype == RHUMB ? polyr.Compute(reverse, sign, perimeter, area) :
206  poly.Compute(reverse, sign, perimeter, area); // geodesic + authalic
207  if (num > 0) {
208  *output << num << " " << Utility::str(perimeter, prec);
209  if (!polyline) {
210  *output << " " << Utility::str(area, std::max(0, prec - 5));
211  }
212  *output << eol;
213  }
214  linetype == EXACT ? polye.Clear() :
215  linetype == RHUMB ? polyr.Clear() : poly.Clear();
216  eol = "\n";
217  } else {
218  linetype == EXACT ? polye.AddPoint(lat, lon) :
219  linetype == RHUMB ? polyr.AddPoint(lat, lon) :
220  poly.AddPoint(linetype == AUTHALIC ? ellip.AuthalicLatitude(lat) :
221  lat, lon);
222  }
223  }
224  num =
225  linetype == EXACT ? polye.Compute(reverse, sign, perimeter, area) :
226  linetype == RHUMB ? polyr.Compute(reverse, sign, perimeter, area) :
227  poly.Compute(reverse, sign, perimeter, area);
228  if (num > 0) {
229  *output << num << " " << Utility::str(perimeter, prec);
230  if (!polyline) {
231  *output << " " << Utility::str(area, std::max(0, prec - 5));
232  }
233  *output << eol;
234  }
235  linetype == EXACT ? polye.Clear() :
236  linetype == RHUMB ? polyr.Clear() : poly.Clear();
237  eol = "\n";
238  return 0;
239  }
240  catch (const std::exception& e) {
241  std::cerr << "Caught exception: " << e.what() << "\n";
242  return 1;
243  }
244  catch (...) {
245  std::cerr << "Caught unknown exception\n";
246  return 1;
247  }
248 }
Header for GeographicLib::DMS class.
Header for GeographicLib::Ellipsoid class.
Header for GeographicLib::GeoCoords class.
GeographicLib::Math::real real
Definition: GeodSolve.cpp:31
int main(int argc, const char *const argv[])
Definition: Planimeter.cpp:29
Header for GeographicLib::PolygonAreaT class.
Header for GeographicLib::Utility class.
static void DecodeLatLon(const std::string &dmsa, const std::string &dmsb, real &lat, real &lon, bool longfirst=false)
Definition: DMS.cpp:368
Properties of an ellipsoid.
Definition: Ellipsoid.hpp:39
Math::real Area() const
Definition: Ellipsoid.cpp:45
Math::real AuthalicLatitude(real phi) const
Definition: Ellipsoid.cpp:77
Conversion between geographic coordinates.
Definition: GeoCoords.hpp:49
void Reset(const std::string &s, bool centerp=true, bool longfirst=false)
Definition: GeoCoords.cpp:19
Math::real Latitude() const
Definition: GeoCoords.hpp:276
Math::real Longitude() const
Definition: GeoCoords.hpp:281
Exact geodesic calculations.
Geodesic calculations
Definition: Geodesic.hpp:172
Exception handling for GeographicLib.
Definition: Constants.hpp:316
static T pi()
Definition: Math.hpp:190
static int extra_digits()
Definition: Math.cpp:51
Solve of the direct and inverse rhumb problems.
Definition: Rhumb.hpp:66
static int set_digits(int ndigits=0)
Definition: Utility.cpp:184
static std::string str(T x, int p=-1)
Definition: Utility.hpp:161
Namespace for GeographicLib.
Definition: Accumulator.cpp:12