Arcgis Engine向量裁剪柵格,呼叫Mask工具相關程式碼
阿新 • • 發佈:2019-01-26
今天試了下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); } } }