Azimuthally average line plot to 2D image
14 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Arthur Moya
le 23 Mar 2021
Commenté : Arthur Moya
le 24 Mar 2021
Hi there,
I am trying to create a 2D image from a line profile.
Effectively, I would like to obtain its radial average and it must be centred in the middle of the image.
I have attached the x and y values of the plot.
I don't know where to begin.
Your help in anyway is appreciated,
Arthur
2 commentaires
Mathieu NOE
le 23 Mar 2021
hello
maybe I'm wrong but is the question about doin a kind of circle fit to your data (finding radial average ? )
Réponse acceptée
Mathieu NOE
le 23 Mar 2021
Modifié(e) : Image Analyst
le 23 Mar 2021
Hello again
Check this if it works for you.
I had to divide y by 1000 so that x and y are in similar ranges (is one in mm and the other one in microns ?) so circfit would work.
Code below + function.
All the best.
a = load('Chip_NPS_fit_x_y.txt');
xs = a(:,1);
ys = a(:,2)/1000;
ind = find(xs>0.15 & xs < 0.2);
[xfit,yfit,Rfit,a] = circfit(xs(ind),ys(ind));
figure(1)
plot(xs,ys,'b',xs(ind),ys(ind),'g')
hold on
rectangle('position',[xfit-Rfit,yfit-Rfit,Rfit*2,Rfit*2],...
'curvature',[1,1],'linestyle','-','edgecolor','r');
title(sprintf('Best fit: R = %0.3f; Ctr = (%0.3f,%0.3f)',...
Rfit,xfit,yfit));
plot(xfit,yfit,'g.')
axis equal
hold off
%%%%%%%%%%%%%%%%%%%%%
function [xc,yc,R,a] = circfit(x,y)
%
% [xc yx R] = circfit(x,y)
%
% fits a circle in x,y plane in a more accurate
% (less prone to ill condition )
% procedure than circfit2 but using more memory
% x,y are column vector where (x(i),y(i)) is a measured point
%
% result is center point (yc,xc) and radius R
% an optional output is the vector of coeficient a
% describing the circle's equation
%
% x^2+y^2+a(1)*x+a(2)*y+a(3)=0
%
% By: Izhak bucher 25/oct /1991,
x=x(:); y=y(:);
a=[x y ones(size(x))]\[-(x.^2+y.^2)];
xc = -.5*a(1);
yc = -.5*a(2);
R = sqrt((a(1)^2+a(2)^2)/4-a(3));
end
6 commentaires
Mathieu NOE
le 24 Mar 2021
hello
this is it
the output is a matrix of size 2048 x 2048 stored in zq
a = load('Chip_NPS_fit_x_y.txt');
r = a(:,1);
zs = movmean(a(:,2),25); % a bit of smoothing (optional)
figure(1);
plot(r,zs);
% create the mexican hat 3D plot
theta=linspace(0,2*pi,length(r));
[x,y]=pol2cart(theta, r);
zsr=zs*ones(1,length(r));
figure(2);
s = surf(x,y,zsr);
s.EdgeColor = 'none';
view([0 90]);
axis('equal');
colormap('gray');
colorbar('vert');
% interpolation on a 2048 x 2048 square grid
% square limit = max(r) / sqrt(2);
ll = max(r) / sqrt(2);
m = linspace(-ll,ll,2048);
[Xq,Yq] = meshgrid(m,m);
zq = griddata(x,y,zsr,Xq,Yq); %the output is a matrix of size 2048 x 2048 stored in zq
figure(3);
s = surf(Xq,Yq,zq);
s.EdgeColor = 'none';
view([0 90]);
axis('equal');
colormap('gray');
colorbar('vert');
xlim([-ll ll]);
ylim([-ll ll]);
Plus de réponses (1)
Image Analyst
le 24 Mar 2021
Here is a simple brute force for-loop way. Simply create an image then scan every pixel location determining the radius from the center for that pixel and assigning the signal value there.
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 short g;
format compact;
fontSize = 18;
fprintf('Beginning to run %s.m ...\n', mfilename);
data = readmatrix('Chip_NPS_fit_x_y.txt');
r = data(:, 1);
signal = data(:, 2);
nexttile
plot(r, signal, 'b-');
title('Radial Profile', 'FontSize', fontSize);
grid on;
drawnow;
numElements = length(signal)
rows = ceil(2 * numElements / sqrt(2))
middlex = rows / 2;
middley = rows / 2;
grayImage = zeros(rows, rows);
for col = 1 : rows
fprintf('Assigning column %d.\n', col);
for row = 1 : rows
% Get the distance from the center of the image.
radiusInPixels = sqrt((col - middlex)^2 + (row - middley)^2);
% Scale it to be in the range 0 - 0.5, like r is.
% Find the signal value there.
radius = r(end) * radiusInPixels / numElements;
[~, index] = min(abs(r - radius));
grayImage(row, col) = signal(index);
end
end
nexttile
imshow(grayImage, []);
title('As a 2-D image', 'FontSize', fontSize);
% Maximize window
g = gcf;
g.WindowState = 'maximized';
Voir également
Catégories
En savoir plus sur Lighting, Transparency, and Shading dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!