GeographicLib
1.21
|
00001 /** 00002 * \file GeoidEval.cpp 00003 * \brief Command line utility for evaluating geoid heights 00004 * 00005 * Copyright (c) Charles Karney (2009-2012) <charles@karney.com> and licensed 00006 * under the MIT/X11 License. For more information, see 00007 * http://geographiclib.sourceforge.net/ 00008 * 00009 * Compile and link with 00010 * g++ -g -O3 -I../include -I../man -o GeoidEval \ 00011 * GeoidEval.cpp \ 00012 * ../src/DMS.cpp \ 00013 * ../src/GeoCoords.cpp \ 00014 * ../src/Geoid.cpp \ 00015 * ../src/MGRS.cpp \ 00016 * ../src/PolarStereographic.cpp \ 00017 * ../src/TransverseMercator.cpp \ 00018 * ../src/UTMUPS.cpp 00019 * 00020 * See the <a href="GeoidEval.1.html">man page</a> for usage 00021 * information. 00022 **********************************************************************/ 00023 00024 #include <iostream> 00025 #include <string> 00026 #include <sstream> 00027 #include <fstream> 00028 #include <GeographicLib/Geoid.hpp> 00029 #include <GeographicLib/DMS.hpp> 00030 #include <GeographicLib/Utility.hpp> 00031 #include <GeographicLib/GeoCoords.hpp> 00032 00033 #include "GeoidEval.usage" 00034 00035 int main(int argc, char* argv[]) { 00036 try { 00037 using namespace GeographicLib; 00038 typedef Math::real real; 00039 bool cacheall = false, cachearea = false, verbose = false, 00040 cubic = true, gradp = false; 00041 real caches, cachew, cachen, cachee; 00042 std::string dir; 00043 std::string geoid = Geoid::DefaultGeoidName(); 00044 Geoid::convertflag heightmult = Geoid::NONE; 00045 std::string istring, ifile, ofile, cdelim; 00046 char lsep = ';'; 00047 bool northp = false; 00048 int zonenum = UTMUPS::INVALID; 00049 00050 for (int m = 1; m < argc; ++m) { 00051 std::string arg(argv[m]); 00052 if (arg == "-a") { 00053 cacheall = true; 00054 cachearea = false; 00055 } 00056 else if (arg == "-c") { 00057 if (m + 4 >= argc) return usage(1, true); 00058 cacheall = false; 00059 cachearea = true; 00060 try { 00061 DMS::DecodeLatLon(std::string(argv[m + 1]), std::string(argv[m + 2]), 00062 caches, cachew); 00063 DMS::DecodeLatLon(std::string(argv[m + 3]), std::string(argv[m + 4]), 00064 cachen, cachee); 00065 } 00066 catch (const std::exception& e) { 00067 std::cerr << "Error decoding argument of -c: " << e.what() << "\n"; 00068 return 1; 00069 } 00070 m += 4; 00071 } else if (arg == "--msltohae") 00072 heightmult = Geoid::GEOIDTOELLIPSOID; 00073 else if (arg == "--haetomsl") 00074 heightmult = Geoid::ELLIPSOIDTOGEOID; 00075 else if (arg == "-z") { 00076 if (++m == argc) return usage(1, true); 00077 std::string zone = argv[m]; 00078 try { 00079 UTMUPS::DecodeZone(zone, zonenum, northp); 00080 } 00081 catch (const std::exception& e) { 00082 std::cerr << "Error decoding zone: " << e.what() << "\n"; 00083 return 1; 00084 } 00085 if (!(zonenum >= UTMUPS::MINZONE && zonenum <= UTMUPS::MAXZONE)) { 00086 std::cerr << "Illegal zone " << zone << "\n"; 00087 return 1; 00088 } 00089 } else if (arg == "-n") { 00090 if (++m == argc) return usage(1, true); 00091 geoid = argv[m]; 00092 } else if (arg == "-d") { 00093 if (++m == argc) return usage(1, true); 00094 dir = argv[m]; 00095 } else if (arg == "-l") 00096 cubic = false; 00097 else if (arg == "-g") 00098 gradp = true; 00099 else if (arg == "-v") 00100 verbose = true; 00101 else if (arg == "--input-string") { 00102 if (++m == argc) return usage(1, true); 00103 istring = argv[m]; 00104 } else if (arg == "--input-file") { 00105 if (++m == argc) return usage(1, true); 00106 ifile = argv[m]; 00107 } else if (arg == "--output-file") { 00108 if (++m == argc) return usage(1, true); 00109 ofile = argv[m]; 00110 } else if (arg == "--line-separator") { 00111 if (++m == argc) return usage(1, true); 00112 if (std::string(argv[m]).size() != 1) { 00113 std::cerr << "Line separator must be a single character\n"; 00114 return 1; 00115 } 00116 lsep = argv[m][0]; 00117 } else if (arg == "--comment-delimiter") { 00118 if (++m == argc) return usage(1, true); 00119 cdelim = argv[m]; 00120 } else if (arg == "--version") { 00121 std::cout 00122 << argv[0] 00123 << ": $Id: 6db1ff0b8309a39d0d9b0250dd73be964c5efb7c $\n" 00124 << "GeographicLib version " << GEOGRAPHICLIB_VERSION_STRING << "\n"; 00125 return 0; 00126 } else { 00127 int retval = usage(!(arg == "-h" || arg == "--help"), arg != "--help"); 00128 if (arg == "-h") 00129 std::cout<< "\nDefault geoid path = \"" << Geoid::DefaultGeoidPath() 00130 << "\"\nDefault geoid name = \"" << Geoid::DefaultGeoidName() 00131 << "\"\n"; 00132 return retval; 00133 } 00134 } 00135 00136 if (!ifile.empty() && !istring.empty()) { 00137 std::cerr << "Cannot specify --input-string and --input-file together\n"; 00138 return 1; 00139 } 00140 if (ifile == "-") ifile.clear(); 00141 std::ifstream infile; 00142 std::istringstream instring; 00143 if (!ifile.empty()) { 00144 infile.open(ifile.c_str()); 00145 if (!infile.is_open()) { 00146 std::cerr << "Cannot open " << ifile << " for reading\n"; 00147 return 1; 00148 } 00149 } else if (!istring.empty()) { 00150 std::string::size_type m = 0; 00151 while (true) { 00152 m = istring.find(lsep, m); 00153 if (m == std::string::npos) 00154 break; 00155 istring[m] = '\n'; 00156 } 00157 instring.str(istring); 00158 } 00159 std::istream* input = !ifile.empty() ? &infile : 00160 (!istring.empty() ? &instring : &std::cin); 00161 00162 std::ofstream outfile; 00163 if (ofile == "-") ofile.clear(); 00164 if (!ofile.empty()) { 00165 outfile.open(ofile.c_str()); 00166 if (!outfile.is_open()) { 00167 std::cerr << "Cannot open " << ofile << " for writing\n"; 00168 return 1; 00169 } 00170 } 00171 std::ostream* output = !ofile.empty() ? &outfile : &std::cout; 00172 00173 int retval = 0; 00174 try { 00175 const Geoid g(geoid, dir, cubic); 00176 try { 00177 if (cacheall) 00178 g.CacheAll(); 00179 else if (cachearea) 00180 g.CacheArea(caches, cachew, cachen, cachee); 00181 } 00182 catch (const std::exception& e) { 00183 std::cerr << "ERROR: " << e.what() << "\nProceeding without a cache\n"; 00184 } 00185 if (verbose) { 00186 std::cerr << "Geoid file: " << g.GeoidFile() << "\n" 00187 << "Description: " << g.Description() << "\n" 00188 << "Interpolation: " << g.Interpolation() << "\n" 00189 << "Date & Time: " << g.DateTime() << "\n" 00190 << "Offset (m): " << g.Offset() << "\n" 00191 << "Scale (m): " << g.Scale() << "\n" 00192 << "Max error (m): " << g.MaxError() << "\n" 00193 << "RMS error (m): " << g.RMSError() << "\n"; 00194 if (g.Cache()) 00195 std::cerr<< "Caching:" 00196 << "\n SW Corner: " << g.CacheSouth() << " " << g.CacheWest() 00197 << "\n NE Corner: " << g.CacheNorth() << " " << g.CacheEast() 00198 << "\n"; 00199 } 00200 00201 GeoCoords p; 00202 std::string s, suff; 00203 const char* spaces = " \t\n\v\f\r,"; // Include comma as space 00204 while (std::getline(*input, s)) { 00205 try { 00206 std::string eol("\n"); 00207 if (!cdelim.empty()) { 00208 std::string::size_type m = s.find(cdelim); 00209 if (m != std::string::npos) { 00210 eol = " " + s.substr(m) + "\n"; 00211 std::string::size_type m1 = 00212 m > 0 ? s.find_last_not_of(spaces, m - 1) : std::string::npos; 00213 s = s.substr(0, m1 != std::string::npos ? m1 + 1 : m); 00214 } 00215 } 00216 real height = 0; 00217 if (zonenum != UTMUPS::INVALID) { 00218 // Expect "easting northing" if heightmult == 0, or 00219 // "easting northing height" if heightmult != 0. 00220 std::string::size_type pa = 0, pb = 0; 00221 real easting = 0, northing = 0; 00222 for (int i = 0; i < (heightmult ? 3 : 2); ++i) { 00223 if (pb == std::string::npos) 00224 throw GeographicErr("Incomplete input: " + s); 00225 // Start of i'th token 00226 pa = s.find_first_not_of(spaces, pb); 00227 if (pa == std::string::npos) 00228 throw GeographicErr("Incomplete input: " + s); 00229 // End of i'th token 00230 pb = s.find_first_of(spaces, pa); 00231 (i == 2 ? height : (i == 0 ? easting : northing)) = 00232 Utility::num<real>(s.substr(pa, (pb == std::string::npos ? 00233 pb : pb - pa))); 00234 } 00235 p.Reset(zonenum, northp, easting, northing); 00236 if (heightmult) { 00237 suff = pb == std::string::npos ? "" : s.substr(pb); 00238 s = s.substr(0, pa); 00239 } 00240 } else { 00241 if (heightmult) { 00242 // Treat last token as height 00243 // pb = last char of last token 00244 // pa = last char preceding white space 00245 // px = last char of 2nd last token 00246 std::string::size_type pb = s.find_last_not_of(spaces); 00247 std::string::size_type pa = s.find_last_of(spaces, pb); 00248 if (pa == std::string::npos || pb == std::string::npos) 00249 throw GeographicErr("Incomplete input: " + s); 00250 height = Utility::num<real>(s.substr(pa + 1, pb - pa)); 00251 s = s.substr(0, pa + 1); 00252 } 00253 p.Reset(s); 00254 } 00255 if (heightmult) { 00256 real h = g(p.Latitude(), p.Longitude()); 00257 *output << s 00258 << Utility::str<real>(height + real(heightmult) * h, 4) 00259 << suff << eol; 00260 } else { 00261 if (gradp) { 00262 real gradn, grade; 00263 real h = g(p.Latitude(), p.Longitude(), gradn, grade); 00264 *output << Utility::str<real>(h, 4) << " " 00265 << Utility::str<real>(gradn * 1e6, 2) 00266 << (Math::isnan(gradn) ? " " : "e-6 ") 00267 << Utility::str<real>(grade * 1e6, 2) 00268 << (Math::isnan(grade) ? "" : "e-6") 00269 << eol; 00270 } else { 00271 real h = g(p.Latitude(), p.Longitude()); 00272 *output << Utility::str<real>(h, 4) << eol; 00273 } 00274 } 00275 } 00276 catch (const std::exception& e) { 00277 *output << "ERROR: " << e.what() << "\n"; 00278 retval = 1; 00279 } 00280 } 00281 } 00282 catch (const std::exception& e) { 00283 std::cerr << "Error reading " << geoid << ": " << e.what() << "\n"; 00284 retval = 1; 00285 } 00286 return retval; 00287 } 00288 catch (const std::exception& e) { 00289 std::cerr << "Caught exception: " << e.what() << "\n"; 00290 return 1; 00291 } 00292 catch (...) { 00293 std::cerr << "Caught unknown exception\n"; 00294 return 1; 00295 } 00296 }