surfer7讲义.docx

1、surfer7讲义第六章 三维绘图软件Surfer76.1 网格文件及其建立6.1.1简介 网格文件规则矩形网格的三维数组X，Y，Z文件。用来画等值线图、网框图、矢量图、粘贴图、立体阴影图、色谱图，计算面积、体积和画剖面图。.GRD文件格式.GRD文件：GS二进制文件(缺省网格文件格式)、GS ASCII文件.DAT文件：ASCII XYZ数据文件 . Surfer7 还可以直接使用USGS .DEM, GTopo30 .HDR和SDTS .DDF 文件。6.1.2 由离散点的XYZ文件创建.GRD 文件为了产生等值线图或线网图, 或执行需要网格文件的操作如网格代数处理、计算体积和面积, 光滑

2、或计算残差，必须预先对数据文件进行网格化处理。由GridData 打开所需的数据文件后，出现散点数据插值法对话框Scattered Data Interpolationd。 该对话框由4个Tab组成。 1 Date数据列组（Data Columns） 设定X、Y和Z的值对应于数据文件内的哪一列。相同数据点处理方式选择。由下拉菜单Duplicates可以选择对有相同的X、Y坐标但有不同的Z值的两个或多个点Z值处理办法，包括取平均值（Average）、保留一个并删除其它值(Delete)、对各点Z值加和(Sum)等。有些网格化方法还可以选择忽略这类数据点(Ignore)。数据筛选（Filter）

3、可以使用布尔表达式对原始数据点进筛选。原始XYZ数据对话框提供数据统计。显示在数据文件内指定列中，相应的X, Y和Z值的数据点数、最小值、最大值、平均值、方差和标准差等。2 Genaral网格线参数（Grid Line Geometry） 选择X,Y方向网格的取值范围；设定产生网格文件时网格间隔和网格线的密度。改变网格线的密度，网格的间隔宽度会自动随着改变。网格线的密度高则较精密，但文件长度大，构筑网格费时间， 网格化方法（Gridding Method）设定采用网格化的方法。有下拉式菜单提供多种可用的网格化方法清单。9种网格化的数学方法见本章第二节。用不同的网格化方法处理数据时会产生不同的结

4、果。选择时，首先应使你的数据最接近实际情况。其次应使表现数据的图形为最好。Option 每种网格化方法都有相应的参数设定。 输出网格文件（Output Grid File）指定生成的网格文件的路径和文件名。缺省名与数据文件同名，扩展名.GRD，可点击浏览按钮以改变其路径或文件名并选择网格文件类型。3 Search搜寻功能决定在内插网格结点时，那些数据点参加插值计算。对于要用到全部数据文件的网格化方法，例如多项式回归法、最小曲率法和线性插值的三角形法等，搜寻功能不起作用，Search键不能激活。有以下4种搜寻类型：全部数据（All Data）：每一个结点Z值的计算都要用到全部数据点。对于较少的数

5、据组，如少于250组数据时，最常用到全部数据搜寻。简单搜寻(Sample)：只用最靠近某结点的一些数据来计算该结点的Z值。数据点的数目由窗口右面Search Rules框内的Data Per Sector来设定。只有在搜寻椭圆（Search ellipse）以内的数据点才被用来计算。象限搜寻（Quadrant）:将一个结点周围分为4个90的象限，软件自动分析每个象限的数据点数目，搜寻4个扇区内最近的数据点，而忽略掉该象限中这些最近点以外的数据点。卦限搜寻（Octant）: 将一个结点周围分为8个45的卦限，软件自动分析每个卦限的数据点数目，搜寻8个扇区内最近的数据点。较象限搜寻慢，但构成比较平

6、滑的网格。当观察点明显成簇状不均匀分布时，推荐使用象限或卦限搜寻方法，以消除或减少简单搜寻所造成的失真。搜寻椭圆（Search Ellipse）：在搜寻方式为简单搜寻、象限搜寻或卦限搜寻情况下，计算每个结点Z值时搜寻数据点的范围。搜寻椭圆以外的数据点，进行插值时不予考虑。缺省为圆形，表示在所有方向搜寻距离相同。搜寻椭圆的参数有：半径：半径1和半径2分别是用数据单位表示的X和Y方向的搜寻距离。角度：+X轴和半径1之间的倾斜角度（-360+360）。结点Z值为0说明在特定位置，根据用户设定的搜寻规律，没有足够的数据点来计算结点Z值。搜寻规律（Search Rules）：在搜寻方式为简单搜寻、象限搜

7、寻或卦限搜寻情况下，设置用来计算每个结点Z值的数据点数。与搜寻椭圆相关。每扇区点数：设置每扇区的数据点数，对于简单搜寻、象限搜寻和卦限搜寻，扇区数分别是1、4、8，缺省每扇区点数分别为24、8、3。最少总点数：对一个网格结点进行插值时，最少所需数据点数。最多空扇区：对一个网格结点进行插值时所容许的最多空扇区数。 4 Faults and Breaklines Faults 用.BLN文件格式保存的X、Y二维多段线或多边形，用来限定网格化的范围。如对多段线一侧的原始数据点进行网格化，忽略另一侧。或对一多边形之内（或之外）的原始数据点进行网格化。File Containing Fault Trac

8、es 框给出Faults 的.BLN文件所在的路径。 有一些网格化模型可以使用Faults。 Breaklines 用.BLN文件格式保存的三维数据集，X、Y、Z给出每个点的坐标。不同于Faults，Breaklines不阻断信息流。当Grid算法遇到Breaklines时，用线形内插的方法联合确定结点的Z值。File Containing Breaklines框给出Breaklines的.BLN文件所在的路径。 Breaklines用于定义流线、山脊和其他斜率的突变。6.1.3 由函数关系生成.GRD文件GridFunction 打开网格函数对话框。在函数编辑框中输入数学函数Z= f(X,Y

9、)。数学函数逐点的计算每一个Z值，并将其写入到输出的网格文件中。其中网格结点的数目取决于所选定的 X和 Y的最小值、 最大值和增量。所采用的数学函数中能使用系统 的任何内部函数。例如, 按给定的函数关系，X、Y的取值范围和增量产生一个网格文件。缺省的输出网格文件名为OUT.GRD，可改变文件的路径和名字。6.1.4 网格代数处理由具有相同X,Y网格尺度的两个网格文件按照一定的函数关系C= f(A,B) 生成一个新的输出网格文件。这里A和B是输入的网格文件，C是输出网格文件。从输入网格文件的相应结点上(具有相同X,Y坐标的网格结点) 计算函数值，计算的结果送到具有相同坐标的输出网格文件内。网格数

10、学也能用于单一的网格文件。例如,使用函数C= log10(A)，则网格文件C的每一个结点的值将为网格A相应的每一网格结点值的以10 为底的对数。也可以进行简单的数学操作,例如 C=A-100。这时输出网格文件C比输入网格文件A低100 单位。由GridMath打开网格代数对话框。输入网格文件A组指定执行网格数学操作的第一个网格文件。输入网格文件B组指定用于网格数学操作的第二个网格文件。如果只想在一个网格文件上执行操作, 可使其为空白。第二个网格文件应包含与第一个网格文件相同的网格行、列数, 和相同的X和Y列。缺省的输出网格文件名为MATHT.GRD，可改变文件的路径和名字。底部的编辑框用于输入

11、函数C= f(A,B)。 6.1.5 网格算法（Grid Calculus）提供解析网格文件的工具，用来解释一些在等值线图和线网图中不明显的网格文件的特点。由GridCalculus打开网格算法对话框。选择网格算法进行处理，生成一个新的网格文件。由网格算法对话框可以选择以下4种算法之一：方向导数、地体模型、微分和积分算符以及Fourier和谱分析1） 方向导数（Directional Derivative）方向导数给出沿着规定的方向网格表面的斜率或斜率变化率的信息。斜率是有方向性的，如果规定方向为东，但倾斜方向为北，则虽然向北有斜率，在那个点上方向导数仍然为零。操作时可分为一阶导数、二阶导数和

12、方向曲率三种方式。一阶导数（First Derivative）一阶导数给出表面沿方向线的斜率. 由一阶导数的网格文件生成的等值线图绘制出沿着指定的方向具有相同斜率的等值线。在各网格结点处, 如向上倾斜则斜率为正，如向下倾斜则斜率为负。若斜率接近正或负的无限大，表示在向上或向下的向上接近垂直。二阶导数（Second Derivative）二阶导数给出沿着方向线斜率的变化率。由二阶导数的网格文件生成的等值线图绘制出穿过表面斜率变化率相同的等值线。其值随表面上升下降也有正值和负值；数值的大小反映表面缓陡变化的快慢。方向曲率（Directional Curvature）方向曲率是曲面沿着规定方向线，在

13、切面倾角的变化率。方向曲率采用变化率的绝对值，因此它是一个正数。方向曲率类似二阶导数。2） 地形模型(Terrain Modeling) 在表面上每一网格结点处,Surfer确定了该点最陡的斜坡和最陡斜率的方向(梯度方向)。地形模型的值是从90到-90 度.最陡斜率零度是水平线, 90 度是垂直上升的，而 -90 是垂直下降的。5种地形模型的操作如下；地形斜率、地形方位、纵剖面曲率、平面曲率和切线曲率。可以选择地形模型的临界值，缺省为一极小值。地形斜率（Terrain Slope）穿过表面，斜率梯度的方向可以改变。由地形斜率的网格文件生成的等值线图或线网图绘制出等最大斜率线。该功能与一阶方向导

14、数相同，用来确定表面任何一点的斜率，在考虑自动确定图中每一个点的的梯度方向时更为有用。梯度方向可由地形方位功能得到。地形方位（Terrain Aspect）给出每一网格结点最大斜率的方向。比如表面的水流方向严格垂直于等值线的方向。其值用方位角来表示，0度为北；90度为东等。可以利用地形方位功能创建一幅反映流动方向的粘贴图。纵剖面曲率（Profile Curvature）纵剖面曲率确定每一个结点延梯度方向斜率变化率的增减情况。由纵剖面曲率的网格文件生成的等值线图显示穿过表面，最大斜率的等变化率线。该功能可与二阶方向导数对比，但更加有用，因为它可以自动地判定表面每一点的下降方向，然后确定该点沿此方

15、向斜率的变化率。平面曲率（Plan Curvature）平面曲率反映了在水平面上地形方位角测量的变化率，描述了表面水流的会聚或发散。由平面曲率的网格文件生成的等值线图显示等会聚或等发散线。在网格文件中的负值表示发散，而正值表示会聚。切线曲率（Tangential Curvature）切线曲率与平面曲率类似，也用来反映会聚或发散，但它测量的是在垂直面上的情况。该功能描述表面水流的会聚或发散，可以想象为沿垂直于等值线方向的垂直面的表面曲率。3） 微分和积分算符(Differential & Integral Operator)包括了4种功能：梯度算符、拉普拉斯算符、双谐算符和积分体积梯度算符（Gr

16、adient Operator）梯度算符在表面的任何点上产生最大斜率的网格，等于地形斜率操作, 但是梯度算符为一数值，而不是角度。水平表面梯度算符为零, 而接近于垂直面时梯度算符接近无限大。拉普拉斯算符（Laplacian）拉普拉斯算符用于提供表面充电和放电的测量。地下水、热量和电荷是三个守恒物理量的例子，其局部流动速率与梯度成正比。拉普拉斯算符对这些物理量的局部体积进行量化，流入时或在充电区 0，而流出时或在放电区 0。双谐（Biharmonic）算符薄板和壳体的弯曲、多孔介质内的粘滞流体和在线性弹力范围内的应力函数是能利用双谐算符数学描述的三个物理量的例子。它等同于两次使用拉普拉斯算符。积

17、分体积（Integral Volume）积分体积操作给出整个网格范围内或子网格的堆积体积。子网格可以是任何外形。例如生成等值线图来表示煤矿内煤层的厚度，积分体积图可用于计算每一个月内开采煤的体积。另一个例子是雷暴雨的降雨密度，如果计算了积分体积并覆盖于分水岭边界图, 就可以直接计算流入每一溪流的水流量。如定义从Xmin到Xmax和从Ymin到Ymax的矩形区域，V(x,y) 是在Z=0的平面和表面覆盖下的区域之间的体积。因此, V(xmin,ymin) = 0,而V(xmax,ymax)等于整个矩形的区域下面的体积，其值由GridVolume命令得到。由积分体积功能也很容易计算矩形子网格的体积

18、。4）傅立叶及谱分析(Fourier and Spectral Analysis)包括2个功能：相关图（Correlogram）相关图用于评估空间模型和空间相关性。相关图显示在网格中潜在的趋势和测量网格各向异性。在右上和左下象限相反的意义上，相关图是对称的，虽然不完全一样，但是非常接近。相关图与父本网格有相同的尺度，但围绕原点有正的和负的单位。网格Peridiogram网格peridiogram是将表面分解成多个两维的正弦曲线的加权和。它显示隐藏的对称性以及等值线图或线网图中不明显的重复图案。Surfer以两维的傅立叶频率计算网格peridiogram，没有任何平滑、细化和拟合。6.2 网格化

19、方法6.2.1 网格化基础知识：网格化是把以XYZ数据文件格式表示的、通常是不规则分布的原始数据点，经过数学处理，构筑一个规则的空间矩形网格的过程。原始数据的不规则分布，造成缺失数据的“空洞”。网格化则用外推或内插的算法填充了这些“空洞”。大多数情况下，采用加权平均插值算法，即所有其它参数相等的条件下，愈靠近结点（计算出的规则点）的数据（原始数据点），对计算该结点的Z值贡献愈大。插值方法分为两种：精确插值(Exact interpolators)和平滑插值(Smoothing interpolators)。视网格化所用的数学模型和设定参数的不同，一种网格化方法可以属于两种插值方法中的一种或另一

20、种。精确插值指当网格结点正好位于原始数据点时，该结点的Z值等于此原始数据点的Z值。对于加权平均内插算法，这就意味着此原始数据点的权重为1，而其它数据点对于该结点的权重为0。增加网格密度，就增大了网格结点正好位于原始数据点的可能性。平滑插值用于并不十分依赖原始数据，只试图了解Z值的总体变化趋势的情况。平滑插值不会给任何数据点以权重1，即使某网格结点正好位于原始数据点。每一种网格化的方法都有自己的一组设置。对于每种方法来说，数据处理和方向性都是类似的。6.2.2距离反比法（Inverse Distance to a Power）一种加权平均插值的网格化方法。在计算一个网格结点的Z值时，一定范围内，

21、所有数据点的权重的和为1，权重与某数据点到该结点距离成反比，愈靠近该结点的原始数据点，其权重愈大。如果网格结点正好位于某原始数据点，该结点的Z值就等于此原始数据点的Z值，即此原始数据点对于该结点的权重为1，而其它数据点对于该结点的权重为0。可见，距离反比法是一种准确插值方法。权重系数（Power）的值反映随着数据点到网格结点距离的增加，其权重降低的程度。当权重系数为0时，所计算出的网格面是一个接近水平的平面，其Z值为所有原始数据点的平均值。随着权重系数的增加，形成最临近点插值，导致表面变成多边形。平滑参数（Smoothing）把不确定性与用户输入的数据联系起来，平滑参数愈大，相邻网格的影响愈小

22、。平滑参数大于0，则没有任何一个数据点对于某个结点的权重为1，即使该数据点正好位于网格结点上。各向异性（Anisotropy）指在用原始数据点计算网格点Z值时，对沿某一个坐标轴方向的数据点比其它方向的数据点给予更多的权重。多数情况下，在构筑网格时，不需要考虑方向性，因为大多数等值线图或线网图的X、Y坐标是同一比例尺。这时，X单位1Y单位。当确实需要考虑各向异性时，Ratio对话框内为1个X单位相对与一个Y单位的比例。如选Ratio为2，则意味X单位2Y单位.下面沿X轴方向伸长的椭圆示意不同方向给予不同的权重。选用不同的网格化方法，各向异性的参数略有不同。距离反比法的特点之一就是在网格区内围绕着

23、某些数据点可能产生牛眼状等值线（Bulls-eyes）。距离反比法是一种快速网格化的方法，在小于500数据点时，可以用全部数据点来生成网格。6.2.3 最小曲率法（Minimum Curvature）这是一种在地学中广泛应用的网格化方法。由最小曲率法构成的插值表面像一个线性弹性薄板，是一个尽可能与原始数据点吻合的最平滑的曲面。最小曲率法不是准确插值，是典型的平滑插值。最小曲率法的对话框包括以下内容：最大残差（Max Reciduals）：单位与数据的相同，比较合适的值是数据精度的10%。缺省的最大残差为 0.001 (Zmax - Z min)最大重复参数（Max lterations）：通常

24、设为网格结点数的1到2倍。例如，对于5050的网格，最大重复参数在2500与5000之间。内部和边缘张性系数（Internal and Boundary Tension）：设定弹性薄板内部和边缘弯曲度的参数。该值愈大，弯曲愈小。缺省值均为0。松弛系数（Relaxation Factor）：算法参数，通常，该值愈大，迭代算法会聚愈快。缺省值为1，一般不用另设定。6.2.4 多项式回归法（Polynomial Regression）多项式回归法用来确定用户数据整体的趋势或构造一种模型。多项式回归法实际上并不是一种插值，因为它并不试图预测未知的Z值。用户可以在表面选定（Surface Difinat

25、ion）框内选择以下4种曲面中的任一种。简单平面（Simple planer surface）z(x,y)=A+Bx+Cy双线性鞍形（Bi-linear saddle）z(x,y)=A+Bx+Cy+Dxy二次表面（Quadratic surface）z(x,y)=A+Bx+Cy+Dx2+Exy+Fy2三次表面（Cubic surface）z(x,y)=A+Bx+Cy+Dx2+Exy+Fy2+Gx3+H x2y+Ix y2+Jy3或由用用户自定义。选定的曲面方程相应显示在下面；右面的参数框内则显示X、Y和总的最高项次。用户也可以利用Parameters框自定义多项式方程。6.2.5三角形线性插值

26、法（Triangulation with Linear Interpolation）三角形线性插值法是一种准确插值，方法是在相邻点之间连线构成三角形，并且保持任一三角形的边都不与其它三角形的边相交。这样在网格范围内由一系列三角形平面构成拼接图。由于数据点平均分布，在通过地形变化显示断层线时，三角形法非常有效。因为每一个三角形都构成一个平面，所有的结点都在三角形中，其坐标被三角形平面方程唯一地确定。对于有200至1000个数据点，且平均地分配在网格区域里时，用三角形线性插值法最好。6.2.6修订了的Shepards法Shepards法是一种距离反比加权的最小二乘法。与距离反比法插值法相似，但由于

27、使用局部最小二乘法，消除或减少了绘制等值线时的“牛眼”效应。Shepards法可以是准确插值或者是平滑插值。可以设置网格化的平滑参数，允许进行平滑插值。随平滑参数值的增加，平滑效果愈明显。通常，该值在0和1 之间最合适。Quadratic Neighbors 指定进行最小二乘法的范围，即计算半径内的原始数据点数。Weighting Neighbors指定进行加权平均的范围，即计算半径内的原始数据点数。缺省值为Renka (1988)推荐。6.2.7 Kriging法Kriging法是用南美采矿工程师D. G. Krige的名字命名的一种地学统计内插方法，原来是试图比较准确地预测矿石储量。已经发

28、现在许多领域非常有用。Kriging法描述数据中隐含的趋势。比如孤立的高值点“牛眼”，在使用Kriging法网格化时可连成“山脊”。变异图模型（Variogram Model）：用来确定插值每一个结点时所用数据点的邻域，以及在计算结点时给予数据点的权重。Surfer提供了多种最常用的变异图模型，它们是指数、高斯模型、线性、对数、矿块效应、幂、二次模型、有理数二次模型、球面模型和波（空洞效应）。如果拿不准用哪一种变异图，可选用线性变异图，大多数情况下，效果较好。每一种模型都有Slope，Scale，Length等参数要求设定。比例系数（变异图方程中的C）用来确定所选择的变异图模型的sill，除了

29、线性变异图以外（没有sill），Sill等于矿块效应加变异图比例。当你没有设定任何矿块效应值时，sill等于比例值。偏移类型（Drift Type）：当对原始数据点分布的“空洞”和边界之外的点进行插值计算时，偏移类型功能将有明显影响。Surfer提供了3种偏移类型：无偏移、线性偏移和二次偏移。拿不准时最好选无偏移，即采用普通Kriging法插值。线性偏移和二次偏移被用于实施普通Ktiging 插值。如果数据的变化趋势围绕着一种线性趋势，则采用线性偏移；如果数据的变化趋势围绕着二次趋势（即抛物线型），则采用二次偏移。矿块效应（Nugget Effect）用于在收集数据时存在潜在错误的情况下。指定

30、矿块效应会导致Kriging 产生更为光滑的插值，即个别数据点吻合较差但反映了全体数据的整体趋势。矿块效应愈高，产生的网格愈光滑。矿块效应有两部分构成：矿块效应误差方差微方差误差方差编辑框允许用户设定测量误差的方差；微方差编辑框允许用户设定小规模结构的方差。当误差方差为时，非的矿块效应具有一般光滑的效果，但产生的网格仍然与每个观察点吻合（可视为一种准确插值）。一个非的误差方差允许网格不同于观察点（是一种平滑插值）。6.2.8 径向基本函数法（Radial Basis Functions）径向基本函数法是一种准确插值的方法。其中的多重二次曲面法被许多人认为是最好的方法。在插值生成一个网格结点时，

31、这些函数确定了使用数据点的最优权重组。径向基本函数法的函数类型包括：反比多重二次曲面法；多重对数；多重二次曲面法；自然三次样条和薄板样条。径向基本函数法类似Kriging法中的变异图。在大多数情况下，多重二次曲面函数是最合乎要求的。R2参数是一个决定锐化或平滑的参数。R2值愈大，山顶愈圆滑，等值线愈平滑。R2合理的实验值是在一个平均样本间距和半个平均样本间距之间。6.2.9 最近临点法（Nearest Neighbor）最近临点法用最临近的数据点来计算每个网格结点的值。这种方法通常用于已有规则网格只需要转换为Sufer网格文件时，或数据点几乎构成网格，只有个别点缺失，该方法可以有效地填充“空洞”。通过设置搜寻椭圆半径的值小于数据点之间距离的方法，给缺少数据点的结点赋值为空白。6.2.10 普通临点法Natural Neighbor Surfer7新增加的算法。一种相临点加权平均的内插算法，权重与borrowed 区成正比。不能外推。6.2.11 推荐的构造网格的数学模型不同的网格化方法可以得到不同的网格文件，用户应当选用最能代表自己数据特点的方法，选择网格化方法时应当考虑原始数据点数量的多寡。10个或10个以下的数据点，除了反映数据的一般趋势外，没有多大意义。这样少的点，三角形法无效，数据点25