################################################################################################### 
# Procedure to detect the cover orientation
################################################################################################### 
import ch.psi.pshell.imaging.Utils.integrateVertically as integrateVertically
img.backgroundEnabled=False

REF = (0,96,125)

line = load_image("{images}/line.png", title="Line")
#line = load_image("{images}/line360.png", title="Line")

line.getProcessor().setBackgroundValue(0.0)


#ip = get_image()
ip = integrate_frames(10)
ip = grayscale(ip, True)


smooth(ip)
#bandpass_filter(ip, 30, 1000)
edges(ip)

#invert(ip)
#auto_threshold(ip, method = "Default")


#auto_threshold(ip, method = "Li")
auto_threshold(ip, method = "MaxEntropy")

"""
for m in AutoThresholder.getMethods():   
    print m
    aux = ip.duplicate()
    auto_threshold(aux, method = m)
    binary_fill_holes(aux, dark_background=False)
    renderer = show_panel(aux.bufferedImage)
    time.sleep(1.0)
"""
#binary_dilate(ip, dark_background=False)
#binary_fill_holes(ip, dark_background=False)
#binary_open(ip, dark_background=Tr)

renderer = show_panel(ip.bufferedImage)








#line = sub_image(line, 325, 325, 512, 512)
#ip = sub_image(ip, 325, 325, 512, 512)
line = sub_image(line, 453, 453, 256, 256)
ip = sub_image(ip, 453, 453, 256, 256)
#op = op_fft(ip, line, "correlate") 


renderer = show_panel(ip.bufferedImage)

#renderer = show_panel(op.bufferedImage)
#line.show()
ydata = []
xdata = range (0,180,1)
for i in xdata:
    l = line.duplicate()
    l.getProcessor().setBackgroundValue(0.0)
    l.getProcessor().rotate(float(i))
    op = op_fft(ip, l, "correlate") 
    bi = op.getBufferedImage()
    p = integrateVertically(bi)
    ydata.append(sum(p))
    
    #renderer = show_panel(op.bufferedImage)
    #time.sleep(0.001)

def moving_average(arr, n) :
    ret = []
    for i in range(len(arr)):
        ret.append(mean(arr[max(i-n,0):min(i+n,len(arr)-1)]))
    return ret    
av = moving_average(ydata, 1)
    
p = plot(ydata, xdata=xdata)[0]

p.addSeries(LinePlotSeries("Moving Average"))
p.getSeries(1).setData(xdata, av)      

peaks = estimate_peak_indexes(ydata, xdata,  (min(ydata) + max(ydata))/2,  25.0)
left, right = min(peaks), max(peaks)
if xdata[left]<5 and xdata[right]>(xdata[-1]-5):
    #del peaks[0 if ydata[right] >  ydata[left]  else -1]
    peaks.remove(right if ydata[right] >  ydata[left]  else left)

peaks = sorted(peaks[:3])
peaks_x = map(lambda x:xdata[x], peaks)
peaks_y = map(lambda x:ydata[x], peaks)

print "Peaks", peaks        
print "Peak indexes: " + str(peaks_x)
print "Peak values:  " + str(peaks_y)


for i in range(len(peaks)):
    peak = xdata[peaks[i]]
    p.addMarker(peak, None, "N="+str(round(peak,2)), Color(80,0,80))
    if ((peaks[i]>160) and (REF[i]<20)):
        peaks[i] = peaks[i] - 180.0

print "Peaks x: " + str(peaks_x)


d = mean(arrabs(arrsub(REF, peaks_x)))

print "Angle = ", d