cleaned up example scripts (#366)

* cleaned up example scripts

* updated testcase to avoid failure on 2022 duplicate extremes in DDL

tests/examples/predictie_2019_frommergedcomp_WSdwarsstroming.py → tests/examples/WSdwarsstroming_predictie.py

@@ -35,7 +35,7 @@
     ts_prediction_HANSWT = hatyan.prediction(comp=COMP_merged_HANSWT, times=times_pred)
     ts_prediction_TERNZN = hatyan.prediction(comp=COMP_merged_TERNZN, times=times_pred)
 
-    ts_prediction_diff = ts_prediction_TERNZN-ts_prediction_HANSWT # TODO: metadata is dropped
+    ts_prediction_diff = ts_prediction_TERNZN-ts_prediction_HANSWT # TODO: metadata is retained from TERNZN
     ts_prediction_diff['values'] = ts_prediction_diff['values'].round(2) #round to cm
     ts_prediction_diff_HWLW = hatyan.calc_HWLW(ts=ts_prediction_diff)
 
@@ -46,19 +46,15 @@
     list_matig.append(ts_prediction_diff_matig)
     list_sterk.append(ts_prediction_diff_sterk)
 
-    #fig, (ax1,ax2) = hatyan.plot_timeseries(ts=ts_prediction_diff)
     ax.plot(ts_prediction_diff.index, ts_prediction_diff['values'],'-',linewidth=0.7,markersize=1, label='verval')
     ax.plot(ts_prediction_diff_matig.index, ts_prediction_diff_matig['values'],'oy',markersize=5, label='matig (>%.2fm)'%(value_matig))
     ax.plot(ts_prediction_diff_sterk.index, ts_prediction_diff_sterk['values'],'or',markersize=5, label='sterk (>%.2fm)'%(value_sterk))
     ax.plot([ts_prediction_diff.index[0],ts_prediction_diff.index[-1]],[value_matig,value_matig],'-y')
     ax.plot([ts_prediction_diff.index[0],ts_prediction_diff.index[-1]],[value_sterk,value_sterk],'-r')
-    #ax1.set_ylim(-1.2,1.7)
     ax.legend(loc=1)
     ax.grid()
     ax.set_xlabel('Tijd')
     ax.set_ylabel('astro verval %d [m] (TERNZN-HANSWT)'%yr)
-    #hatyan.write_dia(ts=ts_ext_prediction_main, station=current_station, vertref='NAP', filename='prediction_HWLW_%s_%s_main.dia'%(times_step_pred, current_station))
-    #hatyan.write_dia(ts=ts_ext_prediction_clean, station=current_station, vertref='NAP', filename='prediction_HWLW_%s_%s_agger345.dia'%(times_step_pred, current_station))
 fig.tight_layout()
 fig.savefig('prediction_WSdwarsstroming')
 

tests/examples/astrog_example.py

@@ -14,7 +14,6 @@
 """
 
 import os
-import numpy as np
 import datetime as dt
 import pandas as pd
 import matplotlib.pyplot as plt
@@ -24,7 +23,6 @@
 #dir_testdata = 'P:\\1209447-kpp-hydraulicaprogrammatuur\\hatyan\\hatyan_data_acceptancetests'
 dir_testdata = 'C:\\DATA\\hatyan_data_acceptancetests'
 
-# script settings
 compare2fortran = True #requires validation data
 
 start_date_utc = pd.Timestamp(2000, 1, 1, tz="UTC")
@@ -38,39 +36,38 @@
 
 dT_fortran = True #True is best comparison to fortran, False is more precise 
 
-pdtocsv_kwargs = dict(index=False, sep=',', date_format='%Y-%m-%d %H:%M:%S %Z', float_format='%9.5f', na_rep='--:--          ')
+pdtocsv_kwargs = dict(index=False, sep=',', date_format='%Y-%m-%d %H:%M:%S %Z', float_format='%9.5f', na_rep='--:--')
 
-#%% calculate astrog arrays
-# lunar culmination times, parallax, declination
+# calculate astrog lunar culmination times, parallax, declination
 culminations_python = hatyan.astrog_culminations(tFirst=start_date_utc, tLast=end_date_utc, dT_fortran=dT_fortran)
 culminations_python.to_csv('moon_culminations.csv',**pdtocsv_kwargs)
 
-# lunar phases
+# calculate astrog lunar phases
 phases_python = hatyan.astrog_phases(tFirst=start_date_met, tLast=end_date_met, dT_fortran=dT_fortran)
 phases_python.to_csv('moon_phases.csv',**pdtocsv_kwargs)
 
-# moonrise and -set
+# calculate astrog moonrise and -set
 moonriseset_python = hatyan.astrog_moonriseset(tFirst=start_date_met, tLast=end_date_met, dT_fortran=dT_fortran)
 moonriseset_python.to_csv('moon_riseset.csv',**pdtocsv_kwargs)
 moonriseset_python_perday = hatyan.convert2perday(moonriseset_python)
 moonriseset_python_perday.to_csv('moon_riseset_perday.csv',**pdtocsv_kwargs)
 
-# sunrise and -set
+# calculate astrog sunrise and -set
 sunriseset_python = hatyan.astrog_sunriseset(tFirst=start_date_met, tLast=end_date_met, dT_fortran=dT_fortran)
 sunriseset_python.to_csv('sun_riseset.csv',**pdtocsv_kwargs)
 sunriseset_python_perday = hatyan.convert2perday(sunriseset_python)
 sunriseset_python_perday.to_csv('sun_riseset_perday.csv',**pdtocsv_kwargs)
 
-# lunar anomalies
+# calculate astrog lunar anomalies
 anomalies_python = hatyan.astrog_anomalies(tFirst=start_date_met, tLast=end_date_met, dT_fortran=dT_fortran)
 anomalies_python.to_csv('anomalies.csv',**pdtocsv_kwargs)
 
-# astronomical seasons
+# calculate astrog astronomical seasons
 seasons_python = hatyan.astrog_seasons(tFirst=start_date_met, tLast=end_date_met, dT_fortran=dT_fortran)
 seasons_python.to_csv('seasons.csv',**pdtocsv_kwargs)
 
 if compare2fortran:
-    #%% load fortran results
+    # load fortran results
     pkl_culm = os.path.join(dir_testdata,'other','astrog20_2000_2011.pkl')
     pkl_phas = os.path.join(dir_testdata,'other','astrog30_phases_2000_2011.pkl')
     txt_phas = os.path.join(dir_testdata,'other','maanfasen.txt')
@@ -108,7 +105,7 @@
     seasons_fortran = pd.read_pickle(pkl_seas).set_index('datetime')
     seasons_fortran = seasons_fortran.loc[start_date_naive:end_date_naive]
 
-    #%% plot results (differences)
+    # plot results (differences)
     fig, (ax1,ax2,ax3) = hatyan.plot_astrog_diff(culminations_python, culminations_fortran, typeLab=['lower','upper'], timeBand=[-.18,.18])
     fig.savefig('culmination_differences.png')
 

tests/examples/compare_foremanschureman_freqs.py

@@ -23,7 +23,6 @@
 freqs_pd['const'] = freqs_pd.index
 freqs_pd['freq_absdiff'] = np.abs(freqs_pd.iloc[:,0]-freqs_pd.iloc[:,1])
 freqs_pd['freq_bigdiff'] = freqs_pd['freq_absdiff']>10e-9
-#freqs_pd['freq_nan'] = (np.isnan(freqs_pd.iloc[:,0]-freqs_pd.iloc[:,1]))
 
 v0_pd = pd.concat([v0_pd_schu[0],v0_pd_for[0]],axis=1)
 v0_pd['v0_absdiff'] = (np.abs(v0_pd.iloc[:,0]-v0_pd.iloc[:,1])+0.5*np.pi)%np.pi-0.5*np.pi

tests/examples/compare_foremanschureman_prediction.py

@@ -18,20 +18,18 @@
 
 for current_station in selected_stations:
     const_list = hatyan.get_const_list_hatyan('year')
-    file_data_comp0_lastyear = os.path.join(dir_testdata,'predictie2019','%s_obs1.txt'%(current_station))
+    file_data_comp0 = os.path.join(dir_testdata,'predictie2019','%s_obs1.txt'%(current_station))
 
-    ts_measurements_group0_lastyear = hatyan.read_dia(filename=file_data_comp0_lastyear, station=current_station)
-    #ts_measurements_group0 = hatyan.read_dia(filename=file_data_comp0, station=current_station)
-    #ts_measurements_group0 = hatyan.crop_timeseries(ts_measurements_group0, times_ext=[dt.datetime(2012,1,1),dt.datetime(2013,1,1)])
+    ts_measurements_group0 = hatyan.read_dia(filename=file_data_comp0, station=current_station)
 
     stats_row = pd.DataFrame(index=[current_station])
     for fu_alltimes in [True,False]:
         xfac = False
         #prediction and comparison to measurements
-        COMP_schu = hatyan.analysis(ts=ts_measurements_group0_lastyear, const_list=const_list, source='schureman', fu_alltimes=fu_alltimes, xfac=xfac)
-        ts_prediction_schu = hatyan.prediction(comp=COMP_schu, times=ts_measurements_group0_lastyear.index)
-        COMP_for = hatyan.analysis(ts=ts_measurements_group0_lastyear, const_list=const_list, source='foreman', fu_alltimes=fu_alltimes, xfac=xfac)
-        ts_prediction_for = hatyan.prediction(comp=COMP_for, times=ts_measurements_group0_lastyear.index)
+        COMP_schu = hatyan.analysis(ts=ts_measurements_group0, const_list=const_list, source='schureman', fu_alltimes=fu_alltimes, xfac=xfac)
+        ts_prediction_schu = hatyan.prediction(comp=COMP_schu, times=ts_measurements_group0.index)
+        COMP_for = hatyan.analysis(ts=ts_measurements_group0, const_list=const_list, source='foreman', fu_alltimes=fu_alltimes, xfac=xfac)
+        ts_prediction_for = hatyan.prediction(comp=COMP_for, times=ts_measurements_group0.index)
 
         fig, (ax1,ax2) = hatyan.plot_timeseries(ts=ts_prediction_schu, ts_validation=ts_prediction_for)
         ax1.set_title(f'{current_station} fualltimes={fu_alltimes}')

tests/examples/numbering_FEWS_PG.py

@@ -17,7 +17,6 @@
 
 analyse_ts_bool = False
 
-#station_name = data_pred.var_stations.loc[0,'node_id']
 station_name = 'HOEKVHLD'
 
 file_meas = os.path.join(dir_testdata,'other','FEWS_202010221200_testdata_S_2.nc')
@@ -68,15 +67,6 @@
 ax2.set_ylim(-1,2)
 fig.savefig(file_ncout.replace('.nc','_nrs.png'))
 
-"""
-hatyan.write_netcdf(ts=ts_prediction, station=station_name, vertref='NAP', filename=file_ncout, ts_ext=ts_ext_prediction_nos, tzone_hr=0, nosidx=False)
-#from dfm_tools.get_nc_helpers import get_ncvardimlist
-#vars_pd, dims_pd = get_ncvardimlist(file_nc=file_ncout)
-#data_nc_checkLWval = get_ncmodeldata(file_nc=file_ncout,varname='time_LW',timestep='all',station=0)
-data_ncout = Dataset(file_ncout)
-data_ncout.variables['waterlevel_astro_LW_numbers']
-data_ncout.variables['waterlevel_astro_HW_numbers']
-"""
 hatyan.write_netcdf(ts=ts_prediction, filename=file_ncout_nosidx, ts_ext=ts_ext_prediction_nos, nosidx=True)
 # add fake extra station to show how this works
 ts_prediction.attrs["station"] = station_name+"_COPY"