I want to connect the dots in the image horizontally

from PIL import Image, ImageDraw

# Image Settings
width, height = 1000, 1000
dot_size = 20
dot_interval = 20
shift_right = 20

# image generation
image = Image.new("RGB", (width, height), "white")
draw = ImageDraw.Draw(image)

# dot description
for i in range(shift_right, width, dot_size + dot_interval):
    for j in range(0, height, dot_size + dot_interval):
        draw.rectangle([i, j, i + dot_size, j + dot_size], fill="black")

# image storage
image.save("shifted_square_dot_pattern.png")


import cv2
import numpy as np

# image loading
image = cv2.imread("shifted_square_dot_pattern.png", cv2.IMREAD_GRAYSCALE)
rows, cols = image.shape

# Mapping Configuration
def gaussian_skew_mapping(x, y, amplitude, sigma):
    return amplitude * np.exp(-((y - rows // 2) ** 2) / (2 * sigma ** 2))

# deformation matrix
transformation_matrix = np.float32([[1, 0, 0], [0, 1, 0]])
skewed_image = np.zeros_like(image)

for y in range(rows):
    for x in range(cols):
        deformation = int(gaussian_skew_mapping(x, y, 20, 100))
        new_x = x + deformation
        if 0 <= new_x < cols:
            skewed_image[y, new_x] = image[y, x]

# Save Deformed Images
cv2.imwrite("gaussian_skewed_shifted_square_dot_pattern.png", skewed_image)

import matplotlib.pyplot as plt
img_h = cv2.imread("gaussian_skewed_shifted_square_dot_pattern.png")
gray_h = cv2.cvtColor(img_h, cv2.COLOR_BGR2GRAY)
# binarization processing
# the process of converting a characteristic part of an image so as to extract the part of interest
threshold=130
ret, gray_h_thresh = cv2.threshold(gray_h,threshold,255,cv2.THRESH_BINARY_INV)
plt.imshow(gray_h_thresh)
#contour extraction
contours, hierarchy = cv2.findContours(gray_h_thresh, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
#Display the number of contours
print(len(contours))

#Remove non-left dot space
count =0
areas = []
list_contours = list(contours)
for i in range(len(contours)):
    area = cv2.contourArea(contours[i], True)
    areas.append(area)
   
    if area < -500:
        print(count)
        del list_contours[i]
    count+=1

#Contours size 2x4 and close coordinates removed using Euclidean distance
for i in range(len(trimed_list_contours)):
    if not trimed_list_contours[i].shape[0] ==4:
        sq_index =[]
        for j in range(trimed_list_contours[i].shape[0]-1):
            euclidean_distance = np.sqrt(sum(sum((trimed_list_contours[i][j]-trimed_list_contours[i][j+1])**2)))
            sq_index.append(euclidean_distance)
    
        for k in range(trimed_list_contours[i].shape[0]-4):
            argmin = np.argmin(sq_index)
            trimed_list_contours[i] =np.delete(trimed_list_contours[i],argmin,0)
            sq_index = np.delete(sq_index,argmin,0)

#connect dots in parallel
count = 0
copy_gray_h=gray_h.copy()
#for i in range(0,26):
for i in range(len(trimed_list_contours)):
    count+=1
    if count == 25:
        count =0
        
        continue
    else:
        
        # # Assuming you have a list or tuple of points from OpenCV, for example:
        opencv_points = [trimed_list_contours[i][0][0], trimed_list_contours[i][1][0], trimed_list_contours[i+1][2][0],trimed_list_contours[i+1][3][0]]

       
        opencv_points =np.vstack(opencv_points)
        cv2.fillPoly(copy_gray_h, [opencv_points], (0, 0, 0))
        
plt.imshow(copy_gray_h)
plt.show()