import numpy as np
import matplotlib.pyplot as plt
from scipy import misc as mc

img = mc.imread('kolka.png')

print img.dtype
print img.shape

R = img[:,:,0]
G = img[:,:,1]
B = img[:,:,2]

img=R*0.2126+G*0.7152+B*0.0722

I1=img



def ffilter(img,mask):
    m=mask.shape[0]
    n=mask.shape[1]
    if np.sum(mask)==0:
        N=1
    else:
        N=np.sum(mask)
    mm=m/2
    if m%2==1:
        mm=mm+1
    nm=n/2
    if n%2==1:
        nm=nm+1
    offsetup=mm-1
    offsetdown=mask.shape[0]-offsetup
    offsetleft=nm-1
    offsetright=mask.shape[1]-offsetleft
    img=powiekszenie(img,mask)
    x=img.shape[0]
    y=img.shape[1]
    for i in range (offsetup,x-offsetdown):
        for j in range (offsetleft,y-offsetright):
            img[i,j]=np.sum(np.multiply(mask,img[i-offsetup+1:i+offsetdown+1,j-offsetleft+1:j+offsetright+1]))
            img[i,j]=img[i,j]/N
    return img


def powiekszenie(img,mask):
    m=mask.shape[0]
    n=mask.shape[1]
    mm=m/2
    if m%2==1:
        mm=mm+1
    nm=n/2
    if n%2==1:
        nm=nm+1
    offsetup=mm-1
    offsetdown=mask.shape[0]-mm
    if m%2==0:
        offsetdown=offsetdown+1
    offsetleft=nm-1
    offsetright=mask.shape[1]-nm
    if n%2==0:
        offsetright=offsetright+1
    A=np.ones((offsetup,offsetleft))*img[0,0]
    B=np.ones((offsetup,img.shape[1]))
    B=np.multiply(B,img[0,:])
    C=np.ones((offsetup,offsetright))*img[0,img.shape[0]-1]
    D=np.ones((img.shape[0],offsetleft))
    D=np.multiply(D,img[:,0:1:1])
    E=np.ones((img.shape[0],offsetright))
    E=np.multiply(E,img[:,img.shape[1]-1:img.shape[1]:1])
    F=np.ones((offsetdown,offsetleft))*img[img.shape[1]-1,0]
    G=np.ones((offsetdown,img.shape[1]))
    G=np.multiply(G,img[img.shape[0]-1,:])
    H=np.ones((offsetdown,offsetright))*img[img.shape[0]-1,img.shape[1]-1]
    img=np.concatenate((B,img),axis=0)
    CE=np.concatenate((C,E),axis=0)
    img=np.concatenate((img,CE),axis=1)
    AD=np.concatenate((A,D),axis=0)
    img=np.concatenate((AD,img),axis=1)
    FGH=np.concatenate((F,G,H),axis=1)
    img=np.concatenate((img,FGH),axis=0)
    return img

def mfilter(img,mask):
    i=mask.shape[0]
    m=mask.shape[0]
    n=mask.shape[1]
    mm=m/2
    if m%2==1:
        mm=mm+1
    nm=n/2
    if n%2==1:
        nm=nm+1
    offsetup=mm-1
    offsetdown=mask.shape[1]-mm
    offsetleft=nm-1
    offsetright=mask.shape[2]-nm
    x=img.shape[0]+offsetup+offsetdown
    y=img.shape[1]+offsetleft+offsetright
    images=np.zeros((x,y))
    for j in range (0,i):
        images=images+abs(ffilter(img,mask[j]))
    return (images/i)


sobel=np.zeros((4,3,3))
sobel[0,:,:]=np.array([[-1,0,1],[-2,0,2],[-1,0,1]])
sobel[1,:,:]=np.array([[-1,-2,-1],[0,0,0],[1,2,1]])
sobel[2,:,:]=np.array([[0,1,2],[-1,0,1],[-2,-1,0]])
sobel[3,:,:]=np.array([[-2,-1,0],[-1,0,1],[0,1,2]])

def kolo(r):
    mask=np.zeros((2*r,2*r))
    m=mask.shape[0]
    n=mask.shape[1]

    mm=r-1
    nm=r-1
    for k in range (1,361,1):
        x=np.round(mm+r*np.sin(k)).astype('int')
        y=np.round(nm+r*np.cos(k)).astype('int')
        mask[x,y]=1

    return mask

def kolofilter(img,mask):
    m=mask.shape[0]
    n=mask.shape[1]
    if np.sum(mask)==0:
        N=1
    else:
        N=np.sum(mask)
    mm=m/2
    if m%2==1:
        mm=mm+1
    nm=n/2
    if n%2==1:
        nm=nm+1
    offsetup=mm-1
    offsetdown=mask.shape[0]-offsetup
    offsetleft=nm-1
    offsetright=mask.shape[1]-offsetleft
    img=kolopowiekszenie(img,mask)
    x=img.shape[0]
    y=img.shape[1]
    for i in range (offsetup,x-offsetdown):
        for j in range (offsetleft,y-offsetright):
            img[i,j]=np.sum(np.multiply(mask,img[i-offsetup+1:i+offsetdown+1,j-offsetleft+1:j+offsetright+1]))
            img[i,j]=img[i,j]/N
    return img


def kolopowiekszenie(img,mask):
    m=mask.shape[0]
    n=mask.shape[1]
    mm=m/2
    if m%2==1:
        mm=mm+1
    nm=n/2
    if n%2==1:
        nm=nm+1
    offsetup=mm-1
    offsetdown=mask.shape[0]-mm
    if m%2==0:
        offsetdown=offsetdown+1
    offsetleft=nm-1
    offsetright=mask.shape[1]-nm
    if n%2==0:
        offsetright=offsetright+1
    A=np.zeros((offsetup,offsetleft))
    B=np.zeros((offsetup,img.shape[1]))
    C=np.zeros((offsetup,offsetright))
    D=np.zeros((img.shape[0],offsetleft))
    E=np.zeros((img.shape[0],offsetright))
    F=np.zeros((offsetdown,offsetleft))
    G=np.zeros((offsetdown,img.shape[1]))
    H=np.zeros((offsetdown,offsetright))
    img=np.concatenate((B,img),axis=0)
    CE=np.concatenate((C,E),axis=0)
    img=np.concatenate((img,CE),axis=1)
    AD=np.concatenate((A,D),axis=0)
    img=np.concatenate((AD,img),axis=1)
    FGH=np.concatenate((F,G,H),axis=1)
    img=np.concatenate((img,FGH),axis=0)
    return img

image=mfilter(I1,sobel)
poProgowaniu = (image[:,:] >= 25).astype('float64')


def ht(img, rmin, rmax, N):
    m=img.shape[0]
    n=img.shape[1]
    acc=np.zeros((m,n))
    licznik=0
    possible=np.zeros((100,3))
    middle=np.zeros((N,3))
    tempak=np.zeros((m,n))
    for ir in range (rmin,rmax,2):
        temp=kolofilter(img,kolo(ir))

        w=np.where(temp*255>=140)
        print w
        print len(w[0][:])
        if len(w[0][:])!=0:
            for x in range (0,len(w[0][:])):
                possible[licznik,0]=w[0][x]
                possible[licznik,1]=w[1][x]
                possible[licznik,2]=temp[w[0][x],w[1][x]]*255
                licznik=licznik+1
                tempak[w[0][x],w[1][x]]=1

        print ir
    plt.imshow(tempak, cmap=plt.cm.gray)
    plt.show()
    #acc=acc/((rmax-rmin)/2)
    for o in range (0,licznik):
        for s in range (0,N):
            if (possible[o,2]>middle[s,2]):
                temp=middle[s,0]
                middle[s,0]=possible[o,0]
                possible[o,0]=temp
                temp=middle[s,1]
                middle[s,1]=possible[o,1]
                possible[o,1]=temp
                temp=middle[s,2]
                middle[s,2]=possible[o,2]
                possible[o,2]=temp

    for d in range (0,N):
        print middle[d,0]
        print middle[d,1]
        print middle[d,2]
        acc[middle[d,0].astype('int'),middle[d,1].astype('int')]=1

    return acc

plt.subplot(1,2,1)
plt.imshow(poProgowaniu, cmap=plt.cm.gray)
plt.subplot(1,2,2)
temp=ht(poProgowaniu,14,58,5)
plt.imshow(temp, cmap=plt.cm.gray)
plt.show()