MATLAB-CATIA导入多段直线段

2022年“盈建科杯”中，需要使用不同的参数生成不同的结构，再导入Midas Civil进行受力分析。由于构件全为直杆，不需要将点阵拟合成样条曲线。提供一种方法：在MATLAB中进行结构设计（方便修改不同的参数），再使用CATIA自带的文件“GSD_PointSplineLoftFromExcel.xls”导入直线段。

示例：
main.m

clear;clc;
%initialize
heigh = [350,350,200];
outdiameter = [350,320,290];
indiameter = [220,190,180];

%% solve geometic

%circle
figure(3)
x = -0.5*indiameter(1):0.01:0.5*indiameter(1);
yptv = sqrt((0.5*indiameter(1))^2-x.^2);
yntv = -sqrt((0.5*indiameter(1))^2-x.^2);
plot(x,yptv)  %%
axis equal
%eight pol
hold on %%
tempa = outdiameter(1)/(1+sqrt(2));
Pts = [tempa/2,tempa*(sqrt(2)+1)/2;
    tempa*(sqrt(2)+1)/2,tempa/2;
    tempa*(sqrt(2)+1)/2,-tempa/2;
    tempa/2,-tempa*(sqrt(2)+1)/2;
    -tempa/2,-tempa*(sqrt(2)+1)/2;
     -tempa*(sqrt(2)+1)/2,-tempa/2;
      -tempa*(sqrt(2)+1)/2,tempa/2;
      -tempa/2,tempa*(sqrt(2)+1)/2;];
scatter(Pts(:,1),Pts(:,2)) %%
% qiedian
qiedian = findqiedian(Pts(8,:),1,indiameter(1),1); % jiaodianzuobiao shangxiabanyuan banjing zuoyouqiedian 
diandata = [Pts(8,:),0,qiedian,0];
scatter(qiedian(1),qiedian(2)) %%

%% get direction vector
C = getrotatemetrix(diandata);
Long = sqrt(sum((diandata(1:3)-diandata(4:6)).^2));
changepos = [90,90,90];
diwidth = Long;%dimiandechangdu
foldtimes = [3,3,3];% zhediedecishu 
connectdir = [1,1,1]; %1 is out to in -1 is in to out

%start
fixedpoint = [0,0,0;
    diwidth,0,0;
    diwidth,0,heigh(1);
    
    indiameter(1)-indiameter(2),0,heigh(1);
    diwidth+indiameter(1)-indiameter(2),0,heigh(1);
    diwidth+indiameter(1)-indiameter(2),0,heigh(1)+heigh(2);

    indiameter(1)-indiameter(3),0,heigh(1)+heigh(2);
    diwidth+indiameter(1)-indiameter(3),0,heigh(1)+heigh(2);
    diwidth+indiameter(1)-indiameter(3),0,heigh(1)+heigh(2)+heigh(3);
    
    0,0,heigh(1);
    0,0,heigh(1)+heigh(2);
    0,0,heigh(1)+heigh(2)+heigh(3)];
figure(1)
plot3(fixedpoint(:,1),fixedpoint(:,2),fixedpoint(:,3));
xlabel("X")
ylabel("Y")
zlabel("Z")
axis equal 
Line = [];
Line(end+1,:) = [fixedpoint(1,:),fixedpoint(2,:)];
Line(end+1,:) = [fixedpoint(2,:),fixedpoint(3,:)];
Line(end+1,:) = [fixedpoint(10,:),fixedpoint(5,:)];
Line(end+1,:) = [fixedpoint(1,:),fixedpoint(10,:)];

Line(end+1,:) = [fixedpoint(4,:),fixedpoint(4,1),0,heigh(1)+heigh(2)];
Line(end+1,:) = [fixedpoint(5,:),fixedpoint(6,:)];
Line(end+1,:) = [fixedpoint(4,1),0,heigh(1)+heigh(2),fixedpoint(8,:)];

Line(end+1,:) = [fixedpoint(8,:),fixedpoint(9,:)];
Line(end+1,:) = [fixedpoint(7,:),fixedpoint(7,1),0,heigh(1)+heigh(2)+heigh(3)];
Line(end+1,:) = [fixedpoint(7,1),0,heigh(1)+heigh(2)+heigh(3),fixedpoint(9,:)];

%midden
for i = 1:3
    if connectdir(i) == 1
    if i == 1
        temp = 0;
    else
        temp = heigh(i-1)+temp;
    end
    temppointout = [repmat(fixedpoint(1+3*(i-1),1),foldtimes(i),1),zeros(foldtimes(i),1),(temp+cumsum(heigh(i)/(foldtimes(i)).*ones(1,foldtimes(i))))'];
    tempa = (temppointout(foldtimes(i),3)-temppointout(foldtimes(i)-1,3))/2;
    temppointin = [repmat(fixedpoint(2+3*(i-1),1),foldtimes(i),1),zeros(foldtimes(i),1),temppointout(:,3)-tempa];
    end
    temppointout = [temppointout(1,:);temppointout];
    temppointout(1,3) = temppointout(1,3)-temppointout(end,3)+temppointout(end-1,3);
    for j = 1:foldtimes(i)
        Line(end+1,:) = [temppointin(j,:),temppointout(j,:)];
        Line(end+1,:) = [temppointin(j,:),temppointout(j+1,:)];
    end
end
figure(2)
hold on
for i = 1:size(Line,1)

    X = [Line(i,1);Line(i,4)];
%     Y = [Line(i,2);Line(i,5)];
    Z = [Line(i,3);Line(i,6)];
    plot(X,Z)
end
axis equal





testLine = Line(:,1:3);
testLine = testLine';
testLineTrans1 = inv(C)*testLine;
testLineTrans1 = testLineTrans1';

testLine = Line(:,4:6);
testLine = testLine';
testLineTrans2 = inv(C)*testLine;
testLineTrans2 = testLineTrans2';

testLineTrans = [testLineTrans1,testLineTrans2];


figure(4)
hold on
for i = 1:size(testLineTrans,1)

    X = [testLineTrans(i,1);testLineTrans(i,4)];
    Y = [testLineTrans(i,2);testLineTrans(i,5)];
    Z = [testLineTrans(i,3);testLineTrans(i,6)];
    plot3(X,Y,Z);

end
xlabel("X")
ylabel("Y")
zlabel("Z")
axis equal


testLineTrans11 = testLineTrans1+Long.*C(:,2)';
testLineTrans22 = testLineTrans2+Long.*C(:,2)';
testLineTrans = [testLineTrans11,testLineTrans22];

figure(5)
hold on
for i = 1:size(testLineTrans,1)

    X = [testLineTrans(i,1);testLineTrans(i,4)];
    Y = [testLineTrans(i,2);testLineTrans(i,5)];
    Z = [testLineTrans(i,3);testLineTrans(i,6)];
    plot3(X,Y,Z);

end
xlabel("X")
ylabel("Y")
zlabel("Z")
axis equal

Line = [Line;testLineTrans];

%% excel
writematrix('StartLoft','GSD.xlsx')
for i = 1:size(Line,1)
    str = ['A',num2str(4*(i-1)+2)];
    writematrix('StartCurve','GSD.xlsx',Range=str)
    str = ['A',num2str(4*(i-1)+3) ':' 'C',num2str(4*(i-1)+4)];
    writematrix(reshape(Line(i,:),3,2)','GSD.xlsx',Range=str);
    str = ['A',num2str(4*(i-1)+5)];
    writematrix('EndCurve','GSD.xlsx',Range=str);
end
str = ['A',num2str(4*(size(Line,1)-1)+6)];
writematrix('EndLoft','GSD.xlsx',Range=str);
str = ['A',num2str(4*(size(Line,1)-1)+7)];
writematrix('End','GSD.xlsx',Range=str);

findqiedian.m

function qiedian = findqiedian(jiaodian,shangxia,diameter,zuoyou)
% jiaodian ---- dingdian
% shangxia ---- 1:shang 2:xia 3:zhong
% diameter ---- banjing
% zouyou   ---- 1:zuo 2:you
syms x0
if shangxia == 1
    y0 = sqrt((0.5*diameter)^2-x0.^2);
    x = solve(jiaodian(2)-y0==-x0/y0*(jiaodian(1)-x0));
    x = double(x);
    x = x(zuoyou);
    qiedian = [x,sqrt((0.5*diameter)^2-x.^2)];
elseif shangxai == 2
    y0 = -sqrt((0.5*diameter(1))^2-x0.^2);
    x = solve(jiaodian(2)-y0==-x0/y0*(jiaodian(1)-x0));
    x = double(x);
    x = x(zuoyou);
    qiedian = [x,-sqrt((0.5*diameter)^2-x.^2)];
elseif shangxia == 3
    if jiaodian(2) > 0
        qiedian = [0,0.5*diameter];
    else
        qiedian = [0,-0.5*diameter];
    end
end

getrotatemetrix.m

function C = getrotatemetrix(diandata)
vecter1 = [0,0,1];
vecter2 = [diandata(4)-diandata(1),diandata(5)-diandata(2),diandata(6)-diandata(3)];
vecter2 = vecter2./sqrt(sum(vecter2.^2));
% satisfy the right hand rule's vector3
A = [vecter1(1:2);vecter2(1:2)];
B = -null(A);
vecter3 = [B',0];
eta = [vecter3',vecter2',vecter1'];
ksi = eye(3);
C = eta\ksi;

其中，可以将母体进行旋转和平移操作（使用坐标转换公式），从而设计不同的结构。

【注】：要将MATLAB生成的"GSD.xlsx"复制进 “GSD_PointSplineLoftFromExcel.xls”，才可使用宏。
【思考】：怎样让MATLAB在生成Excel时可以在里面附带宏？