def makedark(files, output):
    """
    Make dark image for NIRC2 data. Makes a calib/ directory
    and stores all output there. All output and temporary files
    will be created in a darks/ subdirectory. 

    files: integer list of the files. Does not require padded zeros.
    output: output file name. Include the .fits extension.
    """
    redDir = os.getcwd() + "/"  # Reduce directory.
    curDir = redDir + "calib/"
    darkDir = util.trimdir(curDir + "darks/")
    rawDir = util.trimdir(os.path.abspath(redDir + "../raw") + "/")

    util.mkdir(curDir)
    util.mkdir(darkDir)

    _out = darkDir + output
    _outlis = darkDir + "dark.lis"
    util.rmall([_out, _outlis])

    darks = [rawDir + "n" + str(i).zfill(4) + ".fits" for i in files]

    f_on = open(_outlis, "w")
    f_on.write("\n".join(darks) + "\n")
    f_on.close()

    ir.unlearn("imcombine")
    ir.imcombine.combine = "median"
    ir.imcombine.reject = "sigclip"
    ir.imcombine.nlow = 1
    ir.imcombine.nhigh = 1
    ir.imcombine("@" + _outlis, _out)

def ImgCombineWithZeroFloating(imglistfname,outputfile,cmethod="median",czero="median",
                               creject="avsigclip",cstatsection='[150:900,150:900]'):
    """ Returns the combined image with actuall average median flux, 
    It does zero scaling only for sigma rejection of stars. This is needed to remove faint 
    stars in rejection algorithm when the background sky itself is varying from frame to frame. """
    iraf.imcombine.unlearn()
    Xmin=float(cstatsection[1:-1].split(',')[0].split(':')[0])  #Everything now in fits coordinates
    Xmax=float(cstatsection[1:-1].split(',')[0].split(':')[1])
    Ymin=float(cstatsection[1:-1].split(',')[1].split(':')[0])
    Ymax=float(cstatsection[1:-1].split(',')[1].split(':')[1])

    if czero == "median" : statfunction = np.median
    elif czero == "average" : statfunction = np.mean
    else : 
        print('Error: czero should be median or average. Unknown option {0}'.format(czero))
        raise

    with open(imglistfname,'r') as imgfile:
        statlist=[]
        for img in imgfile:
            img = img.rstrip()
            statlist.append(statfunction(fits.getdata(img)[Ymin-1:Ymax,Xmin-1:Xmax]))
    print('{0} of images: {1}'.format(czero,str(statlist)))
    statAvg=np.mean(statlist)
    Zeroshifts= statAvg - np.array(statlist)
    print('Zeroshifts of images: {0} :: ImgAvg ={1}'.format(str(Zeroshifts),statAvg))
    with open(outputfile+'_zeroshifts.txt','w') as zeroshiftFILE:
        for shift in Zeroshifts: 
            zeroshiftFILE.write('{0} \n'.format(shift))
    # Now call iraf imcombine with zero scaling
    iraf.imcombine(input='@'+imglistfname, output=outputfile, combine=cmethod, 
                   reject=creject, statsec=cstatsection, zero='@'+outputfile+'_zeroshifts.txt')

def cos_clear2(filenames):
    outname = 'fake_' + filenames[0]
    if os.path.isfile(outname):
        print outname, 'is already exist'
    else:
        inname = filenames[0] + ',' + filenames[1]
        print 'runing cos_clear2, imcombine...'
        print 'make file', outname
        iraf.imcombine(input = inname
                , output = outname, headers = '', bpmasks = ''
                , rejmasks = '', nrejmasks = '', expmasks = ''
                , sigmas = '', imcmb = '$I', logfile = 'STDOUT'
                , combine = 'average', reject = 'minmax', project = False
                , outtype = 'real', outlimits = 'none', offsets = 'none'
                , masktype = 'none', maskvalue = 0, blank = 0.0
                , scale = 'exposure', zero = 'none', weight = 'none'
                , statsec = '', expname = 'EXPTIME', lthreshold = 'INDEF'
                , hthreshold = 'INDEF', nlow = 0, nhigh = 1
                , nkeep = 1, mclip = True, lsigma = 3.0
                , hsigma = 10.0, rdnoise = 'RDNOISE', gain = 'GAIN'
                , snoise = 0.0, sigscale = 0.1, pclip = -0.5, grow = 0.0)
        iraf.hedit(images = outname
                , fields = 'ncombine', value = '0', add = True
                , addonly = False, delete = False, verify = False
                , show = True, update = True)
        iraf.hedit(images = outname
                , fields = 'EXPTIME', value = '0', add = True
                , addonly = False, delete = False, verify = False
                , show = True, update = True)
    print 'display %s 1' % outname
    iraf.display(image = outname, frame = 1)
    filenames.append(outname)
    return cos_clear3(filenames)

def combine_subset(filter='G141', idx=np.array([0]), root='set1', use_scaled=True):
    """
    Subset, get index array of objects to use from the "show_profile" function above
    """
    
    from pyraf import iraf
    
    bg_flt, bg_field, bg_val = np.loadtxt('background.%s.dat' %(filter), dtype=np.str, unpack=True)
    weights = np.cast[float](bg_val)**2
           
    fp = open('%s.%s.list' %(filter, root),'w')
    fpw = open('%s.%s.weight' %(filter, root),'w')
    for msk, wht in zip(bg_flt[idx], weights[idx]):
        if os.path.exists(msk):
            if use_scaled:
                img = msk.replace('msk','msk.s')
            else:
                img = msk
            fp.write('%s\n' %(img))
            fpw.write('%.4f\n' %(wht))
    #
    fp.close()
    fpw.close()
        
    iraf.imcombine ( input = '@%s.%s.list' %(filter, root), output = 'combine.%s.%s' %(filter, root), 
       headers = '', bpmasks = '', rejmasks = '', nrejmasks = '', 
       expmasks = '', sigmas = '', logfile = 'STDOUT', combine = 'average', 
       reject = 'minmax', project = iraf.no, outtype = 'real', 
       outlimits = '', offsets = 'none', masktype = 'none', 
       maskvalue = '0', blank = 0.0, scale = 'none', zero = 'none', 
       weight = '@%s.%s.weight' %(filter, root), statsec = '', expname = '', lthreshold = 1e-04, 
       hthreshold = 100.0, nlow = 2, nhigh = 2, nkeep = 1, 
       mclip = iraf.yes, lsigma = 3.0, hsigma = 3.0, rdnoise = '0.', 
       gain = '1.', snoise = '0.', sigscale = 0.1, pclip = -0.5)

    def gen_dark(self,exptime):
            iraf.imcombine.reject = 'none'
            iraf.imcombine.combine = self.darkcombine

            name = self.dark_dir+'/Combined_'+str(exptime)+'_DARK.fits'
            iraf.imcombine(','.join(self.darks[exptime]['raw']),name)
            self.darks[exptime]['combined'] = name

def createFlatFiles():
	import re
	print "CreateFlatFiles start"
	for f in FILTERS:
		if(os.listdir(os.path.join(OUTPUTDIR, "flat", f))):
			iraf.imcombine.setParam("input", os.path.join(OUTPUTDIR, "flat", f) + "/*.fits")
			flatFile = os.path.join(OUTPUTDIR, "flat", f , "Flat.fits")
			if os.path.exists(flatFile):
				print("flatFile %s alreday exists deleting"  % flatFile)
				os.remove(flatFile)
			iraf.imcombine.setParam("output", flatFile)
			#from doc:	
			#http://www.iac.es/sieinvens/siepedia/pmwiki.php?n=HOWTOs.PythonianIRAF
			#--> iraf.listpix(mode='ql')     # confirms parameter
			#--> iraf.listpix(mode='h')     # doesn't ask for [email protected]@
			iraf.imcombine(mode="h")
			#NORMALIZE
			#imstat
			res = iraf.imstat(flatFile, Stdout=1)
			print(res[0].strip()) 
			print(res[1].strip()) 
			resArray = re.split("\s+", res[1].strip())
			#max value
			#toDivValue = float(resArray[5])
			#meanValue
			toDivValue = float(resArray[2])
			flatNormFile = os.path.join(OUTPUTDIR, "flat", f , "FlatNorm.fits")
			if os.path.exists(flatNormFile):
				print("flatNormFile %s alreday exists deleting"  % flatNormFile)
				os.remove(flatNormFile)
			#divide by max value
			iraf.imarith(flatFile, '/', toDivValue, flatNormFile)
		else:
			print("NO FLAT FILES for filter %s PRESENT" %f)
	print "CreateFlatFiles end"

def combine_cube_frames(cubefile, range_pairs, target_dir):
    dirname, filename = os.path.split(cubefile)
    filebase = filename.rsplit('.', 1)[0]
    if len(filebase) > 8:
        warn("IRAF doesn't like long filenames. "
             "Consider shortening the cube filename ({0})".format(filebase))
    
    outfiles = []
    
    for fromidx, toidx in range_pairs:
        inlst = os.path.join(target_dir, 'in.lst')
        infiles = []
        for i in range(fromidx, toidx+1):
            infiles.append(cubefile + "[*,*,{0}]".format(i))
        
        f = open(inlst, 'w')
        f.writelines(infiles)
        f.write('\n')
        f.close()
        
        outfile = '{0}/{1}_{2}-{3}.fit'.format(target_dir, filebase, fromidx, toidx)
        debug("imcombine input={input} output={output} combine=sum reject=none".format(
            input="@{0}".format(inlst), #','.join(infiles),
            output=outfile,
        ))
        outfiles.append(outfile)
        iraf.imcombine(
            input=','.join(infiles),
            output=outfile,
            combine="sum",
            reject="none",
            # project='no', # IRAF wants bools specified / default is nonboolean?
            # mclip='no',
        )
    return outfiles

def makedark(files, output):
    """
    Make dark image for NIRC2 data. Makes a calib/ directory
    and stores all output there. All output and temporary files
    will be created in a darks/ subdirectory. 

    files: integer list of the files. Does not require padded zeros.
    output: output file name. Include the .fits extension.
    """
    redDir = os.getcwd() + '/'  # Reduce directory.
    curDir = redDir + 'calib/'
    darkDir = util.trimdir(curDir + 'darks/')
    rawDir = util.trimdir(os.path.abspath(redDir + '../raw') + '/')

    util.mkdir(curDir)
    util.mkdir(darkDir)
    
    _out = darkDir + output
    _outlis = darkDir + 'dark.lis'
    util.rmall([_out, _outlis])

    darks = [rawDir + 'n' + str(i).zfill(4) + '.fits' for i in files]

    f_on = open(_outlis, 'w')
    f_on.write('\n'.join(darks) + '\n')
    f_on.close()
    
    ir.unlearn('imcombine')
    ir.imcombine.combine = 'median'
    ir.imcombine.reject = 'sigclip'
    ir.imcombine.nlow = 1
    ir.imcombine.nhigh = 1
    ir.imcombine('@' + _outlis, _out)

def align_combine(fitsdir, myfilter, examine=True):
    from pyraf import iraf 
    
    iraf.noao(_doprint=0)
    iraf.digiphot(_doprint=0)
    iraf.apphot(_doprint=0)
    
    os.chdir(fitsdir)
    listfiles = glob.glob(myfilter)
    listfiles.sort()
    
    if (examine):
        print "Opening ",listfiles[0]," to examine."
        iraf.imexamine(input=listfiles[0], \
                    logfile="coords.dat", \
                    keeplog="yes")
        
        with open("align.list",'w') as f:
            for i in listfiles:
                f.write(i+"\n")
    
    print "Aligning with reference:",listfiles[0]
    iraf.imalign( input   =  "@align.list", referenc= listfiles[0], coords  =  "coords.dat", output  = "[email protected]")  
    
    listfiles = glob.glob("a_"+myfilter)
    listfiles.sort()
    with open("comb.list",'w') as f:
        for i in listfiles:
            f.write(i+"\n")
            
    print "Combining"        
    iraf.imcombine(input = "@comb.list",\
                   output = "out.fits",\
                   combine= "median")

def cos_clear3(filenames):
    outname = 'c3' + filenames[0]
    if os.path.isfile(outname):
        print outname, 'is already exist'
    else:
        inname = ''
        for i in filenames:
            inname = inname + ',' + i
        inname = inname[1:]
        print 'runing cos_clear3, imcombine...'
        print 'make file', outname
        iraf.imcombine(input = inname
                , output = outname, headers = '', bpmasks = ''
                , rejmasks = '', nrejmasks = '', expmasks = ''
                , sigmas = '', imcmb = '$I', logfile = 'STDOUT'
                , combine = 'average', reject = 'avsigclip', project = False
                , outtype = 'real', outlimits = 'none', offsets = 'none'
                , masktype = '', maskvalue = 0.0, blank = 0.0
                , scale = 'none', zero = 'none', weight = 'exposure'
                , statsec = '', expname = 'EXPTIME', lthreshold = 'INDEF'
                , hthreshold = 'INDEF', nlow = 0, nhigh = 1
                , nkeep = 1, mclip = True, lsigma = 3.0
                , hsigma = 8.0, rdnoise = 'RDNOISE', gain = 'GAIN'
                , snoise = 0.0, sigscale = 0.1, pclip = -0.5, grow = 0.0)
    print 'display %s 1' % outname
    iraf.display(image = outname, frame = 1)
    valget = raw_input('Are you need to run crmedian?(y or n): ')
    if valget == 'Y' or valget == 'y':
        return cos_clear1([outname])
    else:
        return outname

def combine(fdict, data_string):

    iraf.imcombine.combine = 'average'
    iraf.imcombine.reject = 'avsigclip'
    iraf.imcombine.lsigma = 3

    name_dict = {}
    if data_string[data_string.find('*')-1] == '_':
        fill = ''
    else:
        fill = '_'

    for f_ratio in fdict.keys():
        name_dict[f_ratio] = []
        
        for filt in fdict[f_ratio].keys():
            name = data_string[:data_string.find('*')]+fill+filt+'F'+\
                str(int(f_ratio*10))
            for ftype in fdict[f_ratio][filt].keys():
                iraf.imcombine(','.join(fdict[f_ratio][filt][ftype]),\
                                   name+ftype[0]+'.fits')
                
                name_dict[f_ratio].append(name+ftype[0]+'.fits')

    return name_dict

def make_dark_avg(images, out):
    """
    Create an average dark frame from a set of images.
    """
    filenames = [im.filename for im in images]
    iraf.imcombine(input=",".join(filenames), output=out, Stdout=0,
                   combine="median")
    return out

def generate_all_masterbias(target_dir, bias_dir_name):
    """ Calculation of all the masterbias files.
    
    This function search for bias files from current directory.
    The bias images are located in specific directories that only
    contains bias images and have a specific denomination, so searching
    for bias files is searching these directories.
    Once a directory for bias had been found a masterbias is calculated
    with an average operation using all the bias files.
    
    Args:
        target_dir: Directory of the files.
        bias_dir_name: Name of the directories that contain bias images.     
    
    """

    logging.info("Generating all masterbias files from %s ..." % target_dir)
    
    # Walk from current directory.
    for path, dirs, files in os.walk(target_dir):
    	
        # Check if current directory is for bias fits.
        for dr in dirs:
            if dr == bias_dir_name:
    				
                # Get the full path of the directory.                
                full_dir = os.path.join(path, dr)
                logging.debug("Found a directory for 'bias': %s" % (full_dir))
                
                # Get the list of files.
                files = glob.glob(os.path.join(full_dir, WILDCARD_FIT_FILE))
                logging.debug("Found %d bias files" % (len(files)))
                
                # Build the masterbias file name.
                masterbias_name = os.path.join(full_dir, MASTERBIAS_FILENAME) 
                
                # Check if masterbias already exists.
                if os.path.exists(masterbias_name) == True:
                    logging.debug("Masterbias file exists '%s', so resume to next directory." % 
                                  (masterbias_name))
                else:                                        
                    # Put the files list in a string.
                    list_of_files = ",".join(files)
                    
                    #show_bias_files_statistics(list_of_files)
                        	
                    # Combine all the bias files.
                    try:
                        logging.debug("Creating bias file: %s" % 
                                    masterbias_name)
                                    
                        iraf.imcombine(list_of_files, masterbias_name, Stdout=1)
                                                
                    except iraf.IrafError as exc:
                        logging.error("Error executing imcombine combining " + \
                                      "bias with: %s" %
                                      (list_of_files))  
                        logging.error("Iraf error is: %s" % (exc))                        

def makeSkyFlat(dfView):
    """
	dfView is a view from either the {J,H}rawFileList data frame
	the view should contain dither positions of the same observation
	this module will calculate the skyflat from the dither positions
	and save it
	"""
    # gather up the images in a python list
    images = dfView.file.values.tolist()
    # they should all have the same date, just grab the first one
    date = str(dfView.Date.values[0])
    # organize the images and output file name in iraf-friendly ways
    inputFiles = joinStrList(images)
    outputSkyFlat = "scratch/" + date + "sky"
    # use iraf imcombine to create the median skyflat
    iraf.imcombine(
        inputFiles,
        outputSkyFlat,
        headers="",
        bpmasks="",
        rejmasks="",
        nrejmasks="",
        expmasks="",
        sigmas="",
        imcmb="$I",
        logfile="STDOUT",
        combine="median",
        reject="none",
        project="no",
        outtype="real",
        outlimits="",
        offsets="none",
        masktype="none",
        maskvalue="0",
        blank=0.0,
        scale="none",
        zero="none",
        weight="none",
        statsec="",
        expname="",
        lthreshold="INDEF",
        hthreshold="INDEF",
        nlow=1,
        nhigh=2,
        nkeep=1,
        mclip="yes",
        lsigma=3.0,
        hsigma=3.0,
        rdnoise="0.",
        gain="1.",
        snoise="0.",
        sigscale=0.1,
        pclip=-0.5,
        grow=0.0,
    )
    return

def reduceBias(bListFn, bMasterFn):
    if not os.path.isfile(bMasterFn + ".fits"):
        iraf.imcombine(
            "@" + bListFn,
            output=bMasterFn,
            logfile=bMasterFn + ".log",
            weight="none",
            zero="none",
            reject="crreject",
            combine="median",
            scale="none",
        )

def createZeroFile():
	print "CreateZeroFile start"
	if(os.listdir(os.path.join(OUTPUTDIR, "bias"))):
		iraf.imcombine.setParam("input", os.path.join(OUTPUTDIR, "bias") + "/*.fits")
		zeroFile = os.path.join(OUTPUTDIR, "bias", "Zero.fits")
		if os.path.exists(zeroFile):
			os.remove(zeroFile)
		iraf.imcombine.setParam("output", zeroFile)
		iraf.imcombine()
	else:
		print("NO BIAS FILES PRESENT")
	print "CreateZeroFile end"

def imcombine(filelist, out, options, bpmask=None, reject="none", nlow=None,
        nhigh=None):

    '''Convenience wrapper around IRAF task imcombine

    Args:
        filelist: The list of files to imcombine
        out: The full path to the output file
        options: Options dictionary
        bpmask: The full path to the bad pixel mask
        reject: none, minmax, sigclip, avsigclip, pclip
        nlow,nhigh: Parameters for minmax rejection, see iraf docs
    
    Returns:
        None

    Side effects:
        Creates the imcombined file at location `out'
    '''

    #TODO: REMOVE Iraf and use python instead. STSCI Python has
    # A builtin routine.
    from pyraf import iraf
    iraf.images()


    filelist = [("%s[0]" % f) for f in filelist]
    pars = iraf.imcombine.getParList()
    iraf.imcombine.unlearn()

    path = "flatcombine.lst"
    f = open(path, "w")
    for file in filelist:
        f.write(file + "\n")
    f.close()

    s = ("%s," * len(filelist))[0:-1]
    s = s % tuple(filelist)

    f = open("flatcombinelog.txt", "w")
    if reject == 'minmax':
        t = iraf.imcombine("@%s" % path, out, Stdout=f,
            reject=reject, nlow=nlow, nhigh=nhigh)
    else:
        t = iraf.imcombine(s, out, Stdin=filelist, Stdout=f,
            reject=reject)
    f.close()
    f=open("flatcombinelog.txt")
    for line in f:
        info(line.rstrip("\n"))
    f.close()

    iraf.imcombine.setParList(pars)

def combineDithers(color):
    """grab the aligned images and sum them up into one"""
    # grab the paths to the sky subtracted flattened aligned images
    aligned = sorted(glob.glob("scratch/s-f-*_gregister.fits"))
    # pick out the YYMMDD date from the file names
    date = aligned[0][17:23]

    # organize input and output in iraf-friendly ways
    inputFiles = joinStrList(aligned)
    outputFile = "reduced/" + date + "." + color + ".4U1543.s-f-a-c"
    # use iraf imcombine to sum them all up
    iraf.imcombine(
        inputFiles,
        outputFile,
        headers="",
        bpmasks="",
        rejmasks="",
        nrejmasks="",
        expmasks="",
        sigmas="",
        imcmb="$I",
        logfile="STDOUT",
        combine="sum",
        reject="none",
        project="no",
        outtype="real",
        outlimits="",
        offsets="none",
        masktype="none",
        maskvalue="0",
        blank=0.0,
        scale="none",
        zero="none",
        weight="none",
        statsec="",
        expname="",
        lthreshold="INDEF",
        hthreshold="INDEF",
        nlow=1,
        nhigh=2,
        nkeep=1,
        mclip="yes",
        lsigma=3.0,
        hsigma=3.0,
        rdnoise="0.",
        gain="1.",
        snoise="0.",
        sigscale=0.1,
        pclip=-0.5,
        grow=0.0,
    )
    return

def combinePyraf(input, output, median, scale):
    """
    Combines given images using IRAF's imcombine.
    """
    from pyraf import iraf
    from pyraf.iraf import imred

    reject = "sigclip"
    combine = ""
    if median: combine = "median"
    else: combine = "average"

    iraf.imcombine(input=input, output=output, combine=combine, reject=reject, scale=scale)

def avg_images(images, out):
    """
    Create an average image from a list of images using imcombine.
    Writes to the file out.
    """
    filelist = "tmp.txt"
    f = open(filelist, "w")
    for image in images:
        f.write("%s\n" % image.filename)
    f.close()
    iraf.imcombine(input="@%s" % filelist, output=out, Stdout=1)
    os.remove(filelist)
    return out