#################################################################################
from pathlib import Path
import os
import sys
from astropy.stats import mad_std
from astropy.io import fits
from astropy.utils.data import get_pkg_data_filename
import glob, os
from astropy import stats
import matplotlib.pyplot as plt
import numpy as np
from datetime import date
today = date.today()
from astropy.convolution import convolve, Box1DKernel
from astropy.convolution import Box2DKernel
date = today.strftime("%b-%d-%Y")
from astropy.wcs import WCS
import warnings
warnings.filterwarnings('ignore')
from matplotlib.backends.backend_pdf import PdfPages
import shutil
#################################################################################

# Read in day number from command line:  

n = len(sys.argv)
if n == 0 or n == 1 or n > 3:
    print('Usage:  python process_movie.py ### $$$' )
    print('   where ### is a 3-digit identifier for the day number')
    print('   and $$$ is an (optional) 3-digit identifier for the observer code')
    # end program:
    sys.exit(0)

if n == 2:  
    daynum = sys.argv[1]
    obscode = '*'
if n == 3:  
    daynum = sys.argv[1]
    obscode = sys.argv[2]

#################################################################################
# Use glob to create list:
current_dir = os.getcwd()
if obscode == '*':
    image_dir = '/home/cannon/Grism_Focus_Test_WMurillo_2424/'
if obscode != '*':
    image_dir = '/home/cannon/Grism_Focus_Test_WMurillo_2424/'
out_dir = current_dir+'/'
os.makedirs(out_dir, exist_ok = True)
image_list = sorted(glob.glob(image_dir+'*.fts',recursive=True))
N_images = len(image_list)
print('-----------------------------------------------------------')
if obscode == '*':
    print('    --> '+str(N_images)+' images found for day '+str(daynum)+' from all observers')
if obscode != '*':
    print('    --> '+str(N_images)+' images found for day '+str(daynum)+' from observer '+str(obscode))
print('-----------------------------------------------------------')
#################################################################################

#################################################################################
# Plot temp vs. time:
templist = []
timelist = []
amblist=[]
sinklist=[]
explist = []
ratio=[]
wndlist=[]
for ii in range(0, N_images):
    im = fits.getdata(image_list[ii])
    hd = fits.getheader(image_list[ii])
    exptime = hd['EXPTIME']
    logexptime= np.log(exptime)/np.log(2)
    DATE=str.split(hd['DATE-OBS'],'T')
    day = DATE[0]
    hour    = float(str.split(DATE[1],':')[0])
    minute  = float(str.split(DATE[1],':')[1])
    second  = float(str.split(DATE[1],':')[2])    
    time=hour + (minute/60.) + (second/3600.)
    temp=hd['CCD-TEMP']
    hsink=hd['HSINKT']
    timelist.append(time)
    templist.append(temp)
    explist.append(logexptime)
    sinklist.append(hsink)

#################################################################################
# Sort the list by time: 
astlist = [x for _,x in sorted(zip(timelist,image_list))]
#################################################################################

#%matplotlib inline
#daynum = str.split(str.split(astlist[i],'/')[-1],'.')[0][3:6]
print('Processing '+str(N_images)+' images from day number '+str(daynum))
with PdfPages(out_dir+str(daynum)+'.H-alpha-grism.images.pdf') as pdf:
      # for i in range(0,4):
      for i in range(len(image_list)):
            name = str.split(str.split(astlist[i],'/')[-1],'.')[0]
            with fits.open(str(astlist[i])) as hdul:
                minval = np.mean(im[20:4075,20:4075])-np.std(im[20:4075,20:4075])*0.125
                maxval = np.mean(im[20:4075,20:4075])+2*np.std(im[20:4075,20:4075])*0.125
                meanval = np.mean(im[20:4075,20:4075])
                stdval = np.std(im[20:4075,20:4075])
                im = hdul[0].data
                hd = hdul[0].header
                fig = plt.figure(figsize=(8,8))
                ax1 = fig.add_subplot(projection = None)
                im0 = ax1.imshow(im,vmin=minval,vmax=maxval,cmap='cividis')
                ax1.annotate(str(name)+'  ('+str(hd['OBJECT'])+', '+str(hd['FILTER'])+ ' filter)', xy=(100,200), fontsize=10, color='white')
                ax1.annotate('UTC  '+str.split(hd['DATE-OBS'],'T')[0]+'  '+str.split(hd['DATE-OBS'],'T')[1],xy=(100,300), fontsize=10, color='white')
                ax1.annotate('Mean = '+str("{:.3f}".format(meanval)), xy=(100,400), fontsize=10, color='white')
                ax1.annotate(r'$\sigma$ = '+str("{:.3f}".format(stdval)), xy=(100,500), fontsize=10, color='white')
                ax1.annotate(r'Exposure Time ='+str(round(hd['EXPTIME'],3))+' s', xy=(100,600), fontsize=10, color='white')
                ax1.annotate(hd['READOUTM'],xy=(100,700), fontsize=10, color='white')
                ax1.annotate(r'min = mean - 0.125 $\sigma$ = '+str(round(minval,3)),xy=(100,800), fontsize=10, color='white')
                ax1.annotate(r'max = mean + 0.75 $\sigma$ = '+str(round(maxval,3)),xy=(100,900), fontsize=10, color='white')
                pdf.savefig(dpi=100)
                plt.close()
                print('Image #'+str(i+1)+' processed, '+str(round(100.*((i+1)/len(astlist)),1))+'% complete')

if obscode == '*':
    print('File '+out_dir+str(daynum)+'.images.pdf successfully created')
if obscode != '*':
    shutil.move(out_dir+str(daynum)+'.images.pdf',out_dir+str(daynum)+'.'+str(obscode)+'.images.pdf')
    print('File '+out_dir+str(daynum)+'.'+str(obscode)+'.images.pdf successfully created')
#plt.close()