[GRASS-SVN] r42745 - grass/trunk/imagery/i.atcorr

svn_grass at osgeo.org svn_grass at osgeo.org
Fri Jul 9 13:49:08 EDT 2010


Author: aghisla
Date: 2010-07-09 17:49:08 +0000 (Fri, 09 Jul 2010)
New Revision: 42745

Added:
   grass/trunk/imagery/i.atcorr/create_iwave.py
Modified:
   grass/trunk/imagery/i.atcorr/
   grass/trunk/imagery/i.atcorr/i.atcorr.html
Log:
Added Python script that converts csv files to cpp functions, small edits on docpage, svn properties setup.


Property changes on: grass/trunk/imagery/i.atcorr
___________________________________________________________________
Modified: svn:ignore
   - 
*OBJ*

   + *OBJ*
*.pyc


Added: grass/trunk/imagery/i.atcorr/create_iwave.py
===================================================================
--- grass/trunk/imagery/i.atcorr/create_iwave.py	                        (rev 0)
+++ grass/trunk/imagery/i.atcorr/create_iwave.py	2010-07-09 17:49:08 UTC (rev 42745)
@@ -0,0 +1,256 @@
+#!/usr/bin/python2.6
+"""
+Created on Sat Mar 27 11:35:32 2010
+
+Program to interpolate filter function to correct 
+step. Should be 2.5 nm
+Then output filter function in a format similar to
+what is needed in the Iwave.cpp file
+
+Needs numpy and scipy
+
+ at author: daniel
+contact: daniel {dot} victoria {at} gmail {dot} com
+
+usage() explains how this is supposed to work
+Basically it needs a .csv file with spectral response for each
+band in a column. First column has to be wavelength (nm)
+First line (and only first) is a header with Wl, followed by band names
+file name is used for sensor name
+
+Self though programer warning: beware of fractured python bellow
+and fractured english also :)
+"""
+import os
+import sys
+import numpy as np
+from scipy import interpolate
+
+def usage():
+    """How to use this..."""
+    print "create_iwave.py <csv file>"
+    print "Generate filter function IWave.cpp template from csv file"
+    print "csv file must have wl response for each band in each column"
+    print "first line must be a header with wl followed by band names"
+    print "following lines will be the data."
+    print "If response is null, leave empty in csv file. Ex.:"
+    print "WL(nm),band 1,band 2,band 3,band 4"
+    print "455,0.93,,,"
+    print "485,0.94,0.00,,"
+    print "545,0.00,0.87,0.00,"
+    print "Program will interpolate filter function to 2.5 nm steps"
+    print "and output a cpp template file in the IWave format"
+
+def read_input(csvfile):
+    """
+    Func to read input file
+    return number of bands and array of values for each band
+    should be a .csv file with the values
+    of the filter function for each band in the sensor
+    one column for band
+    first line must have a header with sensor band name
+    first column is wavelength
+    """
+    infile = open(csvfile, 'r')
+    
+    # get number of bands and band names
+    bands = infile.readline().split(',')
+    bands.remove(bands[0])
+    bands[-1] = bands[-1].strip()
+    
+    infile.close()
+    
+    # create converter dictionary for import
+    # fix nodata or \n
+    conv = {}
+    for b in range(len(bands)):
+        conv[b+1] = lambda s: float(s or -99)
+    
+    
+    values = np.loadtxt(csvfile, delimiter=',', skiprows=1, converters = conv)
+    
+    return(bands, values)
+
+def interpolate_band(values):
+    """
+    Receive wavelength and response for one band
+    interpolate at 2.5 nm steps
+    return interpolated filter func
+    and min, max wl values
+    values must be numpy array with 2 columns
+    """
+    
+    # remove nodata (-99) lines in values array
+    # where response is nodata?
+    w = values[:,1] >= 0
+    response = values[w]
+    
+    # interpolating
+    f = interpolate.interp1d(response[:,0],response[:,1])
+    
+    filter_f = f(np.arange(response[0,0], response[-1,0], 2.5))
+    
+    return(filter_f, (response[0,0], response[-1,0]))
+
+def plot_filter(values):
+    """Plot wl response values and interpolated
+    filter function. This is just for checking...
+    value is a 2 column numpy array
+    function has to be used inside Spyder python environment
+    """
+    filter_f, limits = interpolate_band(values)
+    
+    # removing nodata
+    w = values[:,1] >= 0
+    response = values[w]
+    
+    plot(response[:,0],response[:,1], 'ro')
+    plot(arange(limits[0], limits[1], 2.5), filter_f)
+    
+    return
+
+def pretty_print(filter_f):
+    """
+    Create pretty string out of filter function
+    8 values per line, with spaces, commas and all the rest
+    """
+    nlines = int(len(filter_f)/8)
+    rest = len(filter_f) - nlines*8
+    pstring = '    '
+    
+    if rest == 0: # multiple of 8; special case
+        for l in range(nlines):
+            for i in range(8):
+                if l == nlines-1 and i == 7: # last value
+                    pstring = pstring+' %.2f };\n' % (filter_f[l*8+i])
+                else:
+                    pstring = pstring+' %.2f,' % (filter_f[l*8+i])
+            pstring = pstring+'\n    '
+    else:
+        for l in range(nlines):
+            for i in range(8):
+                pstring = pstring+'    %.2f,' % (filter_f[l*8+i])
+            pstring = pstring+'\n    '
+        # adding rest of values in filter_f
+        for i in range(rest):
+            if i == rest-1: # last value
+                pstring = pstring+'    %.2f };\n' % (filter_f[nlines*7+i])
+            else:
+                pstring = pstring+'    %.2f,' % (filter_f[nlines*7+i])
+    
+    return pstring
+
+def write_cpp(bands, values, sensor):
+    """
+    from bands, values and sensor name
+    create output file in cpp style
+    needs other functions: interpolate_bands, pretty_print
+    """
+    
+    # getting necessary data
+    # single or multiple bands?
+    if len(bands) == 1:
+        filter_f, limits = interpolate_band(values)
+    else:
+        filter_f = []
+        limits = []
+        for b in range(len(bands)):
+            fi, li = interpolate_band(values[:,[0,b+1]])
+            filter_f.append(fi)
+            limits.append(li)
+    
+    # writing...
+    outfile = open(sensor+"_cpp_template.txt", 'w')
+    outfile.write('/* Following filter function created using create_iwave.py */\n\n')
+    
+    if len(bands) == 1:
+        outfile.write('void IWave::%s()\n{\n\n' % (sensor))
+    else:
+        outfile.write('void IWave::%s(int iwa)\n{\n\n' % (sensor))
+        
+    # single band case
+    if len(bands) == 1:
+        outfile.write('    /* %s of %s */\n' % (bands[0], sensor))
+        outfile.write('    static const float sr[%i] = {\n' % (len(filter_f)))
+        filter_text = pretty_print(filter_f)
+        outfile.write(filter_text)
+        
+        # calculate wl slot for band start
+        # slots range from 250 to 4000 at 2.5 increments (total 1500)
+        s_start = int((limits[0] - 250)/2.5)
+        
+        outfile.write('\n')
+        outfile.write('    ffu.wlinf = %.4ff;\n' % (limits[0]/1000.0))
+        outfile.write('    ffu.wlsup = %.4ff;\n' % (limits[1]/1000.0))
+        outfile.write('    int i = 0;\n')
+        outfile.write('    for(i = 0; i < %i; i++)\tffu.s[i] = 0;\n' % (s_start))
+        outfile.write('    for(i = 0; i < %i; i++)\tffu.s[%i+i] = sr[i];\n' % (len(filter_f), s_start))
+        outfile.write('    for(i = %i; i < 1501; i++)\tffu.s[i] = 0;\n' % (s_start + len(filter_f)))
+        outfile.write('}\n')
+        
+    else: # more than 1 band
+        # writing bands
+        for b in range(len(bands)):
+            outfile.write('    /* %s of %s */\n' % (bands[b], sensor))
+            outfile.write('    static const float sr%i[%i] = {\n' % (b+1,len(filter_f[b])))
+            filter_text = pretty_print(filter_f[b])
+            outfile.write(filter_text+'\n')
+        
+        # writing band limits
+        for b in range(len(bands)):
+            inf = str([i[0] for i in limits])
+            sup = str([i[1] for i in limits])
+            
+            for j in zip(['[',']'],['{','}']):
+                inf = inf.replace(j[0],j[1])
+                sup = sup.replace(j[0],j[1])
+        
+        outfile.write('    static const float wli[%i] = %s;\n' % (len(bands), inf))
+        outfile.write('    static const float wls[%i] = %s;\n' % (len(bands), sup))
+        
+        outfile.write('\n')
+        outfile.write('    ffu.wlinf = (float)wli[iwa-1];\n')
+        outfile.write('    ffu.wlsup = (float)wls[iwa-1];\n\n')
+        
+        outfile.write('    int i;\n')
+        outfile.write('    for(i = 0; i < 1501; i++) ffu.s[i] = 0;\n\n')
+        
+        outfile.write('    switch(iwa)\n    {\n')
+        
+        # now start of case part...
+        for b in range(len(bands)):
+            s_start = int((limits[b][0] - 250)/2.5)
+            outfile.write('    case %i: for(i = 0; i < %i; i++)  ffu.s[%i+i] = sr%i[i];\n' % ((b+1), len(filter_f[b]), s_start, (b+1)))
+            outfile.write('        break;\n')
+        outfile.write('    }\n}\n')
+        
+    return
+
+def main():
+    """ control function """
+    
+    inputfile = sys.argv[1]
+    
+    # getting sensor name from full csv file name
+    sensor = os.path.splitext(os.path.basename(inputfile))[0]
+    
+    print "Getting sensor name from csv file: %s" % (sensor)
+    
+    # getting data from file
+    bands, values = read_input(inputfile)
+    
+    # writing file
+    write_cpp(bands, values, sensor)
+    
+    print "Filter function written to %s" % (sensor+"_cpp_template.txt")
+    print "Please check file for possible errors before inserting into IWave.cpp"
+    print "Don't forget to add necessary data to IWave.h"
+    
+    return
+
+if __name__ == '__main__':
+    if len(sys.argv) == 1:
+        usage()
+        sys.exit()
+    else:
+        main()


Property changes on: grass/trunk/imagery/i.atcorr/create_iwave.py
___________________________________________________________________
Added: svn:mime-type
   + text/python
Added: svn:keywords
   + Author Date Id
Added: svn:eol-style
   + native

Modified: grass/trunk/imagery/i.atcorr/i.atcorr.html
===================================================================
--- grass/trunk/imagery/i.atcorr/i.atcorr.html	2010-07-09 16:09:25 UTC (rev 42744)
+++ grass/trunk/imagery/i.atcorr/i.atcorr.html	2010-07-09 17:49:08 UTC (rev 42745)
@@ -141,7 +141,7 @@
 </tr>
 
 <tr>
-<td>10</td>
+<td>11</td>
 <td><b>avnir</b></td>
 <td>enter month,day,hh.ddd,long.,lat. *</td>
 </tr>



More information about the grass-commit mailing list