#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Fri Apr 13 18:04:48 2018

Versión para la web de EMIR.

@author: fgl
"""
import numpy as np
from astropy.io import fits
import os
import sys


tread = 1.055289     # lectura productiva [s]
tread_nd = 1.048046 # lectura no productiva [s]
trst = 1.056230       # reset [s]

if (len(sys.argv) - 1)<7:
    print('This script needs 7 arguments: ')
    print('\t* the point from which to fit the slope.')
    print('\t* the folder to get the raw data.')
    print('\t* the folder to write the output data.')
    print('\t* the suffix of the raw data.')
    print('\t* the suffix of the processed.')
    print('\t* the first id number of the raw data.')
    print('\t* the last id number of the raw data.')
    sys.exit()


#frseci=5       # primer punto a considerar de la rampa. -1: 2ª mitad
#
#dir_in='/emir/mos_sciver/28may18/raw/'
#dir_out='/emir/mos_sciver/28may18/goya_mask1/'
#suffix='-20180528-EMIR-STARE_SPECTRA_raw.fits'
#suffix_proc='-20180528-EMIR-STARE_SPECTRA_mod.fits'
#regini=1567287
#regfin=1567637

frseci=int(sys.argv[1])
dir_in=sys.argv[2]
dir_out=sys.argv[3]
suffix=sys.argv[4]
suffix_proc=sys.argv[5]
regini=int(sys.argv[6])
regfin=int(sys.argv[7])

ini=True
for i in range(regini,regfin+1):
    fname_fits=dir_in+str("%010d" % i)+suffix
#    print (fname_fits)
    if not os.path.exists(fname_fits):
#      fname_fits=dir_in+str("%010d" % i)+suffix_pp
#      if not os.path.exists(fname_fits):
        continue
    hdur=fits.open(fname_fits)
    frsec=hdur[0].header['frsec']
    ts=hdur[0].header['tsutc2']
    if ini:  # comienzo
        ini=False
        nfrsec=hdur[0].header['nfrsec']
        toper=hdur[0].header['toper']
        nrdil=hdur[0].header['nrdil']
        nrdil_nd=hdur[0].header['nrdil_nd']
        nfrsec=hdur[0].header['nfrsec']
        if frseci < 0 : frseci=int(np.ceil(nfrsec/2))
        nfr=nfrsec-frseci+1
        exptime=hdur[0].header['exptime']
        bzero=hdur[0].header['bzero']
        bscale=hdur[0].header['bscale']
        dt = (toper/1000. + nrdil_nd*tread_nd + nrdil*tread)/exptime
        fac = 12 / (nfr*nfr - 1)/ nfr/dt
        print ('\t\t',bzero,bscale,exptime,dt*exptime,nfr,fac)
#    if frsec==0:
#        data=hdur[0].data
#        d1/=data[900,1000]
#        d2/=data[1000,900]
#        print (i,data[900,1000],data[1000,900],1/d1,1/d2)
    if frsec<frseci:
        hdur.close()
        continue
    elif frsec==frseci:
        reglist=np.zeros(nfr,dtype=int)
        sumy=np.zeros((2048,2048))
        sumxy=np.zeros((2048,2048))
#        print ('\t\t',tread_nd,tread,exptime,nfrsec,dt)
        imgobbl=hdur[0].header['imgobbl']
    data=hdur[0].data
    hdur.close()
    j=frsec-frseci
    sumxy+=np.float_(data)*(j+1)
    sumy+=np.float_(data)
#    sumxy+=data*(j+1)
#    sumy+=data
#    print (j+1,data[900,1000],data[900,1000]*(j+1),sumy[900,1000],sumxy[900,1000])
    reglist[j]=i
    if frsec==nfrsec:
        print ('Procesados ficheros: ',*reglist)
        head=hdur[0].header
        imgobbl=head['imgobbl']
#        data=np.round((sumxy - (nfr+1)*sumy/2.)*fac)
        data=(sumxy - (nfr+1)*sumy/2)*fac
        data[data<0]=0
#        d1=data[900,1000]
#        d2=data[1000,900]
#        d1=(sumxy[900,1000] - (nfr+1)*sumy[900,1000]/2.)*fac
#        print (i+1,sumy[900,1000])
#        print (i+1,sumxy[900,1000])
#        print (i+1,data[900,1000],d1)
#        hduw = fits.PrimaryHDU(np.round(data))
        hduw = fits.PrimaryHDU(data)
        hduw.scale('int16', bzero=bzero,bscale=bscale)
        hduw.header=head
        hduw.header['bscale']=bscale
        hduw.header['frsec']=0
        fname=dir_out+str("%010d" % (i+1))+suffix_proc
        hduw.writeto(fname, overwrite=True)
        print('*'*5,'Escrito ',fname.split('/')[-1])
        
#        hdur=fits.open(fname)
#        data=hdur[0].data
#        print (i+1,data[900,1000],data[1000,900])
#        hdur.close()