如何封装成一个自动生成山顶点的程序？

weyawu

网上找的，本人不会AE，如何封装成一个自动生成山顶点的程序？
ArcGIS Engine实现提取山顶点         (2014-04-10 22:53:02)[url=]转载▼[/url]

分类： 编程技术

前两天帮别人的毕业论文写了一个小程序，包含一个提取山顶点的功能。山顶点在ArcGIS里不能直接提取，需要经过几个步骤，所以AE里就没有直接的方法。不过不需要额外的算法，基本上都是SpatialAnalyst下面的接口，算是比较简单的一个应用。
AE提取山顶点主要包括以下几个步骤：
1.邻域统计
主要使用INeighborhoodOp中的焦点统计（FocalStatistice），用IRasterNeighborhood设置邻域为矩形（Rectangle）。因为窗口的长和宽能够直接影响最后提取山顶点的结果，所以这两个参数由用户自己指定。设置统计类型为最大值（Maximum），这样操作后得到的结果就使得每一个像元的值都是周围矩形范围内所有像元的最大值。
// 输入：geoDataset、width、height
// 输出：geoDataset_result
INeighborhoodOp pNeighborhoodOP = new RasterNeighborhoodOpClass();
IRasterNeighborhood rasterNeighborhood = new RasterNeighborhoodClass();
// 设置邻域窗口
// width和height为邻域窗口的宽和高
rasterNeighborhood.SetRectangle(width, height,ESRI.ArcGIS.GeoAnalyst.esriGeoAnalysisUnitsEnum.esriUnitsCells);
// 邻域统计
// geoDataset为邻域窗口的宽和高
IGeoDataset geoDataset_result = pNeighborhoodOP.FocalStatistics(geoDataset,esriGeoAnalysisStatisticsEnum.esriGeoAnalysisStatsMaximum, rasterNeighborhood,true);


2.栅格计算
用邻域统计得到的最大值栅格，减去原本的DEM栅格数据，所得结果中所有像元值为0的点就是山顶点。这里我没有使用IMathOp或IConditionalOp，而是用了条件表达式IMapAlgebraOp。首先用IMapAlgebraOp中的BindRaster方法设置栅格数据的别名，比如“1”和“2”，然后用Execute方法求取表达式“"[2] - [1] == 0"”。注意操作符号前后要加空格，别名要用中括号括起来。这样就得到了一个二值图像，符合表达式的值为1，即为山顶点。
// 输入：geoDataset、geoDataset2
// 输出：geoDataset_result
IMapAlgebraOp pAlgebra = new RasterMapAlgebraOpClass();
// 设置别名
// geoDataset为原始栅格
// geoDataset2为上一步得到的栅格
pAlgebra.BindRaster(geoDataset, "1");
pAlgebra.BindRaster(geoDataset2, "2");
// 条件表达式
IGeoDataset geoDataset_result = pAlgebra.Execute("[2] - [1] == 0")



3.重分类
用IReclassOp，将值为1的像元还设为1，将值为0的像元设为NoData。
// 输入：geoDataset
// 输出：geoDataset_result
IReclassOp pReclassOp = new RasterReclassOpClass();
INumberRemap pNumRemap = new NumberRemapClass();
// 0设置为NoData

pNumRemap.MapValueToNoData(0);
// 1设置为1
NumRemap.MapValue(1, 1);

IRemap pRemap = pNumRemap as IRemap;


// 重分类
// geoDataset为上一步得到的栅格
IGeoDataset geoDataset_result = pReclassOp.ReclassByRemap(geoDataset, pRemap,true)



4.栅格转点
用IConversionOp中的RasterDataToPointFeatureData方法，需要制定数据保存路径，以创建工作空间（IWorkspace）。
// 输入：geoDataset，filePath，fileName
// 输出：geoDataset_result
IWorkspaceFactory pWorkspaceFactoryShp = new RasterNeighborhoodOpClass();
// filePath是数据存储的文件夹
IWorkspace pWorkspace = pWorkspaceFactoryShp.OpenFromFile(filePath, 0);


IConversionOp pConversionOp = new RasterNeighborhoodOpClass();
// 栅格转点
// geoDataset为上一步得到的栅格
// fileName是文件名
IGeoDataset geoDataset_result = pConversionOp.RasterDataToPointFeatureData(geoDataset, pWorkspace, fileName);



最后把生成的要素加载到地图上：
IDataset pDataset = geoDataset as IDataset;
IFeatureClass pFeatureClass = pDataset as IFeatureClass;
IFeatureLayer pFeatureLayer = new FeatureLayerClass();
pFeatureLayer.FeatureClass = pFeatureClass;
pFeatureLayer.Name = fileName;
AxMapControl.Map.AddLayer(pFeatureLayer as ILayer);

huangzhi1116

模型编辑器里面建一个工作流

weyawu

谢谢huangzhi1116版主，研究了一下，还是没整明白。