i.albedo: add Sentinel-2 broadband albedo support

imagery/i.albedo/bb_alb_sentinel2.c

@@ -0,0 +1,19 @@
+/* Broadband albedo Sentinel-2 MSI
+ *
+ * Narrow-to-broadband conversion coefficients from:
+ * Bonafoni, S. and Sekertekin, A. (2020). Albedo Retrieval From Sentinel-2
+ * by New Narrow-to-Broadband Conversion Coefficients. IEEE Geoscience and
+ * Remote Sensing Letters, 17(9), 1618-1622.
+ *
+ * Input bands are surface reflectance in [0, 1], in Sentinel-2 band order
+ * B2 (blue), B3 (green), B4 (red), B8 (NIR), B11 (SWIR1), B12 (SWIR2).
+ */
+double bb_alb_sentinel2(double b2chan, double b3chan, double b4chan,
+                        double b8chan, double b11chan, double b12chan)
+{
+    double result;
+
+    result = (0.2266 * b2chan + 0.1236 * b3chan + 0.1573 * b4chan +
+              0.3417 * b8chan + 0.1170 * b11chan + 0.0338 * b12chan);
+    return result;
+}

imagery/i.albedo/functions.h

@@ -8,6 +8,8 @@ double bb_alb_aster(double greenchan, double redchan, double nirchan,
                     double swirchan4, double swirchan5, double swirchan6);
 double bb_alb_landsat(double bluechan, double greenchan, double redchan,
                       double nirchan, double chan5, double chan7);
+double bb_alb_sentinel2(double b2chan, double b3chan, double b4chan,
+                        double b8chan, double b11chan, double b12chan);
 double bb_alb_noaa(double redchan, double nirchan);
 
 double bb_alb_modis(double redchan, double nirchan, double chan3, double chan4,

imagery/i.albedo/i.albedo.md

@@ -3,8 +3,9 @@
 *i.albedo* calculates the albedo, that is the Shortwave surface
 reflectance in the range of 0.3-3 micro-meters. It takes as input
 individual bands of surface reflectance originating from MODIS, AVHRR,
-Landsat or Aster satellite sensors and calculates the albedo for those.
-This is a precursor to *r.sun* and any energy-balance processing.
+Landsat, Sentinel-2 or Aster satellite sensors and calculates the albedo
+for those. This is a precursor to *r.sun* and any energy-balance
+processing.
 
 ## NOTES
 
@@ -13,6 +14,11 @@ solution until an algorithm is found).
 
 It assumes MODIS product surface reflectance in \[0;10000\].
 
+The Sentinel-2 (**-s**) mode expects exactly 6 surface reflectance bands
+in \[0,1\], in the order B2, B3, B4, B8, B11, B12 (e.g. as produced by
+Sentinel-2 L2A surface reflectance, scaled from its native \[0,10000\]
+digital numbers).
+
 ## EXAMPLE
 
 The following example creates the raster map "albedo_lsat7_1987" from
@@ -37,6 +43,14 @@ i.albedo -l input=lsat7_2000_10,lsat7_2000_20,lsat7_2000_30,lsat7_2000_40,lsat7_
 ![i.albedo LANDSAT-TM7 2000 example](i_albedo_landsat7.png)  
 *Figure: Resulting albedo map from LANDSAT 7*
 
+The following example creates a Sentinel-2 broadband albedo map from six
+surface reflectance bands (B2, B3, B4, B8, B11, B12) already scaled to
+\[0,1\]:
+
+```sh
+i.albedo -s input=s2_B2,s2_B3,s2_B4,s2_B8,s2_B11,s2_B12 output=albedo_sentinel2
+```
+
 ## TODO
 
 Maybe change input requirement of MODIS to \[0.0-1.0\]?
@@ -49,6 +63,13 @@ Salleh and Chan, 2014. Land Surface Albedo Determination: Remote Sensing
 and Statistical Validation. in proceedings of FIG 2014
 ([PDF](https://www.fig.net/resources/proceedings/fig_proceedings/fig2014/papers/ts05g/TS05G_salleh_chan_6910.pdf))
 
+For the Sentinel-2 narrow-to-broadband conversion coefficients:
+
+Bonafoni, S. and Sekertekin, A., 2020. Albedo Retrieval From Sentinel-2 by
+New Narrow-to-Broadband Conversion Coefficients. *IEEE Geoscience and
+Remote Sensing Letters*, 17(9), 1618-1622.
+([DOI](https://doi.org/10.1109/LGRS.2020.2967085))
+
 ## SEE ALSO
 
 *[r.sun](r.sun.md), [i.vi](i.vi.md)*

imagery/i.albedo/main.c

@@ -29,6 +29,9 @@ double bb_alb_landsat8(double shortbluechan, double bluechan, double greenchan,
                        double redchan, double nirchan, double chan5,
                        double chan7);
 
+double bb_alb_sentinel2(double b2chan, double b3chan, double b4chan,
+                        double b8chan, double b11chan, double b12chan);
+
 double bb_alb_noaa(double redchan, double nirchan);
 
 double bb_alb_modis(double redchan, double nirchan, double chan3, double chan4,
@@ -43,7 +46,7 @@ int main(int argc, char *argv[])
     struct Option *input, *output;
     struct Flag *flag1, *flag2, *flag3;
     struct Flag *flag4, *flag5, *flag6;
-    struct Flag *flag7;
+    struct Flag *flag7, *flag8;
     struct History history; /*metadata */
     struct Colors colors;   /*Color rules */
 
@@ -58,7 +61,7 @@ int main(int argc, char *argv[])
     char **ptr;
     int i = 0;
     int modis = 0, aster = 0, avhrr = 0;
-    int landsat = 0, landsat8 = 0;
+    int landsat = 0, landsat8 = 0, sentinel2 = 0;
     void *inrast[MAXFILES];
     unsigned char *outrast;
 
@@ -94,6 +97,7 @@ int main(int argc, char *argv[])
     G_add_keyword(_("reflectance"));
     G_add_keyword(_("satellite"));
     G_add_keyword(_("Landsat"));
+    G_add_keyword(_("Sentinel"));
     G_add_keyword(_("ASTER"));
     G_add_keyword(_("AVHRR"));
     G_add_keyword(_("MODIS"));
@@ -125,6 +129,12 @@ int main(int argc, char *argv[])
     flag4->key = '8';
     flag4->description = _("Landsat 8 (7 input bands:1,2,3,4,5,6,7)");
 
+    flag8 = G_define_flag();
+    flag8->key = 's';
+    flag8->description =
+        _("Sentinel-2 (6 input bands: B2,B3,B4,B8,B11,B12 surface "
+          "reflectance in [0,1]; Bonafoni & Sekertekin 2020 coefficients)");
+
     flag5 = G_define_flag();
     flag5->key = 'a';
     flag5->description = _("ASTER (6 input bands:1,3,5,6,8,9)");
@@ -158,6 +168,7 @@ int main(int argc, char *argv[])
     landsat = (flag3->answer);
     landsat8 = (flag4->answer);
     aster = (flag5->answer);
+    sentinel2 = (flag8->answer);
 
     for (; *ptr != NULL; ptr++) {
         if (nfiles >= MAXFILES)
@@ -233,6 +244,9 @@ int main(int argc, char *argv[])
                 else if (aster) {
                     de = bb_alb_aster(d[1], d[2], d[3], d[4], d[5], d[6]);
                 }
+                else if (sentinel2) {
+                    de = bb_alb_sentinel2(d[1], d[2], d[3], d[4], d[5], d[6]);
+                }
                 if (Rast_is_d_null_value(&de)) {
                     /*Do nothing */
                 }
@@ -390,6 +404,9 @@ int main(int argc, char *argv[])
             else if (aster) {
                 de = bb_alb_aster(d[1], d[2], d[3], d[4], d[5], d[6]);
             }
+            else if (sentinel2) {
+                de = bb_alb_sentinel2(d[1], d[2], d[3], d[4], d[5], d[6]);
+            }
             if (flag6->answer || flag7->answer) {
                 /* Post-Process Albedo */
                 de = a * de + b;