

In [ ]:

    
import numpy as np
import matplotlib.pyplot as plt
import scipy.signal as signal
import skimage as ski
from skimage import io

Building edge for tests



In [ ]:

    
x=np.linspace(0,10,100)
arg=(x-5.3)/2;
A=10
y=A*np.exp(-arg*arg);
plt.plot(x,y)
y[y<1e-3]=0



In [ ]:

    
y.shape



In [ ]:

    
cy=np.cumsum(y)
plt.plot(x,cy)



In [ ]:

    
dcy=np.diff(cy)



In [ ]:

    
plt.plot(x,y,'.',x,0.5<y,x,np.ones(x.size)*0.5)



In [ ]:

    
(0.5<y).sum()*20/1000



In [ ]:

    
2*np.sqrt(np.log(2))*2



In [ ]:

    
np.log(np.exp(1))

Look at unit test results



In [ ]:

    
a=io.imread('../../UnitTests/data/slantededge.tif')
plt.figure(figsize=[12,12])
plt.subplot(221)
plt.imshow(a,cmap='gray')
plt.xticks(np.arange(0,a.shape[1],2),labels=[])
plt.yticks(np.arange(0,a.shape[0],2),labels=[])
plt.grid(color='orange', linestyle='-', linewidth=0.5)
plt.grid(True)
m= 24.6932 ; k= -0.0522558
print(np.arctan(k)*180/np.pi)
ybx=np.arange(0,a.shape[0])
xbx=m+k*ybx
plt.plot(xbx,ybx,'r',label='Fitted edge line')
plt.title('Image of slanted edge')
plt.legend()

plt.subplot(222)
plt.imshow(np.diff(a,axis=1),cmap='gray')
m= 25.6932 ; k= -0.0522558
print(np.arctan(k)*180/np.pi)
ybx=np.arange(0,a.shape[0])
xbx=m+k*ybx
plt.plot(xbx,ybx,'r',label='Fitted edge line')
plt.xticks(np.arange(0,a.shape[1],2),labels=[])
plt.yticks(np.arange(0,a.shape[0],2),labels=[])
plt.grid(color='orange', linestyle='-', linewidth=0.5)
plt.grid(True)
plt.title('Image of edge derivative')
plt.legend()


plt.subplot(2,2,3)
d=a
plt.step(np.arange(0,d.shape[1]),d[0,:], label='First line')
plt.step(np.arange(0,d.shape[1]),d[25,:], label='Middle line')
plt.step(np.arange(0,d.shape[1]),d[-1,:], label='Last line')
plt.title('Line profiles edge')
plt.legend()
plt.subplot(2,2,4)
d=np.diff(a,axis=1)
plt.step(np.arange(0,d.shape[1]),d[0,:], label='First line')
plt.step(np.arange(0,d.shape[1]),d[25,:], label='Middle line')
plt.step(np.arange(0,d.shape[1]),d[-1,:], label='Last line')

plt.title('Line profiles edge derivative')
plt.legend()

plt.tight_layout()
plt.savefig('edgedemo.pdf')
plt.savefig('edgedemo.png',dpi=150)



In [ ]:

    
plt.subplot(1,2,1)
d=a
plt.step(np.arange(0,d.shape[1]),d[0,:], label='First line')
plt.step(np.arange(0,d.shape[1]),d[25,:], label='Middle line')
plt.step(np.arange(0,d.shape[1]),d[-1,:], label='Last line')
plt.legend()
plt.subplot(1,2,2)
d=np.diff(a,axis=1)
plt.step(np.arange(0,d.shape[1]),d[0,:], label='First line')
plt.step(np.arange(0,d.shape[1]),d[25,:], label='Middle line')
plt.step(np.arange(0,d.shape[1]),d[-1,:], label='Last line')
plt.legend()



In [ ]:

    
b=io.imread('../../UnitTests/data/raw_edge.tif').astype(float)
plt.figure(figsize=[15,6])
plt.subplot(241)
plt.imshow(b)
m= 32.8755 ;k= -0.0226082
yax=np.arange(0,b.shape[0])
xax=m+k*yax
plt.plot(xax,yax,'r')

plt.subplot(245)
plt.imshow(np.diff(b,axis=1),clim=[-1600,200])
plt.plot(xax,yax,'r')


plt.subplot(1,4,2)
d = np.loadtxt('../../UnitTests/data/edgeprofile.txt', delimiter="\t")
plt.plot(d[:,0],(d[:,1]))

plt.subplot(1,4,3)
dd=np.diff(d[:,1])
plt.plot(d[0:-1,0],dd)

plt.subplot(1,4,4)
ftdd=np.abs(np.fft.fft(dd))
plt.plot(ftdd[0:len(ftdd)//2]/ftdd[0])
plt.plot([0,len(ftdd)//2-1],[0.1,0.1],'--')
plt.plot([0,len(ftdd)//2-1],[0.2,0.2],'--')



In [ ]:

    
len(ftdd)

Weighted CoG



In [ ]:

    
y2=y+np.random.normal(0,0.1*A,size=y.size)
c1=(x*y2).sum()/y2.sum()
plt.plot(x,y2)
y3=y2
y3[y2<0]=0
plt.plot(x,y3)

c2=(x*y3).sum()/y3.sum()
print(y2.mean(), y2.std())

y4=y2;
y4=y4-(y4.mean()+y4.std())
y4[y4<0]=0
c3=(x*y4).sum()/y4.sum()
print(c1,c2,c3)



In [ ]:

    
c1=[]
c2=[]
c3=[]
for i in np.arange(0,100) :
    y2=y+np.random.normal(0,0.1*A,size=y.size)
    c1.append((x*y2).sum()/y2.sum())
    
    y3=y2
    y3[y2<0]=0
    c2.append((x*y3).sum()/y3.sum())
    y4=y2;
    y4=y4-(y4.mean()+y4.std())
    y4[y4<0]=0
    c3.append((x*y4).sum()/y4.sum())



In [ ]:

    
plt.plot(c1,label='Plain')
plt.plot(c2,label='Non-zero')
plt.plot(c3,label='Unbias')
plt.legend()



In [ ]:

    
print(np.mean(c1),np.mean(c2),np.mean(c3))



In [ ]:



In [ ]: