一种计算太阳辐照度的遮挡检测方法及装置与流程

本发明属于地表物体之间遮挡检测技术领域，尤其涉及一种基于不规则三角网(triangularirregularnetwork，tin)模型计算太阳辐照度的遮挡检测方法及装置。

背景技术：

近年来，许多研究学者利用激光点云、建筑三维模型和摄影测量的方法对城市地区建筑物表面接收到的太阳辐照度进行估算。在进行太阳辐照度估算的过程中，遮挡检测分析是其中一个重要的环节。如果建筑物表面受到周围其它物体的遮挡，那么其接收到的太阳辐照度会大大减少，所以遮挡检测分析结果的精确与否关系到最后太阳辐照度的计算结果。由于激光点云数据具有精度高、易获取、采集过程中受周围环境限制小的优点，在城市建模、阴影检测和太阳辐照度估算等领域有一定的应用优势。但是激光点云数据也有其不足，例如建筑物立面上的少量点云会导致信息不够详细、数据量冗余，这为遮挡检测计算带来较大的困难。相比之下，建筑物三维模型在某些方面与激光点云数据正好呈现相反的特点，它数据量小、模型准确，但较难获取。在估算太阳辐照度之前，往往需要花费大量时间精力进行建模工作，因此建筑物三维模型不适合在城市区域进行大规模太阳辐照度的快速估算。

随着倾斜摄影测量技术的发展，利用倾斜摄影测量技术可以较容易地获取城市表面的tin模型。除此之外，tin模型对建筑物立面信息、树木等信息的记录较为详细。因此可以考虑使用tin模型对建筑物表面接收到的太阳辐照度进行估算。

目前，大多数激光点云数据采用的遮挡检测分析方法是天空角系数(skyviewfactor,svf)算法。svf算法的关键是计算每个方向的最大遮挡高度角，而每个方向的最大遮挡高度角与水平方向的划分有关。若划分角度太小则影响计算速度，若划分角度太大则影响计算精度。因此svf算法具有一定的缺陷。由于射线法对遮挡检测的计算结果精确，建筑三维模型多采用射线法进行遮挡检测。但对大量tin模型而言，应用射线法则意味着需要逐点计算，这么做显然会影响计算速度，因此射线法对大量tin模型的遮挡检测计算并不适用。

到目前为止，针对大量tin模型进行遮挡检测的问题，还没有一种方法可以同时实现高检测精度和高计算效率。

技术实现要素：

针对现有技术的不足，本发明提供一种计算太阳辐照度的遮挡检测方法及装置，以提高检测精度和计算效率。

本发明是通过如下的技术方案来解决上述技术问题的：一种计算太阳辐照度的遮挡检测方法，包括：

步骤1：定义正交的三维坐标系，其中，所述坐标系中的y轴指向正北方，z轴指向铅垂方向，x轴与y轴和z轴构成右手坐标系；

步骤2：读入待检测建筑物表面重建后的tin模型，获取模型中每个三角形对应的外接圆圆心的坐标以及所有三角形中最大外接圆的半径；以对所述tin模型重采样所获得的点作为待计算点；

步骤3：根据待检测建筑物所处区域的经纬度和计算时间得到太阳高度角β和太阳方位角α，再根据所述太阳高度角β和太阳方位角α得到太阳光线方向向量v；

步骤4：将所述tin模型中的所有待计算点、所有三角形顶点以及所有三角形外接圆圆心同时旋转，使太阳光线方向向量v与z轴平行，然后记录旋转后所有待计算点的坐标以及所有三角形外接圆圆心的坐标；

步骤5：根据所述步骤4中旋转后每个待计算点的坐标，找出该待计算点对应的三角形集合，所述三角形集合中的每个三角形均包含在以该待计算点为中心，以最大外接圆半径为半径的范围内；并确定该待计算点对应的三角形集合中每个三角形外接圆圆心的坐标，由三角形集合中每个三角形外接圆圆心的坐标构成每个待计算点的外接圆圆心坐标集合；

步骤6：在每个待计算点的外接圆圆心坐标集合中，剔除外接圆圆心的z轴坐标值小于该待计算点的z轴坐标值的圆心，即剔除不对该待计算点造成遮挡的三角形，剩余对该待计算点造成遮挡的外接圆圆心坐标集合；

步骤7：在xy平面上，判断每个待计算点是否位于对该待计算点造成遮挡的外接圆圆心坐标集合所对应的三角形内，如果是，则该待计算点受到遮挡，否则该待计算点不受遮挡。

进一步地，所述步骤3中，太阳高度角β和太阳方位角α的计算公式分别为：

β＝sin^-1(sinbsinδ+cosbcosδcosh)，

设太阳光线方向向量则：

其中，b为当地的纬度，δ为赤纬角，h为时角。

进一步地，所述步骤4中，所有待计算点、所有三角形顶点以及所有三角形外接圆圆心的旋转方式为：同时绕z轴顺时针旋转α，再绕x轴顺时针旋转90°－β。

进一步地，所述步骤4中，旋转后每个待计算点以及每个三角形外接圆圆心的坐标的计算公式均为：

其中，α为太阳方位角，β为太阳高度角，(x,y,z)为旋转前的坐标，(x”,y”,z”)为旋转后的坐标。

进一步地，所述步骤5中，先根据旋转后每个三角形外接圆圆心的平面坐标构建kd(k-dimensionaltree，kd-tree)树，再利用kd树找出每个待计算点对应的三角形集合。

进一步地，所述步骤7中，判断每个待计算点是否位于对该待计算点造成遮挡的外接圆圆心坐标集合所对应的三角形内的具体步骤包括：

步骤7.1：判断每个待计算点是否位于对该待计算点造成遮挡的外接圆圆心坐标集合所对应的三角形的外接矩形内，如果是，则转入步骤7.2，否则该待计算点不受遮挡；

步骤7.2：采用面积法判断该待计算点是否位于当前三角形内：计算当前三角形的三个顶点与该待计算点分别构成的三个三角形的面积之和，如果该面积之和大于当前三角形的面积，则认为该待计算点在当前三角形的外部，即该待计算点不受当前三角形遮挡，否则该待计算点受当前三角形遮挡。

本发明还提供一种计算太阳辐照度的遮挡检测装置，包括：

获取模块，用于读入待检测建筑物表面重建后的tin模型，获取模型中每个三角形对应外接圆圆心的坐标以及所有三角形中最大外接圆的半径；以对所述tin模型重采样所获得的点作为待计算点；

角度计算模块，用于根据待检测建筑物所处区域的经纬度和计算时间得到太阳高度角β和太阳方位角α，再根据所述太阳高度角β和太阳方位角α得到太阳光线方向向量v；

坐标旋转模块，用于将所述tin模型中所有待计算点、所有三角形顶点以及所有三角形外接圆圆心同时旋转，使太阳光线方向向量v与z轴平行，然后记录旋转后所有三角形待计算点的坐标以及所有三角形外接圆圆心的坐标；

集合生成模块，用于根据旋转后每个待计算点的坐标，找出该待计算点对应的三角形集合，所述三角形集合中的每个三角形均包含在以该待计算点为中心，以最大外接圆半径为半径的范围内；并确定该待计算点对应的三角形集合中每个三角形外接圆圆心的坐标，由所述三角形集合中每个三角形外接圆圆心的坐标构成该待计算点的外接圆圆心坐标集合；

剔除模块，用于在每个待计算点的外接圆圆心坐标集合中，剔除外接圆圆心的z轴坐标值小于该待计算点的z轴坐标值的圆心，即剔除不对该待计算点造成遮挡的三角形，剩余对该待计算点造成遮挡的外接圆圆心坐标集合；

判断模块，用于在xy平面上，判断每个待计算点是否位于对该待计算点造成遮挡的外接圆圆心坐标集合所对应的三角形内，如果是，则该待计算点受到遮挡，否则该待计算点不受遮挡。

进一步地，所述角度计算模块中，太阳高度角β和太阳方位角α的具体计算公式为：

β＝sin^-1(sinbsinδ+cosbcosδcosh)，

其中，b为当地纬度，δ为赤纬角，h为时角。

进一步地，所述坐标旋转模块中，所有待计算点、所有三角形顶点以及所有三角形外接圆圆心的旋转方式为：同时绕z轴顺时针旋转α，再绕x轴顺时针旋转90°－β。

进一步地，所述集合生成模块中，采用kd树找出每个待计算点对应的三角形集合。

有益效果

与现有技术相比，本发明所提供的一种计算太阳辐照度的遮挡检测方法及装置，通过对tin模型和待计算点的旋转，使任意时刻的太阳光线方向与z轴平行，然后比较待计算点与对应三角形集合中三角形顶点之间的z轴坐标值的大小来简化遮挡检测的计算过程，最后仅判断待计算点是否在对应的三角形内即可判断待计算点是否受到遮挡，该方法大大的简化了遮挡检测的计算过程，提高了遮挡检测的计算效率；该方法采用kd树搜索每个待计算点对应的三角形集合，有助于解决大量待计算点数据带来计算速度缓慢的问题，进一步提高了计算效率。

该方法采用tin模型进行遮挡检测计算，tin模型能够较准确地反映地表物体的三维信息，基于tin模型的遮挡检测方法有效提高了遮挡检测的精度。

附图说明

为了更清楚地说明本发明的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一个实施例，对于本领域普通技术人员来说，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明实施例中遮挡检测的流程图；

图2是本发明实施例中点绕z轴旋转的过程图；

图3是本发明实施例中点绕x轴旋转的过程图；

图4是本发明实施例中遮挡检测的计算结果；

图5是本发明实施例中相同时刻拍摄的倾斜摄影图片，其中，1-阴影区域。

具体实施方式

下面结合本发明实施例中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1所示，本发明所提供的一种计算太阳辐照度的遮挡检测方法，包括：

1、定义正交的三维坐标系，其中，坐标系中的y轴指向正北方，z轴指向铅垂方向，x轴与y轴和z轴构成右手坐标系。

2、读入待检测建筑物表面重建后的tin模型，获取模型中每个三角形对应外接圆圆心的坐标以及所有三角形中最大外接圆的半径；以对tin模型重采样所获得的点作为待计算点。

利用倾斜摄影测量、激光点云等技术手段获得大量的点云数据，经过表面重建后得到大量的不规则三角形(或三角网)，由大量不规则三角形构成待检测建筑物的tin模型，tin模型的构建属于现有技术。tin模型的读入实际上就是读取一系列三角形顶点的坐标，根据每个三角形的顶点坐标可以获取该三角形对应的外接圆圆心oj的坐标(xj,yj,zj)，以及外接圆的半径rj，找出外接圆半径的最大值rmax＝max(rj)，该最大值rmax即为所有三角形中最大外接圆的半径。

采用软件cloudcompare对tin模型进行重采样，重采样所获得的点作为待计算点，待计算点包括多个，对哪个区域进行遮挡检测分析就选取tin模型的哪个区域进行重采样以获取待计算点。待计算点可能是三角形的顶点，也可能是三角形面上的点，再使用软件的hougenormalcomputation功能计算每个待计算点的法向量并导出，第i个待计算点pi的坐标及其法向量表示为(xi,yi,zi,nx,ny,nz)。

3、根据待检测建筑物所处区域的经纬度和计算时间得到太阳高度角β和太阳方位角α，再根据太阳高度角β和太阳方位角α得到太阳光线方向向量v。

太阳高度角β和太阳方位角α的计算公式分别为：

β＝sin^-1(sinbsinδ+cosbcosδcosh)

设太阳光线方向向量则：

其中，b为当地的纬度(例如建筑物所处区域的经纬度为北纬28.169591°，东经112.922825°)，δ为赤纬角(是指地球赤道平面与太阳和地球中心的连线之间的夹角)，h为时角。根据计算时间(例如2019年4月16日16时33分)可以计算出赤纬角和时角，赤纬角和时角的计算为现有技术，可参考《太阳能应用技术》中关于太阳辐射的内容。

4、将tin模型中的所有待计算点pi、所有三角形顶点以及所有三角形外接圆圆心oj同时旋转，使太阳光线方向向量v与z轴平行，然后记录旋转后每个待计算点pi″的坐标以及每个三角形外接圆圆心o″j的坐标。

所有待计算点pi、所有三角形顶点以及所有三角形外接圆圆心oj的旋转方式为：同时绕z轴顺时针旋转α(如图2所示，第一次旋转后的坐标为(x',y',z'))，再绕x轴顺时针旋转90°－β(如图3所示，第二次旋转后的坐标为(x”,y”,z”))，此时太阳光线方向向量与z轴平行。点旋转后的坐标可以通过以下通用公式来求取(此处的点包括待计算点、三角形顶点以及外接圆圆心)：

其中，α为太阳方位角，β为太阳高度角，(x,y,z)为旋转前的坐标，(x”,y”,z”)为旋转后的坐标。旋转后外接圆圆心o″j的坐标记为(x″j,y″j,z″j)，旋转后待计算点pi”的坐标记为(x″i,y″i,z″i)。本申请中，三角形外接圆圆心均指每个三角形对应的外接圆圆心。所有待计算点pi″、三角形顶点以及三角形外接圆圆心o″j的旋转即为tin模型的旋转，通过旋转使任意时刻的太阳光线方向与z轴平行，旋转处理后使对待计算点pi″的遮挡检测计算转变成比较待计算点pi″与对应三角形集合中三角形顶点之间的z轴坐标值(又称高程值)的大小，以及只需考虑在待计算点pi″一定范围(xy平面)内的三角形即可，在不影响计算精确度的情况下，大大简化了遮挡检测计算过程，提高了遮挡检测的计算效率。

5、根据旋转后每个待计算点pi″的坐标，找出该待计算点pi″对应的三角形集合，该三角形集合中的每个三角形均包含在以该待计算点pi″为中心，以最大外接圆半径rmax为半径的范围内；并确定该待计算点pi″对应的三角形集合中每个三角形外接圆圆心的坐标，由三角形集合中每个三角形外接圆圆心的坐标构成该待计算点pi″的外接圆圆心坐标集合。

根据旋转后外接圆圆心o″j的平面坐标(x″j,y″j)构建kd树，构建的kd树能够在输入一个点的xy坐标和查询范围r时，返回满足条件(即外接圆圆心所对应的三角形包含在以待计算点为中心，以最大外接圆半径rmax为半径的范围内)的外接圆圆心的索引号，通过索引号得到满足条件的外接圆圆心坐标，由这些满足条件的外接圆圆心的坐标构成每个待计算点pi″的外接圆圆心坐标集合ω，即

三角形对应的外接圆与该三角形是一一对应的，每个待计算点pi″的外接圆圆心坐标集合ω即为每个待计算点pi″对应的三角形集合(三角形集合中的三角形为可能对pi″造成遮挡的三角形)，通过对外接圆圆心坐标集合ω进行分析和判断，即可找出对该集合ω对应的待计算点pi″造成遮挡的三角形。本申请中，查询过程以最大外接圆半径rmax作为查询范围，确保能够找到所有对pi″造成遮挡的三角形，提高了遮挡检测的精确度。

kd树是一种对k维空间中的实例点进行存储以便对其进行快速检索或查询的树形数据结构，kd树的构建方法为现有技术。利用kd树找出每个待计算点pi″对应的三角形集合和外接圆圆心坐标集合ω，有助于解决大量待计算点造成的计算速度缓慢的问题，大大提高了对待计算点pi″近邻范围查询的效率，进一步提高了遮挡检测的计算效率。

6、在每个待计算点pi″的外接圆圆心坐标集合ω中，剔除外接圆圆心的z轴坐标值小于该待计算点的z轴坐标值的圆心，即剔除不对该待计算点造成遮挡的三角形，剩余对该待计算点造成遮挡的外接圆圆心坐标集合φ。

经过步骤3的旋转变换，太阳光线方向与坐标轴z平行，因此，集合ω中只有z轴坐标值大于pi″的z轴坐标值的外接圆所对应的三角形才有可能对pi″造成遮挡。为了简化后续的判断过程，在每个待计算点pi″的外接圆圆心坐标集合ω中，剔除不可能对pi″造成遮挡的三角形，即剔除集合ω中z轴坐标值小于pi″的z轴坐标值的外接圆圆心，剩余对该待计算点pi″造成遮挡的外接圆圆心坐标集合φ，

7、在xy平面上，判断每个待计算点pi″是否位于对该待计算点造成遮挡的外接圆圆心坐标集合φ所对应的三角形内(包括边界)，如果是，则该待计算点受到周围物体的遮挡，否则该待计算点不受遮挡。

具体的判断步骤包括：

步骤7.1：判断每个待计算点pi″是否位于对该待计算点造成遮挡的外接圆圆心坐标集合φ所对应的三角形的外接矩形内，如果是，则转入步骤7.2，否则该待计算点不受遮挡。先判断pi″是否在对应三角形的外接矩形内(只需要判断待计算点是否在外接矩形的坐标范围内即可)，只有在外接矩形内才需要进行下一步的判断，大大加快了判断速度，进一步提高了遮挡检测的计算效率。

步骤7.2：采用面积法判断该待计算点pi″是否位于当前三角形内：计算当前三角形的三个顶点与该待计算点分别构成的三个三角形的面积之和，如果该面积之和大于当前三角形的面积，则认为该待计算点在当前三角形的外部，即该待计算点不受当前三角形遮挡，否则该待计算点受当前三角形遮挡。

例如，当前三角形的三个顶点分别为(a,b,c)，则当前三角形的两个顶点与该待计算点构成的三个三角形分别为(a,b,pi″)、(b,c,pi″)、(a,c,pi″)。

图4为采用本发明的方法进行遮挡检测的计算结果，黑色区域表示阴影区域，图5是相同时刻拍摄的倾斜摄影照片，带有阴影区域(图5中的标号1)。通过图4和图5的对比可知，采用本发明的方法进行遮挡检测后的阴影区域结果与同时刻照片记录的实际阴影区域高度相似，由此可知，本发明的遮挡检测方法具有很高的检测精度。

本发明还提供一种计算太阳辐照度的遮挡检测装置，包括：

获取模块，用于读入待检测建筑物表面重建后的tin模型，获取模型中每个三角形对应外接圆圆心的坐标以及所有三角形中最大外接圆的半径；以对所述tin模型重采样所获得的点作为待计算点；

角度计算模块，用于根据待检测建筑物所处区域的经纬度和计算时间得到太阳高度角β和太阳方位角α，再根据所述太阳高度角β和太阳方位角α得到太阳光线方向向量v；

坐标旋转模块，用于将tin模型中所有待计算点、所有三角形顶点以及所有三角形外接圆圆心同时旋转，使太阳光线方向向量v与z轴平行，然后记录旋转后所有待计算点的坐标以及所有三角形外接圆圆心的坐标；

集合生成模块，用于根据旋转后每个待计算点的坐标，找出该待计算点对应的三角形集合，所述三角形集合中的每个三角形均包含在以该待计算点为中心，以最大外接圆半径为半径的范围内；并确定该待计算点对应的三角形集合中每个三角形外接圆圆心的坐标，由三角形集合中每个三角形外接圆圆心的坐标构成每个待计算点的外接圆圆心坐标集合；

剔除模块，用于在每个待计算点的外接圆圆心坐标集合中，剔除外接圆圆心的z轴坐标值小于该待计算点的z轴坐标值的圆心，即剔除不对该待计算点造成遮挡的三角形，剩余对该待计算点造成遮挡的外接圆圆心坐标集合；

判断模块，用于在xy平面上，判断每个待计算点是否位于对该待计算点造成遮挡的外接圆圆心坐标集合所对应的三角形内，如果是，则该待计算点受到遮挡，否则该待计算点不受遮挡。

在角度计算模块中，太阳高度角β和太阳方位角α的具体计算公式为：

β＝sin^-1(sinbsinδ+cosbcosδcosh)，

其中，b为当地纬度，δ为赤纬角，h为时角。

在坐标旋转模块中，所有待计算点、所有三角形顶点以及所有三角形外接圆圆心的旋转方式为：同时绕z轴顺时针旋转α，再绕x轴顺时针旋转90°－β。

在集合生成模块中，采用kd树找出每个待计算点对应的三角形集合。

以上所揭露的仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到变化或变型，都应涵盖在本发明的保护范围之内。