[发明专利]一种基于GPU加速的多分辨率大规模森林演替过程仿真方法

摘要

一种基于GPU加速的多分辨率大规模森林演替过程仿真方法，包括以下步骤：1)，建立树木数据组织与样地数据组织；根据GPU架构对树木数据进行组织，形成利于GPU计算的数数据组织结构，对于整个样地单元建立一与样地网格大小相同的二维数组结构M；2)将每棵树木的遍历转换成cuda内核中单位线程，计算plot[k]的种子分布情况，Y纬度值代表样地，对线程网格的纬度Y＝k处，所有X纬度上的线程进行对样地K的种子分布影响计算；3)种子分布内核方法，成年树NCI方法GPU的内核方法以及成年树光照GLI内核方法。本发明提供一种正确性良好、提升运行效率的基于GPU加速的多分辨率大规模森林演替过程仿真方法。

主权项

一种基于GPU加速的多分辨率大规模森林演替过程仿真方法，其特征在于：所述仿真方法包括以下步骤：1)，建立树木数据组织与样地数据组织；根据GPU架构对树木数据进行组织，形成利于GPU计算的数数据组织结构：所有计算的树木数据结构将被设计成数组形式；在多分辨率数据组成中聚类结点第k层与第k+1层的数据表示，对于第k层单位结点数据由m个k‑1层结点数据累加而成，对于每个单位结点数据，由s种树种对应单位成年树结点数据构成；每棵成年树或者幼年树x的树冠被建模成圆柱体，它们在地表的投影都是一个圆C(x)，设样地单元边长为r，则样地单元网格的面积即r2，对于投影完全包含k个样地单元网格的树木，则其暴露面积即为kr2，为正确计算不同树木发生遮挡时暴露面积，设树木树冠上顶面离地面距离为H(x)，所有样地网格单元标记为Mij，在理想不想交投影区域，则如果某样地网格Mij完全包含在树木投影C(x)中，则将Mij分配给树木x，否者不对其分配；在投影相交区域，对于样地网格Mij被多个投影所覆盖，该样地网格将会被分配给能够完全包含它的H(x)最大的那棵树木x；对于整个样地单元建立一与样地网格大小相同的二维数组结构M，所有元素所记录的是所被分配的成年树ID或者未被分配标记；对于具有n棵树木的森林场景，首先遍历树木，遍历中对树木x投影所包含的所有样地网格单元M(i，j)进行分配，如果该网格单元未分配给任何树木，则将其分配给x；如果该网格单元已经被分配给树木t，则通过x树木ID(x)获取树高H(x)，通过记录数组获取t的树木ID(t)，利用ID(t)获取树高H(t)，比较H(t)与H(x)大小，若H(x)大于H(t)，则将该样地单元网格M(i，j)重新分配给树木x，否则不做任何 变化；该方法的最大时间复杂度为O(kn)，k为单棵树木最大可能覆盖的样地网格，在1m*1m的样地网格；2)线程分配设计：根据模型公式(1)： R i = 1 η STR Σ k = 1 N ( ( DBH k 30 ) β e - Dd ik θ ) - - - ( 1 ) 其中，Ri表示样地i的种子数量，μ，STR，D，θ，β都是与树种相关的常数参数，DBH为成年树树干胸径大小，由成年树胸径的大小及树种参数所决定，dik表示第k棵成年树距离样地i的实际物理距离，N为所有具有繁殖能力的成年树总数，通过对树种真实参数的代入可知，随着母树与样地的距离的dik不断增大，母树所对样地产生的种子数量Ri成单调递减趋势；将每棵树木的遍历转换成cuda内核中单位线程，计算plot[k]的种子分布情况，Y纬度值代表样地，对线程网格的纬度Y＝k处，所有X纬度上的线程进行对样地K的种子分布影响计算；对于样地plot(i，j)的种子分布计算，线程网格Y纬度为i+j*PlotWidth的计算其样地的种子分布情况，其中PlotWidth为样地的宽度，对于计算同一样地种子分布的线程划分，是根据邻域树木与样地的距离的不同，线程所计算的结点层次也不同；3)种子分布内核方法具体流程如下：a获取线程所在块的Y纬度k，，获取线程在X维度上的值；b判断该线程所计算结点与线程所计算目标样地单元的距离d，根据不同距离尺度标准从而进行进行判断是否进行细分结点，第二层结点距离尺度标准为ρ2，如果d大于ρ2，则不进行划分，如果小于ρ2，则进行细分结点；c如果结点为叶子结点或不需要细分结点，则直接获取线程所对应母结点聚类数据；通过上一步所获得的母结点与目标样地距离进行种子密度的计算；d如果细分结点，则该线程将获取母聚类结点所对应子节聚类点数据进行计算，重复上面步骤b；成年树NCI方法GPU的内核方法，根据模型公式(6)： NCI = Σ i = 1 S Σ j = 1 N λ i ( DBH ij ) α ( dis tan ce ij ) β - - - ( 6 ) 计算公式中，i＝1...S为树种数，j＝1...N领域内大于某一固定DBH的树木；α是由目标树树种决定的NCIAlpha参数值；β是由目标树树种决定的NCIBeta参数值；DBHjk是领域树j的DBH值；λik是树种j相对于目标树种的NCILambda参数；distanceik是目标树与领域树之间的距离；成年树NCI方法GPU的内核方法具体流程如下：a获取线程所在块的Y纬度k，获取线程在X维度上的值(即线程所计算多分辨率数据根结点层的索引)；b判断该线程所计算结点与线程所计算目标树木的距离d，根据距离判断该结点是否在目标树木的邻域影响圈内；c如果在影响圈内，则直接获取最底层叶子结点数据，精确判断叶子结点是否在目标树木的影响圈，如果是，则根据模型计算其NCI值，如果不是，则不计算，d如果线程所对应结点超出影响圈，则该线程任务结束，成年树光照GLI内核方法具体流程如下：a数据初始化：分配线程，将线程划分为k批进行处理，每批m个线程，线程总量为k*m，为成年树与幼年树木数量，动态建立实际分辨率所需样地单元网格二维数组记录结构assign_tree，对其所有元素初始化为‑1(未分配任何树木ID)，初始化所有成年树暴露面积A为0，获取整个森林场景中成年树与幼年树DBH，位置坐标，树种等属性数组，b获取线程号thread_id与块号block_id，通过thread_id与block_id索引树木 索引tree_id，对于树木tree_id为i的树木，获取其DBH属性DBH(i)，树木所在位置L(i)，并通过模型公式由DBH(i)计算其高度H(i)，树冠半径R(i)，c通过L(i)获取树木i所在样地单元PID(i)，以所在样地单元坐标PLocation(PID)为中心，遍历(x，y)坐标在x∈(L(PID)‑R(i)，L(PID)+R(i))，y∈(L‑R(i)，L(PID)+R(i))内的所有样地网格，若该范围内的样地单元网格k距离L(i)最远顶角坐标小于R(i)，则该样地网格完全包含于树木i的投影圆内，进入步骤4，否者不作任何处理，d样地单元网格k完全包含于投影圆内，通过记录数组获取网格k所被分配的成年树ID：assign_tree(k)，若该样地网格未被分配，则ID为‑1，此时只需将assign_tree(k)赋值为i，并对A(i)累加该样地单元面积r2，若该样地网格已被分配，则获取该样地被分配成年树ID的树高H(assign_tree(k))，比较H(assign_tree(k))与H(i)大小，若H(assign_tree(k))大于等于H(i)，则不做任何处理，否者，将assign_tree(k)赋值为i，并对A(i)累加该样地单元面积r2，并对A(assign_tree(k))减去该样地单元面积r2，e继续进入步骤c遍历阴影范围内的样地网格，若遍历完，则进入步骤b进行下一棵成年树，直至结束该方法。