<commit_before>
"""
specpolflux

in current directory, save fluxcal table for flux standards in ec* files
if an appropriate fluxcal file now exists, apply to listed _stokes.fits files 

"""

import os
import sys
import glob
import shutil
import inspect

import numpy as np
from astropy.io import fits as pyfits
from astropy.table import Table,unique
from scipy.interpolate import interp1d

from pyraf import iraf
from iraf import pysalt
from saltobslog import obslog
from saltsafelog import logging
from scrunch1d import scrunch1d

import reddir
# from . import reddir
DATADIR = os.path.dirname(inspect.getfile(reddir)) + "/data/"
np.set_printoptions(threshold=np.nan)
debug = True

def specpolflux(infilelist, logfile='salt.log', debug=False):
    """Finds/ produces fluxcal tables, and applies to listed _stokes.fits

    Parameters
    ----------
    infilelist: list
        filename or list of finalstokes filenames 

    logfile: str
        Name of file for logging

    """
  # Info on CAL_SPST files:
    calspstname_s,calspstfile_s=np.loadtxt(DATADIR+"spst_filenames.txt",    \
        dtype=str,usecols=(0,1),unpack=True)
    namelen = max(map(len,calspstname_s))

    confitemlist = ['GRATING','GR-ANGLE','CAMANG']

  # Find fluxdb files already in this directory
    fluxdbtab = Table(names=['no','OBJECT']+confitemlist,        \
               dtype=[int,'S'+str(namelen),'S6',float,float])
    fluxdblist = sorted(glob.glob('*fluxdb.txt'))
    olddbentries = len(fluxdblist)
    for e,dbfile in enumerate(fluxdblist):        
        confdat_d = np.genfromtxt(dbfile,usecols=2,comments='?', \
                max_rows=4,dtype=str)
        fluxdbtab.add_row(np.insert(confdat_d,0,e+1))
    if olddbentries: 
        printstdlog('\n    Existing Fluxdb:\n'+str(fluxdbtab), logfile)

  # Create new fluxdb files if new data present
    eclist = sorted(glob.glob('ec*.fits'))
    obs_i,config_i,obstab,configtab = configmap(eclist,confitemlist)
    obss = len(obstab)

    for obs in range(obss):
        object,config = obstab[obs]

      # Find spst standards
        i_j = np.where(obs_i == obs)[0]
        if object not in calspstname_s: continue

        newfluxdbtab = Table( \
            names=['no','OBJECT']+confitemlist,dtype=[int,'S'+str(namelen),'S6',float,float])
        newfluxdbtab.add_row([len(fluxdbtab)+1,object]+list(configtab[obs]))       
        if Table(newfluxdbtab.columns[1:]) in Table(fluxdbtab.columns[1:]): continue

      # It's a new flux standard, process it
        s = np.where(object==calspstname_s)[0][0]

        spstfile=iraf.osfn("pysalt$data/standards/spectroscopic/"+calspstfile_s[s])
        wav_f,ABmag_f = np.loadtxt(spstfile,usecols=(0,1),unpack=True)
        flam_f = 10.**(-0.4*(ABmag_f+2.402))/(wav_f)**2

        wbinedge_f = (wav_f[1:] + wav_f[:-1])/2.
        wbinedge_f = np.insert(wbinedge_f,0,2.*wav_f[0]-wbinedge_f[0])
        wbinedge_f = np.append(wbinedge_f,2.*wav_f[-1]-wbinedge_f[-1])
      
      # average multiple samples,E,O
        hdul = pyfits.open(eclist[i_j[0]])
        grating,grang,artic = configtab[config]
        wav0 = hdul['SCI'].header['CRVAL1']
        dwav = hdul['SCI'].header['CDELT1']
        wavs = hdul['SCI'].data.shape[-1]
        phot_w = np.zeros(wavs)
        count_w = np.zeros(wavs)
        exptime = 0.
        samples = i_j.shape[0]
        for j in range(samples):
            hdul = pyfits.open(eclist[i_j[j]])
            phot_w += hdul['SCI'].data.reshape((2,-1)).sum(axis=0)                
            count_w += (hdul['BPM'].data.reshape((2,-1))==0).sum(axis=0)
            exptime +=  hdul['SCI'].header['EXPTIME']  
        int_w = phot_w/(2*samples*exptime)
        ok_w = (count_w == 2*samples)
      # scrunch onto flux star grid
        wav_w = np.arange(wav0,wav0+wavs*dwav,dwav)
        binedge_f = (wbinedge_f - (wav_w[0] - dwav/2.))/dwav
        int_f = scrunch1d(int_w,binedge_f)
        ok_f = (scrunch1d((~ok_w).astype(int),binedge_f) == 0)   # good flux bins have no bad wavs
        ok_f &= ((wav_f > wav_w[0]) & (wav_f < wav_w[-1]))

      # save flux/intensity mean, extrapolate to edge
        fluxcal_F = flam_f[ok_f]/int_f[ok_f]
        wav_F = wav_f[ok_f]
        fluxcalslope_F = (fluxcal_F[1:]-fluxcal_F[:-1])/(wav_F[1:]-wav_F[:-1])
        wav_F = np.insert(wav_F,0,wav_w[0])
        wav_F = np.append(wav_F,wav_w[-1])
        fluxcal_F = np.insert(fluxcal_F,0,fluxcal_F[0]-fluxcalslope_F[0]*(wav_F[1]-wav_F[0]))
        fluxcal_F = np.append(fluxcal_F,fluxcal_F[-1]+fluxcalslope_F[-1]*(wav_F[-1]-wav_F[-2]))
        fluxdbfile = calspstname_s[s]+'_c'+str(config)+'_fluxdb.txt'
        hdr = ("OBJECT: "+object+"\nGRATING: %s \nARTIC: %s \nGRANG: %s" \
            % (grating,grang,artic))
        np.savetxt(fluxdbfile, np.vstack((wav_F,fluxcal_F)).T,fmt="%8.2f %12.3e ",header=hdr)
        fluxdbtab.add_row(list(newfluxdbtab[0]))
        fluxdblist.append(fluxdbfile)
    dbentries = len(fluxdbtab)<commit_msg>knsep17_4 scripts/specpolflux.py      Changed fluxdb file naming to be consistent with arcdb<commit_after>
"""
specpolflux

in current directory, save fluxcal table for flux standards in ec* files
if an appropriate fluxcal file now exists, apply to listed _stokes.fits files 

"""

import os
import sys
import glob
import shutil
import inspect

import numpy as np
from astropy.io import fits as pyfits
from astropy.table import Table,unique
from scipy.interpolate import interp1d

from pyraf import iraf
from iraf import pysalt
from saltobslog import obslog
from saltsafelog import logging
from scrunch1d import scrunch1d

import reddir
# from . import reddir
DATADIR = os.path.dirname(inspect.getfile(reddir)) + "/data/"
np.set_printoptions(threshold=np.nan)
debug = True

def specpolflux(infilelist, logfile='salt.log', debug=False):
    """Finds/ produces fluxcal tables, and applies to listed _stokes.fits

    Parameters
    ----------
    infilelist: list
        filename or list of finalstokes filenames 

    logfile: str
        Name of file for logging

    """
  # Info on CAL_SPST files:
    calspstname_s,calspstfile_s=np.loadtxt(DATADIR+"spst_filenames.txt",    \
        dtype=str,usecols=(0,1),unpack=True)
    namelen = max(map(len,calspstname_s))

    confitemlist = ['GRATING','GR-ANGLE','CAMANG']

  # Find fluxdb files already in this directory
    fluxdbtab = Table(names=['no','OBJECT']+confitemlist,        \
               dtype=[int,'S'+str(namelen),'S6',float,float])
    fluxdblist = sorted(glob.glob('fluxdb*.txt'))
    olddbentries = len(fluxdblist)
    for e,dbfile in enumerate(fluxdblist):        
        confdat_d = np.genfromtxt(dbfile,usecols=2,comments='?', \
                max_rows=4,dtype=str)
        fluxdbtab.add_row(np.insert(confdat_d,0,e+1))
    if olddbentries: 
        printstdlog('\n    Existing Fluxdb:\n'+str(fluxdbtab), logfile)

  # Create new fluxdb files if new data present
    eclist = sorted(glob.glob('ec*.fits'))
    obs_i,config_i,obstab,configtab = configmap(eclist,confitemlist)
    obss = len(obstab)

    for obs in range(obss):
        object,config = obstab[obs]

      # Find spst standards
        i_j = np.where(obs_i == obs)[0]
        if object not in calspstname_s: continue

        newfluxdbtab = Table( \
            names=['no','OBJECT']+confitemlist,dtype=[int,'S'+str(namelen),'S6',float,float])
        newfluxdbtab.add_row([len(fluxdbtab)+1,object]+list(configtab[obs]))       
        if Table(newfluxdbtab.columns[1:]) in Table(fluxdbtab.columns[1:]): continue

      # It's a new flux standard, process it
        s = np.where(object==calspstname_s)[0][0]

        spstfile=iraf.osfn("pysalt$data/standards/spectroscopic/"+calspstfile_s[s])
        wav_f,ABmag_f = np.loadtxt(spstfile,usecols=(0,1),unpack=True)
        flam_f = 10.**(-0.4*(ABmag_f+2.402))/(wav_f)**2

        wbinedge_f = (wav_f[1:] + wav_f[:-1])/2.
        wbinedge_f = np.insert(wbinedge_f,0,2.*wav_f[0]-wbinedge_f[0])
        wbinedge_f = np.append(wbinedge_f,2.*wav_f[-1]-wbinedge_f[-1])
      
      # average multiple samples,E,O
        hdul = pyfits.open(eclist[i_j[0]])
        grating,grang,artic = configtab[config]
        wav0 = hdul['SCI'].header['CRVAL1']
        dwav = hdul['SCI'].header['CDELT1']
        wavs = hdul['SCI'].data.shape[-1]
        phot_w = np.zeros(wavs)
        count_w = np.zeros(wavs)
        exptime = 0.
        samples = i_j.shape[0]
        for j in range(samples):
            hdul = pyfits.open(eclist[i_j[j]])
            phot_w += hdul['SCI'].data.reshape((2,-1)).sum(axis=0)                
            count_w += (hdul['BPM'].data.reshape((2,-1))==0).sum(axis=0)
            exptime +=  hdul['SCI'].header['EXPTIME']  
        int_w = phot_w/(2*samples*exptime)
        ok_w = (count_w == 2*samples)
      # scrunch onto flux star grid
        wav_w = np.arange(wav0,wav0+wavs*dwav,dwav)
        binedge_f = (wbinedge_f - (wav_w[0] - dwav/2.))/dwav
        int_f = scrunch1d(int_w,binedge_f)
        ok_f = (scrunch1d((~ok_w).astype(int),binedge_f) == 0)   # good flux bins have no bad wavs
        ok_f &= ((wav_f > wav_w[0]) & (wav_f < wav_w[-1]))

      # save flux/intensity mean, extrapolate to edge
        fluxcal_F = flam_f[ok_f]/int_f[ok_f]
        wav_F = wav_f[ok_f]
        fluxcalslope_F = (fluxcal_F[1:]-fluxcal_F[:-1])/(wav_F[1:]-wav_F[:-1])
        wav_F = np.insert(wav_F,0,wav_w[0])
        wav_F = np.append(wav_F,wav_w[-1])
        fluxcal_F = np.insert(fluxcal_F,0,fluxcal_F[0]-fluxcalslope_F[0]*(wav_F[1]-wav_F[0]))
        fluxcal_F = np.append(fluxcal_F,fluxcal_F[-1]+fluxcalslope_F[-1]*(wav_F[-1]-wav_F[-2]))
        fluxdbfile = 'fluxdb_'+calspstname_s[s]+'_c'+str(config)+'.txt'
        hdr = ("OBJECT: "+object+"\nGRATING: %s \nARTIC: %s \nGRANG: %s" \
            % (grating,grang,artic))
        np.savetxt(fluxdbfile, np.vstack((wav_F,fluxcal_F)).T,fmt="%8.2f %12.3e ",header=hdr)
        fluxdbtab.add_row(list(newfluxdbtab[0]))
        fluxdblist.append(fluxdbfile)
    dbentries = len(fluxdbtab)
{'fraction_non_alphanumeric': 0.14736525746343418, 'fraction_numerical': 0.017331196153075535, 'mean_word_length': 3.1150041220115416, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2189152', 'n_tokens_mistral': 4029, 'n_tokens_neox': 3894, 'n_words': 647}