# Path: ugali/utils/constants.py
# MAGLIMS = dict(
#     ### ADW: Need to make sure everything matches
#     # 95% completeness from ERD
#     sdss = {
#         'u': 22.0,
#         'g': 22.2,
#         'r': 22.2,
#         'i': 21.3,
#         'z': 20.5
#     },
#     # 10 sigma mag-limit
#     dr10 = {
#         'u': 22.12,
#         'g': 22.60,
#         'r': 22.29,
#         'i': 21.85,
#         'z': 20.35
#     },
#     des = {
#         'g': 24.6,
#         'r': 24.1,
#         'i': 24.3,
#         'z': 23.8,
#         'Y': 21.5
#     },
#     sva1_gold = {
#         'g': 22.5,
#         'r': 22.5,
#         'i': 22.5,
#         'z': 22.0,
#         'Y': 21.0
#     },
#     y1a1 = {
#         #'g': 23.7, # Probably closer to 23.9 (magerr) or 24.2 (mangle)
#         #'r': 23.7, # Probably closer to 23.9 (mangle)
#         'g': 23.0, # Truncate at 23.0 (star-galaxy separation)
#         'r': 23.0, # Truncate at 23.0
#         },
#     y2n = {
#         'g': 23.0,
#         'r': 23.0,
#         },
#     y2u1 = {
#         'g': 24.0,
#         'r': 24.0,
#         },
#     y2q1 = {
#         'g': 23.5,
#         'r': 23.5,
#         },
#     y17v1 = {
#         'g': 23.5,
#         'r': 23.5,
#         },
#     y3a2 = {
#         'g': 24.0, 
#         'r': 24.0,
#         },
#     ps1_dr1 = {
#         'g': 23.0,
#         'r': 23.0,
#         },
# )


import numpy as np
import astropy.io.fits as pyfits
import pylab as plt
import glob
import scipy.optimize
from ugali.utils.constants import MAGLIMS 

# For DR10
#p0 = [-1.0,0.01,0.02]
# For Y1A1
p0 = [-1.0,-1.0, 0.02]

# DR10
#bins = np.linspace(16,24,100)
# Y1A1
bins = np.linspace(17,26,100)

centers = (bins[:-1]+bins[1:])/2.
sample=100
nfiles = 1000
medians = {b:[] for b in bands}

def expf(p,x,b='r'):
    x0 = maglims[b]
    # x0 really not necessary
    return np.exp(p[0]*(x0-x)+p[1])+p[2]

def errf(p,x,y,b='r'):
    ypred = expf(p,x,b)
    return np.fabs(((y-ypred)/ypred)).sum()

def relerr(p,x,y,b='r'):
    ypred = expf(p,x,b)
    return np.fabs(((y-ypred)/ypred)).sum()

def abserr(p,x,y,b='r'):

# Path: ugali/utils/constants.py
# MAGLIMS = dict(
#     ### ADW: Need to make sure everything matches
#     # 95% completeness from ERD
#     sdss = {
#         'u': 22.0,
#         'g': 22.2,
#         'r': 22.2,
#         'i': 21.3,
#         'z': 20.5
#     },
#     # 10 sigma mag-limit
#     dr10 = {
#         'u': 22.12,
#         'g': 22.60,
#         'r': 22.29,
#         'i': 21.85,
#         'z': 20.35
#     },
#     des = {
#         'g': 24.6,
#         'r': 24.1,
#         'i': 24.3,
#         'z': 23.8,
#         'Y': 21.5
#     },
#     sva1_gold = {
#         'g': 22.5,
#         'r': 22.5,
#         'i': 22.5,
#         'z': 22.0,
#         'Y': 21.0
#     },
#     y1a1 = {
#         #'g': 23.7, # Probably closer to 23.9 (magerr) or 24.2 (mangle)
#         #'r': 23.7, # Probably closer to 23.9 (mangle)
#         'g': 23.0, # Truncate at 23.0 (star-galaxy separation)
#         'r': 23.0, # Truncate at 23.0
#         },
#     y2n = {
#         'g': 23.0,
#         'r': 23.0,
#         },
#     y2u1 = {
#         'g': 24.0,
#         'r': 24.0,
#         },
#     y2q1 = {
#         'g': 23.5,
#         'r': 23.5,
#         },
#     y17v1 = {
#         'g': 23.5,
#         'r': 23.5,
#         },
#     y3a2 = {
#         'g': 24.0, 
#         'r': 24.0,
#         },
#     ps1_dr1 = {
#         'g': 23.0,
#         'r': 23.0,
#         },
# )


# Path: ugali/scratch/maglims.py
import numpy as np
import astropy.io.fits as pyfits
import pylab as plt
import glob
import scipy.optimize
from ugali.utils.constants import MAGLIMS 


# For DR10
#p0 = [-1.0,0.01,0.02]
# For Y1A1
p0 = [-1.0,-1.0, 0.02]

# DR10
#bins = np.linspace(16,24,100)
# Y1A1
bins = np.linspace(17,26,100)

centers = (bins[:-1]+bins[1:])/2.
sample=100
nfiles = 1000
medians = {b:[] for b in bands}

def expf(p,x,b='r'):
    x0 = maglims[b]
    # x0 really not necessary
    return np.exp(p[0]*(x0-x)+p[1])+p[2]

def errf(p,x,y,b='r'):
    ypred = expf(p,x,b)
    return np.fabs(((y-ypred)/ypred)).sum()

def relerr(p,x,y,b='r'):
    ypred = expf(p,x,b)
    return np.fabs(((y-ypred)/ypred)).sum()

def abserr(p,x,y,b='r'):

    ypred = expf(p,x,b)