#Sentinel 5 - CH4 import
#Sentinel 5 data from Sentinel 5 Pre-Operations Data Hub, available at https://s5phub.copernicus.eu/dhus/#/home


from tkinter import *
from tkinter.filedialog import askdirectory
from tkinter.filedialog import askopenfilename
from tkinter.ttk import Progressbar, Style
import os
import sys
import numpy as np
from netCDF4 import Dataset
from osgeo import gdal
import ilwis
import shutil
from threading import Thread

#passing parameters from sub menu
def execute(SentDir, InFileName, WorkDir, OutFileName, exportFormat):
    print(exportFormat)
    if len(SentDir) > 0 and len(InFileName) > 0 and len(WorkDir) > 0 and len(OutFileName) > 0:

        progress['value'] = 0

        os.chdir(WorkDir)

        #process geodetic information
        #using geodetic_in.nc

        nc_f = os.path.join(SentDir,InFileName)
        shutil.copy(SentDir+'\\'+InFileName, WorkDir)

        nc_fid = Dataset(nc_f, 'r')
        lats = nc_fid.groups['PRODUCT']['latitude'][:][0]


        height, width = lats.shape
        offset = 0
		#multiplication_factor_to_convert_to_molecules_per_cm2 = 6.02214E19f
		#scaling used here is 46000 / conversion from mol/m2 to mg/m2
		#46 (g/mol) * 1000 (mg/g) = 46000
        scale = 46000

        #create gdal vrt file
		#check assignment of coordinate origin - here centre of pixel is assumed and NOT the corner
        vrt_file = open(OutFileName + '.vrt', "w")
        vrt_file.write('<VRTDataset rasterXSize="' + str(width) + '" rasterYSize="' + str(height) + '">\n')
        vrt_file.write('  <Metadata domain="GEOLOCATION">\n')
        vrt_file.write('    <MDI key="spatial_resolution">7x3.5km2</MDI>\n')
        vrt_file.write('    <MDI key="title">TROPOMI/S5P NO2 1-Orbit L2 Swath 7x3.5km</MDI>\n')
        vrt_file.write('    <MDI key="SRS">WGS84</MDI>\n')
        vrt_file.write('    <MDI key="X_DATASET">HDF5:"'+(InFileName)+'"://PRODUCT/longitude</MDI>\n')
        vrt_file.write('    <MDI key="X_BAND">1</MDI>\n')
        vrt_file.write('    <MDI key="Y_DATASET">HDF5:"'+(InFileName)+'"://PRODUCT/latitude</MDI>\n')
        vrt_file.write('    <MDI key="Y_BAND">1</MDI>\n')
        vrt_file.write('    <MDI key="PIXEL_OFFSET">0.5</MDI>\n')
        vrt_file.write('    <MDI key="PIXEL_STEP">1</MDI>\n')
        vrt_file.write('    <MDI key="LINE_OFFSET">0.5</MDI>\n')
        vrt_file.write('    <MDI key="LINE_STEP">1</MDI>\n')
        vrt_file.write('  </Metadata>\n')
        vrt_file.write('  <VRTRasterBand dataType="Float32" band="1">\n')
        vrt_file.write('    <Metadata>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_long_name">methane_mixing_ratio</MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_add_offset">' + str(offset) + ' </MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_scale_factor">' + str(scale) + ' </MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_multiplication_factor_to_convert_to_molecules_percm2">6.02214e+19</MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_standard_name">methane_mixing_ratio</MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'_units">mol m-2</MDI>\n')
        vrt_file.write('      <MDI key="'+(InFileName[:-3])+'__FillValue">9.96921e+36</MDI>\n')
        vrt_file.write('    </Metadata>\n')
        vrt_file.write('    <SimpleSource>\n')
        vrt_file.write('      <SourceFilename relativeToVRT="1">HDF5:"'+(InFileName)+'"://PRODUCT/methane_mixing_ratio</SourceFilename>\n')
        vrt_file.write('      <SourceBand>1</SourceBand>\n')
        vrt_file.write('      <SourceProperties RasterXSize="' + str(width) + '" RasterYSize="' + str(height) + '" DataType="Float32" BlockXSize="' + str(width) + '" BlockYSize="' + str(height) + '" />\n')
        vrt_file.write('      <SrcRect xOff="0" yOff="0" xSize="' + str(width) + '" ySize="' + str(height) + '" />\n')
        vrt_file.write('      <DstRect xOff="0" yOff="0" xSize="' + str(width) + '" ySize="' + str(height) + '" />\n')
        vrt_file.write('    </SimpleSource>\n')
        vrt_file.write('  </VRTRasterBand>\n')
        vrt_file.write('</VRTDataset>\n')
        vrt_file.close()

        progress['value'] = 40

        #conduct gdal warp operation using vrt file created and save file as ilwis mpr
        input_raster = gdal.Open(OutFileName + '.vrt')
        output_raster = (OutFileName+'.mpr')
        gdal.Warp(output_raster,input_raster,format='ILWIS',srcSRS='WGS84',resampleAlg='cubic',geoloc=True, srcNodata='9.96921e+36' ,dstNodata='-999')
        os.remove ((WorkDir.replace('\\','/')) + '/' + OutFileName +'.vrt')

        progress['value'] = 65 

        #final radiometric correction applying ilwis-objects - loaded as plugin under python
        ilwis.Engine.setWorkingCatalog('file:///'+ WorkDir)
        rc = ilwis.RasterCoverage(output_raster)
        print (rc.name())
        if exportFormat == 'ilwis':
            rc2=ilwis.Engine.do('mapcalc','iff((@1 != -999),@1*'+str(scale)+'+'+str(offset)+',?)',rc)
            rc2.store((OutFileName) +'.mpr', 'map', 'ilwis3')
            os.remove ((WorkDir.replace('\\','/')) + '/' + OutFileName +'.mpr.aux.xml')

        elif exportFormat == 'qgis':
        #print ('Output file in QGIS-Geotiff format created = '+ OutFileName +'.tif')
        #set nodata (NaN) (of surrounding area) to -999 for display in QGIS
            rc2=ilwis.Engine.do('mapcalc','iff((@1 != -999),@1*'+str(scale)+'+'+str(offset)+',-999)',rc)
            rc2.store((OutFileName) +'.tif', "GTiff", "gdal")
            os.remove ((WorkDir.replace('\\','/')) + '/' + OutFileName +'.mpr.aux.xml')
            os.remove((WorkDir.replace('\\', '/'))+'/'+ OutFileName +'.csy')
            os.remove((WorkDir.replace('\\', '/'))+'/'+ OutFileName +'.grf')
            os.remove((WorkDir.replace('\\', '/'))+'/'+ OutFileName +'.mpr')
            os.remove((WorkDir.replace('\\', '/'))+'/'+ OutFileName +'.mp#')
            
        progress['value'] = 100
    button0.config(state='normal')
    button1.config(state='normal')

#redundant files deleted in main menu script

def executeInThread(SentDir, InFileName, WorkDir, OutFileName, exportFormat):
    button0.config(state='disabled')
    button1.config(state='disabled')
    thread = Thread(target = execute, args = (SentDir, InFileName, WorkDir, OutFileName, exportFormat))
    thread.start()

#main()

sent_path = os.path.realpath('.')
print (sent_path)
os.chdir(sent_path)
os.environ["sent_lib"] = sent_path

app = Tk()
app.title("Import S5P L2-CH4")
app.iconbitmap(sent_path +'\pics\esa.ico')

w = 650 # width for the Tk app
h = 210 # height for the Tk app

app.geometry('%dx%d' % (w, h))
app.resizable(0,0)

labelText1 = StringVar()
labelText1.set ("""S5P Sentinel-5 precursor/TROPOMI Level 2 Methane (mixing ratio)
Typical filename: S5P_NRTI_L2__CH4____Start of granule in UTC_End of granule in UTC_Orbit number_*_*_*.nc""")

label1 = Label(app, textvariable=labelText1, anchor='w', justify = LEFT)
label1.pack(side = LEFT)
label1.place(x = 10, y = 5, width=650, height=35)

L1 = Label(app, text="Input Directory:", anchor='w')
L1.pack(side = LEFT)
L1.place(x = 10, y = 50, width=150, height=25)

inputDir = StringVar()
E1 = Entry(app, textvariable=inputDir, bd =2)
E1.pack()
E1.place(x = 165, y = 50, width=445, height=25)

button_dirin = Button(app, text="....", width=2, command = lambda: inputDir.set(askdirectory(initialdir = inputDir.get())))
button_dirin.pack()
button_dirin.place(x = 615, y = 49)

L2 = Label(app, text="Input Filename:", anchor='w')
L2.pack(side = LEFT)
L2.place(x = 10, y = 80, width=150, height=25)

inputFile = StringVar()
E2 = Entry(app, textvariable=inputFile, bd =2)
E2.pack()
E2.place(x = 165, y = 80, width=445, height=25)

button_filein = Button(app, text="....", width=2, command = lambda: inputFile.set(askopenfilename(initialdir = inputDir.get()).split('/')[-1]))
button_filein.pack()
button_filein.place(x = 615, y = 79)

L3 = Label(app, text="Output Directory:", anchor='w')
L3.pack(side = LEFT)
L3.place(x = 10, y = 110, width=150, height=25)

outputDir = StringVar()
E3 = Entry(app, textvariable=outputDir, bd =2)
E3.pack()
E3.place(x = 165, y = 110, width=445, height=25)

button_dirout = Button(app, text="....", width=2, command = lambda: outputDir.set(askdirectory(initialdir = outputDir.get())))
button_dirout.pack()
button_dirout.place(x = 615, y = 109) 

L4 = Label(app, text="Output Filename:", anchor='w')
L4.pack(side = LEFT)
L4.place(x = 10, y = 140, width=150, height=25)

outputFile = StringVar()
E4 = Entry(app, textvariable=outputFile, bd =2)
E4.pack()
E4.place(x = 165, y = 140, width=445, height=25)

button0 = Button(app, text="Quit", width=15, command=app.destroy)
button0.pack()
button0.place(x = 525, y = 180)

exportFormat = StringVar()
   
R1 = Radiobutton(app, text="ILWIS", variable=exportFormat, value='ilwis')
R1.pack()
R1.place(x = 165, y = 180)

R2 = Radiobutton(app, text="QGIS", variable=exportFormat, value='qgis')
R2.pack()
R2.place(x = 285, y = 180)

R1.select()
R2.deselect()

button1 = Button(app, text="Execute", width=15, command=lambda: executeInThread(inputDir.get(), inputFile.get(), outputDir.get(), outputFile.get(), exportFormat.get()))
button1.pack()
button1.place(x = 405, y = 180)

style = Style()
style.theme_use('clam')
style.configure('blue.Horizontal.TProgressbar', foreground='blue', background='blue')

progress = Progressbar(app, style='blue.Horizontal.TProgressbar', orient=HORIZONTAL, length=100, mode='determinate')
progress.pack()
progress.place(x = 10, y = 180)
app.mainloop()

print("inputFile=" + inputFile.get())
print("outputDir=" + outputDir.get())
print("outputFile=" + outputFile.get())

