From 161f348e432e123e81561d1d9ebd31ea3187f09c Mon Sep 17 00:00:00 2001
From: Lionel GUEZ <guez@lmd.ens.fr>
Date: Thu, 14 Apr 2022 16:00:12 +0200
Subject: [PATCH] Do not store in memory all the shapefiles extremum

---
 gt_segmented_cf.py | 42 +++++++++++++++++++++++-------------------
 1 file changed, 23 insertions(+), 19 deletions(-)

diff --git a/gt_segmented_cf.py b/gt_segmented_cf.py
index 07810487..dd644218 100755
--- a/gt_segmented_cf.py
+++ b/gt_segmented_cf.py
@@ -38,7 +38,8 @@ def comp_ishape(d_init, ishape_last, days_1950, eddy_index):
 # Necessary functions #
 #######################
 
-def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, ext_shp_rcd, max_rcd, ishape_last):
+def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, extremum,
+                               max_rcd, ishape_last):
     
     radii = 0 #[m]
     rossby = 0 #[1/s]
@@ -51,17 +52,19 @@ def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, ext_shp_rc
         year = date.year
         
         # calculate the location in the shapefile
-        d_init_1 = ext_shp_rcd[year][0].record[1]
+        d_init_1 = extremum[year].record(0)[1]
         location = comp_ishape(d_init_1, ishape_last[year], 
                                current_eddy['date_index'], current_eddy['eddy_index'])
         
         # now that we have the location in the shapefiles, we need to get the radius and the rossby number
         
+        shapeRec = extremum[year].shapeRecord(location)
         # rossby:
-        lat_in_deg = ext_shp_rcd[year][location].shape.points[0][1] #[deg]
+        lat_in_deg = shapeRec.shape.points[0][1]
+        #[deg]
         f = 2*2*math.pi/(24*3600)*math.sin(math.radians(lat_in_deg)) # [1/s]
         
-        V_max = ext_shp_rcd[year][location].record[4] #[m/s]
+        V_max = shapeRec.record[4] #[m/s]
         R_Vmax = max_rcd[year][location]['r_eq_area'] * 1000 #[m]
         
         
@@ -151,19 +154,15 @@ for yr in range(1993, 2019):
     # load the shapefiles and the ishape_last files
     print('Loading extrema and ishape_last files.')
 
-    ext_shp_rcd = {}
     max_rcd = {}
     
-    ext_shp_rcd[yr] = extremum[yr].shapeRecords()
     max_rcd[yr] = max_speed[yr].records()
         
     # check if we are reaching the last two years or the last year
     if (yr <= 2017):
-        ext_shp_rcd[yr+1] = extremum[yr + 1].shapeRecords()
         max_rcd[yr+1] = max_speed[yr + 1].records()
         
         if (yr <= 2016):
-            ext_shp_rcd[yr+2] = extremum[yr + 2].shapeRecords()
             max_rcd[yr+2] = max_speed[yr + 2].records()    
 
     print("Finished reading shapefiles.")
@@ -205,16 +204,16 @@ for yr in range(1993, 2019):
         
             # calculate the location in the shapefile
             
-            d_init_1 = ext_shp_rcd[first_year][0].record[1]
-            d_init_2 = ext_shp_rcd[last_year][0].record[1]
+            d_init_1 = extremum[first_year].record(0)[1]
+            d_init_2 = extremum[last_year].record(0)[1]
             
             first_loc = comp_ishape(d_init_1, ishape_last[first_year], first['date_index'], first['eddy_index'])
             last_loc = comp_ishape(d_init_2, ishape_last[last_year], last['date_index'], last['eddy_index'])
             
             # grab the centers
             
-            first_pos = ext_shp_rcd[first_year][first_loc].shape.points[0]
-            last_pos = ext_shp_rcd[last_year][last_loc].shape.points[0]
+            first_pos = extremum[first_year].shape(first_loc).points[0]
+            last_pos = extremum[last_year].shape(last_loc).points[0]
             
             ##### STORE POSITIONS IN THE VPS ######
             g.vp.pos_first[n] = first_pos # [deg, deg]
@@ -230,8 +229,10 @@ for yr in range(1993, 2019):
                 last_rossby = 0 # []
                 
                 # First 7 days calculation
-                first_res = calculate_radii_and_rossby(0, num_of_days_to_avg, 1, segment, e_overestim, 
-                                                 ext_shp_rcd, max_rcd, ishape_last)
+                first_res = calculate_radii_and_rossby(0, num_of_days_to_avg,
+                                                       1, segment, e_overestim,
+                                                       extremum, max_rcd,
+                                                       ishape_last)
                 
                 # average and assign radii
                 first_radii = first_res['radii'] / num_of_days_to_avg
@@ -250,9 +251,11 @@ for yr in range(1993, 2019):
                     
                 # Last 7 days calculation
                 last_res = calculate_radii_and_rossby(len(segment) - 1, 
-                                                      len(segment) - (num_of_days_to_avg + 1), -1,
-                                                      segment, e_overestim, ext_shp_rcd,
-                                                      max_rcd, ishape_last)                  
+                                                      len(segment) - (num_of_days_to_avg + 1),
+                                                      -1,
+                                                      segment, e_overestim,
+                                                      extremum,
+                                                      max_rcd, ishape_last)
                     
                 # Average and assign the last radii
                 last_radii = last_res['radii'] / num_of_days_to_avg
@@ -273,8 +276,9 @@ for yr in range(1993, 2019):
             # days over which to average, the values will be the same except
             # for the positions
             else:
-                res = calculate_radii_and_rossby(0, num_of_days, 1, segment, e_overestim, 
-                                                 ext_shp_rcd, max_rcd, ishape_last)
+                res = calculate_radii_and_rossby(0, num_of_days, 1, segment,
+                                                 e_overestim,
+                                                 extremum, max_rcd, ishape_last)
                 
                 # grab the days modifier
                 modifier = res['days_modifier']
-- 
GitLab