diff --git a/modules/amv/intercompare.py b/modules/amv/intercompare.py index eaf15ead1a9cdfe2c219b770152098f76981f2bd..b876474c3399879431d461693ec0af76910a1772 100644 --- a/modules/amv/intercompare.py +++ b/modules/amv/intercompare.py @@ -419,25 +419,22 @@ def run_best_fit(amvs, raob_dct, dist_threshold=200, min_num_levs=20): for j in range(num_amvs): tup = raob_nbor_dct[j] if tup is not None: + cnt += 1 nbor_keys = tup[0] nbor_dsts = tup[1] + s_i = np.argsort(nbor_dsts).tolist() - last_dist = dist_threshold - keep_raob = None - for k, ll_tup in enumerate(nbor_keys): - dist = nbor_dsts[k] + for k in range(len(nbor_keys)): + dist = nbor_dsts[s_i[k]] + ll_tup = nbor_keys[s_i[k]] raob = raob_dct.get(ll_tup) if raob.shape[0] < min_num_levs: continue - if dist < last_dist: - last_dist = dist - keep_raob = raob - if keep_raob is not None: - rspd = keep_raob[:, raob_spd_idx] - rdir = keep_raob[:, raob_dir_idx] - rprs = keep_raob[:, raob_pres_idx] + rspd = raob[:, raob_spd_idx] + rdir = raob[:, raob_dir_idx] + rprs = raob[:, raob_pres_idx] aspd = amvs[j, amv_spd_idx] adir = amvs[j, amv_dir_idx] @@ -445,9 +442,11 @@ def run_best_fit(amvs, raob_dct, dist_threshold=200, min_num_levs=20): alat = amvs[j, amv_lat_idx] alon = amvs[j, amv_lon_idx] - bf = best_fit(aspd, adir, aprs, alat, alon, rspd, rdir, rprs) + bf = best_fit(aspd, adir, aprs, alat, alon, rspd, rdir, rprs, verbose=True) bf = (j, bf[0], bf[1], bf[2], bf[3]) - bestfits.append(bf) + if bf[4] == 0: + bestfits.append(bf) + break return bestfits @@ -954,8 +953,10 @@ def best_fit(amv_spd, amv_dir, amv_prs, amv_lat, amv_lon, fcst_spd, fcst_dir, fc SatwindGoodConstraint = 1 flag = 0 # Check for overall good agreement - vdiff = np.sqrt((amv_uwind - SatwindBestFitU) ** 2 + (amv_vwind - SatwindBestFitV) ** 2) - if vdiff > 4.0: + #vdiff = np.sqrt((amv_uwind - SatwindBestFitU) ** 2 + (amv_vwind - SatwindBestFitV) ** 2) + vdiff = np.sqrt((fcst_uwind - SatwindBestFitU) ** 2 + (fcst_vwind - SatwindBestFitV) ** 2) + min_vdiff = min(vdiff) + if min_vdiff > 4.0: SatwindGoodConstraint = 0 flag = 2 @@ -963,7 +964,8 @@ def best_fit(amv_spd, amv_dir, amv_prs, amv_lat, amv_lon, fcst_spd, fcst_dir, fc if SatwindGoodConstraint == 1: mm = np.where(fcst_prs > (SatwindBestFitPress + 100))[0] nn = np.where(fcst_prs < (SatwindBestFitPress - 100))[0] - if (np.sum(VecDiff[mm] < (vdiff + 2.0)) + np.sum(VecDiff[nn] < (vdiff + 2.0))) > 0: + #if (np.sum(VecDiff[mm] < (vdiff + 2.0)) + np.sum(VecDiff[nn] < (vdiff + 2.0))) > 0: + if (np.sum(vdiff[mm] < (min_vdiff + 2.0)) + np.sum(vdiff[nn] < (min_vdiff + 2.0))) > 0: SatwindGoodConstraint = 0 flag = 1 @@ -977,7 +979,7 @@ def best_fit(amv_spd, amv_dir, amv_prs, amv_lat, amv_lon, fcst_spd, fcst_dir, fc bfit_v = undef if verbose: - print('*** AMV best-fit ***') + print('*** AMV best-fit *********************************') print('AMV -> p/minspd/maxspd: {0} {1} {2}'.format(amv_prs,SatwindMinSpeed,SatwindMaxSpeed)) # print('Bestfit -> p1,p2,p3,v1,v2,v3: {0} {1} {2} {3} {4} {5}'.format(p1,p2,p3,v1,v2,v3)) print('Bestfit -> pbest,bfu,bfv,amvu,amvv,bgu,bgv: {0} {1} {2} {3} {4} {5} {6}'.format( @@ -985,7 +987,7 @@ def best_fit(amv_spd, amv_dir, amv_prs, amv_lat, amv_lon, fcst_spd, fcst_dir, fc print('Good Constraint: {0}'.format(SatwindGoodConstraint)) print('Minimum Vector Difference: {0}'.format(VecDiff[imin])) print('Vector Difference Profile: ') - print(VecDiff) + print(vdiff) print('Pressure Profile: ') print(fcst_prs) @@ -1090,8 +1092,9 @@ def best_fit_altitude(amv_spd, amv_dir, amv_alt, amv_lat, amv_lon, fcst_spd, fcs good_constraint = 1 flag = 0 # Check for overall good agreement: want vdiff less than 4 m/s - vdiff = np.sqrt((amv_uwind - sat_wind_best_fit_u) ** 2 + (amv_vwind - sat_wind_best_fit_v) ** 2) - if vdiff > 4.0: + vdiff = np.sqrt((fcst_uwind - sat_wind_best_fit_u) ** 2 + (fcst_vwind - sat_wind_best_fit_v) ** 2) + min_vdiff = min(vdiff) + if min_vdiff > 4.0: good_constraint = 0 flag = 2 @@ -1099,7 +1102,7 @@ def best_fit_altitude(amv_spd, amv_dir, amv_alt, amv_lat, amv_lon, fcst_spd, fcs if good_constraint == 1: mm = np.where(fcst_alt > (sat_wind_best_fit_alt + constraint_half_width))[0] nn = np.where(fcst_alt < (sat_wind_best_fit_alt - constraint_half_width))[0] - if (np.sum(vec_diff[mm] < (vdiff + 2.0)) + np.sum(vec_diff[nn] < (vdiff + 2.0))) > 0: + if (np.sum(vdiff[mm] < (min_vdiff + 2.0)) + np.sum(vdiff[nn] < (min_vdiff + 2.0))) > 0: good_constraint = 0 flag = 1 @@ -1117,7 +1120,7 @@ def best_fit_altitude(amv_spd, amv_dir, amv_alt, amv_lat, amv_lon, fcst_spd, fcs print('Good Constraint: {0}'.format(good_constraint)) print('Minimum Vector Difference: {0}'.format(vec_diff[imin])) print('Vector Difference Profile: ') - print(vec_diff) + print(vdiff) print('Altitude Profile: ') print(fcst_alt)