%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import os

##  load (and compile if necessary) the ridge_directed_ring_detector
import pyximport; pyximport.install()
import ridge_directed_ring_detector as ring_detector

#ring_detector.__file__

## a function to plot the rings
pi = np.pi

colourmap = plt.jet
txt_colour = 'r'
ring_colour = 'y'

def plot_rings(xyr, minrad=1, ring_colour = 'g', ring_style=[None,'dotted','dashed'][0], lw=1, txt_colour='b', title=None):
    '''
    given a 3-column array of (x,y) for ring centres & a radius, plots the
    rings
    '''
    #import matplotlib.pyplot as plt
    for c in xyr:
        if c[-1]>minrad:
            circ = plt.Circle(c[:2][::-1],c[-1],ec=ring_colour,fill=False, lw=lw, ls=ring_style) 
            plt.gca().add_patch(circ)
            i,j,r = c
            if txt_colour:
                plt.text(j,i,'.', fontsize=20, color=txt_colour)
    plt.draw()

ridge_hough = ring_detector.RidgeHoughTransform()

default_params = ridge_hough.params.copy()

ridge_hough.params['sigma'] = 0.8
ridge_hough.params['curv_thresh'] = -50
ridge_hough.params['circle_thresh'] = 2 * np.pi * .33
ridge_hough.params['Rmin'] = 4
ridge_hough.params['Rmax'] = 50
ridge_hough.params['vote_thresh'] = 4
ridge_hough.params['dr'] = 3
ridge_hough.params['eccentricity'] = 0
experimental_params = ridge_hough.params.copy()

# print 'default img analysis parameters:'
# for key in sorted(default_params.keys()):
#     print '\t %s = %s'  % (key, default_params[key])
    
print '\nimg analysis parameters in the experiment:'
for key in sorted(experimental_params.keys()):
    print '\t %s = %s'  % (key, experimental_params[key])

## Load picture
imgName = '../data/unprocessed_fig_46.png'
img = plt.imread(imgName )

if img.max()<=1.:
    img*=255./img.max()
    
IMG = img.astype(np.float32)

%%time

ridge_hough = ring_detector.RidgeHoughTransform(IMG)
ridge_hough.img_preprocess()
#ridge_hough.debugging_rings_detection()
ridge_hough.rings_detection()
#detected_rings = ridge_hough.output['rings']
rings = ridge_hough.output['rings_subpxl'] # keep for later (sub-pxl mask example)

print 'NB: first run may be several times slower and not indicative.\n'

fig,ax = plt.subplots(figsize=(15,10))
im = ax.matshow(IMG, cmap=plt.cm.bone_r, vmax=80)
# fig.colorbar(im)
plot_rings(rings, ring_colour='r')
plt.grid(True)
print 'found {} rings'.format(rings.shape[0])

if True:
    import pandas as pd
    #pd.set_option('display.max_rows', 10)
    rings_df = pd.DataFrame(rings[rings.T[-1].argsort()[::-1]], 
                            columns=['row','col','rad'], 
                            index=range(1,rings.shape[0]+1),
                            )
rings_df

%%timeit
ridge_hough = ring_detector.RidgeHoughTransform(IMG)
ridge_hough.img_preprocess()
# ridge_hough.debugging_rings_detection()
# ridge_hough.rings_detection()
# # detected_rings = ridge_hough.output['rings']
# rings = ridge_hough.output['rings_subpxl'] # keep for later (sub-pxl mask example)

from skimage import data
from scipy import ndimage, misc

image = data.coins()[0:95, 70:370]

%%time

## as the algorithm is currently meant for ridges, convert the image to a gradient magnitude image
IMG = image.astype(np.float32)
sigma = .5
d_mag = ndimage.gaussian_gradient_magnitude(IMG, sigma)

IMG = d_mag

ridge_hough = ring_detector.RidgeHoughTransform(IMG)
ridge_hough.params['sigma'] = 3
ridge_hough.params['Rmin'] = 15
ridge_hough.params['Rmax'] = 30
ridge_hough.params['curv_thresh'] = -35
ridge_hough.params['circle_thresh'] = .4*2*pi
ridge_hough.params['vote_thresh'] = 3
ridge_hough.params['dr'] = 3
ridge_hough.img_preprocess()
#ridge_hough.debugging_rings_detection()
ridge_hough.rings_detection()
#detected_rings = ridge_hough.output['rings']
CHT_rings = ridge_hough.output['rings_subpxl'] # keep for later (sub-pxl mask example)

print 'NB: first run may be several times slower and not indicative.\n'

plt.matshow(image, cmap=plt.cm.gray)
plot_rings(CHT_rings, ring_colour='r')
plt.grid(True)
print CHT_rings

if False:
    image = data.coins()[0:95, 70:370].astype(np.float32)
    sigma = .5
    di,dj = [ndimage.gaussian_filter(image,sigma,order=der) for der in ((0,1),(1,0))]
    d_mag = np.sqrt(di**2+dj**2)
    
plt.matshow(d_mag)
plt.colorbar()
print 'gradient magnitude: \n\tmean = {} \n\tstd = {}'.format(d_mag.mean(), d_mag.std())

plt.matshow(np.asarray(ridge_hough.deriv['principal_curv']))
plt.colorbar()
#A = [plot(circ[1],circ[0],'o') for circ in CHT_rings]
plot_rings(CHT_rings)
plt.grid(True)
#print CHT_rings