Arcgis Engine矢量裁剪栅格，调用Mask工具相关代码

今天试了下arcgis 矢量裁剪栅格，首先在arcmap上面操作，发现有好几种方式都能实现：

方法1，矢量转栅格，把像元值设置1(方法是添加一个字段，将字段值全部赋值1）然后调用toolbox的时候，选择该字段；这样就生成了一个所有cell的值为1的栅格RasterA。

     接下来，将RasterA与我们要裁剪的栅格对象RasterB进行栅格计算(栅格计算也用toolbox里面的，不要用菜单里面的，因为我们最后代码里面调用的是toolbox)，将RasterA*RasterB就是最终结果啦。

方法2，和方法1的区别就是，不用添加字段，而是重分类一下，Reclassify重分类，将某个字段的值重分类为1。

方法3，也是最屌的方法，直接用Extract by Polygon工具。

3种方法区别：

方法1和方法2算是同一种，他们比方法3稳定，因为他们只是单纯栅格计算，在大数据量的时候比较靠谱。方法3虽然简单，但面对大数据量处理的时候，容易出错，不太稳定。当然小数据的话，还是木有问题的。

下边直接贴代码啦，简单的封装了一个类

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using ESRI.ArcGIS.Geodatabase;
using ESRI.ArcGIS.Geometry;
using ESRI.ArcGIS.Geoprocessor;
using ESRI.ArcGIS.Geoprocessing;
using ESRI.ArcGIS.esriSystem;
using ESRI.ArcGIS.SpatialAnalystTools;
using ESRI.ArcGIS.SpatialAnalyst;
using ESRI.ArcGIS.DataSourcesRaster;
using ESRI.ArcGIS.Carto;
using ESRI.ArcGIS.DataSourcesGDB;

namespace GPOperator
{
    public class ClipHelper : GPTool
    {
        /// <summary>
        /// 栅格裁剪
        /// </summary>
        /// <param name="inputRaster">要裁剪的栅格(支持对象和栅格文件路径)</param>
        /// <param name="maskData">掩膜数据(支持对象和shp文件路径)</param>
        /// <returns></returns>
        public IGeoDataset ClipRaster2(object inputRaster, object maskData, out string resultInfo)
        {
            string outRaster = System.IO.Path.Combine(System.IO.Path.GetTempPath(), "ClipRaster_" +DateTime.Now.ToString("yyyy-MM-dd HH-mm-ss")+ Guid.NewGuid().ToString

("D").Substring(0,4) + ".tif");
            //参数：对象或路径都支持
            ESRI.ArcGIS.SpatialAnalystTools.ExtractByMask mask = new ESRI.ArcGIS.SpatialAnalystTools.ExtractByMask(inputRaster, maskData, outRaster);//@"C:\tmp\新建文件夹

\ret.tif"
            bool isok = this.Excute(mask, out resultInfo);
            //IWorkspaceFactory wsf = new RasterWorkspaceFactory();
            //IRasterWorkspace rasterWS = wsf.OpenFromFile(System.IO.Path.GetDirectoryName(outRaster), 0) as IRasterWorkspace;
            //IRasterDataset rds = rasterWS.OpenRasterDataset(System.IO.Path.GetFileName(outRaster));
            IRasterDataset rasterDataset = new RasterDataset();
            rasterDataset.OpenFromFile(outRaster);
            return rasterDataset as IGeoDataset;
        }
        /// <summary>
        /// 栅格裁剪
        /// </summary>
        /// <param name="inputRaster">要裁剪的栅格(支持对象和栅格文件路径)</param>
        /// <param name="maskData">裁剪多边形</param>
        /// <returns>裁剪得到的栅格</returns>
        public IGeoDataset ClipRaster(object inputRaster, IPolygon polygon, out string resultInfo)
        {
            var maskData = CreateFeatureClass(polygon);
            string outRaster = System.IO.Path.Combine(System.IO.Path.GetTempPath(), "ClipRaster_" + DateTime.Now.ToString("yyyy-MM-dd HH-mm-ss") + Guid.NewGuid().ToString

("D").Substring(0, 4) + ".tif");
            ESRI.ArcGIS.SpatialAnalystTools.ExtractByMask mask = new ESRI.ArcGIS.SpatialAnalystTools.ExtractByMask(inputRaster, maskData, outRaster);

            if (this.Excute(mask, out resultInfo))
            {
                return null;
            }
            IRasterDataset rasterDataset = new RasterDataset();
            try
            {
                rasterDataset.OpenFromFile(outRaster);
            }
            catch (Exception ex)
            {
                resultInfo = ex.ToString();
                return null;
            }
            finally
            {
                try
                {
                    ESRI.ArcGIS.ADF.ComReleaser.ReleaseCOMObject(maskData);
                }
                catch { }
            }
            return rasterDataset as IGeoDataset;
        }

        /// <summary>
        /// 栅格裁剪
        /// Op方法，直接返回对象
        /// </summary>
        /// <param name="inputRaster">要裁剪的栅格</param>
        /// <param name="maskData">掩膜数据</param>
        /// <returns></returns>
        public IGeoDataset ClipRaster3(IGeoDataset inputRaster, IGeoDataset maskData,out string resultInfo)
        {
            resultInfo = "执行成功";            
            IExtractionOp2 op = new RasterExtractionOpClass();
            try
            {
                return op.Raster(inputRaster, maskData);
            }
            catch (Exception ex)
            {
                resultInfo = ex.ToString();
                return null;
            }
        }
        public IFeatureClass CreateFeatureClass(IPolygon polygon)
        {
            IFeatureClassDescription fcDescription = new FeatureClassDescriptionClass();
            IObjectClassDescription ocDescription = (IObjectClassDescription)fcDescription;
            IFeatureClass fc= createFeatureClassInmemeory("temp", "temp", polygon.SpatialReference, esriGeometryType.esriGeometryPolygon, ocDescription.RequiredFields);
            if (null == fc)
            {
                return null;
            }
            IFeature f= fc.CreateFeature();
            f.Shape = polygon;
            f.Store();

            return fc;
        }

        //IFeatureCache
        /// <summary>
        /// 在内存中创建临时要素类
        /// </summary>
        /// <param name="DataSetName">数据集名称</param>
        /// <param name="AliaseName">别名</param>
        /// <param name="SpatialRef">空间参考</param>
        /// <param name="GeometryType">几何类型</param>
        /// <param name="PropertyFields">属性字段集合</param>
        /// <returns>IfeatureLayer</returns>
        private IFeatureClass createFeatureClassInmemeory(string DataSetName, string AliaseName, ISpatialReference SpatialRef, esriGeometryType GeometryType, IFields 

PropertyFields)
        {
            IWorkspaceFactory workspaceFactory = new InMemoryWorkspaceFactoryClass();
            ESRI.ArcGIS.Geodatabase.IWorkspaceName workspaceName = workspaceFactory.Create("", "MyWorkspace", null, 0);
            ESRI.ArcGIS.esriSystem.IName name = (IName)workspaceName;
            ESRI.ArcGIS.Geodatabase.IWorkspace inmemWor = (IWorkspace)name.Open();
            IField oField = new FieldClass();
            IFields oFields = new FieldsClass();
            IFieldsEdit oFieldsEdit = null;
            IFieldEdit oFieldEdit = null;
            IFeatureClass oFeatureClass = null;
            try
            {
                oFieldsEdit = oFields as IFieldsEdit;
                oFieldEdit = oField as IFieldEdit;
                for (int i = 0; i < PropertyFields.FieldCount; i++)
                {
                    IField field = PropertyFields.get_Field(i);
                    if (field.Type != esriFieldType.esriFieldTypeGeometry)
                    {
                        oFieldsEdit.AddField(field);
                    }
                }
                IGeometryDef geometryDef = new GeometryDefClass();
                IGeometryDefEdit geometryDefEdit = (IGeometryDefEdit)geometryDef;
                geometryDefEdit.AvgNumPoints_2 = 5;
                geometryDefEdit.GeometryType_2 = GeometryType;
                geometryDefEdit.GridCount_2 = 1;
                geometryDefEdit.HasM_2 = false;
                geometryDefEdit.HasZ_2 = false;
                geometryDefEdit.SpatialReference_2 = SpatialRef;
                oFieldEdit.Name_2 = "SHAPE";
                oFieldEdit.Type_2 = esriFieldType.esriFieldTypeGeometry;
                oFieldEdit.GeometryDef_2 = geometryDef;
                oFieldEdit.IsNullable_2 = true;
                oFieldEdit.Required_2 = true;
                oFieldsEdit.AddField(oField);
                oFeatureClass = (inmemWor as IFeatureWorkspace).CreateFeatureClass(DataSetName, oFields, null, null, esriFeatureType.esriFTSimple, "SHAPE", "");
                (oFeatureClass as IDataset).BrowseName = DataSetName;

                return oFeatureClass;
            }
            catch(Exception ex)
            {
                return null;
            }
        }
    }

    /// <summary>
    /// GP工具基类，实现执行工具公共方法和释放资源方法的定义。
    /// </summary>
    public class GPTool
    {
        private Geoprocessor gp = null;
        private IGeoProcessorResult gpResult = null;

        public GPTool()
        {
            gp = new Geoprocessor();
            gp.OverwriteOutput = true;
        }

        /// <summary>
        /// 释放资源
        /// </summary>
        public void Release()
        {
            if (gp != null)
            {
                ESRI.ArcGIS.ADF.ComReleaser.ReleaseCOMObject(gp);

                gp = null;
            }
            if (gpResult != null)
            {
                ESRI.ArcGIS.ADF.ComReleaser.ReleaseCOMObject(gpResult);
                gpResult = null;
            }
        }

        /// <summary>
        /// 执行已定义参数的GP工具。
        /// </summary>
        protected bool Excute(IGPProcess gpProcess, out string resultInfo)
        {
            resultInfo = "执行成功";
            gpResult = (IGeoProcessorResult)gp.Execute(gpProcess, null);

            if (gpResult == null || gpResult.Status != esriJobStatus.esriJobSucceeded)
            {
                // 执行错误，输入错误报告。
                StringBuilder message = new StringBuilder(gpProcess.ToolName);

                for (int i = 0; i < gp.MessageCount; i++)
                {
                    message.AppendLine(gp.GetMessage(i));
                }
                resultInfo = message.ToString();
                return false;
            }

            ESRI.ArcGIS.ADF.ComReleaser.ReleaseCOMObject(gpProcess);
            return true;
        }
    }
}

下边是调用代码：

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using ESRI.ArcGIS.esriSystem;
using GPOperator;
using ESRI.ArcGIS.DataSourcesRaster;
using ESRI.ArcGIS.Geodatabase;
using ESRI.ArcGIS.DataSourcesFile;
using ESRI.ArcGIS.Geometry;

namespace ConsoleApplication1
{
    class Program
    {
        private static LicenseInitializer m_AOLicenseInitializer = new ConsoleApplication1.LicenseInitializer();
    
        static void Main(string[] args)
        {
            //ESRI License Initializer generated code.
            m_AOLicenseInitializer.InitializeApplication(new esriLicenseProductCode[] { esriLicenseProductCode.esriLicenseProductCodeBasic, 

esriLicenseProductCode.esriLicenseProductCodeAdvanced },
            new esriLicenseExtensionCode[] { esriLicenseExtensionCode.esriLicenseExtensionCodeNetwork, esriLicenseExtensionCode.esriLicenseExtensionCodeSpatialAnalyst, 

esriLicenseExtensionCode.esriLicenseExtensionCodeSchematics, esriLicenseExtensionCode.esriLicenseExtensionCodeMLE, 

esriLicenseExtensionCode.esriLicenseExtensionCodeDataInteroperability, esriLicenseExtensionCode.esriLicenseExtensionCodeTracking, 

esriLicenseExtensionCode.esriLicenseExtensionCodeArcScan, esriLicenseExtensionCode.esriLicenseExtensionCodeBusiness, esriLicenseExtensionCode.esriLicenseExtensionCodeCOGO, 

esriLicenseExtensionCode.esriLicenseExtensionCodeGeoStats, esriLicenseExtensionCode.esriLicenseExtensionCodePublisher });
            //ESRI License Initializer generated code.

            test();               

            //Do not make any call to ArcObjects after ShutDownApplication()
            m_AOLicenseInitializer.ShutdownApplication();
        }

        public static void test()
        {
            //准备数据
            string polygon = @"D:\quickly\其它\矢量裁剪栅格\shp\LD.shp";
            string inputRastger = @"D:\quickly\其它\矢量裁剪栅格\temp\risk11.tif";

            ClipHelper c = new ClipHelper();//必须要在打开栅格之前调用，不然打开栅格会出错，还不知道原因。

            string errInfo = "";            
            IRasterDataset rasterDataset = new RasterDataset();
            rasterDataset.OpenFromFile(inputRastger);
            IWorkspaceFactory wsf = new ShapefileWorkspaceFactory();
            IFeatureWorkspace fws = wsf.OpenFromFile(System.IO.Path.GetDirectoryName(polygon), 0) as IFeatureWorkspace;
            IFeatureClass f = fws.OpenFeatureClass("LD");
            var feat = f.Search(null, false).NextFeature();

            var ret = c.ClipRaster(rasterDataset, feat.Shape as IPolygon, out errInfo);
        }
    }
}