Python scipy.histogram() Examples

def __findKnuth(x, minv, maxv):
    """
	Implement Knuth method for histogram bin selection
		
	"""
    N = ((x >= minv) & (x <= maxv)).sum()

    def funcer(M):
        hh, loc = scipy.histogram(x, bins=M, range=[minv, maxv])
        return np.log(M) + 1. / N * (scipy.special.gammaln(M / 2.) -
                                     M * scipy.special.gammaln(0.5) -
                                     scipy.special.gammaln(N + M / 2.) +
                                     scipy.special.gammaln(hh + 0.5).sum())

    maxN = 1000
    ns = np.arange(1, maxN + 1)
    vals = ns * 0.
    for i in range(len(ns)):
        vals[i] = funcer(ns[i])
    bestn = ns[np.argmax(vals)]
    if bestn == maxN:
        print('WARNING the best number of bins is > maxbin(%d)' % (maxn))
    return bestn 

def _auto_color(self, url:str, ranks):
        phrases = ["Calculating colors..."] # in case I want more
        #try:
        await self.bot.say("**{}**".format(random.choice(phrases)))
        clusters = 10

        async with aiohttp.get(url) as r:
            image = await r.content.read()
        with open('data/leveler/temp_auto.png','wb') as f:
            f.write(image)

        im = Image.open('data/leveler/temp_auto.png').convert('RGBA')
        im = im.resize((290, 290)) # resized to reduce time
        ar = scipy.misc.fromimage(im)
        shape = ar.shape
        ar = ar.reshape(scipy.product(shape[:2]), shape[2])

        codes, dist = scipy.cluster.vq.kmeans(ar.astype(float), clusters)
        vecs, dist = scipy.cluster.vq.vq(ar, codes)         # assign codes
        counts, bins = scipy.histogram(vecs, len(codes))    # count occurrences

        # sort counts
        freq_index = []
        index = 0
        for count in counts:
            freq_index.append((index, count))
            index += 1
        sorted_list = sorted(freq_index, key=operator.itemgetter(1), reverse=True)

        colors = []
        for rank in ranks:
            color_index = min(rank, len(codes))
            peak = codes[sorted_list[color_index][0]] # gets the original index
            peak = peak.astype(int)

            colors.append(''.join(format(c, '02x') for c in peak))
        return colors # returns array
        #except:
            #await self.bot.say("```Error or no scipy. Install scipy doing 'pip3 install numpy' and 'pip3 install scipy' or read here: https://github.com/AznStevy/Maybe-Useful-Cogs/blob/master/README.md```")

    # converts hex to rgb 

def window_func(x, y, func, xmin=None, xmax=None, nbin=100, empty=False,
			xlog=False):
	"""This function does compute a user-supplied function 
	on the subsets of y grouped by values of x
	E.g. imagine that you have x from 0 to 100 and you want
	to know the mean values of y for 0<x<10, 10<x<2- etc..
	In that case you need to do
	xbin,funcy,nperbin=window_func(x,y,lambda x: x.mean(),0,100,nbin=10)
	where xbin is the array with the centers of bins, 
	funcy is the func -- evaluated for y where x is within the appropriate bin
	and nperbin is the number of points in the appropriate bin
	empty keyword is needed if you want to retrieve the function evaluation 
	in empty bins too, otherwise they aren't returned at all
	"""
	
	if xmin is None:
		xmin = x.min()
	if xmax is None:
		xmax = x.max()
	if xlog:
		xmin,xmax,x=[numpy.log10(tmp) for tmp in [xmin,xmax,x]]
	if (len(x)!=len(y) and x.ndim==1) or (x.shape!=y.shape and x.ndim>1):
		raise ValueError('Input arrays must have the same size')
	#hh,loc=scipy.histogram(x,range=(xmin,xmax),bins=nbin)
	inds = ((x-xmin)/float(xmax-xmin)*nbin).astype(int)
	mask = numpy.zeros(nbin, bool)
	retv = numpy.zeros(nbin)
	hh = numpy.zeros(nbin,int)

	for i in range(nbin):
		cury=y[inds==i]
		hh[i]=len(cury)
		mask[i]=len(cury)>0
		if len(cury)>0 or empty:
			retv[i]=func(cury)
	retx = xmin+(xmax-xmin)*1./nbin*(0.5+numpy.arange(nbin))
	if xlog:
		retx = 10**retx
	if empty:
		mask |= True
	return retx[mask], retv[mask],hh[mask] 

def discrete_kldiv(X_baseline: np.array, X_sampled: np.array) -> float:
    P, bins = histogram(X_baseline, bins=10, density=True)
    P += 1e-10
    Q, _ = histogram(X_sampled, bins=bins, density=True)
    Q += 1e-10

    return entropy(P, Q) 

def  extract_dsift(cls, rawfilename, winsize, nbins):
    """read_ppm(rawfilename, filename)

    Read in raw pixel data from rawfilename (.ppm).
    Create a histogram around each pixel to become
    the feature vector for that obsevation (pixel).
    Pickle the result and save it to filename.
    Note: does NOT update object fields.
    Follow this with a call to readin().
    """
    im  = Image.open(rawfilename)
    (width, height) = im.size

    # To be removed in the future
    # Pick up all windows, stepping by half of the window size
    labels  = []
    halfwin = int(winsize/2)
    for y in range(halfwin, height-halfwin, int(halfwin/2)):
      for x in range(halfwin, width-halfwin, int(halfwin/2)):
        labels    += ['(%d,%d)' % (y,x)]

    mlab.bb_dsift(N.array(im), winsize, 'temp.mat')
    sift_features = scipy.io.loadmat('temp.mat')
    sift_features = sift_features['d_']

    return (sift_features, labels, width, height) 

def extract_hist_subimg(sub_image):
    hist_bins = range(0,260,1)
    hist_features = N.histogram(sub_image.ravel(), hist_bins)[0]
    return hist_features 

def get_dominant_image_colors(image, num_clusters=4):
    """
    Returns the dominant image color that isn't pure white or black.  Uses
    kmeans on the colors.  Returns the result as RGB hex strings in the format
    ['#rrggbb', '#rrggbb', ...].

    :param image: PIL image or path
    """

    if isinstance(image, basestring):
        image = Image.open(image)

    # downsample for speed
    im = image.resize((512, 512), Image.ANTIALIAS)

    # reshape
    ar0 = scipy.misc.fromimage(im)
    shape = ar0.shape
    npixels = scipy.product(shape[:2])
    ar0 = ar0.reshape(npixels, shape[2])

    # keep only nontransparent elements
    ar = ar0[ar0[:, 3] == 255][:, 0:3]

    try:
        # kmeans clustering
        codes, dist = scipy.cluster.vq.kmeans(ar, num_clusters)
    except:
        # kmeans sometimes fails -- if that is the case, use the mean color and
        # nothing else.
        arf = ar.astype(float)
        clamp = lambda p: max(0, min(255, int(p)))
        return ['#' + ''.join(['%0.2x' % clamp(arf[:, i].sum() / float(arf.shape[1])) for i in (0, 1, 2)])]

    vecs, dist = scipy.cluster.vq.vq(ar, codes)         # assign codes
    counts, bins = scipy.histogram(vecs, len(codes))    # count occurrences

    # sort by count frequency
    indices = [i[0] for i in
               sorted(enumerate(counts), key=lambda x:x[1], reverse=True)]

    # convert to hex strings
    colors = [''.join(chr(c) for c in code).encode('hex') for code in codes]

    results = []
    for idx in indices:
        color = colors[idx]
        if color != 'ffffff' and color != '000000':
            results.append('#' + color)

    return results 

def  extract_sift(cls, rawfilename, winsize, nbins):
    """read_ppm(rawfilename, filename)

    Read in raw pixel data from rawfilename (.ppm).
    Create a histogram around each pixel to become
    the feature vector for that obsevation (pixel).
    Pickle the result and save it to filename.
    Note: does NOT update object fields.
    Follow this with a call to readin().
    """
    if cls._VL_SIFT_:
      # VLSIFT matlab 

      im  = Image.open(rawfilename)
      (width, height) = im.size

      mlab.bb_sift(N.array(im), 'temp.mat')
      sift_features = scipy.io.loadmat('temp.mat')
      kp = sift_features['f_']
      sift_features = sift_features['d_']
      sift_features  = scipy.concatenate((sift_features.transpose(), kp[2:4].transpose()), 1).transpose()

      labels = [];
      for ikp in kp.transpose():
        (x,y) = ikp[0:2]
        labels    += ['(%d,%d)' % (y,x)]
    else:
      #Opencv SIFT 
      img = cv2.imread(rawfilename)
      gray= cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
      height, width = gray.shape

      # Computing SIFT
      sift = cv2.SIFT(edgeThreshold = 3)
      kp, des = sift.detectAndCompute(gray,None)

      labels  = []
      sift_features = N.transpose(des)
      scale_angle = []

      for ikp in kp:
        (x,y) = ikp.pt
        scale_angle.append([ikp.size/12, ikp.angle])
        labels    += ['(%d,%d)' % (y,x)]
    
      scale_angle = N.array(scale_angle)
      sift_features  = scipy.concatenate((sift_features.transpose(), scale_angle), 1).transpose()

    return (sift_features, labels, width, height) 

def  extract_hist(cls, rawfilename, winsize, nbins):
    # This function extracts the histogram features from the image

    im  = Image.open(rawfilename)
    
    (width, height) = im.size
    npixels = width * height
    pix = scipy.array(im)

    # Generate one feature vector (histogram) per pixel
    #winsize = 20  # for test.pgm
    #winsize = 0  # for RGB
    halfwin = int(winsize/2)

    bins    = scipy.linspace(0, 255, nbins)

    # Only use windows that are fully populated
    mywidth  = width-winsize
    myheight = height-winsize
    #data     = scipy.zeros((nbins-1, mywidth * myheight))
    #data     = scipy.zeros((3*winsize*winsize, mywidth * myheight))
    data    = []
    labels  = []

    # Pick up all windows, stepping by half of the window size
    for y in range(halfwin, height-halfwin, int(halfwin/2)):
      for x in range(halfwin, width-halfwin, int(halfwin/2)):
        # Read in data in row-major order
        ind = (y-halfwin)*mywidth + (x-halfwin)
        #data[:,ind] = \
        #    scipy.histogram(pix[y-halfwin:y+halfwin,
        #                        x-halfwin:x+halfwin],
        #                        bins)[0]
        # Just RGB
        #data[:,ind] = pix[y,x]
        # RGB window
        #data[:,ind] = pix[y-halfwin:y+halfwin,x-halfwin:x+halfwin].flat
        hist_features = TCData.extract_hist_subimg(pix[y-halfwin:y+halfwin,x-halfwin:x+halfwin])
        if data == []:
          data = hist_features.reshape(-1,1)
        else:
          data = scipy.concatenate((data, hist_features.reshape(-1,1)),1)
        labels    += ['(%d,%d)' % (y,x)]

    return (data, labels, width, height) 

def read_ppm(self, rawfilename, filename):
    # This function reads the ppm/jpg file and extracts the features if the 
    # features pkl file doesn't exist. It is also compatible for extension 
    # of the feauture vector and doesn't compute the already computed features

    new_feature_string = []
    updated_feature = 0
    data = N.array([], dtype=int)
    if os.path.exists(filename):
      pkl_f = open(filename, 'r')
      (data, labels, feature_string, width, height, winsize, nbins)= pickle.load(pkl_f)
      self.winsize = winsize
      self.nbins = nbins
      new_feature_string = list(feature_string)
      pkl_f.close()      

    if not new_feature_string.count('dsift'):
      updated_feature = 1
      (sift_features, labels, width, height) = self.extract_dsift(rawfilename, self.winsize, self.nbins)
      if data.size:
        data = scipy.concatenate((data.transpose(), sift_features.transpose()), 1).transpose()
      else:
        data = sift_features
      new_feature_string.append('dsift')

    if not new_feature_string.count('histogram'):
      updated_feature = 1 
      (hist_features, labels, width, height) = self.extract_hist(rawfilename, self.winsize, self.nbins)
      hist_features = hist_features/(self.winsize)
      if data.size:
        data = scipy.concatenate((data.transpose(), hist_features.transpose()), 1).transpose()
      else:
        data = hist_features
      new_feature_string.append('histogram')

    '''
    if not new_feature_string.count('position'):
      updated_feature = 1 
      
      position_features = []
      for label in labels:
        (y,x) = map(int, label.strip('()').split(','))
        position_features.append([x,y]) 
      position_features = N.array(position_features)
    
      if data.size:
        data = scipy.concatenate((data.transpose(), position_features), 1).transpose()
      else:
        data = position_features
      new_feature_string.append('position')
    '''
    if updated_feature:
      outf = open(filename, 'w')
      pickle.dump((data, labels, new_feature_string, width, height, self.winsize, self.nbins),outf)
      outf.close()
      print 'Saved data to %s.' % filename
    
    return (data, labels, new_feature_string, width, height, self.winsize, self.nbins)