From f23c7c89e613d37c217f648e91566b006a90ad95 Mon Sep 17 00:00:00 2001 From: Bo Dong postdoc Date: Sat, 13 Jan 2024 23:02:30 -0800 Subject: [PATCH] Changes to be committed: (use "git reset HEAD ..." to unstage) modified: driver_monsoon_sperber.py Untracked files: (use "git add ..." to include in what will be committed) ../../pcmdi_utils/ --- .../monsoon_sperber/driver_monsoon_sperber.py | 71 ++++++++++++------- 1 file changed, 47 insertions(+), 24 deletions(-) diff --git a/pcmdi_metrics/monsoon_sperber/driver_monsoon_sperber.py b/pcmdi_metrics/monsoon_sperber/driver_monsoon_sperber.py index e2fca5304..499a4bcc3 100644 --- a/pcmdi_metrics/monsoon_sperber/driver_monsoon_sperber.py +++ b/pcmdi_metrics/monsoon_sperber/driver_monsoon_sperber.py @@ -50,6 +50,7 @@ import numpy as np import xarray as xr import xcdat as xc +import pandas as pd import pcmdi_metrics from pcmdi_metrics import resources @@ -346,13 +347,14 @@ def tree(): mip, model, exp, run, "monsoon_sperber", startYear, endYear ) - fout = cdms2.open( - os.path.join( + file_path = os.path.join( outdir(output_type="diagnostic_results"), output_filename + ".nc", - ), - "w", - ) + ) + try: + fout = xr.open_dataset(file_path, mode='a') # 'a' stands for append mode + except FileNotFoundError: + fout = xr.Dataset() # Plotting setup if plot: @@ -519,6 +521,7 @@ def tree(): ) pentad_time_series = [] + time_coords = np.array([], dtype='datetime64') for d_sub_aave_chunk in list_d_sub_aave_chunks: # ignore when chunk length is shorter than defined @@ -529,6 +532,10 @@ def tree(): axis=0, skipna=True ).compute() pentad_time_series.append(float(ave_chunk)) + datetime_str = str(d_sub_aave_chunk['time'][0].values) + datetime = pd.to_datetime([datetime_str[:10]]) + time_coords = np.concatenate([time_coords, datetime]) + time_coords = pd.to_datetime(time_coords) if debug: print( @@ -543,17 +550,20 @@ def tree(): pentad_time_series, ref_length, debug=debug ) - pentad_time_series_cumsum = np.cumsum(pentad_time_series) - pentad_time_series = xr.DataArray(pentad_time_series) - pentad_time_series.attrs["units"] = d.units + pentad_time_series_cumsum = np.cumsum(pentad_time_series) + pentad_time_series = xr.DataArray(pentad_time_series,dims='time', name=region + "_" + str(year)) + pentad_time_series.attrs['units'] = str(d.units.values) + pentad_time_series.coords['time'] = time_coords + + pentad_time_series_cumsum = xr.DataArray(pentad_time_series_cumsum,dims='time',\ + name=region + "_" + str(year) + "_cumsum") + pentad_time_series_cumsum.attrs['units'] = str(d.units.values) + pentad_time_series_cumsum.coords['time'] = time_coords if nc_out: # Archive individual year time series in netCDF file - fout.write(pentad_time_series, id=region + "_" + str(year)) - fout.write( - pentad_time_series_cumsum, - id=region + "_" + str(year) + "_cumsum", - ) + pentad_time_series.to_netcdf(file_path, mode='a') + pentad_time_series_cumsum.to_netcdf(file_path, mode='a') """ if plot: @@ -610,6 +620,23 @@ def tree(): composite_pentad_time_series_cumsum_normalized = metrics_result[ "frac_accum" ] + + composite_pentad_time_series = xr.DataArray(composite_pentad_time_series, \ + dims='time', name=region + "_comp") + composite_pentad_time_series.attrs['units'] = str(d.units) + composite_pentad_time_series.coords['time'] = time_coords + + composite_pentad_time_series_cumsum = xr.DataArray(composite_pentad_time_series_cumsum, \ + dims='time', name=region + "_comp_cumsum") + composite_pentad_time_series_cumsum.attrs['units'] = str(d.units) + composite_pentad_time_series_cumsum.coords['time'] = time_coords + + composite_pentad_time_series_cumsum_normalized = xr.DataArray(\ + composite_pentad_time_series_cumsum_normalized,\ + dims='time', name=region + "_comp_cumsum_fraction") + composite_pentad_time_series_cumsum_normalized.attrs['units'] = str(d.units) + composite_pentad_time_series_cumsum_normalized.coords['time'] = time_coords + if model == "obs": dict_obs_composite[reference_data_name][region] = {} dict_obs_composite[reference_data_name][ @@ -632,15 +659,11 @@ def tree(): # Archice in netCDF file if nc_out: - fout.write(composite_pentad_time_series, id=region + "_comp") - fout.write( - composite_pentad_time_series_cumsum, - id=region + "_comp_cumsum", - ) - fout.write( - composite_pentad_time_series_cumsum_normalized, - id=region + "_comp_cumsum_fraction", - ) + + composite_pentad_time_series.to_netcdf(file_path, mode='a') + composite_pentad_time_series_cumsum.to_netcdf(file_path, mode='a') + composite_pentad_time_series_cumsum_normalized.to_netcdf(file_path, mode='a') + if region == list_monsoon_regions[-1]: fout.close() @@ -664,7 +687,7 @@ def tree(): ymin=0, ymax=composite_pentad_time_series_cumsum_normalized[ idx - ], + ].item(), c="red", ls="--", ) @@ -687,7 +710,7 @@ def tree(): ymin=0, ymax=dict_obs_composite[reference_data_name][region][ idx - ], + ].item(), c="blue", ls="--", )