i.oif: now it can work with several bands

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@59350 15284696-431f-4ddb-bdfa-cd5b030d7da7

scripts/i.oif/i.oif.html

@@ -17,7 +17,7 @@ required for <em>d.rgb</em> or <em>r.composite</em>.
 
 <h2>NOTES</h2>
 
-Landsat TM:
+Landsat 1-7 TM:
 Colour Composites in BGR order as important Landsat TM band combinations
 (example: 234 in BGR order means: B=2, G=3, R=4):
 
@@ -55,8 +55,7 @@ limited to that number of jobs.
 North Carolina sample dataset:
 <div class="code"><pre>
 g.region rast=lsat7_2002_10 -p
-i.oif image1=lsat7_2002_10 image2=lsat7_2002_20 image3=lsat7_2002_30 \
-      image4=lsat7_2002_40 image5=lsat7_2002_50 image7=lsat7_2002_70
+i.oif input=lsat7_2002_10,lsat7_2002_20,lsat7_2002_30,lsat7_2002_40,lsat7_2002_50,lsat7_2002_70
 </pre></div>
 
 

scripts/i.oif/i.oif.py

@@ -3,8 +3,11 @@
 ############################################################################
 #
 # MODULE:       i.oif
-# AUTHOR(S):    Markus Neteler
+# AUTHOR(S):    Markus Neteler 21.July 1998
+#               Updated for GRASS 5.7 by Michael Barton 2004/04/05
 #               Converted to Python by Glynn Clements
+#               Customised by Nikos Alexandris 04. August 2013
+#               Customised by Luca Delucchi Vienna Code Sprint 2014
 # PURPOSE:      calculates the Optimum Index factor of all band combinations
 #               for LANDSAT TM 1,2,3,4,5,7
 # COPYRIGHT:    (C) 1999,2008 by the GRASS Development Team
@@ -13,45 +16,20 @@
 #               License (>=v2). Read the file COPYING that comes with GRASS
 #               for details.
 #
-#############################################################################
 # Ref.: Jensen: Introductory digital image processing 1996, p.98
-#
-# Input: tm1 - tm5, tm7 (not tm6)
-#
-# written by Markus Neteler 21.July 1998
-#            neteler geog.uni-hannover.de
-# updated for GRASS 5.7 by Michael Barton 2004/04/05
-
+#############################################################################
 
 #% Module
-#% description: Calculates Optimum-Index-Factor table for LANDSAT TM bands 1-5, & 7
+#% description: Calculates Optimum-Index-Factor table for spectral bands
 #% keywords: imagery
-#% keywords: landsat
+#% keywords: multispectral
 #% keywords: statistics
 #% End
-#% option G_OPT_R_INPUT
-#% key: image1
-#% description: LANDSAT TM band 1.
+#% option G_OPT_R_INPUTS
 #% end
-#% option G_OPT_R_INPUT
-#% key: image2
-#% description: LANDSAT TM band 2.
-#% end
-#% option G_OPT_R_INPUT
-#% key: image3
-#% description: LANDSAT TM band 3.
-#% end
-#% option G_OPT_R_INPUT
-#% key: image4
-#% description: LANDSAT TM band 4.
-#% end
-#% option G_OPT_R_INPUT
-#% key: image5
-#% description: LANDSAT TM band 5.
-#% end
-#% option G_OPT_R_INPUT
-#% key: image7
-#% description: LANDSAT TM band 7.
+#% option G_OPT_F_OUTPUT
+#% description: Name for output file (if omitted or "-" output to stdout)
+#% required: no
 #% end
 #% Flag
 #% key: g
@@ -66,106 +44,111 @@ import sys
 import os
 from grass.script import core as grass
 
-bands = [1,2,3,4,5,7]
 
 def oifcalc(sdev, corr, k1, k2, k3):
-
+    grass.debug(_("Calculating OIF for combination: %s, %s, %s" % (k1, k2,
+                                                                   k3)), 1)
     # calculate SUM of Stddeviations:
     ssdev = [sdev[k1], sdev[k2], sdev[k3]]
     numer = sum(ssdev)
 
     # calculate SUM of absolute(Correlation values):
-    scorr = [corr[k1,k2], corr[k1,k3], corr[k2,k3]]
+    scorr = [corr[k1, k2], corr[k1, k3], corr[k2, k3]]
     denom = sum(map(abs, scorr))
 
     # Calculate OIF index:
     #     Divide (SUM of Stddeviations) and (SUM of Correlation)
     return numer / denom
 
-def perms():
-    for i in range(0,4):
-	for j in range(i+1,5):
-	    for k in range(j+1,6):
-		yield (bands[i], bands[j], bands[k])
+
+def perms(bands):
+    n = len(bands)
+    for i in range(0, n - 2):
+        for j in range(i + 1, n - 1):
+            for k in range(j + 1, n):
+                yield (bands[i], bands[j], bands[k])
+
 
 def main():
     shell = flags['g']
     serial = flags['s']
-    image = {}
-    for band in bands:
-	image[band] = options['image%d' % band]
-
+    bands = options['input'].split(',')
+    output = options['output']
     # calculate the Stddev for TM bands
     grass.message(_("Calculating Standard deviations for all bands..."))
     stddev = {}
 
     if serial:
         for band in bands:
-	    grass.verbose("band %d" % band)
-	    s = grass.read_command('r.univar', flags = 'g', map = image[band])
-	    kv = grass.parse_key_val(s)
-	    stddev[band] = float(kv['stddev'])
+            grass.verbose("band %d" % band)
+            s = grass.read_command('r.univar', flags='g', map=band)
+            kv = grass.parse_key_val(s)
+            stddev[band] = float(kv['stddev'])
     else:
-	# run all bands in parallel
-	if "WORKERS" in os.environ:
+        # run all bands in parallel
+        if "WORKERS" in os.environ:
             workers = int(os.environ["WORKERS"])
-	else:
-	    workers = 6
-
-	proc = {}
-	pout = {}
-
-	# spawn jobs in the background
-	for band in bands:
-	    grass.debug("band %d, <%s>  %% %d" % (band, image[band], band % workers))
-	    proc[band] = grass.pipe_command('r.univar', flags = 'g', map = image[band])
-	    if band % workers is 0:
-		# wait for the ones launched so far to finish
-		for bandp in bands[:band]:
-		    if not proc[bandp].stdout.closed:
-			pout[bandp] = proc[bandp].communicate()[0]
-		    proc[bandp].wait()
-
-	# wait for jobs to finish, collect the output
-	for band in bands:
-	    if not proc[band].stdout.closed:
-		pout[band] = proc[band].communicate()[0]
-	    proc[band].wait()
-
-	# parse the results
+        else:
+            workers = len(bands)
+        proc = {}
+        pout = {}
+
+        # spawn jobs in the background
+        n = 0
+        for band in bands:
+            proc[band] = grass.pipe_command('r.univar', flags='g', map=band)
+            if n % workers is 0:
+                # wait for the ones launched so far to finish
+                for bandp in bands[:n]:
+                    if not proc[bandp].stdout.closed:
+                        pout[bandp] = proc[bandp].communicate()[0]
+                    proc[bandp].wait()
+            n = n + 1
+
+        # wait for jobs to finish, collect the output
         for band in bands:
-	    kv = grass.parse_key_val(pout[band])
-	    stddev[band] = float(kv['stddev'])
+            if not proc[band].stdout.closed:
+                pout[band] = proc[band].communicate()[0]
+            proc[band].wait()
+
+    # parse the results
+        for band in bands:
+            kv = grass.parse_key_val(pout[band])
+            stddev[band] = float(kv['stddev'])
 
 
     grass.message(_("Calculating Correlation Matrix..."))
     correlation = {}
-    s = grass.read_command('r.covar', flags = 'r', map = [image[band] for band in bands])
+    s = grass.read_command('r.covar', flags='r', map=[band for band in bands],
+                           quiet=True)
     for i, row in zip(bands, s.splitlines()):
-	for j, cell in zip(bands, row.split(' ')):
-	    correlation[i,j] = float(cell)
+        for j, cell in zip(bands, row.split(' ')):
+            correlation[i, j] = float(cell)
 
     # Calculate all combinations
-    grass.message(_("Calculating OIF for the 20 band combinations..."))
+    grass.message(_("Calculating OIF for all band combinations..."))
 
     oif = []
-    for p in perms():
-	oif.append((oifcalc(stddev, correlation, *p), p))
-    oif.sort(reverse = True)
+    for p in perms(bands):
+        oif.append((oifcalc(stddev, correlation, *p), p))
+    oif.sort(reverse=True)
 
-    grass.verbose(_("The Optimum Index Factor analysis result "
+    grass.verbose(_("The Optimum Index Factor analysis result " \
                     "(Best combination comes first):"))
-    
+
     if shell:
-	fmt = "%d%d%d:%.4f\n"
+        fmt = "%s%s%s:%.4f\n"
     else:
-	fmt = "%d%d%d:  %.4f\n"
+        fmt = "%s%s%s:  %.4f\n"
 
-    outf = file('i.oif.result', 'w')
-    for v, p in oif:
-	sys.stdout.write(fmt % (p + (v,)))
-	outf.write(fmt % (p + (v,)))
-    outf.close()
+    if not output or output == '-':
+        for v, p in oif:
+            sys.stdout.write(fmt % (p + (v,)))
+    else:
+        outf = file(output, 'w')
+        for v, p in oif:
+            outf.write(fmt % (p + (v,)))
+        outf.close()
 
 if __name__ == "__main__":
     options, flags = grass.parser()