git.itanic.dy.fi Git - maemo-mapper/blob - src/util.c

   1 /*
   2  * Copyright (C) 2006, 2007 John Costigan.
   3  *
   4  * POI and GPS-Info code originally written by Cezary Jackiewicz.
   5  *
   6  * Default map data provided by http://www.openstreetmap.org/
   7  *
   8  * This file is part of Maemo Mapper.
   9  *
  10  * Maemo Mapper is free software: you can redistribute it and/or modify
  11  * it under the terms of the GNU General Public License as published by
  12  * the Free Software Foundation, either version 3 of the License, or
  13  * (at your option) any later version.
  14  *
  15  * Maemo Mapper is distributed in the hope that it will be useful,
  16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  18  * GNU General Public License for more details.
  19  *
  20  * You should have received a copy of the GNU General Public License
  21  * along with Maemo Mapper.  If not, see <http://www.gnu.org/licenses/>.
  22  */
  23
  24 #ifdef HAVE_CONFIG_H
  25 #    include "config.h"
  26 #endif
  27
  28 #define _GNU_SOURCE
  29
  30 #include <ctype.h>
  31 #include <string.h>
  32 #include <math.h>
  33
  34 #ifndef LEGACY
  35 #    include <hildon/hildon-note.h>
  36 #else
  37 #    include <hildon-widgets/hildon-note.h>
  38 #endif
  39
  40 #include "types.h"
  41 #include "data.h"
  42 #include "defines.h"
  43
  44 #include "gpx.h"
  45 #include "util.h"
  46
  47 gboolean convert_iaru_loc_to_lat_lon(const gchar* txt_lon, gdouble* lat, gdouble* lon);
  48
  49 /**
  50  * Pop up a modal dialog box with simple error information in it.
  51  */
  52 void
  53 popup_error(GtkWidget *window, const gchar *error)
  54 {
  55     GtkWidget *dialog;
  56     printf("%s(\"%s\")\n", __PRETTY_FUNCTION__, error);
  57
  58     dialog = hildon_note_new_information(GTK_WINDOW(window), error);
  59
  60     gtk_dialog_run(GTK_DIALOG(dialog));
  61     gtk_widget_destroy(dialog);
  62
  63     vprintf("%s(): return\n", __PRETTY_FUNCTION__);
  64 }
  65
  66 void
  67 deg_format(gdouble coor, gchar *scoor, gchar neg_char, gchar pos_char)
  68 {
  69     gdouble min;
  70     gdouble acoor = fabs(coor);
  71     printf("%s()\n", __PRETTY_FUNCTION__);
  72
  73     switch(_degformat)
  74     {
  75             case IARU_LOC:
  76         case UK_OSGB:
  77         case UK_NGR:
  78         case UK_NGR6:
  79                 // These formats should not be formatted in the same way
  80                 // - they need to be converted first, therefore if we reach
  81                 // this bit of code use the first available format - drop through.
  82         case DDPDDDDD:
  83             sprintf(scoor, "%.5f°", coor);
  84             break;
  85
  86         case DDPDDDDD_NSEW:
  87             sprintf(scoor, "%.5f° %c", acoor,
  88                     coor < 0.0 ? neg_char : pos_char);
  89             break;
  90         case NSEW_DDPDDDDD:
  91             sprintf(scoor, "%c %.5f°",
  92                     coor < 0.0 ? neg_char : pos_char,
  93                     acoor);
  94             break;
  95
  96         case DD_MMPMMM:
  97             sprintf(scoor, "%d°%06.3f'",
  98                     (int)coor, (acoor - (int)acoor)*60.0);
  99             break;
 100         case DD_MMPMMM_NSEW:
 101             sprintf(scoor, "%d°%06.3f' %c",
 102                     (int)acoor, (acoor - (int)acoor)*60.0,
 103                     coor < 0.0 ? neg_char : pos_char);
 104             break;
 105         case NSEW_DD_MMPMMM:
 106             sprintf(scoor, "%c %d° %06.3f'",
 107                     coor < 0.0 ? neg_char : pos_char,
 108                     (int)acoor, (acoor - (int)acoor)*60.0);
 109             break;
 110
 111         case DD_MM_SSPS:
 112             min = (acoor - (int)acoor)*60.0;
 113             sprintf(scoor, "%d°%02d'%04.1f\"", (int)coor, (int)min,
 114                     ((min - (int)min)*60.0));
 115             break;
 116         case DD_MM_SSPS_NSEW:
 117             min = (acoor - (int)acoor)*60.0;
 118             sprintf(scoor, "%d°%02d'%04.1f\" %c", (int)acoor, (int)min,
 119                     ((min - (int)min)*60.0),
 120                     coor < 0.0 ? neg_char : pos_char);
 121             break;
 122         case NSEW_DD_MM_SSPS:
 123             min = (acoor - (int)acoor)*60.0;
 124             sprintf(scoor, "%c %d° %02d' %04.1f\"",
 125                     coor < 0.0 ? neg_char : pos_char,
 126                     (int)acoor, (int)min,
 127                     ((min - (int)min)*60.0));
 128             break;
 129     }
 130     vprintf("%s(): return\n", __PRETTY_FUNCTION__);
 131 }
 132
 133 /** Return the location (in units) of the given string address.  This function
 134   * makes a call to the internet, so it may take some time. */
 135 Point locate_address(GtkWidget *parent, const gchar *addr)
 136 {
 137     Path temp;
 138     Point retval = _point_null;
 139     GnomeVFSResult vfs_result;
 140     gint size;
 141     gchar *bytes = NULL;
 142     gchar *addr_escaped;
 143     gchar *buffer;
 144     printf("%s(%s)\n", __PRETTY_FUNCTION__, addr);
 145
 146     addr_escaped = gnome_vfs_escape_string(addr);
 147     buffer = g_strdup_printf(_route_dl_url, addr_escaped, addr_escaped);
 148     g_free(addr_escaped);
 149
 150     /* Attempt to download the route from the server. */
 151     if(GNOME_VFS_OK != (vfs_result = gnome_vfs_read_entire_file(
 152                 buffer, &size, &bytes)))
 153     {
 154         g_free(buffer);
 155         g_free(bytes);
 156         popup_error(parent, _("Failed to connect to GPX Directions server"));
 157         return _point_null;
 158     }
 159
 160     g_free(buffer);
 161
 162     MACRO_PATH_INIT(temp);
 163
 164     if(strncmp(bytes, "<?xml", strlen("<?xml")))
 165     {
 166         /* Not an XML document - must be bad locations. */
 167         popup_error(parent, _("Invalid address."));
 168     }
 169     /* Else, if GPS is enabled, replace the route, otherwise append it. */
 170     else if(!gpx_path_parse(&temp, bytes, size, 0) || !temp.head[1].unity)
 171     {
 172         popup_error(parent, _("Unknown error while locating address."));
 173     }
 174     else
 175     {
 176         /* Save Destination in Route Locations list. */
 177         GtkTreeIter iter;
 178         if(!g_slist_find_custom(_loc_list, addr, (GCompareFunc)strcmp))
 179         {
 180             _loc_list = g_slist_prepend(_loc_list, g_strdup(addr));
 181             gtk_list_store_insert_with_values(_loc_model, &iter,
 182                     INT_MAX, 0, addr, -1);
 183         }
 184
 185         retval = temp.head[1];
 186     }
 187
 188     MACRO_PATH_FREE(temp);
 189     g_free(bytes);
 190
 191     vprintf("%s(): return (%d, %d)\n", __PRETTY_FUNCTION__,
 192             retval.unitx, retval.unity);
 193     return retval;
 194 }
 195
 196 /**
 197  * Calculate the distance between two lat/lon pairs.  The distance is returned
 198  * in nautical miles and should be converted using UNITS_CONVERT[_units].
 199  */
 200 gdouble
 201 calculate_distance(gdouble lat1, gdouble lon1, gdouble lat2, gdouble lon2)
 202 {
 203     gdouble dlat, dlon, slat, slon, a;
 204
 205     /* Convert to radians. */
 206     lat1 = deg2rad(lat1);
 207     lon1 = deg2rad(lon1);
 208     lat2 = deg2rad(lat2);
 209     lon2 = deg2rad(lon2);
 210
 211     dlat = lat2 - lat1;
 212     dlon = lon2 - lon1;
 213
 214     slat = sin(dlat / 2.0);
 215     slon = sin(dlon / 2.0);
 216     a = (slat * slat) + (cos(lat1) * cos(lat2) * slon * slon);
 217
 218     vprintf("%s(): return\n", __PRETTY_FUNCTION__);
 219     return ((2.0 * atan2(sqrt(a), sqrt(1.0 - a))) * EARTH_RADIUS);
 220 }
 221
 222 /**
 223  * Calculate the bearing between two lat/lon pairs.  The bearing is returned
 224  * as the angle from lat1/lon1 to lat2/lon2.
 225  */
 226 gdouble
 227 calculate_bearing(gdouble lat1, gdouble lon1, gdouble lat2, gdouble lon2)
 228 {
 229     gdouble x, y;
 230     gdouble dlon = deg2rad(lon2 - lon1);
 231     lat1 = deg2rad(lat1);
 232     lat2 = deg2rad(lat2);
 233
 234     y = sin(dlon) * cos(lat2);
 235     x = (cos(lat1) * sin(lat2)) - (sin(lat1) * cos(lat2) * cos(dlon));
 236
 237     dlon = rad2deg(atan2(y, x));
 238     if(dlon < 0.0)
 239         dlon += 360.0;
 240     return dlon;
 241 }
 242
 243
 244
 245 void
 246 force_min_visible_bars(HildonControlbar *control_bar, gint num_bars)
 247 {
 248 #ifdef LEGACY
 249     GValue val;
 250     printf("%s()\n", __PRETTY_FUNCTION__);
 251     memset(&val, 0, sizeof(val));
 252     g_value_init(&val, G_TYPE_INT);
 253     g_value_set_int(&val, num_bars);
 254     g_object_set_property(G_OBJECT(control_bar), "minimum-visible-bars", &val);
 255     vprintf("%s(): return\n", __PRETTY_FUNCTION__);
 256 #endif
 257 }
 258
 259 gboolean
 260 banner_reset()
 261 {
 262     printf("%s()\n", __PRETTY_FUNCTION__);
 263
 264     /* Re-enable any banners that might have been up. */
 265     {
 266         if(_connect_banner)
 267         {
 268             gtk_widget_hide(_connect_banner);
 269             gtk_widget_unrealize(_connect_banner);
 270             gtk_widget_realize(_connect_banner);
 271             gtk_widget_show(_connect_banner);
 272         }
 273
 274         if(_fix_banner)
 275         {
 276             gtk_widget_hide(_fix_banner);
 277             gtk_widget_unrealize(_fix_banner);
 278             gtk_widget_realize(_fix_banner);
 279             gtk_widget_show(_fix_banner);
 280         }
 281
 282         if(_waypoint_banner)
 283         {
 284             gtk_widget_hide(_waypoint_banner);
 285             gtk_widget_unrealize(_waypoint_banner);
 286             gtk_widget_realize(_waypoint_banner);
 287             gtk_widget_show(_waypoint_banner);
 288         }
 289
 290         if(_download_banner)
 291         {
 292             gtk_widget_hide(_download_banner);
 293             gtk_widget_unrealize(_download_banner);
 294             gtk_widget_realize(_download_banner);
 295             gtk_widget_show(_download_banner);
 296         }
 297
 298         /*
 299         ConnState old_state = _gps_state;
 300         set_conn_state(RCVR_OFF);
 301         set_conn_state(old_state);
 302         */
 303     }
 304
 305     vprintf("%s(): return FALSE\n", __PRETTY_FUNCTION__);
 306     return FALSE;
 307 }
 308 /*
 309  * Get one numeric token off the string, fractional part separator
 310  * is ',' or '.' and number may follow any token from sep.
 311  * Return value found is in d.
 312  *
 313  * If utf8_deg is set then accept also Unicode degree symbol U+00B0
 314  * encoded as UTF-8 octets 0xc2 0xb0.
 315  *
 316  * @returns
 317  *    0 : on error.
 318  *    1 : when nothing or just white space was seen.
 319  *    2 : when fractional number was seen.
 320  *    3 : when whole number was seen.
 321  *  In case 0, endptr points to the beginning of the input string nptr,
 322  *  d is undefined.
 323  *
 324  *  In case 1 endptr points past the whitespace, d is undefined.
 325  *
 326  *  In cases 2 and 3 the found number is stored in d and endptr points
 327  *  to first char not part of the number.
 328  *
 329  */
 330 static gint
 331 strdmstod_1(gdouble *d, gchar *nptr, gchar **endptr, gchar *sep, gint utf8_deg)
 332 {
 333     guchar *p;
 334     gint v_flag, f_flag;
 335     gdouble v;
 336
 337     p = nptr;
 338
 339     while (*p && isspace(*p)) {
 340         p++;
 341     }
 342
 343     v_flag = 0;
 344     f_flag = 0;
 345
 346     /* whole part */
 347     v = 0.0;
 348     while (*p && isdigit(*p)) {
 349         v *= 10;
 350         v += *p++ - '0';
 351         v_flag = 1;
 352     }
 353
 354     if (v_flag) {
 355         if (*p && (*p == '.' || *p == ',')) {
 356             gdouble f;
 357             gint n;
 358
 359             p++;
 360
 361             /* fractional part */
 362             f = 0.0;
 363             n = 1;
 364             while (*p && isdigit(*p)) {
 365                 f *= 10.0;
 366                 f += *p++ - '0';
 367                 n *= 10;
 368             }
 369
 370             if (n > 1) {
 371                 f_flag = 1;
 372                 v += f/n;
 373             }
 374         }
 375
 376         /* allow for extra sep char at the end */
 377         if (*p && strchr(sep, *p)) {
 378             p++;
 379         } else if (utf8_deg && *p == 0xc2 && *(p+1) == 0xb0) {
 380             p += 2;
 381         }
 382
 383         *d = v;
 384         if (endptr) *endptr = p;
 385         return f_flag ? 2 : 3;
 386     }
 387
 388     if (endptr) *endptr = p;
 389     return *p == 0;
 390 }
 391
 392 static gdouble
 393 strdmstod_2(gchar *nptr, gchar **endptr)
 394 {
 395     gint ret;
 396     guchar *p;
 397     gdouble d, m, s;
 398
 399     p = nptr;
 400
 401     /* degrees */
 402     ret = strdmstod_1(&d, p, endptr, "dD@", 1);
 403     switch (ret) {
 404         case 0: return 0.0;
 405         case 1:
 406                 if (endptr) *endptr = (char *)nptr;
 407                 return 0.0;
 408         case 2: return d;
 409     }
 410
 411     /* minutes */
 412     p = *endptr;
 413     m = 0.0;
 414     ret = strdmstod_1(&m, p, endptr, "mM'", 0);
 415     switch (ret) {
 416         case 0: return 0.0;
 417         case 1: return d;
 418         case 2: return d + m/60.0;
 419     }
 420
 421     /* seconds */
 422     p = *endptr;
 423     ret = strdmstod_1(&s, p, endptr, "sS\"", 0);
 424     switch (ret) {
 425         case 0: return 0.0;
 426         case 1: return d + m/60.0;
 427         case 2:
 428         case 3: return d + m/60.0 + s/3600.0;
 429     }
 430
 431     /* can't be here */
 432     return 0.0;
 433 }
 434
 435 /*
 436  * format: / \s* [+-NSWE]? \s* ( d | D \s+ m | D \s+ M \s+ s ) [NSWE]? /x
 437  * where D := / \d+[@d°]? /ix
 438  *       M := / \d+['m]? /ix
 439  *       d := / D | \d+[,.]\d+[@d°]? /ix
 440  *       m := / M | \d+[,.]\d+['m]? /ix
 441  *       s := / S | \d+[,.]\d+["s]? /ix
 442  *
 443  *   and N or W are treated as positive, S or E are treated as negative,
 444  *   they may occur only once.
 445  */
 446 gdouble
 447 strdmstod(const gchar *nptr, gchar **endptr)
 448 {
 449     gint s;
 450     gchar *p, *end;
 451     gchar *sign = 0;
 452     gdouble d;
 453
 454     p = (char *)nptr;
 455
 456     while(*p && isspace(*p))
 457         p++;
 458
 459     if(!*p) {
 460         if(endptr)
 461             *endptr = (char *)nptr;
 462         return 0.0;
 463     }
 464
 465     if(strchr("nwseNWSE-+", *p)) {
 466         sign = p;
 467         p++;
 468     }
 469
 470     d = strdmstod_2(p, &end);
 471     if(p == end && d == 0.0) {
 472         if(endptr) *endptr = end;
 473         return d;
 474     }
 475
 476     p = end;
 477     while(*p && isspace(*p))
 478         p++;
 479
 480     s = 1;
 481     if(sign == 0) {
 482         if(*p && strchr("nwseNWSE", *p)) {
 483             if(tolower(*p) == 's' || tolower(*p) == 'w')
 484                 s = -1;
 485             p++;
 486         }
 487     } else {
 488         if(tolower(*sign) == 's' || tolower(*sign) == 'w' || *sign == '-')
 489             s = -1;
 490         printf("s: %d\n", s);
 491     }
 492
 493     if(endptr) *endptr = p;
 494     return s * d;
 495 }
 496
 497 double marc(double bf0, double n, double phi0, double phi)
 498 {
 499         return bf0 * (((1 + n + ((5 / 4) * (n * n)) + ((5 / 4) * (n * n * n))) * (phi - phi0))
 500             - (((3 * n) + (3 * (n * n)) + ((21 / 8) * (n * n * n))) * (sin(phi - phi0)) * (cos(phi + phi0)))
 501             + ((((15 / 8) * (n * n)) + ((15 / 8) * (n * n * n))) * (sin(2 * (phi - phi0))) * (cos(2 * (phi + phi0))))
 502             - (((35 / 24) * (n * n * n)) * (sin(3 * (phi - phi0))) * (cos(3 * (phi + phi0)))));
 503 }
 504
 505 gboolean os_grid_check_lat_lon(double lat, double lon)
 506 {
 507         // TODO - Check exact OS Grid range
 508         if(lat < 50.0 || lat > 62 || lon < -7.5 || lon > 2.2 )
 509         {
 510                 return FALSE;
 511         }
 512         else
 513         {
 514                 return TRUE;
 515         }
 516 }
 517
 518 gboolean coord_system_check_lat_lon (gdouble lat, gdouble lon, gint *fallback_deg_format)
 519 {
 520         // Is the current coordinate system applicable to the provided lat and lon?
 521         gboolean valid = FALSE;
 522
 523         switch(_degformat)
 524         {
 525         case UK_OSGB:
 526         case UK_NGR:
 527         case UK_NGR6:
 528                 valid = os_grid_check_lat_lon(lat, lon);
 529                 if(fallback_deg_format != NULL) *fallback_deg_format = DDPDDDDD;
 530                 break;
 531         default:
 532                 valid = TRUE;
 533                 break;
 534         }
 535
 536         return valid;
 537 }
 538
 539 gboolean convert_lat_long_to_os_grid(double lat, double lon, int *easting, int *northing)
 540 {
 541         if(!os_grid_check_lat_lon(lat, lon))
 542         {
 543                 return FALSE;
 544         }
 545
 546         const double deg2rad_multi = (2 * PI / 360);
 547
 548         const double phi = lat * deg2rad_multi;          // convert latitude to radians
 549         const double lam = lon * deg2rad_multi;          // convert longitude to radians
 550         const double a = 6377563.396;              // OSGB semi-major axis
 551         const double b = 6356256.91;               // OSGB semi-minor axis
 552         const double e0 = 400000;                  // easting of false origin
 553         const double n0 = -100000;                 // northing of false origin
 554         const double f0 = 0.9996012717;            // OSGB scale factor on central meridian
 555         const double e2 = 0.0066705397616;         // OSGB eccentricity squared
 556         const double lam0 = -0.034906585039886591; // OSGB false east
 557         const double phi0 = 0.85521133347722145;   // OSGB false north
 558         const double af0 = a * f0;
 559         const double bf0 = b * f0;
 560
 561         // easting
 562         double slat2 = sin(phi) * sin(phi);
 563         double nu = af0 / (sqrt(1 - (e2 * (slat2))));
 564         double rho = (nu * (1 - e2)) / (1 - (e2 * slat2));
 565         double eta2 = (nu / rho) - 1;
 566         double p = lam - lam0;
 567         double IV = nu * cos(phi);
 568         double clat3 = pow(cos(phi), 3);
 569         double tlat2 = tan(phi) * tan(phi);
 570         double V = (nu / 6) * clat3 * ((nu / rho) - tlat2);
 571         double clat5 = pow(cos(phi), 5);
 572         double tlat4 = pow(tan(phi), 4);
 573         double VI = (nu / 120) * clat5 * ((5 - (18 * tlat2)) + tlat4 + (14 * eta2) - (58 * tlat2 * eta2));
 574         double east = e0 + (p * IV) + (pow(p, 3) * V) + (pow(p, 5) * VI);
 575
 576         // northing
 577         double n = (af0 - bf0) / (af0 + bf0);
 578         double M = marc(bf0, n, phi0, phi);
 579         double I = M + (n0);
 580         double II = (nu / 2) * sin(phi) * cos(phi);
 581         double III = ((nu / 24) * sin(phi) * pow(cos(phi), 3)) * (5 - pow(tan(phi), 2) + (9 * eta2));
 582         double IIIA = ((nu / 720) * sin(phi) * clat5) * (61 - (58 * tlat2) + tlat4);
 583         double north = I + ((p * p) * II) + (pow(p, 4) * III) + (pow(p, 6) * IIIA);
 584
 585          // make whole number values
 586          *easting = round(east);   // round to whole number
 587          *northing = round(north); // round to whole number
 588
 589          return TRUE;
 590 }
 591
 592 gboolean convert_os_grid_to_bng(gint easting, gint northing, gchar* bng)
 593 {
 594         gdouble eX = (gdouble)easting / 500000.0;
 595         gdouble nX = (gdouble)northing / 500000.0;
 596         gdouble tmp = floor(eX) - 5.0 * floor(nX) + 17.0;
 597         gchar eing[12];
 598         gchar ning[12];
 599
 600         nX = 5.0 * (nX - floor(nX));
 601         eX = 20.0 - 5.0 * floor(nX) + floor(5.0 * (eX - floor(eX)));
 602         if (eX > 7.5) eX = eX + 1;    // I is not used
 603         if (tmp > 7.5) tmp = tmp + 1; // I is not used
 604
 605         snprintf(eing, 12, "%u", easting);
 606         snprintf(ning, 12, "%u", northing);
 607
 608         snprintf(eing, (_degformat == UK_NGR ? 5 : 4), "%s", eing+1);
 609         snprintf(ning, (_degformat == UK_NGR ? 5 : 4), "%s", ning+1);
 610
 611         sprintf(bng, "%c%c%s%s",
 612                         (char)(tmp + 65),
 613                         (char)(eX + 65),
 614                         eing, ning
 615                 );
 616
 617
 618         return TRUE;
 619 }
 620
 621 gint convert_str_to_int(const gchar *str)
 622 {
 623         gint i=0;
 624         gint result = 0;
 625         while(str[i] >= '0' && str[i] <= '9')
 626         {
 627                 gint v = str[i] - '0';
 628
 629                 result = (result * 10) + v;
 630
 631                 i++;
 632         }
 633
 634         return result;
 635 }
 636 gboolean convert_os_xy_to_latlon(const gchar *easting, const gchar *northing, gdouble *d_lat, gdouble *d_lon)
 637 {
 638     const double deg2rad_multi = (2 * PI / 360);
 639
 640         const gdouble N = (gdouble)convert_str_to_int(northing);
 641         const gdouble E = (gdouble)convert_str_to_int(easting);
 642
 643         const gdouble a = 6377563.396, b = 6356256.910;              // Airy 1830 major & minor semi-axes
 644         const gdouble F0 = 0.9996012717;                             // NatGrid scale factor on central meridian
 645         const gdouble lat0 = 49*PI/180, lon0 = -2*PI/180;  // NatGrid true origin
 646         const gdouble N0 = -100000, E0 = 400000;                     // northing & easting of true origin, metres
 647         const gdouble e2 = 1 - (b*b)/(a*a);                          // eccentricity squared
 648         const gdouble n = (a-b)/(a+b), n2 = n*n, n3 = n*n*n;
 649
 650         gdouble lat=lat0, M=0;
 651         do {
 652                 lat = (N-N0-M)/(a*F0) + lat;
 653
 654                 const gdouble Ma = (1 + n + (5/4)*n2 + (5/4)*n3) * (lat-lat0);
 655                 const gdouble Mb = (3*n + 3*n*n + (21/8)*n3) * sin(lat-lat0) * cos(lat+lat0);
 656                 const gdouble Mc = ((15/8)*n2 + (15/8)*n3) * sin(2*(lat-lat0)) * cos(2*(lat+lat0));
 657                 const gdouble Md = (35/24)*n3 * sin(3*(lat-lat0)) * cos(3*(lat+lat0));
 658             M = b * F0 * (Ma - Mb + Mc - Md);                // meridional arc
 659
 660         } while (N-N0-M >= 0.00001);  // ie until < 0.01mm
 661
 662         const gdouble cosLat = cos(lat), sinLat = sin(lat);
 663         const gdouble nu = a*F0/sqrt(1-e2*sinLat*sinLat);              // transverse radius of curvature
 664         const gdouble rho = a*F0*(1-e2)/pow(1-e2*sinLat*sinLat, 1.5);  // meridional radius of curvature
 665         const gdouble eta2 = nu/rho-1;
 666
 667         const gdouble tanLat = tan(lat);
 668         const gdouble tan2lat = tanLat*tanLat, tan4lat = tan2lat*tan2lat, tan6lat = tan4lat*tan2lat;
 669         const gdouble secLat = 1/cosLat;
 670         const gdouble nu3 = nu*nu*nu, nu5 = nu3*nu*nu, nu7 = nu5*nu*nu;
 671         const gdouble VII = tanLat/(2*rho*nu);
 672         const gdouble VIII = tanLat/(24*rho*nu3)*(5+3*tan2lat+eta2-9*tan2lat*eta2);
 673         const gdouble IX = tanLat/(720*rho*nu5)*(61+90*tan2lat+45*tan4lat);
 674         const gdouble X = secLat/nu;
 675         const gdouble XI = secLat/(6*nu3)*(nu/rho+2*tan2lat);
 676         const gdouble XII = secLat/(120*nu5)*(5+28*tan2lat+24*tan4lat);
 677         const gdouble XIIA = secLat/(5040*nu7)*(61+662*tan2lat+1320*tan4lat+720*tan6lat);
 678
 679         const gdouble dE = (E-E0), dE2 = dE*dE, dE3 = dE2*dE, dE4 = dE2*dE2, dE5 = dE3*dE2;
 680         const gdouble dE6 = dE4*dE2, dE7 = dE5*dE2;
 681         lat = lat - VII*dE2 + VIII*dE4 - IX*dE6;
 682         const gdouble lon = lon0 + X*dE - XI*dE3 + XII*dE5 - XIIA*dE7;
 683
 684         *d_lon = lon / deg2rad_multi;
 685         *d_lat = lat / deg2rad_multi;
 686
 687         return TRUE;
 688 }
 689
 690 gboolean convert_os_ngr_to_latlon(const gchar *text, gdouble *d_lat, gdouble *d_lon)
 691 {
 692         // get numeric values of letter references, mapping A->0, B->1, C->2, etc:
 693         gint l1;
 694         gint l2;
 695
 696         gchar easting[7], northing[7];
 697
 698         if( ((gchar)text[0])>='a' && ((gchar)text[0]) <= 'z' )
 699                 l1 = text[0] - (gint)'a'; // lower case
 700         else if( ((gchar)text[0])>='A' && ((gchar)text[0]) <= 'Z' )
 701                 l1 = text[0] - (gint)'A'; // upper case
 702         else
 703                 return FALSE; // Not a letter - invalid grid ref
 704
 705         if( ((gchar)text[1])>='a' && ((gchar)text[1]) <= 'z' )
 706                 l2 = text[1] - (gint)'a'; // lower case
 707         else if( ((gchar)text[1])>='A' && ((gchar)text[1]) <= 'Z' )
 708                 l2 = text[1] - (gint)'A'; // upper case
 709         else
 710                 return FALSE; // Not a letter - invalid grid ref
 711
 712
 713         // shuffle down letters after 'I' since 'I' is not used in grid:
 714         if (l1 > 7) l1--;
 715         if (l2 > 7) l2--;
 716
 717         // convert grid letters into 100km-square indexes from false origin (grid square SV):
 718         gdouble e = ((l1-2)%5)*5 + (l2%5);
 719         gdouble n = (19-floor(l1/5)*5) - floor(l2/5);
 720
 721         // skip grid letters to get numeric part of ref, stripping any spaces:
 722         gchar *gridref = (gchar*)(text+2);
 723
 724         // user may have entered a space, so remove any spaces
 725         while(gridref[0] == ' ')
 726         {
 727                 gridref = (gchar*)(text+1);
 728         }
 729
 730
 731         // floor the length incase a space has been added
 732         const gint len = (gint)floor((gdouble)strlen(gridref)/2.00); // normally this will be 4, often 3
 733         if(len>5 || len <3) return FALSE;
 734
 735         snprintf(easting, 2+len, "%u%s", (gint)e, gridref);
 736         snprintf(northing, 2+len, "%u%s", (gint)n, gridref+len);
 737
 738         fprintf(stderr, "easting = %s, northing = %s\n", easting, northing);
 739
 740         switch (len)
 741         {
 742         case 3:
 743                 easting[4] = '0';
 744                 easting[5] = '0';
 745                 northing[4] = '0';
 746                 northing[5] = '0';
 747                 break;
 748         case 4:
 749                 easting[5] = '0';
 750                 northing[5] = '0';
 751                 break;
 752         // 10-digit refs are already 1m
 753         }
 754
 755         easting[6] = '\0';
 756         northing[6] = '\0';
 757
 758
 759         convert_os_xy_to_latlon(easting, northing, d_lat, d_lon);
 760
 761         return TRUE;
 762 }
 763
 764
 765 // Attempt to convert any user entered grid reference to a double lat/lon
 766 // return TRUE on valid
 767 gboolean parse_coords(const gchar* txt_lat, const gchar* txt_lon, gdouble* lat, gdouble* lon)
 768 {
 769         gboolean valid = FALSE;
 770
 771         // UK_NGR starts with two letters, and then all numbers - it may contain spaces - no lon will be entered
 772         if( _degformat == UK_NGR || _degformat == UK_NGR6 )
 773         {
 774                 valid = convert_os_ngr_to_latlon(txt_lat, lat, lon);
 775
 776         if(!valid || *lat < -90. || *lat > 90.) { valid = FALSE; }
 777         else   if(*lon < -180. || *lon > 180.)  { valid = FALSE; }
 778
 779         }
 780         // UK_OSGB contains two 6 digit integers
 781         else if( _degformat == UK_OSGB)
 782         {
 783                 valid = convert_os_xy_to_latlon(txt_lat, txt_lon, lat, lon);
 784
 785         if(!valid || *lat < -90. || *lat > 90.) { valid = FALSE; }
 786         else   if(*lon < -180. || *lon > 180.)  { valid = FALSE; }
 787
 788         }
 789         else if( _degformat == IARU_LOC)
 790         {
 791                 valid = convert_iaru_loc_to_lat_lon(txt_lat, lat, lon);
 792
 793         if(!valid || *lat < -90. || *lat > 90.) { valid = FALSE; }
 794         else   if(*lon < -180. || *lon > 180.)  { valid = FALSE; }
 795         }
 796         // It must either be invalid, or a lat/lon format
 797         else
 798         {
 799                 gchar* error_check;
 800                 *lat = strdmstod(txt_lat, &error_check);
 801
 802         if(txt_lat == error_check || *lat < -90. || *lat > 90.) { valid = FALSE; }
 803         else { valid = TRUE; }
 804
 805         if(valid == TRUE)
 806         {
 807                 *lon = strdmstod(txt_lon, &error_check);
 808
 809             if(txt_lon == error_check || *lon < -180. || *lon > 180.)  { valid = FALSE; }
 810         }
 811         }
 812
 813
 814
 815         return valid;
 816 }
 817
 818 gboolean convert_iaru_loc_to_lat_lon(const gchar* txt_lon, gdouble* lat, gdouble* lon)
 819 {
 820         gint u_first = 0;
 821         gint u_second = 0;
 822         gint u_third = 0;
 823         gint u_fourth = 0;
 824         gint u_fifth = 0;
 825         gint u_sixth = 0;
 826         gint u_seventh = 0;
 827         gint u_eighth = 0;
 828
 829         if(strlen(txt_lon) >= 1)
 830         {
 831                 if( ((gchar)txt_lon[0])>='a' && ((gchar)txt_lon[0]) <= 'z' )
 832                         u_first = txt_lon[0] - (gint)'a'; // lower case
 833                 else if( ((gchar)txt_lon[0])>='A' && ((gchar)txt_lon[0]) <= 'Z' )
 834                         u_first = txt_lon[0] - (gint)'A'; // upper case
 835         }
 836
 837         if(strlen(txt_lon) >= 2)
 838         {
 839                 if( ((gchar)txt_lon[1])>='a' && ((gchar)txt_lon[1]) <= 'z' )
 840                         u_second = txt_lon[1] - (gint)'a'; // lower case
 841                 else if( ((gchar)txt_lon[1])>='A' && ((gchar)txt_lon[1]) <= 'Z' )
 842                         u_second = txt_lon[1] - (gint)'A'; // upper case
 843         }
 844
 845         if(strlen(txt_lon) >= 3)
 846                 u_third = txt_lon[2] - (gint)'0';
 847
 848         if(strlen(txt_lon) >= 4)
 849                 u_fourth = txt_lon[3] - (gint)'0';
 850
 851         if(strlen(txt_lon) >= 5)
 852         {
 853                 if( ((gchar)txt_lon[4])>='a' && ((gchar)txt_lon[4]) <= 'z' )
 854                         u_fifth = txt_lon[4] - (gint)'a'; // lower case
 855                 else if( ((gchar)txt_lon[4])>='A' && ((gchar)txt_lon[4]) <= 'Z' )
 856                         u_fifth = txt_lon[4] - (gint)'A'; // upper case
 857         }
 858
 859         if(strlen(txt_lon) >= 6)
 860         {
 861                 if( ((gchar)txt_lon[5])>='a' && ((gchar)txt_lon[5]) <= 'z' )
 862                         u_sixth = txt_lon[5] - (gint)'a'; // lower case
 863                 else if( ((gchar)txt_lon[5])>='A' && ((gchar)txt_lon[5]) <= 'Z' )
 864                         u_sixth = txt_lon[5] - (gint)'A'; // upper case
 865         }
 866
 867
 868         if(strlen(txt_lon) >= 7)
 869                 u_seventh= txt_lon[6] - (gint)'0';
 870
 871         if(strlen(txt_lon) >= 8)
 872                 u_eighth = txt_lon[7] - (gint)'0';
 873
 874         *lat = ((gdouble)u_first * 20.0) + ((gdouble)u_third * 2.0) + ((gdouble)u_fifth * (2.0/24.0)) + ((gdouble)u_seventh * (2.0/240.0)) - 90.0;
 875         *lon = ((gdouble)u_second * 10.0) + ((gdouble)u_fourth) + ((gdouble)u_sixth * (1.0/24.0)) + ((gdouble)u_eighth * (1.0/240.0)) - 180.0;
 876
 877         return TRUE;
 878 }
 879
 880 void convert_lat_lon_to_iaru_loc(gdouble d_lat, gdouble d_lon, gchar *loc)
 881 {
 882         const gdouble d_a_lat = (d_lat+90.0);
 883         const gdouble d_a_lon = (d_lon+180.0);
 884
 885         const gint i_first  = (gint)floor(d_a_lon/20.0);
 886         const gint i_second = (gint)floor(d_a_lat/10.0);
 887         const gint i_third  = (gint)floor((d_a_lon - (20.0*i_first))/2.0);
 888         const gint i_fourth = (gint)floor((d_a_lat - (10.0*i_second)));
 889         const gint i_fifth  = (gint)floor((d_a_lon - (20.0*i_first) - (2.0*i_third))/(5.0/60.0));
 890         const gint i_sixth  = (gint)floor((d_a_lat - (10.0*i_second) - (i_fourth))/(2.5/60.0));
 891
 892         const gint i_seventh  = (gint)floor((d_a_lon - (20.0*i_first) - (2.0*i_third)
 893                         - ((2.0/24.0)*i_fifth))/(2.0/240.0)  );
 894         const gint i_eighth = (gint)floor((d_a_lat - (10.0*i_second) - (i_fourth)
 895                         - ((1.0/24.0)*i_sixth))/(1.0/240.0));
 896
 897
 898         sprintf(loc, "%c%c%u%u%c%c%u%u",
 899                         'A'+i_first,
 900                         'A'+i_second,
 901                         i_third,
 902                         i_fourth,
 903                         'a' + i_fifth,
 904                         'a' + i_sixth,
 905                         i_seventh,
 906                         i_eighth);
 907 }
 908
 909 gboolean convert_lat_lon_to_bng(gdouble lat, gdouble lon, gchar* bng)
 910 {
 911         gint easting, northing;
 912
 913         if( convert_lat_long_to_os_grid(lat, lon, &easting, &northing) )
 914         {
 915                 if( convert_os_grid_to_bng(easting, northing, bng) )
 916                 {
 917                         return TRUE;
 918                 }
 919         }
 920
 921         return FALSE;
 922 }
 923
 924
 925 void format_lat_lon(gdouble d_lat, gdouble d_lon, gchar* lat, gchar* lon)
 926 {
 927         gint east = 0;
 928         gint north = 0;
 929
 930         switch (_degformat)
 931         {
 932         case UK_OSGB:
 933
 934                 if(convert_lat_long_to_os_grid(d_lat, d_lon, &east, &north))
 935                 {
 936                         sprintf(lat, "%06d", east);
 937                         sprintf(lon, "%06d", north);
 938                 }
 939                 else
 940                 {
 941                         // Failed (possibly out of range), so use defaults
 942                         lat_format(d_lat, lat);
 943                         lon_format(d_lon, lon);
 944                 }
 945                 break;
 946         case UK_NGR:
 947         case UK_NGR6:
 948
 949                 if(convert_lat_lon_to_bng(d_lat, d_lon, lat))
 950                 {
 951                         lon[0] = 0;
 952                 }
 953                 else
 954                 {
 955                         // Failed (possibly out of range), so use defaults
 956                         lat_format(d_lat, lat);
 957                         lat_format(d_lon, lon);
 958                 }
 959                 break;
 960
 961         case IARU_LOC:
 962                 convert_lat_lon_to_iaru_loc(d_lat, d_lon, lat);
 963
 964                 break;
 965
 966         default:
 967                 lat_format(d_lat, lat);
 968                 lon_format(d_lon, lon);
 969
 970                 break;
 971         }
 972 }
 973
 974 /* Custom version of g_ascii_strtoll, since Gregale does not support
 975  * GLIB 2.12. */
 976 gint64
 977 g_ascii_strtoll(const gchar *nptr, gchar **endptr, guint base)
 978 {
 979     gchar *minus = g_strstr_len(nptr, -1, "");
 980     if(minus)
 981         return -g_ascii_strtoull(minus + 1, endptr, base);
 982     else
 983         return g_ascii_strtoull(nptr, endptr, base);
 984 }
 985
 986 #if 0
 987 struct t_case {
 988     gchar *fmt;
 989     gdouble value;
 990 } t_cases[] = {
 991     { "12°", 12 },
 992     { "+12d", 12 },
 993     { "-12.345d", -12.345 },
 994     { "12.345 E", -12.345 },
 995     { "12d34m", 12.5666667 },
 996     { "N 12 34", 12.5666667 },
 997     { "S 12d34.56m", -12.576 },
 998     { "W 12 34.56", 12.576 },
 999     { "E 12d34m56s", -12.582222 },
1000     { "12 34 56 S", -12.582222 },
1001     { "12 34 56", 12.582222 },
1002     { "12d34m56.7s E", -12.582417 },
1003     { "12 34 56.7 W", 12.582417 },
1004     { "12° 34 56.7 W", 12.582417 },
1005
1006     { 0, 0 }
1007 };
1008
1009 int
1010 strdmstod_test()
1011 {
1012     struct t_case *t;
1013     gint fail, ok;
1014
1015     fail = ok = 0;
1016
1017     for (t = t_cases; t->fmt; t++) {
1018         gdouble v;
1019         gchar *endp;
1020
1021         v = strdmstod(t->fmt, &endp);
1022         if (endp == t->fmt || *endp != 0) {
1023             fprintf(stderr, "FAIL syntax %s\n", t->fmt);
1024             fail++;
1025         } else if (fabs(v - t->value) > 0.000001) {
1026             fprintf(stderr, "FAIL value %s -> %.10g (%.10g)\n",
1027                     t->fmt, v, t->value);
1028             fail++;
1029         } else {
1030             ok++;
1031         }
1032     }
1033
1034     if (fail == 0) {
1035         fprintf(stderr, "ALL TESTS OK\n");
1036     } else {
1037         fprintf(stderr, "FAIL %d, OK %d\n", fail, ok);
1038     }
1039 }
1040 #endif