How to use the mahotas.features function in mahotas

def colors(image):
    image = image // 64
    r,g,b = image.transpose((2,0,1))
    pixels = 1 * r + 4 * b + 16 * g
    hist = np.bincount(pixels.ravel(), minlength=64)
    hist = hist.astype(float)
    hist = np.log1p(hist)
    return hist

for i in range(len(images)):
    print "processing image %i of %i" % (i+1, len(images)) 
    labels.append(images[i][:-len('00.jpg')])
    im = mh.imread(images[i])
    imgrey = mh.colors.rgb2gray(im, dtype=np.uint8)
    features.append(np.concatenate([mh.features.haralick(im).ravel(), mh.features.lbp(imgrey, 30, 10).ravel(), colors(im)]))
    surfim = mh.imread(images[i], as_grey=True)
    surfim = surfim.astype(np.uint8)
    alldescriptors.append(surf.dense(surfim, spacing=16))

concatenated = np.concatenate(alldescriptors)
print "fitting k mean clusters for surf descriptors"
km = KMeans(15)
km.fit(concatenated)
print "creating surf features"
sfeatures = []
for d in alldescriptors:
    c = km.predict(d)
    sfeatures.append(np.array([np.sum(c == ci) for ci in range(15)]))

features = np.array(features) 
sfeatures = np.array(sfeatures, dtype=float)

def fd_haralick(image):
    # convert the image to grayscale
    gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)
    # compute the haralick texture feature vector
    haralick = mahotas.features.haralick(gray).mean(axis=0)
    # return the result
    return haralick

def describe(self, image):
		# return the Zernike moments for the image
		return mahotas.features.zernike_moments(image, self.radius)

def calc_haralick(roi):

    feature_vec = []

    texture_features = mt.features.haralick(roi)
    mean_ht = texture_features.mean(axis=0)

    [feature_vec.append(i) for i in mean_ht[0:9]]

    return np.array(feature_vec)

mask[~m2] = False
                else:
                    mask = m1

        pixel_data[~mask] = 0
        # mahotas.features.haralick bricks itself when provided a dtype larger than uint8 (version 1.4.3)
        pixel_data = skimage.util.img_as_ubyte(pixel_data)

        features = numpy.empty((n_directions, 13, len(unique_labels)))

        for index, label in enumerate(unique_labels):
            label_data = numpy.zeros_like(pixel_data)
            label_data[labels == label] = pixel_data[labels == label]

            try:
                features[:, :, index] = mahotas.features.haralick(
                    label_data, distance=scale, ignore_zeros=True
                )
            except ValueError:
                features[:, :, index] = numpy.nan

        for direction, direction_features in enumerate(features):
            for feature_name, feature in zip(F_HARALICK, direction_features):
                statistics += self.record_measurement(
                    image=image_name,
                    feature=feature_name,
                    obj=object_name,
                    result=feature,
                    scale="{:d}_{:02d}".format(scale, direction),
                    workspace=workspace,
                )

def describe(self, image):
		# return the Zernike moments for the image
		return mahotas.features.zernike_moments(image, self.radius)

hlabels = mh.features.texture.haralick_labels[:-1]
        colnames.extend(["mean:{}".format(ell) for ell in hlabels])
        colnames.extend(["ptp:{}".format(ell) for ell in hlabels])
    if args.lbp:
        from mahotas.features.lbp import lbp_names
        colnames.extend(lbp_names(args.lbp_radius, args.lbp_points))
    _write_row(output, colnames)

    for fname in args.fnames:
        cur = []
        im = read_bw(fname, args)
        if args.haralick:
            har = mh.features.haralick(im, return_mean_ptp=True)
            cur.append(har)
        if args.lbp:
            cur.append(mh.features.lbp(im, args.lbp_radius, args.lbp_points))

        _write_row(output, chain.from_iterable(cur), fname)
    output.close()

(cnts, _) = cv2.findContours(thresh.copy(), cv2.RETR_EXTERNAL,
		cv2.CHAIN_APPROX_SIMPLE)

	# loop over the contours
	for c in cnts:
		# create an empty mask for the contour and draw it
		mask = np.zeros(image.shape[:2], dtype="uint8")
		cv2.drawContours(mask, [c], -1, 255, -1)

		# extract the bounding box ROI from the mask
		(x, y, w, h) = cv2.boundingRect(c)
		roi = mask[y:y + h, x:x + w]

		# compute Zernike Moments for the ROI and update the list
		# of shape features
		features = mahotas.features.zernike_moments(roi, cv2.minEnclosingCircle(c)[1], degree=8)
		shapeFeatures.append(features)

	# return a tuple of the contours and shapes
	return (cnts, shapeFeatures)

output = open(args.output, 'w')
    colnames = []
    if args.haralick:
        hlabels = mh.features.texture.haralick_labels[:-1]
        colnames.extend(["mean:{}".format(ell) for ell in hlabels])
        colnames.extend(["ptp:{}".format(ell) for ell in hlabels])
    if args.lbp:
        from mahotas.features.lbp import lbp_names
        colnames.extend(lbp_names(args.lbp_radius, args.lbp_points))
    _write_row(output, colnames)

    for fname in args.fnames:
        cur = []
        im = read_bw(fname, args)
        if args.haralick:
            har = mh.features.haralick(im, return_mean_ptp=True)
            cur.append(har)
        if args.lbp:
            cur.append(mh.features.lbp(im, args.lbp_radius, args.lbp_points))

        _write_row(output, chain.from_iterable(cur), fname)
    output.close()