find circles with Hough transform

Question

the code receives an img with a circle in greyscale and need to output an

Returns: img_c: RGB image array in float format (range: 0..1) which contains the original image and the marked circles which were found. circles: a list - each element represent a circle and contains a list of 3 values: origin_x, origin_y, radius and I get some error and I don't know why, I already implement edge_detect function using canny edge detector, and now I need to implement this, I can't use cv2 just numpy, and from scipy.signal import convolved thanks!!

def find_circles(img: np.array, t_low: np.float32, t_high: np.float32) -> (np.array, list):
    edges, tg_theta = edge_detect(img, t_low, t_high)
    r_min = 14
    ncols, nrows = edges.shape
    r_max = round(min(ncols, nrows) / 2)
    M = np.zeros((ncols, nrows, r_max))
    for x in range(ncols):
      for y in range(nrows):
        if edges[x,y] == 0: 
          continue
        if abs(tg_theta[x,y]) <= 1:
          for cx in range(max(x - r_max, 1), min(x + r_max, ncols)):
            cy = round(tg_theta[x, y] * (cx - x) + y)
            r = round(np.sqrt((cx - x)**2 +(cy - y)**2))
            if 0 < cy < nrows and 0 < r < r_max:
              M[cx, cy, r] +=1
          else:
            for cy in range(max(y - r_max, 1), min(y + r_max, nrows)):
              cx = round((cy-y)/tg_theta[x, y] + x)
              r = round(np.sqrt((cx - x)**2 +(cy - y)**2))
              if 0 < cx < ncols and 0 < r < r_max:
                M[cx, cy, r] +=1
    kernel = np.ones((5, 5, 5)) / 125.0
    M_smooth = convolve(M, kernel, mode='same')
    img_c = np.zeros_like(edges)
    circles=[]
    for cx in range(ncols):
      for cy in range(nrows):
        for r in range(r_min,r_max):
          if M_smooth[cx, cy, r] <= r/25: 
            continue
          D=15
          neigbours = M_smooth[ max(cx-D, 1):min(cx+D, ncols), max(cy-D, 1):min(cy+D, nrows), max(r-D, r_min):min(r+D, r_max)]
          if M_smooth[cx, cy, r] < np.max(neigbours):
            continue
          circles.append([cx, cy, r])
          theta = 0
          while theta <= 2*np.pi:
              pt = np.round(np.array([cx,cy]) + np.array([r*np.cos(theta),r*np.sin(theta)]))
              irow = int(pt[1])
              icol = int(pt[0])
              if 0 < irow < nrows and 0 < icol < ncols:
                  img_c[irow, icol] = 1
              theta = theta + 1/r
    return img_c, circles

edges is array in int format (values: 0, 1) - binary image with edges pixels set to 1. I get this error

TypeError: loop of ufunc does not support argument 0 of type int which has no callable sqrt method

in line

r = round(np.sqrt((cx - x)**2 +(cy - y)**2))

EDIT: now the code run, I change it to float like u said

def find_circles(img: np.array, t_low: np.float32, t_high: np.float32) -> (np.array, list):
    edges, tg_theta = edge_detect(img, t_low, t_high)
    r_min = 14
    ncols, nrows = edges.shape
    r_max = round(min(ncols, nrows) / 2)
    M = np.zeros((ncols, nrows, r_max))
    for x in range(ncols):
      for y in range(nrows):
        if edges[x,y] == 0: 
          continue
        if abs(tg_theta[x,y]) <= 1:
          for cx in range(max(x - r_max, 1), min(x + r_max, ncols)):
            cy = round(tg_theta[x, y] * (cx - x) + y)
            r = round(np.sqrt(float((cx - x)**2) +float((cy - y)**2)))
            if 0 < cy < nrows and 0 < r < r_max:
              M[cx, cy, r] +=1
          else:
            for cy in range(max(y - r_max, 1), min(y + r_max, nrows)):
              cx = round((cy-y)/tg_theta[x, y] + x)
              r = round(np.sqrt(float((cx - x)**2) +float((cy - y)**2)))
              if 0 < cx < ncols and 0 < r < r_max:
                M[cx, cy, r] +=1
    kernel = np.ones((5, 5, 5)) / 125.0
    M_smooth = convolve(M, kernel, mode='same')
    img_c = np.zeros_like(img)
    circles=[]
    for cx in range(ncols):
      for cy in range(nrows):
        for r in range(r_min,r_max):
          if M_smooth[cx, cy, r] <= r/25: 
            continue
          D=15
          neigbours = M_smooth[ max(cx-D, 1):min(cx+D, ncols), max(cy-D, 1):min(cy+D, nrows), max(r-D, r_min):min(r+D, r_max)]
          if M_smooth[cx, cy, r] < np.max(neigbours):
            continue
          circles.append([cx, cy, r])
          theta = 0
          while theta <= 2*np.pi:
              pt = np.round(np.array([cx,cy]) + np.array([r*np.cos(theta),r*np.sin(theta)]))
              irow = int(pt[1])
              icol = int(pt[0])
              if 0 < irow < nrows and 0 < icol < ncols:
                  img_c[icol, irow] = 1
              theta = theta + 1/r
    return img_c, circles

and my main code is:

img = read_file(file_name='balls5.tiff') / 255
img_c, circles = find_circles(img=img, t_low=9/255, t_high=10/255)

img_rgb = np.dstack((img, img, img))
img_rgb[:,:,2] += img_c
img_rgb = img_rgb.clip(0,1)
plt.imshow(img_rgb);
print(circles)

but the problem is I get the output as the same image as the input

I get:

the correct answer is:

EDIT again: now I change to this:

def find_circles(img: np.array, t_low: np.float32, t_high: np.float32) -> (np.array, list):
    edges, tg_theta = edge_detect(img, t_low, t_high)
    r_min = 14
    ncols, nrows = img.shape
    r_max = round(min(ncols, nrows) / 2)
    M = np.zeros((ncols, nrows, r_max))
    e_col, e_row = edges.shape
    for x in range(e_col):
        for y in range(e_row):
            if edges[x, y] == 0:
                continue
            if abs(tg_theta[x, y]) <= 1:
                for cx in range(max(x - r_max, 0), min(x + r_max, ncols)):
                    cy = round(tg_theta[x, y] * (cx - x) + y)
                    r = round(np.sqrt(float((cx - x) ** 2) + float((cy - y) ** 2)))
                    if 0 < cy < nrows and 0 < r < r_max:
                        M[cx, cy, r] += 1
            else:
                for cy in range(max(y - r_max, 0), min(y + r_max, nrows)):
                    cx = round((cy - y) / tg_theta[x, y] + x)
                    r = round(np.sqrt(float((cx - x) ** 2) + float((cy - y) ** 2)))
                    if 0 < cx < ncols and 0 < r < r_max:
                        M[cx, cy, r] += 1
    kernel = np.ones((5, 5, 5)) / 125.0
    M_smooth = convolve(M, kernel, mode='same')
    img_c = np.zeros_like(img)
    circles = []
    for cx in range(ncols):
        for cy in range(nrows):
            for r in range(r_min, r_max):
                if M_smooth[cx, cy, r] <= r / 25:
                    continue
                D = 15
                neighbours = M_smooth[max(cx - D, 1):min(cx + D, ncols), max(cy - D, 1):min(cy + D, nrows),
                             max(r - D, r_min):min(r + D, r_max)]
                if M_smooth[cx, cy, r] < np.max(neighbours):
                    continue
                circles.append([cx, cy, r])
                theta = 0
                while theta <= 2 * np.pi:
                    pt = np.round(np.array([cx, cy]) + np.array([r * np.cos(theta), r * np.sin(theta)]))
                    irow = int(pt[1])
                    icol = int(pt[0])
                    if 0 < irow < nrows and 0 < icol < ncols:
                        img_c[icol, irow] = 1
                    theta = theta + 1 / r
    return img_c, circles

and I get this output: