利用克立格插值法计算土方的实例研究.docxVIP

利用克立格插值法计算土方的实例研究 土壤处理是土壤处理的主要项目之一。其主要目标是通过填充土壤，使平坦区域的平面位置和高度平整，以满足设计要求。计算土方的常用方法有格网法、断面法或等高线法等。格网法模型简单,易于实现,比较适于手算。但计算烦琐,工作量大且精度不高;利用地形图进行土方计算,势必会受到地形图精度及其现势性的影响,当精度要求较高时,则需现场实测格网点高程、断面图或更大比例尺的地形图,野外工作量较大。基于此,尝试运用地统计学中插值的理论方法,结合Esri公司的Arcview软件来进行土方量的计算,并对其精度进行分析。 1 空间自然现象的科学保证 地统计学是以区域化变量理论为基础,以变异函数为主要工具,研究那些在空间分布上既有随机性又有结构性,或空间相关和依赖性的自然现象的科学。简单地说,它主要完成两个方面的任务:利用半方差函数生成半方差图来量化研究对象的空间结构;通过插值方法利用半方差图中拟合模型和研究对象周围的实测值来对未知值进行预测。 1.1 半方差图的模拟 半方差函数是用来描述区域化变量结构性和随机性并存这一空间特征而提出的。在满足假设的条件下,半方差函数公式为: r(h)=12Ν(h)Ν(h)Σi=1[Ζ(xi)-Ζ(xi+h)]2(1)r(h)=12N(h)Σi=1N(h)[Z(xi)?Z(xi+h)]2(1) 式中:r(h)——所有空间相距h的点队的平均方差;N(h)——在空间上具有相同间隔距离h的离散点对数目;Z(xi)、Z(xi+h)——xi、与xi距h的点的某一属性观测值。将以r(h)为纵轴,h为横轴绘制出的r(h)随h的增加而变化的曲线,称为半方差图,用来表示变量在研究区域内的空间变异和相关程度。图1即为高程关于滞后距离h的半方差图。图中对半方差值的模拟采用指数模型(见公式2),从中可以看出半方差函数的几个重要参数:块金系数C0=0.0059,它是h=0时的半方差函数值;基台值C+C0=0.0139,它是拟合的曲线趋于平稳时的半方差值,C为拱高。块金系数与基台值的比值C0/C+C0=42%,说明这些高程值有中等程度的空间相关性,半方差值呈现抛物线形状表示这些高程值在空间上有非常好的连续性;变程a=138,是空间内插的极限距离,对整理区高程的插值要在该范围内进行才具有意义。指数模型的一般公式为: r(h)={0,h=0C0+C(l-eha)?h＜0(2)r(h)={0,h=0C0+C(l?eha)?h＜0(2) 1.2 实测高程点权重的确定 有大量的方法供地统计进行空间插值计算,采用普通克立格法(0rdinarykriging),根据待插值点与临近实测高程点的空间位置,对待插值点的高程值进行线性无偏最优估计,通过生成一个关于高程的克立格插值图来表达研究区域的原始地形。总的公式是: Ζ(x0)=ΝΣi=1λiΖ(xi)(3)Z(x0)=Σi=1NλiZ(xi)(3) 式中:Z(x0)——待预测的高程值;Z(si)——在60个点处采集的实测高程值;N=60;λi是分配给每个实测高程点的权重。它的确定是通过半方差图分析获取的,根据统计学上无偏和最优的要求,利用拉格朗日极小化原理,可推导出权重值和半方差之间的公式: A?[λθ]=B(4) 其中 A=[r11r12?r1Ν1r21r22?r2Ν1?????rΝ1rΝ2?rΝΝ111?10]?B=[r10r20?rΝ01]?[λθ]=[λ1λ2?λΝθ] 矩阵A中的r11……rNN是实测高程值之间的半方插值,距阵B中的r10……rN0为实测高程点xi和内插点x0之间的半方差值,θ是拉格郎日乘子。由(4)式解出λi,然后代入(3)式就可获得点xo上的内插估计高程值Z(xo)。 由于克立格法不仅考虑了落在该样点的数据,还考虑了临近样点的数据;不仅考虑了待估样点与临近实测样点的空间关系,而且还考虑了各临近实测样点彼此间的位置关系,因此它在地理统计意义上优于其他插值方法(象逆距离加权法,多项式插值法,趋势面法,样条插值等)。 一般来说,对参与克立格插值计算的实测样点的数目要有所限制,因该值的大小选择会影响插值的计算结果。若太大,则内插结果过于平滑;太小,则无法反映地表的变化;距离预测点较远的实测点可能与待估样点已经不存在自相关关系,也不能参与插值计算。采取以插补点为圆心,以R为半径的圆来确定取样的范围和参加计算的实测样点数目(如果存在各向异性,则可考虑划定一椭圆作为研究区域)。为了避免方向上的偏差,将圆平均地分为4个扇区,每个扇区内实测点数目在2到5之间,这样总共参与每个待估点预测的实测点数目平均达到8个。 2 计算地表高、填方最小二乘 在施工中,为了估计土方工程的工程量和费用,须根据设计平面位置对土方工程量进行计算,土方工程量实际上就是原始地表与设计地表之间的体积值