Calculating effective area in Soccer

I am currently trying to interpret how to code the attached algorithm (pictured).

Beginning information states "In order to create a polygon on the planar dimension, at least three points are necessary (i.e., triangle). Therefore, three players need to be considered to build triangles as the combinations of N players, in which N is the total number of players within a team."

The algorithm is then given.

I am wondering how this translates into MATLAB.

Answer 1

[David Goodmanson]

Hi William,

The algorithm draws a lot of triangles with players at the vertices. I appears that the final result for P is the polygon enclosed by the red boundary in the code below (the convex hull). After that, what the last line is doing is not so clear. Are you just looking for the area inside the boundary? If it's no more complicated than that, then the following code should get it done.

Let the players be at locations (x,y) and let vx be the vector of x coordinates and vy be the vector of corresponding y coordinates. Then

% make up some data
vx = 100*rand(1,11);
vy = 50*rand(1,11);

ind = convhull(vx,vy)             % indices of players on the boundary
plot(vx,vy,'o',vx(ind),vy(ind));
xlim([0 100]);
ylim([0 50]);
A = polyarea(vx(ind),vy(ind))

[William Sheehan]

Thanks David, that looks awesome!

Is there away to formulate the code in a way that measures the change in area over time, i.e. on a second to second basis over a 30 second window?

I would like to apply this code to different 30 second periods throughout the game using (x,y) coordinate data.

[David Goodmanson]

Hi William,

that would not be hard to do, you just need vx and vy for each of those one-second intervals. How would you be storing/accessing the vx and vy data?

Answer 2

[William Sheehan]

Data is collected and exported to excel with a separate spreadsheet for each individual (as pictured). Only values that will be relevant will be the latitude (vx) and longitude (vy) columns. Is it best to create two new spreadsheets - one with all individuals vx values over time side by side and another separate one with all the vy values?

[David Goodmanson]

Hi William,

Here's an example. First, in the answer I posted, vx and vy sounded too much like velocity, so in the example below px and py are now position.

Suppose that with some combination of Excel and Matlab you concatenate side-by-side the x position columns for the players into a matrix Mpx, and the same for y into a matrix Mpy. So players are columns of the matrix, and position varies down the columns. (It also works well to concatenate the pair of columns F&G into a matrix and separate out x and y later). Suppose also there is a column vector t of times which is Excel column B, and I assume that B is the same for all players. The first part of the code below sets all that up (with more active players than can be on the field at one time, just for demo purposes).

In the second part 11 players are selected to plot, as well as the start time and the duration.

I used pause(.2) to slow the plots down enough that I could see them. Pause time will be machine dependent.

What is not included is that you show lat/lon rather than position in meters. Maybe you have other Excel columns for that. Anyway it's easy to import lat/lon into Matlab and do the conversion to meters, assuming lat/lon is provided for one corner of the field.

% -----------------------------
% make up some data. 20 players, 48 minutes, field size 100 x 50.
n = 20;      
p0x = 100*rand(1,n);  p0y = 50*rand(1,n);
v0x = 10*rand(1,n)-5;  v0y = 10*rand(1,n)-5;
t = (0:.1:2880)';  onz = ones(size(t));
Mpx = onz*p0x + t*v0x;  Mpx = mod(Mpx,200);  Mpx(Mpx>100) = 200 - Mpx(Mpx>100);
Mpy = onz*p0y + t*v0y;  Mpy = mod(Mpy,100);  Mpy(Mpy>50) = 100 - Mpy(Mpy>50);
% ----------------------------------

% select players and times
players = [1 2 3 4 6 7 9 10 14 15 17];
tstart = 1600;    % seconds
tduration = 30;

ind0 = max(find(t<= tstart));
ind = (ind0:10:ind0 + tduration*10);
Area = zeros(size(ind));
count = 0;
for k = ind
  px = Mpx(k,players);
  py = Mpy(k,players);
  istart = convhull(px,py);     % indices of players on the boundary
  plot(px,py,'o',px(istart),py(istart))
  hold on
  plot([0 0 100 100 0],[0 50 50 0 0],'k-');  % field boundary
  xlim([-10 110]);  ylim([-10 80]);
  hold off
  pause(.2)  
  count = count+1;
  Area(count) = polyarea(px(istart),py(istart));
end
Area'

[William Sheehan]

This is awesome David!

I'm just trying to get my head around the variables you have included. Am I correct in saying that 'ind0' finds the start point of the specified point selected at 'tstart' and the 'ind' variable takes every tenth value and thereby represents the change in convex hull/area at every second?

Thank you again this is perfect

[David Goodmanson]

That's correct. It would have been possible to use ind0 = find(t==tstart) except that the times are not strictly in .1 second intervals and sometimes have an extra .01 sec as in 3.01 etc.

This is cool stuff but kind of scary in a way. Back in a different era it would just be a coach on the sidelines yelling "Get your butt in gear, Best, you're doggin' it!" Now it's a guy telling a guy that GPS says he only ran 6.37 miles for the game and what's up with that.

[William Sheehan]

Hi David, thank you for your feedback and answers as they have proven to be of great use! Touching on something you mentioned earlier, to convert these variables to metres, how would I go about doing that?

[William Sheehan]

Is there a way to calculate convhull without removing NaNs from the data?

[David Goodmanson]

Hi William,

First, it's easier than I thought to just use px and py as lon and lat values in degrees, respectively. After calculating the area in square degrees, you can convert to square meters plenty accurately enough by multiplying by

(pi/180)^2*(6370e3)^2*cosd(-33.7698)
ans =  1.0275e+10

Here is a function that does convhull with nans. 'Nan players' are ignored, and in case you want to keep track of the players at the corners of the convex hull, for each player the convhullnan index stays the same as it was for the full list in convhull.

I had not realized before that the second output of convhull supplies the area so there is no need to use polyarea.

% example (in no particular units)

px = rand(1,20);
py = rand(1,20);
px([10 2 6]) = nan;
py([7 2 11]) = nan;
[ind A] = convhullnan(px,py)

function [ind A] = convhullnan(px,py)

% convex hull when position row vectors px,py might contain nans 

n = (1:length(px))';
indnan = find(any(isnan([px;py])));
px(indnan) = [];
py(indnan) = [];
n(indnan) = [];
[ind1 A] = convhull(px,py);
ind = n(ind1);
end

[William Sheehan]

Hi David,

This appears to work well! Instead of using the area from convhull and the conversion you provided, I am using areaint function with "earthradius" which provides the area in metres squared. Both the results were very similar. Which would be more reliable ?

[David Goodmanson]

Hi William, probably areaint, since it is part of a Matlab toolbox. According to wikipedia, the earth's radius varies from 6353 km to 6384 km depending on location. For all I know, areaint has a lot of data and uses the radius appropriate to Sydney. earthradius gives 6371. I used 6370 so I would be interested to know how the answers compare if you changed the 6370 I used to 6371.