# -*- coding: utf-8 -*-
"""
#!/bin/python
# Filename em27_fix_sunrun.py

Created on Wed Sep  2 17:21:40 2015

@author: Jacob Hedelius

version Apr-01-2020

"""

## import some sub-functions
import numpy as np
import re
import os
import pandas as pd
curr_dir = os.getcwd()
os.chdir('../python_subs')
import egi_subs as egs
import readusbgpsfiles as rgps
import read_timesig_fb as rtime
import findclosest_n as fclsn
import readzenofiles as rzeno
import read_hdr
os.chdir(curr_dir)





def em27_fix_sunrun2020(fname_pre, sunrundat=None):
    """
    Fixes sunrun files for EGI 2020. I've modified things so that the met and 
    location matching is done before in metmatch_sub.py.
    I'm assuming filenames based on an input path name.
    All that's left to do is just read in the existing sunrun and fix it
    where needed.
    
    example input
    fname_pre = '/home/osboxes/research/tccon/egi2020/results/cn_spec/daily/20140609/matchup/cn20140609'
    sunrundat = '/home/osboxes/research/tccon/ggg_devel/tccon/cn_sunrun.dat'
    """
    
    #### Set paths to input files
    mm_fname = fname_pre + '_mmtch.dat'
    day_dir = os.path.split(os.path.split(fname_pre)[0])[0] 
    basesunrun = day_dir + '/runlog/' + fname_pre[-10:] + '_unfixed.gop' 
    outsunrun = day_dir + '/runlog/' + fname_pre[-10:] + '.gop' 
    
    ####### Read in files...twice for sunruns...
    with open(basesunrun, 'r') as fobj: # gets file object for reading
      origtext = fobj.readlines() #reads in entire file all at once, line by line
        
    sr_file = egs.read_sunrun(basesunrun) #sunrun read into dictionary of lists
    
    hdr_raw = read_hdr.read_hdr(fname_pre + '_hdr1.dat')
    if 'TIM1' in hdr_raw:
      ifg_utc_off = float(hdr_raw['TIM1'].split('GMT')[1].split(')')[0]) #E.g., Harvard=-4
    else: #likely I2S files, NOTE: changed "GMT" to "UTC" as well
      ifg_utc_off = float(hdr_raw['TIM'].split('UTC')[1].split(')')[0]) #E.g., Harvard=-4
    time_corr_sec = np.array(sr_file['tcorr'], dtype=np.float) - ifg_utc_off*3600
    
    
    # Find the column header line by looking for line beginning with "Iflag"
    with open(mm_fname, 'r') as fobj: # gets file object for reading
      mm_raw = fobj.readlines() 
    n_mm_hdr = next(i1 for i1, val in enumerate(mm_raw) if val[0:5] == 'Iflag')-1
    
    mm_df = pd.read_csv(mm_fname, delim_whitespace=True, header=n_mm_hdr) 
    mm_noflag = mm_df.loc[mm_df.loc[:,'Iflag']==0,:] #get rid of flagged rows
    
    ## read in xx_sunrun.dat file, get pressure correction
    pcorr = 0
    if sunrundat is not None:
        f_sr = open(sunrundat, 'r') # gets file object for reading
        sunrundat_cell = f_sr.readlines() #reads in entire file all at once, line by line
        f_sr.close()
        pline = re.split('[!]+', sunrundat_cell[12]) #splits line
        pcorr = float(pline[0])
    
    
    
    
    
    ########  format output lines except header
    nSR_hdr = 3 #yes hardcoded, I'm being lazy
    
    outtext5 = (len(origtext)-nSR_hdr)*[None]
    for i1 in range(0, len(origtext) - nSR_hdr): #writing line-by-line into a temporary list
        outline_i1 = origtext[i1 + nSR_hdr][0:62]
        outline_i1 = outline_i1 + format(time_corr_sec[i1], '6.0f')
        
        outline_i1 = outline_i1 + '. ' #py col 67-68 
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'lat'], '8.4f') #lat, 8 char
        outline_i1 = outline_i1 + ' ' #py col 77
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'lon'], '9.4f') #long, 9 char
        outline_i1 = outline_i1 + '  ' #col 87-88
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'alt']/1000, '5.3f') #alt, 5 char
        outline_i1 = outline_i1 + ' ' #col 94
        
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'Tins'], '5.1f') #temperature, 5 characters
        outline_i1 = outline_i1 + ' ' #col 100
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'Pout'], '7.2f') #pressure, 7 characters
        outline_i1 = outline_i1 + ' ' #col 108   
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'RH'], '5.1f') #RH, 5 char
        outline_i1 = outline_i1 + ' ' #col 114
        
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'Tout'], '5.1f') #tout
        outline_i1 = outline_i1 + ' ' #col 120
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'Pout'] + pcorr, '7.2f') #pout
        outline_i1 = outline_i1 + ' ' #col 128
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'RH'], '5.1f') #hout
        outline_i1 = outline_i1 + ' ' #col 134
        
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'sia'], '6.1f') #SIA, 6 char
        
        outline_i1 = outline_i1 + origtext[i1 + nSR_hdr][142:150] #actual columns 143-150. FVSI
        
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'WSPD'], '5.1f') #WSPD, 5 char
        outline_i1 = outline_i1 + ' ' #col 154
        outline_i1 = outline_i1 + format(mm_noflag.loc[i1, 'WDIR'], '4.0f') #WDIR, 4 char
        
        outline_i1 = outline_i1 + origtext[i1 + nSR_hdr][160:]
        
        outtext5[i1] = outline_i1
        
    
    
    
    ########### Actual file writing 
    with open(outsunrun, 'w') as fileID: #open file for writing
        #write header lines
        i1 = 0
        while i1<nSR_hdr:
            fileID.write(origtext[i1])
            i1 += 1
    
        #write all other lines (data lines)
        for outline_i1 in outtext5: #over each ifg name
            fileID.write(outline_i1)
    
    
    
    out = outsunrun
    return out
    
    



def em27_fix_sunrun(basesunrun, spec_timesigfile, zenofilename, location_file, sunrundat, time_corr):
        
    #This is script is meant to fix sunrun files, it adds in the necessary
    #weather data where needed.  
    #
    #Pre-operations needed include: (1) get a Zeno weather file
    #(2) Extract time headers using the read_timesig_fb script. 
    #(3) Create the base sunrun without the weather data.  You will need a
    #.gnd file first to make the sunrun.  Go to the TCCON-wiki page under
    #sunrun for instructions on creating this list of spectra.
    #Note: by running the shell script egiopusrldaily, all this should be
    #automatically done
    #
    #NOTE: sunruns are particular to the number of spaces, this file keeps
    #them the exact same at the expense of some versatility.  If there is ever
    #a change, this file will need to be changed manually.
    
    ## Test values and code
    
#    basesunrun = 'C:\\Users\\Jacob Hedelius\\Documents\\Research\\python\\example_data\\cn20141005uf.gop'
#    spec_timesigfile = 'C:\\Users\\Jacob Hedelius\\Documents\\Research\\python\\example_data\\cnspectra_sigtime20141005.dat' 
#    zenofilename = 'C:\\Users\\Jacob Hedelius\\Documents\\Research\\python\\example_data\\cn_cmb_zeno.zno'
#    sunrundat = 'C:\\Users\\Jacob Hedelius\\Documents\\Research\\python\\example_data\\cn_sunrun.dat';
#    location_file = [34.2, -118.1, 237]
#    time_corr=0
        
    
    
    
    max_timediff = 300 # maximum time difference for location/weather files (seconds)
    
    
    ## Open, read in outline file
    fobj = open(basesunrun, 'r') # gets file object for reading
    origtext = fobj.readlines() #reads in entire file all at once, line by line
    fobj.close()
    
    specname_lines = np.ones((len(origtext), 1), dtype=bool)
    rl_head = 3
    specname_lines[0:rl_head] = np.zeros((rl_head,1,), dtype=bool) #First 3 lines should always be false if the format remains the same
    
    outtext = origtext # copy the original text
    
    time_fix2utc = time_corr
    

    ## set longitude
    if not(location_file == 0 or egs.is_number(location_file[0])): # first character not a number, filepath given
        gps_loc_raw = rgps.readusbgpsfiles(location_file)
        if 'CompSrlDate' in gps_loc_raw:
            gps_time = gps_loc_raw['CompSrlDate'] #serial date (matlab version) exists.  Same time as ifgs. 
        elif gps_loc_raw['utc_lcl'] == 2: #local time
            gps_time_p = egs.dt_from6vec(gps_loc_raw['Cyear'],gps_loc_raw['Cmonth'],gps_loc_raw['Cday'],gps_loc_raw['Chour'],gps_loc_raw['Cmin'],gps_loc_raw['Csec']) #matches spectrometer/computer time
            gps_time = gps_time_p #+ 719529 #converts to matlab type
        elif gps_loc_raw['utc_lcl'] == 1: #utc time
            gps_time_p = egs.dt_from6vec(gps_loc_raw['Uyear'],gps_loc_raw['Umonth'],gps_loc_raw['Uday'],gps_loc_raw['Uhour'],gps_loc_raw['Umin'],gps_loc_raw['Usec']) - time_fix2utc/24 #uses UTC time
            gps_time = gps_time_p # with new function already in to matlab type 
        
        
        long_in = np.nanmedian(gps_loc_raw['Long']) #manually sets first longitude
        
    elif location_file == 0: #Default to Caltech, already formatted
        long_in = -118.1 
    
    elif len(location_file) == 3: #lat, long, alt_m given
        long_in = location_file[1]
    
    else: #all others default to Caltech
        long_in = -118.1 
    
    ## read in time/sig file
    time_sig = rtime.read_timesig_fb(spec_timesigfile, long_in, time_fix2utc) #MATLAB ~4.4 s, python version now ~1 s updated egi_subs to make 150x faster :)
    n_intfgs = len(time_sig['intfg_name'])
    
    ## replace location data with correct location
    if not(location_file == 0 or egs.is_number(location_file[0])): # first character not a number, filepath given
        pass        
        #gps_avg_corr_data[0,1] = gps_loc_raw.Long[0] # may be unused
        
    elif location_file == 0: #Default to Caltech, already formatted
        gps_avg_corr_data = np.outer( np.ones([n_intfgs,1]), [34.13625, -118.12697, 237])
    
    elif len(location_file) == 3: #lat, long, alt_m given
        gps_avg_corr_data = np.outer( np.ones([n_intfgs,1]), location_file)
    
    else: #all others default to Caltech
        gps_avg_corr_data = np.outer( np.ones([n_intfgs,1]), [34.13625, -118.12697, 237])
    
    ## finish GPS cooridnation if needed
    filedates = time_sig['raw_time'] # + 719529 #converts to matlab type 
    
    if not(location_file == 0 or egs.is_number(location_file[0])): 
        gps_closedate_p = fclsn.findclosest_n(filedates, gps_time, 3) #up to 3 points
        gps_closedate = gps_closedate_p[0]
        gpsind_win_min = abs(gps_closedate['closest_indx_diff']) <= float(1)/(24*60) # %indexes where data was taken within one minute
        
        # gps_avg_corr_data_orig = gps_avg_corr_data; %original values, without nans - UNUSED
        
        if 'masl' in gps_loc_raw:
            alt_all = gps_loc_raw['masl'] 
        else:
            alt_all = gps_loc_raw['masl_inacc'] + gps_loc_raw['geoidal_WGS84_sep'] + 18; #the +18 seems like it would make the measurement more accurate (from an empirical observation at Caltech)
        
        i1 = 0
        for fdate_i1 in filedates:
            gps_data_indx_i1_p = gps_closedate['closest_indx'][i1, :]
            gps_data_indx_i1 = gps_data_indx_i1_p[gpsind_win_min[i1,:]]
            # alt_i1 = gps_loc_raw.data(gps_data_indx_i1, 8) + gps_loc_raw.data(gps_data_indx_i1, 9);
             
            gps_avg_corr_data[i1,0] = np.nanmean( gps_loc_raw['Lat'][gps_data_indx_i1])
            gps_avg_corr_data[i1,1] = np.nanmean(gps_loc_raw['Long'][gps_data_indx_i1])
            gps_avg_corr_data[i1,2] = np.nanmean(alt_all[gps_data_indx_i1])
            
            i1 += 1 #iterates by 1
    
    
    
    if np.mean(gps_avg_corr_data[:,2])>=10:
        gps_avg_corr_data[:,2] = gps_avg_corr_data[:,2]/1000 #km to meters
    
    ## read in spectra time files
    spec_names = time_sig['intfg_name']
    
    ## read in xx_sunrun.dat file, get pressure correction
    f_sr = open(sunrundat, 'r') # gets file object for reading
    sunrundat_cell = f_sr.readlines() #reads in entire file all at once, line by line
    f_sr.close()
    pline = re.split('[!]+', sunrundat_cell[12]) #splits line
    pcorr = float(pline[0])
    
        
    ## Zeno data read in and average and coordinate with times interferograms were taken
    zenofraw = rzeno.readzenofiles(zenofilename)
    
    if 'CompSrlDate' in zenofraw:
        zeno_date = zenofraw['CompSrlDate'] #serial date (matlab version) exists.  Same time as ifgs. 
    elif zenofraw['utc_lcl'] == 2: #local time
        zeno_date_p = egs.dt_from6vec(zenofraw['Cyear'],zenofraw['Cmonth'],zenofraw['Cday'],zenofraw['Chour'],zenofraw['Cmin'],zenofraw['Csec']) #matches spectrometer/computer time
        zeno_date = zeno_date_p #+ 719529 #converts to matlab type
    elif zenofraw['utc_lcl'] == 1: #utc time
        zeno_date_p = egs.dt_from6vec(zenofraw['Uyear'],zenofraw['Umonth'],zenofraw['Uday'],zenofraw['Uhour'],zenofraw['Umin'],zenofraw['Usec']) - time_fix2utc/24 #matches spectrometer/computer time
        zeno_date = zeno_date_p #+ 719529 #converts to matlab type
    
        
    # add in -9 if no WSPD and WDIR data
    if not('WSPD' in zenofraw):
        zenofraw['WSPD'] = np.ones([len(zeno_date),1])*0
    
    if not('WDIR' in zenofraw):
        zenofraw['WDIR'] = np.ones([len(zeno_date),1])*0
    
    if not('RH' in zenofraw): #if no RH data in the met files
        zenofraw['RH'] = np.ones([len(zeno_date),1])*0
    
    zeno_closedate_p = fclsn.findclosest_n(filedates, zeno_date, 5) #up to 5 points
    zeno_closedate = zeno_closedate_p[0]
    zenoind_win_m = abs(zeno_closedate['closest_ind_diff']) <= float(max_timediff)/(24*60*60) # %indexes where data was taken within accepted time. Needs float or Python 2 rounds.
    
    zeno_avg_corr_data = -999*np.ones([len(filedates),6]) #data from zeno that goes to input file includes temp, pres, RH, winddir, windspd - pre-allocating to see what's been written
    
    if 'Temp' in zenofraw:
        zeno_temp = zenofraw['Temp']
    elif 'Tout' in zenofraw:
        zeno_temp = zenofraw['Tout']
    else: #no Temperature data
        zeno_temp = np.ones([len(zeno_date),1])*0
    
    
    i1 = 0
    for fdate_i1 in filedates: #over each ifg/timepoint
        zeno_data_indx_i1_p2 = zeno_closedate['closest_ind'][i1, :]
        zeno_data_indx_i1_p1 = zeno_data_indx_i1_p2[zenoind_win_m[i1,:]]
        zeno_data_indx_i1 = zeno_data_indx_i1_p1.astype(int)
         
        zeno_avg_corr_data[i1,0] = egs.nanmean(zeno_temp[zeno_data_indx_i1])
        zeno_avg_corr_data[i1,1] = egs.nanmean(zenofraw['RH'][zeno_data_indx_i1])
        zeno_avg_corr_data[i1,2] = egs.nanmean(zenofraw['Pout'][zeno_data_indx_i1])
        wdirspd_avg = egs.avg_wdir(zenofraw['WDIR'][zeno_data_indx_i1], zenofraw['WSPD'][zeno_data_indx_i1])
        zeno_avg_corr_data[i1,3] = wdirspd_avg[1] #
        zeno_avg_corr_data[i1,4] = wdirspd_avg[0]
        
        i1 += 1 #iterates by 1
    
    zeno_avg_corr_data[:,0] = np.round(zeno_avg_corr_data[:,0],1) #Tout
    zeno_avg_corr_data[:,1] = np.round(zeno_avg_corr_data[:,1]) #RH
    zeno_avg_corr_data[:,2] = np.round(zeno_avg_corr_data[:,2],1) #Pins
    zeno_avg_corr_data[:,3] = np.round(zeno_avg_corr_data[:,3],1) #WSPD
    zeno_avg_corr_data[:,4] = np.round(zeno_avg_corr_data[:,4]) #WDIR
    #zeno_avg_corr_data[:,5] = egs.nan_ary([len(zeno_avg_corr_data[:,4]),1])
    
    n_a = np.isnan(zeno_avg_corr_data)
    zeno_avg_corr_data[n_a] = -999
    
    
    ## TODO: This is where the internal EM27 data on pressure, humidity, etc. could be read in
    #NONE IS AVALIABLE
    
    ## Replace false weather data with Zeno data
    sr_file = egs.read_sunrun(basesunrun) #sunrun read into dictionary of lists
    
    zeno_avg_corr_data[zeno_avg_corr_data[:,3]>100,3] = 99 # Replaces wind speeds that are over 100 mps with 99
    zeno_avg_corr_data[:,5] = zeno_avg_corr_data[:,2] + pcorr #includes pressure correction for pout
    sr_file['tins'] = list(time_sig['TSC']) #scanner temperature for tins
    sr_file['pins'] = zeno_avg_corr_data[:,2]
    sr_file['hins'] = zeno_avg_corr_data[:,1]
    sr_file['tout'] = zeno_avg_corr_data[:,0]
    sr_file['pout'] = zeno_avg_corr_data[:,5]
    sr_file['hout'] = zeno_avg_corr_data[:,1]
    sr_file['wspd'] = zeno_avg_corr_data[:,3]
    sr_file['wdir'] = zeno_avg_corr_data[:,4]
    
    ## replace GPS
    sr_file['oblat'] = gps_avg_corr_data[:,0]
    sr_file['oblon'] = gps_avg_corr_data[:,1]
    sr_file['obalt'] = gps_avg_corr_data[:,2]
    
    ## Arbitrary SIA
    sr_file['sia'] = list(time_sig['arb_sia']) 
    
    ## Time Correction to UTC (besides what's in xx_sunrun.dat)
    time_corr_sec = np.array(sr_file['tcorr'],dtype=np.float) + time_corr*3600
    sr_file['tcorr'] = time_corr_sec
    
    # NOTE: the values here are spacing sensitive, and manually input
    # need to be in appropriate format of strings, in terms of leading spaces 
    # and trailing zeros
    
    outtext5 = origtext[rl_head:]
    i1 = 0
    for inline_i1 in outtext5: #writing line-by-line
        outline_i1 = inline_i1[0:62]
        outline_i1 = outline_i1 + format(sr_file['tcorr'][i1], '6.0f')
        
        outline_i1 = outline_i1 + '. ' #py col 67-68 
        outline_i1 = outline_i1 + format(sr_file['oblat'][i1], '8.4f') #lat, 8 char
        outline_i1 = outline_i1 + ' ' #py col 77
        outline_i1 = outline_i1 + format(sr_file['oblon'][i1], '9.4f') #long, 9 char
        outline_i1 = outline_i1 + '  ' #col 87-88
        outline_i1 = outline_i1 + format(sr_file['obalt'][i1], '5.3f') #alt, 5 char
        outline_i1 = outline_i1 + ' ' #col 94
        
        outline_i1 = outline_i1 + format(sr_file['tins'][i1], '5.1f') #temperature, 5 characters
        outline_i1 = outline_i1 + ' ' #col 100
        outline_i1 = outline_i1 + format(sr_file['pins'][i1], '7.2f') #pressure, 7 characters
        outline_i1 = outline_i1 + ' ' #col 108   
        outline_i1 = outline_i1 + format(sr_file['hins'][i1], '5.1f') #RH, 5 char
        outline_i1 = outline_i1 + ' ' #col 114
        
        outline_i1 = outline_i1 + format(sr_file['tout'][i1], '5.1f') #tout
        outline_i1 = outline_i1 + ' ' #col 120
        outline_i1 = outline_i1 + format(sr_file['pout'][i1], '7.2f') #pout
        outline_i1 = outline_i1 + ' ' #col 128
        outline_i1 = outline_i1 + format(sr_file['hout'][i1], '5.1f') #hout
        outline_i1 = outline_i1 + ' ' #col 134
        
        outline_i1 = outline_i1 + format(sr_file['sia'][i1], '6.1f') #SIA, 6 char
        
        outline_i1 = outline_i1 + inline_i1[142:150] #actual columns 143-150. FVSI
        
        outline_i1 = outline_i1 + format(sr_file['wspd'][i1], '5.1f') #WSPD, 5 char
        outline_i1 = outline_i1 + ' ' #col 154
        outline_i1 = outline_i1 + format(sr_file['wdir'][i1], '4.0f') #WDIR, 4 char
        
        outline_i1 = outline_i1 + inline_i1[160:]
        
        outtext5[i1] = outline_i1
        i1 +=1 #advance the loop
    
    
    ## write the file
    filenameout = basesunrun[0:-6] + basesunrun[-4:] #removes the "uf" from the name
    
    fileID = open(filenameout, 'w') #open file for writing
    
    #write header lines
    i1 = 0
    while i1<rl_head:
        fileID.write(origtext[i1])
        i1 += 1
    
    #write all other lines (data lines)
    for outline_i1 in outtext5: #over each ifg name
        fileID.write(outline_i1)
    
    fileID.close()
    
    out = {'filename': filenameout}
    return out

