From 6fd428329a1f1c2dfbb681fe1fd4041e25c17807 Mon Sep 17 00:00:00 2001 From: Lionel GUEZ <guez@lmd.ens.fr> Date: Thu, 14 Apr 2022 17:22:41 +0200 Subject: [PATCH] Use `util_eddies` Use `comp_ishape` and `open_shpc` from `util_eddies`. --- gt_segmented_cf.py | 71 +++++++++++++++------------------------------- 1 file changed, 23 insertions(+), 48 deletions(-) diff --git a/gt_segmented_cf.py b/gt_segmented_cf.py index 99e3068f..4a383ed8 100755 --- a/gt_segmented_cf.py +++ b/gt_segmented_cf.py @@ -19,27 +19,9 @@ import shapefile import datetime from numpy import loadtxt import report_graph +import util_eddies -###################### -# Utility functions: # -###################### - -def comp_ishape(d_init, ishape_last, days_1950, eddy_index): - k = days_1950 - d_init - - if k == 0: - ishape = eddy_index - 1 - else: - ishape = ishape_last[k - 1] + eddy_index - - return ishape - -####################### -# Necessary functions # -####################### - -def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, extremum, - max_speed, ishape_last): +def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, handlers): radii = 0 #[m] rossby = 0 #[1/s] @@ -53,22 +35,21 @@ def calculate_radii_and_rossby(start, end, inc, segment, e_overestim, extremum, year = date.year # calculate the location in the shapefile - 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']) + location = util_eddies.comp_ishape(handlers[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) + shapeRec = handlers[year]["readers"]["extremum"].shapeRecord(location) # rossby: 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 = shapeRec.record[4] #[m/s] - R_Vmax = max_speed[year].record(location)['r_eq_area'] * 1000 #[m] - + R_Vmax = handlers[year]["readers"]["max_speed_contour"]\ + .record(location)['r_eq_area'] * 1000 #[m] if (V_max < 100): # calculate Ro and Delta_Ro @@ -125,18 +106,13 @@ g.set_fast_edge_removal() ################################## # setup the dicts for the shapefiles and the ishape_last files -extremum = {} -max_speed = {} -ishape_last = {} shp_tr_dir = input('Directory containing %Y/SHPC_(anti|cyclo): ') +handlers = {} for year in range(1993,2019): shp_dir = path.join(shp_tr_dir, f'{year}/SHPC_{orientation}') - extremum[year] = shapefile.Reader(path.join(shp_dir, 'extremum')) - max_speed[year] = shapefile.Reader(path.join(shp_dir, 'max_speed_contour')) - ishape_last[year] = loadtxt(path.join(shp_dir, 'ishape_last.txt'), - dtype = int, delimiter = "\n", unpack = False) + handlers[year] = util_eddies.open_shpc(shp_dir) # change if there is a change over the number of days to average num_of_days_to_avg = 7 @@ -178,17 +154,19 @@ for yr in range(1993, 2019): last_year = last_date.year # calculate the location in the shapefile - - 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']) + first_loc = util_eddies.comp_ishape(handlers[first_year], + first['date_index'], + first['eddy_index']) + last_loc = util_eddies.comp_ishape(handlers[last_year], + last['date_index'], + last['eddy_index']) # grab the centers - first_pos = extremum[first_year].shape(first_loc).points[0] - last_pos = extremum[last_year].shape(last_loc).points[0] + first_pos = handlers[first_year]["readers"]["extremum"]\ + .shape(first_loc).points[0] + last_pos = handlers[last_year]["readers"]["extremum"]\ + .shape(last_loc).points[0] ##### STORE POSITIONS IN THE VPS ###### g.vp.pos_first[n] = first_pos # [deg, deg] @@ -206,8 +184,7 @@ for yr in range(1993, 2019): # First 7 days calculation first_res = calculate_radii_and_rossby(0, num_of_days_to_avg, 1, segment, e_overestim, - extremum, max_speed, - ishape_last) + handlers) # average and assign radii first_radii = first_res['radii'] / num_of_days_to_avg @@ -229,8 +206,7 @@ for yr in range(1993, 2019): len(segment) - (num_of_days_to_avg + 1), -1, segment, e_overestim, - extremum, - max_speed, ishape_last) + handlers) # Average and assign the last radii last_radii = last_res['radii'] / num_of_days_to_avg @@ -252,8 +228,7 @@ for yr in range(1993, 2019): # for the positions else: res = calculate_radii_and_rossby(0, num_of_days, 1, segment, - e_overestim, - extremum, max_speed, ishape_last) + e_overestim, handlers) # grab the days modifier modifier = res['days_modifier'] -- GitLab