ArcEngine+C# TIN相关三维功能模块介绍(四)
功能综合演示
作者:刘志远
本篇打算对TIN的介绍做个小结,把前面将的功能都集中到一个小程序中。除了前面将的功能外,本篇还添加了如下一些功能。只要有数据,这些功能都是可以单独使用的,方便读者根据具体的需要直接参考相应代码。
Ø 由矢量点或线数据,通过核函数密度制图方法生成栅格数据(类DEM);
Ø 按给定间距,从栅格数据(DEM)中提取矢量等高线数据;
Ø 由栅格数据(DEM)生成TIN模型;
Ø 动态变化DEM、TIN的渲染样式;
Ø 矢量图层叠加到TIN模型;
本文主要对这几个功能进行简单的讲解,希望对接触TIN的朋友有所帮助。具体的实现,可以参看程序中相应的代码。同时,在压缩包中有个简短的“使用说明”,可以通过该程序尝试讲过的功能。
2.核函数密度制度生成“类DEM”
通常情况下,我们通过对高程点数据进行插值操作(如反距离权插值、样条插值、克里金插值等)生成DEM(数字高程模型)栅格数据,生成的栅格数据的像元灰度值即代表对应的高程值。
本文介绍的是密度制图,生成的栅格灰度值代表密度大小,所以这里成为“类DEM”。密度制图有核函数密度制图(Kernal)和简单密度制图(Simple)两种方法,本文用的是Kernal方法。对于二次开发来说,由矢量数据生成栅格数据,不论是DEM也好,还是密度制图,只是调用相应接口的某个方法而已,内部的具体操作都是封装好的,程序员只需根据实际需要进行选择即可。这也体现了面向对象的巨大优势。主要代码如下:
代码
/// <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
try41
{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
catch50
{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
代码
//***************生成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
else25
{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图层表面数据。具体代码如下:
代码
/// <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

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


