本申请涉及遥感技术领域，特别是涉及一种地下空间资源量分析方法。
背景技术：
在城市中开展地下空间资源规划和建设，需要对已利用的地下空间资源进行了解，做到合理布局与开发。
目前，关于已利用地下空间资源量研究方面，大多数学者均以资源质量评价为研究目标；也有部分学者开展具体的评价手段、技术方法研究；还有部分学者开展评价体系的研究；但对于已利用资源量的调查研究很少，现在主要还是利用人工走访调研的形式，实地统计已利用资源量，存在调查成本高的问题。
技术实现要素：
基于此，有必要针对上述技术问题，提供一种能够解决调查成本高的问题的地下空间资源量分析方法。
一种地下空间资源量分析方法，所述方法包括：
获取待分析区域的卫星影像；
对所述卫星影像进行预处理，获得处理后的卫星影像；
采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果；
根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果；
根据太阳方位线和所述阴影提取结果进行分析，获得所述阴影提取结果中各阴影的阴影长度；
根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度；
将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
在其中一个实施例中，所述采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果的步骤，包括：
采用形态学建筑指数从所述处理后的卫星影像中进行建筑物的提取，获得第一提取结果；
采用深度学习从所述处理后的卫星影像中进行建筑物的提取，获得第二提取结果；
将所述第一提取结果和所述第二提取结果进行融合，获得融合后的提取结果；
对所述融合后的提取结果进行楼顶形状优化，获得建筑物提取结果。
在其中一个实施例中，所述根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果的步骤，包括：
将各特征指数的阴影表现力的强度进行对比，确定阴影表现力的强度最高的特征指数作为阴影提取的特征指数；
采用所述阴影提取的特征指数进行阴影提取工作，获得初步阴影提取结果；
对所述初步阴影提取结果采用聚类处理进行优化，获得阴影提取结果。
在其中一个实施例中，所述根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度的步骤，包括：
根据所述建筑物提取结果，确定待进行建筑物高度反演的建筑物；
基于太阳、卫星、所述建筑物三者之间的空间位置关系，以及所述建筑物的建筑物高度和阴影的关系，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度。
在其中一个实施例中，当太阳和卫星的方位相同，而太阳的高度角小于卫星的高度角，同时太阳射线在地面的投影与建筑物走向垂直时，所述建筑物高度反演模型的高度反演公式为：
h＝s·(tanα·tanβ)/(tanα-tanβ)
其中，α是卫星的高度角，β是太阳的高度角，s为建筑物在阴影提取结果中的阴影长度，h是建筑物的建筑物高度。
在其中一个实施例中，当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线在地面的投影和卫星视线在地面的投影位于建筑物的两侧时，所述建筑物高度反演模型的高度反演公式为：
h＝s·tanβ/sinγ
其中，γ为建筑物走向与太阳射线投影方向的夹角。
在其中一个实施例中，当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线投影和卫星视线投影位于建筑物同侧时，所述建筑物高度反演模型的高度反演公式为：
h＝s·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
其中，δ为建筑物走向与卫星视线投影方向的夹角。
在其中一个实施中，所述将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括：
根据建筑物高度与建筑物影响深度的转换关系表，将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度；
根据各所述建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
在其中一个实施中，所述根据各所述建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括:
根据各所述建筑物的建筑物影响深度，分别以建筑物影响深度为0m-10m、10m-30m、30m-50m和50m-100m为计算单元；
对各所述计算单元进行地下空间资源量利用情况叠加分析，获得地下空间资源量分析结果。
上述地下空间资源量分析方法，通过获取待分析区域的卫星影像；对所述卫星影像进行预处理，获得处理后的卫星影像；采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果；根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果；根据太阳方位线和所述阴影提取结果进行分析，获得所述阴影提取结果中各阴影的阴影长度；根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度；将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果，解决了实地统计已利用资源量，调查成本高的问题。
附图说明
图1为一个实施例中地下空间资源量分析方法的流程示意图；
图2为一个实施例中部分处理后的卫星影像的建筑物提取结果示意图；
图3为一个实施例中部分处理后的卫星影像示意图；
图4为一个实施例中部分处理后的卫星影像的阴影提取结果示意图；
图5为一个实施例中地下空间资源量分析方法中理想情况下星日同侧时几何关系示意图；
图6为一个实施例中地下空间资源量分析方法中理想情况下星日异侧时几何关系示意图；
图7为一个实施例中地下空间资源量分析方法中非理想情况下星日同侧时几何关系示意图；
图8为一个实施例中地下空间资源量分析方法中非理想情况下星日异侧时几何关系示意图。
具体实施方式
为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处描述的具体实施例仅仅用以解释本申请，并不用于限定本申请。
在一个实施例中，如图1所示，提供了一种地下空间资源量分析方法，包括以下步骤：
步骤s220，获取待分析区域的卫星影像。
其中，待分析区域是指需要进行地下空间资源量分析的地理区域。卫星影像是通过卫星获取的地面影像，具有分辨率高、覆盖区域广、更新周期短、光谱特征丰富等特点，被称为太空中的“千里眼”“透视眼”。
步骤s240，对卫星影像进行预处理，获得处理后的卫星影像。
其中，对卫星影像进行预处理的方式，包括辐射定标、大气校正、正射校正、影像融合、影像裁剪中的至少一种以上。
步骤s260，采用形态学建筑指数和深度学习从处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果。
其中，采用形态学建筑指数和深度学习分别从处理后的卫星影像进行建筑物的提取工作，并对两种提取的结果进行融合，再对融合后的结果进行楼顶形状优化，最终得到明显可识、清楚完整、边界清晰的建筑物单体提取结果，该建筑物单体提取结果即为建筑物提取结果，如图2所示的部分处理后的卫星影像的建筑物提取结果示意图。
在一个实施例中，采用形态学建筑指数和深度学习从处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果的步骤，包括：采用形态学建筑指数从处理后的卫星影像中进行建筑物的提取，获得第一提取结果；采用深度学习从处理后的卫星影像中进行建筑物的提取，获得第二提取结果；将第一提取结果和第二提取结果进行融合，获得融合后的提取结果；对融合后的提取结果进行楼顶形状优化，获得建筑物提取结果。
其中，形态学建筑指数(mbi)的基本思想是在建筑物的隐含特征(例如亮度、大小和对比度)和形态学操作符(例如通过重构、粒度测定和方向性获得高帽)之间建立关系，根据建筑物的隐含特征和建筑物的形态学操作符之间建立的关系，从处理后的卫星影像中进行建筑物的提取，获得第一提取结果。
形态学建筑指数建立在建筑物的边缘光谱变化较大，而建筑物的内部光谱变化较小的基础上。该形态学建筑指数的构建考虑到了建筑物的形状、方向、亮度和对比度等特征，具体实现步骤如下：
第一步，计算处理后的卫星影像的亮度b，其中，处理后的卫星影像的亮度包括处理后的卫星影像中各像素x的亮度值，像素x的亮度值计算：像素x的多光谱波段的最大光谱值被记录为像素x的亮度值：
b(x)＝max1≤k≤k(bandk(x))(1)
其中，b(x)表示为像素x的亮度值，k为可见光光谱波段数目，bandk(x)表示第k个光谱带像素x的光谱值，因为多光谱波段中的可见光波段对建筑物的光谱信息影响较大，所以这里选择可见光波段每个像素的最大值做为该像素x的亮度值。
第二步，mbi构建：利用top-hat的dmp通过一系列线性se重构来表示建筑的光谱结构特征(如对比度、大小、方向性)。用于构造mbi的形态学算子如下：
1、形态学白帽重建(w-th)
其中，w-th(d,s)是处理后的卫星影像重建后的形态学白帽，
2、计算形态学剖面mp：
其中，mpw-th(d,s)是处理后的卫星影像的形态学剖面，mpw-th(d,0)是在线性结构元素的尺度为0时，处理后的卫星影像的形态学剖面。
3、计算微分形态学剖面(dmp)：
dmpw-th(d,s)＝|mpw-th(d,(s δs))-mpw-th(d,s)|(4)
其中，dmpw-th(d,s)表示处理后的卫星影像的微分形态学剖面，δs是轮廓间隔，范围是smin≤s≤smax。
4、计算形态建筑指数mbi：
其中，mbi是处理后的卫星影像的形态建筑指数，d表示轮廓的方向数，s表示轮廓的尺度数。由于轮廓的方向数d的增加并不能提高建筑检测的精度，所以本研究考虑了四个方向(0°、45°、90°、135°)。结构元素smin、smax和δs大小应根据建筑的空间特征和空间分辨率来确定，在一个实施例中smax＝52、smin＝2、δs＝5。轮廓的尺度s的计算公式为：
mbi的构建是基于建筑结构在白色顶帽dmp直方图的大部分方向上都有较大的值，因为它们在这些方向上具有较高的局部对比度，因此，经常引用较大的mbi值，从处理后的卫星影像中进行建筑物的提取，获得第一提取结果。
采用深度学习从处理后的卫星影像中进行建筑物的提取，获得第二提取结果，具体地：利用包含建筑物5488个，总面积为84.0095km2的区域的卫星影像建立了训练集，训练集又经过旋转、镜像、模糊、光照调整等操作进行增强，最终获得了数量为10万张256*256大小图片的数据集，该数据集中只包括建筑物和非建筑两种类别，该卫星影像覆盖了不同的密集程度的建筑区域，既包括建筑稠密的地区又包括建筑物稀疏分布的区域，然后利用deeplab-v3网络，根据该训练集进行训练获得建筑物提取模型，将处理后的卫星影像输入建筑物提取模型，根据输入的建筑物提取模型根据输入的处理后的卫星影像提取出对应的建筑物，即第二提取结果。
deeplab是谷歌在fcn的基础上设计出来的一种算法。fcn为了得到一个更加密集的位置敏感图，将一张500x500的输入图像，直接在第一个卷积层上加了一个100的填充，最终在第七卷积层得到一个16x16的位置敏感图。deeplab这里使用了一种更加巧妙的做法：将vgg网络的14和15池化层的步幅由原来的2改为了1，再加上1个补白。就是这样一个改动，使得vgg网络总的步幅由原来的32变成8，进而使得在输入图像为514x514时，第七卷积层能得到67x67的位置敏感图，与fcn相比更加密集。
deeplab-v3仍然主要在空洞卷积进行探究，设计采用多比例的带孔卷积级联或并行来捕获多尺度背景，利用deeplab-v3网络训练获得的建筑物提取模型，通过多个卷积传递图像，将卷积均设置步长为1，使用相同填充，而非使用池化和全连接层，通过这种处理，使每一次卷积都保留了输入的空间维度，通过多个卷积堆叠，并最终得到一个用于进行建筑物的提取的分割模型(即建筑物提取模型)。
步骤s280，根据各特征指数从处理后的卫星影像中进行阴影提取，获得阴影提取结果。
其中，各特征指数是阴影提取指数、绿通道植被指数(gndvi)、影像第四波段等指数。通过将各特征指数的阴影表现力做对比，选用阴影表现力最为清晰的特征指数进行阴影提取工作，在本实施例中，阴影表现力最为明显的特征指数是影像第四波段，因此采用影像第四波段进行阴影提取工作，简单易行实用性高。但由于提取阴影的方法是基于像元的，所以在提取结果当中存在一定的椒盐现象，导致提取到的初步阴影提取结果存在小图斑，与此同时，由于楼顶并不平整，有部分诸如楼梯间的小型建筑，所以会产生一部分小型阴影，因此要对小图斑和小型阴阴影进行剔除或重新分类，本申请采用聚类处理来进行阴影提取结果分类后处理。如图3所示的部分处理后的卫星影像示意图，处理后的卫星影像中有被太阳照射的阴影，根据各特征指数从处理后的卫星影像中进行阴影提取，获得如图4所示的部分阴影提取结果示意图。
在一个实施例中，根据各特征指数从处理后的卫星影像中进行阴影提取，获得阴影提取结果的步骤，包括：
将各特征指数的阴影表现力的强度进行对比，确定阴影表现力的强度最高的特征指数作为阴影提取的特征指数；采用阴影提取的特征指数进行阴影提取工作，获得初步阴影提取结果；对初步阴影提取结果采用聚类处理进行优化，获得阴影提取结果。
其中，阴影表现力的强度最高的特征指数是阴影表现力最为清晰明显的特征指数。采用聚类处理(clump)对提取后的阴影进行优化。聚类处理(clump)是运用数学形态学算子(腐蚀和膨胀)，将临近的类似分类区域聚类并进行合并。分类图像经常缺少空间连续性(分类区域中斑点或洞的存在)。低通滤波虽然可以用来平滑这些图像，但是类别信息常常会被临近类别的编码干扰，而聚类处理解决了这个问题，聚类处理首先将被选的分类用一个膨胀操作合并到一块，然后用变换核对分类图像进行腐蚀操作，从而避免了类别信息常常会被临近类别的编码干扰的问题，因此，采用聚类处理(clump)对提取后的阴影进行优化可以提高阴影提取结果的准确性。
步骤s300，根据太阳方位线和阴影提取结果进行分析，获得阴影提取结果中各阴影的阴影长度。
其中，根据卫星影像数据，绘制太阳方位线，根据太阳方位线与阴影提取结果中的各个阴影的交集，获得阴影提取结果中各阴影的阴影长度。
步骤s320，根据建筑物提取结果和阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定建筑物提取结果中各建筑物的建筑物高度。
其中，根据建筑物提取结果和阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，针对太阳、卫星、建筑物三者之间不同的位置关系分四种情况建立相应模型进行模拟，并利用太阳方位角、太阳高度角、卫星高度角、卫星方位角、建筑物走向线等数据，按照建筑物高度反演模型进行计算完成高度反演，并对建筑物高度对应赋予到建筑物提取结果的各建筑物上，确定建筑物提取结果中各建筑物的建筑物高度。
在一个实施例中，根据建筑物提取结果和阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定建筑物提取结果中各建筑物的建筑物高度的步骤，包括：
根据建筑物提取结果，确定待进行建筑物高度反演的建筑物；基于太阳、卫星、建筑物三者之间的空间位置关系，以及建筑物的建筑物高度和阴影的关系，利用建筑物高度反演模型进行建筑物高度的反演工作，确定建筑物提取结果中各建筑物的建筑物高度。
其中，利用建筑物高度反演模型进行建筑物高度的反演工作，是基于太阳、卫星、建筑物三者之间的空间位置关系，建筑物高度和阴影的关系，该建筑物高度和阴影的关系可以分为太阳与卫星方位相同与不同、太阳与卫星在建筑物的异侧与同侧共四种情况。理想情况下太阳与卫星在垂直平面内位于同一平面时，当太阳与卫星的方位和高度角都相同，或太阳与卫星的方位相同而太阳高度角大于或等于卫星高度角时，这种情况下没有阴影。
根据不同情况建立了不同的高度反演模型来实现建筑物的精准反演，共包含以下高度反演公式：
在一个实施例中，当太阳和卫星的方位相同，而太阳的高度角小于卫星的高度角，同时太阳射线在地面的投影与建筑物走向垂直时，建筑物高度反演模型的高度反演公式为：
h＝s·(tanα·tanβ)/(tanα-tanβ)
其中，α是卫星的高度角，β是太阳的高度角，s为建筑物在阴影提取结果中的阴影长度，h是建筑物的建筑物高度。
如图5所示，其中，当太阳和卫星的方位相同，而太阳的高度角小于卫星的高度角，同时太阳射线在地面的投影与建筑物走向垂直时，太阳、卫星、建筑物和阴影的关系如下：
s＝l2-l1＝h/tanβ-h/tanα(7)
式中，α是卫星的高度角，β是太阳的高度角，l2是建筑物实际的阴影长度，l1为被建筑物遮挡的阴影长度，s为建筑物在阴影提取结果中的阴影长度，h是建筑物的建筑物高度。一旦获取到了建筑物在阴影提取结果中的阴影长度s，就可以依照下式推算出建筑物的建筑物高度：
h＝s·(tanα·tanβ)/(tanα-tanβ)(8)
若卫星拍摄时为正视姿态时，卫星将可以观测到建筑物的整个阴影，不存在有卫星观测不到的阴影部分，因此l1＝0，l2＝s，此时
h＝s·tanβ(9)
如图6所示，当太阳和卫星位于建筑物异侧，同时太阳射线在地面的投影与建筑物走向垂直，此时，可以看到建筑物的全部落影而且建筑物有本影，则：
建筑物阴影长：
l2＝s＝h/tanβ(10)
建筑物高度：
h＝s·tanβ(11)
以上两种情况只是考虑到卫星、建筑物、太阳三者位于同一竖直平面内时，卫星和太阳位于建筑物异侧或者同侧两种特殊的情况，实际上二者是有一定的夹角的，而且建筑物走向与太阳射线在地面的投影也不一定垂直。考虑到太阳和卫星方位以及二者与建筑物走向的关系，建筑物高度的计算可以利用下面的非理想情况：
具体地：当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线在地面的投影和卫星视线在地面的投影位于建筑物的两侧时，建筑物高度反演模型的高度反演公式为：
h＝s·tanβ/sinγ
其中，γ为建筑物走向与太阳射线投影方向的夹角。
如图7所示，当太阳射线在地面的投影和卫星视线在地面的投影位于建筑物的两侧时，一般情况下绝大部分阴影(包括本影和落影)能够清晰地显示在遥感影像上，此时，用于建筑物高度计算的阴影只是全部的落影，几何关系如下：
sinγ＝l2/s(式拾贰)
h＝l2·tanβ(式拾叁)
式中s为建筑物在阴影提取结果中阴影长度，l2为建筑物实际的阴影长度，h为建筑物高度，γ为建筑物走向与太阳射线投影方向的夹角。
h＝s·tanβ/sinγ(式拾肆)
当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线投影和卫星视线投影位于建筑物同侧时，建筑物高度反演模型的高度反演公式为：
h＝s·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
其中，δ为建筑物走向与卫星视线投影方向的夹角。
如图8所示，当太阳射线投影和卫星视线投影位于建筑物同侧，建筑物的向光面成像时会使部分阴影被遮挡，从而使得卫星影像上的可视阴影只是建筑物本影的一部分，几何关系如下：
h＝s1·tanβ＝s2·tanα(15)
l2＝s1·sinγ(式拾陆)
l1＝s2·sinδ(式拾柒)
其中，δ为建筑物走向与卫星视线投影方向的夹角，s1为太阳投射方向阴影长度，s2为卫星投射方向阴影长度，又因为
s＝l2-l1(式拾捌)
所以
h＝s·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)(式拾玖)
式中：s是建筑物在阴影提取结果中的阴影长度(即卫星影像上垂直于建筑物走向线的可视阴影长度)。
式中：a1是太阳方位角，a2是卫星方位角，(x1，y1)、(x2，y2)为建筑物走向线或其平行线上两点的坐标。
步骤s340，将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
以基于中心经纬度为(120.3e，36.2n)、(120.5e，36.1n)、(120.2e，36.0n)、(120.5e，35.9n)，(120.3e，36.4n)、(120.6e、36.3n)的六幅青岛市2016年8月底9月初附近的高分二号影像为待分析区域的卫星影像，对卫星影像进行预处理，获得处理后的卫星影像，对卫星影像进行预处理，获得处理后的卫星影像，根据各特征指数从处理后的卫星影像中进行阴影提取，获得阴影提取结果，根据太阳方位线和阴影提取结果进行分析，获得阴影提取结果中各阴影的阴影长度，根据建筑物提取结果和阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定出的建筑物提取结果中各建筑物的建筑物高度，与各建筑物的实际建筑物高度进行对比，得出建筑物高度反演精度评价表1.1。
由此可以看出，本申请在可以较为精准的完成建筑物高度反演工作，得出较为精准的建筑物高度。
在一个实施例中，将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括：根据建筑物高度与建筑物影响深度的转换关系表，将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度；根据各建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
在一个实施例中，根据各建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括:根据各建筑物的建筑物影响深度，分别以建筑物影响深度为0m-10m、10m-30m、30m-50m和50m-100m为计算单元；对各计算单元进行地下空间资源量利用情况叠加分析，获得地下空间资源量分析结果。
其中，按照建筑物高度与建筑物影响深度的转换关系表1.2，将各建筑物的建筑物高度转换为建筑物影响深度，然后利用建筑物影响深度和面积以0m-10m、10m-30m、30m-50m、50m-100m共4层为计算单元，开展叠加分析进行分层的地下空间资源量利用情况统计计算，获得地下空间资源量分析结果。
表1.2：建筑物高度与建筑物影响深度的转换关系表
具体地，对于地下空间可利用资源量进行叠加计算，城市地下空间资源量在不同深度条件下，有不同的用途。在可利用资源量叠加分析时，需要考虑到这一点。拟以0m-10m、10m-30m、30m-50m、50m-100m共4层为计算单元，开展叠加分析，计算可利用资源量。以0m-10m地下空间可利用资源量叠加计算方法为例：计算方法为0m-10m总资源量扣除10m已利用资源量以及深度超过10m地下空间资源已利用的0m-10m部分，其他部分以此类推。经过上述步骤，得到0m-10m、10m-30m、30m-50m、50m-100m共4层计算单元的空间资源利用量分布图，最终获得地下空间资源量分析结果。
0m-10m地下空间可利用资源量叠加计算方法为0m-10m总资源量扣除影响深度10m已利用资源量以及影响深度超过10m地下空间资源已利用的0m-10m部分，计算公式如下：
k0-10＝(s10 s30 s50 s100)*10(式贰拾贰)
v(0-10m可利用资源量)＝v(0-10m总资源量)-v(0-10m已利用资源量)-k0-10(式贰拾叁)
式中：k0-10是深度超过10m地下空间资源已利用的0m-10m部分，s10是影响深度10m的地块面积，s30是影响深度30m的地块面积，s50是影响深度50m的地块面积，s100是影响深度100m的地块面积，v(0-10m可利用资源量)是0m-10m的地下空间可利用资源量，v(0-10m总资源量)是0m-10m的地下空间总资源量，v(0-10m已利用资源量)是0m-10m的地下空间已利用资源量。
10m-30m地下空间可利用资源量叠加计算方法为：10m-30m总资源量扣除影响深度超过10m地下空间资源已利用的10m-30m部分，计算公式如下：
k10-30＝(s30 s50 s100)*20(式贰拾肆)
v(10m-30m可利用资源量)＝v(10m-30m总资源量)-k10-30(式贰拾伍)
式中：k10-30为深度超过10m地下空间资源已利用的10m-30m部分，v(10m-30m可利用资源量)是10m-30m的地下空间可利用资源量，v(10m-30m总资源量)是10m-30m的地下空间总资源量。
30m-50m地下空间可利用资源量叠加计算方法为：30m-50m总资源量扣除影响深度超过30m地下空间资源已利用的30m-50m部分，计算公式如下：
k30-50＝(s50 s100)*20(式贰拾陆)
v(30m-50m可利用资源量)＝v(30m-50m总资源量)-k30-50(式贰拾柒)
式中：k30-50为深度超过30m地下空间资源已利用的30m-50m部分，v(30m-50m可利用资源量)是30m-50m的地下空间可利用资源量，v(30m-50m总资源量)是30m-50m的地下空间总资源量。
50m-100m以下地下空间可利用资源量叠加计算方法为：50m-100m总资源量扣除影响深度超过50m地下空间资源已利用的50m-100m部分，计算公式如下：
k50-100＝s100*50(式贰拾捌)
v(50-100m可利用资源量)＝v(50-100m总资源量)-k50-100(式贰拾玖)
式中：k50-100为深度超过50m地下空间资源已利用的50m-100m部分，v(50-100m可利用资源量)是50m-100m的地下空间可利用资源量，v(50-100m总资源量)是50m-100m的地下空间总资源量。
上述地下空间资源量分析方法，通过获取待分析区域的卫星影像；对所述卫星影像进行预处理，获得处理后的卫星影像；采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果；根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果；根据太阳方位线和所述阴影提取结果进行分析，获得所述阴影提取结果中各阴影的阴影长度；根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度；将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果，解决了实地统计已利用资源量，调查成本高的问题。
利用高分2号遥感影像，采用建筑物指数和深度学习相结合的办法进行了待分析区域的建筑物提取；利用高分影像第四波段阈值分割方法，能够简单有效的将大部分建筑物阴影提取出来，非常适用于实际的工作生产当中；建筑物高度反演方面仔细分析了不同情况下的高度反演模型，能够具体准确的完成高度反演工作；城市地下空间资源量在不同深度条件下，有不同的用途，考虑到这一点对地下空间可利用资源量叠加分析时进行了分层分析。最后通过将建筑物提取结果与建筑物影像深度结果相结合完成了待分析区域的地下空间资源量计算工作。
应该理解的是，虽然图1的流程图中的各个步骤按照箭头的指示依次显示，但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明，这些步骤的执行并没有严格的顺序限制，这些步骤可以以其它的顺序执行。而且，图1中的至少一部分步骤可以包括多个子步骤或者多个阶段，这些子步骤或者阶段并不必然是在同一时刻执行完成，而是可以在不同的时刻执行，这些子步骤或者阶段的执行顺序也不必然是依次进行，而是可以与其它步骤或者其它步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
以上实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。
以上所述实施例仅表达了本申请的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。因此，本申请专利的保护范围应以所附权利要求为准。
1.一种地下空间资源量分析方法，其特征在于，所述方法包括：
获取待分析区域的卫星影像；
对所述卫星影像进行预处理，获得处理后的卫星影像；
采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果；
根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果；
根据太阳方位线和所述阴影提取结果进行分析，获得所述阴影提取结果中各阴影的阴影长度；
根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度；
将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
2.根据权利要求1所述的方法，其特征在于，所述采用形态学建筑指数和深度学习从所述处理后的卫星影像中进行建筑物的提取，获得建筑物提取结果的步骤，包括：
采用形态学建筑指数从所述处理后的卫星影像中进行建筑物的提取，获得第一提取结果；
采用深度学习从所述处理后的卫星影像中进行建筑物的提取，获得第二提取结果；
将所述第一提取结果和所述第二提取结果进行融合，获得融合后的提取结果；
对所述融合后的提取结果进行楼顶形状优化，获得建筑物提取结果。
3.根据权利要求1所述的方法，其特征在于，所述根据各特征指数从所述处理后的卫星影像中进行阴影提取，获得阴影提取结果的步骤，包括：
将各特征指数的阴影表现力的强度进行对比，确定阴影表现力的强度最高的特征指数作为阴影提取的特征指数；
采用所述阴影提取的特征指数进行阴影提取工作，获得初步阴影提取结果；
对所述初步阴影提取结果采用聚类处理进行优化，获得阴影提取结果。
4.根据权利要求1所述的方法，其特征在于，所述根据所述建筑物提取结果和所述阴影提取结果中各阴影的阴影长度，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度的步骤，包括：
根据所述建筑物提取结果，确定待进行建筑物高度反演的建筑物；
基于太阳、卫星、所述建筑物三者之间的空间位置关系，以及所述建筑物的建筑物高度和阴影的关系，利用建筑物高度反演模型进行建筑物高度的反演工作，确定所述建筑物提取结果中各建筑物的建筑物高度。
5.根据权利要求4所述的方法，其特征在于，当太阳和卫星的方位相同，而太阳的高度角小于卫星的高度角，同时太阳射线在地面的投影与建筑物走向垂直时，所述建筑物高度反演模型的高度反演公式为：
h＝s·(tanα·tanβ)/(tanα-tanβ)
其中，α是卫星的高度角，β是太阳的高度角，s为建筑物在阴影提取结果中的阴影长度，h是建筑物的建筑物高度。
6.根据权利要求4所述的方法，其特征在于，当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线在地面的投影和卫星视线在地面的投影位于建筑物的两侧时，所述建筑物高度反演模型的高度反演公式为：
h＝s·tanβ/sinγ
其中，γ为建筑物走向与太阳射线投影方向的夹角。
7.根据权利要求4所述的方法，其特征在于，当建筑物走向与太阳射线在地面的投影不垂直，且太阳射线投影和卫星视线投影位于建筑物同侧时，所述建筑物高度反演模型的高度反演公式为：
h＝s·sinα·sinβ/(sinα·cosβ·sinγ-cosα·sinβ·sinδ)
其中，δ为建筑物走向与卫星视线投影方向的夹角。
8.根据权利要求1所述的方法，其特征在于，所述将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度，利用各建筑物的建筑物影响深度和面积进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括：
根据建筑物高度与建筑物影响深度的转换关系表，将各建筑物的建筑物高度转换为各建筑物的建筑物影响深度；
根据各所述建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果。
9.根据权利要求8所述的方法，其特征在于，所述根据各所述建筑物的建筑物影响深度和面积，进行分层的地下空间资源量利用情况分析，获得地下空间资源量分析结果的步骤，包括:
根据各所述建筑物的建筑物影响深度，分别以建筑物影响深度为0m-10m、10m-30m、30m-50m和50m-100m为计算单元；
对各所述计算单元进行地下空间资源量利用情况叠加分析，获得地下空间资源量分析结果。
技术总结