# PROGRAMA para fazer filtragem linear
##################################################################
import numpy as np
import matplotlib.pyplot as plt

def mat2vec(M): 
    tl,tc = M.shape
    M1=np.zeros(( tl*tc),dtype = float)
    c=0
    for i in range(tl):      # varrer template e criar vetor
        for j in range(tc):
            M1[c]=M[i,j]
            c=c+1
#    print(M1)
    return(M1)

print('######################################')
print(' correlacao .......................')

# ler imagem de entrada e recuperar dimensoes
X1= plt.imread('recorte.tif')
nl,nc = X1.shape
X=np.array(X1, dtype=float)   #transforma em float para facilitar multiplicacoes
Y=np.zeros((nl,nc),dtype = float)        # replica imagem para gerar uma saida
Y=np.uint8(Y)

# Ler template
T1= plt.imread('template.tif')
tl,tc = T1.shape
T=np.array(T1, dtype=float)
print('Image size: ',nl,nc, "  Template size: ", tl, tc)


# ################# cria janela de filtro tamanho DIM  ################
dim=tl              # dimensao do filtro exe 3x3 dim=3, 5x5 dim=5...
lado=(dim-1)/2     # numero de vizinhos antes ou depois do central
lado=np.uint8(lado)
print('dim=',dim, 'vizinho=',lado)

TT=mat2vec(T)

## correlacao propriamente dita
for L in range(lado,nl-lado):          # varrer em linhas lado+1, para vizinhos+central
    for C in range (lado, nc-lado):    # varrer em colunas
        # posicao inicial
        p=(L-lado)             # determina o elemento do filtro a ser usado
        q=(C-lado)
        R=X[p:p+tl, q:q+tc]
        RR=mat2vec(R)
        rm=np.corrcoef(RR,TT)
        r=rm[0,1]
        if r<0:
            r=0
        rp=r*255
        v=np.uint8( np.round( rp ) )   # arredonda e muda a uint8
        Y[L,C]=v                      # salva na posicao do central

plt.imsave('Matriz.png',Y,cmap='gray')

## buscar maximo local que supere limiar
Limiar=250
maximo=0
for L in range(lado,nl-lado):          # varrer em linhas lado+1, para vizinhos+central
    for C in range (lado, nc-lado):    # varrer em colunas
        if Y[L,C]>Limiar:              #supera limiar, vejamos se é maximo local
	    # se é maximo local, é maior que os quatro vizinhos
            if Y[L,C]>Y[L+1,C] and Y[L,C]>Y[L-1,C] and Y[L,C]>Y[L,C+1] and Y[L,C]>Y[L,C-1]:
                print('max=', L,C, Y[L,C], maximo)
                if Y[L,C]>maximo:     # e se é maior que o máximo anteriormente regstrado
                    maxL=L            # coordenadas e valor do máximo
                    maxC=C
                    maximo=Y[L,C]

print('localizado em: ', maxL, maxC, '| ',maximo)

# visualizar imagem com o local marcado
imgplot = plt.imshow(X, cmap='gray', vmin=0, vmax=255)
plt.plot(maxC, maxL, 'o', color='red');
plt.show()