#!/dmsp/reference/bin/ruby

#
# ara's lib
#
  require 'alib'
  require 'mmap'
  require 'narray'
  require 'envi'
#
# main program class
#
  class Main < ALib::SimpleMain 
#--{{{
    version '0.0.0'
    author 'ara.t.howard@noaa.gov'

    usage(
#--{{{
      <<-usage
        NAME
          #{ prognam } v#{ version }

        SYNOPSIS
          #{ prognam } [options]+ red_grid green_grid blue_grid rgb_grid

        DESCRIPTTION

        ENVIRONMENT
          none

        DIAGNOSTICS
          success -> $? == 0
          failure -> $? != 0

        AUTHOR
          #{ author }

        BUGS
          1 < bugno && bugno < 42
      usage
#--}}}
    )

    examples(
#--{{{
      <<-examples
        EXAMPLES

          0) generate an rgb image from vis an tir channels

            #{ prognam } F15200509050235.katrina.OIS.vis F15200509050235.katrina.OIS.tir F15200509050235.katrina.OIS.rgb
      examples
#--}}}
    )

    argument 'vis_grid', 'tir_grid', 'rgb_grid'

    def main
#--{{{
      load_grids
      gen_rgb
      write_rgb
      write_rgb_header
      exit EXIT_SUCCESS
#--}}}
    end
    def load_grids
#--{{{
      grid_opts = { 
        'cast' => 'int',
        'mmap' => {'mode' => 'r', 'flags' => Mmap::MAP_SHARED}
      }
      @rg = ENVI::Grid::new vis_grid, grid_opts
      @gg = ENVI::Grid::new vis_grid, grid_opts
      @bg = ENVI::Grid::new tir_grid, grid_opts
#--}}}
    end
    def gen_rgb
#--{{{
      r = @rg.narray
      g = @gg.narray
      b = @bg.narray

      r[ (r.eq 255).where ] = 0 # hide missing
      g[ (g.eq 255).where ] = 0 # hide missing

      linear_stretch r, 7 .. 30
      linear_stretch g, 7 .. 30
      linear_stretch b, 175 .. 230

      r[] = r * (1 << 24) 
      g[] = g * (1 << 16) 
      b[] = (255 - b) * (1 << 8)

      @rgb = r | g | b
      @rgb = @rgb.hton
      r = g = b = nil
#--}}}
    end
    def linear_stretch na, irange, orange = (0..255)
#--{{{
      fna = NArray::float(*na.shape)
      fna[] = na

      imax, imin = irange.max, irange.min
      isize = imax - imin + 1

      omax, omin = orange.max, orange.min
      osize = omax - omin + 1

      # floor input
      lte_imin, gt_imin = (fna <= imin).where2
      fna[ lte_imin ] = 0
      fna[ gt_imin ] = fna[ gt_imin ] - imin + 1

      # ceil input
      gt_imax = (fna > isize).where
      fna[ gt_imax ] = isize + 1

      # percentage of in range 
      fna[] = fna / (isize + 1)

      # scale into out range
      fna[] = fna * (osize - 1)

      na[] = fna
#--}}}
    end
    def write_rgb
#--{{{
      open(rgb_grid, 'w'){|f| f.write @rgb.to_s}
#--}}}
    end
    def write_rgb_header
#--{{{
      h = @rg.header
      open(rgb_grid + ".hdr", 'w') do |f|
        hdr = <<-hdr
          ENVI
          description = { R (#{ vis_grid }), G (#{ vis_grid }), B (#{ tir_grid }), A (alpha) }
          samples = #{ h.samples } 
          lines   = #{ h.lines } 
          bands   = 4
          header offset = 0
          file type = ENVI Standard
          data type = 1
          interleave = bip
          sensor type = Unknown
          byte order = 1
          map info = #{ h.map_info } 
          wavelength units = Unknown
          default stretch = 0.000000 255.000000 linear
          data origin = NGDC
          band names = { R (Band 1:#{ vis_grid }), G (Band 2:#{ vis_grid }), B (Band 3:#{ tir_grid }), A (Band 4:alpha) }
        hdr
        hdr.each do |line|
          f.puts line.strip
        end
      end
#--}}}
    end
#--}}}
  end
#
# run unless included as a lib 
#
  Main::new.run if $0 == __FILE__