Merge pull request #268 from inbo/split_abiotic

Split abiotic: allow specifying only nutrient_level or acidity

LICENSE

@@ -1,6 +1,6 @@
 MIT License
 
-Copyright (c) 2017 Instituut voor Natuur en Bosonderzoek (INBO)
+Copyright (c) 2017-2020 Instituut voor Natuur en Bosonderzoek (INBO)
 
 Permission is hereby granted, free of charge, to any person obtaining a copy
 of this software and associated documentation files (the "Software"), to deal

docs/advanced_usage.ipynb

@@ -201,7 +201,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "Next we will create a new niche model using the same options as our previous full models, but we will also add the previously calculated acidity and nutrient level values as input, and run with the `abiotic=True` option. Note that we use the `read_config_file` method (and not  `run_config_file`) because we still want to edit the configuration before running."
+    "Next we will create a new niche model using the same options as our previous full models, but we will add the previously calculated acidity as input: in this case, this acidity grid will be used (and not calculated from the other grids). Note that we use the `read_config_file` method (and not  `run_config_file`) because we still want to edit the configuration before running."
    ]
   },
   {
@@ -212,10 +212,9 @@
    "source": [
     "adjusted = nv.Niche()\n",
     "adjusted.read_config_file(\"full.yml\")\n",
-    "adjusted.set_input(\"acidity\", \"output_abiotic/full_acidity.tif\")\n",
     "adjusted.set_input(\"nutrient_level\", \"output_abiotic/adjusted_nutrient.tif\")\n",
-    "adjusted.name = \"adjusted\"\n",
-    "adjusted.run(abiotic=True)\n",
+    "adjusted.name = \"adjusted_nutrient_level\"\n",
+    "adjusted.run()\n",
     "\n",
     "adjusted.plot(7)\n",
     "full.plot(7)\n",
@@ -283,9 +282,9 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Python 3",
+   "display_name": "niche",
    "language": "python",
-   "name": "python3"
+   "name": "niche"
   },
   "language_info": {
    "codemirror_mode": {
@@ -297,7 +296,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.5.3"
+   "version": "3.7.3"
   }
  },
  "nbformat": 4,

docs/cli.rst

@@ -69,8 +69,7 @@ Abiotic and/or deviation
 Using a configuration file, it is also possible to use abiotic values, like
 previously demonstrated in `Using abiotic grids`_.
 
-To do this the values ``acidity`` and ``nutrient_level`` must be specified,
-together with the ``abiotic`` model option.
+To do this the values ``acidity`` or ``nutrient_level`` must be specified.
 
 The option ``deviation`` creates deviation maps, which show the difference between
 the borders specified in the niche table and the actual values of mhw and mlw for

docs/conf.py

@@ -58,7 +58,7 @@
 
 # General information about the project.
 project = 'niche_vlaanderen'
-copyright = '2017 - 2018, Research Institute for Nature and Forest (INBO)'
+copyright = '2017 - 2020, Research Institute for Nature and Forest (INBO)'
 author = u'Johan Van de Wauw, Cécile Herr and Dries Adriaens'
 
 # The version info for the project you're documenting, acts as replacement for

docs/full.yml

@@ -1,5 +1,4 @@
 model_options:
-  abiotic: false
   deviation: false
   full_model: true
   name: 'full'

niche_vlaanderen/acidity.py

@@ -7,7 +7,10 @@
 
 
 class Acidity(object):
-    '''
+    ''' Class to calculate the Acidity
+
+     The used codetables can be overwritten by using the corresponding ct_*
+     arguments.
     '''
 
     nodata = 255  # uint8 data type
@@ -124,6 +127,29 @@ def _get_seepage(self, seepage):
 
     def calculate(self, soil_class, mlw, inundation, seepage, minerality,
                   rainwater):
+        """
+
+
+        Parameters:
+        ==========
+        soil_class: numpy.array
+            Array containing the soil codes. Values must be present
+            in the soil_code table. -99 is used as no data value.
+        mlw: numpy.array
+            Array containing the mean lowest waterlevel.
+        inundation:  numpy.array
+            Array containing rate of inundation.
+            https://inbo.github.io/niche_vlaanderen/invoer.html#overstroming-trofie-inundation-nutrient
+        seepage:  numpy.array
+            Array containing the flux of groundwater at the location
+            https://inbo.github.io/niche_vlaanderen/invoer.html#overstroming-trofie-inundation-nutrient
+        minerality: numpy.array
+            Array noting whether the groundwater is rich(1) or poor in minerals
+            https://inbo.github.io/niche_vlaanderen/invoer.html#mineraalrijkdom-minerality
+        rainwater: numpy.array (optional)
+            Array denoting whether rainwater lenses occur.
+            https://inbo.github.io/niche_vlaanderen/invoer.html#regenlens-rainwater
+        """
         soil_mlw = self._calculate_soil_mlw(soil_class, mlw)
         seepage = self._get_seepage(seepage)
         acidity = self._get_acidity(rainwater, minerality, inundation,

niche_vlaanderen/niche.py

@@ -34,30 +34,23 @@
     "rainwater", "inundation_vegetation", "management_vegetation", "acidity",
     "nutrient_level"}
 
-_minimal_input = {
-    "soil_code", "mlw", "msw", "mhw", "seepage", "inundation_acidity",
-    "nitrogen_atmospheric", "nitrogen_animal", "nitrogen_fertilizer",
-    "management", "minerality", "rainwater",
-    "inundation_nutrient"}
-
-
-def minimal_input():  # pragma: no cover
-    """Get a set with all required input codes for a full Niche model
-    Returns
-    =======
-       :set
-       set of input_codes required for a full Niche model
-    """
-    return _minimal_input
+_minimal_input = {"mlw", "mhw", "soil_code"}
+
+_input_nutrient_level = {"msw", "soil_code", "nitrogen_atmospheric",
+                         "nitrogen_fertilizer", "nitrogen_animal",
+                         "management", "inundation_nutrient"}
 
+_input_acidity = {"soil_code", "mlw", "minerality", "seepage", "rainwater",
+                  "inundation_acidity"}
 
 _abiotic_keys = {"nutrient_level", "acidity"}
 
-_code_tables = ["ct_acidity", "ct_soil_mlw_class", "ct_soil_codes",
+_code_tables = {"ct_acidity", "ct_soil_mlw_class", "ct_soil_codes",
                 "lnk_acidity", "ct_seepage", "ct_vegetation", "ct_management",
-                "ct_nutrient_level", "ct_mineralisation"]
-_code_tables_fp = ["duration",
-                   "frequency", "lnk_potential", "potential"]
+                "ct_nutrient_level", "ct_mineralisation"}
+
+_code_tables_fp = {"duration",
+                   "frequency", "lnk_potential", "potential"}
 
 logging.basicConfig()
 
@@ -266,9 +259,7 @@ def read_config_file(self, config, overwrite_ct=False):
             if "output_dir" in config_loaded["model_options"].keys():
                 self._options["output_dir"] = \
                     config_loaded["model_options"]["output_dir"]
-            if "abiotic" in config_loaded["model_options"].keys():
-                self._options["abiotic"] = \
-                    config_loaded["model_options"]["abiotic"]
+
 
         if "flooding" in config_loaded.keys():
             self._options["flooding"] = []
@@ -413,10 +404,11 @@ def _check_input_files(self, full_model):
 
         self._check_all_lower(inputarray, "mhw", "mlw")
 
-        if full_model:
+        if 'msw' in inputarray.keys():
             self._check_all_lower(inputarray, "msw", "mlw")
             self._check_all_lower(inputarray, "mhw", "msw")
 
+        if full_model and not 'nutrient_level' in inputarray.keys():
             with np.errstate(invalid='ignore'):  # ignore NaN comparison errors
                 if np.any((inputarray['nitrogen_animal'] < 0)
                           | (inputarray['nitrogen_animal'] > 10000)
@@ -430,8 +422,7 @@ def _check_input_files(self, full_model):
         # if all is successful:
         self._inputarray = inputarray
 
-    def run(self, full_model=True, deviation=False, abiotic=False,
-            strict_checks=True):
+    def run(self, full_model=True, deviation=False, strict_checks=True):
         """Run the niche model
 
         Runs niche Vlaanderen model. Requires that the necessary input values
@@ -447,9 +438,6 @@ def run(self, full_model=True, deviation=False, abiotic=False,
                     Create the maps with the difference between the needed MHW
                     and MLW and the actual MHW/MLW for a vegetation type.
 
-        abiotic:  bool
-                Specify the abiotic grids as input rather than calculating
-                them. See tutorial at: `Using abiotic grids`_
         strict_checks: bool
                 By default running a model will fail if impossible combinations
                 of MxW exist somewhere in the input files. By disabling strict
@@ -460,29 +448,21 @@ def run(self, full_model=True, deviation=False, abiotic=False,
 
         self._options["full_model"] = full_model
         self._options["deviation"] = deviation
-        self._options["abiotic"] = abiotic
         self._options["strict_checks"] = strict_checks
 
-        if abiotic:
-            missing_keys = (_abiotic_keys
-                            - set(self._inputfiles.keys())
-                            - set(self._inputvalues.keys()))
-            if len(missing_keys) > 0:
-                print("Abiotic input are missing: (abiotic=True)")
-                print(missing_keys)
-                raise NicheException(
-                    "Error, abiotic keys are missing")
+        if full_model:
+            required_input = set(_minimal_input)
 
-        if not abiotic and (
-                (_abiotic_keys & set(self._inputfiles.keys()))
-                or (_abiotic_keys & set(self._inputvalues.keys()))):
-            warnings.warn(
-                "abiotic inputs specified but not specified in model options\n"
-                "abiotic inputs will not be used")
+            given_input = set(self._inputfiles.keys()) | set(self._inputvalues.keys())
+
+            if 'nutrient_level' not in given_input:
+                required_input |= set(_input_nutrient_level)
+
+            if 'acidity' not in given_input:
+                required_input |= set(_input_acidity)
+
+            missing_keys = required_input - given_input
 
-        if full_model:
-            missing_keys = _minimal_input - set(self._inputfiles.keys()) \
-                           - set(self._inputvalues.keys())
             if len(missing_keys) > 0:
                 print("Different keys are missing: ")
                 print(missing_keys)
@@ -491,41 +471,43 @@ def run(self, full_model=True, deviation=False, abiotic=False,
 
         self._check_input_files(full_model)
 
-        if full_model and not abiotic:
-            ct_nl = dict()
+        if full_model:
+            if 'nutrient_level' not in self._inputarray:
+                ct_nl = dict()
 
-            keys = set(NutrientLevel.__init__.__code__.co_varnames) \
-                & set(self._code_tables)
+                keys = set(NutrientLevel.__init__.__code__.co_varnames) \
+                    & set(self._code_tables)
 
-            for k in keys:
-                ct_nl[k] = self._code_tables[k]
+                for k in keys:
+                    ct_nl[k] = self._code_tables[k]
 
-            nl = NutrientLevel(**ct_nl)
+                nl = NutrientLevel(**ct_nl)
 
-            self._abiotic["nutrient_level"] = nl.calculate(
-                soil_code=self._inputarray["soil_code"],
-                msw=self._inputarray["msw"],
-                nitrogen_atmospheric=self._inputarray["nitrogen_atmospheric"],
-                nitrogen_animal=self._inputarray["nitrogen_animal"],
-                nitrogen_fertilizer=self._inputarray["nitrogen_fertilizer"],
-                management=self._inputarray["management"],
-                inundation=self._inputarray["inundation_nutrient"])
+                self._abiotic["nutrient_level"] = nl.calculate(
+                    soil_code=self._inputarray["soil_code"],
+                    msw=self._inputarray["msw"],
+                    nitrogen_atmospheric=self._inputarray["nitrogen_atmospheric"],
+                    nitrogen_animal=self._inputarray["nitrogen_animal"],
+                    nitrogen_fertilizer=self._inputarray["nitrogen_fertilizer"],
+                    management=self._inputarray["management"],
+                    inundation=self._inputarray["inundation_nutrient"])
 
-            ct_acidity = dict()
+            if  'acidity' not in self._inputarray:
+                ct_acidity = dict()
 
-            keys = set(Acidity.__init__.__code__.co_varnames) \
-                & set(self._code_tables)
+                keys = set(Acidity.__init__.__code__.co_varnames) \
+                    & set(self._code_tables)
 
-            for k in keys:
-                ct_acidity[k] = self._code_tables[k]
+                for k in keys:
+                    ct_acidity[k] = self._code_tables[k]
 
-            acidity = Acidity(**ct_acidity)
-            self._abiotic["acidity"] = acidity.calculate(
-                self._inputarray["soil_code"], self._inputarray["mlw"],
-                self._inputarray["inundation_acidity"],
-                self._inputarray["seepage"],
-                self._inputarray["minerality"],
-                self._inputarray["rainwater"])
+                acidity = Acidity(**ct_acidity)
+                self._abiotic["acidity"] = acidity.calculate(
+                    self._inputarray["soil_code"], self._inputarray["mlw"],
+                    self._inputarray["inundation_acidity"],
+                    self._inputarray["seepage"],
+                    self._inputarray["minerality"],
+                    self._inputarray["rainwater"])
 
         ct_veg = dict()
 
@@ -551,14 +533,18 @@ def run(self, full_model=True, deviation=False, abiotic=False,
                 inundation=self._inputarray["inundation_vegetation"],
                 management=self._inputarray["management_vegetation"]
             )
-            if not abiotic:
-                veg_arguments.update(
-                    nutrient_level=self._abiotic["nutrient_level"],
-                    acidity=self._abiotic["acidity"])
+
+            if 'nutrient_level' in self._inputarray:
+                veg_arguments['nutrient_level'] = self._inputarray['nutrient_level']
+            else:
+                veg_arguments['nutrient_level'] = self._abiotic[
+                    'nutrient_level']
+
+            if 'acidity' in self._inputarray:
+                veg_arguments['acidity'] = self._inputarray['acidity']
             else:
-                veg_arguments.update(
-                    nutrient_level=self._inputarray["nutrient_level"],
-                    acidity=self._inputarray["acidity"])
+                veg_arguments['acidity'] = self._abiotic[
+                    'acidity']
 
         self._vegetation, self.occurrence = vegetation.calculate(
             full_model=full_model, **veg_arguments)
@@ -575,7 +561,7 @@ def write(self, folder, overwrite_files=False):
         Saves the model results to a folder. Files will be written as geotiff.
         Vegetation files have names V1 ... V28
         Abiotic files are exported as well (nutrient_level.tif and
-        acidity.tif).
+        acidity.tif) if they were not input files.
 
         Parameters
         ----------

tests/small_abiotic.yaml

@@ -4,20 +4,16 @@ input_layers:
   mlw: data/small/mlw.asc
   soil_code: data/small/soil_code.asc
   seepage: 0
-  inundation_acidity: 0
   inundation_nutrient: 0
   nitrogen_atmospheric: 30
   nitrogen_animal: 0
   nitrogen_fertilizer: 0
   management: 1
-  minerality: 1
-  rainwater: 0
   acidity: data/small/acidity.asc
   nutrient_level: data/small/nutrient_level.asc
 
 model_options:
   simple_model: False
-  abiotic: True
   deviation: True
   output_dir: _output
   overwrite_files: True

tests/small_abiotic_extra.yaml

@@ -0,0 +1,12 @@
+input_layers:
+  mhw: data/small/mhw.asc
+  mlw: data/small/mlw.asc
+  soil_code: data/small/soil_code.asc
+  acidity: data/small/acidity.asc
+  nutrient_level: data/small/nutrient_level.asc
+
+model_options:
+  simple_model: False
+  deviation: True
+  output_dir: _output
+  overwrite_files: True