diff --git a/rainforest/performance/eval_calculate.py b/rainforest/performance/eval_calculate.py index 0119f4a..ff736e0 100644 --- a/rainforest/performance/eval_calculate.py +++ b/rainforest/performance/eval_calculate.py @@ -56,9 +56,11 @@ def calcScoresStations(precip_ref : pd.DataFrame, th_ref, th_est = threshold if isinstance(threshold, list) else [threshold] * 2 - perf_scores = pd.DataFrame(columns=['YESref_YESpred','NOref_NOpred','NOref_YESpred','YESref_NOpred', + score_names = ['YESref_YESpred','NOref_NOpred','NOref_YESpred','YESref_NOpred', 'RMSE','corr_p','scatter','logBias', 'n_values', - 'n_events_db','sum_ref_db','sum_pred_db'], + 'n_events_db','sum_ref_db','sum_pred_db'] + score_db_names = ['n_events_db','sum_ref_db','sum_pred_db', 'RMSE','corr_p','scatter','logBias'] + perf_scores = pd.DataFrame(columns=score_names, index=precip_ref.columns.unique()) stations = METSTATIONS.copy() @@ -94,6 +96,7 @@ def calcScoresStations(precip_ref : pd.DataFrame, np.round(10*np.log10(sum_pred_db/sum_ref_db),decimals=4) else: logging.info('No measurements for station {}'.format(station_id)) + perf_scores.loc[station_id, ['logBias','RMSE','corr_p','scatter']] = np.nan continue try: scores = det_cont_fct(precip_est[station_id][double_cond.index].to_numpy(dtype=float), @@ -104,6 +107,7 @@ def calcScoresStations(precip_ref : pd.DataFrame, perf_scores.at[station_id, 'scatter']= np.round(scores['scatter'],decimals=4) except: logging.info('Could not calculate scores for station {}'.format(station_id)) + perf_scores.loc[station_id, ['RMSE','corr_p','scatter']] = np.nan return perf_scores diff --git a/rainforest/performance/eval_get_estimates.py b/rainforest/performance/eval_get_estimates.py index 9e38d13..f0a5ac2 100644 --- a/rainforest/performance/eval_get_estimates.py +++ b/rainforest/performance/eval_get_estimates.py @@ -28,7 +28,8 @@ from ..common.io_data import read_cart -def getGaugeObservations(gaugefolder, t0=None, t1=None, slf_stations=False): +def getGaugeObservations(gaugefolder, t0=None, t1=None, slf_stations=False, + missing2nan=False): """_summary_ Args: @@ -50,7 +51,9 @@ def getGaugeObservations(gaugefolder, t0=None, t1=None, slf_stations=False): dtype = {'TIMESTAMP':int, 'STATION': str}) gauge_all = gauge_all.compute().drop_duplicates() - gauge_all = gauge_all.replace(-9999,np.nan) + + if missing2nan: + gauge_all = gauge_all.replace(-9999,np.nan) # Assure that datetime object is in UTC if t0 != None : @@ -291,7 +294,8 @@ def extractEstimatesFromMaps(self, slf_stations=False, tagg_hourly=True, save_ou #----------------------------------- logging.info('Get gauge observations') gauge_all = getGaugeObservations(gaugefolder=self.gaugefolder, - t0 = self.tstart, t1= self.tend, slf_stations = False) + t0 = self.tstart, t1= self.tend, slf_stations = False, + missing2nan=False) stations = np.unique(gauge_all['STATION']) # Get all model files @@ -360,6 +364,8 @@ def extractEstimatesFromMaps(self, slf_stations=False, tagg_hourly=True, save_ou col = gauge_all.loc[gauge_all['STATION'] == ss, 'RRE150Z0'].to_frame(name=ss) col.set_index(gauge_all.loc[(gauge_all['STATION'] == ss),'TIME'], inplace=True) precip_qpe['GAUGE'][ss] = col + precip_qpe['GAUGE'][ss].fillna(0, inplace=True) + precip_qpe['GAUGE'][ss].loc[precip_qpe['GAUGE'][ss] < 0] = np.nan for m in list_models: precip_qpe[m] = pd.DataFrame(precip_qpe[m], index=tstamps_10min, columns=precip_qpe['GAUGE'].columns)