# Extracting profiles from a point cloud surface along the mazimum dip of the plane

44 Ansichten (letzte 30 Tage)
Elisa am 19 Jul. 2024 um 9:40
Kommentiert: Elisa am 29 Jul. 2024 um 13:59
Dear all,
Based on the result of this code (blue cloud rotated), I need to extract a series of profiles (with a fixed diastance) to get the Z values of all the points constituting the profile. These profiles should be extracted from vertical planes oriented according to the maximum slope of the cloud.
you have other ideas, they are accepted!! I have included a sketch of what I would like.
clear all
close all
clc
X = XYZ(:,1);
Y = XYZ(:,2);
Z = XYZ(:,3);
xyz0=mean(XYZ,1);
A=XYZ-xyz0; % center the data at zero
% Find the direction of most variance using SVD and rotate the data to make
% that the x-axis
[~,~,V]=svd(A,0);
a=cross(V(:,3),[0;0;1]);
T=makehgtform('axisrotate', a, -atan2(norm(a),V(3,3)));
R=T(1:3,1:3);
A_rot = A*R;
A_rot = A_rot + [xyz0(1) xyz0(2) 0]; % move so the centers are aligned in z-diection
% Plot the raw data
close all
scatter3(X,Y,Z,0.1,"magenta")
hold on
scatter3(A_rot(:,1),A_rot(:,2),A_rot(:,3),0.1,'blue');
xlabel('X-Axis','FontSize',14,'FontWeight','bold')
ylabel('Y-Axis','FontSize',14,'FontWeight','bold')
zlabel('Z-Axis','FontSize',14,'FontWeight','bold')
axis equal
##### 3 Kommentare1 älteren Kommentar anzeigen1 älteren Kommentar ausblenden
Elisa am 22 Jul. 2024 um 8:21
dear @Umar thank you so much!! i get an error: "At least one END is missing. The statement beginning here does not have a matching end. "
Umar am 22 Jul. 2024 um 11:15
Bearbeitet: Umar am 23 Jul. 2024 um 14:14
Hi Elisa,
Please see Avni response below. She has already provided the solution you are seeking.

Melden Sie sich an, um zu kommentieren.

### Antworten (1)

Avni Agrawal am 23 Jul. 2024 um 11:50
Hi @Elisa,
I understand that you are trying to extract the profiles oriented according to the maximum slope of the cloud. To achieve this:
• Compute the maximum slope direction and rotate the point cloud:
% code already mentioned in the question
[~,~,V]=svd(A,0);
a=cross(V(:,3),[0;0;1]);
T=makehgtform('axisrotate', a, -atan2(norm(a),V(3,3)));
R=T(1:3,1:3);
A_rot = A*R;
A_rot = A_rot + [xyz0(1) xyz0(2) 0]; % move so the centers are aligned in z-diection
• Extract Vertical planes: Create vertical planes at fixed intervals and extract points that lie on these planes
• Extract Z values: For each plane, extract the Z values of the points that lie on that
clear all
close all
clc
% Load the data
X = XYZ(:,1);
Y = XYZ(:,2);
Z = XYZ(:,3);
xyz0 = mean(XYZ, 1);
A = XYZ - xyz0; % Center the data at zero
% Find the direction of most variance using SVD and rotate the data to make that the x-axis
[~,~,V] = svd(A, 0);
a = cross(V(:,3), [0;0;1]);
T = makehgtform('axisrotate', a, -atan2(norm(a), V(3,3)));
R = T(1:3,1:3);
A_rot = A * R;
A_rot = A_rot + [xyz0(1) xyz0(2) 0]; % Move so the centers are aligned in z-direction
% Plot the raw data
close all
scatter3(X, Y, Z, 0.1, "magenta")
hold on
scatter3(A_rot(:,1), A_rot(:,2), A_rot(:,3), 0.1, 'blue');
xlabel('X-Axis', 'FontSize', 14, 'FontWeight', 'bold')
ylabel('Y-Axis', 'FontSize', 14, 'FontWeight', 'bold')
zlabel('Z-Axis', 'FontSize', 14, 'FontWeight', 'bold')
axis equal
% Extract profiles
fixed_distance = 1; % Define the fixed distance between profiles
x_min = min(A_rot(:,1));
x_max = max(A_rot(:,1));
profile_planes = x_min:fixed_distance:x_max;
for i = 1:length(profile_planes)
plane_x = profile_planes(i);
% Find points within a small tolerance of the plane
tolerance = 0.05; % Adjust the tolerance as needed
profile_points = A_rot(abs(A_rot(:,1) - plane_x) < tolerance, :);
if ~isempty(profile_points)
% Plot the profile points
scatter3(profile_points(:,1), profile_points(:,2), profile_points(:,3), 10, 'filled');
% Extract Z values
Z_values = profile_points(:,3);
% Display Z values (or process further as needed)
disp(['Profile at X = ', num2str(plane_x), ': Z values = ', num2str(Z_values')]);
end
end
hold off
I hope this helps!
##### 1 Kommentar-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden
Elisa am 29 Jul. 2024 um 13:59
dear @Avni Agrawal thank you so much for the help! I tried my self to implement the code to do that but I am myssing the way to determine the maximum slope and thus the equation of the line (that becomes a plane then) parallel to it. If you see in my code I inserted an equation that I manually calculated based on the coordinates of selected points on the point cloud plot. Can you help me implementing this code to automatically find the planes aligned on the maxium slope of the projected plane? thank you so much
clear all
close all
clc
X = XYZ(:,1);
Y = XYZ(:,2);
Z = XYZ(:,3);
xyz0=mean(XYZ,1);
A=XYZ-xyz0; % center the data at zero
% Find the direction of most variance using SVD and rotate the data to make
% that the x-axis
[~,~,V]=svd(A,0);
a=cross(V(:,3),[0;0;1]);
T=makehgtform('axisrotate', a, -atan2(norm(a),V(3,3)));
R=T(1:3,1:3);
A_rot = A*R;
A_rot = A_rot + [xyz0(1) xyz0(2) 0]; % move so the centers are aligned in z-diection
% Plot the original point cloud
close all
scatter3(X,Y,Z,0.1,"magenta")
hold on
scatter3(A_rot(:,1),A_rot(:,2),A_rot(:,3),0.1,'blue');
xlabel('X-Axis','FontSize',14,'FontWeight','bold')
ylabel('Y-Axis','FontSize',14,'FontWeight','bold')
zlabel('Z-Axis','FontSize',14,'FontWeight','bold')
axis equal
% line equation parameters
slope = -69 / 41;
intercept = 45999 / 41;
% Define the coordinate ranges
x_min = min(A_rot(1,:));
x_max = max(A_rot(1,:));
z_min = 0;
z_max = max(A_rot(3,:));
% Number of planes
num_planes = 5;
% Extend the length of the line by increasing the range of x values
extended_x_min = x_min - 50;
extended_x_max = x_max + 50;
% Generate points on the initial line for plotting
x_line = linspace(extended_x_min, extended_x_max, 100);
y_line = slope * x_line + intercept;
z_line = zeros(size(x_line)); % Line in the XY plane at z = 0
% Define the spacing between planes
intercept_spacing = linspace(intercept, intercept + (num_planes - 1) * 100, num_planes);
% Generate the grid for the planes
[x_plane, z_plane] = meshgrid(linspace(extended_x_min, extended_x_max, 100), linspace(-100, 100, 100));
% Plot the line
hold on;
plot3(x_line, y_line, z_line, 'r-', 'LineWidth', 2);
min_z_values = zeros(num_planes, 1);
max_z_values = zeros(num_planes, 1);
% Plot the vertical planes and calculate min/max z values
for i = 1:num_planes
current_intercept = intercept_spacing(i);
y_plane = slope * x_plane + current_intercept;
surf(x_plane, y_plane, z_plane, 'FaceAlpha', 0.5, 'EdgeColor', 'none');
% Calculate the distance of each point to the current plane
y_plane_point = slope * A_rot(:,2) + current_intercept;
distances = abs(A_rot(:,2) - y_plane_point);
% Select points within a mm range of the plane
within_range = distances <= 150;
% Extract the z values of the points within the range
% pippo=A_rot(:,3)
z_values_within_range = squeeze(A_rot(within_range,3));
% Calculate min and max z values for the selected points
if ~isempty(z_values_within_range)
min_z_values(i) = min(z_values_within_range); %,'omitnan');
max_z_values(i) = max(z_values_within_range); %,'omitnan');
else
min_z_values(i) = NaN; % No points within range
max_z_values(i) = NaN; % No points within range
end
end
% Add labels and title
xlabel('X');
ylabel('Y');
zlabel('Z');
title('Parallel Vertical Planes and Point Cloud Analysis');
% Add grid and axis equal
grid on;
axis equal;
legend('Line', 'Vertical Planes');
hold off;
% Display the min and max z values
% disp('Min Z values for each plane:');
disp(min_z_values);
% disp('Max Z values for each plane:');
disp(max_z_values);
Max_z_values_mm = 70.69/1000*max_z_values
Min_z_values_mm = 70.69/1000*min_z_values
Rz = Max_z_values_mm - Min_z_values_mm
E5 = 1.5715*Rz+4.0318
mean_E5 = mean(E5)
E6 = 4.4192*Rz.^0.6482
mean_E6 = mean(E6)

Melden Sie sich an, um zu kommentieren.

### Kategorien

Mehr zu Point Cloud Processing finden Sie in Help Center und File Exchange

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by