From 14ef373f3e5996a4e6e33e551301f31654324f8b Mon Sep 17 00:00:00 2001 From: Lionel GUEZ <guez@lmd.ens.fr> Date: Tue, 19 Apr 2022 17:48:26 +0200 Subject: [PATCH] Polish --- cost_function.py | 30 +++++++++++++++--------------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/cost_function.py b/cost_function.py index 8fbe63f3..859cac62 100755 --- a/cost_function.py +++ b/cost_function.py @@ -244,13 +244,13 @@ for n in g.vertices(): for edge in g.edges(): source_node = edge.source() target_node = edge.target() - + cf = -10000 - lat_for_conv = (g.vp.pos_last[source_node][1] + + lat_for_conv = (g.vp.pos_last[source_node][1] + g.vp.pos_first[target_node][1]) / 2 # latitude needed for conversion of degrees to kilometers lat_for_conv = math.radians(lat_for_conv) # need to convert to radians - + # because of the wrapping issue (360° wrapping incorrectly to 0°), we check for that here lon_diff = abs(g.vp.pos_last[source_node][0] - g.vp.pos_first[target_node][0]) if (lon_diff > 300): @@ -258,38 +258,38 @@ for edge in g.edges(): # calculate Delta_cent: numbers used for conversion obtained from: # https://stackoverflow.com/questions/1253499/simple-calculations-for-working-with-lat-lon-and-km-distance - Delta_Cent = math.sqrt(( lon_diff * 111.32 * math.cos(lat_for_conv) )**2 + + Delta_Cent = math.sqrt(( lon_diff * 111.32 * math.cos(lat_for_conv) )**2 + ( (g.vp.pos_last[source_node][1] - g.vp.pos_first[target_node][1]) * 110.574 )**2) - + # calculate the first term first_term = ((Delta_Cent - delta_cent_mean)/delta_cent_std) ** 2 - - # Rossbies: + + # Rossbies: if (g.vp.first_av_ros[target_node] and g.vp.last_av_ros[source_node]): Delta_Ro = g.vp.last_av_ros[source_node] - g.vp.first_av_ros[target_node] else: print("At least one of the rossbies is invalid.") #Delta_Ro = delta_ro_mean Delta_Ro = 0 - + # Calculate the second term second_term = ((Delta_Ro - delta_ro_mean)/delta_ro_std ) ** 2 - + # R_Vmax 1 and 2 already exist, just get the delta - + Delta_R_Vmax = g.vp.last_av_rad[source_node] - g.vp.first_av_rad[target_node] - + # Calculate the third term third_term = ((Delta_R_Vmax - delta_r_mean)/delta_r_std) ** 2 - + ############################# # calculate the cost function ############################# - + cf = math.sqrt(first_term + second_term + third_term) - + # assign as weight to the edge - g.ep.nl_cost_function[edge] = cf + g.ep.nl_cost_function[edge] = cf ################################ -- GitLab