Linear fit for a intensity plot (surface plot)
显示 更早的评论
Hi,
I want to find the rotation angle of this plot. (Angle with respect to the horizontal axis.)
In order to do that I'm thinking of doing a linear fit and get the gradient. Something like this:
How can I do a linear fit to this surface plot?
Thank you.
5 个评论
Mathieu NOE
2021-3-8
hello
do you have the data or do you have only this picture ?
S
2021-3-8
Mathieu NOE
2021-3-8
So first attempt here .... but my linear fit seem a bit off
I smoothed first a bit the data and search for peak location at each row.
the red line is the line of maximum values , but visually it seems not truly aligned with our expectation
data = readmatrix('Rotation.txt');
[m,n] = size(data);
% [val,ind] = max(data,[],'all');
dataS = smooth2a(data,1,7);
x=0;
y=0;
k=0;
for ci = 15:40
k = k+1;
[val,ind] = max(dataS(ci,:));
x(k) = ind;
% zz= cumtrapz(dataS(ci,:));
% x(k) = interp1(zz,1:109,max(zz/2));
y(k) = ci;
end
figure(1),imagesc(dataS)
hold on
plot(x,y,'r')
hold off
axis('xy');
% axis('equal');
slope = mean(diff(y))./mean(diff(x));
angl = atan(slope) % in radians
Image Analyst
2021-3-8
Unrecognized function or variable 'smooth2a'.
Error in test5 (line 14)
dataS = smooth2a(data,1,7);
S
2021-3-8
采纳的回答
Image Analyst
2021-3-8
regionprops() fits the image to an ellipse and returns the center and angle in degrees. Try this:
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;
fprintf('Beginning to run %s.m ...\n', mfilename);
fontSize = 15;
data = readmatrix('Rotation.txt');
[rows, columns] = size(data);
subplot(2, 2, 1);
imshow(data, []);
impixelinfo
axis('on', 'xy');
title('Original Data', 'FontSize', fontSize);
subplot(2, 2, 2);
histogram(data);
grid on;
title('Histogram of Data', 'FontSize', fontSize);
mask = data > 30;
% Take largest blob only.
mask = bwareafilt(mask, 1);
subplot(2, 2, 3:4);
imshow(mask);
axis('on', 'xy');
title('Mask', 'FontSize', fontSize);
% Get the angle (Orientation) using regionprops
props = regionprops(mask, 'Orientation', 'Centroid')
angle = -props.Orientation;
xCenter = props.Centroid(1)
yCenter = props.Centroid(2)
fprintf('The angle is %.2f degrees.\n', angle);
x = 1 : columns;
slope = tand(angle)
y = slope * (x - xCenter) + yCenter;
hold on;
plot(xCenter, yCenter, 'r+', 'LineWidth', 2, 'MarkerSize', 30);
plot(x, y, 'r-', 'LineWidth', 2);
caption = sprintf('Angle = %.2f degrees. Slope = %.2f. xCenter = %.1f. yCenter = %.1f', ...
angle, slope, xCenter, yCenter);
title(caption, 'FontSize', fontSize);
g = gcf;
g.WindowState = 'maximized'
props =
struct with fields:
Centroid: [50.5731822474032 28.7478753541076]
Orientation: -17.0388178322049
xCenter =
50.5731822474032
yCenter =
28.7478753541076
The angle is 17.04 degrees.
slope =
0.30647166070947
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Geographic Plots 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
