Source code for tshirt.pipeline.analysis

#import phot_pipeline
import numpy as np
from scipy.stats import binned_statistic
import matplotlib.pyplot as plt
from astropy.table import Table
from import fits, ascii
from scipy import signal
from scipy.interpolate import interp1d
from scipy import ndimage
import matplotlib.pyplot as plt
import pdb
from .utils import robust_poly
from .phot_pipeline import phot
from .spec_pipeline import spec
from copy import deepcopy
import os
import tqdm

def bin_examine():
    t ='tser_data/refcor_phot/refcor_phot_S02illum1miAll_GLrun104.fits')
    allStd, allBinSize = [], []

    binsTry = [500,100,50,30,10]
    for oneBin in binsTry:
        yBin, xEdges, binNum = binned_statistic(t['Time'],t['Y Corrected'],
        stdBinned = np.std(yBin)
        binSize = (xEdges[1] - xEdges[0]) * 24. * 60.

    print("{} bin ({} min) std = {} ppm".format(binsTry[-1],allBinSize[-1],allStd[-1]))

[docs] def adjust_aperture_set(phot_obj,param,srcSize,backStart,backEnd, showPlot=False,plotHeight=None): """ Takes a photometry or spectroscopy object and sets the aperture parameters Parameters ------------ phot_obj: a tshirt phot or spec object Input object to be adjusted param: dictionary A parameter dictionary for a phot or spec object showPlot: bool Show a plot of the aperture set? plotHeight: float A size of the photometry plot (only works for photometry) Returns -------- new_phot: a tshirt phot or spec object Spec or phot object w/ modified parameters """ if phot_obj.pipeType == 'photometry': if phot_obj.param['scaleAperture'] == True: param['apRadius'] = 1.0 param['apScale'] = srcSize else: param['apRadius'] = srcSize param['backStart'] = backStart param['backEnd'] = backEnd new_phot = phot(directParam=param) if showPlot == True: new_phot.showStamps(boxsize=plotHeight) elif phot_obj.pipeType == 'spectroscopy': param['apWidth'] = srcSize backLoc_1_start = param['starPositions'][0] - backEnd backLoc_1_end = param['starPositions'][0] - backStart backLoc_2_start = param['starPositions'][0] + backStart backLoc_2_end = param['starPositions'][0] + backEnd backLocList = [[backLoc_1_start,backLoc_1_end], [backLoc_2_start,backLoc_2_end]] if param['dispDirection'] == 'x': param['bkgRegionsY'] = backLocList else: param['bkgRegionsX'] = backLocList new_phot = spec(directParam=param) if showPlot == True: new_phot.showStarChoices() else: raise Exception("Unrecognized pipeType") return new_phot
[docs] def aperture_size_sweep(phot_obj,stepSize=5,srcRange=[5,20],backRange=[5,28], minBackground=2,stepSizeSrc=None,stepSizeBack=None, shorten=False): """ Calculate the Noise Statistics for a "Sweep" of Aperture Sizes Loops through a series of source sizes and background sizes in a grid search Parameters ----------- stepSize: float The step size. It will be superseded by stepSize_bck or stepSizeSrc if used. srcRange: two element list The minimum and maximum src radii to explore backRange: two element list The minimum and maximum aperture radii to explore (both for the inner & outer) minBackground: float The minimum thickness for the background annulus or rectangle width stepSizeSrc: float (optional) Specify the step size for the source that will supersed the general stepSize stepSizeBack: float (optional) Specify the step size for the background that will supersed the general stepSize shorten: bool Shorten the time series? (This is passed to print_phot_statistics) """ if stepSizeSrc == None: stepSizeSrc = stepSize if stepSizeBack == None: stepSizeBack = stepSize ## change the short name to avoid over-writing previous files origParam = deepcopy(phot_obj.param) origName = origParam['srcNameShort'] param = deepcopy(origParam) ## show the most compact and most expanded configurations srcPlots = srcRange backStartPlots = [np.max([srcRange[0],backRange[0]]), backRange[1] - minBackground] backEndPlots = [backStartPlots[0] + minBackground, backRange[1]] plotHeights = [backRange[1]+5,backRange[1] + 5] for i, oneConfig in enumerate(['compact','expanded']): param['srcNameShort'] = origName + '_' + oneConfig new_phot = adjust_aperture_set(phot_obj,param,srcRange[i], backStartPlots[i],backEndPlots[i], showPlot=True,plotHeight=plotHeights[i]) apertureSets = [] t = Table(names=['src','back_st','back_end']) for srcSize in np.arange(srcRange[0],srcRange[1],stepSizeSrc): ## start from the backRange min or the src, whichever is bigger back_st_minimum = np.max([srcSize,backRange[0]]) ## finish at the backRange max, but allow thickness back_st_maximum = backRange[1] - minBackground for back_st in np.arange(back_st_minimum,back_st_maximum,stepSizeBack): ## start the outer background annulus, at least minBackground away back_end_minimum = back_st + minBackground back_end_maximum = backRange[1] for back_end in np.arange(back_end_minimum,back_end_maximum,stepSize): apertureSets.append([srcSize,back_st,back_end]) t.add_row([srcSize,back_st,back_end]) #print("src: {}, back st: {}, back end: {}".format(srcSize,back_st,back_end)) ## for the rest, it will save a common file name param['srcNameShort'] = origName + '_aperture_sizing' stdevArr, theo_err, mad_arr = [], [], [] for i,apSet in enumerate(apertureSets): print("src: {}, back st: {}, back end: {}".format(apSet[0],apSet[1],apSet[2])) new_phot = adjust_aperture_set(phot_obj,param,apSet[0],apSet[1],apSet[2], showPlot=False) if phot_obj.pipeType == 'photometry': new_phot.do_phot(useMultiprocessing=True) noiseTable = new_phot.print_phot_statistics(refCorrect=True,returnOnly=True,shorten=shorten) elif phot_obj.pipeType == 'spectroscopy': new_phot.do_extraction(useMultiprocessing=True) noiseTable = new_phot.print_noise_wavebin(shorten=shorten,nbins=1,recalculate=True) else: raise Exception("Unrecognized pipeType") stdevArr.append(noiseTable['Stdev (%)'][0]) theo_err.append(noiseTable['Theo Err (%)'][0]) mad_arr.append(noiseTable['MAD (%)'][0]) t['stdev'] = stdevArr t['theo_err'] = theo_err t['mad_arr'] = mad_arr outTable_name = 'aperture_opt_{}_src_{}_{}_step_{}_back_{}_{}_step_{}.csv'.format(new_phot.dataFileDescrip, srcRange[0],srcRange[1],stepSizeSrc, backRange[0],backRange[1], stepSizeBack) if phot_obj.pipeType == 'photometry': tableDir = 'phot_aperture_optimization' elif phot_obj.pipeType == 'spectroscopy': tableDir = 'spec_aperture_optimization' else: raise Exception("Unrecognized pipeType") outTable_path = os.path.join(new_phot.baseDir,'tser_data',tableDir,outTable_name) t.write(outTable_path,overwrite=True) print('Writing table to {}'.format(outTable_path)) ind = np.argmin(t['stdev']) print("Min Stdev results:") print(t[ind]) return t
[docs] def plot_apsizes(apertureSweepFile,showPlot=True): """ Plot the aperture sizes calculated from :any:`aperture_size_sweep` Parameters ---------- apertureSweepFile: str A .csv file created by aperture_size_sweep showPlot: bool Show the plot w/ matplotlib? otherwise, it saves to file """ dat = fig, axArr2D = plt.subplots(3,3,sharey=True) ## share axes along columns for oneColumn in [0,1,2]: axTop = axArr2D[0,oneColumn] axMid = axArr2D[1,oneColumn] axBot = axArr2D[2,oneColumn] axTop.get_shared_x_axes().join(axMid, axBot) labels = ['Source Radius','Back Start','Back End'] keys = ['src','back_st','back_end'] statistics = ['stdev','theo_err','mad_arr'] for statInd,statistic in enumerate(statistics): axArr1D = axArr2D[statInd] for ind, ax in enumerate(axArr1D): ax.semilogy(dat[keys[ind]],dat[statistic],'.') ax.set_xlabel(labels[ind]) if ind==0: ax.set_ylabel(statistic) if showPlot == True: else: raise NotImplementedError
[docs] def backsub_list(specObj,outDirectory='tmp'): """ Save all the background-subtracted images Parameters ----------- specObj: a `any::tshirt.pipeline.spec_pipeline.spec` object """ for i in tqdm.tqdm(range(len(specObj.fileL))): oneFile = specObj.fileL[i] # read in img, head = specObj.getImg(oneFile) # do backsub subImg, bkgModel, subHead = specObj.do_backsub(img,head,i,saveFits=False, directions=specObj.param['bkgSubDirections']) ## save baseName = os.path.splitext(os.path.basename(oneFile))[0] outName = baseName + '_backsub.fits' outPath = os.path.join(outDirectory,outName) HDUList = fits.PrimaryHDU(subImg,subHead) HDUList.writeto(outPath)