forked from netstim/leaddbs
-
Notifications
You must be signed in to change notification settings - Fork 0
/
ea_cylinder.m
134 lines (105 loc) · 3.71 KB
/
ea_cylinder.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
function [Cylinder, EndPlate1, EndPlate2] = ea_cylinder(X1,X2,r,n,cyl_color,closed,lines)
%
% This function constructs a cylinder connecting two center points
%
% Usage :
% [Cylinder EndPlate1 EndPlate2] = Cylinder(X1+20,X2,r,n,'r',closed,lines)
%
% Cylinder-------Handle of the cylinder
% EndPlate1------Handle of the Starting End plate
% EndPlate2------Handle of the Ending End plate
% X1 and X2 are the 3x1 vectors of the two points
% r is the radius of the cylinder
% n is the no. of elements on the cylinder circumference (more--> refined)
% cyl_color is the color definition like 'r','b',[0.52 0.52 0.52]
% closed=1 for closed cylinder or 0 for hollow open cylinder
% lines=1 for displaying the line segments on the cylinder 0 for only
% surface
%
% Typical Inputs
% X1=[10 10 10];
% X2=[35 20 40];
% r=1;
% n=20;
% cyl_color='b';
% closed=1;
%
% Calculating the length of the cylinder
length_cyl=norm(abs(X2-X1));
% Creating a circle in the YZ plane
t=linspace(0,2*pi,n)';
x2=r*cos(t);
sintheta = sin(t); sintheta(n) = 0;
x3=r*sintheta;
% Creating the points in the X-Direction
x1=[length_cyl 0];
% Creating (Extruding) the cylinder points in the X-Directions
xx3=repmat(x1,length(x2),1);
xx2=repmat(x2,1,2);
xx1=repmat(x3,1,2);
% Drawing two filled cirlces to close the cylinder
if closed==1
hold on
%EndPlate1=fill3(xx1(:,1),xx2(:,1),xx3(:,1),'r');
f.vertices=[xx1(:,1),xx2(:,1),xx3(:,1)];
f.vertices=[mean(f.vertices);f.vertices];
f.faces=zeros(length(xx1)-1,3);
for face=1:length(f.faces)
f.faces(face,:)=[1,face+1,face+2];
end
EndPlate1=patch('Vertices',f.vertices,'Faces',f.faces);
f.vertices=[xx1(:,2),xx2(:,2),xx3(:,2)];
f.vertices=[[0,0,0];f.vertices];
f.faces=zeros(length(xx1)-1,3);
for face=1:length(f.faces)
f.faces(face,:)=[1,face+1,face+2];
end
EndPlate2=patch('Vertices',f.vertices,'Faces',f.faces);
end
% Plotting the cylinder along the X-Direction with required length starting
% from Origin
Cylinder=mesh(xx1,xx2,xx3);
% Defining Unit vector along the Z-direction
unit_Vx=[0 0 1];
% Calulating the angle between the x direction and the required direction
% of cylinder through dot product
angle_X1X2=acos( dot( unit_Vx,(X2-X1) )/( norm(unit_Vx)*norm(X2-X1)) )*180/pi;
% Finding the axis of rotation (single rotation) to roate the cylinder in
% X-direction to the required arbitrary direction through cross product
axis_rot=cross([1 0 0],(X2-X1) );
% Rotating the plotted cylinder and the end plate circles to the required
% angles
% if angle_X1X2~=0 % Rotation is not needed if required direction is along X
% rotate(Cylinder,axis_rot,angle_X1X2,[0 0 0])
% if closed==1
%
% rotate(EndPlate1,axis_rot,angle_X1X2,[0 0 0])
% rotate(EndPlate2,axis_rot,angle_X1X2,[0 0 0])
% end
% end
% Till now cylinder has only been aligned with the required direction, but
% position starts from the origin. so it will now be shifted to the right
% position
if closed==1
vx=get(EndPlate1,'Vertices');
vx(:,1)=vx(:,1)+X2(1); vx(:,2)=vx(:,2)+X2(2); vx(:,3)=vx(:,3)+X2(3);
set(EndPlate1,'Vertices',vx);
vx=get(EndPlate2,'Vertices');
vx(:,1)=vx(:,1)+X2(1); vx(:,2)=vx(:,2)+X2(2); vx(:,3)=vx(:,3)+X2(3);
set(EndPlate2,'Vertices',vx);
end
set(Cylinder,'XData',get(Cylinder,'XData')+X2(1))
set(Cylinder,'YData',get(Cylinder,'YData')+X2(2))
set(Cylinder,'ZData',get(Cylinder,'ZData')+X2(3))
% Setting the color to the cylinder and the end plates
set(Cylinder,'FaceColor',cyl_color)
if closed==1
set([EndPlate1 EndPlate2],'FaceColor',cyl_color,'EdgeColor','none')
else
EndPlate1=[];
EndPlate2=[];
end
% If lines are not needed making it disapear
if lines==0
set(Cylinder,'EdgeAlpha',0)
end