ArcEngine+C# TIN相关三维功能模块介绍（四）

ArcEngine+C#  TIN相关三维功能模块介绍（四）

功能综合演示

 

作者：刘志远

PDF浏览  示例程序下载

                                                  

1.简介

       本篇打算对TIN的介绍做个小结，把前面将的功能都集中到一个小程序中。除了前面将的功能外，本篇还添加了如下一些功能。只要有数据，这些功能都是可以单独使用的，方便读者根据具体的需要直接参考相应代码。

Ø         由矢量点或线数据，通过核函数密度制图方法生成栅格数据（类DEM）；

Ø         按给定间距，从栅格数据（DEM）中提取矢量等高线数据；

Ø         由栅格数据（DEM）生成TIN模型；

Ø         动态变化DEM、TIN的渲染样式；

Ø         矢量图层叠加到TIN模型；

       本文主要对这几个功能进行简单的讲解，希望对接触TIN的朋友有所帮助。具体的实现，可以参看程序中相应的代码。同时，在压缩包中有个简短的“使用说明”，可以通过该程序尝试讲过的功能。

 

2.核函数密度制度生成“类DEM”

       通常情况下，我们通过对高程点数据进行插值操作（如反距离权插值、样条插值、克里金插值等）生成DEM（数字高程模型）栅格数据，生成的栅格数据的像元灰度值即代表对应的高程值。

       本文介绍的是密度制图，生成的栅格灰度值代表密度大小，所以这里成为“类DEM”。密度制图有核函数密度制图（Kernal）和简单密度制图（Simple）两种方法，本文用的是Kernal方法。对于二次开发来说，由矢量数据生成栅格数据，不论是DEM也好，还是密度制图，只是调用相应接口的某个方法而已，内部的具体操作都是封装好的，程序员只需根据实际需要进行选择即可。这也体现了面向对象的巨大优势。主要代码如下：

 

 1    /**//// <summary>
 2        /// 核函数密度制图
 3         /// </summary>
 4        /// <param name="pFeatureLayer">进行密度分析的矢量数据图层</param>
 5        /// <param name="cellSize">输出栅格格网大小</param>
 6        /// <param name="radius">制图搜索半径</param>
 7        /// <param name="path">栅格文件保存路径  可以为空，则作为临时数据处理</param>
 8        /// <returns></returns>
 9        public IRaster DensityAnalyst(IFeatureLayer pFeatureLayer, double cellSize, double radius, string path)
10        {
11            IFeatureClass pFClass01 = pFeatureLayer.FeatureClass;
12            IDensityOp pDensityOp = new RasterDensityOpClass();
13            IRasterAnalysisEnvironment pEnv = pDensityOp as IRasterAnalysisEnvironment;
14
15            Double double_cellSize = cellSize;
16            object object_cellSize = (System.Object)double_cellSize;
17            pEnv.SetCellSize(esriRasterEnvSettingEnum.esriRasterEnvValue, ref object_cellSize);
18            IFeatureClassDescriptor pFDescr = new FeatureClassDescriptorClass();
19            pFDescr.Create(pFClass01, null, "NONE");
20
21            System.Double double_radio_dis = radius;
22            object object_radio_dis = (System.Object)double_radio_dis;
23            object Missing = Type.Missing;
24            //调用核函数密度制图方法生成栅格数据
25            IRaster pRasOut1 = pDensityOp.KernelDensity(pFDescr as IGeoDataset, ref object_radio_dis, ref Missing) as IRaster;
26            
27            //栅格值重计算，将栅格值缩放到合理范围（主要是纠正由于坐标系统不同而造成的数据夸张变大）
28            //如数据的坐标系等没有问题则可以不用重计算
29            IMapAlgebraOp pMapAlgebraOp = (IMapAlgebraOp)new RasterMapAlgebraOp();
30            IRasterAnalysisEnvironment pRasAnaEnv = (IRasterAnalysisEnvironment)pMapAlgebraOp;
31            pMapAlgebraOp.BindRaster(pRasOut1 as IGeoDataset, "R");
32            string expresion = "[R] * 0.0000915"; //坐标系统转换参数
33            IRaster pRasOut = pMapAlgebraOp.Execute(expresion) as IRaster;
34
35            //输出栅格数据，否则存放为临时数据
36            if (path != "")
37            {
38                string foldPath = System.IO.Path.GetDirectoryName(path);
39                string rasterName = System.IO.Path.GetFileName(path);
40                try
41                {
42                    IRasterBandCollection bandCollection = (IRasterBandCollection)pRasOut;
43                    IRasterDataset dataset = bandCollection.Item(0).RasterDataset;
44                    ITemporaryDataset temp = (ITemporaryDataset)dataset;
45                    RasterWorkspaceFactory factory = new RasterWorkspaceFactory();
46                    IWorkspace workSpace = factory.OpenFromFile(foldPath, 0);
47                    temp.MakePermanentAs(rasterName, workSpace, "IMAGINE Image");
48                }
49                catch
50                {
51                    return pRasOut;
52                }
53            }
54            return pRasOut;
55        }
56

 

 

3.DEM提取等高线

       根据指定间隔，从DEM栅格数据中提取等高线数据。主要用了ISurfaceOp接口下的Contour（）方法。主要代码如下：

 

    ISurfaceOp pSurfaceOP = new RasterSurfaceOpClass();
2    IFeatureLayer pFL = new FeatureLayerClass(); 
3    object Missing = Type.Missing;
4    pFL.FeatureClass = pSurfaceOP.Contour(iRaster as IGeoDataset, interval, ref Missing) as IFeatureClass;
5

 

 

4.由DEM生成TIN

 

 

  1               //***************生成TIN模型*********************************************
  2                IGeoDataset pGeoData = iRaster as IGeoDataset;
  3                IEnvelope pExtent = pGeoData.Extent;
  4                IRasterBandCollection pRasBC = iRaster as IRasterBandCollection;
  5                IRasterBand pRasBand = pRasBC.Item(0);
  6                IRawPixels pRawPixels = pRasBand as IRawPixels;
  7                IRasterProps pProps = pRawPixels as IRasterProps;
  8
  9                int iWid = pProps.Width;
 10                int iHei = pProps.Height;
 11
 12                double w = iWid / 1000.0f;
 13                double h = iHei / 1000.0f;
 14
 15                IPnt pBlockSize = new DblPntClass();
 16                bool IterationFlag;
 17
 18                if (w < 1 && h < 1) //横纵都小于1000个像素
 19                {
 20                    pBlockSize.X = iWid;
 21                    pBlockSize.Y = iHei;
 22                    IterationFlag = false;
 23                }
 24                else
 25                {
 26                    pBlockSize.X = 1001.0f;
 27                    pBlockSize.Y = 1001.0f;
 28                    IterationFlag = true;
 29                }
 30
 31                double cellsize = 0.0f;     //栅格大小
 32                IPnt pPnt1 = pProps.MeanCellSize(); //栅格平均大小
 33                cellsize = pPnt1.X;
 34
 35                ITinEdit pTinEdit = new TinClass() as ITinEdit;
 36                pTinEdit.InitNew(pExtent);
 37
 38                ISpatialReference pSpatial = pGeoData.SpatialReference;
 39                pExtent.SpatialReference = pSpatial;
 40
 41                IPnt pOrigin = new DblPntClass();
 42                IPnt pPixelBlockOrigin = new DblPntClass();
 43
 44                //栅格左上角像素中心坐标
 45                double bX = pBlockSize.X;
 46                double bY = pBlockSize.Y;
 47
 48                pBlockSize.SetCoords(bX, bY);
 49                IPixelBlock pPixelBlock = pRawPixels.CreatePixelBlock(pBlockSize);
 50
 51                object nodata = pProps.NoDataValue;     //无值标记
 52                ITinAdvanced2 pTinNodeCount = pTinEdit as ITinAdvanced2;
 53                int nodeCount = pTinNodeCount.NodeCount;
 54
 55                object vtMissing = Type.Missing;
 56
 57                object vPixels = null;      //格子
 58                if (IterationFlag)   //当为一个处理单元格子时
 59                {
 60                    pPixelBlockOrigin.SetCoords(0.0f, 0.0f);
 61                    pRawPixels.Read(pPixelBlockOrigin, pPixelBlock);
 62
 63                    vPixels = pPixelBlock.get_SafeArray(0);
 64                    double xMin = pExtent.XMin;
 65                    double yMax = pExtent.YMax;
 66                    pOrigin.X = xMin + cellsize / 2;
 67                    pOrigin.Y = yMax - cellsize / 2;
 68                    bX = pOrigin.X;
 69                    bY = pOrigin.Y;
 70
 71                    pTinEdit.AddFromPixelBlock(bX, bY, cellsize, cellsize, nodata, vPixels, m_zTolerance, ref vtMissing, out vtMissing);
 72                }
 73                else  //当有多个处理单元格时，依次循环处理每个单元格
 74                {
 75                    int i = 0, j = 0, count = 0;
 76                    int FirstGoNodeCount = 0;
 77                    while (nodeCount != FirstGoNodeCount)
 78                    {
 79                        count++;
 80                        nodeCount = pTinNodeCount.NodeCount;
 81                        //依次循环处理
 82                        for (i = 0; i < h + 1; i++)
 83                        {
 84                            for (j = 0; j < w + 1; j++)
 85                            {
 86                                double bX1, bY1, xMin1, yMax1;
 87                                bX1 = pBlockSize.X;
 88                                bY1 = pBlockSize.Y;
 89
 90                                pPixelBlockOrigin.SetCoords(j * bX1, i * bY1);
 91                                pRawPixels.Read(pPixelBlockOrigin, pPixelBlock);
 92                                vPixels = pPixelBlock.get_SafeArray(0);
 93
 94                                xMin1 = pExtent.XMin;
 95                                yMax1 = pExtent.YMax;
 96
 97                                bX1 = pBlockSize.X;
 98                                bY1 = pBlockSize.Y;
 99
100                                pOrigin.X = xMin1 + j * bX1 * cellsize + cellsize / 2.0f;
101                                pOrigin.Y = yMax1 + i * bY1 * cellsize - cellsize / 2.0f;
102
103                                bX1 = pOrigin.X;
104                                bY1 = pOrigin.Y;
105
106                                pTinEdit.AddFromPixelBlock(bX1, bY1, cellsize, cellsize, nodata, vPixels, m_zTolerance, ref vtMissing, out vtMissing);
107
108                                FirstGoNodeCount = pTinNodeCount.NodeCount;
109                            }
110                        }
111                    }
112                }
113
114                //保存TIN文件
115                pTinEdit.SaveAs(tinFileName, ref vtMissing);
116              pTinEdit.StopEditing(true);
117

 

 

5.动态变化DEM、TIN的渲染样式

       动态变化DEM、TIN的渲染样式其实是根据用户的设定重新生成对应的渲染样式而已，这里就不贴代码了，具体的可以参考程序。

 

6.矢量图层叠加到TIN模型

       矢量图层叠加到TIN模型，其实是改变图层的属性，将图层中对应的ILayerExtensions接口变为I3DProperties接口，然后将I3DProperties. BaseSurface属性设为TIN图层表面数据。具体代码如下：

 

 1        /**//// <summary>
 2        /// 将矢量图层叠加到TIN模型上
 3        /// </summary>
 4        /// <param name="pLayer">要叠加的图层</param>
 5        /// <param name="pSceneControl">三维控件</param>
 6        public void setLayerToTIN(ILayer pLayer, AxSceneControl pSceneControl)
 7        {
 8            I3DProperties p3DProperties;
 9            p3DProperties = get3DProps(pLayer);
10            p3DProperties.BaseOption = esriBaseOption.esriBaseSurface;
11            ISurface pSurface = getTinlayer(pSceneControl).Dataset as ISurface;
12            p3DProperties.BaseSurface = pSurface;
13            p3DProperties.Apply3DProperties(pLayer);
14        }
15
16        /**//// <summary>
17        /// 恢复到平常样式
18        /// </summary>
19        /// <param name="pLayer"></param>
20        /// <param name="pSceneControl"></param>
21        public void setLayerOutTIN(ILayer pLayer, AxSceneControl pSceneControl)
22        {
23            I3DProperties p3DProperties;
24            p3DProperties = get3DProps(pLayer);
25            p3DProperties.BaseOption = esriBaseOption.esriBaseExpression;
26            p3DProperties.Apply3DProperties(pLayer);
27        }
28
29        /**//// <summary>
30        /// 获得普通要素图层的3D属性
31        /// </summary>
32        /// <param name="pLayer">普通矢量图层</param>
33        /// <returns></returns>
34        public I3DProperties get3DProps(ILayer pLayer)
35        {
36            ILayerExtensions pLayerExts = pLayer as ILayerExtensions;
37            if (pLayerExts != null)
38            {
39                for (int i = 0; i < pLayerExts.ExtensionCount; i++)
40                {
41                    I3DProperties p3DProps = pLayerExts.get_Extension(i) as I3DProperties;
42                    if (p3DProps != null)
43                        return p3DProps;
44                }
45            }
46            return null;
47        }
48
49       /**//// <summary>
50        /// 返回TIN数据图层
51        /// </summary>
52        /// <param name="pSceneControl">三维控件</param>
53        /// <returns></returns>
54        public ITinLayer getTinlayer(AxSceneControl pSceneControl)
55        {
56            ITinLayer pTinlayer = null;
57            IScene map = pSceneControl.Scene;
58            if (map == null)
59                return null;
60            for (int i = 0; i < map.LayerCount; i++)
61            {
62                ILayer lyr = map.get_Layer(i);
63                if (lyr is ITinLayer)
64                {
65                    pTinlayer = lyr as ITinLayer;
66                    break;
67                }
68            }
69            return pTinlayer;
70        }
71

 

 

7.总结

       至此，关于TIN的一些介绍就到这告一段落，文中没有阐述清楚的部分可以参看程序中的代码。当然，其中不足之处还请各位留言批评指正。