pyJoinSpotCentres

import cv2
import math
def getDist(coord1, coord2):
    return ((coord1[0] - coord2[0]) ** 2 + (coord1[1] - coord2[1]) ** 2) ** 0.5

def getAngle(coord1, coord2):
    return math.atan2(coord1[0] - coord2[0], coord1[1] - coord2[1])

img = cv2.imread("FourSpotsV2.jpg")
grayimg = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
foundSpots = {}             # Dictionary of dark spots.
spotCount = 0
for row, aline in enumerate(grayimg):           # aline contains one horizontal line of pixels
    for col, pixel in enumerate(aline):         # pixel contains whiteness level of one pixel
        if pixel < 85:                          # Find a dark pixel
            for spotKey in foundSpots:          # Compare with all spots found so far
                pData = foundSpots[spotKey]     # [coords 1st pixel], [minx, miny], [maxx, maxy]
                spot = pData[0]                 # get coords of first pixel of spot
                if getDist((col, row), spot) < 100:         # Part of same spot?
                    pData[1][0] = min(pData[1][0], col)         # keep minimum x coord
                    pData[1][1] = min(pData[1][1], row)         # keep minimum y coord
                    pData[2][0] = max(pData[2][0], col)         # keep maximum x coord
                    pData[2][1] = max(pData[2][1], row)         # keep maximum y coord
                    foundSpots[spotKey] = pData
                    break                       # pixel is close to an already known spot
            else:
                spotCount += 1                  # start of a new spot
                foundSpots[spotCount] = [[col, row], [col, row], [col, row]]            # Start a new spot

spotCentres = []
for spotKey in foundSpots:              # # Find the centre pixel of each dark spot
    pData = foundSpots[spotKey]
    avgCol = (pData[1][0] + pData[2][0]) // 2       # mean of min and max column
    avgRow = (pData[1][1] + pData[2][1]) // 2       # mean of min and max row
    spotCentres.append((avgCol, avgRow))

# Find the coordinates of the centre of the whole group of spots
totCol = sum([spot[0] for spot in spotCentres])     # Add up all the spot columns
totRow = sum([spot[1] for spot in spotCentres])     # Add up all the spot rows
figCentre = (totCol // len(spotCentres), totRow // len(spotCentres))

# Order the found spots in rotation order as viewed from the centre point
spotsByAngle = sorted([(getAngle(spot, figCentre), spot) for spot in spotCentres])
last = spotsByAngle[-1][1]
for item in spotsByAngle:       # draw line joining spot pairs, in rotation sequence
    cv2.line(img, last, item[1], (0, 0, 255), 1)
    last = item[1]

cv2.imshow("img", img)
cv2.waitKey(0)