Best Way to find density of [x,y] coordinates within set window
5 views (last 30 days)
Show older comments
I am using Matlab 2015. What is the best way to code a set window and have it move through some [x,y] coordinate data (from start to end), count the plots inside the window for the current iteration and save those values in a seperate array? The size of the window can remain the same as it moves throughout, and each window can be centered around the coordinate point.
I am doing this to detect clusters of coordinates. Below is an illustration which hopefully aids in explaining what I am after.
8 Comments
Adam Danz
on 16 Aug 2019
To add to Walter's comment above, if the point is to detect clusters, histcounts2() or histrogram2() might be more suitable (as I mentioned in your first question).
Walter Roberson
on 16 Aug 2019
To handle the moving window, you can do histcounts2() using bins that are the width of a constant increment in x and y, after which you would do conv2() to calculate the total number of entries that are in each larger box -- like a moving sum but in 2D. This can easily be vectorized -- but it does depend upon the size of the bins being known ahead of time.
Accepted Answer
Adam Danz
on 16 Aug 2019
Edited: Adam Danz
on 18 Aug 2020
Here's an alternative that uses inpolygon() to determine whether coordinates are within a defined rectangular window centered at each of your coordinates. The size of the window can be adjusted in the variable "windowSize".
'inCount' is the number of points within the window defined by windowSize, centered on each coodinate.
'onCount' is the number of points that fall exactly on the window (which will be rare).
It would be very easy to also store the coordinates of all points that fall into each window if you needd that information, too.
load('xy_coordinates.mat'); % better to use a full path
windowSize = [20,20]; % window size, [width,height];
inCount = nan(size(xy_coordinates,1),1);
ouCount = inCount;
onCount = 0;
% Loop through each data point
for i = 1:size(xy_coordinates,1)
% Determine square window coordinates, centered at i-th datapoint
windCoord = xy_coordinates(i,:) + windowSize/2.*[1,1;-1,-1]; %[UpperLeft, UpperRight; LowerLeft, LowerRight] coordinates
% Alternative: windCoord = xy_coordinates([i,i],:) + [windowSize/2;-windowSize/2]
% Determine how many points are within window
[in,on] = inpolygon(xy_coordinates(:,1),xy_coordinates(:,2),windCoord([2,1,1,2,2]),windCoord([4,4,3,3,4]));
% to plot the window : hold on; plot(windCoord([2,1,1,2,2]),windCoord([4,4,3,3,4]), 'r-')
% Count the number of points IN and ON the window
inCount(i) = sum(in);
onCount(i) = sum(on);
end
With animation....
In this version below, you can watch the animated window move through your data and there is text in the upper left corner of the plot showing the point counts within the window.
load('xy_coordinates.mat'); % better to use a full path
figure();
plot(xy_coordinates(:,1), xy_coordinates(:,2), '-ob')
axis equal
hold on
windowSize = [20,20]; % window size, [width,height];
inCount = nan(size(xy_coordinates,1),1);
ouCount = inCount;
onCount = 0;
wh = gobjects(1);
th = text(min(xlim()),max(ylim()),'Count: 0', 'VerticalAlignment', 'top','FontSize',14);
% Loop through each data point
for i = 1:size(xy_coordinates,1)
% Determine square window coordinates, centered at i-th datapoint
windCoord = xy_coordinates(i,:) + windowSize/2.*[1,1;-1,-1]; %[UpperLeft, UpperRight; LowerLeft, LowerRight] coordinates
% Alternative: windCoord = xy_coordinates([i,i],:) + [windowSize/2;-windowSize/2]
% Determine how many points are within window
[in,on] = inpolygon(xy_coordinates(:,1),xy_coordinates(:,2),windCoord([2,1,1,2,2]),windCoord([4,4,3,3,4]));
% Count the number of points IN and ON the window
inCount(i) = sum(in);
onCount(i) = sum(on);
% Draw window and label number of memebers
delete(wh)
wh = plot(windCoord([2,1,1,2,2]),windCoord([4,4,3,3,4]), 'r-');
th.String = sprintf('Count: %d',inCount(i));
drawnow
pause(.1)
end
4 Comments
More Answers (2)
Image Analyst
on 16 Aug 2019
Edited: Image Analyst
on 16 Aug 2019
James:
Attached is how you do it using a rectangular window. Choose the window width to be square if you want. It counts the number of other data points ANYWHERE in the data. If you want it just considering the past N elements (like recent history), you'd need to change the deltax and deltay computations to just include x and y from k-N to k ONLY instead of the whole vector like this: x(max([1, k-N]) : k). Same for y. Let me know how it goes.
clc; % Clear the command window.
close all; %Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;
% Get original data.
s = load('xy_coordinates.mat')
x = s.xy_coordinates(:, 1);
y = s.xy_coordinates(:, 2);
% Plot original data.
subplot(2, 1, 1);
plot(x, y, 'r.-', 'LineWidth', 2, 'MarkerSize', 20);
grid on;
xlabel('x', 'FontSize', fontSize);
ylabel('y', 'FontSize', fontSize);
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0, 0.04, 1, 0.96]);
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst for James Morris', 'NumberTitle', 'Off')
% Make both x and y axes have the same scale so it's not distorted and
% we get a true idea of what the window looks like.
axis equal;
xticks(0:50:1000);
% Now find out how many are within a sliding window centered about each data point.
windowWidthX = 30; % Adjust as needed.
windowWidthY = 50; % Adjust as needed.
% Plot the window over the first point so we can see its size.
hold on;
rectangle('Position', [x(1)-windowWidthX/2, y(1)-windowWidthY/2, windowWidthX, windowWidthY]);
caption = sprintf('Box size: %.1f tall by %.1f wide. One box is shown over the very first point.', windowWidthY, windowWidthX);
title(caption, 'FontSize', fontSize);
withinWindow = zeros(1, length(x)); % Initialize for speed.
for k = 1 : length(x)
thisX = x(k);
thisY = y(k);
% Compute distance from ALL other data points to this data point
% (not just to those within a certain distance in the past along the trace).
deltax = abs(thisX - x);
deltay = abs(thisY - y);
indexesInsideWindow = abs(deltax <= windowWidthX/2) & abs(deltay < windowWidthY/2);
withinWindow(k) = sum(indexesInsideWindow);
end
% Plot point density as a function of index number.
subplot(2, 1, 2);
plot(1:length(withinWindow), withinWindow, 'b.-', 'LineWidth', 2, 'MarkerSize', 20);
grid on;
xlabel('Element number (index)', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
title('Count within window.', 'FontSize', fontSize);
4 Comments
Image Analyst
on 17 Aug 2019
Yes, that's a newer function and it's not crucial. You can take it out and it will just automatically decide on the tick mark values.
Since it worked perfectly for you, maybe you can "Accept this Answer" and vote for any other answers that also worked (if you tried them).
Image Analyst
on 16 Aug 2019
Can you supply your x and y in a .mat file
save('answers.mat', 'x', 'y');
Use the paper clip icon.
I'll probably just scan along and for each point, find the distance to each point, then count the number less than some window (radius) that you want to consider.
for k = 1 : length(x)
thisX = x(k);
thisY = y(k);
allDistances = sqrt((thisX - x).^2 + (thisY-y).^2);
withinWindow(k) = sum(allDistances < someValue);
end
withinWindow is a vector that will contain the count of all other points that are closer than "someValue" to the point at that index.
2 Comments
See Also
Categories
Find more on Startup and Shutdown 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!