import os
import time
import numpy
import scipy
import base64
from io import BytesIO
import scipy.interpolate
from matplotlib.figure import Figure

from utils import TimeTools

def skyNoise(jd, ut=-5.0, longitude=-76.87, filename='/app/utils/galaxy.txt'):
    """
    hydrapos returns RA and Dec provided by Bill Coles (Oct 2003).
    
    Parameters
    ----------
    jd = The julian date of the day (and time), scalar or vector.
    
    dec_cut = A scalar giving  the declination  to get a cut of the skynoise over Jica-
        marca. The default value is -11.95.
    filename = A string to specify name the skynoise file. The default value is skynoi-
        se_jro.dat
    
    Return
    ------
    maxra = The maximum right ascension to the declination used to get a cut.
    ra = The right ascension.
    
    Examples
        --------
    >> [maxra,ra] = skynoise_jro()
    >> print maxra, ra
    139.45 -12.0833333333
    
    Modification history
    --------------------             
    Converted to Python by Freddy R. Galindo, ROJ, 06 October 2009.
    """

    # Defining date to compute SkyNoise.
    [year, month, dom, hour, mis, secs] = TimeTools.Julian(jd).change2time()
    is_dom = (month==9) & (dom==21)
    if is_dom:
        tmp = jd
        jd = TimeTools.Time(year,9,22).change2julian()
        dom = 22

    # Reading SkyNoise        
    
    g = os.getcwd()

    f = open(filename)
    
    lines = f.read()
    f.close()
    
    nlines = 99
    lines = lines.split('\n')
    data = numpy.zeros((2,nlines))*numpy.float32(0.)
    for ii in numpy.arange(nlines):
        line = numpy.array([lines[ii][0:6],lines[ii][6:]])
        data[:,ii] = numpy.float32(line)

    # Getting SkyNoise to the date desired.
    otime = data[0,:]*60.0
    opowr = data[1,:]
    
    hour = numpy.array([0,23])
    mins = numpy.array([0,59])
    secs = numpy.array([0,59])
    LTrange = TimeTools.Time(year,month,dom,hour,mins,secs).change2julday()
    LTtime  = LTrange[0] + numpy.arange(1440)*((LTrange[1] - LTrange[0])/(1440.-1))
    lst = TimeTools.Julian(LTtime + (-3600.*ut/86400.)).change2lst()

    ipowr = lst*0.0
    # Interpolating using scipy (inside max and min "x")
    otime = otime/3600.
    val = numpy.where((lst>numpy.min(otime)) & (lst<numpy.max(otime))); val = val[0]
    tck = scipy.interpolate.interp1d(otime,opowr)
    ipowr[val] = tck(lst[val])

    # Extrapolating above maximum time data (23.75).
    uval = numpy.where(lst>numpy.max(otime))
    if uval[0].size>0:
        ii = numpy.min(uval[0])
        m = (ipowr[ii-1] - ipowr[ii-2])/(lst[ii-1] - lst[ii-2])
        b = ipowr[ii-1] - m*lst[ii-1]            
        ipowr[uval] = m*lst[uval] + b

    if is_dom:
        lst = numpy.roll(lst,4)        
        ipowr = numpy.roll(ipowr,4)

    new_lst = numpy.int32(lst*3600.)
    new_pow = ipowr

    return ipowr, LTtime, lst

def skynoise_plot(year, month, day):
    """
    getPlot method creates a skynoise map over Jicamarca for a desired day. Additionally
    save a PNG file of this plot.
    
    Examples
    --------
    >> SkyNoisePlot(skypower,skytime,skytime_lst).getPlot()
    
    Modification history
    --------------------             
    Created by Freddy R. Galindo, ROJ, 18 October 2009.
    """

    # Working with the time before to plot the SkyNoise

    julian = TimeTools.Time(year, month, day).change2julday()
    power, tm, time_lst = skyNoise(julian)

    secs = TimeTools.Julian(tm).change2secs()    
    secs_lst = time_lst*1.
    if time_lst.argmin()>0:
        secs_lst[time_lst.argmin():] = secs_lst[time_lst.argmin():] + 24.
        secs_lst = secs_lst*3600. 

    label_secs = time.localtime(secs[power.argmax()] + time.timezone)
    label_secs_lst = time.localtime(secs_lst[power.argmax()] + time.timezone)
    
    # Creating labels for main x-labelticks (Local Time):
    snow = TimeTools.Time(year, month, day, 0, 0, 0).change2secs()
    today = secs - snow
    xticks_dpos = []
    xticks_dval = []
    for ix in [0,120,240,360,480,600,720,840,960,1080,1200,1320,1439]:             
        xticks_dpos.append(today[ix])
        time_tuple = time.localtime(today[ix] + time.timezone)
        xticks_dval.append(time.strftime('%H:%M',time_tuple))
        if ix==1439:xticks_dval[12] = ''

    # Creating labels for secondary x-labelticks (Sidereal Time):
    xticks_upos = [secs_lst[0],secs_lst[-1]]
    xticks_uval = ['','']
    for ix in [0,120,240,360,480,600,720,840,960,1080,1200,1320]:               
        index = numpy.argmin(numpy.abs(time_lst - (ix/60. + 1.)))
        xticks_upos.append(secs_lst[index])
        time_tuple = time.localtime(secs_lst[index] + time.timezone)
        if time_tuple[4]==59:
            tmp_list = list(time_tuple)
            tmp_list[4] = 0
            tmp_list[3] = tmp_list[3] + 1
            time_tuple = tuple(tmp_list)            
        xticks_uval.append(time.strftime('%H:%M',time_tuple))

    # Creating SkyNoise Map.
    fig = Figure()
    fig.subplots_adjust(top=0.80)
    
    main_mesg = "SKY BRIGHTNESS AT 50Mhz - Date: " 
    main_mesg = main_mesg + time.strftime("%d-%b-%Y (%j)",label_secs)
    main_mesg = main_mesg + "\nGalaxy Pass at " + time.strftime("%H:%M:%S LT",label_secs)
    main_mesg = main_mesg + time.strftime("  (%H:%M:%S LST)",label_secs_lst)
    
    fig.suptitle(main_mesg, fontsize=12)
    axl = fig.add_subplot(1,1,1)
    
    axl.plot(today,power,'b-')       
    axr = axl.twinx()        
    axl.set_xlabel('Local Time (Hours)')
    axl.set_ylabel('Signal Strengh (mA)')
    axr.set_ylabel('Temperature [K]x10^3')
    axl.set_xlim(0,24)
    axl.set_ylim(0,115)
    axr.set_ylim(0,50)
    axl.set_xticks(xticks_dpos)
    axl.set_xticklabels(xticks_dval, fontsize=8)
    axl.grid()

    axt = axl.twiny()
    axt.set_xlabel('Local Sidereal Time (Hours)')
    axt.set_xlim(secs_lst[0],secs_lst[-1])
    axt.set_xticks(xticks_upos)
    axt.set_xticklabels(xticks_uval, fontsize=8)
    
    buf = BytesIO()
    fig.savefig(buf, format="png")
    
    return buf