in multifocal/code/multifocal_decomposition_from_focalstack/train/util.py [0:0]
def circ_filter(w): # reimplement from matlab function fspecial('disk', rad)
rad = w/2
rad = np.maximum(rad, 1e-6) # prevent zero-divide
rad2 = rad**2
crad = np.ceil(rad-0.5)
l = np.linspace(-crad,crad,2*crad+1)
[x,y] = np.meshgrid(l,l)
maxxy = np.maximum(np.abs(x), np.abs(y))
minxy = np.minimum(np.abs(x), np.abs(y))
m1 = (rad2 < (maxxy+0.5)**2 + (minxy-0.5)**2)*(minxy-0.5) + (rad2 >= (maxxy+0.5)**2 + (minxy-0.5)**2)*np.sqrt(np.abs(rad2 - (maxxy + 0.5)**2))
m2 = (rad2 > (maxxy-0.5)**2 + (minxy+0.5)**2)*(minxy+0.5) + (rad2 <= (maxxy-0.5)**2 + (minxy+0.5)**2)*np.sqrt(np.abs(rad2 - (maxxy - 0.5)**2))
sgrid = (rad2*(0.5*(arcsin_hack(m2/rad) - arcsin_hack(m1/rad)) + 0.25*(np.sin(2*arcsin_hack(m2/rad)) - np.sin(2*arcsin_hack(m1/rad)))) -
(maxxy-0.5)*(m2-m1) + (m1-minxy+0.5))*((np.logical_or(np.logical_and((rad2 < (maxxy+0.5)**2 + (minxy+0.5)**2) , (rad2 > (maxxy-0.5)**2 + (minxy-0.5)**2)),
(np.logical_and(np.logical_and((minxy==0),(maxxy-0.5 < rad)),(maxxy+0.5>=rad))))))
sgrid = sgrid + ((maxxy+0.5)**2 + (minxy+0.5)**2 < rad2)
cradint = np.int(crad)
sgrid[cradint,cradint] = np.minimum(np.pi*rad2,np.pi/2)
if ((crad>0) and (rad > crad-0.5) and (rad2 < (crad-0.5)**2+0.25)):
m1 = np.sqrt(rad2 - (crad - 0.5)**2)
m1n = m1/rad
sg0 = 2*(rad2*(0.5*arcsin_hack(m1n) + 0.25*np.sin(2*arcsin_hack(m1n)))-m1*(crad-0.5))
sgrid[2*cradint,cradint] = sg0
sgrid[cradint,2*cradint] = sg0
sgrid[cradint,0] = sg0
sgrid[0,cradint] = sg0
sgrid[2*cradint-1,cradint] = sgrid[2*cradint-1,cradint] - sg0
sgrid[cradint,2*cradint-1] = sgrid[cradint,2*cradint-1] - sg0
sgrid[cradint,1] = sgrid[cradint,1] - sg0
sgrid[1,cradint] = sgrid[1,cradint] - sg0
sgrid[cradint,cradint] = np.minimum(sgrid[cradint,cradint],1)
kernel = (sgrid/np.sum(sgrid)).astype(dtype=np.float32)
return kernel.shape[0], kernel