def demo_histeq():
    """
    Gray histogram equalization, make every gray value have same distribution.
    This process could imporve contrast of image.
    """
    im = np.array(Image.open('./computer_vision/scenery.jpg').convert('L'))
    imhist, bins = np.histogram(im.flatten(), bins = 256, normed = True)
    cdf = imhist.cumsum()
    """Normalization"""
    cdf = 255 * cdf / cdf[-1]
    im_n = np.interp(im.flatten(),bins[:-1], cdf)
    """
    Image of equalization.
    I don't know why reshape() function can't modify shape attribute, 
    it just change the arrangement of output if you print it out.
    I've tried when you create a array with reshape(), it works.
    """
    #im_n.reshape(im.shape)
    im_n.shape = im.shape
    
    pl.figure('Histogram equalization')
    pl.subplot(1,2,1)
    pl.gray()
    pl.imshow(im)
    pl.title('Orignal gray image')
    pl.subplot(1,2,2)
    pl.gray()
    pl.imshow(im_n)
    pl.title('equalized image')
    pl.show()

def explore_data(data, images, target):

	# try to determine the type of data...
	print "data_type belonging to key data:"
	try: 
		print np.dtype(data)
	except TypeError as err:
		print err
	
	print "It has dimension", np.shape(data)

	# plot a 3
	
	# get indices of all threes in target
	threes = np.where(target == 3)
	#assert threes is not empty
	assert(len(threes) > 0)
	# choose the first 3
	three_indx = threes[0]
	# get the image
	img = images[three_indx][0]

	#plot it
	plot.figure()
	plot.gray()
	plot.imshow(img, interpolation = "nearest")
	plot.show()
	plot.close()

def display_head(set_x, set_y, n = 5):
    '''
    show some figures based on gray image matrixs
    
    @type set_x: TensorSharedVariable, 
    @param set_x: gray level value matrix of the
    
    @type set_y: TensorVariable, 
    @param set_y: label of the figures    
    
    @type n: int, 
    @param n: numbers of figure to be display, less than 10, default 5
    '''
    import pylab
    
    if n > 10: n = 10
    img_x = set_x.get_value()[0:n].reshape(n, 28, 28)
    img_y = set_y.eval()[0:n]
    
    for i in range(n): 
        pylab.subplot(1, n, i+1); 
        pylab.axis('off'); 
        pylab.title(' %d' % img_y[i])
        pylab.gray()
        pylab.imshow(img_x[i])

def compare_keypoints(im1, im2, pos1, pos2, filename = None, separation = 0) :
    """ Show two images next to each other with the keypoints marked
    """

    # Construct unified image
    im3 = append_images(im1,im2, separation)

    # Find the offset and add it
    offset = im1.shape[1]
    pos2_o = [(x+offset + separation,y) for (x,y) in pos2]

    # Create figure
    fig = pylab.figure(frameon=False, figsize=(12.0, 8.0))
    #ax = pylab.Axes(fig, [0., 0., 1., 1.])

    # Show images
    pylab.gray()
    pylab.imshow(im3)
    pylab.plot([x for x,y in pos1], [y for x,y in pos1], marker='o', color = '#00aaff', lw=0)
    pylab.plot([x for x,y in pos2_o], [y for x,y in pos2_o], marker='o', color = '#00aaff', lw=0)
    pylab.axis('off')

    pylab.xlim(0,im3.shape[1])
    pylab.ylim(im3.shape[0],0)

    if filename != None :
        fig.savefig(filename, bbox_inches='tight', dpi=300)

    def updateColorTable(self, cItem):
        print "now viz!"+str(cItem.row())+","+str(cItem.column())

        row = cItem.row()
        col = cItem.column()

        pl.clf()
        #pl.ion()
        x = pl.arange(self.dataDimen+1)
        y = pl.arange(self.dataDimen+1)
        X, Y = pl.meshgrid(x, y)
        pl.subplot(1,2,1)
        pl.pcolor(X, Y, self.mWx[row*self.dataMaxRange+col])
        pl.gca().set_aspect('equal')
        pl.colorbar()
        pl.gray()
        pl.title("user 1")

        pl.subplot(1,2,2)
        pl.pcolor(X, Y, self.mWy[row*self.dataMaxRange+col])
        pl.gca().set_aspect('equal')
        pl.colorbar()
        pl.gray()
        pl.title("user 2")
        #pl.tight_layout()

        pl.draw()
        #pl.show()
        pl.show(block=False) 

def apply_conv(name):
    import pylab
    from PIL import Image

    # open random image of dimensions 639x516
    img = Image.open(open(name))

    # capture image height and width
    width,height = img.size

    # dimensions are (height, width, channel)
    img = numpy.asarray(img, dtype='float64') / 256.

    # put image in 4D tensor of shape (1, 3, height, width)
    img_ = img.transpose(2, 0, 1).reshape(1, 3, height, width)
    filtered_img = f(img_)

    # plot original image and first and second components of output
    pylab.subplot(1, 3, 1); pylab.axis('off'); pylab.imshow(img)
    pylab.gray();
    # recall that the convOp output (filtered image) is actually a "minibatch",
    # of size 1 here, so we take index 0 in the first dimension:
    pylab.subplot(1, 3, 2); pylab.axis('off'); pylab.imshow(filtered_img[0, 0, :, :])
    pylab.subplot(1, 3, 3); pylab.axis('off'); pylab.imshow(filtered_img[0, 1, :, :])
    pylab.show()

def main():
    """Compute the SSIM index on two input images specified on the cmd line."""
    import pylab
    argv = sys.argv
    if len(argv) != 3:
        print('usage: python -m sp.ssim image1.tif image2.tif', file=sys.stderr)
        sys.exit(2)

    try:
        from PIL import Image
        img1 = numpy.asarray(Image.open(argv[1]))
        img2 = numpy.asarray(Image.open(argv[2]))
    except Exception as e:
        e = 'Cannot load images' + str(e)
        print(e, file=sys.stderr)

    ssim_map = ssim(img1, img2)
    ms_ssim = msssim(img1, img2)

    pylab.figure()
    pylab.subplot(131)
    pylab.title('Image1')
    pylab.imshow(img1, interpolation='nearest', cmap=pylab.gray())
    pylab.subplot(132)
    pylab.title('Image2')
    pylab.imshow(img2, interpolation='nearest', cmap=pylab.gray())
    pylab.subplot(133)
    pylab.title('SSIM Map\n SSIM: %f\n MSSSIM: %f' % (ssim_map.mean(), ms_ssim))
    pylab.imshow(ssim_map, interpolation='nearest', cmap=pylab.gray())
    pylab.show()
    
    return 0

def test_rot180():
    import cudamat
    import numpy as np
    num_images = 100
    img_size = 50
    img_tot_size = 50*50


    inputs = np.random.randn(num_images, img_tot_size)
    inputs[1:] = (np.random.rand(*inputs[1:].shape)<0.5)
    inputs[0] = (np.random.rand(*inputs[1].shape)<0.005)

    targets = np.random.randn(num_images, img_tot_size)

    cu_inputs = cudamat.CUDAMatrix(inputs.T)
    cu_targets = cudamat.CUDAMatrix(targets.T)

    cudamat._cudamat.rot180(cu_inputs.p_mat, cu_targets.p_mat, 0)

    cua_targets = cu_targets.asarray().T

    targets = np.array([x[::-1,::-1]
                        for x in inputs.reshape(num_images, img_size, img_size)]).reshape(num_images, img_tot_size)

    print abs(targets - cua_targets).max()

    from pylab import imshow, subplot, gray
    gray()
    subplot(221)
    imshow(inputs[0].reshape(img_size, img_size), interpolation='nearest')
    subplot(222)
    imshow(targets[0].reshape(img_size, img_size), interpolation='nearest')

    subplot(223)
    imshow(cua_targets[0].reshape(img_size, img_size), interpolation='nearest')

def test_copyOutOf():
    import cudamat
    import numpy as np
    num_images = 100
    img_size = 50
    target_size = 72
    img_tot_size = img_size**2
    target_tot_size = target_size**2

    targets = np.random.randn(target_tot_size, num_images)<-2
    inputs = np.zeros((img_tot_size, num_images))

    cu_inputs = cudamat.CUDAMatrix(inputs)
    cu_targets = cudamat.CUDAMatrix(targets)

    assert (target_size - img_size) % 2 == 0
    padding = (target_size - img_size)/2
    cudamat._cudamat.copy_out_of_center(cu_targets.p_mat, cu_inputs.p_mat, padding, 0)

    cua_inputs = cu_inputs.asarray()

    #print abs(targets - cua_targets).max()

    from pylab import imshow, subplot, gray
    gray()

    #subplot(221)
    #imshow(inputs[0].reshape(img_size, img_size), interpolation='nearest')
    subplot(222)
    imshow(targets[:,1].reshape(target_size, target_size), interpolation='nearest')

    subplot(223)
    imshow(cua_inputs[:,1].reshape(img_size, img_size), interpolation='nearest')

def makeTestPair(paths, homography, collection, location=".", size=(250,250), scale = 1.0) :
    """ Given a pair of paths to two images and a homography between them,
        this function creates two crops and calculates a new homography.
        input: paths [strings] (paths to images)
               homography [numpy.ndarray] (3 by 3 array homography)
               collection [string] (The name of the testset)
               location [string] (The location (path) of the testset
               size [(int, int)] (The size of an image crop in pixels)
               scale [double] (The scale by which we resize the crops after they've been cropped)
        out:   nothing
    """
    
    # Get width and height
    width, height = size
    
    # Load images in black/white
    images = map(loadImage, paths)
    
    # Crop part of first image and part of second image:
    (top_o, left_o) = (random.randint(0, images[0].shape[0]-height), random.randint(0, images[0].shape[1]-width))
    (top_n, left_n) = (random.randint(0, images[1].shape[0]-height), random.randint(0, images[1].shape[1]-width))
    
    # Get two file names
    c_path = getRandPath("%s/%s/" % (location, collection))
    if not exists(dirname(c_path)) : makedirs(dirname(c_path))
        
    # Make sure we save as gray
    pylab.gray()
    
    im1 = images[0][top_o: top_o + height, left_o: left_o + width]
    im2 = images[1][top_n: top_n + height, left_n: left_n + width]
    im1_scaled = imresize(im1, size=float(scale), interp='bicubic')
    im2_scaled = imresize(im2, size=float(scale), interp='bicubic')
    pylab.imsave(c_path + "_1.jpg", im1_scaled)
    pylab.imsave(c_path + "_2.jpg", im2_scaled)
    
    # Homography for transpose
    T1 = numpy.identity(3)
    T1[0,2] = left_o
    T1[1,2] = top_o
    
    # Homography for transpose back
    T2 = numpy.identity(3)
    T2[0,2] = -1*left_n
    T2[1,2] = -1*top_n
    
    # Homography for scale
    Ts = numpy.identity(3)
    Ts[0,0] = scale
    Ts[1,1] = scale
    
    # Homography for scale back
    Tsinv = numpy.identity(3)
    Tsinv[0,0] = 1.0/scale
    Tsinv[1,1] = 1.0/scale
    
    # Combine homographyies and save
    hom = Ts.dot(T2).dot(homography).dot(T1).dot(Tsinv)
    hom = hom / hom[2,2]
    numpy.savetxt(c_path, hom)

def plotSources(im, sources, schema):
    plt.clf()
    plt.imshow(im.getArray(), origin='lower', interpolation='nearest',
               vmin=-100, vmax=500)
    plt.gray()
    shapekey = schema.find('shape.sdss').key
    xykey = schema.find('centroid.sdss').key

    flagkeys = [schema.find(x).key for x in [
        'shape.sdss.flags.maxiter', 'shape.sdss.flags.shift',
        'shape.sdss.flags.unweighted', 'shape.sdss.flags.unweightedbad']]
    flagabbr = ['M','S','U','B']

    for source in sources:
        quad = source.get(shapekey)
        ixx,iyy = quad.getIxx(), quad.getIyy()
        x,y = source.get(xykey)
        sx,sy = sqrt(ixx),sqrt(iyy)
        plt.gca().add_artist(Ellipse([x,y], 2.*sx, 2.*sy, angle=0.,
                                     ec='r', fc='none', lw=2, alpha=0.5))
        fs = ''
        for j,key in enumerate(flagkeys):
            val = source.get(key)
            if val:
                fs += flagabbr[j]
        if len(fs):
            plt.text(x+5, y, fs, ha='left', va='center')

def show_weights(layer):        
    while isinstance(layer, lasagne.layers.InputLayer) == False:
        if isinstance(layer, dnn.Conv2DDNNLayer):
            first_conv_layer =  layer
        layer = layer.input_layer
    
    weights = first_conv_layer.get_params()[0].get_value()
    weights_no = weights.shape[0]
    
    display_size = int(math.sqrt(weights_no)) + 1
    
    print 'display_size : %s' % display_size
    
    pylab.gray() 
    for i in range(display_size):
        for j in range(display_size):
            index = i * display_size + j + 1
            
            if index >= weights_no:
                break

            print 'index : %s' % index
    
            one_weight = weights[index][0]
            pylab.subplot(display_size, display_size, index) 
            pylab.axis('off') 
            pylab.imshow(one_weight)
    
    pylab.show()

def _test():
    # make a unit circle going counter-clockwise
    radius = 60
    theta = np.linspace(0, 2*np.pi, 10)
    circle_ccw = np.array([radius*np.cos(theta), radius*np.sin(theta)])
    area = ContourArea(circle_ccw)
    assert area > 0
    circle_cw = MakeClockwise(circle_ccw)
    area = ContourArea(circle_cw)
    assert area < 0
    assert (circle_cw == circle_ccw[:,::-1]).all() # it actually got reversed
    p = circle_cw + np.array([[280],[430]])

    plb.ion()
    plb.figure(0)
    plb.gray()
    i = plb.imread('mri2.png')
    i = np.mean(i, axis=2)
    plb.imshow(i)
    global _contour
    _contour, = plb.plot(np.append(p[0,-1], p[0]),np.append(p[1,-1], p[1]))
    plb.draw()
    Snake2D(i, p, iterations=500)
    print 'done'
    plb.ioff()
    plb.savefig('mri-result.png')
    plb.show()

def main():
    butils = measDeblend.BaselineUtilsF

    foot = buildExample2()

    fbb = foot.getBBox()
    mask1 = afwImg.MaskU(fbb.getWidth(), fbb.getHeight())
    mask1.setXY0(fbb.getMinX(), fbb.getMinY())
    afwDet.setMaskFromFootprint(mask1, foot, 1)

    if plt:
        plt.clf()
        plt.imshow(mask1.getArray(), origin='lower', interpolation='nearest',
                   extent=(fbb.getMinX(), fbb.getMaxX(), fbb.getMinY(), fbb.getMaxY()))
        plt.gray()
        plt.savefig('foot2.png')

        sfoot = butils.symmetrizeFootprint(foot, 355, 227)

        mask2 = afwImg.MaskU(fbb.getWidth(), fbb.getHeight())
        mask2.setXY0(fbb.getMinX(), fbb.getMinY())
        afwDet.setMaskFromFootprint(mask2, sfoot, 1)

        plt.clf()
        plt.imshow(mask2.getArray(), origin='lower', interpolation='nearest',
                   extent=(fbb.getMinX(), fbb.getMaxX(), fbb.getMinY(), fbb.getMaxY()))
        plt.gray()
        plt.savefig('sfoot3.png')

def test():

    from PIL import Image
    import numpy
    import pylab

    import os
    p = r'C:\repository\research_code\gwb_cropped\gwb_cropped'

    imlist = [p+'/'+f for f in os.listdir(p) if 'jpg' in f]

    im = numpy.array(Image.open(imlist[0])) #open one image to get the size
    m,n = im.shape[0:2] #get the size of the images
    imnbr = len(imlist) #get the number of images

#create matrix to store all flattened images
    immatrix = numpy.array([numpy.array(Image.open(imlist[i])).flatten() for i in range(imnbr)],'f')

#perform PCA
    V,S,immean = pca(immatrix)

#mean image and first mode of variation
    immean = immean.reshape(m,n)
    mode = V[0].reshape(m,n)

#show the images
    pylab.figure()
    pylab.gray()
    pylab.imshow(immean)

    pylab.figure()
    pylab.gray()
    pylab.imshow(mode)

    pylab.show()

    def plot(self):
        """
        
        .. plot::
            :include-source:
            :width: 80%
            
            from cellnopt.simulate import *
            from cellnopt.core import *
            pkn = cnodata("PKN-ToyPB.sif")
            midas = cnodata("MD-ToyPB.csv")
            s = boolean.BooleanSimulator(CNOGraph(pkn, midas))
            s.simulate(30)
            s.plot()
        """
        
        pylab.clf()

        data = numpy.array([self.data[x] for x in self.species if x in self.species])
        data = data.transpose()
        data = 1 - pylab.flipud(data)

        pylab.pcolor(data, vmin=0, vmax=1, edgecolors="k")
        pylab.xlabel("species"); 
        pylab.ylabel("Time (tick)");
        pylab.gray()

        pylab.xticks([0.5+x for x in range(0,30)], self.species, rotation=90)

        pylab.ylim([0, self.tick])
        pylab.xlim([0, len(self.species)])

def plotAVTable(experiment):
    pylab.figure()
    pylab.gray()
    pylab.pcolor(experiment.agent.module.params.reshape(81,4).max(1).reshape(9,9),
            shading='faceted')
    pylab.title("Action-Value table, %s, Run %d" %
            (experiment.agent.learner.__class__.__name__, experiment.stepid))

def main():
    A = pl.imread(IMAGE_FILE)

    i = 1
    pc_values = (1, 5, 10, 20, 30, 40)
    for num_pcs in pc_values:

        # perform (truncated) pca
        egvecs, proj, egvals = pca(A, num_pcs)

        # reconstruct image
        A_rec = np.dot(egvecs, proj).T + np.mean(A, axis=0)
         
        # create sublplot
        ax = pl.subplot(2, 3, i, frame_on=False)
        ax.xaxis.set_major_locator(pl.NullLocator())
        ax.yaxis.set_major_locator(pl.NullLocator())

        # draw
        pl.imshow(A_rec)
        pl.title("{} pc's".format(num_pcs))
        pl.gray()

        i += 1

    pl.show()

 def plot_matches(self, name = "harris_match.jpg", show_below = True, match_maximum = None):
     """ 対応点を線で結んで画像を表示する
     入力： 
     show_below（対応の下に画像を表示するならTrue）"""
     if self._append_image is None:
         self.appendimages()
     im1   = self._image_1.get_array_image()
     im2   = self._image_2.get_array_image()
     im3   = self._append_image
     if self._image_1.get_harris_point() is None:
         self._image_1.make_harris_points()
     if self._image_2.get_harris_point() is None:
         self._image_2.make_harris_points()
     locs1 = self._image_1.get_harris_point()
     locs2 = self._image_2.get_harris_point()
     if show_below:
         im3 = numpy.vstack((im3,im3))
     pylab.figure(dpi=160)
     pylab.gray()
     pylab.imshow(im3, aspect = "auto")
     
     cols1 = im1.shape[1]
     if self._match_score is None:
         self.match()
     if match_maximum is not None:
         self._match_score = self._match_score[:match_maximum]
     for i,m in enumerate(self._match_score):
         if m>0: 
             pylab.plot([locs1[i][1],locs2[m][1]+cols1],[locs1[i][0],locs2[m][0]],'c')
     pylab.axis('off')
     pylab.savefig(name, dpi=160)

	def interference_maker(self):
		wavelength = input("input wavelength: ")
		k = 2*pi/wavelength
		amplitude = input("input amplitude: ")
		distance_between_centers = input("input distance between centers: ")
		length = 100.0
		num_points = 500
		spacing = length/num_points

		first_x = length/2 + distance_between_centers/2
		second_x = length/2 - distance_between_centers/2

		first_y = length/2
		second_y = first_y

		points_array = empty([num_points, num_points], float)

		for i in range(num_points):
			y = spacing*i
			for j in range(num_points):
				x = spacing*j
				r1 = sqrt((x-first_x)**2 + (y - first_y)**2)
				r2 = sqrt((x-second_x)**2 + (y - second_y)**2)
				points_array[i, j] = amplitude*(sin(k*r1)+sin(k*r2))

		imshow(points_array, origin="lower", extent=[0, length, 0, length])
		gray()
		show()