Repository: inbo/niche_vlaanderen

Issue:
Abiotic: use custom acidity OR nutrient_level
* niche_vlaanderen version: 1.0
* Python version: 3
* Operating System: Win 

## Niche Python Package question

### Description
If you use  mymodel.run( abiotic= True), you need to provide both an input layer for acidity and nutrient_level. 

It happens regularly that we need only one custom layer (in my case: most often a custom layer for acidity while I use the standard nutrient_level), so I have to run the model twice (first run to produce the nutrient_level, and another with the custom acidity and standard nutrient_level) . 
It would be very nice to be able to use only a custom acidity or only a custom nutrient_level without having to run the model twice.




Discussion:

Patch:
diff --git a/LICENSE b/LICENSE
index 03bf2b2..c70f54f 100644
--- a/LICENSE
+++ b/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
diff --git a/docs/advanced_usage.ipynb b/docs/advanced_usage.ipynb
index 4a56af2..f68eb23 100644
--- a/docs/advanced_usage.ipynb
+++ b/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,
diff --git a/docs/cli.rst b/docs/cli.rst
index 308e0fa..e39b2b8 100644
--- a/docs/cli.rst
+++ b/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
diff --git a/docs/conf.py b/docs/conf.py
index cfdef3e..63f286d 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -58,7 +58,7 @@ master_doc = 'index'
 
 # 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
diff --git a/docs/full.yml b/docs/full.yml
index 21a5a80..9c5afd9 100644
--- a/docs/full.yml
+++ b/docs/full.yml
@@ -1,5 +1,4 @@
 model_options:
-  abiotic: false
   deviation: false
   full_model: true
   name: 'full'
diff --git a/niche_vlaanderen/acidity.py b/niche_vlaanderen/acidity.py
index 84d07d3..047c116 100644
--- a/niche_vlaanderen/acidity.py
+++ b/niche_vlaanderen/acidity.py
@@ -7,7 +7,10 @@ from .codetables import validate_tables_acidity, check_codes_used
 
 
 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 @@ class Acidity(object):
 
     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,
diff --git a/niche_vlaanderen/niche.py b/niche_vlaanderen/niche.py
index 657f4a7..e8bd7af 100644
--- a/niche_vlaanderen/niche.py
+++ b/niche_vlaanderen/niche.py
@@ -34,30 +34,23 @@ _allowed_input = {
     "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 @@ class Niche(object):
             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 @@ class Niche(object):
 
         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 @@ class Niche(object):
         # 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 @@ class Niche(object):
                     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 @@ class Niche(object):
 
         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 @@ class Niche(object):
 
         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 @@ class Niche(object):
                 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 @@ class Niche(object):
         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
         ----------


Test patch:
diff --git a/tests/small_abiotic.yaml b/tests/small_abiotic.yaml
index d58d3ec..4af98e0 100644
--- a/tests/small_abiotic.yaml
+++ b/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
diff --git a/tests/small_abiotic_extra.yaml b/tests/small_abiotic_extra.yaml
new file mode 100644
index 0000000..693d206
--- /dev/null
+++ b/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
diff --git a/tests/small_abiotic_invalid.yaml b/tests/small_abiotic_invalid.yaml
new file mode 100644
index 0000000..6baca3b
--- /dev/null
+++ b/tests/small_abiotic_invalid.yaml
@@ -0,0 +1,18 @@
+input_layers:
+  mhw: data/small/mhw.asc
+  msw: data/small/msw.asc
+  mlw: data/small/mlw.asc
+  soil_code: data/small/soil_code.asc
+  seepage: 0
+  inundation_nutrient: 0
+  nitrogen_animal: 0
+  nitrogen_fertilizer: 0
+  management: 1
+  acidity: data/small/acidity.asc
+  # invalid as nutrient_level is not defined and nitrogen_environment is missing
+
+model_options:
+  simple_model: False
+  deviation: True
+  output_dir: _output
+  overwrite_files: True
diff --git a/tests/test_niche.py b/tests/test_niche.py
index 6856562..79a807c 100644
--- a/tests/test_niche.py
+++ b/tests/test_niche.py
@@ -103,12 +103,31 @@ class TestNiche(TestCase):
 
         dir = os.listdir(tmpdir)
 
+        if sys.version_info < (3, 2):
+            self.assertItemsEqual(expected_files, dir)
+        else:
+            self.assertCountEqual(expected_files, dir)
+
+        # run check after predefined acidity
+        myniche = self.create_zwarte_beek_niche()
+        myniche.set_input("acidity", tmpdir+'/acidity.tif')
+        myniche.run()
+
+        tmpdir_acidity = tempfile.mkdtemp()
+        myniche.write(tmpdir_acidity)
+
+        dir = os.listdir(tmpdir_acidity)
+
+        expected_files.remove('acidity.tif')
+
         if sys.version_info < (3, 2):
             self.assertItemsEqual(expected_files, dir)
         else:
             self.assertCountEqual(expected_files, dir)
 
         shutil.rmtree(tmpdir)
+        shutil.rmtree(tmpdir_acidity)
+
 
     def test_zwarte_beek_constant_values(self):
         myniche = self.create_zwarte_beek_niche()
@@ -139,7 +158,7 @@ class TestNiche(TestCase):
         myniche.run(full_model=False)
         tmpdir = tempfile.mkdtemp()
         myniche.write(tmpdir)
-        # check tempdir contains the vegation and the abiotic files
+        # check tempdir contains the vegetation files
         expected_files = [
              'V01.tif', 'V02.tif', 'V03.tif', 'V04.tif', 'V05.tif', 'V06.tif',
              'V07.tif', 'V08.tif', 'V09.tif', 'V10.tif', 'V11.tif', 'V12.tif',
@@ -294,20 +313,17 @@ class TestNiche(TestCase):
         myniche = niche_vlaanderen.Niche()
         myniche.run_config_file(config)
 
-    def test_run_abiotic_error(self):
-        myniche = self.create_small()
-        with pytest.raises(NicheException):
-            myniche.run(abiotic=True)
+        config = 'tests/small_abiotic_extra.yaml'
+        myniche = niche_vlaanderen.Niche()
+        myniche.run_config_file(config)
 
-        myniche.set_input("nutrient_level", 1)
-        myniche.set_input("acidity", 1)
+        # rerun with a file with missing abiotic values
+        config = 'tests/small_abiotic_invalid.yaml'
+        myniche = niche_vlaanderen.Niche()
+        with pytest.raises(NicheException):
+            myniche.run_config_file(config)
 
-        myniche.run(abiotic=True)
 
-        # The next run should raise a warning as abiotic grids are added but
-        # not used
-        with pytest.warns(UserWarning):
-            myniche.run(abiotic=False)
 
     def test_rereadoutput(self):
         """