/** * \addtogroup utilities * @{ */ /** * \file src/locator.c * \brief locator and bearing conversion interface * \author Stephane Fillod and the Hamlib Group * \date 2000-2010 * * Hamlib Interface - locator, bearing, and conversion calls */ /* * Hamlib Interface - locator and bearing conversion calls * Copyright (c) 2001-2010 by Stephane Fillod * Copyright (c) 2003 by Nate Bargmann * Copyright (c) 2003 by Dave Hines * * * Code to determine bearing and range was taken from the Great Circle, * by S. R. Sampson, N5OWK. * Ref: "Air Navigation", Air Force Manual 51-40, 1 February 1987 * Ref: "ARRL Satellite Experimenters Handbook", August 1990 * * Code to calculate distance and azimuth between two Maidenhead locators, * taken from wwl, by IK0ZSN Mirko Caserta. * * New bearing code added by N0NB was found at: * http://williams.best.vwh.net/avform.htm#Crs * * * This library is free software; you can redistribute it and/or * modify it under the terms of the GNU Lesser General Public * License as published by the Free Software Foundation; either * version 2.1 of the License, or (at your option) any later version. * * This library is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * Lesser General Public License for more details. * * You should have received a copy of the GNU Lesser General Public * License along with this library; if not, write to the Free Software * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA * */ /*! \page hamlib Hamlib general purpose API * * Here are grouped some often used functions, like locator conversion * routines. */ #ifdef HAVE_CONFIG_H # include "config.h" #endif #include #include #include #include #include #include #ifndef DOC_HIDDEN #define RADIAN (180.0 / M_PI) /* arc length for 1 degree, 60 Nautical Miles */ #define ARC_IN_KM 111.2 /* The following is contributed by Dave Hines M1CXW * * begin dph */ /* * These are the constants used when converting between Maidenhead grid * locators and longitude/latitude values. MAX_LOCATOR_PAIRS is the maximum * number of locator character pairs to convert. This number MUST NOT exceed * the number of pairs of values in loc_char_range[]. * Setting MAX_LOCATOR_PAIRS to 3 will convert the currently defined 6 * character locators. A value of 4 will convert the extended 8 character * locators described in section 3L of "The IARU region 1 VHF managers * handbook". Values of 5 and 6 will extent the format even more, to the * longest definition I have seen for locators, see * http://www.btinternet.com/~g8yoa/geog/non-ra.html * Beware that there seems to be no universally accepted standard for 10 & 12 * character locators. * * The ranges of characters which will be accepted by locator2longlat, and * generated by longlat2locator, are specified by the loc_char_range[] array. * This array may be changed without requiring any other code changes. * * For the fifth pair to range from aa to xx use: * const static int loc_char_range[] = { 18, 10, 24, 10, 24, 10 }; * * For the fifth pair to range from aa to yy use: * const static int loc_char_range[] = { 18, 10, 24, 10, 25, 10 }; * * MAX_LOCATOR_PAIRS now sets the limit locator2longlat() will convert and * sets the maximum length longlat2locator() will generate. Each function * properly handles any value from 1 to 6 so MAX_LOCATOR_PAIRS should be * left at 6. MIN_LOCATOR_PAIRS sets a floor on the shortest locator that * should be handled. -N0NB */ const static int loc_char_range[] = { 18, 10, 24, 10, 24, 10 }; #define MAX_LOCATOR_PAIRS 6 #define MIN_LOCATOR_PAIRS 1 /* end dph */ #endif /* !DOC_HIDDEN */ /** * \brief Convert DMS to decimal degrees * \param degrees Degrees, whole degrees * \param minutes Minutes, whole minutes * \param seconds Seconds, decimal seconds * \param sw South or West * * Convert degree/minute/second angle to decimal degrees angle. * \a degrees >360, \a minutes > 60, and \a seconds > 60.0 are allowed, * but resulting angle won't be normalized. * * When the variable sw is passed a value of 1, the returned decimal * degrees value will be negative (south or west). When passed a * value of 0 the returned decimal degrees value will be positive * (north or east). * * \return The angle in decimal degrees. * * \sa dec2dms() */ double HAMLIB_API dms2dec(int degrees, int minutes, double seconds, int sw) { double st; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); if (degrees < 0) { degrees = abs(degrees); } if (minutes < 0) { minutes = abs(minutes); } if (seconds < 0) { seconds = fabs(seconds); } st = (double)degrees + (double)minutes / 60. + seconds / 3600.; if (sw == 1) { return -st; } else { return st; } } /** * \brief Convert D M.MMM notation to decimal degrees * \param degrees Degrees, whole degrees * \param minutes Minutes, decimal minutes * \param sw South or West * * Convert a degrees, decimal minutes notation common on * many GPS units to its decimal degrees value. * * \a degrees > 360, \a minutes > 60.0 are allowed, but * resulting angle won't be normalized. * * When the variable sw is passed a value of 1, the returned decimal * degrees value will be negative (south or west). When passed a * value of 0 the returned decimal degrees value will be positive * (north or east). * * \return The angle in decimal degrees. * * \sa dec2dmmm() */ double HAMLIB_API dmmm2dec(int degrees, double minutes, int sw) { double st; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); if (degrees < 0) { degrees = abs(degrees); } if (minutes < 0) { minutes = fabs(minutes); } st = (double)degrees + minutes / 60.; if (sw == 1) { return -st; } else { return st; } } /** * \brief Convert decimal degrees angle into DMS notation * \param dec Decimal degrees * \param degrees Pointer for the calculated whole Degrees * \param minutes Pointer for the calculated whole Minutes * \param seconds Pointer for the calculated decimal Seconds * \param sw Pointer for the calculated SW flag * * Convert decimal degrees angle into its degree/minute/second * notation. * * When \a dec < -180 or \a dec > 180, the angle will be normalized * within these limits and the sign set appropriately. * * Upon return dec2dms guarantees 0 >= \a degrees <= 180, * 0 >= \a minutes < 60, and 0.0 >= \a seconds < 60.0. * * When \a dec is < 0.0 \a sw will be set to 1. When \a dec is * >= 0.0 \a sw will be set to 0. This flag allows the application * to determine whether the DMS angle should be treated as negative * (south or west). * * \retval -RIG_EINVAL if any of the pointers are NULL. * \retval RIG_OK if conversion went OK. * * \sa dms2dec() */ int HAMLIB_API dec2dms(double dec, int *degrees, int *minutes, double *seconds, int *sw) { int deg, min; double st; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); /* bail if NULL pointers passed */ if (!degrees || !minutes || !seconds || !sw) { return -RIG_EINVAL; } /* reverse the sign if dec has a magnitude greater * than 180 and factor out multiples of 360. * e.g. when passed 270 st will be set to -90 * and when passed -270 st will be set to 90. If * passed 361 st will be set to 1, etc. If passed * a value > -180 || < 180, value will be unchanged. */ if (dec >= 0.0) { st = fmod(dec + 180, 360) - 180; } else { st = fmod(dec - 180, 360) + 180; } /* if after all of that st is negative, we want deg * to be negative as well except for 180 which we want * to be positive. */ if (st < 0.0 && st != -180) { *sw = 1; } else { *sw = 0; } /* work on st as a positive value to remove a * bug introduced by the effect of floor() when * passed a negative value. e.g. when passed * -96.8333 floor() returns -95! Also avoids * a rounding error introduced on negative values. */ st = fabs(st); deg = (int)floor(st); st = 60. * (st - (double)deg); min = (int)floor(st); st = 60. * (st - (double)min); *degrees = deg; *minutes = min; *seconds = st; return RIG_OK; } /** * \brief Convert a decimal angle into D M.MMM notation * \param dec Decimal degrees * \param degrees Pointer for the calculated whole Degrees * \param minutes Pointer for the calculated decimal Minutes * \param sw Pointer for the calculated SW flag * * Convert a decimal angle into its degree, decimal minute * notation common on many GPS units. * * When passed a value < -180 or > 180, the value will be normalized * within these limits and the sign set apropriately. * * Upon return dec2dmmm guarantees 0 >= \a degrees <= 180, * 0.0 >= \a minutes < 60.0. * * When \a dec is < 0.0 \a sw will be set to 1. When \a dec is * >= 0.0 \a sw will be set to 0. This flag allows the application * to determine whether the D M.MMM angle should be treated as negative * (south or west). * * \retval -RIG_EINVAL if any of the pointers are NULL. * \retval RIG_OK if conversion went OK. * * \sa dmmm2dec() */ int HAMLIB_API dec2dmmm(double dec, int *degrees, double *minutes, int *sw) { int r, min; double sec; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); /* bail if NULL pointers passed */ if (!degrees || !minutes || !sw) { return -RIG_EINVAL; } r = dec2dms(dec, degrees, &min, &sec, sw); if (r != RIG_OK) { return r; } *minutes = (double)min + sec / 60; return RIG_OK; } /** * \brief Convert Maidenhead grid locator to Longitude/Latitude * \param longitude Pointer for the calculated Longitude * \param latitude Pointer for the calculated Latitude * \param locator The Maidenhead grid locator--2 through 12 char + nul string * * Convert Maidenhead grid locator to Longitude/Latitude (decimal degrees). * The locator should be in 2 through 12 chars long format. * \a locator2longlat is case insensitive, however it checks for * locator validity. * * Decimal long/lat is computed to center of grid square, i.e. given * EM19 will return coordinates equivalent to the southwest corner * of EM19mm. * * \retval -RIG_EINVAL if locator exceeds RR99xx99xx99 or exceeds length * limit--currently 1 to 6 lon/lat pairs. * \retval RIG_OK if conversion went OK. * * \bug The fifth pair ranges from aa to xx, there is another convention * that ranges from aa to yy. At some point both conventions should be * supported. * * \sa longlat2locator() */ /* begin dph */ int HAMLIB_API locator2longlat(double *longitude, double *latitude, const char *locator) { int x_or_y, paircount; int locvalue, pair; int divisions; double xy[2], ordinate; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); /* bail if NULL pointers passed */ if (!longitude || !latitude) { return -RIG_EINVAL; } paircount = strlen(locator) / 2; /* verify paircount is within limits */ if (paircount > MAX_LOCATOR_PAIRS) { paircount = MAX_LOCATOR_PAIRS; } else if (paircount < MIN_LOCATOR_PAIRS) { return -RIG_EINVAL; } /* For x(=longitude) and y(=latitude) */ for (x_or_y = 0; x_or_y < 2; ++x_or_y) { ordinate = -90.0; divisions = 1; for (pair = 0; pair < paircount; ++pair) { locvalue = locator[pair * 2 + x_or_y]; /* Value of digit or letter */ locvalue -= (loc_char_range[pair] == 10) ? '0' : (isupper(locvalue)) ? 'A' : 'a'; /* Check range for non-letter/digit or out of range */ if ((locvalue < 0) || (locvalue >= loc_char_range[pair])) { return -RIG_EINVAL; } divisions *= loc_char_range[pair]; ordinate += locvalue * 180.0 / divisions; } /* Center ordinate in the Maidenhead "square" or "subsquare" */ ordinate += 90.0 / divisions; xy[x_or_y] = ordinate; } *longitude = xy[0] * 2.0; *latitude = xy[1]; return RIG_OK; } /* end dph */ /** * \brief Convert longitude/latitude to Maidenhead grid locator * \param longitude Longitude, decimal degrees * \param latitude Latitude, decimal degrees * \param locator Pointer for the Maidenhead Locator * \param pair_count Precision expressed as lon/lat pairs in the locator * * Convert longitude/latitude (decimal degrees) to Maidenhead grid locator. * \a locator must point to an array at least \a pair_count * 2 char + '\\0'. * * \retval -RIG_EINVAL if \a locator is NULL or \a pair_count exceeds * length limit. Currently 1 to 6 lon/lat pairs. * \retval RIG_OK if conversion went OK. * * \bug \a locator is not tested for overflow. * \bug The fifth pair ranges from aa to yy, there is another convention * that ranges from aa to xx. At some point both conventions should be * supported. * * \sa locator2longlat() */ /* begin dph */ int HAMLIB_API longlat2locator(double longitude, double latitude, char *locator, int pair_count) { int x_or_y, pair, locvalue, divisions; double square_size, ordinate; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); if (!locator) { return -RIG_EINVAL; } if (pair_count < MIN_LOCATOR_PAIRS || pair_count > MAX_LOCATOR_PAIRS) { return -RIG_EINVAL; } for (x_or_y = 0; x_or_y < 2; ++x_or_y) { ordinate = (x_or_y == 0) ? longitude / 2.0 : latitude; divisions = 1; /* The 1e-6 here guards against floating point rounding errors */ ordinate = fmod(ordinate + 270.000001, 180.0); for (pair = 0; pair < pair_count; ++pair) { divisions *= loc_char_range[pair]; square_size = 180.0 / divisions; locvalue = (int)(ordinate / square_size); ordinate -= square_size * locvalue; locvalue += (loc_char_range[pair] == 10) ? '0' : 'A'; locator[pair * 2 + x_or_y] = locvalue; } } locator[pair_count * 2] = '\0'; return RIG_OK; } /* end dph */ /** * \brief Calculate the distance and bearing between two points. * \param lon1 The local Longitude, decimal degrees * \param lat1 The local Latitude, decimal degrees * \param lon2 The remote Longitude, decimal degrees * \param lat2 The remote Latitude, decimal degrees * \param distance Pointer for the distance, km * \param azimuth Pointer for the bearing, decimal degrees * * Calculate the QRB between \a lon1, \a lat1 and \a lon2, \a lat2. * * This version will calculate the QRB to a precision sufficient * for 12 character locators. Antipodal points, which are easily * calculated, are considered equidistant and the bearing is * simply resolved to be true north (0.0°). * * \retval -RIG_EINVAL if NULL pointer passed or lat and lon values * exceed -90 to 90 or -180 to 180. * \retval RIG_OK if calculations are successful. * * \return The distance in kilometers and azimuth in decimal degrees * for the short path are stored in \a distance and \a azimuth. * * \sa distance_long_path(), azimuth_long_path() */ int HAMLIB_API qrb(double lon1, double lat1, double lon2, double lat2, double *distance, double *azimuth) { double delta_long, tmp, arc, az; rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); /* bail if NULL pointers passed */ if (!distance || !azimuth) { return -RIG_EINVAL; } if ((lat1 > 90.0 || lat1 < -90.0) || (lat2 > 90.0 || lat2 < -90.0)) { return -RIG_EINVAL; } if ((lon1 > 180.0 || lon1 < -180.0) || (lon2 > 180.0 || lon2 < -180.0)) { return -RIG_EINVAL; } /* Prevent ACOS() Domain Error */ if (lat1 == 90.0) { lat1 = 89.999999999; } else if (lat1 == -90.0) { lat1 = -89.999999999; } if (lat2 == 90.0) { lat2 = 89.999999999; } else if (lat2 == -90.0) { lat2 = -89.999999999; } /* Convert variables to Radians */ lat1 /= RADIAN; lon1 /= RADIAN; lat2 /= RADIAN; lon2 /= RADIAN; delta_long = lon2 - lon1; tmp = sin(lat1) * sin(lat2) + cos(lat1) * cos(lat2) * cos(delta_long); if (tmp > .999999999999999) { /* Station points coincide, use an Omni! */ *distance = 0.0; *azimuth = 0.0; return RIG_OK; } if (tmp < -.999999) { /* * points are antipodal, it's straight down. * Station is equal distance in all Azimuths. * So take 180 Degrees of arc times 60 nm, * and you get 10800 nm, or whatever units... */ *distance = 180.0 * ARC_IN_KM; *azimuth = 0.0; return RIG_OK; } arc = acos(tmp); /* * One degree of arc is 60 Nautical miles * at the surface of the earth, 111.2 km, or 69.1 sm * This method is easier than the one in the handbook */ *distance = ARC_IN_KM * RADIAN * arc; /* Short Path */ /* Change to azimuth computation by Dave Freese, W1HKJ */ az = RADIAN * atan2(sin(lon2 - lon1) * cos(lat2), (cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(lon2 - lon1))); az = fmod(360.0 + az, 360.0); if (az < 0.0) { az += 360.0; } else if (az >= 360.0) { az -= 360.0; } *azimuth = floor(az + 0.5); return RIG_OK; } /** * \brief Calculate the long path distance between two points. * \param distance The shortpath distance * * Calculate the long path (respective of the short path) * of a given distance. * * \return the distance in kilometers for the opposite path. * * \sa qrb() */ double HAMLIB_API distance_long_path(double distance) { rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); return (ARC_IN_KM * 360.0) - distance; } /** * \brief Calculate the long path bearing between two points. * \param azimuth The shortpath bearing--0.0 to 360.0 degrees * * Calculate the long path (respective of the short path) * of a given bearing. * * \return the azimuth in decimal degrees for the opposite path or * -RIG_EINVAL upon input error (outside the range of 0.0 to 360.0). * * \sa qrb() */ double HAMLIB_API azimuth_long_path(double azimuth) { rot_debug(RIG_DEBUG_VERBOSE, "%s called\n", __func__); if (azimuth == 0.0 || azimuth == 360.0) { return 180.0; } else if (azimuth > 0.0 && azimuth < 180.0) { return 180.0 + azimuth; } else if (azimuth == 180.0) { return 0.0; } else if (azimuth > 180.0 && azimuth < 360.0) { return (180.0 - azimuth) * -1.0; } else { return -RIG_EINVAL; } } /*! @} */