%matplotlib inline
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np
import sys,os
ia898path = os.path.abspath('/etc/jupyterhub/ia898_1s2017/')
if ia898path not in sys.path:
    sys.path.append(ia898path)
import ia898.src as ia

nb = ia.nbshow(2)
f = mpimg.imread('../data/cameraman.tif')
nb.nbshow(f,'Imagem original')
fsort = np.sort(f.ravel()).reshape(f.shape)
nb.nbshow(fsort, 'Imagem pixels ordenados')
nb.nbshow()

"""
Explanation: Equalização de Histograma
A transformação de contraste que procura distribuir a ocorrência dos níveis de cinza igualmente
na faixa de tons de cinza é denominada equalização de histograma.
O objetivo do exemplo a seguir é o de fazer a equalização de histograma de uma imagem. Ele é
dado pela aplicação de uma transformação de contraste T(r) que é calculado pelo histograma
acumulado normalizado da imagem.
Equação da Transformação de Intensidade para Equalizar Histograma
A equação da transformação de intensidade T(r) a partir do histograma h(i) da imagem que equaliza
a imagem, é dada por:
$$   T(r) = \frac{(L-1)}{n} \sum_{i = 0}^{r} h(i), \quad r = 0, 1,..., L-1 $$

Os valores de intensidade dos pixels da imagem (indices i e r) variam de 0 a L-1;
n é o número total de pixels na imagem;
$h(i)$ é o histograma da imagem, i. e., é o número de vezes que o nível de cinza i aparece na imagem;    
$T(r)$ é a função de mapeamento que equaliza a imagem dada pelo histograma $h$.
$L$ é o número de níveis na imagem final, usualmente 256 para imagens uint8.

Existem dois pontos que devem ser considerados no uso da equalização de histograma:

A formulação acima é exata no caso contínuo. Para o caso discreto, o histograma equalizado
   não fica totalmente planar pois não é possível no modelo acima distribuir um nível de cinza
   com muita ocorrência em mais de um valor de cinza. Normalmente mostra-se o histograma acumulado
   que é normalmente uma reta de 45 graus, como será visto no exemplo abaixo. Veremos num outro
   tutorial uma forma de contornar este problema.

Ilustração da equalização de histograma
A imagem original é mostrada a seguir, e para fins ilustrativos, a imagem com seus pixels ordenados também é mostrada ao
lado mostrando por exemplo que a imagem tem uma predominância de níveis de cinza escuros, bem 
destacados.
End of explanation
"""


h = ia.histogram(f)
plt.plot(h),plt.title('Histograma de f');

n = f.size
T = 255./n * np.cumsum(h)
T = T.astype('uint8')
plt.plot(T),plt.title('Transformação de intensidade para equalizar');

"""
Explanation: Plotamos seu histograma e calculamos a transformação de contraste que equaliza o histograma baseado
na equação vista anteriormente. O somatório da equação é eficientemente calculado com a função
np.cumsum que calcula a soma acumulada de um vetor. Visualizamos a transformação T[r] pelo gráfico:
End of explanation
"""


nb = ia.nbshow(3)
nb.nbshow(f,'imagem original, média=%d' % (f.mean()))
g = T[f]
nb.nbshow(g, 'imagem equalizada, média=%d' % (g.mean()))
gsort = np.sort(g.ravel()).reshape(g.shape)
nb.nbshow(gsort, 'imagem equalizada ordenada')
nb.nbshow()

"""
Explanation: A aplicação da transformação T em f: fazendo-se g = T[f], resulta na imagem g equalizada.
Para fins ilustrativos, colocamos ao lado a imagem equalizada com seus pixels ordenados.
Observa-se que a distribuição dos níveis de cinza ficou uniforme:
End of explanation
"""


plt.figure(0)
hg = ia.histogram(g)
plt.plot(hg),plt.title('Histograma da imagem equalizada')
plt.figure(1)
hgc = np.cumsum(hg)
plt.plot(hgc),plt.title('Histograma acumulado da imagem equalizada');

"""
Explanation: Finalmente, plotamos o histograma da imagem equalizada. Note o efeito mencionado acima em que o
histograma equalizado fica espalhado. Quando se calcula o histograma acumulado, nota-se daí que
o histograma de fato está normalizado.
End of explanation
"""


f = mpimg.imread('../data/angiogr.tif')
f = np.clip(f,75,255)
ia.adshow(f)
h = ia.histogram(f)
plt.plot(h);
#print('info:',ia.iaimginfo(f))

"""
Explanation: Um problema da formulação simplificada acima é que no caso da imagem original não ter nenhum
pixel igual a zero, a equalização da imagem usando esta formulação não irá fazer com que o
menor pixel seja zero. Veja o exemplo a seguir, onde o menor valor do pixel na imagem original é 75:
End of explanation
"""


n = f.size
T = 255./n * np.cumsum(h)
T = T.astype('uint8')
print('T:',T)
plt.plot(T),plt.title('Transformação de intensidade para equalizar')
g = T[f]
#print('info:', ia.iaimginfo(g))
ia.adshow(g, 'imagem equalizada')

"""
Explanation: Observe que a Transformação que equaliza a imagem, o seu primeiro valor
não zero é 8 (pois T[75] = 8). Isto faz com que o menor valor da imagem equalizada resultante seja 8 e não zero como
desejado:
End of explanation
"""


gn = ia.normalize(g)
#print 'info:',ia.iaimginfo(gn)
ia.adshow(gn, 'imagem equalizada e normalizada')
hgn = ia.histogram(gn)
plt.plot(hgn),plt.title('histograma');

"""
Explanation: Para fazer com que o valor do menor pixel da imagem equalizada seja zero, temos duas opções básicas:

Após a equalização, normalizar a imagem entre 0 e 255. Isto pode ser feito com  ia898:normalize.
Outra solução é alterar a equação do início desta página já incorporando esta normalização final. Esta
   equação é a mesma que aparece na descrição da Wikipedia, no link no final da página:
End of explanation
"""


wiki=np.array([[52,55,61,66,70,61,64,73],
              [63,59,55,90,109,85,69,72],
              [62,59,68,113,144,104,66,73],
              [63,58,71,122,154,106,70,69],
              [67,61,68,104,126,88,68,70],
              [79,65,60,70,77,68,58,75],
              [85,71,64,59,55,61,65,83],
              [87,79,69,68,65,76,78,94]])
print('wiki=\n',wiki)
h = ia.histogram(wiki)
n = wiki.size
T = 255./n * np.cumsum(h)
T = np.floor(T).astype('uint8')
g = T[wiki]
print('g=\n',g)
gn = ia.normalize(g)
print('gn=\n',gn)

"""
Explanation: Verificando a equação da wikipedia
End of explanation
"""


faux = g.ravel().astype('float')
minimum = min(faux)
maximum = max(faux)
lower = 0
upper = 255
gnn = np.round((faux - minimum) * (upper - lower) / (maximum - minimum) + lower,0)
gnn = gnn.reshape(g.shape).astype(np.int)
print('gnn=\n',gnn)

"""
Explanation: Comparando-se o resultado (gn) com o valor da Wikipedia, percebemos que existe uma
pequena diferença nos valores de alguns pixels. Esta diferença é devido ao fato que
na equação da Wikipedia, é usado um arredondamento (round) enquanto que na
função ia898:normalize, é usado um truncamento. A seguir foi feita uma outra
função similar à ianormalize, porém utilizando a função (round). Note que
neste caso o resultado confere com a Wikipedia.
End of explanation
"""


