Can someone help us modify our current code in order to store all of the desired output variables as columns of a single file?
1 次查看(过去 30 天)
显示 更早的评论
Below is the code that we need to modify in order to get all the desired output variables for each curved-edge polygon/arc-polygon (i.e., vertex numbers/coordinates for each arc-polygon, number of edges for each arc-polygon, total perimeter length for each arc-polygon, area of each arc-polygon, or in other words, areas_by_algebra for each arc-polygon), and store all of these outputs in columns within a single file for future analysis:
% close all figures and clear all variables
close all;
clearvars;
% load data
load('fpep.mat');
% create the object
obj=TripletGraph2(fpep);
% draw the figure
% obj.plotspline;
obj.plotcirc;
% hold the figure and draw upon the previous figure
axis equal;
hold on;
% get polygons from the object
[polygons, vertices_index] = getPolyshapes(obj);
% get the arc point information
arcpoint = obj.ArcPoints;
% interpolation length
sq=linspace(0,1,100)-0.5;
% for each polygon, calculate its area and draw it
for i= 1:1:length(polygons)
% get the vertices of current polygon
vertices = polygons(i).Vertices;
% calculate its approximation area
area(i) = polyarea(vertices (:,1), vertices (:,2));
% plot the current polygon
plot(polygons(i));
axis equal;
hold on;
% for each edge of this polygon, build the circular segment from the arc
for edge = 1:1:size(vertices,1)
% use start and end points C1, C2 of this edge as a key, find the two
% bridge points V1 V2 between them
start_point = vertices(edge, :);
if edge + 1 > size(vertices,1)
end_point = vertices(mod(edge+1, size(vertices,1)), :);
else
end_point = vertices(edge+1, :);
end
temp = arcpoint(:,:,sum(arcpoint(:,1,:)==start_point')==2);
res = temp(:,:,sum(end_point'==temp(:,4,:))==2);
if isempty(res)
temp = arcpoint(:,:,sum(arcpoint(:,1,:)==end_point')==2);
res = temp(:,:,sum(start_point'==temp(:,4,:))==2);
end
if isempty(res)
error('can not find V1 and V2!!')
end
% C1 start point, C2 end point, V1 V2 bridge points
C1 = res(:,1,:);
C2 = res(:,4,:);
V1 = res(:,2,:);
V2 = res(:,3,:);
% calcuate by polyarea
L=norm(C2-C1);
U=(C2-C1)/L;
dV1=(V1-C1)/norm(V1-C1);
dV2=(V2-C2)/norm(V2-C2);
theta1=acosd(dot( dV1, U) );
theta2=acosd(dot( dV2, -U) );
theta=(theta1+theta2)/2;
arcradius=norm(C1-C2)/(2*sind(theta));
W=cross( cross([U;0],[dV1;0]), [U;0]);
W=W(1:2)/norm(W);
D=L/2;
t=2*D*sq;
Dcot=D*cotd(theta);
s=sqrt(Dcot.^2-(t.^2-D.^2))-Dcot;
XY=(C1+C2)/2+U*t+W(1:2)*s;
X = XY(1,:);
Y = XY(2,:);
% estimate arc segment as polygon
circle_segement = polyshape(X, Y);
plot(circle_segement);
area_by_polygon = polyarea(X, Y);
% end here
% calculate arc segment w/ algebra
ca = 2* theta;
r = arcradius;
circle = pi * r ^ 2;
triangle = r ^ 2 * sind(ca) / 2;
area_by_algebra = circle * ca / 360 - triangle;
% algebra ends here
% if V1 V2 live inside the polygonal area, we should substract the circle_segement
if inpolygon([V1(1),V2(1)], [V1(2),V2(2)], vertices (:,1), vertices (:,2))
area(i) = area(i) - area_by_algebra;
% if V1 V2 live outside the polygonal area, we should add the circle_segement
else
area(i) = area(i) + area_by_algebra;
end
end
end
% display the area values
disp(area);
Thanks in advance for your help!
0 个评论
回答(1 个)
Pravin Jagtap
2019-12-23
编辑:Pravin Jagtap
2019-12-23
Hello Steve,
The objective mentioned in the question can be easily achieved by using the ‘save’ command. Refer following example for saving the variables the variables
p = rand(1,10); % column vector
q = rand(1,10); % column vector
% make sure that the variable you want save are in require format
data = [p' q'] % collect the variables you want to save in one matrix
save('pqfile.mat','data') % 'pqfile.mat' file will get created in the same directory
% load it for further analysis
For more details, refer following documentation
~Pravin
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Computational Geometry 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!