首页 > 编程语言 >《空间三角面片对相交判断算法》的matlab实现

《空间三角面片对相交判断算法》的matlab实现

时间:2023-11-11 12:55:56浏览次数:47  
标签:面片 T2 flag T1 A1 算法 matlab dot B1

function [flag] = InsectTriPatch(T1,T2)
% 判断两个空间三角形面片是否相交    
% T1=[0 0 0;
    %     2 0 0;
    %     0 1.5 0;
    %     0 0 1];
    % T2=[0 0 -1;
    %     2 0 -1;
    %     0 2 -1;
    %     0 0 1];
    
    % 出自:《空间三角面片对相交判断算法》
    %   T=[A;B;C;n]
    %   A=[x y z]
    A1=T1(1,:);  % A1点的坐标:[x y z]
    B1=T1(2,:);
    C1=T1(3,:);
    n1=T1(4,:);
    
    A2=T2(1,:);
    B2=T2(2,:);
    C2=T2(3,:);
    n2=T2(4,:);
    
    %判断T2的3个顶点是否在平面π1 同侧
    A1A2=A2-A1;
    A1B2=B2-A1;
    A1C2=C2-A1;
    temp1=[dot(A1A2,n1) dot(A1B2,n1) dot(A1C2,n1)];
    %判断T1的3个顶点是否在平面π2同侧
    A2A1=A1-A2;
    A2B1=B1-A2;
    A2C1=C1-A2;
    temp2=[dot(A2A1,n2) dot(A2B1,n2) dot(A2C1,n2) ];
    
    % 判断temp1 或 temp2 是否同号
    if any([all(temp1<0) all(temp1>0) all(temp2<0) all(temp2>0)])
        % T1,T2 不相交,return 0
        flag=0;
        return
    else
        if all([isequal(cross(n1,n2),[0 0 0]),isequal(temp1,[0 0 0])])
            % T1和 T2共面
            if  DetectVertexInTri(T1,T2)
                % 存在T1的顶点在T2内部或者存在T2的顶点在T1内部,判定T1、T2相交, 返回真
                flag=1;
                return
            else
                if DetectVertex(T1,T2)
                    flag=1;
                    return
                else
                    flag=0;
                    return
                end
            end
        else
            % T1 T2异面
            if any([EdgeInsectTri(T1,T2) EdgeInsectTri(T2,T1)])
                %存在 T2 边与 T1 相交或存在 T1 的边与 T2相交
                flag=1;
                return
            else
                flag=0;
                return
            end
        end
    end
end

function is_in= PointInTri(T1,P)
% 函数实现: 判断一点是否在三角形T1内
    %   T=[A;B;C;n]
    %   A=[x y z]
    A1=T1(1,:);  % A1点的坐标:[x y z]
    B1=T1(2,:);
    C1=T1(3,:);
    n1=T1(4,:);
    
    A1B1=B1-A1;
    A1C1=C1-A1;
    B1C1=C1-B1;
    B1A1=A1-B1;
    C1A1=A1-C1;
    C1B1=B1-C1;
    
    A1P=P-A1;
    B1P=p-B1;
    C1P=P-C1;
    
    k1=dot(dot(A1P,cross(A1B1,n1)),dot(A1C1,cross(A1B1,n1)));
    k2=dot(dot(B1P,cross(B1C1,n1)),dot(B1A1,cross(B1C1,n1)));
    k3=dot(dot(C1P,cross(C1A1,n1)),dot(C1B1,cross(C1A1,n1)));
    
    if any([k1 k2 k3]<0)
        % k1 k2 k3 存在负值,则P不在T1中
        is_in=0;
        return
    else
        is_in=1;
    end
end

function flag = DetectVertexInTri(T1,T2)
% 比较:存在T1 的顶点在 T2 内部或者存在 T2 的顶点在T1 内部?
% loop T1
    for i=1:3 
        if PointInTri(T2,T1(i,:)) 
            % 说明T1的顶点在T2内部,T1,T2 相交,return 1
            flag=1;
            return
        end
    end
    % loop T2
    for i=1:3 
        if PointInTri(T1,T2(i,:)) 
            % 说明T2的顶点在T1内部,T1,T2 相交,return 1
            flag=1;
            return
        end
    end
end

function flag=fun2(T2,p1,p2,n)

    A2=T2(1,:);
    B2=T2(2,:);
    C2=T2(3,:);

    h1=dot(A2-p1,cross(p2-p1,n));
    h2=dot(B2-p1,cross(p2-p1,n));
    h3=dot(C2-p1,cross(p2-p1,n));
    temp=[h1 h2 h3];
    if any([all(temp<0) all(temp>0)])
        flag=0;
    else
        flag=1;
    end
end

function flag=DetectVertex(T1,T2)
    % 共面-步骤二
    %   T=[A;B;C;n]
    %   A=[x y z]
    % A1点的坐标:[x y z]
    n1=T1(4,:);
    n2=T2(4,:);

    index=[1 2;2 3;3 1];
    % T2的3个顶点
    for i=1:3
        if fun2(T2,T1(index(i,1),:),T1(index(i,2),:),n1)
            flag=1;
        else
            %有一条边不符合要求,直接返回0,不相交。
            flag=0;
            return
        end
    end
    % T1的3个顶点
    for i=1:3
        if fun2(T1,T2(index(i,1),:),T2(index(i,2),:),n2)
            flag=1;
        else
            %有一条边不符合要求,直接返回0,不相交。
            flag=0;
            return
        end
    end
    flag=1;
end

function flag= EdgeInsectTri(T1,T2)
% 函数实现:判断T2的线段A2B2 B2C2 C2A2是否与三角面片T1相交
    index=[1 2;2 3;3 1];
    tf=0;
    for i=1:3
        if fun3(T2(index(i,1),:),T2(index(i,2),:),T1)
            % 存在线段和面片相交
            flag=1;
            return
        else
            tf=0;
        end
    end
    flag=tf;
end

function flag=fun3(p1,p2,T1)
% 实现:判断边p1p2是否与T1相交
    %   T=[A;B;C;n]
    %   A=[x y z]
    A1=T1(1,:);  % A1点的坐标:[x y z]
    B1=T1(2,:);
    C1=T1(3,:);
    n1=T1(4,:);
    
    s1=cross(p1-A1,n1);
    s2=cross(p2-A1,n1);

    if all([s1==[0 0 0] s2==[0 0 0]])
        % p1p2和T1共面 
        if any([PointInTri(T1,p1) PointInTri(T1,p2)])
            % p1 p2,至少一个点在三角形T1中,则p1p2和T1相交
            flag=1;
            return
        else
            flag=0;
            return
        end
    else
        % p1p2和T1异面
        if dot(s1,s2)>0
            % P1P2在T1同侧, 判定不相交;
            flag=0;
            return
        else
            g1=dot(dot(p2-A1,cross(B1-A1,p1-A1)),dot(C1-A1,cross(B1-A1,p1-A1)));
            g2=dot(dot(p2-B1,cross(C1-B1,p1-B1)),dot(A1-B1,cross(C1-B1,p1-B1)));
            g3=dot(dot(p2-C1,cross(A1-C1,p1-C1)),dot(B1-C1,cross(A1-C1,p1-C1)));

            if all([g1>0 g2>0 g3>0])
                % g1 g2 g3 的计算结果均为非负,p1p2和T1相交
                flag=1;
                return
            elseif any([g1<0 g2<0 g3<0])
                % g1 g2 g3 的计算结果存在负值,p1p2和T1不相交
                flag=0;
                return
            end
        end
    end
end

标签:面片,T2,flag,T1,A1,算法,matlab,dot,B1
From: https://www.cnblogs.com/aksoam/p/17825798.html

相关文章

  • 11.11算法
    题目填充每个节点的下一个右侧节点指针给定一个 完美二叉树 ,其所有叶子节点都在同一层,每个父节点都有两个子节点。二叉树定义如下:structNode{intval;Node*left;Node*right;Node*next;}填充它的每个next指针,让这个指针指向其下一个右侧节点。如果找不到下一个右......
  • 线性代数 · 矩阵 · Matlab | Moore-Penrose 伪逆矩阵代码实现
    背景-Moore-Penrose伪逆矩阵:对任意矩阵\(A\in\mathbbC^{m\timesn}\),其Moore-Penrose逆矩阵\(A^+\in\mathbbC^{n\timesm}\)存在且唯一。定义:若矩阵G满足\(AGA=A,~GAG=G,~(AG)^H=AG,~(GA)^H=GA\),则G是Moore-Penrose逆矩阵,可以记作\(A^+\)。性质:\(A^......
  • matlab 字符串处理函数
    ​ %字符串处理a=' a';b='b b';c='cccc';m=''%获取字符串长度length(a)    %连接两个字符串,每个字符串最右边的空格被裁切d=strcat(a,c) length(d)%连接多行字符串,每行长度可不等,自动把非最长字符串最右边补空格%使与最长字符串相等,会忽略空字符串e=s......
  • MATLAB快捷键
    ​ 一、索引混排版备注:删除了如F1(帮助)等类型的常见快捷命令 ​编辑SHIFT+DELETE永久删除DELETE删除ALT+ENTER属性ALT+F4关闭CTRL+F4关闭ALT+TAB切换ALT+ESC切换ALT+空格键窗口菜单CTRL+ESC开始菜单拖动某一项时按CTRL复制所选项目拖动某一项时按CTRL......
  • MATLAB寻找最大值和最小值
     最大值C=max(A)最小值C=min(A)如果A是一个向量,max(A)返回A中的最大/最小元素。如果A是一个矩阵,max(A)将A中的每一列作为一个向量,并返回一个行向量,这个行向量包含了每一列的最大/小元素。比如:a=[1,7,10];b=min(a);得到:b=1而:a=[1,7,10;6,5,4;6,5,5]......
  • MATLAB对矩阵按照某一列排序
    转载:matlab对矩阵按照某一列排序_matlab对矩阵按列升序排列-CSDN博客升序排列:命令:data=[1,2,3;7,8,9;4,5,6];a1=sortrows(data,1);%按照第一列排序(升序),其他列与排序结果一一对应。a2=sortrows(data,2);%按照第二列排序(升序),其他列与排序结果一一对应......
  • MATLAB中的disp函数
    disp函数会直接将内容输出在Matlab命令窗口中,比如:par1=csvread('front_surface_pressure_005.csv',1,0);disp(par1);运行之后,会在命令窗口输出 front_surface_pressure_005.csv 文件中的数据内容:-0.014801.38680.053941.3740109.3200-0......
  • Matlab绘制信号包络线
    ✅作者简介:热爱科研的算法开发者,Python、Matlab项目可交流、沟通、学习。......
  • Matlab代码优化之道
    ​ 一、遵守PerformanceAcceleration的规则关于什么是“PerformanceAcceleration”请参阅matlab的帮助文件。1、只有使用以下数据类型,matlab才会对其加速:logical,char,int8,uint8,int16,uint16,int32,uint32,double而语句中如果使用了非以上的数据类型则不会加速,如numeric......
  • 基于OFDM的水下图像传输通信系统matlab仿真
    1.算法运行效果图预览  2.算法运行软件版本matlab2022a 3.算法理论概述      基于OFDM的水下图像传输通信系统是一种用于在水下环境中传输图像数据的通信系统。它采用了OFDM(OrthogonalFrequencyDivisionMultiplexing)技术,这种技术在水下通信中具有一些优......