diff --git a/modules/amv/intercompare.py b/modules/amv/intercompare.py index c8064b39a9e38f2284d4b459b1701648acf2c845..14fdd06c638aad7d069f639750c88cdbc99298e9 100644 --- a/modules/amv/intercompare.py +++ b/modules/amv/intercompare.py @@ -1039,6 +1039,79 @@ def best_fit(amv_spd, amv_dir, amv_prs, amv_lat, amv_lon, fcst_spd, fcst_dir, fc return bf_data +def best_fit_altitude(amv_spd, amv_dir, amv_alt, amv_lat, amv_lon, fcst_spd, fcst_dir, fcst_alt, test_a=True, test_b=True): + verbose = False + + fcst_num_levels = fcst_spd.shape[0] + + alt_diff = 1000.0 + alt_top = 20000.0 + alt_bot = 0 + + flag = 3 + bf_tup = (None, None, None, flag) + if amv_alt > alt_top: + if verbose: + print('AMV location lat,lon,lat ({0},{1},{2}) is higher than max altitude {3}'.format(amv_lat, amv_lon, amv_alt, alt_top)) + return bf_tup + + # Calculate the height +/- alt_diff from the AMV height + alt_max = min((amv_alt + alt_diff), alt_top) + alt_min = max((amv_lat - alt_diff), alt_bot) + + # 1d array of indices to consider for best fit height + kk = np.where((fcst_alt > alt_min) & (fcst_alt < alt_max)) + if len(kk[0]) == 0: + if verbose: + print('AMV location lat,lon,alt ({0},{1},{2}) failed to find fcst alt around AMV'.format(amv_lat,amv_lon,amv_alt)) + return bf_tup + + # Diagnostic field: Find the model minimum speed and maximum speed within PressDiff of the AMV. + if verbose: + sat_wind_min_spd = min(fcst_spd[kk]) + sat_wind_max_spd = max(fcst_spd[kk]) + + # Compute U anv V for both AMVs and forecast + amv_uwind = -amv_spd * np.sin((np.pi/180.0)*amv_dir) + amv_vwind = -amv_spd * np.cos((np.pi/180.0)*amv_dir) + fcst_uwind = -fcst_spd * np.sin((np.pi/180.0)*fcst_dir) + fcst_vwind = -fcst_spd * np.cos((np.pi/180.0)*fcst_dir) + + # Calculate the vector difference between the AMV and model background at all levels. + vec_diff = np.sqrt((amv_uwind - fcst_uwind) ** 2 + (amv_vwind - fcst_vwind) ** 2) + + # Find the model level of best-fit pressure, from the minimum vector difference. + min_vec_diff = min(vec_diff[kk]) + imin = np.where(vec_diff == min_vec_diff) + if len(imin[0]) == 0: + if verbose: + print('AMV location lat,lon,alt ({0},{1},{2}) failed to find min vector difference in layers around AMV'.format(amv_lat,amv_lon,amv_alt)) + return bf_tup + + imin = imin[0][0] + imin = kk[0][imin] + + if imin == fcst_num_levels - 1: + if verbose: + print('AMV location lat,lon,alt ({0},{1},{2}) failed to find min vector difference in layers around AMV'.format(amv_lat,amv_lon,amv_alt)) + return bf_tup + + a2 = fcst_alt[imin] + v2 = vec_diff[imin] + + if imin == 0: + return bf_tup + + a2 = fcst_alt[imin] + v2 = vec_diff[imin] + + a1 = fcst_alt[imin-1] + v2 = vec_diff[imin-1] + a3 = fcst_alt[imin+1] + v3 = vec_diff[imin+1] + + return bf_tup + def get_press_bin_ranges(lop, hip, bin_size=100): bin_ranges = [] delp = hip - lop