Commit 885cdad8 authored by Fabian Kovac's avatar Fabian Kovac
Browse files

[i] general improvements

parent d0943f43
......@@ -2,8 +2,8 @@
# Title: Data Preparation for LINK Configs and Transmissions
# Author: Fabian Kovac <ds191008@fhstp.ac.at>
# Team: University of Applied Sciences St. Pölten
# Version: 1.4
# Last changed: 2021-06-24
# Version: 1.5
# Last changed: 2021-06-28
#
import sys
......@@ -41,7 +41,7 @@ def parse_arguments() -> argparse.Namespace:
# add transmissions parameter to parser
required_args = parser.add_argument_group('required arguments')
required_args.add_argument('-c', '--config', type = str, required = True, help = 'Path to Config-File')
required_args.add_argument('-t', '--transmissions', type = str, required = True, help = 'Path to Config-File')
required_args.add_argument('-t', '--transmissions', type = str, required = True, help = 'Path to Transmissions-File')
required_args.add_argument('-i', '--inca', type = str, required = True, help = 'Path to Inca-Dir')
# parse arguments
......@@ -71,10 +71,10 @@ def _log(msg: str) -> None:
def load_inca_file(file_inca: pathlib.Path) -> np.array:
"""Loads loads gzipped INCA data to a given file
Parameters:
file_inca (pathlib.Path): Path to INCA file
Returns:
x (np.array): Matrix with shape (401, 701)
"""
......@@ -90,17 +90,17 @@ def load_inca_file(file_inca: pathlib.Path) -> np.array:
def load_inca_data(dir_inca: pathlib.Path) -> np.array:
"""Loads inca files to a given inca dir
Parameters:
dir_inca (pathlib.Path): Directory for INCA files
Returns:
inca_Data (np.array): Tensor with shape (96, 401, 701)
"""
# initialize tensor with 96 15min intervals
inca_data = np.zeros((96, 401, 701))
# load inca dates from inca dir
for i, file_inca in enumerate(sorted([file for file in dir_inca.iterdir() if file.is_file()])):
# load zipped ascii data
......@@ -114,13 +114,13 @@ def load_inca_data(dir_inca: pathlib.Path) -> np.array:
def get_midpoint(lon_a: np.array, lat_a: np.array, lon_b: np.array, lat_b: np.array) -> np.array:
"""Calculcates the midpoint between two utm coordinate pairs
Parameters:
lon_a (np.array): Longitudes of point A
lat_a (np.array): Latitudes of point A
lon_b (np.array): Longitudes of point B
lat_b (np.array): Latitudes of point B
Returns:
lon, lat (np.array): Vectors with longitudes and latitudes of the midpoints
"""
......@@ -143,56 +143,59 @@ def get_midpoint(lon_a: np.array, lat_a: np.array, lon_b: np.array, lat_b: np.ar
def get_distance(lon_a: np.array, lat_a: np.array, lon_b: np.array, lat_b: np.array) -> np.array:
"""Calculcates distance between two coordinates in km
using a rotation-ellipsoid in cartesian coordinates out of polar coordiantes
"""Calculcates distance between two coordinates in km using a WGS84 rotation-ellipsoid
Parameters:
lon_a (np.array): Longitudes of point A
lat_a (np.array): Latitudes of point A
lon_b (np.array): Longitudes of point B
lat_b (np.array): Latitudes of point B
Returns:
length (np.array): Vector with distances in km (can directly be assigned a pandas column)
length (np.array): Vector with distances in km
"""
# constants (euqator radius and pole radius in km)
# constants (equator radius and pole radius in km)
# r_equator is the 'semi-major axis' and r_pole the 'semi-minor axis' on a WGS84 ellipsoid
r_equator = 6378.137
r_pole = 6356.7523142
# calculate rotation-ellipsoid in cartesian coordinates out of polar coordiantes
# calculate points on rotation-ellipsoid
# point A
za = np.sin(np.radians(lat_a)) * r_pole
ra = np.cos(np.radians(lat_a)) * r_equator
xa = np.sin(np.radians(lon_a)) * ra
ya = np.cos(np.radians(lon_a)) * ra
# point B
zb = np.sin(np.radians(lat_b)) * r_pole
rb = np.cos(np.radians(lat_b)) * r_equator
xb = np.sin(np.radians(lon_b)) * rb
yb = np.cos(np.radians(lon_b)) * rb
# calculate distances between point a and point b
# calculate differences between A and B
dx = xa - xb
dy = ya - yb
dz = za - zb
length = np.sqrt(np.square(dx) + np.square(dy) + np.square(dz))
return length
# calculate distances between points using the square root of the sum of squares
distance = np.sqrt(np.square(dx) + np.square(dy) + np.square(dz))
return distance
def utm_to_lambert(lon: np.array, lat: np.array) -> tuple:
"""Convert WGS 84 UTM coordinates to a Lambert Conic Conformal Projection
"""Convert UTM coordinates to a Lambert Conic Conformal Projection
Parameters:
lon (np.array): Vector containing longitudes
lat (np.array): Vector containing latitudes
Returns:
x, y (tuple(np.array, np.array)): tuple of two vectors with coordinates
x, y (np.array): tuple of two vectors with coordinates
"""
# convert utm coordinates (angles of degree) to radians
# convert utm coordinates to radians
lon = np.radians(lon)
lat = np.radians(lat)
......@@ -208,7 +211,7 @@ def utm_to_lambert(lon: np.array, lat: np.array) -> tuple:
lat0 = np.radians(47.5)
# INCA grid: 701x401 km
# compensate for point of reference is in the middle of the grid
# compensate for point of reference (middle of the grid)
# false easting: half of 701km = 350500m
# false northing: half of 401km = 200500m
x0 = 350500
......@@ -224,7 +227,7 @@ def utm_to_lambert(lon: np.array, lat: np.array) -> tuple:
p = R * F * (1 / np.tan(np.pi/4 + lat/2))**n
p0 = R * F * (1 / np.tan(np.pi/4 + lat0/2))**n
# calculate lambert conic conformal x and y
# calculate x and y
x = p * np.sin(n * (lon - lon0)) + x0
y = p0 - p * np.cos(n * (lon - lon0)) + y0
......@@ -233,17 +236,24 @@ def utm_to_lambert(lon: np.array, lat: np.array) -> tuple:
def lambert_to_inca_idx(x: np.array, y: np.array) -> tuple:
"""Convert x and y of Lambert Conic Conformal Projection to INCA index
(rounded lambert coordinates)
Parameters:
x (np.array): Vector containing x values in meter of Lambert Conic Conformal Projection
y (np.array): Vector containing y values in meter of in Lambert Conic Conformal Projection
x (np.array): Vector containing x values in meter of Lambert Conformal Conic Projection
y (np.array): Vector containing y values in meter of in Lambert Conformal Conic Projection
Returns:
ix, iy (tuple(np.array, np.array)): tuple of indices for INCA data
ix, iy (np.array): indices for INCA data
"""
return np.round(x/1000, decimals = 0).astype(int), np.round(y/1000, decimals = 0).astype(int)
# convert from m to km
x /= 1000
y /= 1000
# convert coordinates to INCA index (coordinates as int)
ix = x.astype(int)
iy = y.astype(int)
return ix, iy
def get_inca_data(datetimes: np.array, lon_a: np.array, lat_a: np.array, lon_b: np.array, lat_b: np.array, lengths: np.array) -> np.array:
......@@ -269,8 +279,9 @@ def get_inca_data(datetimes: np.array, lon_a: np.array, lat_a: np.array, lon_b:
window = 15
inca_times = sorted([str(i * datetime.timedelta(seconds = window))[2:] for i in range(24*60//window)])
# generate times of LINK data
link_times = datetimes.map(lambda x: f'{x[-4:-2]}:{x[-2:]}')
# generate times of LINK data (last 4 chars of datetime vector)
time_func = np.vectorize(lambda x: f'{x[-4:-2]}:{x[-2:]}')
link_times = time_func(datetimes)
# get INCA indices of LINK times
idx_times = np.searchsorted(inca_times, link_times)
......@@ -311,13 +322,7 @@ def get_inca_data(datetimes: np.array, lon_a: np.array, lat_a: np.array, lon_b:
def prep() -> None:
"""Data preparation for LINK config and transmissions
Parameters:
file_config (pathlib.Path): Config File
file_trans (pathlib.Path): Transmissions File
dir_inca (pathlib.Path): INCA Directory
"""
"""Data preparation for LINK config and transmissions"""
_log('\n******************************** READ FILES ********************************')
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment