SURF Surface-Plot

Hi there,

I plotted something in a regualr 3D plot. Now my task is to code a surface-plot of the exact same thing.

So that the it's not just connected points, as in my actual plot, but a whole surface around the lines.

I hope that somebody can understand my problem and help me with it.

n = 3;
n1 = n-1;
a = 20;             
b = 10;             
P = [0 b;0 0;a 0];  
T = 15;             
H = 8;              
R = 2;              
h = 6;              
syms t s(t)
B = bernsteinMatrix(n1,t);
bezierCurve = B*P;
s(t) = int(norm(diff(bezierCurve)),0,t);
snum = linspace(0,s(1),T);
for i = 1:T
    tnum(i) = vpasolve(snum(i)==s(t),t);   
end
px = double(subs(bezierCurve(:,1),t,tnum)).';
py = zeros(T,1);
pz = double(subs(bezierCurve(:,2),t,tnum)).';
normalToCurve = diff(bezierCurve)*[0 1;-1 0];
normalToCurve = normalToCurve/norm(normalToCurve);
newNormalToCurve = [double(subs(normalToCurve(1),t,tnum))' double(subs(normalToCurve(2),t,tnum))'];
newNormalToCurvex = newNormalToCurve(:,1);
newNormalToCurvez = newNormalToCurve(:,2);

%%%%%%%%%%%Obere Linie
for i = 1:T
    S = double(s((i-1)/(T-1)));
    d = R*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    delta(i) = d;
end
delta = delta';
pxnew1 = px+newNormalToCurvex.*delta;
pznew1 = pz+newNormalToCurvez.*delta;
for i = 1:T
    S = double(s((i-1)/(T-1)));
    rhilfe = (h/2)*S/double(s(1));
    r(i) = rhilfe;
end
r = r';
pynew1a = - r;
pynew1b = r;

%%%%%%%%%%Untere Linie
delta = - delta;
pxnew2 = px+newNormalToCurvex.*delta;
pznew2 = pz+newNormalToCurvez.*delta;
pynew2a = pynew1a;
pynew2b = pynew1b;

%%%%%%%%%%seitliche Linien (schwarz)
for i = 1:T
    S = double(s((i-1)/(T-1)));
    chilfe = R*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    c(i) = chilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    lhilfe = H/2*S/double(s(1));
    l(i) = lhilfe;
end
l = l';
c = c';
pynew3 = ones(T,1).*c;
pxnew3a = px-newNormalToCurvex.*l;
pxnew3b = px+newNormalToCurvex.*l;
pznew3a = pz-newNormalToCurvez.*l;
pznew3b = pz+newNormalToCurvez.*l;
pynew3a = pynew3;
pynew3b = -pynew3;

%%%%%%%%%%%%4Ecken
%%%vorne links/rechts (magenta)
Rnewx = R*cos(pi/4);
Rnewy = R*cos(pi/4);
for i = 1:T
    S = double(s((i-1)/(T-1)));
    jhilfe = Rnewx*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    j(i) = jhilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    ehilfe = Rnewy*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    e(i) = ehilfe;
end
j = j';
e = e';
pxnew4 = px+newNormalToCurvex.*j;
pznew4 = pz+newNormalToCurvez.*j;
pynew4a = -ones(T,1).*e;
pynew4b = ones(T,1).*e;

%%%%%vorne links/recht(grün)
pxnew5 = px-newNormalToCurvex.*j;
pznew5 = pz-newNormalToCurvez.*j;
pynew5a = pynew4a;
pynew5b = pynew4b;

%%%%%%%%%%PLOTS
plot3(px,py,pz,'b-')
axis equal
axis tight
grid on
hold on
plot3(pxnew1,pynew1a,pznew1,'r-')
plot3(pxnew1,pynew1b,pznew1,'r-')
plot3(pxnew2,pynew1a,pznew2,'r-')
plot3(pxnew2,pynew1b,pznew2,'r-')
plot3(pxnew3a,pynew3a,pznew3a,'k-')
plot3(pxnew3b,pynew3a,pznew3b,'k-')
plot3(pxnew3a,pynew3b,pznew3a,'k-')
plot3(pxnew3b,pynew3b,pznew3b,'k-')
plot3(pxnew4,pynew4a,pznew4,'m-')
plot3(pxnew4,pynew4b,pznew4,'m-')
plot3(pxnew5,pynew5a,pznew5,'g-')
plot3(pxnew5,pynew5b,pznew5,'g-')
hold off

pxgesamt = [px;pxnew1; pxnew1; pxnew2; pxnew2; pxnew3a;pxnew3b;pxnew3a;pxnew3b;pxnew4; pxnew4; pxnew5; pxnew5];
pygesamt = [py;pynew1a;pynew1b;pynew2a;pynew2b;pynew3a;pynew3a;pynew3b;pynew3b;pynew4a;pynew4b;pynew5a;pynew5b];
pzgesamt = [pz;pznew1; pznew1; pznew2; pznew2; pznew3a;pznew3b;pznew3a;pznew3b;pznew4; pznew4; pznew5; pznew5];
pgesamt = [pxgesamt pygesamt pzgesamt];

surf([pxgesamt pygesamt pzgesamt]) %%%% that obviously does not work

Answer 1

[Star Strider]

Concatenate the column vectors horizontally (to create matrices) instead of vertically (to create column vectors) and it sort of works.

There are still some ptoblems. I leave the resolutions of those to you.

I also used the three column vectors with griddata to create matrices from them using interpolation. You can experiment with that as well.

I’m not certain how to resolve the inconsistiencies in creating matrices from the vectors so that the surf provides a more pleasing and representative appearance. Therea are several interpolation techniques to experiment with.

n = 3;
n1 = n-1;
a = 20;             
b = 10;             
P = [0 b;0 0;a 0];  
T = 15;             
H = 8;              
R = 2;              
h = 6;              
syms t s(t)
B = bernsteinMatrix(n1,t);
bezierCurve = B*P;
s(t) = int(norm(diff(bezierCurve)),0,t);
snum = linspace(0,s(1),T);
for i = 1:T
    tnum(i) = vpasolve(snum(i)==s(t),t);   
end
px = double(subs(bezierCurve(:,1),t,tnum)).';
py = zeros(T,1);
pz = double(subs(bezierCurve(:,2),t,tnum)).';
normalToCurve = diff(bezierCurve)*[0 1;-1 0];
normalToCurve = normalToCurve/norm(normalToCurve);
newNormalToCurve = [double(subs(normalToCurve(1),t,tnum))' double(subs(normalToCurve(2),t,tnum))'];
newNormalToCurvex = newNormalToCurve(:,1);
newNormalToCurvez = newNormalToCurve(:,2);

%%%%%%%%%%%Obere Linie
for i = 1:T
    S = double(s((i-1)/(T-1)));
    d = R*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    delta(i) = d;
end
delta = delta';
pxnew1 = px+newNormalToCurvex.*delta;
pznew1 = pz+newNormalToCurvez.*delta;
for i = 1:T
    S = double(s((i-1)/(T-1)));
    rhilfe = (h/2)*S/double(s(1));
    r(i) = rhilfe;
end
r = r';
pynew1a = - r;
pynew1b = r;

%%%%%%%%%%Untere Linie
delta = - delta;
pxnew2 = px+newNormalToCurvex.*delta;
pznew2 = pz+newNormalToCurvez.*delta;
pynew2a = pynew1a;
pynew2b = pynew1b;

%%%%%%%%%%seitliche Linien (schwarz)
for i = 1:T
    S = double(s((i-1)/(T-1)));
    chilfe = R*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    c(i) = chilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    lhilfe = H/2*S/double(s(1));
    l(i) = lhilfe;
end
l = l';
c = c';
pynew3 = ones(T,1).*c;
pxnew3a = px-newNormalToCurvex.*l;
pxnew3b = px+newNormalToCurvex.*l;
pznew3a = pz-newNormalToCurvez.*l;
pznew3b = pz+newNormalToCurvez.*l;
pynew3a = pynew3;
pynew3b = -pynew3;

%%%%%%%%%%%%4Ecken
%%%vorne links/rechts (magenta)
Rnewx = R*cos(pi/4);
Rnewy = R*cos(pi/4);
for i = 1:T
    S = double(s((i-1)/(T-1)));
    jhilfe = Rnewx*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    j(i) = jhilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    ehilfe = Rnewy*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    e(i) = ehilfe;
end
j = j';
e = e';
pxnew4 = px+newNormalToCurvex.*j;
pznew4 = pz+newNormalToCurvez.*j;
pynew4a = -ones(T,1).*e;
pynew4b = ones(T,1).*e;

%%%%%vorne links/recht(grün)
pxnew5 = px-newNormalToCurvex.*j;
pznew5 = pz-newNormalToCurvez.*j;
pynew5a = pynew4a;
pynew5b = pynew4b;

%%%%%%%%%%PLOTS
plot3(px,py,pz,'b-')
axis equal
axis tight
grid on
hold on
plot3(pxnew1,pynew1a,pznew1,'r-')
plot3(pxnew1,pynew1b,pznew1,'r-')
plot3(pxnew2,pynew1a,pznew2,'r-')
plot3(pxnew2,pynew1b,pznew2,'r-')
plot3(pxnew3a,pynew3a,pznew3a,'k-')
plot3(pxnew3b,pynew3a,pznew3b,'k-')
plot3(pxnew3a,pynew3b,pznew3a,'k-')
plot3(pxnew3b,pynew3b,pznew3b,'k-')
plot3(pxnew4,pynew4a,pznew4,'m-')
plot3(pxnew4,pynew4b,pznew4,'m-')
plot3(pxnew5,pynew5a,pznew5,'g-')
plot3(pxnew5,pynew5b,pznew5,'g-')
hold off

% Matrices — 
pxgesamt = [px, pxnew1, pxnew1, pxnew2, pxnew2, pxnew3a, pxnew3b, pxnew3a, pxnew3b, pxnew4, pxnew4, pxnew5, pxnew5];
pygesamt = [py, pynew1a, pynew1b, pynew2a, pynew2b, pynew3a, pynew3a, pynew3b, pynew3b, pynew4a, pynew4b, pynew5a, pynew5b];
pzgesamt = [pz, pznew1,  pznew1,  pznew2,  pznew2,  pznew3a, pznew3b, pznew3a, pznew3b, pznew4,  pznew4,  pznew5,  pznew5];

% Column Vectors
pxgesamtv = [px;pxnew1; pxnew1; pxnew2; pxnew2; pxnew3a;pxnew3b;pxnew3a;pxnew3b;pxnew4; pxnew4; pxnew5; pxnew5];
pygesamtv = [py;pynew1a;pynew1b;pynew2a;pynew2b;pynew3a;pynew3a;pynew3b;pynew3b;pynew4a;pynew4b;pynew5a;pynew5b];
pzgesamtv = [pz;pznew1; pznew1; pznew2; pznew2; pznew3a;pznew3b;pznew3a;pznew3b;pznew4; pznew4; pznew5; pznew5];
% pgesamt = [pxgesamt pygesamt pzgesamt];
% 
figure
surfc(pxgesamt, pygesamt, pzgesamt) %%%% that obviously does not work
colormap(turbo)
xlabel('pxgesamt')
ylabel('pygesamt')
zlabel('pzgesamt')
axis('equal')

N = numel(pxgesamtv);
N = ceil(N/5)

N = 39

xv = linspace(min(pxgesamtv), max(pxgesamtv), N);
yv = linspace(min(pygesamtv), max(pygesamtv), N);
[X, Y] = ndgrid(xv, yv);
Z = griddata(pxgesamtv,pygesamtv, pzgesamtv, X, Y);

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

figure
surfc(X, Y, Z, 'EdgeAlpha',0.25)
grid on
colormap(turbo)
xlabel('pxgesamt')
ylabel('pygesamt')
zlabel('pzgesamt')
axis('equal')

.

[Tim Schaller]

Thanks a lot!!

As you said there are still some problems in the first plot, but it really gets coles to what I am looking for.

Know I need to figure out, how to do it properly.

But to be honest, I don't konw at all how to do that.

[Star Strider]

My pleasure!

Those vectors are not going to lend themselves easily to any interpolation method. There appears to be more to the first surf plot I posted that is immediately obvious.

That can be seen by changing the 'FaceAlpha' value —

n = 3;
n1 = n-1;
a = 20;             
b = 10;             
P = [0 b;0 0;a 0];  
T = 15;             
H = 8;              
R = 2;              
h = 6;              
syms t s(t)
B = bernsteinMatrix(n1,t);
bezierCurve = B*P;
s(t) = int(norm(diff(bezierCurve)),0,t);
snum = linspace(0,s(1),T);
for i = 1:T
    tnum(i) = vpasolve(snum(i)==s(t),t);   
end
px = double(subs(bezierCurve(:,1),t,tnum)).';
py = zeros(T,1);
pz = double(subs(bezierCurve(:,2),t,tnum)).';
normalToCurve = diff(bezierCurve)*[0 1;-1 0];
normalToCurve = normalToCurve/norm(normalToCurve);
newNormalToCurve = [double(subs(normalToCurve(1),t,tnum))' double(subs(normalToCurve(2),t,tnum))'];
newNormalToCurvex = newNormalToCurve(:,1);
newNormalToCurvez = newNormalToCurve(:,2);

%%%%%%%%%%%Obere Linie
for i = 1:T
    S = double(s((i-1)/(T-1)));
    d = R*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    delta(i) = d;
end
delta = delta';
pxnew1 = px+newNormalToCurvex.*delta;
pznew1 = pz+newNormalToCurvez.*delta;
for i = 1:T
    S = double(s((i-1)/(T-1)));
    rhilfe = (h/2)*S/double(s(1));
    r(i) = rhilfe;
end
r = r';
pynew1a = - r;
pynew1b = r;

%%%%%%%%%%Untere Linie
delta = - delta;
pxnew2 = px+newNormalToCurvex.*delta;
pznew2 = pz+newNormalToCurvez.*delta;
pynew2a = pynew1a;
pynew2b = pynew1b;

%%%%%%%%%%seitliche Linien (schwarz)
for i = 1:T
    S = double(s((i-1)/(T-1)));
    chilfe = R*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    c(i) = chilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    lhilfe = H/2*S/double(s(1));
    l(i) = lhilfe;
end
l = l';
c = c';
pynew3 = ones(T,1).*c;
pxnew3a = px-newNormalToCurvex.*l;
pxnew3b = px+newNormalToCurvex.*l;
pznew3a = pz-newNormalToCurvez.*l;
pznew3b = pz+newNormalToCurvez.*l;
pynew3a = pynew3;
pynew3b = -pynew3;

%%%%%%%%%%%%4Ecken
%%%vorne links/rechts (magenta)
Rnewx = R*cos(pi/4);
Rnewy = R*cos(pi/4);
for i = 1:T
    S = double(s((i-1)/(T-1)));
    jhilfe = Rnewx*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    j(i) = jhilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    ehilfe = Rnewy*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    e(i) = ehilfe;
end
j = j';
e = e';
pxnew4 = px+newNormalToCurvex.*j;
pznew4 = pz+newNormalToCurvez.*j;
pynew4a = -ones(T,1).*e;
pynew4b = ones(T,1).*e;

%%%%%vorne links/recht(grün)
pxnew5 = px-newNormalToCurvex.*j;
pznew5 = pz-newNormalToCurvez.*j;
pynew5a = pynew4a;
pynew5b = pynew4b;

%%%%%%%%%%PLOTS
plot3(px,py,pz,'b-')
axis equal
axis tight
grid on
hold on
plot3(pxnew1,pynew1a,pznew1,'r-')
plot3(pxnew1,pynew1b,pznew1,'r-')
plot3(pxnew2,pynew1a,pznew2,'r-')
plot3(pxnew2,pynew1b,pznew2,'r-')
plot3(pxnew3a,pynew3a,pznew3a,'k-')
plot3(pxnew3b,pynew3a,pznew3b,'k-')
plot3(pxnew3a,pynew3b,pznew3a,'k-')
plot3(pxnew3b,pynew3b,pznew3b,'k-')
plot3(pxnew4,pynew4a,pznew4,'m-')
plot3(pxnew4,pynew4b,pznew4,'m-')
plot3(pxnew5,pynew5a,pznew5,'g-')
plot3(pxnew5,pynew5b,pznew5,'g-')
hold off

% Matrices — 
pxgesamt = [px, pxnew1, pxnew1, pxnew2, pxnew2, pxnew3a, pxnew3b, pxnew3a, pxnew3b, pxnew4, pxnew4, pxnew5, pxnew5];
pygesamt = [py, pynew1a, pynew1b, pynew2a, pynew2b, pynew3a, pynew3a, pynew3b, pynew3b, pynew4a, pynew4b, pynew5a, pynew5b];
pzgesamt = [pz, pznew1,  pznew1,  pznew2,  pznew2,  pznew3a, pznew3b, pznew3a, pznew3b, pznew4,  pznew4,  pznew5,  pznew5];

% Column Vectors
pxgesamtv = [px;pxnew1; pxnew1; pxnew2; pxnew2; pxnew3a;pxnew3b;pxnew3a;pxnew3b;pxnew4; pxnew4; pxnew5; pxnew5];
pygesamtv = [py;pynew1a;pynew1b;pynew2a;pynew2b;pynew3a;pynew3a;pynew3b;pynew3b;pynew4a;pynew4b;pynew5a;pynew5b];
pzgesamtv = [pz;pznew1; pznew1; pznew2; pznew2; pznew3a;pznew3b;pznew3a;pznew3b;pznew4; pznew4; pznew5; pznew5];
% pgesamt = [pxgesamt pygesamt pzgesamt];
% 
figure
surfc(pxgesamt, pygesamt, pzgesamt, 'FaceAlpha',0.25) %%%% that obviously does not work
colormap(turbo)
xlabel('pxgesamt')
ylabel('pygesamt')
zlabel('pzgesamt')
axis('equal')

N = numel(pxgesamtv);
N = ceil(N/5)

N = 39

xv = linspace(min(pxgesamtv), max(pxgesamtv), N);
yv = linspace(min(pygesamtv), max(pygesamtv), N);
[X, Y] = ndgrid(xv, yv);
Z = griddata(pxgesamtv,pygesamtv, pzgesamtv, X, Y);

Warning: Duplicate data points have been detected and removed - corresponding values have been averaged.

figure
surfc(X, Y, Z, 'EdgeAlpha',0.25)
grid on
colormap(turbo)
xlabel('pxgesamt')
ylabel('pygesamt')
zlabel('pzgesamt')
axis('equal')

% Matrices — 
pxgesamt = [px, pxnew1, pxnew1, pxnew2, pxnew2, pxnew3a, pxnew3b, pxnew3a, pxnew3b, pxnew4, pxnew4, pxnew5, pxnew5];
pygesamt = [py, pynew1a, pynew1b, pynew2a, pynew2b, pynew3a, pynew3a, pynew3b, pynew3b, pynew4a, pynew4b, pynew5a, pynew5b];
pzgesamt = [pz, pznew1,  pznew1,  pznew2,  pznew2,  pznew3a, pznew3b, pznew3a, pznew3b, pznew4,  pznew4,  pznew5,  pznew5];

[pygesamtsrt,idx] = sortrows(pygesamt.',1);
pygesamtsrt = pygesamtsrt.';
pxgesamtsrt = pxgesamt(:,idx);
pzgesamtsrt = pzgesamt(:,idx);

figure
surfc(pxgesamtsrt, pygesamtsrt, pzgesamtsrt, 'FaceAlpha',0.25, 'EdgeAlpha',0.5) %%%% that obviously does not work
colormap(turbo)
xlabel('pxgesamt')
ylabel('pygesamt')
zlabel('pzgesamt')
axis('equal')

[az,el] = view

az = -37.5000

el = 30

% figure
% hold on
% for k = 1:numel(idx)-1
%     surfc(pxgesamtsrt(:,k:k+1), pygesamtsrt(:,k:k+1), pzgesamtsrt(:,k:k+1), 'FaceAlpha',0.25, 'EdgeAlpha',0.5) %%%% that obviously does not work
% end
% hold off
% grid on
% colormap(turbo)
% xlabel('pxgesamt')
% ylabel('pygesamt')
% zlabel('pzgesamt')
% axis('equal')
% view(az,el)

It might be worthwhile to experiment with the matrices in the first surf plot, for example sorting ‘pygesamt’ column-wise (using sortrows on the temporarily transposed matrix, perhaps sorting by the first column) and then using that index vector to sort the columns of the others. Doing that experiment appears to be slightly better, however it did not completely solve the problem.

.

[Tim Schaller]

I got it now, thank you

n = 3;
n1 = n-1;
a = 20;             
b = 10;             
P = [0 b;0 0;a 0];  
T = 15;             
H = 8;              
R = 2;              
h = 6;              
syms t s(t)
B = bernsteinMatrix(n1,t);
bezierCurve = B*P;
s(t) = int(norm(diff(bezierCurve)),0,t);
snum = linspace(0,s(1),T);
for i = 1:T
    tnum(i) = vpasolve(snum(i)==s(t),t);   
end
px = double(subs(bezierCurve(:,1),t,tnum)).';
py = zeros(T,1);
pz = double(subs(bezierCurve(:,2),t,tnum)).';
normalToCurve = diff(bezierCurve)*[0 1;-1 0];
normalToCurve = normalToCurve/norm(normalToCurve);
newNormalToCurve = [double(subs(normalToCurve(1),t,tnum))' double(subs(normalToCurve(2),t,tnum))'];
newNormalToCurvex = newNormalToCurve(:,1);
newNormalToCurvez = newNormalToCurve(:,2);
%%%%%%%%%%%Obere Linie
for i = 1:T
    S = double(s((i-1)/(T-1)));
    d = R*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    delta(i) = d;
end
delta = delta';
pxnew1 = px+newNormalToCurvex.*delta;
pznew1 = pz+newNormalToCurvez.*delta;
for i = 1:T
    S = double(s((i-1)/(T-1)));
    rhilfe = (h/2)*S/double(s(1));
    r(i) = rhilfe;
end
r = r';
pynew1a = - r;
pynew1b = r;
%%%%%%%%%%Untere Linie
delta = - delta;
pxnew2 = px+newNormalToCurvex.*delta;
pznew2 = pz+newNormalToCurvez.*delta;
pynew2a = pynew1a;
pynew2b = pynew1b;
%%%%%%%%%%seitliche Linien (schwarz)
for i = 1:T
    S = double(s((i-1)/(T-1)));
    chilfe = R*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    c(i) = chilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    lhilfe = H/2*S/double(s(1));
    l(i) = lhilfe;
end
l = l';
c = c';
pynew3 = ones(T,1).*c;
pxnew3a = px-newNormalToCurvex.*l;
pxnew3b = px+newNormalToCurvex.*l;
pznew3a = pz-newNormalToCurvez.*l;
pznew3b = pz+newNormalToCurvez.*l;
pynew3a = pynew3;
pynew3b = -pynew3;
%%%%%%%%%%%%4Ecken
%%%vorne links/rechts (magenta)
Rnewx = R*cos(pi/4);
Rnewy = R*cos(pi/4);
for i = 1:T
    S = double(s((i-1)/(T-1)));
    jhilfe = Rnewx*(1-S/double(s(1)))+(H/2)*S/double(s(1));
    j(i) = jhilfe;
end
for i = 1:T
    S = double(s((i-1)/(T-1)));
    ehilfe = Rnewy*(1-S/double(s(1)))+(h/2)*S/double(s(1));
    e(i) = ehilfe;
end
j = j';
e = e';
pxnew4 = px+newNormalToCurvex.*j;
pznew4 = pz+newNormalToCurvez.*j;
pynew4a = -ones(T,1).*e;
pynew4b = ones(T,1).*e;
%%%%%vorne links/recht(grün)
pxnew5 = px-newNormalToCurvex.*j;
pznew5 = pz-newNormalToCurvez.*j;
pynew5a = pynew4a;
pynew5b = pynew4b;
%%%%%%%%%%PLOTS
plot3(px,py,pz,'b-')
axis equal
axis tight
grid on
hold on
plot3(pxnew1,pynew1a,pznew1,'r-')
plot3(pxnew1,pynew1b,pznew1,'r-')
plot3(pxnew2,pynew1a,pznew2,'r-')
plot3(pxnew2,pynew1b,pznew2,'r-')
plot3(pxnew3a,pynew3a,pznew3a,'k-')
plot3(pxnew3b,pynew3a,pznew3b,'k-')
plot3(pxnew3a,pynew3b,pznew3a,'k-')
plot3(pxnew3b,pynew3b,pznew3b,'k-')
plot3(pxnew4,pynew4a,pznew4,'m-')
plot3(pxnew4,pynew4b,pznew4,'m-')
plot3(pxnew5,pynew5a,pznew5,'g-')
plot3(pxnew5,pynew5b,pznew5,'g-')
hold off

pxgesamt = [pxnew2,  pxnew5,  pxnew3b, pxnew3b, pxnew4,  pxnew1,  pxnew1,  pxnew4,  pxnew3b, pxnew3b, pxnew5,  pxnew2 , pxnew2];
pygesamt = [pynew2b, pynew5b, pynew3,  pynew3,  pynew4b, pynew1b, pynew1a, pynew4a, pynew3b  pynew3b, pynew5a, pynew2a, pynew2b];
pzgesamt = [pznew2,  pznew5,  pznew3a, pznew3b, pznew4,  pznew1,  pznew1,  pznew4,  pznew3b, pznew3a, pznew5,  pznew2,  pznew2 ];


surf(pxgesamt, pygesamt, pzgesamt)
colormap(turbo)
axis equal
axis tight

[Star Strider]

As always, my pleasure!

I wasn’t certain what it was supposed to look like.

.