import numpy as np


def correlation(image, angrange, brange, sigma):
    xsize, ysize = image.shape
    
    subimage = image.astype(np.float32)
    
    medianint = np.median(image[image > 0])
    angles = np.arange(angrange[0], angrange[2], angrange[1])
    trans = np.arange(brange[0], brange[2], brange[1])
    corrmap = np.zeros((len(angles), len(trans)))
    
    for thetaind in range(len(angles)):
        for bind in range(len(trans)):
            theta = angles[thetaind]
            b = trans[bind]
            refline = reference_line_gen(xsize, ysize, theta, b, sigma)
            temp = refline * subimage
            temp2 = np.sum(temp)
            temp3 = temp2 / medianint
            corrmap[thetaind, bind] = temp3
            
    return corrmap


def reference_line_gen(xsize, ysize, theta, b, sigma):
    Xmesh, Ymesh = np.meshgrid(np.arange(1, xsize+1), np.arange(1, ysize+1))
    Xmesh = Xmesh - xsize / 2
    Ymesh = Ymesh - ysize / 2
    
    angle = np.deg2rad(theta)
    dist = b - (Xmesh * np.sin(angle)) - (Ymesh * np.cos(angle))
    distance = dist ** 2 / (sigma ** 2)
    result = np.exp(-distance)
    
    return result


# Procedure 3: Peak Extraction Algorithm
def max_intensity_finding(correlation_map, dist_input=None):
    if dist_input is not None:
        distance_threshold = dist_input
    else:
        distance_threshold = 10
    
    size_corr = np.array(correlation_map.shape)
    label_map = np.zeros(size_corr)
    xgrid, ygrid = np.meshgrid(np.arange(1, size_corr[1]+1), np.arange(1, size_corr[0]+1))
    
    for indx in range(size_corr[1]):
        for indy in range(size_corr[0]):
            dist = np.sqrt((xgrid - indx) ** 2 + (ygrid - indy) ** 2)
            pixuseful = dist < distance_threshold
            if correlation_map[indy, indx] >= np.max(correlation_map[pixuseful]):
                label_map[indy, indx] = 1
    
    label_map[correlation_map == 0] = 0
    
    return label_map


def sub_window_analysis(correlation_map, dist_input, image, subimgsz, sigma):
    imgsize = np.array(image.shape)
    sgimgsize = imgsize.astype(np.float32)
    sgsubimgsz = np.array(subimgsz).astype(np.float32)
    sgtestimage = image.astype(np.float32)
    
    if imgsize[0] % subimgsz[0] != 0 or imgsize[1] % subimgsz[1] != 0:
        print("The subimage size does not match the whole image size")
        return
    
    angles = np.arange(angrange[0], angrange[2], angrange[1])
    trans = np.arange(brange[0], brange[2], brange[1])
    corr = np.zeros((len(angles), len(trans), imgsize[0] // subimgsz[0], imgsize[1] // subimgsz[1]))
    
    xsize = sgimgsize[1]
    ysize = sgimgsize[0]
    Xmesh, Ymesh = np.meshgrid(np.arange(xsize), np.arange(ysize))
    Xmesh = np.mod(Xmesh, subimgsz[1])
    Ymesh = np.mod(Ymesh, subimgsz[0])
    Xmesh = Xmesh - sgsubimgsz[1] / 2
    Ymesh = Ymesh - sgsubimgsz[0] / 2
    
    for thetaind in range(len(angles)):
        for bind in range(len(trans)):
            theta = angles[thetaind]
            b = trans[bind]
            angle = np.deg2rad(theta)
            dist = b - (Xmesh * np.sin(angle)) - (Ymesh * np.cos(angle))
            distpw = (dist ** 2) / (sigma ** 2)
            refline = np.exp(-distpw)
            temp = refline * sgtestimage
            temp2 = np.sum(temp, axis=1)
            temp3 = np.squeeze(np.sum(temp2, axis=2))
            corr[thetaind, bind, :, :] = temp3
    
    corrsz = np.array(corr.shape)
    spcorr = np.reshape(np.transpose(corr, (0, 2, 1, 3)), (corrsz[0] * corrsz[2], corrsz[1] * corrsz[3]))
    colsubimgs = np.reshape(sgtestimage, (subimgsz[0], imgsize[0] // subimgsz[0], subimgsz[1], imgsize[1] // subimgsz[1]))
    
    return spcorr, colsubimgs

# Example Usage
image = cv2.imread("0.0/01.pgm", cv2.COLOR_BGR2GRAY)
angrange = (0, 10, 180)   # Range of angles for correlation analysis
brange = (0, 1, 10)       # Range of translation values for correlation analysis
sigma = 1.0               # Sigma value for reference line generation

# Perform whole image correlation analysis
corr_map = correlation(image, angrange, brange, sigma)
print("Correlation map shape:", corr_map.shape)

# Extract maximum intensity peaks from the correlation map
label_map = max_intensity_finding(corr_map)
print("Label map shape:", label_map.shape)

# Plot the correlation map and label map
plt.figure(figsize=(10, 4))
plt.subplot(1, 2, 1)
plt.imshow(corr_map, cmap="gray")
plt.title("Correlation Map")
plt.colorbar()
plt.subplot(1, 2, 2)
plt.imshow(label_map, cmap="gray")
plt.title("Label Map")
plt.colorbar()
plt.show()

# Perform sub-window analysis
subimgsz = (10, 10)  # Sub-window size for analysis
spcorr, colsubimgs = sub_window_analysis(corr_map, brange, image, subimgsz, sigma)
print("Spaced correlation shape:", spcorr.shape)
print("Column sub-images shape:", colsubimgs.shape)

Вот мой код, получаю ошибку:
---------------------------------------------------------------------------
AxisError                                 Traceback (most recent call last)
<ipython-input-11-5beb91802097> in <cell line: 132>()
    130 # Perform sub-window analysis
    131 subimgsz = (10, 10)  # Sub-window size for analysis
--> 132 spcorr, colsubimgs = sub_window_analysis(corr_map, brange, image, subimgsz, sigma)
    133 print("Spaced correlation shape:", spcorr.shape)
    134 print("Column sub-images shape:", colsubimgs.shape)

3 frames
/usr/local/lib/python3.10/dist-packages/numpy/core/fromnumeric.py in _wrapreduction(obj, ufunc, method, axis, dtype, out, **kwargs)
     84                 return reduction(axis=axis, out=out, **passkwargs)
     85 
---> 86     return ufunc.reduce(obj, axis, dtype, out, **passkwargs)
     87 
     88 

AxisError: axis 2 is out of bounds for array of dimension 1