diff --git a/nlmod/data/geotop/REF_GTP_LITHO_CLASS.csv b/nlmod/data/geotop/REF_GTP_LITHO_CLASS.csv new file mode 100644 index 00000000..195bcf7e --- /dev/null +++ b/nlmod/data/geotop/REF_GTP_LITHO_CLASS.csv @@ -0,0 +1,12 @@ +LITHO_CLASS_CD,DESCRIPTION,VOXEL_NR,SEQ_NR,RED_DEC,GREEN_DEC,BLUE_DEC +a,antropogeen,0,0,200,200,200 +v,organisch materiaal (veen),1,1,157,78,64 +k,klei,2,2,0,146,0 +kz,"kleiig zand, zandige klei en leem",3,3,194,207,92 +zf,zand fijn,5,5,255,255,0 +zm,zand midden,6,6,243,225,6 +zg,zand grof,7,7,231,195,22 +g,grind,8,8,216,163,32 +she,schelpen,9,9,95,95,255 +z,zand niet gedifferentieerd,Null,Null,243,225,6 +o,overig,Null,Null,144,144,144 diff --git a/nlmod/data/geotop/REF_GTP_LITHO_CLASS.xlsx b/nlmod/data/geotop/REF_GTP_LITHO_CLASS.xlsx new file mode 100644 index 00000000..65e4f8c0 Binary files /dev/null and b/nlmod/data/geotop/REF_GTP_LITHO_CLASS.xlsx differ diff --git a/nlmod/data/geotop/REF_GTP_STR_UNIT.csv b/nlmod/data/geotop/REF_GTP_STR_UNIT.csv new file mode 100644 index 00000000..702cea8c --- /dev/null +++ b/nlmod/data/geotop/REF_GTP_STR_UNIT.csv @@ -0,0 +1,107 @@ +STR_UNIT_CD,DESCRIPTION,VOXEL_NR,SEQ_NR,RED_DEC,GREEN_DEC,BLUE_DEC +AAOP,"Antropogene afzettingen, opgebrachte grond",1000,1,200,200,200 +AAES,"Antropogene afzettingen, esdekken",1005,2,110,110,110 +NIGR,"Formatie van Nieuwkoop, Laagpakket van Griendtsveen",1010,3,132,88,44 +NINB,"Formatie van Nieuwkoop, Laag van Nij Beets",1045,4,180,90,20 +NASC,"Formatie van Naaldwijk, Laagpakket van Schoorl",1020,5,255,255,160 +ONAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (gedeelte boven NAZA)",1030,6,152,220,117 +NAZA,"Formatie van Naaldwijk, Laagpakket van Zandvoort",1040,7,254,183,56 +NAWAZU,"Formatie van Naaldwijk, Laagpakket van Walcheren, Zuiderzee Laag",1048,14,177,160,199 +NAWAAL,"Formatie van Naaldwijk, Laagpakket van Walcheren, Almere Laag",1049,15,96,73,122 +NAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren, gelegen onder Formatie van Naaldwijk, Laagpakket van Zandvoort",1050,16,152,220,117 +BHEC,Formatie van Echteld (gedeelte buiten NIHO),1060,,170,255,245 +OEC,Formatie van Echteld (gedeelte boven NIHO),1070,13,170,255,245 +NAWOBE,"Formatie van Naaldwijk, Laagpakket van Wormer, Laag van Bergen",1080,18,0,157,155 +KK1,Kreekrak Formatie (gedeelte boven NIHO),1085,17,68,138,112 +NIFL,"Formatie van Nieuwkoop, Flevomeer Laag",1089,19,252,213,180 +NIHO,"Formatie van Nieuwkoop, Hollandveen Laagpakket",1090,20,208,130,40 +NAZA2,"Formatie van Naaldwijk, Laagpakket van Zandvoort (gedeelte onder NIHO)",1095,21,254,183,56 +NAWO,"Formatie van Naaldwijk, Laagpakket van Wormer",1100,27,63,165,76 +NWNZ,"Formatie van Naaldwijk, laagpakketten van Wormer en Zandvoort",1110,,63,165,76 +NAWOVE,"Formatie van Naaldwijk, Laagpakket van Wormer, Laag van Velsen",1120,29,189,183,107 +KK2,Kreekrak Formatie (gedeelte onder NIHO),1125,30,68,138,112 +NIBA,"Formatie van Nieuwkoop, Basisveen Laag",1130,31,152,47,10 +NA,Formatie van Naaldwijk,2000,12,105,200,100 +EC,Formatie van Echteld,2010,28,170,255,245 +NI,Formatie van Nieuwkoop,2020,,208,130,40 +KK,Kreekrak Formatie,2030,,68,138,112 +BXKO,"Formatie van Boxtel, Laagpakket van Kootwijk",3000,34,255,255,190 +BXSI,"Formatie van Boxtel, Laagpakket van Singraven",3010,,190,190,0 +BXSI1,Formatie van Boxtel Laagpakket van Singraven (bovenste deel),3011,35,190,190,0 +BXWI,"Formatie van Boxtel, Laagpakket van Wierden",3020,36,255,255,80 +BXSI2,Formatie van Boxtel Laagpakket van Singraven (onderste deel),3012,37,255,255,130 +BXWIKO,"Formatie van Boxtel, laagpakketten van Wierden en Kootwijk",3025,38,255,255,80 +BXWISIKO,"Formatie van Boxtel, laagpakketten van Wierden, Singraven en Kootwijk",3030,39,255,255,80 +BXDE,"Formatie van Boxtel, Laagpakket van Delwijnen",3040,40,225,225,130 +BXDEKO,"Formatie van Boxtel, laagpakketten van Delwijnen en Kootwijk",3045,41,225,225,130 +BXSC,"Formatie van Boxtel, Laagpakket van Schimmert",3050,42,255,205,0 +BXLM,"Formatie van Boxtel, Laagpakket van Liempde",3060,43,225,190,0 +BXBS,"Formatie van Boxtel, Laagpakket van Best",3090,46,235,205,0 +BX,Formatie van Boxtel,3100,47,255,235,0 +KRWY,"Formatie van Kreftenheye, Laag van Wijchen",4000,48,86,0,0 +KRBXDE,"Formatie van Kreftenheye en Formatie van Boxtel, Laagpakket van Delwijnen",4010,49,176,48,96 +KRZU,"Formatie van Kreftenheye, Laagpakket van Zutphen",4020,50,136,0,0 +KROE,"Formatie van Kreftenheye, gelegen onder de Eem Formatie",4030,61,176,48,96 +KRTW,"Formatie van Kreftenheye, Laagpakket van Twello",4040,51,156,18,46 +KR,Formatie van Kreftenheye,4050,52,176,48,96 +BEOM,"Formatie van Beegden, Laagpakket van Oost-Maarland",4055,32,129,105,123 +BEWY,"Formatie van Beegden, Laag van Wijchen",4060,53,170,155,180 +BERO,"Formatie van Beegden, Laag van Rosmalen",4070,54,160,140,155 +BE,Formatie van Beegden,4080,55,200,200,255 +KW1,Formatie van Koewacht (kleiige top),4085,56,152,139,0 +KW,Formatie van Koewacht,4090,57,172,169,43 +WB,Formatie van Woudenberg,4100,58,137,67,30 +EE,Eem Formatie,4110,59,210,255,115 +EEWB,Formatie van Woudenberg en Eem Formatie,4120,60,210,255,115 +DR,Formatie van Drente,5000,62,255,127,80 +DRGI,"Formatie van Drente, Laagpakket van Gieten",5010,63,235,97,30 +GE,Door landijs gestuwde afzettingen,5020,64,156,156,156 +DN,Formatie van Drachten,5030,65,250,250,210 +URTY,"Formatie van Urk, Laagpakket van Tynje",5040,66,169,163,87 +PE,Formatie van Peelo,5050,67,238,130,238 +UR,Formatie van Urk,5060,68,189,183,107 +ST,Formatie van Sterksel,5070,69,205,92,92 +AP,Formatie van Appelscha,5080,70,218,165,32 +SY,Formatie van Stramproy,5090,71,255,228,181 +PZ,Formatie van Peize,5100,72,255,255,0 +WA,Formatie van Waalre,5110,73,255,165,0 +PZWA,Formatie van Peize en Formatie van Waalre,5120,74,255,204,0 +MS,Formatie van Maassluis,5130,75,135,206,235 +KI,Kiezeloöliet Formatie,5140,76,188,143,143 +OO,Formatie van Oosterhout,5150,77,118,157,39 +IE,Formatie van Inden,5160,78,236,121,193 +VI,Formatie van Ville,5170,79,153,102,0 +BR,Formatie van Breda,5180,80,108,188,150 +VE,Formatie van Veldhoven,5185,,102,100,16 +RUBO,"Rupel Formatie, Laagpakket van Boom",5190,81,154,78,163 +RU,Rupel Formatie,5200,82,184,123,238 +TOZEWA,"Formatie van Tongeren, Laagpakket van Zelzate, Laag van Watervliet",5210,83,80,144,194 +TOGO,"Formatie van Tongeren, Laagpakket van Goudsberg",5220,84,70,129,169 +TO,Formatie van Tongeren,5230,85,90,159,219 +DOAS,"Formatie van Dongen, Laagpakket van Asse",5240,86,186,146,141 +DOIE,"Formatie van Dongen, Laagpakket van Ieper",5250,87,206,176,191 +DO,Formatie van Dongen,5260,88,216,191,216 +LA,Formatie van Landen,5270,89,208,32,144 +HT,Formatie van Heijenrath,5280,90,178,34,34 +HO,Formatie van Houthem,5290,91,210,105,30 +MT,Formatie van Maastricht,5300,92,255,160,102 +GU,Formatie van Gulpen,5310,93,245,222,179 +VA,Formatie van Vaals,5320,94,21,153,79 +AK,Formatie van Aken,5330,95,152,231,205 +AEC,Formatie van Echteld (geulafzettingen generatie A),6000,9,118,147,60 +ABEOM,"Formatie van Beegden, Laagpakket van Oost-Maarland (geulafzettingen generatie A)",6005,33,118,147,60 +ANAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (geulafzettingen generatie A)",6010,8,118,147,60 +ANAWO,"Formatie van Naaldwijk, Laagpakket van Wormer (geulafzettingen generatie A)",6020,,118,147,60 +BEC,Formatie van Echteld (geulafzettingen generatie B),6100,11,102,205,171 +BNAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (geulafzettingen generatie B)",6110,10,102,205,171 +BNAWO,"Formatie van Naaldwijk, Laagpakket van Wormer (geulafzettingen generatie B)",6120,,102,205,171 +CEC,Formatie van Echteld (geulafzettingen generatie C),6200,22,170,196,255 +CNAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (geulafzettingen generatie C)",6210,,170,196,255 +CNAWO,"Formatie van Naaldwijk, Laagpakket van Wormer (geulafzettingen generatie C)",6220,,170,196,255 +DEC,Formatie van Echteld (geulafzettingen generatie D),6300,24,102,153,205 +DNAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (geulafzettingen generatie D)",6310,,102,153,205 +DNAWO,"Formatie van Naaldwijk, Laagpakket van Wormer (geulafzettingen generatie D)",6320,23,102,153,205 +EEC,Formatie van Echteld (geulafzettingen generatie E),6400,26,27,101,175 +ENAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren (geulafzettingen generatie E)",6410,,27,101,175 +ENAWO,"Formatie van Naaldwijk, Laagpakket van Wormer (geulafzettingen generatie E)",6420,25,27,101,175 +NN,Niet formeel ingedeelde afzettingen of onbekend,0,,85,85,85 diff --git a/nlmod/data/geotop/REF_GTP_STR_UNIT.xlsx b/nlmod/data/geotop/REF_GTP_STR_UNIT.xlsx new file mode 100644 index 00000000..15447611 Binary files /dev/null and b/nlmod/data/geotop/REF_GTP_STR_UNIT.xlsx differ diff --git a/nlmod/data/geotop/geo_eenheden.csv b/nlmod/data/geotop/geo_eenheden.csv index 4a114bf8..25d84280 100644 --- a/nlmod/data/geotop/geo_eenheden.csv +++ b/nlmod/data/geotop/geo_eenheden.csv @@ -6,13 +6,14 @@ strat,code,name 1030,ONAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren, gelegen boven Formatie van Naaldwijk, Laagpakket van Zandvoort" 1040,NAZA,"Formatie van Naaldwijk, Laagpakket van Zandvoort" +1045,NINB,"Formatie van Nieuwkoop, Laag van Nij Beets" 1050,NAWA,"Formatie van Naaldwijk, Laagpakket van Walcheren, gelegen onder Formatie van Naaldwijk, Laagpakket van Zandvoort" 1060,BHEC,"Formatie van Echteld, gelegen buiten de verbreiding van de Formatie van Nieuwkoop, Hollandveen Laagpakket" 1070,OEC,"Formatie van Echteld, gelegen boven de Formatie van Nieuwkoop, Hollandveen Laagpakket" 1080,NAWOBE,"Formatie van Naaldwijk, Laagpakket van Wormer, Laag van Bergen" 1090,NIHO,"Formatie van Nieuwkoop, Hollandveen Laagpakket" -1095,NINB,"Formatie van Nieuwkoop, Laagpakket van Nij Beets" +1095,NAZA2,"Formatie van Naaldwijk, Laagpakket van Zandvoort (gedeelte onder NIHO)" 1100,NAWO,"Formatie van Naaldwijk, Laagpakket van Wormer" 1110,NWNZ,"Formatie van Naaldwijk, Laagpakketten van Wormer en Zandvoort (gecombineerde eenheid in modelgebied Zeeland)" diff --git a/nlmod/dims/base.py b/nlmod/dims/base.py index 5c32aacd..a53fa7c4 100644 --- a/nlmod/dims/base.py +++ b/nlmod/dims/base.py @@ -206,7 +206,7 @@ def to_model_ds( return ds -def extrapolate_ds(ds, mask=None, layer="layer"): +def extrapolate_ds(ds, mask=None, layer="layer", mask_values=None): """Fill missing data in layermodel based on nearest interpolation. Used for ensuring layer model contains data everywhere. Useful for @@ -223,6 +223,10 @@ def extrapolate_ds(ds, mask=None, layer="layer"): variable. The default is None. layer : str, optional The name of the layer dimension. The default is 'layer'. + mask_values : np.ndarray, optional + Boolean mask for each cell, with a value of True if its value is used to fill + data in mask. When mask_values is None, it is determined from mask. The default + is None. Returns ------- @@ -234,6 +238,8 @@ def extrapolate_ds(ds, mask=None, layer="layer"): if not mask.any(): # all of the model cells are is inside the known area return ds + if mask_values is None: + mask_values = ~mask if mask.all(): raise (ValueError("The model only contains NaNs")) if "gridtype" in ds.attrs and ds.gridtype == "vertex": @@ -243,7 +249,7 @@ def extrapolate_ds(ds, mask=None, layer="layer"): else: x, y = np.meshgrid(ds.x, ds.y) dims = ("y", "x") - points = np.stack((x[~mask], y[~mask]), axis=1) + points = np.stack((x[mask_values], y[mask_values]), axis=1) xi = np.stack((x[mask], y[mask]), axis=1) # geneterate the tree only once, to increase speed tree = cKDTree(points) @@ -254,11 +260,11 @@ def extrapolate_ds(ds, mask=None, layer="layer"): data = ds[key].data if ds[key].dims == dims: if np.isnan(data[mask]).sum() > 0: # do not update if no NaNs - data[mask] = data[~mask][i] + data[mask] = data[mask_values][i] elif ds[key].dims == (layer,) + dims: for lay in range(len(ds[layer])): if np.isnan(data[lay, mask]).sum() > 0: # do not update if no NaNs - data[lay, mask] = data[lay, ~mask][i] + data[lay, mask] = data[lay, mask_values][i] else: logger.warning( f"Data variable '{key}' not extrapolated because " diff --git a/nlmod/dims/grid.py b/nlmod/dims/grid.py index df7bf1b9..de2fd6aa 100644 --- a/nlmod/dims/grid.py +++ b/nlmod/dims/grid.py @@ -702,8 +702,6 @@ def update_ds_from_layer_ds(ds, layer_ds, method="nearest", **kwargs): Dataset with variables from layer_ds. """ if not layer_ds.layer.equals(ds.layer): - # do not change the original Dataset - layer_ds = layer_ds.copy() # update layers in ds drop_vars = [] for var in ds.data_vars: @@ -1201,7 +1199,7 @@ def gdf_to_data_array_struc( DESCRIPTION. """ warnings.warn( - "The method gdf_to_data_array_struc is deprecated. Please use gdf_to_da instead", + "The method gdf_to_data_array_struc is deprecated. Please use gdf_to_da instead.", DeprecationWarning, ) @@ -1770,7 +1768,7 @@ def get_thickness_from_topbot(top, bot): or (layer, icell2d). """ warnings.warn( - "The method get_thickness_from_topbot is deprecated. Please use nlmod.layers.calculate_thickness instead", + "The method get_thickness_from_topbot is deprecated. Please use nlmod.layers.calculate_thickness instead.", DeprecationWarning, ) diff --git a/nlmod/dims/layers.py b/nlmod/dims/layers.py index c468582f..31d2d77e 100644 --- a/nlmod/dims/layers.py +++ b/nlmod/dims/layers.py @@ -472,6 +472,7 @@ def kheq_combined_layers(kh, thickness, reindexer): kheq = np.nansum( thickness.data[v, :, :] * kh.data[v, :, :], axis=0 ) / np.nansum(thickness.data[v, :, :], axis=0) + kheq[np.isinf(kheq)] = np.nan else: kheq = kh.data[v] da_kh.data[k] = kheq @@ -512,6 +513,7 @@ def kveq_combined_layers(kv, thickness, reindexer): kveq = np.nansum(thickness.data[v, :, :], axis=0) / np.nansum( thickness.data[v, :, :] / kv.data[v, :, :], axis=0 ) + kveq[np.isinf(kveq)] = np.nan else: kveq = kv.data[v] da_kv.data[k] = kveq @@ -536,7 +538,7 @@ def combine_layers_ds( ds : xarray.Dataset xarray Dataset containing information about layers (layers, top and bot) - combine_layers : list of tuple of ints + combine_layers : list of tuple of ints, or dict of layer names list of tuples, with each tuple containing integers indicating layer indices to combine into one layer. E.g. [(0, 1), (2, 3)] will combine layers 0 and 1 into a single layer (with index 0) and layers @@ -586,6 +588,17 @@ def combine_layers_ds( ds = ds.copy() ds["top"] = ds["botm"] + calculate_thickness(ds) + if isinstance(combine_layers, dict): + new_layer_names = combine_layers.keys() + combine_layers = [ + tuple(np.where(ds.layer.isin(x))[0]) for x in combine_layers.values() + ] + # make sure there are no layers in between + assert np.all([(np.diff(x) == 1).all() for x in combine_layers]) + else: + new_layer_names = [ds.layer.data[x[0]] for x in combine_layers] + new_layer_names = dict(zip(combine_layers, new_layer_names)) + da_dict = {} new_top, new_bot, reindexer = layer_combine_top_bot( @@ -604,10 +617,10 @@ def combine_layers_ds( if kv is not None: logger.info(f"Calculate equivalent '{kv}' for combined layers.") da_dict[kv] = kveq_combined_layers(ds[kv], thickness, reindexer) - if kD is not None: + if kD is not None and kD in ds: logger.info(f"Calculate value '{kD}' for combined layers with sum.") da_dict[kD] = sum_param_combined_layers(ds[kD], reindexer) - if c is not None: + if c is not None and c in ds: logger.info(f"Calculate value '{c}' for combined layers with sum.") da_dict[c] = sum_param_combined_layers(ds[c], reindexer) @@ -615,8 +628,9 @@ def combine_layers_ds( layer_names = [] for _, i in reindexer.items(): if isinstance(i, tuple): - i = i[0] - layercode = ds[layer].data[i] + layercode = new_layer_names[i] + else: + layercode = ds[layer].data[i] layer_names.append(layercode) # assign new layer names diff --git a/nlmod/dims/resample.py b/nlmod/dims/resample.py index dc5594a5..e14b5ca5 100644 --- a/nlmod/dims/resample.py +++ b/nlmod/dims/resample.py @@ -583,7 +583,7 @@ def structured_da_to_ds(da, ds, method="average", nodata=np.nan): def extent_to_polygon(extent): logger.warning( "nlmod.resample.extent_to_polygon is deprecated. " - "Use nlmod.util.extent_to_polygon instead" + "Use nlmod.util.extent_to_polygon instead." ) from ..util import extent_to_polygon @@ -594,7 +594,7 @@ def get_extent_polygon(ds, rotated=True): """Get the model extent, as a shapely Polygon.""" logger.warning( "nlmod.resample.get_extent_polygon is deprecated. " - "Use nlmod.grid.get_extent_polygon instead" + "Use nlmod.grid.get_extent_polygon instead." ) from .grid import get_extent_polygon @@ -605,7 +605,7 @@ def affine_transform_gdf(gdf, affine): """Apply an affine transformation to a geopandas GeoDataFrame.""" logger.warning( "nlmod.resample.affine_transform_gdf is deprecated. " - "Use nlmod.grid.affine_transform_gdf instead" + "Use nlmod.grid.affine_transform_gdf instead." ) from .grid import affine_transform_gdf @@ -615,7 +615,7 @@ def affine_transform_gdf(gdf, affine): def get_extent(ds, rotated=True): """Get the model extent, corrected for angrot if necessary.""" logger.warning( - "nlmod.resample.get_extent is deprecated. Use nlmod.grid.get_extent instead" + "nlmod.resample.get_extent is deprecated. Use nlmod.grid.get_extent instead." ) from .grid import get_extent @@ -626,7 +626,7 @@ def get_affine_mod_to_world(ds): """Get the affine-transformation from model to real-world coordinates.""" logger.warning( "nlmod.resample.get_affine_mod_to_world is deprecated. " - "Use nlmod.grid.get_affine_mod_to_world instead" + "Use nlmod.grid.get_affine_mod_to_world instead." ) from .grid import get_affine_mod_to_world @@ -637,7 +637,7 @@ def get_affine_world_to_mod(ds): """Get the affine-transformation from real-world to model coordinates.""" logger.warning( "nlmod.resample.get_affine_world_to_mod is deprecated. " - "Use nlmod.grid.get_affine_world_to_mod instead" + "Use nlmod.grid.get_affine_world_to_mod instead." ) from .grid import get_affine_world_to_mod @@ -647,7 +647,7 @@ def get_affine_world_to_mod(ds): def get_affine(ds, sx=None, sy=None): """Get the affine-transformation, from pixel to real-world coordinates.""" logger.warning( - "nlmod.resample.get_affine is deprecated. Use nlmod.grid.get_affine instead" + "nlmod.resample.get_affine is deprecated. Use nlmod.grid.get_affine instead." ) from .grid import get_affine diff --git a/nlmod/read/geotop.py b/nlmod/read/geotop.py index 9457e53a..772851d2 100644 --- a/nlmod/read/geotop.py +++ b/nlmod/read/geotop.py @@ -41,8 +41,24 @@ def get_lithok_colors(): def get_strat_props(): - fname = os.path.join(NLMOD_DATADIR, "geotop", "geo_eenheden.csv") - df = pd.read_csv(fname, index_col=0, keep_default_na=False) + fname = os.path.join(NLMOD_DATADIR, "geotop", "REF_GTP_STR_UNIT.csv") + df = pd.read_csv(fname) + # rename the columns to previously used values + # so existing nlmod-code will keep working + df = df.rename( + columns={"STR_UNIT_CD": "code", "VOXEL_NR": "strat", "DESCRIPTION": "name"} + ) + # calculate color from red, green and blue columns + color = {} + for index in df.index: + color[index] = ( + df.at[index, "RED_DEC"] / 255, + df.at[index, "GREEN_DEC"] / 255, + df.at[index, "BLUE_DEC"] / 255, + ) + df["color"] = color + df = df.drop(columns=["RED_DEC", "GREEN_DEC", "BLUE_DEC"]).set_index("strat") + return df @@ -116,9 +132,13 @@ def to_model_layers( units = units[~np.isnan(units)].astype(int) shape = (len(units), len(geotop_ds.y), len(geotop_ds.x)) - # stratigraphy unit (geo eenheid) 2000 is above 1130 - if (2000 in units) and (1130 in units): - units[(units == 2000) + (units == 1130)] = [2000, 1130] + if "SEQ_NR" in strat_props.columns: + # sort units based on SEQ_NR in strat_props + units = strat_props.loc[units, "SEQ_NR"].sort_values().index.values + else: + # stratigraphy unit (geo eenheid) 2000 is above 1130 + if (2000 in units) and (1130 in units): + units[(units == 2000) + (units == 1130)] = [2000, 1130] # fill top and bot top = np.full(shape, np.nan) @@ -620,3 +640,21 @@ def aggregate_to_ds( ds[kh] = xr.concat(kh_ar, ds.layer) ds[kv] = xr.concat(kv_ar, ds.layer) return ds + + +def _save_excel_files_as_csv(): + """ + This method takes the files REF_GTP_STR_UNIT.xlsx and REF_GTP_LITHO_CLASS.xlsx that + are taken from the GeoTOP 1.6 zipfile downloaded from DINOloket, and saves them as + csv-files. In this way version-control can better process the changes in future + versions of GeoTOP. + + Returns + ------- + None. + + """ + for name in ["REF_GTP_STR_UNIT.xlsx", "REF_GTP_LITHO_CLASS.xlsx"]: + fname = os.path.join(NLMOD_DATADIR, "geotop", name) + df = pd.read_excel(fname, keep_default_na=False) + df.to_csv(fname.replace(".xlsx", ".csv"), index=False)