#! /usr/bin/env python
"""
@author: nick.bearson@ssec.wisc.edu

Example usage:
  python bin2hdf.py -f myNpseh_lats.3305x3305x1.float myNpseh_lons.3305x3305x1.float NPe1.bin NPe1.hdf

Dependencies: 
  + Python > 2.5 (tested on 2.7)
  + numpy (http://numpy.scipy.org/)
  + pyhdf (http://pysclint.sourceforge.net/pyhdf/)
"""

# External packages:
import pyhdf.SD as h4
import numpy as np

# Python stdlib:
from optparse import OptionParser
import os
from struct import pack, unpack


DEFAULT_ELEMENTS = 3305
DEFAULT_DTYPE = "h"


DTYPE_SIZES={}
DTYPE_SIZES['c']=1
DTYPE_SIZES['b']=1
DTYPE_SIZES['B']=1
DTYPE_SIZES['?']=1
DTYPE_SIZES['h']=2
DTYPE_SIZES['H']=2
DTYPE_SIZES['i']=4
DTYPE_SIZES['I']=4
DTYPE_SIZES['l']=4
DTYPE_SIZES['L']=4
DTYPE_SIZES['q']=8
DTYPE_SIZES['Q']=8
DTYPE_SIZES['f']=4
DTYPE_SIZES['D']=8


SDC_DTYPES={}
SDC_DTYPES['c']=h4.SDC.CHAR8
SDC_DTYPES['b']=h4.SDC.CHAR8
SDC_DTYPES['B']=h4.SDC.UCHAR8
SDC_DTYPES['h']=h4.SDC.INT16
SDC_DTYPES['H']=h4.SDC.UINT16
SDC_DTYPES['i']=h4.SDC.INT32
SDC_DTYPES['I']=h4.SDC.UINT32
SDC_DTYPES['l']=h4.SDC.INT32
SDC_DTYPES['L']=h4.SDC.UINT32
SDC_DTYPES['f']=h4.SDC.FLOAT32
SDC_DTYPES['D']=h4.SDC.FLOAT64
# The rest don't have any equivalent in the SDC package
# using same-sized data types but the output is likely to be WRONG.
SDC_DTYPES['?']=h4.SDC.INT8
SDC_DTYPES['q']=h4.SDC.FLOAT64
SDC_DTYPES['Q']=h4.SDC.FLOAT64


NP_DTYPES={}
NP_DTYPES['c']=np.str
NP_DTYPES['b']=np.int8
NP_DTYPES['B']=np.uint8
NP_DTYPES['h']=np.int16
NP_DTYPES['H']=np.uint16
NP_DTYPES['i']=np.int32
NP_DTYPES['I']=np.uint32
NP_DTYPES['l']=np.int32
NP_DTYPES['L']=np.uint32
NP_DTYPES['f']=np.float32
NP_DTYPES['D']=np.float64
NP_DTYPES['?']=np.bool
NP_DTYPES['q']=np.int64
NP_DTYPES['Q']=np.uint64


# For valid values of dtype, see: http://docs.python.org/library/struct.html#format-characters
def convert_binaries(latfilename, lonfilename, datafilename, outfilename, dtype=DEFAULT_DTYPE, nele=DEFAULT_ELEMENTS):
  lats = read_binary(latfilename, dtype="f")
  lons = read_binary(lonfilename, dtype="f")
  data = read_binary(datafilename)

  outfile = h4.SD(outfilename, h4.SDC.WRITE|h4.SDC.CREATE)
  
  if lats.shape != lons.shape or lons.shape != data.shape:
    print "WARNING: DIMENSIONS DO NOT MATCH!"

  
  dslat = outfile.create("latitude", SDC_DTYPES["f"], lats.shape)
  dim1 = dslat.dim(0)
  dim2 = dslat.dim(1)
  dim1.setname('lines')
  dim2.setname('elements')
  dslat[:] = lats[:]
  dslat.endaccess()

  dslon = outfile.create("longitude", SDC_DTYPES["f"], lons.shape)
  dim1 = dslon.dim(0)
  dim2 = dslon.dim(1)
  dim1.setname('lines')
  dim2.setname('elements')
  dslon[:] = lons[:]
  dslon.endaccess()

  dsdata = outfile.create("data", SDC_DTYPES[dtype], data.shape)
  dim1 = dsdata.dim(0)
  dim2 = dsdata.dim(1)
  dim1.setname('lines')
  dim2.setname('elements')
  dsdata[:] = data[:]
  dsdata.endaccess()

  outfile.end()



def read_binary(bfilename, dtype=DEFAULT_DTYPE, nele=DEFAULT_ELEMENTS):
  blockstr = str(nele) + dtype
  unpacked = []
  with open(bfilename, 'rb') as f:
    raw = f.read(nele * DTYPE_SIZES[dtype])
    while raw != "":
      unpacked.append(unpack(blockstr, raw))
      raw = f.read(nele * DTYPE_SIZES[dtype])
  return np.array(unpacked, dtype=NP_DTYPES[dtype])


def usage():
  print ""


def handle_args():
  usage = "usage: %prog [options] lat.bin lon.bin data.bin output.hdf\n"
  usage+= "  lat.bin and lon.bin assumed to be of type float32\n"
  parser = OptionParser(usage=usage)
  parser.add_option("-d", "--dtype", action="store", type="string", dest="dtype", help="data type of binary data file\n(see: http://docs.python.org/library/struct.html#format-characters)\ndefault: 'h'") 
  parser.add_option("-n", "--nele", action="store", type="int", dest="nele", help="number of elements in binary data file") 
  parser.add_option("-f", "--force", action="store_true", dest="force", help="force completion, doing things like overwriting the output file if it already exists")
  return parser.parse_args()



if __name__ == '__main__':
  (options, args) = handle_args()

  if len(args) < 4:
    usage()
    exit(-1)

  for file in args[0:2]:
    if not os.path.exists(file):
      print "File doesn't exist: " + file
      usage()
      exit(-1)

  if os.path.exists(args[3]):
    if options.force:
      print "Output file already exists, force set (removing): " + args[3]
      os.remove(args[3])
    else:
      print "Output file already exists, exiting: " + args[3]
      exit(-1)

  latfilename = args[0]
  lonfilename = args[1]
  datafilename = args[2]
  outfilename = args[3]

  convert_binaries(latfilename, lonfilename, datafilename, outfilename)