# ============================ # (R)etrieve (A)ngles from (C)ontour sampling (P)oints # TRY TO RETRIEVE ANGLES AT WHICH # CONTOUR SAMPLING TAKES PLACE # AROUND THE OLIVA CONTOUR # ---------------------------------------------------------- import argparse import imutils import cv2 import numpy as np import matplotlib.pyplot as plt import math from PIL import Image, ImageDraw from numpy import asarray # ============================================ # Set general defaults for figure size and DPI # ============================================ plt.rcParams["figure.figsize"] = (3,4) plt.rcParams["figure.dpi"] = 200 # ============================================ # Define a function that returns the angle # of a line betwween two points # ============================================ def get_angle(x1, y1, x2, y2): # Calculate change in y and x dy = y2 - y1 dx = x2 - x1 # atan2 handles all quadrants correctly angle_rad = math.atan2(dy, dx) return angle_rad # ==================================================== # Define a function for resizing images maintaining by # default the same aspect ratio of the original image # https://stackoverflow.com/questions/44650888/resize-an-image-without-distortion-opencv # ==================================================== def image_resize(image, width = None, height = None, inter = cv2.INTER_AREA): # Initialize the dimensions of the image to be resized and # grab the image size dim = None (h, w) = image.shape[:2] # If both the width and height are None, then return the # original image if width is None and height is None: return image # Check to see if the width is None if width is None: # Calculate the ratio of the height and construct the # dimensions r = height / float(h) dim = (int(w * r), height) # Otherwise, the height is None else: # Calculate the ratio of the width and construct the # dimensions r = width / float(w) dim = (width, int(h * r)) # Resize the image resized = cv2.resize(image, dim, interpolation = inter) # Return the resized image return resized # ============================================= # FILL A PYTHON LIST WITH ANGLES IN RADIANS # TO THE SAMPLED POINTS ON THE CONTOUR # ============================================= ANG_RAD_TO_SAMPLED_PTS = [] ANG_RAD_TO_SAMPLED_PTS = [0 for i in range(100)] # ============================================= # FILL A PYTHON LIST WITH SPECIMEN UNIQUE ID's # ============================================= SPECIMENS = [] SPECIMENS.append("VI0642") # ===================================== # FILL A PYTHON LIST WITH SPECIES NAMES # ===================================== SPECIES = [] SPECIES.append("vidua") # =========================================================== # PREPARE A 2-D ARRAY (=LIST) # FOR THE 100 EUCLIDEAN DISTANCES (as many columns + 2, first # column will contain the species and second column # the specimen ID) FOR THE 199 SILHOUETTES (as many rows) # =========================================================== print("INITIALIZATION - CREATION OF EUCLIDEAN DISTANCES ARRAY (=LIST)") rows_euc_dist, cols_euc_dist = (199, 102) arr_euc_dist = [[0 for i in range(cols_euc_dist)] for j in range(rows_euc_dist)] howmany = len(SPECIMENS) for num_spec in range(0,howmany): print("MAIN LOOP - ITEM " + str(num_spec) + " OF " + str(howmany)) #input("Press Enter to continue...") print("MAIN LOOP - INITIALIZING ARRAY ROWS") # --> Euclidean distances array : put current specimen id # --> in column 0 and current species in column 1 of row num_spec arr_euc_dist[num_spec][0] = SPECIES[num_spec] arr_euc_dist[num_spec][1] = SPECIMENS[num_spec] # --> Prepare image file name curr_img_name ="INPUT_SILHOUETTES/"+SPECIMENS[num_spec]+".jpg" print("MAIN LOOP - Current silhouette being processed= " + curr_img_name) # --> Read original image curr_image = cv2.imread(curr_img_name) #cv2.imshow("1 - Original image from INPUT_SILHOUETTES folder ", curr_image) #cv2.waitKey(0) # --> Transform to greyscale. gray = cv2.cvtColor(curr_image, cv2.COLOR_BGR2GRAY) #cv2.imshow("2 - Graytones Image", gray) #cv2.waitKey(0) # --> Blur to reduce noise blurred = cv2.GaussianBlur(gray, (5, 5), 1) #cv2.imshow("3 - Graytones, Blurred Image", blurred) #cv2.waitKey(0) # ===================================== # Apply threshold to binarize the image # ===================================== # --> The following parameters leave "holes" # --> (=generate more than one contour): # --> blurred, 60, 255, cv2.THRESH_BINARY leaves "holes" # --> also tested cv2.THRESH_TOZERO thresh = cv2.threshold(blurred,200,250,cv2.THRESH_BINARY)[1] #cv2.imshow("4 - Graytones, Blurred and Thresholded Image", thresh) #cv2.waitKey(0) # ========================== # Apply CANNY edge detection # ========================== # --> Canny edge detection has the purpose # --> of cropping the original image at the contour # --> minimum bounding box, so that it will contain # --> just the silhouette. # --> Variable name "mid" refers to a mid-range # --> setting of the parameters (initially, also # --> "low" and "high" were tested). # --> Aperture_size is set at the default value of 3 # --> and it could have been omiotted. It's explicitly # --> declared because initially also different values # --> were tested aperture_size = 3 mid = cv2.Canny(thresh, 30, 150, apertureSize=aperture_size) #cv2.imshow("5 - AFTER CANNY DETECTION WITH 30, 150 Edge Map", mid) #cv2.waitKey(0) # --> Find contour. The operation will be repeated later, # --> after the image will have been resized to 1000ph # --> of height. contoursM, hierarchy = cv2.findContours(mid, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_NONE) # --> Debug: one print line # --> We need ONE and only one continuous contour print("FULL-SIZE IMAGE: Number of Contours Found (MID) = " + str(len(contoursM))) # --> ContoursM is an array with all the contours (in this case, 1) # --> found in the mid edge map. To compute the bounding box, we must # --> use a for loop. for c in contoursM: # --> Get the bounding rectangle ("minimum bounding box") x, y, w, h = cv2.boundingRect(c) # --> Crop at the minimum bounding box, # --> WITH A LITTLE SLACK of 3 pixel in # --> each direction mid_cropped = mid[y-3:y+h+3,x-3:x+w+3] #cv2.imshow("6 - Cropped Mid Edge Map Image", mid_cropped) #cv2.waitKey(0) # --> Normalize image size at height = 1000 px mid_crop_res = image_resize(mid_cropped, height = 1000) #cv2.imshow("7 - Cropped Mid Edge Map Image, resized h=1000px", mid_crop_res) #cv2.imwrite("mid_crop_res.jpg",mid_crop_res) #cv2.waitKey(0) # ======================================================== # Re-compute the contours on the cropped and resized image # ======================================================== # --> ContoursMCR refers to the Mid setting (M), # --> cropped (C) and resized (R) image contoursMCR, hierarchy = cv2.findContours(mid_crop_res, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_NONE) # --> Debug: one print line # --> We need ONE and only one continuous contour print("RESIZED IMAGE: Number of Contours Found (contoursMCR) = " + str(len(contoursMCR))) # ===================================================== # ===================================================== # ===================================================== # ===================================================== # EUCLIDEAN DISTANCES CURVE GENERATION # ===================================================== # ===================================================== # ===================================================== # ===================================================== # ===================================================== # For the whole image, we meed to sample 100 points of # the contour, at equidistant intervals. We start by # computing the centroid of the cropped and resized image # ===================================================== # --> ContoursMCR is an array of contours, in our case containing # --> just one contour: anyway, a for loop is needed for c in contoursMCR: # --> Compute the center of each contour M = cv2.moments(c) cX = int(M["m10"] / M["m00"]) cY = int(M["m01"] / M["m00"]) # --> Debug: two print lines print("EUCLIDEAN DISTANCES CURVE GENERATION - Contour center X = " + str(cX)) print("EUCLIDEAN DISTANCES CURVE GENERATION - Contour center Y = " + str(cY)) # --> Draw the contour cv2.drawContours(mid_crop_res, [c], -1, (255, 255, 255), 3) #cv2.drawContours(mid_crop_res, [c], -1, (0, 255, 0), 20) # --> Show the center cv2.circle(mid_crop_res, (cX, cY), 3, (255, 0, 255), -1) cv2.putText(mid_crop_res, "Center", (cX - 20, cY - 20), cv2.FONT_HERSHEY_SIMPLEX, 0.5, (255, 0, 255), 2) #cv2.imshow('EUCLIDEAN DISTANCES CURVE GENERATION - Contours Crop Res Mid Edge Map', mid_crop_res) #cv2.waitKey(0) # ================================== # EQUIDISTANT CONTOUR POINT SAMPLING # ================================== # Source - https://stackoverflow.com/a/77946738 # Posted by Tino D, modified by community. See post 'Timeline' for change history # Retrieved 2025-11-23, License - CC BY-SA 4.0 # -------------------------------------------- # --> To sample 100 points on the contour evenly, we # --> need to unpack contour to X and Y coordinates X,Y = np.array(contoursMCR).T # --> Debug: two print lines print("EUCLIDEAN DISTANCES CURVE GENERATION - Array X shape = " + str(X.shape)) print("EUCLIDEAN DISTANCES CURVE GENERATION - Array Y shape = " + str(Y.shape)) X = X[0] # ==//== Y = Y[0] # ==//== # --> Debug: two print lines print("EUCLIDEAN DISTANCES CURVE GENERATION - len(Y) = " + str(len(Y))) print("EUCLIDEAN DISTANCES CURVE GENERATION - len(X) = " + str(len(X))) # --> We then must specify how many # --> points we want to sample numPoints = 100 # --> Generate uniformally distant indices until # --> len(X) - len(X)//numPoints to account # --> for the end of the contour resampleIndices = np.linspace(0, len(Y) - len(Y)//numPoints, numPoints, dtype=int) resampledX = X[resampleIndices] # ==//== resampledY = Y[resampleIndices] # ==//== # --> Debug: two print lines print("EUCLIDEAN DISTANCES CURVE GENERATION - len(resampledY) = " + str(len(resampledY))) print("EUCLIDEAN DISTANCES CURVE GENERATION - len(resampledX) = " + str(len(resampledX))) # ===================================================== # Fill an array of euclidean distances between each # evenly spaced contour points and the centroid cx, cy # ===================================================== # IMPORTANT OBSERVATION: by checking the output of the # following loop, it was ascertained that the points # on the contour are sampled at increasing (sub?)equal # angular distances from the topmost, in a clockwise # 360° rotation. This predictable and uniform progression # is exactly what's needed to compare the resulting # monodimensional array of euclidean distances that # represent the shape. # ------------------------------------------------------ # Initialize the local vector for the current contour Eucl_Distances = [] # --> Loop thru the 100 uniformally distant points for ind in range(0,numPoints): # --> Debug: two print lines #print("EUCLIDEAN DISTANCES CURVE GENERATION - resampledX["+str(ind)+"]="+str(resampledX[ind])) #print("EUCLIDEAN DISTANCES CURVE GENERATION - resampledY["+str(ind)+"]="+str(resampledY[ind])) # --> Calculate Euclidean distance to the centroid currEucDist = np.sqrt((resampledX[ind]-cX)**2 + (resampledY[ind]-cY)**2) # --> Debug: one print line #print("EUCLIDEAN DISTANCES CURVE GENERATION - Eucl_Distances["+str(ind)+"]="+str(currEucDist)) # --> Append the Eucliedan distance to the vector for the current contour Eucl_Distances.append(currEucDist) # --> Save euclidean distances in the global array # --> as string with three decimal places # --> Append starting from column three # --> (one is species name, two is specimen ID) arr_euc_dist[num_spec][ind+2] = "%.3f" % currEucDist ANG_RAD_TO_SAMPLED_PTS[ind]=get_angle(cX,cY,resampledX[ind],resampledY[ind]) print("From centroid X ="+str(cX)+", Y ="+str(cY)) print("to contour point X ="+str(resampledX[ind])+", Y ="+str(resampledY[ind])) print("my_angle_rad["+str(ind)+"] = "+ str(ANG_RAD_TO_SAMPLED_PTS[ind])) input("Press Enter to continue...") # ================================= # Plot an image of the silhouette # surrounded by the contour # ================================= # --> Re-do contour based only on the 100 sampled points modifiedContour = np.dstack((resampledX, resampledY)) # --> Draw the raw contour around a copy of mid_crop_res image imModifiedContour = cv2.drawContours(mid_crop_res.copy(), contoursMCR, -1, (0,0,255), 1) # --> Draw the modified contour based only on the 100 sampled points imModifiedContour = cv2.drawContours(imModifiedContour, modifiedContour, -1, (255,0,0), 1) # --> The plt. instructions are for plotting plt.figure() plt.imshow(imModifiedContour) plt.scatter(resampledX, resampledY, color = "g", edgecolor = 'none') plt.axis("off") # ==//== #plt.title("EUCLIDEAN DISTANCES CURVE GENERATION - Approximated with "+str(numPoints)+" points") # ==//== # --> Define the unique name of the plot image CURR_SAVE_PLOT ="CONTOURS/CONTOUR_"+SPECIMENS[num_spec]+".jpg" # --> Save plotted image plt.savefig(CURR_SAVE_PLOT) # --> Show plot plt.show() plt.close() # --> Wait for keypress input("Press Enter to continue...") # =================================== # HERE WHEN THE LOOP IS OVER! # DEBUG - PRINT ALL THE GLOBAL ARRAYS # =================================== print('ARRAY OF EUCLIDEAN DISTANCES') print(arr_euc_dist) print('ARRAY OF ANGLES TO SAMPLING POINTS') print(ANG_RAD_TO_SAMPLED_PTS) input("Press Enter to continue...") # =============================== # Save arr_euc_dist as a csv file # =============================== np.savetxt('EUC_DIST_CURVES/ANG_RADS_TO_SAMPLED_PTS.csv', ANG_RAD_TO_SAMPLED_PTS, delimiter = ',', fmt='%s') cv2.destroyAllWindows()