Radial average of 2D matrix

Dear Community

I have a filter with a 256 by 256 double dimension. I would like to take the radial average (starting at the center of the image, 1 pixel at a time) in order to convert this into a 1D filter to apply to my data.

I know that

ceil((n+1)/2)

directs me to the center of my matrix, but how would I take radial averages outward from this point?

Thank you so much!

[Sergey Loginov]

For 256x256 it should not matter much, but consider using accumarray based solutions!

function [Tics,Average]=radial_profile(data,radial_step)

%main axii cpecified:

x=(1:size(data,2))-size(data,2)/2;

y=(1:size(data,1))-size(data,1)/2;

% coordinate grid:

[X,Y]=meshgrid(x,y);

% creating circular layers

Z_integer=round(abs(X+1i*Y)/radial_step)+1;

% % illustrating the principle:

% % figure;imagesc(Z_integer.*data)

% very fast MatLab calculations:

Tics=accumarray(Z_integer(:),abs(X(:)+1i*Y(:)),[],@mean);

Average=accumarray(Z_integer(:),data(:),[],@mean);

end

Sergey Loginov (2021). radial_profile (https://www.mathworks.com/matlabcentral/fileexchange/101480-radial_profile), MATLAB Central File Exchange. Retrieved November 3, 2021.

Answer 1

[Tushar Sinha]

Hi Brian

I understand that you wish to know how to compute the radial average for a 2D matrix. Please refer to the MATLAB script in the File Exchange submission link below which calculates the radial average of a 2D matrix:

Radial Average

I hope this helps resolve the issue.

Thanks,

Tushar

[Brian]

This is fantastic! Thank you so much, Tushar.

[Sergey Loginov]

Dear Tushar,

There is a much faster approach based on the matrix features of MatLab itself. I have used accumarray function and get a very substantial speed-up compared to the for loops!

Sergey Loginov (2021). radial_profile (https://www.mathworks.com/matlabcentral/fileexchange/101480-radial_profile), MATLAB Central File Exchange. Retrieved November 3, 2021.

Answer 2

[Sergey Loginov]

This is very sad! For many years people ask this question and get the same for loop answer. The code below works 50 times faster at size data 2048x1792 .

function [Tics,Average]=radial_profile(data,radial_step)

%main axii cpecified:

x=(1:size(data,2))-size(data,2)/2;

y=(1:size(data,1))-size(data,1)/2;

% coordinate grid:

[X,Y]=meshgrid(x,y);

% creating circular layers

Z_integer=round(abs(X+1i*Y)/radial_step)+1;

% % illustrating the principle:

% % figure;imagesc(Z_integer.*data)

% very fast MatLab calculations:

Tics=accumarray(Z_integer(:),abs(X(:)+1i*Y(:)),[],@mean);

Average=accumarray(Z_integer(:),data(:),[],@mean);

end

[tamp]

Hi Sergey,

How could I extend this to the entire diameter, instead of just radial? I'd like a full cross-section profile.

Thank you,

Tamas

[Sergey Loginov]

Hi Tamas,

Sorry for the late reply.

Unfortunately, I do not quite get your question. Have checked Radon transform already realized in MatLab?

Regards,

Sergey

Answer 3

[Hugo Trentesaux]

You can program a function like this :

function profile = radialAverage(IMG, cx, cy, w)
    % computes the radial average of the image IMG around the cx,cy point
    % w is the vector of radii starting from zero
    [a,b] = size(IMG);
    [X, Y] = meshgrid( (1:a)-cx, (1:b)-cy);
    R = sqrt(X.^2 + Y.^2);
    profile = [];
    for i = w % radius of the circle
        mask = (i-1<R & R<i+1); % smooth 1 px around the radius
        values = (1-abs(R(mask)-i)) .* double(IMG(mask)); % smooth based on distance to ring
        % values = IMG(mask); % without smooth
        profile(end+1) = mean( values(:) );
    end
end

[Nathaniel Conrad]

Hi, as a note, this code will only work if the image of interest is centered properly. i.e., if cx and cy are not exactly the center of the picture, then the mask will not index the IMG matrix properly. Do you have a solution to this problem?

For example, say you have a 20x20 image, but you want to radial average about the (5,5) point in the image. Choosing w = 20, cx = 5, and cy = 5 leads to a mismatch in indexing when applying IMG(mask). In particular, looking at just i =1 in the for loop (i.e., < 2 pixel out in radii from the cx and cy), IMG(mask) will give you intensity values from the following points: (14,2), (15,2), (16,2), (4,3), (4,6), (14,3), (15,3), (16,3). It should be giving you intensity values from (4,4), (5,4), (6,4), (4,5), (6,4), (4,6), (5,6), (6,6).

Anyone know why that is happening and how to fix it?

[Hugo Trentesaux]

You mean subpixel position of the center? In this case, you can do a subpixel translation of the mask with the imtranslate function.

Answer 4

[David Fischer]

https://www.mathworks.com/matlabcentral/fileexchange/46468-radialavg-zip

I updated my radial average function from the recent recommendations to use accumarray. It's faster now, but I'm not tested it extensively to see how much faster. There is still some improvement to be done I'm sure in smarter handing of the matrix data. And it no longer safely handles NaN in the input data set. This thread referenced in update in acknowledgement.