NDWI index added (#49)

* NDWI index added

NDWI (Normalized Difference Water Index) after McFeeters (1996), https://doi.org/10.3390/rs5073544

* NDWI added to testsuite/test_vi.py

* NDWI example: set color table

imagery/i.vi/i.vi.html

@@ -15,6 +15,7 @@ parameters.
   <li>MSAVI2: second Modified Soil Adjusted Vegetation Index</li>
   <li>MSAVI: Modified Soil Adjusted Vegetation Index</li>
   <li>NDVI: Normalized Difference Vegetation Index</li>
+  <li>NDWI: Normalized Difference Water Index</li>
   <li>PVI: Perpendicular Vegetation Index</li>
   <li>RVI: ratio vegetation index</li>
   <li>SAVI: Soil Adjusted Vegetation Index</li>
@@ -216,6 +217,17 @@ NDVI = (NIR - Red) / (NIR + Red)
 </pre></div>
 
 <p>
+<b>NDWI: Normalized Difference Water Index</b>
+
+(after McFeeters, 1996)
+
+<div class="code"><pre>
+ndwi( greenchan, nirchan )
+
+NDWI = (green - NIR) / (green + NIR)
+</pre></div>
+
+<p>
 <b>PVI: Perpendicular Vegetation Index</b>
 
 <div class="code"><pre>
@@ -380,6 +392,16 @@ i.vi red=band.3 nir=band.4 viname=ndvi output=ndvi
 r.univar -e ndvi
 </pre></div>
 
+<h3>Calculation of NDWI</h3>
+
+The calculation of NDWI from the reflectance values is done as follows:
+
+<div class="code"><pre>
+g.region raster=band.2 -p
+i.vi green=band.2 nir=band.4 viname=ndwi output=ndwi
+r.colors ndwi color=byg -n
+r.univar -e ndwi
+</pre></div>
 
 <h3>Calculation of PVI</h3>
 

imagery/i.vi/main.c

@@ -1,29 +1,30 @@
 
 /****************************************************************************
  *
- * MODULE:       i.vi
- * AUTHOR(S):    Baburao Kamble baburaokamble@gmail.com
- *		 Yann Chemin - yann.chemin@gmail.com
- *		 Nikos Alexandris - nik@nikosalexandris.net
- * PURPOSE:      Calculates 15 vegetation indices
- * 		 based on biophysical parameters.
+ * MODULE:      i.vi
+ * AUTHOR(S):   Baburao Kamble baburaokamble@gmail.com
+ *              Yann Chemin - yann.chemin@gmail.com
+ *              Nikos Alexandris - nik@nikosalexandris.net
+ * PURPOSE:     Calculates 16 vegetation and related indices
+ *              based on biophysical parameters.
  *
- * COPYRIGHT:    (C) 2002-2013 by the GRASS Development Team
+ * COPYRIGHT:    (C) 2002-2019 by the GRASS Development Team
  *
  *               This program is free software under the GNU General Public
  *   	    	 License (>=v2). Read the file COPYING that comes with GRASS
  *   	    	 for details.
  *
- * Remark:
- *		 These are generic indices that use red and nir for most of them.
- *               Those can be any use by standard satellite having V and IR.
- *		 However arvi uses red, nir and blue;
- *		 GVI uses B,G,R,NIR, chan5 and chan 7 of landsat;
- *		 and GARI uses B,G,R and NIR.
+ * Remarks:
+ *           These are generic indices that use red and nir for most of them.
+ *           Those can be any use by standard satellite having V and IR.
+ *           However, arvi uses red, nir and blue;
+ *           GVI uses B,G,R,NIR, chan5 and chan 7 of landsat;
+ *           and GARI uses B,G,R and NIR.
  *
- * Changelog:	 Added EVI on 20080718 (Yann)
- * 		 Added VARI on 20081014 (Yann)
- * 		 Added EVI2 on 20130208 (NikosA)
+ * Changelog:   Added EVI on 20080718 (Yann)
+ *              Added VARI on 20081014 (Yann)
+ *              Added EVI2 on 20130208 (NikosA)
+ *              Added NDWI on 20190711 (neteler)
  *
  *****************************************************************************/
 
@@ -37,6 +38,7 @@
 
 double s_r(double redchan, double nirchan);
 double nd_vi(double redchan, double nirchan);
+double nd_wi(double greenchan, double nirchan);
 double ip_vi(double redchan, double nirchan);
 double d_vi(double redchan, double nirchan);
 double e_vi(double bluechan, double redchan, double nirchan);
@@ -112,7 +114,7 @@ int main(int argc, char *argv[])
     desc = NULL;
     G_asprintf(&desc,
 	       "arvi;%s;dvi;%s;evi;%s;evi2;%s;gvi;%s;gari;%s;gemi;%s;ipvi;%s;msavi;%s;"
-	       "msavi2;%s;ndvi;%s;pvi;%s;savi;%s;sr;%s;vari;%s;wdvi;%s",
+	       "msavi2;%s;ndvi;%s;ndwi;%s;pvi;%s;savi;%s;sr;%s;vari;%s;wdvi;%s",
 	       _("Atmospherically Resistant Vegetation Index"),
 	       _("Difference Vegetation Index"),
 	       _("Enhanced Vegetation Index"),
@@ -124,13 +126,14 @@ int main(int argc, char *argv[])
 	       _("Modified Soil Adjusted Vegetation Index"),
 	       _("second Modified Soil Adjusted Vegetation Index"),
 	       _("Normalized Difference Vegetation Index"),
+	       _("Normalized Difference Water Index"),
 	       _("Perpendicular Vegetation Index"),
 	       _("Soil Adjusted Vegetation Index"),
 	       _("Simple Ratio"),
 	       _("Visible Atmospherically Resistant Index"),
 	       _("Weighted Difference Vegetation Index"));
     opt.viname->descriptions = desc;
-    opt.viname->options = "arvi,dvi,evi,evi2,gvi,gari,gemi,ipvi,msavi,msavi2,ndvi,pvi,savi,sr,vari,wdvi";
+    opt.viname->options = "arvi,dvi,evi,evi2,gvi,gari,gemi,ipvi,msavi,msavi2,ndvi,ndwi,pvi,savi,sr,vari,wdvi";
     opt.viname->answer = "ndvi";
     opt.viname->key_desc = _("type");
 
@@ -231,6 +234,9 @@ int main(int argc, char *argv[])
     if (!strcasecmp(viflag, "ndvi") && (!(opt.red->answer) || !(opt.nir->answer)) )
 	G_fatal_error(_("ndvi index requires red and nir maps"));
 
+    if (!strcasecmp(viflag, "ndwi") && (!(opt.green->answer) || !(opt.nir->answer)) )
+	G_fatal_error(_("ndwi index requires green and nir maps"));
+
     if (!strcasecmp(viflag, "ipvi") && (!(opt.red->answer) || !(opt.nir->answer)) )
 	G_fatal_error(_("ipvi index requires red and nir maps"));
 
@@ -458,12 +464,15 @@ int main(int argc, char *argv[])
 
 		/* calculate ndvi                    */
 		if (!strcasecmp(viflag, "ndvi")) {
-		    if (d_redchan + d_nirchan < 0.001)
+		    if (d_redchan + d_nirchan < 0.001)  /* TODO: why this? */
 			Rast_set_d_null_value(&outrast[col], 1);
 		    else
 			outrast[col] = nd_vi(d_redchan, d_nirchan);
 		}
 
+		if (!strcasecmp(viflag, "ndwi"))
+		    outrast[col] = nd_wi(d_greenchan, d_nirchan);
+
 		if (!strcasecmp(viflag, "ipvi"))
 		    outrast[col] = ip_vi(d_redchan, d_nirchan);
 

imagery/i.vi/ndwi.c

@@ -0,0 +1,20 @@
+#include<stdio.h>
+#include<math.h>
+#include<stdlib.h>
+
+    /* Normalized Difference Water Index
+     * after McFeeters (1996), https://doi.org/10.3390/rs5073544 */
+double nd_wi(double greenchan, double nirchan)
+{
+    double result;
+
+    if ((greenchan + nirchan) == 0.0) {
+	result = -1.0; /* TODO: -1 or 0 */
+    }
+    else {
+	result = (greenchan - nirchan) / (greenchan + nirchan);
+    }
+    return result;
+}
+
+

imagery/i.vi/testsuite/test_vi.py

@@ -26,6 +26,7 @@ class TestReport(TestCase):
     @classmethod
     def tearDownClass(cls):
         cls.runModule("g.remove", flags='f', type="raster", name="ipvi")
+        cls.runModule("g.remove", flags='f', type="raster", name="ndwi")
         cls.runModule("g.remove", flags='f', type="raster", name="dvi")
         cls.runModule("g.remove", flags='f', type="raster", name="sr")
         cls.runModule("g.remove", flags='f', type="raster", name="evi")
@@ -43,6 +44,14 @@ class TestReport(TestCase):
                                 refmax=0.906666666667,
                                 msg="ipvi in degrees must be between 0.0454545454545  and 0.906666666667")
 
+    def test_vinamendwi(self):
+        """Testing viname ndwi"""
+        map_output = 'ndwi'
+        self.assertModule('i.vi', red=self.red, green=self.green, nir=self.nir, viname='ndwi',
+                          output=map_output)
+        self.assertRasterMinMax(map=map_output, refmin=-0.71, refmax=0.93,
+                                msg="ndwi in percent must be between -0.71 and 0.93")
+
     def test_vinamedvi(self):
         """Testing viname dvi"""
         map_output = 'dvi'