"""

SNAPPY TOPSAR Coreg Interferogram script from Alaska Satellite Facility (ASF) Distributed Active Archive Center (DAAC)

TO USE:
    (1.) Install SNAP, Snappy, and S1Toolbox from ESA (http://step.esa.int/main/download/)
    (2.) search file for the words <your path> and replace with appropriate path on your system
    (3.) Download the sample granules from ASF (https://vertex.daac.asf.alaska.edu/) by searching for these granule names:
            S1A_IW_SLC__1SSV_20150417T001852_20150417T001922_005516_0070C1_460B.zip
            S1A_IW_SLC__1SDV_20150429T001842_20150429T001909_005691_0074DC_1CA1.zip
    (O.) optionally search for <OPTIONAL> and do what's required

"""

import os, sys
sys.path.append('/<your path>/snappy')
sys.path.append('/<your path>/snap/')
sys.path.append('/<your path>/s1tbx')

from snappy import GPF
from snappy import ProductIO
from snappy import HashMap
from Queue import *
from thread import *
from time import *


iDebug = True
iWriteEachStep = False

GPF.getDefaultInstance().getOperatorSpiRegistry().loadOperatorSpis()
prods = []
file1 = ""
file2 = ""
strFile1 = ""
strFile2 = ""
strInSAR = ""
target1 = ""
target2 = ""
targetTOPSplit1 = ""
targetTOPSplit2 = ""
targetbackGeo = ""
targetInterferogram = ""
targetTOPSARDeburst = ""


q = Queue()



def readFiles(filename1, filename2):
    global file1
    global file2
    global strFile1
    global strFile2
    global strInSAR

    # get my 2 files
    if iDebug: print filename1
    if iDebug: print filename2
    
    file1 = ProductIO.readProduct(filename1)
    if iDebug: print "Read file: " + file1.getName()
    strFile1 = file1.getName()
    file2 = ProductIO.readProduct(filename2)
    if iDebug: print "Read file: " + file2.getName() + ('\n')
    strFile2 = file2.getName()

    # InSAR Stack name: grab the leading common bits, and remove the dates et al
    strFind = "20"
    strInSAR = file1.getName()
    i = file1.getName().find(strFind)
    if i: strInSAR = file1.getName()[:i]
    if iDebug: print strInSAR


def writeFiles(target, filename, strStep):
    ofile = '/<your path>/' + filename + strStep
    if iDebug: print ofile
    if target:
        ProductIO.writeProduct(target, ofile, 'BEAM-DIMAP')   # TO DO:  allow other types here: BEAM-DIMAP, GeoTIFF-BigTIFF, HDF5
    else:
        print "writeFiles Error: No target given"


def tOPSARSplit():
    ### TOPSAR-Split
    global file1
    global file2
    global targetTOPSplit1
    global targetTOPSplit2

    parameters = ""
    parameters = HashMap()
    parameters.put('subswath', 'IW1')
    parameters.put('selectedPolarisations', 'VV')
    if iDebug: print "Entering: TOPSAR-Split"

    targetTOPSplit1 = GPF.createProduct("TOPSAR-Split", parameters, file1)

    targetTOPSplit2 = GPF.createProduct("TOPSAR-Split", parameters, file2)
    if iDebug: print "TOPSAR-Split working"


def applyOrbitFile():
    global targetTOPSplit1
    global targetTOPSplit2
    global target1
    global target2
    global prods

    if iDebug: print "Entering: ApplyOrbit"

    parameters = HashMap()
    parameters.put("Orbit State Vectors", "Sentinel Precise (Auto Download)")
    parameters.put("Polynomial Degree", 3)
    target1 = GPF.createProduct("Apply-Orbit-File", parameters, targetTOPSplit1)
    if iDebug: print "ApplyOrbit working: 1"


    target2 = GPF.createProduct("Apply-Orbit-File", parameters, targetTOPSplit2)
    if iDebug: print "ApplyOrbit working: 2"
    prods.append(target2)
    prods.append(target1)
    if iDebug: print "Added files in reverse order" # pick a master?


def backGeocoding():
    ### backGeocoding
    global prods
    global targetbackGeo


    parameters = ""
    parameters = HashMap()
    parameters.put("Digital Elevation Model", "SRTM 3Sec (Auto Download)")
    parameters.put("DEM Resampling Method", "BICUBIC_INTERPOLATION")
    parameters.put("Resampling Type", "BISINC_5_POINT_INTERPOLATION")
    parameters.put("Mask out areas with no elevation", True)
    parameters.put("Output Deramp and Demod Phase", False)

    if iDebug:
        for u in prods:
            print u

    targetbackGeo = GPF.createProduct("Back-Geocoding", parameters, prods)
    if iDebug: print "Back-Geocoding working"


def interferogram():
    ### Interferogram
    # print GPF.getDefaultInstance().getOperatorSpiRegistry().getOperatorSpi('Interferogram')
    # Interferogram = jpy.get_type('org.esa.s1tbx.insar.gpf.CreateInterferogramOp')
    global targetbackGeo
    global targetInterferogram

    if iDebug: print "Entering: Interferogram"

    parameters = ""
    parameters = HashMap()
    parameters.put("Subtract flat-earth phase", True)
    parameters.put("Degree of \"Flat Earth\" polynomial", 5)
    parameters.put("Number of \"Flat Earth\" estimation points", 501)
    parameters.put("Orbit interpolation degree", 3)
    parameters.put("Include coherence estimation", True)
    parameters.put("Square Pixel", False)
    parameters.put("Independent Window Sizes", False)   
    parameters.put("Coherence Azimuth Window Size", 10)
    parameters.put("Coherence Range Window Size", 10)
    targetInterferogram = GPF.createProduct("Interferogram", parameters, targetbackGeo)
    if iDebug: print "Interferogram working"


def tOPSARDeburst():
    ### TOPSAR-Deburst
    global targetInterferogram
    global targetTOPSARDeburst

    if iDebug: print "Entering: TOPSARDeburst"

    parameters = ""
    parameters = HashMap()
    parameters.put("Polarisations", "VV")
    targetTOPSARDeburst = GPF.createProduct("TOPSAR-Deburst", parameters, targetInterferogram)
    if iDebug: print "TOPSARDeburst working"


def DestinationThread() :
    while True :
        f, args = q.get()
        f(*args)


### Main ###
#
# TO DO:  input filenames
#
#

# NOTE:  add them in date order, oldest to newest
filename1 = "/<your path>/S1A_IW_SLC__1SSV_20150417T001852_20150417T001922_005516_0070C1_460B.zip"
filename2 = "/<your path>/S1A_IW_SLC__1SDV_20150429T001842_20150429T001909_005691_0074DC_1CA1.zip"


start_new_thread( DestinationThread, tuple() )

print "Snappy TOPSAR InSAR: start"

### Read in the files
sleep( 1 )
q.put( (readFiles, [filename1, filename2]), True)
sleep( 1 )


### tOPSAR-Split step
q.put( (tOPSARSplit, ""), True)
sleep( 1 )
if iWriteEachStep: 
    q.put( (writeFiles, [targetTOPSplit1, strFile1, "_tOPSARSplit"]), True)
    sleep( 1 )
    q.put( (writeFiles, [targetTOPSplit2, strFile2, "_tOPSARSplit"]), True)
    sleep( 1 )


### ApplyOrbitFile step
q.put( (applyOrbitFile, ""), True)
sleep( 1 )
if iWriteEachStep: 
    q.put( (writeFiles, [target1, strFile1, "_AppOrb"]), True)
    sleep( 1 )
    q.put( (writeFiles, [target2, strFile2, "_AppOrb"]), True)
    sleep( 1 )


### backGeocoding step
q.put( (backGeocoding, ""), True)
sleep( 1 )
if iWriteEachStep:
    q.put( (writeFiles, [targetbackGeo, strInSAR, "_CFCoReg"]), True)
    sleep( 1 )


### Interferogram step
q.put( (interferogram, ""), True)
sleep( 1 )
if iWriteEachStep:
    q.put( (writeFiles, [targetInterferogram, strInSAR, "_Interf"]), True)
    sleep( 1 )


### TOPSARDeburst step  (LAST step so naming this _InSAR)
q.put( (tOPSARDeburst, ""), True)
sleep( 1 )
q.put( (writeFiles, [targetTOPSARDeburst, strInSAR, "_srp_ifgs"]), True)
sleep( 1 )



print "Snappy TOPSAR InSAR: end"


if iDebug: print "Cleanup tmp files"
myStr = "image"
path = "/<your path>/.snap/var/cache/temp"
dirs = os.listdir( path )

# list files
for file in dirs:
    if iDebug: print file
    # make sure it's an image file
    if file.startswith(myStr):
        # delete them
        os.remove(path + "/" + file)
if iDebug: print "Cleanup complete"