From 26cdc248e2ebe9142a0fb1c95151b71d55164c55 Mon Sep 17 00:00:00 2001 From: Carlo Russo Date: Tue, 19 Nov 2024 16:48:34 +0100 Subject: [PATCH] Use of fill algorithm to merge disputed territories in water demand historic maps. Extended MSWX dataset to 1975 --- .../water-demand-historic/step2_domestic_demand.py | 5 ++++- .../water-demand-historic/step3a_industrial_demand.py | 5 ++++- src/lisfloodutilities/water-demand-historic/tools.py | 4 ++-- 3 files changed, 10 insertions(+), 4 deletions(-) diff --git a/src/lisfloodutilities/water-demand-historic/step2_domestic_demand.py b/src/lisfloodutilities/water-demand-historic/step2_domestic_demand.py index a38859a..8da7f00 100644 --- a/src/lisfloodutilities/water-demand-historic/step2_domestic_demand.py +++ b/src/lisfloodutilities/water-demand-historic/step2_domestic_demand.py @@ -340,7 +340,10 @@ def main(): # Load MSWX air temperature data mswx_monthly = np.zeros((mapsize_global[0],mapsize_global[1],12),dtype=np.single)*np.NaN for month in np.arange(1,13): - dset = Dataset(os.path.join(config['mswx_folder'],'Past','Temp','Monthly',str(year)+str(month).zfill(2)+'.nc')) + if year<1979: #MSWX dataset starts from 1979 + dset = Dataset(os.path.join(config['mswx_folder'],'Past','Temp','Monthly',str(1979)+str(month).zfill(2)+'.nc')) + else: + dset = Dataset(os.path.join(config['mswx_folder'],'Past','Temp','Monthly',str(year)+str(month).zfill(2)+'.nc')) data = np.squeeze(np.array(dset.variables['air_temperature'])) mswx_monthly[:,:,month-1] = imresize_mean(data,mapsize_global) mswx_avg = np.mean(mswx_monthly,axis=2) diff --git a/src/lisfloodutilities/water-demand-historic/step3a_industrial_demand.py b/src/lisfloodutilities/water-demand-historic/step3a_industrial_demand.py index 2455de0..70fa833 100644 --- a/src/lisfloodutilities/water-demand-historic/step3a_industrial_demand.py +++ b/src/lisfloodutilities/water-demand-historic/step3a_industrial_demand.py @@ -731,7 +731,10 @@ def main(): ndays = monthrange(year,month)[1] mswx_month = np.zeros((1800,3600,ndays),dtype=np.single)*np.NaN for day in np.arange(1,ndays+1): - filepath = os.path.join(config['mswx_folder'],'Past','Temp','Daily',datetime(year,month,day).strftime('%Y%j')+'.nc') + if year<1979: #MSWX dataset starts from 1979 + filepath = os.path.join(config['mswx_folder'],'Past','Temp','Daily',datetime(1979,month,day).strftime('%Y%j')+'.nc') + else: + filepath = os.path.join(config['mswx_folder'],'Past','Temp','Daily',datetime(year,month,day).strftime('%Y%j')+'.nc') if os.path.isfile(filepath): dset = Dataset(filepath) mswx_month[:,:,day-1] = np.squeeze(np.array(dset.variables['air_temperature'])) diff --git a/src/lisfloodutilities/water-demand-historic/tools.py b/src/lisfloodutilities/water-demand-historic/tools.py index f3ae8fc..dfc744c 100644 --- a/src/lisfloodutilities/water-demand-historic/tools.py +++ b/src/lisfloodutilities/water-demand-historic/tools.py @@ -35,8 +35,8 @@ def load_country_code_map(filepath,mapsize): assert (country_code_map_refmatrix[0]==-180) & (country_code_map_refmatrix[3]==90) country_code_map = resize(country_code_map,mapsize,order=0,mode='edge',anti_aliasing=False) - country_code_map[country_code_map==158] = 156 # Add Taiwan to China - country_code_map[country_code_map==736] = 729 # South Sudan missing from map + # apply the fill algorithm to fill the gaps of the Disputed Territories (code>1000), so that the fill algorithm will set each pixel value to their nearest neighbour + country_code_map = fill(country_code_map, country_code_map>1000) # CR: check here: the fill algorithm includes also zero as a neighbour country code, so that isles and some land area will be replaced by ocean country_code_map[country_code_map==0] = np.NaN return country_code_map