[f] added link midpoint coordinates and RR data

FHSTP/prep.py

@@ -112,6 +112,36 @@ def load_inca_data(dir_inca: pathlib.Path) -> np.array:
     return inca_data
 
 
+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
+    """
+    
+    # convert degrees to radians
+    lon_a, lat_a = np.radians(lon_a), np.radians(lat_a)
+    lon_b, lat_b = np.radians(lon_b), np.radians(lat_b)
+    
+    # calculate midpoint
+    Bx = np.cos(lat_b) * np.cos(lon_b - lon_a)
+    By = np.cos(lat_b) * np.sin(lon_b - lon_a)
+    
+    lon_mid = lon_a + np.arctan2(By, np.cos(lat_a) + Bx)
+    lat_mid = np.arctan2(
+        np.sin(lat_a) + np.sin(lat_b),
+        np.sqrt(np.square(np.cos(lat_a) + Bx) + np.square(By))
+    )
+    
+    return np.degrees(lon_mid), np.degrees(lat_mid)
+
+
 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
@@ -348,6 +378,10 @@ def prep() -> None:
     df_trans = df_trans.drop(['TEMP_LOC_TUPLE'], axis = 1)
     
     
+    # calculate midpoint of links
+    df_config['LONGITUDE_MID'], df_config['LATITUDE_MID'] = get_midpoint(df_config['LONGITUDE_A'], df_config['LATITUDE_A'], df_config['LONGITUDE_B'], df_config['LATITUDE_B'])
+    _log('Calculated midpoint of links')
+    
     # calculate LENGTH in km between links
     df_config['LENGTH'] = get_distance(df_config['LONGITUDE_A'], df_config['LATITUDE_A'], df_config['LONGITUDE_B'], df_config['LATITUDE_B'])
     _log('Calculated distances between sites using a WGS84 ellipsoid')
@@ -385,29 +419,27 @@ def prep() -> None:
     df_link['BEGINTIME'] = df_link['BEGINTIME'].dt.tz_localize('Europe/Vienna').dt.tz_convert('UTC').dt.tz_localize(None)
     _log('Converted BEGINTIME to UTC')
     
-    
     # only use transmissions with begintime matching INCA date
     date_inca = np.datetime64(f'{str(dir_inca.stem)[0:4]}-{str(dir_inca.stem)[-4:-2]}-{str(dir_inca.stem)[-2:]}')
     df_link = df_link[df_link['BEGINTIME'].dt.date == date_inca]
     _log('Filtered transmissions to match BEGINTIME with INCA date')
     
+    
     # copy REMOTERXLEVEL to PMIN and PMAX (for aggregation in 15min window conversion)
     df_link['PMIN'] = df_link['REMOTERXLEVEL']
     df_link['PMAX'] = df_link['REMOTERXLEVEL']
     _log('Created PMIN and PMAX of REMOTERXLEVEL')
     
-    
     # convert 3min windows to 15min windows
     group_cols = [df_link['BEGINTIME'].dt.floor('15Min'), 'RADIOLINKID']
     agg_cols = {'TXLEVEL' : 'mean', 'REMOTERXLEVEL' : 'mean', 'PMIN' : 'min', 'PMAX' : 'max'}
     df_link = df_link.groupby(group_cols).agg(agg_cols).reset_index()
     _log('Converted 3min windows to 15min windows')
     
-    
     # convert BEGINTIME to RAINLINK format
     df_link['BEGINTIME'] = df_link['BEGINTIME'].dt.strftime('%Y%m%d%H%M')
     _log('Converted BEGINTIME to RAINLINK format "%Y%m%d%H%M"')
-    
+
     
     # build df with differences of sending and receiving levels
     df_diff = df_link[['RADIOLINKID', 'TXLEVEL', 'REMOTERXLEVEL']].copy()
@@ -443,6 +475,8 @@ def prep() -> None:
         'LATITUDE_A' : 'YStart',
         'LONGITUDE_B' : 'XEnd',
         'LATITUDE_B' : 'YEnd',
+        'LONGITUDE_MID' : 'XMid',
+        'LATITUDE_MID' : 'YMid',
         'LENGTH' : 'PathLength',
         'FREQUENCY' : 'Frequency',
     }
@@ -467,6 +501,7 @@ def prep() -> None:
     # set INCA RR data based on datetime and coordinates
     df_link['RRStart'] = get_inca_data(inca_data, df_link['DateTime'], df_link['XStart'], df_link['YStart'])
     df_link['RREnd'] = get_inca_data(inca_data, df_link['DateTime'], df_link['XEnd'], df_link['YEnd'])
+    df_link['RRMid'] = get_inca_data(inca_data, df_link['DateTime'], df_link['XMid'], df_link['YMid'])
     _log('Merged INCA RR data to LINK dataframe')