A=imread( '1.jpg');                       % Original Image
A1=rgb2gray(A);                           % Convert to grayscale
ABi=imbinarize(A1);                       % Convert the grayscale image to binary
AF=imfill(ABi,'holes');                   % Fill the holes of the binary image
figure(1);
imshow(AF)
AC=imcomplement(ABi);                     % Invert the filled image
C=intersect(AF,AC);                       % Intersection of the filled image and the inverted image
intersection = (double(AF) - double(AC)) == 0;
figure(2)
imshow(intersection);
CC = bwconncomp(intersection);            % Connected components of the intersect image

stat=regionprops(CC, 'Area','Centroid');  % Calculate properties of the intersect image
area = cat(1, stat.Area);
area_sorted=sort(area,'descend');         % sort areas of the connected components in descending order
centroids = cat(1, stat.Centroid);
centroids_sorted=sort(centroids,'descend');% sort centroids of the connected components in descending order
desired_points= centroids_sorted(1:3,:);