Extracting spot coordinates from an given line in an image for an almost regular hex pattern of spots but when there is a slight tilt on Y coordinate values.
Show older comments
Hello, I have greyscale image that I am trying to get the positons of the centroids within a given row (i.e the red line below)
The problem is that there is a slight "twist" in my image meaning the y-centroids for a given row actually drift as can be seen below.
The solid CYAN line is a horizontal line which passes thru' the centroid of the most central spot, and you can see at the extreme right and left, the required spots are above or below this line due to the "image twist". But I need to extract the centroids of all the objects on this (tilted) line. The green dashed lines are the centre of the image
This is the most central part, and its easy to extract the centroids of all the objects. The most centralised spot has cooridinates xc,yc (magenta circle). All my other centroids are in vectors x1,y1
But as I only need the centroids on the red line, I filter out all spots above a y-value threshold of say 6 pixels above or below (cyan dashed lines) as I don't want and spots from the row above or below
But its inevitable at the extremes, the spots i want are above or below this threshold and thus not being considered, but instead the spots from one of these upper or lower rows now are counted which I don't want
Im not sure how to proceed, this is my code so far:
ax=app.UIAxes2; IM=getimage(ax); % My function to get the image from the axes [sy,sx]=size(IM); yc=round(sy/2); xc=round(sx/2); % Centre of Image
% Get centroids that have already been found
x1=app.x1; y1=app.y1; % Get saved centroids
% Find most central ones in y, need ROI dimensions
yline(ax,yc,'g--');
xline(ax,xc,'g--');
%Get distance of all x,y's from centre (yc)
n=numel(x1);
data=[];
for i=1:n
dabs = abs(y1(i)-yc); % Y distance of each spot from centre y value
%Below gets the most central spot coordinate
X = [x1(i),y1(i);xc,yc];
dcen = pdist(X,'euclidean');
data(i,1)=i; data(i,2)=x1(i); data(i,3)=y1(i); data(i,4)=dabs; data(i,6)=dcen;
% get intensity at each x,y location (so can remove
% lowest as spurious spots if required
Intensity=IM(round(y1(i)),round(x1(i)));
data(i,5)=Intensity;
end
% Find spot nearest to centre of image (xc,yc)
A = sortrows(data , 6,'ascend'); % sort on distance to centre of image
hold (ax,"on");
XX=A(1,2); YY=A(1,3); % ******* Most Spot Central coordinate ******
yline(ax,YY,'c');
plot(ax,XX,YY,'mo','MarkerSize',20);
disp('distance from cen, all spots')
% This is where I remove all spots with a y value above a distance
% from the most centralised spot
% Remove all rows with y diffences above a ythresh, say 6 pixels
% Effectively ignore all spots except on central line
ythresh=6;
A(A(:,4) > ythresh, :)=[];
plot(ax,A(:,2),A(:,3),'r+','MarkerSize',20);
yline(ax,yc+ythresh,'c--'); yline(ax,yc-ythresh,'c--'); % Show lines constraints on y coordinate
%Remove Spurious spots, i.e. any low intensity ones if any;
Icol=A(:,5); % Get intensity Column
med1=median(Icol); % get median of this column
thr=med1*0.3; % define a threshold (0.3) of median value to remove
A(A(:,5) < thr, :)=[]; % Remove weak spots
A= sortrows(A , 5,"ascend");
% MY AIM IS TO HAVE A VECTOR OF X VALUES OF THE RELEVANT SPOTS ALL
% REFERENCE TO THE MOST CENTRAL SPOT.
% Extract X coords and extract off 1st (=central location spot) X value
% so all relative to centralised spot
C = sortrows(A , 6,'ascend'); % sort on dcen
xvalues=C(:,2);
disp('Central Spot X')
XX
xvalueCen=C(1,2)
xref=xvalues-xvalueCen;
%This is where Im running into difficulty as Im picking up the
%spot centroids from the wrong row due to the tilt
1 Comment
Mathieu NOE
on 10 Jul 2025
hello
you should share also the data if you want someone to help you
Accepted Answer
hello again
maybe you could implement this method :
%% create some dummy pattern of centroids
% center of image (central spot coordinates) is 0,0
xc = 0;
yc = 0;
N = 10;
dx = 1;
x = -N:dx:N;
for k=1:7
y = (k-4)*ones(size(x)) - 0.5*(k)*(cos(pi*x/N)-1);
y_all(k,:) = y; % storage
end
%% selection based on threshold (as you did) - it will fail
threshold = 0.5;
ydist = abs(y_all - yc);
ind = (ydist<threshold);
x_select = x(any(ind,1));
y_select = y_all(ind);
figure
hold on
plot(x,y_all,'*')
plot(x_select,y_select,'dk','markersize',15)
hold off
%% "smart" selection based on closest neighbour search
nx = numel(x);
threshold = 0.5;
indxc = find(x==xc);
% init
y_select2 = zeros(1,nx);
x_select2(indxc) = xc;
y_select2(indxc) = yc;
% right half closest neighbour search (for x > xc)
for ii=indxc+1:nx
ydist = abs(y_all(:,ii) - y_select2(ii-1));
[val,ind] = min(ydist);
y_select2(ii) = y_all(ind,ii);
end
% left half closest neighbour search (for x < xc)
for ii=indxc-1:-1:1
ydist = abs(y_all(:,ii) - y_select2(ii+1));
[val,ind] = min(ydist);
y_select2(ii) = y_all(ind,ii);
end
figure
hold on
plot(x,y_all,'*')
plot(x,y_select2,'dk','markersize',15)
hold off
15 Comments
For some reason I am getting an error
Here is smaller Image Im just trying to apply your code to
and the x and y coordinates of the centroids are in the attached mat file (variables are x1 & y1)
Index in position 2 exceeds array bounds. Index must not exceed 1.
Error in HTS_TestSoftware/DistortionButtonPushed (line 15628)
ydist = abs(y_all(:,ii) - y_select2(ii-1));
Note, the most central spot centroid isn't always exactly in the centre of the image
I also asigned the following
x=x1;
y_all=y1;
Correction, I had to transpose the x1 and y1 like this
x=x1';
y_all=y1';
and the error went away, but its picking out all of the spots, not just the central Y spots...
I don't see threshold being used in your "smart selection" code
%% "smart" selection based on closest neighbour search
nx = numel(x);
threshold = 0.5;
indxc = find(x==xc);
Ive been looking at this all day to work out why my data doesn't work (like your simulated data). Could it because my central spot coordinate isn't at (0,0)?
If I look at just one half of the data, and draw some lines to show me where its found the next spot based on min Y distance,
plot(ax,xc,yc,'c+','Markersize',20,'LineWidth',4)
% right half closest neighbour search (for x > xc)
for ii=indxc+1:nx
ydist = abs(y_all(:,ii) - y_select2(ii-1));
[val,ind] = min(ydist)
y_select2(ii) = y_all(ind,ii); %ii is the index of the next spot required
xline(ax,x(ii),'r--'); pause(1);
But as you can see its including the next spot that shouldn't be included (shown by the red arrows). Furthermore, it actually counts these twice (you can see the red dashed line is thicker), presumably because there are two spots equidistance (i.e. above the line and below the line)
1: The green dashed lines are the centre of the image
2: Cyan Cross is central spot identified
3: Cyan solid line is a horizontal lne thru the y-centroid coordinate of the central spot
4: Cyan dashed lines just represent my ythresh boundaries
hello again
sorry for my late answer , was busy those last days ...
the last data you provided show a nice horizontal layout without too much distorsion or tilt
here the simpler approach with threshold is working fine .
as the vertical separation is about 10 units , so I picked a threshold around 3 .
Of course , with large distorsion or tilt you may have to switch to my suggestion (and BTW, there is no threshold there, sorry for the confusion)
the second code must be adapted according how the data is organized...
load('coords.mat')
x=x1';
y_all=y1';
% central spot coordinates is xc,yc
d = (x-mean(x)).^2 + (y_all-mean(y_all)).^2;
[~,indc] = min(d);
xc = x(indc);
yc = y_all(indc);
%% selection based on threshold (as you did)
threshold = 3;
ydist = abs(y_all - yc);
ind = (ydist<threshold);
x_select = x(any(ind,1));
y_select = y_all(ind);
figure
hold on
plot(x,y_all,'*')
plot(x_select,y_select,'dr','markersize',10)
plot(xc,yc,'dk','markersize',15)
hold off
Mathieu NOE
on 16 Jul 2025
are the coordinates data the same as in your initial question ? seems as your code should have worked fine from the beginning...
Jason
on 16 Jul 2025
hello again
seems that it still work even with the larger data set , butit's probably the max tilt we can afford with this approach
the threshold is now computed based on what is the height gap between the lines
I'll continue to find out what could be a better approach for higher tilt situations - do you have data with more tilt ?
load('coords2.mat')
x=x1';
y_all=y1';
% central spot coordinates is xc,yc
d = (x-mean(x)).^2 + (y_all-mean(y_all)).^2;
[~,indc] = min(d);
xc = x(indc);
yc = y_all(indc);
%% selection based on threshold (as you did)
mdy = min(abs(diff(y_all)));
threshold = ceil(mdy/2);
ydist = abs(y_all - yc);
ind = (ydist<threshold);
x_select = x(any(ind,1));
y_select = y_all(ind);
figure
hold on
plot(x,y_all,'.')
plot(x_select,y_select,'*r')
plot(xc,yc,'dg','markersize',15)
hold off
with more tilt , you would probably use the next suggestion (below)
as the lines are approx straight and parallel , the simplest trick I can think of is to rotate the entire data set based on a first small batch of selected points around the center
so start with a smaller threshold (like 1), use the small batch of selected points to compute the required angle of rotation (around the center point) to have the lines as horizontal as can be and repeat the selection process , still using the simple threshold approach that now works on the entire width of data as they are (almost) horizontal.
the index of the rotated points is the same as the raw data points so you don't even need to do a backward rotation (is also doable)
load('coords2.mat')
x=x1';
y_all=y1';
% central spot coordinates is xc,yc
d = (x-mean(x)).^2 + (y_all-mean(y_all)).^2;
[~,indc] = min(d);
xc = x(indc);
yc = y_all(indc);
%% "smart" selection
threshold = 1;
ydist = abs(y_all - yc);
ind = (ydist<threshold);
x_select = x(ind);
y_select = y_all(ind);
slope_angle = atan(mean(diff(y_select)./diff(x_select))); % in rad
figure
hold on
plot(x,y_all,'.')
plot(x_select,y_select,'*r')
plot(xc,yc,'dg','markersize',15)
hold off
% Define the points as a 2xN matrix (each column is a point)
points = [x-xc; y_all-yc]; % NB : remove center point coordinates
% Create the rotation matrix
theta = -slope_angle;
R = [cos(theta), -sin(theta); sin(theta), cos(theta)];
% Rotate the points
rotated_points = R * points;
xr = rotated_points(1,:) + xc;
yr = rotated_points(2,:) + yc;
% select rotated points
threshold = 1;
ydist = abs(yr - yc);
ind = (ydist<threshold);
xr_select = xr(ind);
yr_select = yr(ind);
% plot rotated points selection
figure
hold on
plot(xr,yr,'.')
plot(xr_select,yr_select,'*r')
plot(xc,yc,'dg','markersize',15)
hold off
% original data selection
x_select = x(ind);
y_select = y_all(ind);
figure
hold on
plot(x,y_all,'.')
plot(x_select,y_select,'*r')
plot(xc,yc,'dg','markersize',15)
hold off
hello again
just for my own fun, I wanted to try another approach , that could be of some help if you data would shows some zig zags , or more complex distorsion so that a simple rotation as shown before would probably not suffice.
here I tried to find the closest neighboors (starting at the center) of the right hand side and left hand side
a threshold criteria is needed to find jumps in the "rearranged" data (see function provided below in the code) so that we keep only the data from the line where the center belongs.
all the best
load('coords2.mat')
x=x1';
y_all=y1';
% central spot coordinates is xc,yc
d = (x-mean(x)).^2 + (y_all-mean(y_all)).^2;
[~,indc] = min(d);
xc = x(indc);
yc = y_all(indc);
%% "smart" selection
p = reorder_data([x' y_all'],indc); % create a "reordered" line from scatter points - find closest distance neighbour)
x_select1 = p(:,1);
y_select1 = p(:,2);
% find trailing segment - "right hand side" of center point
threshold = 0.5;
id1 = find(diff(y_select1)>threshold,1,'first');
x_select1a = x_select1(1:id1);
y_select1a = y_select1(1:id1);
% find leading segment ("left hand side" from center point)
[m,indm] = min(diff(y_select1));
id3 = find(diff(y_select1(indm+1:end))<-threshold,1,'first');
id3 = id3 + indm;
x_select1b = x_select1(indm+1:id3);
y_select1b = y_select1(indm+1:id3);
% finally combine segments , remove duplicates (may be some) and sort in x ascending order
xfinal = [xc;x_select1a;x_select1b];
yfinal = [yc;y_select1a;y_select1b];
[xfinal,ia,ib] = unique(xfinal);
yfinal = yfinal(ia);
figure
hold on
plot(x,y_all,'.')
% plot(x_select1a,y_select1a,'r')
% plot(x_select1b,y_select1b,'m')
plot(xfinal,yfinal,'*r')
plot(xc,yc,'dg','markersize',15)
hold off
%%%%%%%%%%%%%%%%%%%%%%%%%
function p = reorder_data(pts,startIdx)
% trying to put back data in correct order : Greedy Nearest Neighbor Loop
% For each next step, compute distances from the current point to all unvisited points, pick the smallest, and mark it:
% startIdx = 1;
% Create a logical array to flag visited points and an order array:
N = size(pts,1);
visited = false(N,1);
order = zeros(N,1);
order(1) = startIdx;
visited(startIdx) = true;
for k = 2:N
% % original code for x,y,z coordinates points
% currentPt = pts(order(k-1), :); % last visited
% diffs = pts - currentPt(ones(N,1), :); % vector differences
% dists = sqrt(sum(diffs.^2, 2)); % Euclidean distances
% dists(visited) = inf; % ignore already visited
% [~, idx] = min(dists); % nearest
% order(k) = idx;
% visited(idx) = true;
% modified code for y distance only
currentPt = pts(order(k-1), :); % last visited
diffs = pts(:,2) - currentPt(ones(N,1), 2); % y differences
dists = sqrt(sum(diffs.^2, 2)); % Euclidean distances
dists(visited) = inf; % ignore already visited
[~, idx] = min(dists); % nearest
order(k) = idx;
visited(idx) = true;
end
p = pts(order, :);
end
Jason
on 17 Jul 2025
Mathieu NOE
on 17 Jul 2025
Edited: Mathieu NOE
on 17 Jul 2025
hello again
maybe we are not exactly talking about the same thing
in y code xc and yc is supposed to be the most central point in your coordinates data but I don't have any idea where is the center of the image as I don't have it's dimensions
so that's why your definition of the central spot may differ from mine
Jason
on 17 Jul 2025
Mathieu NOE
on 18 Jul 2025
we are trying to find the point which is closest to the centroid of your data (centroid = spatial mean or center of gravity of your data)
the centroid is given by (mean(x),mean(y_all)) , assuming here all dots have same weight (importance)
d = (x-mean(x)).^2 + (y_all-mean(y_all)).^2; is simply the euclidian distance of the dots to this centroid
the index corresponding of the min of d gives us the dot (coordinates) which is the most "central" in your data
Jason
on 18 Jul 2025
Mathieu NOE
on 21 Jul 2025
as always, my pleasure !
More Answers (0)
Categories
Find more on Region and Image Properties in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
