import matplotlib
# Allow running headless from the command line
matplotlib.use("agg")

from pylab import *

import numpy as np
import pandas as pd
import mutils as m
import os
import io
from PIL import Image
import psutil
import sys
import datetime
import json
import time

elevThres = 7

def savePlot(directory,name):
  # Save the data as a PNG file
  if not os.path.exists(directory):
    os.makedirs(directory)

  # Save the png to memory
  ram = io.BytesIO()
  savefig(ram,format='png',dpi=150)
  ram.seek(0)
  # Compress the png data before saving to the file system, the
  # matplotlib png's are not well compressed
  im = Image.open(ram)
  im2 = im.convert('RGB').convert('P', palette=Image.ADAPTIVE)
  im2.save(directory + "/" + name,format='PNG')

# Override pylab find() to avoid deprecation warning
# pylint: disable=function-redefined
def find(x):
  return where(x)[0]

sysData = [ {"sysName":"GPS", "sys":0, "code":0, "freq":0, "sigStr":"L1 C/A"},
            {"sysName":"SBAS", "sys":1, "code":0, "freq":0, "sigStr":"L1 C/A"},
            {"sysName":"GLONASS", "sys":2, "code":0, "freq":0, "sigStr":"L1 C/A"},
            {"sysName":"GALILEO", "sys":3, "code":23, "freq":0, "sigStr":"L1 MBOC(1,1)"},
            {"sysName":"QZSS", "sys":4, "code":0, "freq":0, "sigStr":"L1 C/A"},
            {"sysName":"BDS", "sys":10, "code":26, "freq":6, "sigStr":"B1I"},
            {"sysName":"IRNSS", "sys":9, "code":0, "freq":2, "sigStr":"L5 C/A"} ]

col = ['#1f77b4', '#ff7f0e', '#2ca02c', '#d62728']
monthStr = ['Jan','Feb','Mar','Apr','May','Jun','Jul','Aug','Sep','Oct','Nov','Dec']

#dataType = [{"station":"Beijing",
#             "baseName":"Alloy",
#             "roverName":"NetR9",
#             "baseFileType":"T04",
#             "roverFileType":"T04",
#             "baseType":dFTP, "basePath":"svqa-auto.trimble.com/testing/MON_2020/0142-Alloy_BJ/",
#             "roverType":dFTP, "roverPath":"svqa-auto.trimble.com/testing/MON_2020/0066-NetR9_BJ/"},
#            {"station":"Tokyo",
#             "baseName":"Alloy",
#             "roverName":"NetR9",
#             "baseFileType":"T04",
#             "roverFileType":"T04",
#             "baseType":dFTP, "basePath":"svqa-auto.trimble.com/testing/MON_2020/0132-Alloy/",
#             "roverType":dLocal, "roverPath":"/net/fermion/mnt/data_drive/T01.QZSS/"},
#            {"station":"Terrasat",
#             "baseName":"Alloy",
#             "roverName":"NetR9",
#             "baseFileType":"T02",
#             "roverFileType":"T02",
#             "baseType":dTerrasat, "basePath":"HKRA",
#             "roverType":dTerrasat, "roverPath":"HKR1"},
#            {"station":"Perth",
#             "baseName":"Alloy",
#             "roverName":"NetR9",
#             "baseFileType":"T02",
#             "roverFileType":"T02",
#             "baseType":dTerrasat, "basePath":"STIR",
#             "roverType":dTerrasat, "roverPath":"STI9"},
#            {"station":"Sunnyvale",
#             "baseName":"BD940(2)",
#             "roverName":"BD970(2)",
#             "baseFileType":"T04",
#             "roverFileType":"T02",
#             "baseType":dNoPi, "basePath":"BD940_2_ROVER_",
#             "roverType":dNoPi, "roverPath":"NOV2_ROVER_" },
#            {"station":"Terrasat2",
#             "baseName":"Alloy2",
#             "roverName":"NetR9",
#             "baseFileType":"T04",
#             "roverFileType":"T02",
#             "baseType":dLocal, "basePath":"/net/donald/mnt/data_drive/wlentz/T01_Munich_Alloy/",
#             "roverType":dTerrasat, "roverPath":"HKR1"} ]
#             

with open('station.json', 'r') as fid:
  stationData = json.load(fid)

print(stationData)

baseName  = stationData['baseName']
roverName = stationData['roverName']

now = datetime.datetime.utcnow()
yesterday = now - datetime.timedelta(days=1)

year  = yesterday.year
month = yesterday.month
day   = yesterday.day

print(year,month,day)
startDate = datetime.datetime(year,1,1)

#for month in range(3,7):
  #for day in range(1,32):
if(True):
  if(True):
    DateStr = str(year) + str(month).zfill(2) + str(day).zfill(2)
 
    # Base
    if(stationData['baseType'] == "dFTP"):
      baseFileFull = str(year) + '_' + str(month).zfill(2) + '_' + str(day).zfill(2) + '.' + stationData['baseFileType']
      os.system('wget ftp://'+ stationData['basePath'] + str(year) + '/' + str(month).zfill(2) + '/' + baseFileFull)
      print(baseFileFull)
      os.system('mv ' + baseFileFull + ' base.T04')
      baseFileFull = 'base.T04'
    elif(stationData['baseType'] == "dLocal"):
      baseFileFull = 'Alloy.T04'
      os.system('rm ' + baseFileFull)
      os.system('cat ' + stationData['basePath'] + '*' + str(year) +  str(month).zfill(2) + str(day).zfill(2) + '*.' + stationData['baseFileType'] + ' >> ' + baseFileFull)
      print(baseFileFull)
    elif(stationData['baseType'] == "dTerrasat"):
      directory = "/home/stuart/det-svr/RefData." + str(year-2000) + "/Month." + monthStr[month-1] + "/Day." + str(day).zfill(2) + '/'
      dayNum = int((datetime.datetime(year,month,day) - startDate).total_seconds()/86400) + 1
      station = stationData['basePath']
      baseFileFull = station + str(dayNum).zfill(3) + "0.T02"

      # The script reported an error message that the file was deleted as it had changed during
      # the copy. Therefore, add a retry loop. Give up after a few trys, code later in the file
      # will prevent a fatal if the data is missing.
      for trys in range(10):
        if(os.path.isfile(baseFileFull)==False):
          if(trys == 9):
            time.sleep(3600) # wait an hour for the system to stabilize
          elif(trys >= 1):
            time.sleep(30)
          os.system('cp ' + directory + baseFileFull + ' .')

      print(baseFileFull)
    elif(stationData['baseType'] == "dNoPi"):
      rootDir = '/net/higgs/mnt/data_drive/GNSSResults/' + str(year) + '/' 
      DateStr = str(year) + str(month).zfill(2) + str(day).zfill(2)
      baseFile = stationData['basePath'] + DateStr + '.' + stationData['baseFileType']
      baseFileFull = rootDir + DateStr + '/' + baseFile
      print(baseFileFull)

    # Rover
    if(stationData['roverType'] == "dFTP"):
      roverFileFull = str(year) + '_' + str(month).zfill(2) + '_' + str(day).zfill(2) + '.' + stationData['roverFileType']
      os.system('wget ftp://'+ stationData['roverPath'] + str(year) + '/' + str(month).zfill(2) + '/' + roverFileFull)
      print(roverFileFull)
      os.system('mv ' + roverFileFull + ' rover.T04')
      roverFileFull = 'rover.T04'
    elif(stationData['roverType'] == "dLocal"):
      roverFileFull = 'NetR9.T04'
      os.system('rm ' + roverFileFull)
      os.system('cat ' + stationData['roverPath'] + '*' + str(year) +  str(month).zfill(2) + str(day).zfill(2) + '*.' + stationData['roverFileType'] + ' >> ' + roverFileFull)
      print(roverFileFull)
    elif(stationData['roverType'] == "dTerrasat"):
      directory = "/home/stuart/det-svr/RefData." + str(year-2000) + "/Month." + monthStr[month-1] + "/Day." + str(day).zfill(2) + '/'
      dayNum = int((datetime.datetime(year,month,day) - startDate).total_seconds()/86400) + 1
      station = stationData['roverPath']
      roverFileFull = station + str(dayNum).zfill(3) + "0.T02"
      # We reprocess the Munich NetR9 twice, see if the T02 exists in the alternate directory
      roverFileAlt = '../Alloy-NetR9-Terrasat/'  + station + str(dayNum).zfill(3) + "0.T02"
      if((os.path.isfile(roverFileFull)==False) and (os.path.isfile(roverFileAlt)==False)):
        # The script reported an error message that the file was deleted as it had changed during
        # the copy. Therefore, add a retry loop. Give up after a few trys, code later in the file
        # will prevent a fatal if the data is missing.
        for trys in range(10):
          if(os.path.isfile(roverFileFull)==False):
            if(trys == 9):
              time.sleep(3600) # wait an hour for the system to stabilize
            elif(trys >= 1):
              time.sleep(30)
            os.system('cp ' + directory + roverFileFull + ' .')
      elif((os.path.isfile(roverFileFull)==False) and (os.path.isfile(roverFileAlt)==True)):
        roverFileFull = roverFileAlt
      print(roverFileFull)
    elif(stationData['roverType'] == "dNoPi"):
      rootDir = '/net/higgs/mnt/data_drive/GNSSResults/' + str(year) + '/' 
      DateStr = str(year) + str(month).zfill(2) + str(day).zfill(2)
      roverFile = stationData['roverPath'] + DateStr + '.' + stationData['roverFileType']
      roverFileFull = rootDir + DateStr + '/' + roverFile
      print(roverFileFull)

    if(    os.path.isfile(baseFileFull) 
       and (os.path.getsize(baseFileFull) > 0)
       and os.path.isfile(roverFileFull) 
       and (os.path.getsize(roverFileFull) > 0) ):
      # Read the first file - the "base"
      if(stationData['baseFileType'] == "T04"):
        (baseTmp,k)  =m.vd2arr(baseFileFull,rec='-d35:19',opt=['--dec=1000'])
        i = np.where(baseTmp[:,k.ANTENNA] == 0)[0]
        baseTmp = baseTmp[i,:]
        # Drop data to reduce the memory footprint
        base = baseTmp[:,[k.TIME,k.SAT_TYPE,k.SV,k.FREQ,k.TRACK,k.CNO,k.EL,k.CSTATE]]
        # Create an updated key
        keyBase = m.dotdict({'TIME':0, 'SAT_TYPE':1, 'SV':2, 'FREQ':3, 'TRACK':4, 'CNO':5, 'EL':6, 'CSTATE':7})
      else:
        (baseTmp,k)  =m.vd2arr(baseFileFull,rec='-d27',opt=['--dec=1000'])
        i = np.where(baseTmp[:,k.ANTENNA] == 0)[0]
        baseTmp = baseTmp[i,:]
        base = baseTmp[:,[k.TIME,k.SAT_TYPE,k.SV,k.FREQ,k.TRACK,k.CNO,k.EL]]
        keyBase = m.dotdict({'TIME':0, 'SAT_TYPE':1, 'SV':2, 'FREQ':3, 'TRACK':4, 'CNO':5, 'EL':6})
      del baseTmp
      
      # Read the second file - the "rover"
      if(stationData['roverFileType'] == "T04"):
        (roverTmp,k) =m.vd2arr(roverFileFull,rec='-d35:19',opt=['--dec=1000'])
        i = np.where(roverTmp[:,k.ANTENNA] == 0)[0]
        roverTmp = roverTmp[i,:]
        rover = roverTmp[:,[k.TIME,k.SAT_TYPE,k.SV,k.FREQ,k.TRACK,k.CNO,k.EL,k.CSTATE]]
        keyRover = m.dotdict({'TIME':0, 'SAT_TYPE':1, 'SV':2, 'FREQ':3, 'TRACK':4, 'CNO':5, 'EL':6, 'CSTATE':7})
      else:
        (roverTmp,k) =m.vd2arr(roverFileFull,rec='-d27',opt=['--dec=1000'])
        i = np.where(roverTmp[:,k.ANTENNA] == 0)[0]
        roverTmp = roverTmp[i,:]
        rover = roverTmp[:,[k.TIME,k.SAT_TYPE,k.SV,k.FREQ,k.TRACK,k.CNO,k.EL]]
        keyRover = m.dotdict({'TIME':0, 'SAT_TYPE':1, 'SV':2, 'FREQ':3, 'TRACK':4, 'CNO':5, 'EL':6})
      del roverTmp
      
      psutil.virtual_memory()
      print(dict(psutil.virtual_memory()._asdict()))

      # Now get the base/rover firmware versions

      baseFW  = m.getFirmwareVersion(baseFileFull)
      roverFW = m.getFirmwareVersion(roverFileFull)

      # Get the number of unique epochs, we output this so we can determine whether we
      # have a full day of data.
      baseEpochs  = np.unique( base[:,keyBase.TIME] )
      roverEpochs = np.unique( rover[:,keyRover.TIME] )

      for sys in range(len(sysData)):
        base_i  = find(base[:,keyBase.SAT_TYPE] == sysData[sys]['sys'])
        rover_i = find(rover[:,keyRover.SAT_TYPE] == sysData[sys]['sys'])

        # Data for the current system only
        baseData  = base[base_i,:]
        roverData = rover[rover_i,:]

        svs_base  = np.unique( baseData[:,keyBase.SV] )
        svs_rover = np.unique( roverData[:,keyRover.SV] )

        # Get a unique list of satellites across the two files
        svs = np.unique(np.concatenate([svs_base,svs_rover]))

        for i in range(len(svs)):
          thisSV = int(svs[i])

          print("SV Details %d %d %d %d" % (thisSV, sysData[sys]['sys'], sysData[sys]['freq'], sysData[sys]['code']) )

          # Base Data
          i = find(   (baseData[:,keyBase.SV] == thisSV)
                    & (baseData[:,keyBase.FREQ] == sysData[sys]['freq']) 
                    & (baseData[:,keyBase.TRACK] == sysData[sys]['code']) ) 
          thisBase = baseData[i,:]

          # Rover Data
          i = find(   (roverData[:,keyRover.SV] == thisSV)
                    & (roverData[:,keyRover.FREQ] == sysData[sys]['freq']) 
                    & (roverData[:,keyRover.TRACK] == sysData[sys]['code']) ) 
          thisRover = roverData[i,:]

          _,i_base,i_rover = np.intersect1d( thisBase[:,keyBase.TIME], thisRover[:,keyRover.TIME], return_indices=True)
           
          ######### Find the FLL only data 
          if(stationData['baseFileType'] == "T04"):
            i = find( (thisBase[:,keyBase.CSTATE] < 14) | (thisBase[:,keyBase.CSTATE] >= 21) )
            if(len(i) > 0):
              thisBase_FLL = thisBase[i,:]
            else:
              thisBase_FLL = None
          else:
            thisBase_FLL = None

          if(stationData['roverFileType'] == "T04"):
            i = find( (thisRover[:,keyRover.CSTATE] < 14) | (thisRover[:,keyRover.CSTATE] >= 21) )
            if(len(i) > 0):
              thisRover_FLL = thisRover[i,:]
            else:
              thisRover_FLL = None
          else:
            thisRover_FLL = None
            
          ######### Find the PLL only data 
          if(stationData['baseFileType'] == "T04"):
            i = find( (thisBase[:,keyBase.CSTATE] >= 14) & (thisBase[:,keyBase.CSTATE] < 21) )
            if(len(i) > 0):
              thisBase_PLL = thisBase[i,:]
            else:
              thisBase_PLL = None
          else: # Assume all data is PLL
            thisBase_PLL = thisBase

          if(stationData['roverFileType'] == "T04"):
            i = find( (thisRover[:,keyRover.CSTATE] >= 14) & (thisRover[:,keyRover.CSTATE] < 21) )
            if(len(i) > 0):
              thisRover_PLL = thisRover[i,:]
            else:
              thisRover_PLL = None
          else: # Assume all data is PLL
            thisRover_PLL = thisRover

          if( (len(i_base) > 0) and (len(i_rover) > 0) ):
            fig=figure()
            ax=fig.add_subplot(311)

            if(thisBase_FLL is not None):
              plot( thisBase_FLL[:,keyBase.TIME], thisBase_FLL[:,keyBase.CNO],'.',label=(baseName+' FLL'),color=col[0])
            if(thisBase_PLL is not None):
              plot( thisBase_PLL[:,keyBase.TIME], thisBase_PLL[:,keyBase.CNO],'.',label=(baseName+' PLL'),color=col[1])

            if(thisRover_FLL is not None):
              plot( thisRover_FLL[:,keyRover.TIME], thisRover_FLL[:,keyRover.CNO],'.',label=(roverName+' FLL'),color=col[2])
            if(thisRover_PLL is not None):
              plot( thisRover_PLL[:,keyRover.TIME], thisRover_PLL[:,keyRover.CNO],'.',label=(roverName+' PLL'),color=col[3])

            fid = open("results.txt","a")
            fid.write("%d,%d,%d,%d,%d," % (year,month,day,thisSV,sysData[sys]['sys']))
            fid.write("%d,%d,%d," % (len(thisBase), len(thisRover), len(i_base)))
            delta = thisBase[i_base,keyBase.CNO] - thisRover[i_rover,keyRover.CNO]
            fid.write("%.3f,%.3f," % (np.mean(delta), np.std(delta)))

            # Now find the C/No difference for matching PLL epochs

            if( (thisBase_PLL is not None) and (thisRover_PLL is not None) ):
              _,i_base_PLL,i_rover_PLL = np.intersect1d( thisBase_PLL[:,keyBase.TIME], thisRover_PLL[:,keyRover.TIME], return_indices=True)
              if( (len(i_base_PLL) > 0) and (len(i_rover_PLL) > 0) ):
                delta_PLL = thisBase_PLL[i_base_PLL,keyBase.CNO] - thisRover_PLL[i_rover_PLL,keyRover.CNO]
                fid.write("%.3f,%.3f," % (np.mean(delta_PLL), np.std(delta_PLL)))
              else:
                fid.write("0.0,0.0,")
            else:
              fid.write("0.0,0.0,")

            # Now get the num epochs >= elevThresh degree elevation
            i = find(thisBase[:,keyBase.EL] >= elevThres)
            filteredBase = thisBase[i,:]
            fid.write("%d," % (len(i)))
            i = find(thisRover[:,keyRover.EL] >= elevThres)
            filteredRover = thisRover[i,:]
            fid.write("%d," % (len(i)))

            # The above provided all observables irrespective of tracking state. Now filter
            # the PLL only data.
            # 14 = PLL_A_WIDE
            # 20 = WAAS PLL
            if(stationData['baseFileType'] == "T04"):
              i = find( (filteredBase[:,keyBase.CSTATE] >= 14) & (filteredBase[:,keyBase.CSTATE] < 21) )
              fid.write("%d," % (len(i)))
            else:
              fid.write("%d," % len(filteredBase) )

            if(stationData['roverFileType'] == "T04"):
              i = find( (filteredRover[:,keyRover.CSTATE] >= 14) & (filteredRover[:,keyRover.CSTATE] < 21) )
              fid.write("%d," % (len(i)))
            else:
              fid.write("%d," % len(filteredRover) )

            fid.write("%d,%d," % (len(baseEpochs), len(roverEpochs)))
            
            # Output the number of gaps > 20 seconds
            fid.write("%d,%d," % (len(find(np.diff(baseEpochs)  > 20)),
                                   len(find(np.diff(roverEpochs) > 20))) )  

            fid.write("%s,%s\n" % (baseFW,roverFW))
            fid.close()

            ylabel('C/No [dBHz]')
            title('All Data - ' + str(year)  + '-' + str(month).zfill(2) + '-' + str(day).zfill(2) + ' ' + sysData[sys]['sysName'] + ' ' + str(thisSV) + ' - ' + sysData[sys]['sigStr'] + ' C/No [dBHz]')
            grid(True)
            # Prevent the axis numers having an offset
            ax.get_xaxis().get_major_formatter().set_useOffset(False)
            ax.get_yaxis().get_major_formatter().set_useOffset(False)
            legend(loc='best')
            tight_layout()

            left, right = xlim()

            ax=fig.add_subplot(312)
            plot(thisBase[:,keyBase.TIME], thisBase[:,keyBase.EL],'r.')
            plot(thisRover[:,keyRover.TIME], thisRover[:,keyRover.EL],'r.')
            ylabel('Elev[deg]')
            title('Elevation Angle [deg]')
            grid(True)
            ax.get_xaxis().get_major_formatter().set_useOffset(False)
            ax.get_yaxis().get_major_formatter().set_useOffset(False)
            xlim(left,right)
            tight_layout()

            ax=fig.add_subplot(313)
            label  = r'$\mu$ = ' + str.format('{0:.2f}', np.mean(delta)) + 'm\n' 
            label += r'$\sigma$ = ' + str.format('{0:.2f}', np.std(delta)) + 'm' 
            plot(thisBase[i_base,keyBase.TIME], delta,'.',label=label)
            xlabel('Time [GPS Seconds]');
            ylabel('C/No Diff [dBHz]')
            title('All Data - ' + sysData[sys]['sigStr'] + ' C/No difference (' + baseName + ' - ' + roverName + ') [dBHz]')
            grid(True)
            ax.get_xaxis().get_major_formatter().set_useOffset(False)
            ax.get_yaxis().get_major_formatter().set_useOffset(False)
            xlim(left,right)
            tight_layout()
            legend(loc='best')


            savePlot(DateStr,"/" + DateStr + sysData[sys]['sysName'] + '-' + str(thisSV) + '.png')
            close() # The figure

            #if( (len(i_base_PLL) > 0) and (len(i_rover_PLL) > 0) ):
            #  fig = figure()
            #  ax = fig.add_subplot(211)

            #  if(thisBase_FLL is not None):
            #    plot( thisBase_FLL[:,keyBase.TIME], thisBase_FLL[:,keyBase.CNO],'.',label=(baseName+' FLL'),color=col[0])
            #  if(thisBase_PLL is not None):
            #    plot( thisBase_PLL[:,keyBase.TIME], thisBase_PLL[:,keyBase.CNO],'.',label=(baseName+' PLL'),color=col[1])

            #  if(thisRover_FLL is not None):
            #    plot( thisRover_FLL[:,keyRover.TIME], thisRover_FLL[:,keyRover.CNO],'.',label=(roverName+' FLL'),color=col[2])
            #  if(thisRover_PLL is not None):
            #    plot( thisRover_PLL[:,keyRover.TIME], thisRover_PLL[:,keyRover.CNO],'.',label=(roverName+' PLL'),color=col[3])
            
            #  ylabel('C/No [dBHz]')
            #  title(str(year)  + '-' + str(month).zfill(2) + '-' + str(day).zfill(2) + ' ' + sysData[sys]['sysName'] + ' ' + str(thisSV) + ' - ' + sysData[sys]['sigStr'] + ' C/No [dBHz]')
            #  grid(True)
            #  ax.get_xaxis().get_major_formatter().set_useOffset(False)
            #  ax.get_yaxis().get_major_formatter().set_useOffset(False)
            #  legend()
            #  tight_layout()
            #  left, right = xlim()
              
            #  ax=fig.add_subplot(212)
            #  plot(thisBase[:,keyBase.TIME], thisBase[:,keyBase.EL],'r.')
            #  plot(thisRover[:,keyRover.TIME], thisRover[:,keyRover.EL],'r.')
            #  ylabel('Elev[deg]')
            #  title('Elevation Angle [deg]')
            #  grid(True)
            #  ax.get_xaxis().get_major_formatter().set_useOffset(False)
            #  ax.get_yaxis().get_major_formatter().set_useOffset(False)
            #  xlim(left,right)
            #  xlabel('Time [GPS Seconds]');
            #  tight_layout()
            #  savePlot(DateStr,"/" + DateStr + sysData[sys]['sysName'] + '-' + str(thisSV) + 'PLLFLL.png')
            #  close() # The figure
    else:
      print("Missing File")

    if(stationData['baseType'] == "dFTP"):
      os.system('rm ' + baseFileFull)
    if(stationData['roverType'] == "dFTP"):
      os.system('rm ' + roverFileFull)

    #os.system('convert -delay 100 -loop 0 */*GLONASS-13.png GLONASS13.gif')
    #os.system('convert -delay 100 -loop 0 */*GLONASS-20.png GLONASS20.gif')