def align_img(image_name, loss_fn=sum_of_squared_diff, max_disp = 15, big=False, name=None, hist_eq = False):
    b,g,r = get_bgr(image_name)

    if hist_eq == True:
        b, g, r = exposure.equalize_hist(b), exposure.equalize_hist(g), exposure.equalize_hist(r)

    print("Aligning green and blue: ")
    ag = align(g, b, loss_fn=loss_fn, max_disp = max_disp, big=big, name=name)
    print("Aligning blue and red: ")
    ar = align(r, b, loss_fn=loss_fn, max_disp = max_disp, big=big, name=name)
    # create a color image
    im_out = np.dstack([ar, ag, b])
    plt.show()

    # save the image
    iname = image_name.split('.')
    iname[-1] = 'jpg'
    image_name = '.'.join(iname)
    fname = 'out_' + image_name
    skio.imsave(fname, im_out)

    # display the image
    skio.imshow(im_out)
    plt.show()
    skio.show()

def cutouts_plot(outpath, base_name, order_num, obj, flat, top_trace, bot_trace, trace):
    
    pl.figure('traces', facecolor='white', figsize=(8, 5))
    pl.cla()
    pl.suptitle('order cutouts, {}, order {}'.format(base_name, order_num), fontsize=14)
    pl.set_cmap('Blues_r')

    obj_plot = pl.subplot(2, 1, 1)
    try:
        obj_plot.imshow(exposure.equalize_hist(obj))
    except:
        obj_plot.imshow(obj)
    obj_plot.plot(np.arange(1024), top_trace, 'y-', linewidth=1.5)
    obj_plot.plot(np.arange(1024), bot_trace, 'y-', linewidth=1.5)
    obj_plot.plot(np.arange(1024), trace, 'y-', linewidth=1.5)
    obj_plot.set_title('object')
#     obj_plot.set_ylim([1023, 0])
    obj_plot.set_xlim([0, 1023])
    
    flat_plot = pl.subplot(2, 1, 2)
    try:
        flat_plot.imshow(exposure.equalize_hist(flat))
    except:
        flat_plot.imshow(flat)
    flat_plot.plot(np.arange(1024), top_trace, 'y-', linewidth=1.5)
    flat_plot.plot(np.arange(1024), bot_trace, 'y-', linewidth=1.5)    
    flat_plot.plot(np.arange(1024), trace, 'y-', linewidth=1.5)    
    flat_plot.set_title('flat')
#     flat_plot.set_ylim([1023, 0])
    flat_plot.set_xlim([0, 1023])
 
    pl.tight_layout()
    pl.savefig(constructFileName(outpath, base_name, order_num, 'cutouts.png'))
    pl.close()

def sparect_plot(outpath, base_name, order_num, obj, flat):

    pl.figure('spatially rectified', facecolor='white', figsize=(8, 5))
    pl.cla()
    pl.suptitle('spatially rectified, {}, order {}'.format(base_name, order_num), fontsize=14)
    pl.set_cmap('Blues_r')

    obj_plot = pl.subplot(2, 1, 1)
    try:
        obj_plot.imshow(exposure.equalize_hist(obj))
    except:
        obj_plot.imshow(obj)
    obj_plot.set_title('object')
#     obj_plot.set_ylim([1023, 0])
    obj_plot.set_xlim([0, 1023])
    
    flat_plot = pl.subplot(2, 1, 2)
    try:
        flat_plot.imshow(exposure.equalize_hist(flat))
    except:
        flat_plot.imshow(flat)
    flat_plot.set_title('flat')
#     flat_plot.set_ylim([1023, 0])
    flat_plot.set_xlim([0, 1023])
 
    pl.tight_layout()
    pl.savefig(constructFileName(outpath, base_name, order_num, 'sparect.png'))
    pl.close()

def traces_plot(outpath, base_name, order_num, obj, flat, top_trace, bot_trace):
    
    pl.figure('traces', facecolor='white', figsize=(8, 5))
    pl.cla()
    pl.suptitle('order edge traces, {}, order {}'.format(base_name, order_num), fontsize=14)
    pl.set_cmap('Blues_r')
    pl.rcParams['ytick.labelsize'] = 8

    obj_plot = pl.subplot(1, 2, 1)
    obj_plot.imshow(exposure.equalize_hist(obj))
    obj_plot.plot(np.arange(1024), top_trace, 'y-', linewidth=1.5)
    obj_plot.plot(np.arange(1024), bot_trace, 'y-', linewidth=1.5)

    obj_plot.set_title('object')
    obj_plot.set_ylim([1023, 0])
    obj_plot.set_xlim([0, 1023])

    
    flat_plot = pl.subplot(1, 2, 2)
    flat_plot.imshow(exposure.equalize_hist(flat))
    flat_plot.plot(np.arange(1024), top_trace, 'y-', linewidth=1.5)
    flat_plot.plot(np.arange(1024), bot_trace, 'y-', linewidth=1.5)    
    flat_plot.set_title('flat')
    flat_plot.set_ylim([1023, 0])
    flat_plot.set_xlim([0, 1023])
 
    pl.tight_layout()
    pl.savefig(constructFileName(outpath, base_name, order_num, 'traces.png'))
    pl.close()

def tops_bots_plot(outpath, base_name, tops, bots):
    
    pl.figure('edges', facecolor='white', figsize=(8, 5))
    pl.cla()
    pl.suptitle('top and bottom order edges, {}'.format(base_name), fontsize=14)
#     pl.set_cmap('Blues_r')
    pl.rcParams['ytick.labelsize'] = 8

    obj_plot = pl.subplot(1, 2, 1)
    try:
        obj_plot.imshow(exposure.equalize_hist(tops))
    except:
        obj_plot.imshow(tops)

    obj_plot.set_title('top edges')
    obj_plot.set_ylim([1023, 0])
    obj_plot.set_xlim([0, 1023])

    
    flat_plot = pl.subplot(1, 2, 2)
    try:
        flat_plot.imshow(exposure.equalize_hist(bots))
    except:
        flat_plot.imshow(bots)
    flat_plot.set_title('bottom edges')
    flat_plot.set_ylim([1023, 0])
    flat_plot.set_xlim([0, 1023])
 
    pl.tight_layout()
    pl.savefig(constructFileName(outpath, base_name, None, 'top_bot_edges.png'))
    pl.close()

def test_equalize_masked():
    img = skimage.img_as_float(test_img)
    mask = np.zeros(test_img.shape)
    mask[50:150, 50:250] = 1
    img_mask_eq = exposure.equalize_hist(img, mask=mask)
    img_eq = exposure.equalize_hist(img)

    cdf, bin_edges = exposure.cumulative_distribution(img_mask_eq)
    check_cdf_slope(cdf)

    assert not (img_eq == img_mask_eq).all()

def slidingWindow(P,inX=3,outX=32,inY=3,outY=64,maxM=50,norm=True):
	""" Enhance the constrast

		Cut off extreme values and demean the image
		Utilize scipy convolve2d to get the mean at a given pixel
		Remove local mean with inner exclusion region

		Args:
			P: 2-d numpy array image
			inX: inner exclusion region in the x-dimension
			outX: length of the window in the x-dimension
			inY: inner exclusion region in the y-dimension
			outY: length of the window in the y-dimension
			maxM: size of the output image in the y-dimension
			norm: boolean to cut off extreme values

		Returns:
			Q: 2-d numpy contrast enhanced
	"""
	Q = P.copy()
	m, n = Q.shape
	
	Q = exposure.equalize_hist(Q.astype('Float32'), nbins = 65)
	Q = detrend(Q.astype('Float32'), axis = 1)
	Q = detrend(Q.astype('Float32'), axis = 0)
	Q = wiener(Q.astype('Float32'), 4)

	return Q[:maxM,:]

def twoDimOrderPlot(outpath, base_name, title, obj_name, base_filename, order_num, data, x):
    pl.figure('2d order image', facecolor='white', figsize=(8, 5))
    pl.cla()
    pl.title(title + ', ' + base_name + ", order " + str(order_num), fontsize=14)
    pl.xlabel('wavelength($\AA$)', fontsize=12)
    pl.ylabel('row (pixel)', fontsize=12)
    #pl.imshow(img, aspect='auto')
    #pl.imshow(data, vmin=0, vmax=1024, aspect='auto')
    
    pl.imshow(exposure.equalize_hist(data), origin='lower', 
                  extent=[x[0], x[-1], 0, data.shape[0]], aspect='auto')      
#     from matplotlib import colors
#     norm = colors.LogNorm(data.mean() + 0.5 * data.std(), data.max(), clip='True')
#     pl.imshow(data, norm=norm, origin='lower',
#                   extent=[x[0], x[-1], 0, data.shape[0]], aspect='auto')               
    pl.colorbar()
    pl.set_cmap('jet')
#     pl.set_cmap('Blues_r')
    fn = constructFileName(outpath, base_name, order_num, base_filename)
    pl.savefig(fn)
    log_fn(fn)
    pl.close()
    
#     np.save(fn[:fn.rfind('.')], data)
    
    return

def mirror_mirror():
    pca, clf = gen_classifier()
    cam = cv2.VideoCapture(0)

    last_emo = []
    while True:
        ret, img = cam.read()
        if not ret:
            continue

        gray = cv2.cvtColor(img, cv.CV_BGR2GRAY)
        faces = detect_and_scale_face(gray)
        if not faces:
            continue

        for scaled in faces:
            test_x = exposure.equalize_hist(scaled.reshape((1, -1)))
            face_pca = pca.transform(test_x)
            pred = int(clf.predict(face_pca)[0])
            last_emo.append(pred)

            if len(last_emo) == 4:
                last_emo.pop(0)
                if last_emo[0] != OTHER_LABEL and \
                   np.all(np.array(last_emo) == last_emo[0]):
                    sb.call(['say', TARGET_NAMES[pred]])

def preprocessing(image):
    # Contrast stretching
    #TODO TEST
    #p2, p98 = np.percentile(image, (2, 98))
    #image_rescale = exposure.rescale_intensity(image, in_range=(p2, p98))
    gray = rgb2gray(image)
    return exposure.equalize_hist(gray)

    def preproc(self, img, size, pixel_spacing, equalize=True, crop=True):
       """crop center and resize"""
        #    TODO: this is stupid, you could crop out the heart
        # But should test this
       if img.shape[0] < img.shape[1]:
           img = img.T
       # Standardize based on pixel spacing
       img = transform.resize(img, (int(img.shape[0]*(1.0/np.float32(pixel_spacing[0]))), int(img.shape[1]*(1.0/np.float32(pixel_spacing[1])))))
       # we crop image from center
       short_egde = min(img.shape[:2])
       yy = int((img.shape[0] - short_egde) / 2)
       xx = int((img.shape[1] - short_egde) / 2)
       if crop:
           crop_img = img[yy : yy + short_egde, xx : xx + short_egde]
       # resize to 64, 64
           resized_img = transform.resize(crop_img, (size, size))
       else:
           resized_img = img
       #resized_img = gaussian_filter(resized_img, sigma=1)
       #resized_img = median_filter(resized_img, size=(3,3))
       if equalize:
           resized_img = equalize_hist(resized_img)
           resized_img = adjust_sigmoid(resized_img)
       resized_img *= 255.

       return resized_img.astype("float32")

def alineacion(image,rgbCopia,escala):
    centerPoint=(140,165)    
    extents=15
    #if escala>1:
    #        extents=extents*escala
    cropImg=np.empty((extents*2,extents*2,3),dtype=float)   
    cropSizes=[]     
    cropSizes.append([0,0,0,0])
    cropImg[:,:,0]=cropImage(image[:,:,0],centerPoint,extents)
    for i in range(image.shape[2]-1):
        cropImg[:,:,i+1]=cropImage(image[:,:,i+1],centerPoint,extents*escala)[::escala,::escala]
        cropImg[:,:,i+1]=cropImg[:,:,i+1]-np.mean(cropImg[:,:,i+1])
    channelR=cropImg[:,:,0]
    for i in range(2):
        channelRGradient=np.gradient(channelR)
        cropGradient=np.gradient(cropImg[:,:,i+1])
        #correlation=correlationMatrix(channelR,cropImg[:,:,i+1])
        correlation=correlationMatrix(channelR,cropImg[:,:,i+1])        
        correlation=correlation-np.mean(correlation)
        print "Filtro size ",correlation.shape
        position=np.where(correlation==np.amax(correlation))
        despVector=[position[0][0]-correlation.shape[0]/2,position[1][0]-correlation.shape[1]/2]
        
        fila=despVector[0]
        columna=despVector[1]
        
        print "fila",fila
        print "columna", columna
        
        generateNewRGB(i,rgbCopia,fila,columna)
    
    rgbCopia=equalize_hist(rgbCopia)
    rgbCopia=rgbCopia*255/np.max(rgbCopia)
        
    return rgbCopia

    def compute(self, src):
        image = img_as_ubyte(src)

        # denoise image
        denoised = denoise_tv_chambolle(image, weight=0.05)
        denoised_equalize= exposure.equalize_hist(denoised)

        # find continuous region (low gradient) --> markers
        markers = rank.gradient(denoised_equalize, disk(5)) < 10
        markers = ndi.label(markers)[0]

        # local gradient
        gradient = rank.gradient(denoised, disk(2))

        # labels
        labels = watershed(gradient, markers)

        # display results
        fig, axes = plt.subplots(2,3)
        axes[0, 0].imshow(image)#, cmap=plt.cm.spectral, interpolation='nearest')
        axes[0, 1].imshow(denoised, cmap=plt.cm.spectral, interpolation='nearest')
        axes[0, 2].imshow(markers, cmap=plt.cm.spectral, interpolation='nearest')
        axes[1, 0].imshow(gradient, cmap=plt.cm.spectral, interpolation='nearest')
        axes[1, 1].imshow(labels, cmap=plt.cm.spectral, interpolation='nearest', alpha=.7)
        plt.show()

def loadDataMontgomery(df, path, im_shape):
    """Function for loading Montgomery dataset"""
    X, y = [], []
    for i, item in df.iterrows():
        img = img_as_float(io.imread(path + item[0]))
        gt = io.imread(path + item[1])
        l, r = np.where(img.sum(0) > 1)[0][[0, -1]]
        t, b = np.where(img.sum(1) > 1)[0][[0, -1]]
        img = img[t:b, l:r]
        mask = gt[t:b, l:r]
        img = transform.resize(img, im_shape)
        img = exposure.equalize_hist(img)
        img = np.expand_dims(img, -1)
        mask = transform.resize(mask, im_shape)
        mask = np.expand_dims(mask, -1)
        X.append(img)
        y.append(mask)
    X = np.array(X)
    y = np.array(y)
    X -= X.mean()
    X /= X.std()

    print '### Data loaded'
    print '\t{}'.format(path)
    print '\t{}\t{}'.format(X.shape, y.shape)
    print '\tX:{:.1f}-{:.1f}\ty:{:.1f}-{:.1f}\n'.format(X.min(), X.max(), y.min(), y.max())
    print '\tX.mean = {}, X.std = {}'.format(X.mean(), X.std())
    return X, y

    def histograms_of_stuff(cls, ref):
        # First we show the image, and its histogram
        neg = 1. - ref
        figs = []
        f1, axs1 = plt.subplots(2, 2, figsize=(8, 8))
        for i, im in enumerate(((ref, 'tire'), (neg, 'negative'))):
            axs1[i,0].imshow(im[0], cmap='gray')
            axs1[i,1].hist(im[0].ravel(), 256, range=(0., 1.), fc='k', ec='k')
            axs1[i,0].set_title(im[1])
            axs1[i,1].set_title(im[1] + ': histogram')

        figs.append(f1)

        # Some gamma correction ...
        gammas = (0.5, 1.3)
        f2, axs2 = plt.subplots(2, 2, figsize=(8,8))
        for i, gamma in enumerate(gammas):
            x = ref ** gamma
            axs2[i,0].imshow(x, cmap='gray')
            axs2[i,1].hist(x.ravel(), 256, range=(0.,1.), fc='k', ec='k')
            for z in axs2[i]:
                z.set_title('gamma = ' + str(gamma))
        figs.append(f2)

        # Histogram Equalization
        f3, axs3 = plt.subplots(1, 2, figsize=(8, 4))
        eq_ref = exposure.equalize_hist(ref)
        axs3[0].imshow(eq_ref, cmap='gray')
        axs3[1].hist(eq_ref.ravel(), 256, range=(0.,1.), fc='k', ec='k')
        for x in axs3:
            x.set_title('Histogram Equalization')
        figs.append(f3)
        return figs

 def transform(self, Xb, yb):
     Xb, yb = super(EqualizeHistBatchIteratorMixin, self).transform(Xb, yb)
     Xb_transformed = np.asarray([
         [equalize_hist(img_ch) for img_ch in img] for img in Xb
     ])
     Xb_transformed = Xb_transformed.astype(Xb.dtype)
     return Xb_transformed, yb

def visualize_depth_image(data):

    data[data == 0.0] = np.nan

    maxdepth = np.nanmax(data)
    mindepth = np.nanmin(data)
    data = data.copy()
    data -= mindepth
    data /= (maxdepth - mindepth)

    gray = np.zeros(list(data.shape) + [3], dtype=data.dtype)
    data = (1.0 - data)
    gray[..., :3] = np.dstack((data, data, data))

    # use a greenish color to visualize missing depth
    gray[np.isnan(data), :] = (97, 160, 123)
    gray[np.isnan(data), :] /= 255

    gray = exposure.equalize_hist(gray)

    # set alpha channel
    gray = np.dstack((gray, np.ones(data.shape[:2])))
    gray[np.isnan(data), -1] = 0.5

    return gray * 255

def _equalizeHistogram(img):
    '''
    histogram equalisation not bounded to int() or an image depth of 8 bit
    works also with negative numbers
    '''
    #to float if int:
    intType = None
    if 'f' not in img.dtype.str:
        TO_FLOAT_TYPES = {  np.dtype('uint8'):np.float16,
                    np.dtype('uint16'):np.float32,
                    np.dtype('uint32'):np.float64,
                    np.dtype('uint64'):np.float64} 

        intType = img.dtype
        img = img.astype(TO_FLOAT_TYPES[intType], copy=False)
    
    #get image deph
    DEPTH_TO_NBINS = {np.dtype('float16'):256, #uint8
                      np.dtype('float32'):32768, #uint16
                      np.dtype('float64'):2147483648} #uint32
    nBins = DEPTH_TO_NBINS[img.dtype]

    #scale to -1 to 1 due to skikit-image restrictions
    mn, mx = np.amin(img), np.amax(img)
    if abs(mn) > abs(mx):
        mx = mn
    img /= mx
    img = exposure.equalize_hist(img, nbins=nBins)
    
    img *= mx
    
    if intType:
        img = img.astype(intType)
    return img

def hsi_equalize_hist():
    image=data.astronaut()
    h=color.rgb2hsv(image)
    h[:,:,2]=exposure.equalize_hist(h[:,:,2])
    image_equal=color.hsv2rgb(h)
    io.imshow(image_equal)
    io.imsave('astronautequal.png',image_equal)

def test_equalize_ubyte():
    with expected_warnings(['precision loss']):
        img = skimage.img_as_ubyte(test_img)
    img_eq = exposure.equalize_hist(img)

    cdf, bin_edges = exposure.cumulative_distribution(img_eq)
    check_cdf_slope(cdf)