#!/usr/bin/env python

import os,sys,numpy,pyfits,scipy.stats


if __name__ == "__main__":

    sky_levels = []
    for fn in sys.argv[1:]:

        hdulist = pyfits.open(fn)
        data = hdulist[0].data

        good = numpy.isfinite(data) ##numpy.ones_like(data, dtype=numpy.bool)
        for iter in range(3):
            
            _s = numpy.percentile(data[good], [16,50,84])
            median = _s[1]
            sigma = 0.5*(_s[2]-_s[0])

            bad = (data < median-3*sigma) | (data > median+3*sigma)
            good[bad] = False

            # print iter+1, median, sigma

        # sky = numpy.min(data)
        print fn, median
        sky_levels.append(median)

    print
    print " ".join(["%.3f" % f for f in sky_levels])
    print ",".join(["%.3f" % f for f in sky_levels])