From 32efc75c64beafa0dffb81875d619f5917386837 Mon Sep 17 00:00:00 2001 From: Blake Harnden <32446120+bharnden@users.noreply.github.com> Date: Tue, 25 Feb 2020 20:40:51 -0800 Subject: [PATCH] removed legacy location translation --- daemon/core/location/corelocation.py | 279 --------------------------- daemon/core/location/utm.py | 259 ------------------------- 2 files changed, 538 deletions(-) delete mode 100644 daemon/core/location/corelocation.py delete mode 100644 daemon/core/location/utm.py diff --git a/daemon/core/location/corelocation.py b/daemon/core/location/corelocation.py deleted file mode 100644 index 6eb7d16d..00000000 --- a/daemon/core/location/corelocation.py +++ /dev/null @@ -1,279 +0,0 @@ -""" -location.py: definition of CoreLocation class that is a member of the -Session object. Provides conversions between Cartesian and geographic coordinate -systems. Depends on utm contributed module, from -https://pypi.python.org/pypi/utm (version 0.3.0). -""" - -import logging -from typing import Optional, Tuple - -from core.emulator.enumerations import RegisterTlvs -from core.location import utm - - -class CoreLocation: - """ - Member of session class for handling global location data. This keeps - track of a latitude/longitude/altitude reference point and scale in - order to convert between X,Y and geo coordinates. - """ - - name = "location" - config_type = RegisterTlvs.UTILITY.value - - def __init__(self) -> None: - """ - Creates a MobilityManager instance. - - :return: nothing - """ - # ConfigurableManager.__init__(self) - self.reset() - self.zonemap = {} - self.refxyz = (0.0, 0.0, 0.0) - self.refscale = 1.0 - self.zoneshifts = {} - self.refgeo = (0.0, 0.0, 0.0) - for n, l in utm.ZONE_LETTERS: - self.zonemap[l] = n - - def reset(self) -> None: - """ - Reset to initial state. - """ - # (x, y, z) coordinates of the point given by self.refgeo - self.refxyz = (0.0, 0.0, 0.0) - # decimal latitude, longitude, and altitude at the point (x, y, z) - self.setrefgeo(0.0, 0.0, 0.0) - # 100 pixels equals this many meters - self.refscale = 1.0 - # cached distance to refpt in other zones - self.zoneshifts = {} - - def px2m(self, val: float) -> float: - """ - Convert the specified value in pixels to meters using the - configured scale. The scale is given as s, where - 100 pixels = s meters. - - :param val: value to use in converting to meters - :return: value converted to meters - """ - return (val / 100.0) * self.refscale - - def m2px(self, val: float) -> float: - """ - Convert the specified value in meters to pixels using the - configured scale. The scale is given as s, where - 100 pixels = s meters. - - :param val: value to convert to pixels - :return: value converted to pixels - """ - if self.refscale == 0.0: - return 0.0 - return 100.0 * (val / self.refscale) - - def setrefgeo(self, lat: float, lon: float, alt: float) -> None: - """ - Record the geographical reference point decimal (lat, lon, alt) - and convert and store its UTM equivalent for later use. - - :param lat: latitude - :param lon: longitude - :param alt: altitude - :return: nothing - """ - self.refgeo = (lat, lon, alt) - # easting, northing, zone - e, n, zonen, zonel = utm.from_latlon(lat, lon) - self.refutm = ((zonen, zonel), e, n, alt) - - def getgeo(self, x: float, y: float, z: float) -> Tuple[float, float, float]: - """ - Given (x, y, z) Cartesian coordinates, convert them to latitude, - longitude, and altitude based on the configured reference point - and scale. - - :param x: x value - :param y: y value - :param z: z value - :return: lat, lon, alt values for provided coordinates - """ - # shift (x,y,z) over to reference point (x,y,z) - x -= self.refxyz[0] - y = -(y - self.refxyz[1]) - if z is None: - z = self.refxyz[2] - else: - z -= self.refxyz[2] - # use UTM coordinates since unit is meters - zone = self.refutm[0] - if zone == "": - raise ValueError("reference point not configured") - e = self.refutm[1] + self.px2m(x) - n = self.refutm[2] + self.px2m(y) - alt = self.refutm[3] + self.px2m(z) - (e, n, zone) = self.getutmzoneshift(e, n) - try: - lat, lon = utm.to_latlon(e, n, zone[0], zone[1]) - except utm.OutOfRangeError: - logging.exception( - "UTM out of range error for n=%s zone=%s xyz=(%s,%s,%s)", - n, - zone, - x, - y, - z, - ) - lat, lon = self.refgeo[:2] - return lat, lon, alt - - def getxyz(self, lat: float, lon: float, alt: float) -> Tuple[float, float, float]: - """ - Given latitude, longitude, and altitude location data, convert them - to (x, y, z) Cartesian coordinates based on the configured - reference point and scale. Lat/lon is converted to UTM meter - coordinates, UTM zones are accounted for, and the scale turns - meters to pixels. - - :param lat: latitude - :param lon: longitude - :param alt: altitude - :return: converted x, y, z coordinates - """ - # convert lat/lon to UTM coordinates in meters - e, n, zonen, zonel = utm.from_latlon(lat, lon) - _rlat, _rlon, ralt = self.refgeo - xshift = self.geteastingshift(zonen, zonel) - if xshift is None: - xm = e - self.refutm[1] - else: - xm = e + xshift - yshift = self.getnorthingshift(zonen, zonel) - if yshift is None: - ym = n - self.refutm[2] - else: - ym = n + yshift - zm = alt - ralt - - # shift (x,y,z) over to reference point (x,y,z) - x = self.m2px(xm) + self.refxyz[0] - y = -(self.m2px(ym) + self.refxyz[1]) - z = self.m2px(zm) + self.refxyz[2] - return x, y, z - - def geteastingshift(self, zonen: float, zonel: float) -> Optional[float]: - """ - If the lat, lon coordinates being converted are located in a - different UTM zone than the canvas reference point, the UTM meters - may need to be shifted. - This picks a reference point in the same longitudinal band - (UTM zone number) as the provided zone, to calculate the shift in - meters for the x coordinate. - - :param zonen: zonen - :param zonel: zone1 - :return: the x shift value - """ - rzonen = int(self.refutm[0][0]) - # same zone number, no x shift required - if zonen == rzonen: - return None - z = (zonen, zonel) - # x shift already calculated, cached - if z in self.zoneshifts and self.zoneshifts[z][0] is not None: - return self.zoneshifts[z][0] - - rlat, rlon, _ralt = self.refgeo - # ea. zone is 6deg band - lon2 = rlon + 6 * (zonen - rzonen) - # ignore northing - e2, _n2, _zonen2, _zonel2 = utm.from_latlon(rlat, lon2) - # NOTE: great circle distance used here, not reference ellipsoid! - xshift = utm.haversine(rlon, rlat, lon2, rlat) - e2 - # cache the return value - yshift = None - if z in self.zoneshifts: - yshift = self.zoneshifts[z][1] - self.zoneshifts[z] = (xshift, yshift) - return xshift - - def getnorthingshift(self, zonen: float, zonel: float) -> Optional[float]: - """ - If the lat, lon coordinates being converted are located in a - different UTM zone than the canvas reference point, the UTM meters - may need to be shifted. - This picks a reference point in the same latitude band (UTM zone letter) - as the provided zone, to calculate the shift in meters for the - y coordinate. - - :param zonen: zonen - :param zonel: zone1 - :return: calculated y shift - """ - rzonel = self.refutm[0][1] - # same zone letter, no y shift required - if zonel == rzonel: - return None - z = (zonen, zonel) - # y shift already calculated, cached - if z in self.zoneshifts and self.zoneshifts[z][1] is not None: - return self.zoneshifts[z][1] - - rlat, rlon, _ralt = self.refgeo - # zonemap is used to calculate degrees difference between zone letters - latshift = self.zonemap[zonel] - self.zonemap[rzonel] - # ea. latitude band is 8deg high - lat2 = rlat + latshift - _e2, n2, _zonen2, _zonel2 = utm.from_latlon(lat2, rlon) - # NOTE: great circle distance used here, not reference ellipsoid - yshift = -(utm.haversine(rlon, rlat, rlon, lat2) + n2) - # cache the return value - xshift = None - if z in self.zoneshifts: - xshift = self.zoneshifts[z][0] - self.zoneshifts[z] = (xshift, yshift) - return yshift - - def getutmzoneshift( - self, e: float, n: float - ) -> Tuple[float, float, Tuple[float, str]]: - """ - Given UTM easting and northing values, check if they fall outside - the reference point's zone boundary. Return the UTM coordinates in a - different zone and the new zone if they do. Zone lettering is only - changed when the reference point is in the opposite hemisphere. - - :param e: easting value - :param n: northing value - :return: modified easting, northing, and zone values - """ - zone = self.refutm[0] - rlat, rlon, _ralt = self.refgeo - if e > 834000 or e < 166000: - num_zones = (int(e) - 166000) / (utm.R / 10) - # estimate number of zones to shift, E (positive) or W (negative) - rlon2 = self.refgeo[1] + (num_zones * 6) - _e2, _n2, zonen2, zonel2 = utm.from_latlon(rlat, rlon2) - xshift = utm.haversine(rlon, rlat, rlon2, rlat) - # after >3 zones away from refpt, the above estimate won't work - # (the above estimate could be improved) - if not 100000 <= (e - xshift) < 1000000: - # move one more zone away - num_zones = (abs(num_zones) + 1) * (abs(num_zones) / num_zones) - rlon2 = self.refgeo[1] + (num_zones * 6) - _e2, _n2, zonen2, zonel2 = utm.from_latlon(rlat, rlon2) - xshift = utm.haversine(rlon, rlat, rlon2, rlat) - e = e - xshift - zone = (zonen2, zonel2) - if n < 0: - # refpt in northern hemisphere and we crossed south of equator - n += 10000000 - zone = (zone[0], "M") - elif n > 10000000: - # refpt in southern hemisphere and we crossed north of equator - n -= 10000000 - zone = (zone[0], "N") - return e, n, zone diff --git a/daemon/core/location/utm.py b/daemon/core/location/utm.py deleted file mode 100644 index b80a7d6d..00000000 --- a/daemon/core/location/utm.py +++ /dev/null @@ -1,259 +0,0 @@ -""" -utm -=== - -.. image:: https://travis-ci.org/Turbo87/utm.png - -Bidirectional UTM-WGS84 converter for python - -Usage ------ - -:: - - import utm - -Convert a (latitude, longitude) tuple into an UTM coordinate:: - - utm.from_latlon(51.2, 7.5) - >>> (395201.3103811303, 5673135.241182375, 32, 'U') - -Convert an UTM coordinate into a (latitude, longitude) tuple:: - - utm.to_latlon(340000, 5710000, 32, 'U') - >>> (51.51852098408468, 6.693872395145327) - -Speed ------ - -The library has been compared to the more generic pyproj library by running the -unit test suite through pyproj instead of utm. These are the results: - -* with pyproj (without projection cache): 4.0 - 4.5 sec -* with pyproj (with projection cache): 0.9 - 1.0 sec -* with utm: 0.4 - 0.5 sec - -Authors -------- - -* Tobias Bieniek - -License -------- - -Copyright (C) 2012 Tobias Bieniek - -Permission is hereby granted, free of charge, to any person obtaining a copy of -this software and associated documentation files (the "Software"), to deal in -the Software without restriction, including without limitation the rights to -use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies -of the Software, and to permit persons to whom the Software is furnished to do -so, subject to the following conditions: - -The above copyright notice and this permission notice shall be included in all -copies or substantial portions of the Software. - -THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE -SOFTWARE. -""" - -import math - -__all__ = ['to_latlon', 'from_latlon'] - - -class OutOfRangeError(ValueError): - pass - - -K0 = 0.9996 - -E = 0.00669438 -E2 = E * E -E3 = E2 * E -E_P2 = E / (1.0 - E) - -SQRT_E = math.sqrt(1 - E) -_E = (1 - SQRT_E) / (1 + SQRT_E) -_E3 = _E * _E * _E -_E4 = _E3 * _E - -M1 = (1 - E / 4 - 3 * E2 / 64 - 5 * E3 / 256) -M2 = (3 * E / 8 + 3 * E2 / 32 + 45 * E3 / 1024) -M3 = (15 * E2 / 256 + 45 * E3 / 1024) -M4 = (35 * E3 / 3072) - -P2 = (3 * _E / 2 - 27 * _E3 / 32) -P3 = (21 * _E3 / 16 - 55 * _E4 / 32) -P4 = (151 * _E3 / 96) - -R = 6378137 - -ZONE_LETTERS = [ - (84, None), (72, 'X'), (64, 'W'), (56, 'V'), (48, 'U'), (40, 'T'), - (32, 'S'), (24, 'R'), (16, 'Q'), (8, 'P'), (0, 'N'), (-8, 'M'), (-16, 'L'), - (-24, 'K'), (-32, 'J'), (-40, 'H'), (-48, 'G'), (-56, 'F'), (-64, 'E'), - (-72, 'D'), (-80, 'C') -] - - -def to_latlon(easting, northing, zone_number, zone_letter): - zone_letter = zone_letter.upper() - - if not 100000 <= easting < 1000000: - raise OutOfRangeError('easting out of range (must be between 100.000 m and 999.999 m)') - if not 0 <= northing <= 10000000: - raise OutOfRangeError('northing out of range (must be between 0 m and 10.000.000 m)') - if not 1 <= zone_number <= 60: - raise OutOfRangeError('zone number out of range (must be between 1 and 60)') - if not 'C' <= zone_letter <= 'X' or zone_letter in ['I', 'O']: - raise OutOfRangeError('zone letter out of range (must be between C and X)') - - x = easting - 500000 - y = northing - - if zone_letter < 'N': - y -= 10000000 - - m = y / K0 - mu = m / (R * M1) - - p_rad = (mu + P2 * math.sin(2 * mu) + P3 * math.sin(4 * mu) + P4 * math.sin(6 * mu)) - - p_sin = math.sin(p_rad) - p_sin2 = p_sin * p_sin - - p_cos = math.cos(p_rad) - - p_tan = p_sin / p_cos - p_tan2 = p_tan * p_tan - p_tan4 = p_tan2 * p_tan2 - - ep_sin = 1 - E * p_sin2 - ep_sin_sqrt = math.sqrt(1 - E * p_sin2) - - n = R / ep_sin_sqrt - r = (1 - E) / ep_sin - - c = _E * p_cos ** 2 - c2 = c * c - - d = x / (n * K0) - d2 = d * d - d3 = d2 * d - d4 = d3 * d - d5 = d4 * d - d6 = d5 * d - - latitude = (p_rad - (p_tan / r) * - (d2 / 2 - - d4 / 24 * (5 + 3 * p_tan2 + 10 * c - 4 * c2 - 9 * E_P2)) + - d6 / 720 * (61 + 90 * p_tan2 + 298 * c + 45 * p_tan4 - 252 * E_P2 - 3 * c2)) - - longitude = (d - - d3 / 6 * (1 + 2 * p_tan2 + c) + - d5 / 120 * (5 - 2 * c + 28 * p_tan2 - 3 * c2 + 8 * E_P2 + 24 * p_tan4)) / p_cos - - return (math.degrees(latitude), - math.degrees(longitude) + zone_number_to_central_longitude(zone_number)) - - -def from_latlon(latitude, longitude): - if not -80.0 <= latitude <= 84.0: - raise OutOfRangeError('latitude out of range (must be between 80 deg S and 84 deg N)') - if not -180.0 <= longitude <= 180.0: - raise OutOfRangeError('northing out of range (must be between 180 deg W and 180 deg E)') - - lat_rad = math.radians(latitude) - lat_sin = math.sin(lat_rad) - lat_cos = math.cos(lat_rad) - - lat_tan = lat_sin / lat_cos - lat_tan2 = lat_tan * lat_tan - lat_tan4 = lat_tan2 * lat_tan2 - - lon_rad = math.radians(longitude) - - zone_number = latlon_to_zone_number(latitude, longitude) - central_lon = zone_number_to_central_longitude(zone_number) - central_lon_rad = math.radians(central_lon) - - zone_letter = latitude_to_zone_letter(latitude) - - n = R / math.sqrt(1 - E * lat_sin ** 2) - c = E_P2 * lat_cos ** 2 - - a = lat_cos * (lon_rad - central_lon_rad) - a2 = a * a - a3 = a2 * a - a4 = a3 * a - a5 = a4 * a - a6 = a5 * a - - m = R * (M1 * lat_rad - - M2 * math.sin(2 * lat_rad) + - M3 * math.sin(4 * lat_rad) - - M4 * math.sin(6 * lat_rad)) - - easting = K0 * n * (a + - a3 / 6 * (1 - lat_tan2 + c) + - a5 / 120 * (5 - 18 * lat_tan2 + lat_tan4 + 72 * c - 58 * E_P2)) + 500000 - - northing = K0 * (m + n * lat_tan * (a2 / 2 + - a4 / 24 * (5 - lat_tan2 + 9 * c + 4 * c ** 2) + - a6 / 720 * (61 - 58 * lat_tan2 + lat_tan4 + 600 * c - 330 * E_P2))) - - if latitude < 0: - northing += 10000000 - - return easting, northing, zone_number, zone_letter - - -def latitude_to_zone_letter(latitude): - for lat_min, zone_letter in ZONE_LETTERS: - if latitude >= lat_min: - return zone_letter - - return None - - -def latlon_to_zone_number(latitude, longitude): - if 56 <= latitude <= 64 and 3 <= longitude <= 12: - return 32 - - if 72 <= latitude <= 84 and longitude >= 0: - if longitude <= 9: - return 31 - elif longitude <= 21: - return 33 - elif longitude <= 33: - return 35 - elif longitude <= 42: - return 37 - - return int((longitude + 180) / 6) + 1 - - -def zone_number_to_central_longitude(zone_number): - return (zone_number - 1) * 6 - 180 + 3 - - -def haversine(lon1, lat1, lon2, lat2): - """ - Calculate the great circle distance between two points - on the earth (specified in decimal degrees) - """ - # convert decimal degrees to radians - lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2]) - # haversine formula - dlon = lon2 - lon1 - dlat = lat2 - lat1 - a = math.sin(dlat / 2) ** 2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon / 2) ** 2 - c = 2 * math.asin(math.sqrt(a)) - m = 6367000 * c - return m