def test_plot_plugin():
    viewer = ImageViewer(data.moon())
    plugin = PlotPlugin(image_filter=lambda x: x)
    viewer += plugin

    assert_equal(viewer.image, data.moon())
    plugin._update_original_image(data.coins())
    assert_equal(viewer.image, data.coins())
    viewer.close()

 def coins(self):
     """Prepare some example frames using images from the skimage
     library."""
     coins = np.array([data.coins() for i in range(0, 3*61)])
     coins = coins.reshape(3, 61, *data.coins().shape)
     # Adjust each frame to mimic an X-ray edge with a sigmoid
     S = 1/(1+np.exp(-(self.K_Es-8353))) + 0.1*np.sin(4*self.K_Es-4*8353)
     coins = (coins * S.reshape(3, 61,1,1))
     # Add some noise otherwise some functions div by zero.
     coins = coins * (0.975 + np.random.rand(*coins.shape)/20)
     coins = coins.astype(np.int32)
     return coins

def scikit_example_plot_label():
    image = data.coins()[50:-50, 50:-50]
    
    # apply threshold
    thresh = threshold_otsu(image)
    bw = closing(image > thresh, square(3))
    
    # remove artifacts connected to image border
    cleared = bw.copy()
    clear_border(cleared)
    
    # label image regions
    label_image = label(cleared)
    borders = np.logical_xor(bw, cleared)
    label_image[borders] = -1
    
    fig, ax = plt.subplots(ncols=1, nrows=1, figsize=(6, 6))
    ax.imshow(label_image, cmap='jet')
    
    for region in regionprops(label_image, ['Area', 'BoundingBox']):
    
        # skip small images
        if region['Area'] < 100:
            continue
    
        # draw rectangle around segmented coins
        minr, minc, maxr, maxc = region['BoundingBox']
        rect = mpatches.Rectangle((minc, minr), maxc - minc, maxr - minr,
                                  fill=False, edgecolor='red', linewidth=2)
        ax.add_patch(rect)
    
    plt.show()

def test_uniform_mode():
    """Verify the computed BRIEF descriptors with expected for uniform mode."""
    img = data.coins()

    keypoints = corner_peaks(corner_harris(img), min_distance=5, threshold_abs=0, threshold_rel=0.1)

    extractor = BRIEF(descriptor_size=8, sigma=2, mode="uniform")

    extractor.extract(img, keypoints[:8])

    expected = np.array(
        [
            [False, False, False, True, True, True, False, False],
            [True, True, True, False, True, False, False, True],
            [True, True, True, False, True, True, False, True],
            [True, True, True, True, False, True, False, True],
            [True, True, True, True, True, True, False, False],
            [True, True, True, True, True, True, True, True],
            [False, False, False, True, True, True, True, True],
            [False, True, False, True, False, True, True, True],
        ],
        dtype=bool,
    )

    assert_array_equal(extractor.descriptors, expected)

def test_isodata_coins_image():
    coins = skimage.img_as_ubyte(data.coins())

    threshold = threshold_isodata(coins)
    assert np.floor((coins[coins <= threshold].mean() + coins[coins > threshold].mean()) / 2.0) == threshold
    assert threshold == 107

    assert threshold_isodata(coins, return_all=True) == [107]

def main():
    """Load image, calculate optimal threshold, binarize, plot."""
    # load image
    img = data.coins()
    height, width = img.shape
    nb_pixels = height * width

    # precalculate some values for speedup
    # average pixel value
    g_avg = np.average(img)
    # P(pixel-value), i.e. #pixels-with-value / #all-pixels
    p_g = [0] * 256
    for g in range(0, 256):
        p_g[g] = np.sum(img == g) / nb_pixels

    # Otsu method
    # calculations are based on standard formulas
    q_best = None
    threshold_best = None
    img_bin_best = None
    # iterate over all possible thresholds
    for t in range(1, 255):
        img_bin = np.zeros(img.shape)
        img_bin[img >= t] = 1

        p1 = np.sum(img_bin) / nb_pixels
        p0 = 1 - p1

        g0 = np.average(img[img_bin == 0]) if np.sum(img[img_bin == 0]) > 0 else 0
        g1 = np.average(img[img_bin == 1]) if np.sum(img[img_bin == 1]) > 0 else 0

        var0 = sum([(g-g0)**2 * p_g[g] for g in range(0, t+1)])
        var1 = sum([(g-g1)**2 * p_g[g] for g in range(t+1, 256)])

        var_between = p0 * (g0 - g_avg)**2 + p1 * (g1 - g_avg)**2
        var_inner = p0 * var0**2 + p1 * var1**2

        # q is the relation of variance between classes and variance within classes
        q = var_between / var_inner if var_inner > 0 else 0

        print(t, p0, p1, g0, g1, g_avg, var_between, var_inner, q)
        if q_best is None or q_best < q:
            q_best = q
            threshold_best = t
            img_bin_best = img <= t

    # ground truth, based on scikit-image
    gt_tresh = skifilters.threshold_otsu(img)
    ground_truth = img <= gt_tresh

    # plot
    util.plot_images_grayscale(
        [img, img_bin_best, ground_truth],
        ["Image", "Otsu", "Otsu (Ground Truth)"]
    )

def only_phase():
    im=data.coins()
    
    imf=np.fft.fft2(im)
    amp=np.abs(imf)
    phase=np.angle(imf)
    
    onlyample=np.uint8(np.abs(np.fft.ifft2(amp)))
    io.imsave('onlyample.png',onlyample)
    onlyphase=np.uint8(np.mean(amp)*np.abs(np.fft.ifft2(np.exp(1j*phase))))
    io.imsave('onlyphase.png',onlyphase)

def test_li_coins_image():
    image = skimage.img_as_ubyte(data.coins())
    threshold = threshold_li(image)
    ce_actual = _cross_entropy(image, threshold)
    assert 94 < threshold_li(image) < 95
    assert ce_actual < _cross_entropy(image, threshold + 1)
    # in the case of the coins image, the minimum cross-entropy is achieved one
    # threshold below that found by the iterative method. Not sure why that is
    # but `threshold_li` does find the stationary point of the function (ie the
    # tolerance can be reduced arbitrarily but the exact same threshold is
    # found), so my guess some kind of histogram binning effect.
    assert ce_actual < _cross_entropy(image, threshold - 2)

 def test_mosiac_reference(self):
     """Check that a repeated mosaic of images is converted to optical
     depth.
     
     """
     # Prepare data
     coins = data.coins()
     mosaic = np.tile(coins, (4, 4))
     ref = np.random.rand(*coins.shape) + 1
     od = -np.log(coins/ref)
     expected = np.tile(od, (4, 4))
     # Call the reference correction function
     result = apply_mosaic_reference(mosaic, ref)
     np.testing.assert_almost_equal(result, expected)

def _cv_main():
    from skimage.data import camera
    from skimage.data import coins
    #U0 = plt.imread("rgb_tile_014_i01_j05 - crop.png")[:,:,0].astype('float')-0.5
    U0 = coins().astype('float')/255.
    print np.max(U0)
    cv = chan_vese(U0,mu=0.8,lambda1=1,lambda2=1,tol=1,maxiter=15,dt=100)
    print ("Chan-Vese algorithm finished after "+str(len(cv[1]))+" iterations.")
    print cv[1]
    plt.imshow(cv[0])
    plt.colorbar()
    plt.show()
    plt.plot(cv[1])
    plt.show()
    return

def test_minsize():
    # single-channel:
    img = data.coins()[20:168, 0:128]
    for min_size in np.arange(10, 100, 10):
        segments = felzenszwalb(img, min_size=min_size, sigma=3)
        counts = np.bincount(segments.ravel())
        # actually want to test greater or equal.
        assert_greater(counts.min() + 1, min_size)
    # multi-channel:
    coffee = data.coffee()[::4, ::4]
    for min_size in np.arange(10, 100, 10):
        segments = felzenszwalb(coffee, min_size=min_size, sigma=3)
        counts = np.bincount(segments.ravel())
        # actually want to test greater or equal.
        assert_greater(counts.min() + 1, min_size)

def test_minsize():
    # single-channel:
    img = data.coins()[20:168,0:128]
    for min_size in np.arange(10, 100, 10):
        segments = felzenszwalb(img, min_size=min_size, sigma=3)
        counts = np.bincount(segments.ravel())
        # actually want to test greater or equal.
        assert_greater(counts.min() + 1, min_size)
    # multi-channel:
    coffee = data.coffee()[::4, ::4]
    for min_size in np.arange(10, 100, 10):
        segments = felzenszwalb(coffee, min_size=min_size, sigma=3)
        counts = np.bincount(segments.ravel())
        # actually want to test greater or equal.
        # the construction doesn't guarantee min_size is respected
        # after intersecting the sementations for the colors
        assert_greater(np.mean(counts) + 1, min_size)

def wiener_filter():
    im=data.coins()
    imf=np.fft.fft2(im)

    kernel=np.ones((1,20))
    kernel=kernel/np.sum(kernel)
    kf=np.fft.fft2(kernel,(im.shape[0],im.shape[1]))
    
    g=imf*kf
    im_g=np.uint8(np.abs(np.fft.ifft2(g)))
    h=np.fft.fft2(im_g)*np.conj(kf)/(0.001+np.abs(kf)**2)
    
    io.imsave('coins-wiener.png',np.uint8(np.abs(np.fft.ifft2(h))))
    
    im_gn=np.uint8(util.noise.random_noise(im_g,var=0.00001)*255)
    h=np.fft.fft2(im_gn)*np.conj(kf)/(0.001+np.abs(kf)**2)
    
    io.imsave('coins-wiener-noise.png',np.uint8(np.abs(np.fft.ifft2(h))))

def test_viewer():
    astro = data.astronaut()
    coins = data.coins()

    view = ImageViewer(astro)
    import tempfile
    _, filename = tempfile.mkstemp(suffix='.png')

    view.show(False)
    view.close()
    view.save_to_file(filename)
    view.open_file(filename)
    assert_equal(view.image, astro)
    view.image = coins
    assert_equal(view.image, coins),
    view.save_to_file(filename),
    view.open_file(filename),
    view.reset_image(),
    assert_equal(view.image, coins)

def test_viewer_with_overlay():
    img = data.coins()
    ov = OverlayPlugin(image_filter=sobel)
    viewer = ImageViewer(img)
    viewer += ov

    import tempfile
    _, filename = tempfile.mkstemp(suffix='.png')

    ov.color = 3
    assert_equal(ov.color, 'yellow')
    viewer.save_to_file(filename)
    ov.display_filtered_image(img)
    assert_equal(ov.overlay, img)
    ov.overlay = None
    assert_equal(ov.overlay, None)
    ov.overlay = img
    assert_equal(ov.overlay, img)
    assert_equal(ov.filtered_image, img)

def test_normal_mode():
    """Verify the computed BRIEF descriptors with expected for normal mode."""
    img = data.coins()

    keypoints = corner_peaks(corner_harris(img), min_distance=5)

    extractor = BRIEF(descriptor_size=8, sigma=2)

    extractor.extract(img, keypoints[:8])

    expected = np.array([[False,  True, False, False,  True, False,  True, False],
                         [ True, False,  True,  True, False,  True, False, False],
                         [ True, False, False,  True, False,  True, False,  True],
                         [ True,  True,  True,  True, False,  True, False,  True],
                         [ True,  True,  True, False, False,  True,  True,  True],
                         [False, False, False, False,  True, False, False, False],
                         [False,  True, False, False,  True, False,  True, False],
                         [False, False, False, False, False, False, False, False]], dtype=bool)

    assert_array_equal(extractor.descriptors, expected)

def geometry_transform():
    im=data.coins()
    imtf=np.zeros(im.shape,dtype=np.uint8)
    x,y=im.shape
    s=y/2
    for i in range(x):
        for j in range(y):
            newj=s+np.sign(j-s)*(np.abs(j-s)/s)**2*s
            if newj>=0 and newj<y:
                imtf[i,j]=im[i,np.int(newj)]
    io.imsave('coins-tf.png',imtf)
    
    iminv=np.zeros(im.shape,dtype=np.uint8)
    for i in range(x):
        for j in range(y):
            newj=s+np.sign(j-s)*np.sqrt(np.abs(j-s)/s)*s
            if newj>=0 and newj<y:
                iminv[i,j]=imtf[i,np.int(newj)]
    io.imsave('coins-inverse-tf.png',iminv)
    
    

def inverse_filter():
    im=data.coins()
    imf=np.fft.fft2(im)

    kernel=np.ones((1,20))
    kernel=kernel/np.sum(kernel)
    kf=np.fft.fft2(kernel,(im.shape[0],im.shape[1]))
    
    g=imf*kf
    im_g=np.uint8(np.abs(np.fft.ifft2(g)))
    h=np.fft.fft2(im_g)/(0.01+kf)
    
    io.imsave('coins.png',im)
    io.imsave('coins-blur.png',im_g)
    io.imsave('coins-deblur.png',np.uint8(np.abs(np.fft.ifft2(h))))
    
    im_gn=np.uint8(util.noise.random_noise(im_g,var=0.00001)*255)
    h=np.fft.fft2(im_gn)/(0.01+kf)
    
    io.imsave('coins-blur-noise.png',im_gn)
    io.imsave('coins-deblur-noise.png',np.uint8(np.abs(np.fft.ifft2(h))))

    def draw(self):
        if not hasattr(self, 'ax'):
            self.axOriginal = self.figure.add_subplot(221)
            self.axGreyScale = self.figure.add_subplot(222)
            self.axFiltered = self.figure.add_subplot(223)
            self.axSegments = self.figure.add_subplot(224)
            self.image = data.coins()
#            self.image = imread(self.image_file)

        self.axOriginal.set_title("Original Image", fontsize=12)
        self.axOriginal.imshow(self.image)
        
        self.axGreyScale.set_title("Greyscale Image", fontsize=12)
        self.grey_image = color.rgb2grey(self.image)
        self.axGreyScale.imshow(self.grey_image, cmap = cm.Greys_r)

        self.filter_image()
#        thresh = threshold_otsu(self.grey_image)
#        self.bw_image = closing(self.grey_image > thresh, square(1))
#        self.axThreshold.imshow(self.bw_image)
        
        self.axSegments.set_title("Segmented Image", fontsize=12)
        
        self.label_image = label(self.filtered)
#        borders = np.logical_xor(self.bw_image, self.cleared)
#        self.label_image[borders] = -1
        
#        fig, ax = plt.subplots(ncols=1, nrows=1, figsize=(6, 6))
        self.axSegments.imshow(self.label_image, cmap='jet')
        
        for region in regionprops(self.label_image, ['Area', 'BoundingBox']):
        
            # skip small images
            if region['Area'] < 100:
                continue
        
            # draw rectangle around segmented coins
            minr, minc, maxr, maxc = region['BoundingBox']
            rect = mpatches.Rectangle((minc, minr), maxc - minc, maxr - minr, fill=False, edgecolor='red', linewidth=2)
            self.axSegments.add_patch(rect)

def profile():
    import time
    from iib.simulation import CLContext
    from skimage import io, data, transform
    gs, wgs = 256, 16

    # Load some test data
    r = transform.resize
    sigs = np.empty((gs, gs, 4), np.float32)
    sigs[:, :, 0] = r(data.coins().astype(np.float32) / 255.0, (gs, gs))
    sigs[:, :, 1] = r(data.camera().astype(np.float32) / 255.0, (gs, gs))
    sigs[:, :, 2] = r(data.text().astype(np.float32) / 255.0, (gs, gs))
    sigs[:, :, 3] = r(data.checkerboard().astype(np.float32) / 255.0, (gs, gs))
    sigs[:, :, 2] = r(io.imread("../scoring/corpus/rds/turing_001.png",
                                as_grey=True), (gs, gs))
    sigs[:, :, 3] = io.imread("../scoring/corpus/synthetic/blobs.png",
                              as_grey=True)
    sigs = sigs.reshape(gs*gs*4)

    # Set up OpenCL
    ctx = cl.create_some_context(interactive=False)
    queue = cl.CommandQueue(ctx)
    mf = cl.mem_flags
    ifmt_f = cl.ImageFormat(cl.channel_order.RGBA, cl.channel_type.FLOAT)
    bufi = cl.Image(ctx, mf.READ_ONLY, ifmt_f, (gs, gs))
    cl.enqueue_copy(queue, bufi, sigs, origin=(0, 0), region=(gs, gs))
    clctx = CLContext(ctx, queue, ifmt_f, gs, wgs)

    # Compile the kernels
    feats = cl.Program(ctx, features_cl()).build()
    rdctn = cl.Program(ctx, reduction.reduction_sum_cl()).build()
    blur2 = cl.Program(ctx, convolution.gaussian_cl([np.sqrt(2.0)]*4)).build()
    blur4 = cl.Program(ctx, convolution.gaussian_cl([np.sqrt(4.0)]*4)).build()

    iters = 500
    t0 = time.time()
    for i in range(iters):
        get_features(clctx, feats, rdctn, blur2, blur4, bufi)
    print((time.time() - t0)/iters)