monthly computationsw

arrival_time_routes.py

@@ -15,6 +15,7 @@ import shapely
 from postpro import ArrivalTimeMap
 import calendar
 from typing import List
+import itertools
 
 # --- USEFUL VARIABLES --- #
 minlon, maxlon, minlat, maxlat = 50, 53, 24, 27.5
@@ -23,11 +24,19 @@ basedir = "/export/miro/students/tanselain/OpenOil/"
 sources = ['ras_laffan','ras_laffan_port', 'abu_fontas', 'umm_al_houl']  
 routes = ['Barhain', 'Doha', 'Iran', 'Ras_laffan_field', 'Ras_laffan']
 path_routes = basedir + f'/Indicators/Shipping_selected/Selected_shipping_route_2019'
+route_names = {
+    "Barhain" : "Bahrain",
+    "Iran" : "Qatar-Iran",
+    "Ras_laffan_field": 'Ras Laffan to North Field',
+    "Ras_laffan": "Ras Laffan LNG export",
+    "Doha" : 'Doha-Mesaieed'
+}
+
 
 # --- USEFUL FUNCTIONS --- #
 
 def compute_arrival_time(
-   df: gpd.GeoDataFrame, sources: List[str], year: int, month: str, day: int
+   df: gpd.GeoDataFrame, sources: List[str], year: int, month: str
 ) -> None:
     atm = ArrivalTimeMap(minlon, maxlon, minlat, maxlat, res)
     for src in sources:
@@ -35,10 +44,11 @@ def compute_arrival_time(
             nc_path = basedir + f'Output_backward_lusail/{year}/{month}'
         else:
             nc_path = basedir + f'Output_backward/{year}/{src}/{month}'
-        ncfile = nc_path+f'/out_{month}_{day}.nc'
-        if not os.path.isfile(ncfile):
-            continue
-        atm.add_arrival_times_from_file(ncfile)
+        for day in range(1,32):
+            ncfile = nc_path+f'/out_{month}_{day}.nc'
+            if not os.path.isfile(ncfile):
+                continue
+            atm.add_arrival_times_from_file(ncfile)
     a = atm.get_arrival_time_days()
     df.loc[:,'atime_days'] = a.reshape(-1)[:]
 
@@ -57,23 +67,23 @@ def main() -> None:
             poly = shapely.geometry.Polygon(np.column_stack([px,py]))
             data.append((0.0, poly))
     df = gpd.GeoDataFrame(data, columns=['atime_days','geometry'], crs=df_routes.crs)
-    with open("atime_by_routes.csv","w") as csv:
-        csv.write("year,month,day,route,arrival_time\n")
-        for yr in range(2016,2021):
-            for mn in range(1,13):
-                month = calendar.month_name[mn].lower()
-                print(f'{month} {yr}')
-                for d in range(1,32):
-                    compute_arrival_time(df, sources, yr, month, d)
-                    df_arrival = df[~np.isnan(df['atime_days'])]
-                    if df_arrival.shape[0] == 0:
-                        continue
-                    df_inter = gpd.sjoin(df_arrival, df_routes, op='intersects')
-                    df_inter['dn_total']= df_inter['atime_days'] * df_inter['int']
-                    for route in routes:
-                        df_inter_route = df_inter.loc[df_inter['Route'] == route]
-                        v = df_inter_route['atime_days'].mean()
-                        csv.write(f"{yr},{month},{d},{route},{v}\n")
+    with open("arrival_time_route.csv","w") as csv:
+        csv.write("year,month,route,arrival_time\n")
+        years = np.arange(2016,2021)
+        months = np.arange(1,13)
+        for yr,mn in itertools.product(years,months):
+            month = calendar.month_name[mn].lower()
+            print(f"{month} {yr}")
+            compute_arrival_time(df, sources, yr, month)
+            df_arrival = df[~np.isnan(df['atime_days'])]
+            if df_arrival.shape[0] == 0:
+                continue
+            df_inter = gpd.sjoin(df_arrival, df_routes, op='intersects')
+            df_inter['dn_total']= df_inter['atime_days'] * df_inter['int']
+            for route in routes:
+                df_inter_route = df_inter.loc[df_inter['Route'] == route]
+                v = df_inter_route['atime_days'].mean()
+                csv.write(f"{yr},{calendar.month_abbr[mn]},{route_names[route]},{v}\n")
                 
 if __name__ == '__main__':
     main()
\ No newline at end of file