ae 生成等值线

2013-04-1719:47:14

---

首先用IInterpolationOp接口生成等值面的Raster，再导出成.jpg（或者指定格式的图片），在导出图片时按照图例进行渲染；
等值线则是在等值面Raster的基础上，通过ISurfaceOp接口生成的。

            DateTime dtbg = DateTime.Now;

            //得到站点要插值的要素类
            DataTable dt = GetLastDataTable(pName);
            IFeatureClass pInPointFClass = ModifySationFields(dt, pName);
            //工作空间路径
            string rasterPath = outPath + "\\TempData";
            //要素类图层
            IFeatureLayer pFeatLayer = UtilityFuntion.ReturnFeatureLayer(pInPointFClass);
            //这个值可以通过配置文件设置
            string pZValueField = string.Empty;
            if (this.pName == "自动站雨量")
            {
                pZValueField = "YULANG";
            }
            else if (this.pName == "自动站温度")
            {
                pZValueField = "WENDU";
            }
            else if (this.pName == "自动站气压")
            {
                pZValueField = "QIYA";
            }
            //pZValueField = "elevation";
            IFeatureClassDescriptor pFeatClsDes = ReturnFeatureClassDescriptor(pInPointFClass, pZValueField);
            //克里金插值
            IRaster pOutRaster = MapSpatialAnalyst.ExcuteRasterKrigie(pFeatClsDes, rasterPath, pFeatLayer);
            //IFeatureLayer featContour = null;
            if (pOutRaster == null)
            {
                return;
            }
            IRasterLayer pRasterLayer = new RasterLayerClass();

            pRasterLayer.CreateFromRaster(pOutRaster);
            pRasterLayer.Name = "yulang";
            string inputPath = pRasterLayer.FilePath;
            string outRasterPath = rasterPath + @"\yulang";

            IWorkspaceFactory pFWSF = new ShapefileWorkspaceFactoryClass();
            IFeatureWorkspace pFWS = pFWSF.OpenFromFile(rasterPath, 0) as IFeatureWorkspace;
            // 栅格重分类
            //ConversionTools pConversionTools = new ConversionTools(axMap.SpatialReference, pFWS);
            ConversionTools pConversionTools = new ConversionTools(axMap.SpatialReference, pFWS, rasterPath);
                    
            string maskPolygon = outPath +@"\Data\CQ.shp";
            //得到等值面要素类
            IFeatureClass pFeatureClass = pConversionTools.RasterToContourSurface(pOutRaster, maskPolygon, outRasterPath);

            string fieldName = "GRIDCODE";
           
            //新要素层
            IFeatureLayer pFeatureLayer = new FeatureLayerClass();
            pFeatureLayer.Name = "计算结果";
            pFeatureLayer.FeatureClass = pFeatureClass;
            MapTheme mapTheme = new MapTheme();
            UniqueValueModel model = new UniqueValueModel();
            model.FeatureLayer = pFeatureLayer;
            model.FieldName = fieldName;
            mapTheme.ExcuteUniqueValueRenderer(model);
            //添加图层
            axMap.AddLayer(pFeatureLayer);
            axMap.Refresh();
            axTOC.Update();
            DateTime dtend = DateTime.Now;
            MessageBox.Show(dtend.Subtract(dtbg).Seconds.ToString());   

/// <summary>
        /// 返回等值面
        /// </summary>
        /// <param name="pInRaster">栅格</param>
        /// <param name="pMaskPolygonPath">掩膜路径</param>
        /// <param name="outRasterPath">输出栅格路径</param>
        /// <returns></returns>
        public IFeatureClass RasterToContourSurface(IRaster pInRaster, string pMaskPolygonPath,string outRasterPath)
        {
            IRaster pOutRas = ReclassifyRaster(pInRaster);
            //string sRasterToPolygonLayerName = "MyRasterToPolygon";
            IFeatureClass pFeatureClass = null;
            if (pOutRas != null)
            {
                RasterByMaskPolygon(pOutRas, pMaskPolygonPath, outRasterPath);
            }
            IRasterLayer pRasterLayer = new RasterLayerClass();

            pRasterLayer.CreateFromFilePath(outRasterPath);

            string sRasterToPolygonLayerName = "MyRasterToPolygon";
            pFeatureClass=RasterToPolygon((pRasterLayer.Raster as IRaster2).RasterDataset,sRasterToPolygonLayerName);
            return pFeatureClass;
        }

/// <summary>
        /// 栅格数据重分类
        /// </summary>
        /// <param name="inRaster"></param>
        /// <returns></returns>
        private IRaster ReclassifyRaster(IRaster inRaster)
        {
            if (inRaster == null)
            {
                return null;
            }
            IReclassOp pReclassOp = new RasterReclassOpClass();
            IRasterAnalysisEnvironment pRAE = pReclassOp as IRasterAnalysisEnvironment;
            pRAE.OutWorkspace = m_pMemoryWS;
            pRAE.OutSpatialReference = m_pSpatialReference;
            IGeoDataset pGeodataset = inRaster as IGeoDataset;

            IRasterBandCollection pRsBandCol = pGeodataset as IRasterBandCollection;
            IRasterBand pRasterBand = pRsBandCol.Item(0);
            pRasterBand.ComputeStatsAndHist();
            //统计最大最小值
            IRasterStatistics pRasterStatistic = pRasterBand.Statistics;
            double dMaxValue = pRasterStatistic.Maximum;
            double dMinValue = pRasterStatistic.Minimum;

            double lInterval = Convert.ToDouble((dMaxValue - dMinValue) / 10);
            double lCurrentContour = dMinValue;
            INumberRemap pNumberRemap = new NumberRemapClass();
            int i = 0;
            while (lCurrentContour <= dMaxValue)
            {
                pNumberRemap.MapRange(lCurrentContour, lCurrentContour + lInterval, i + 1);
                lCurrentContour += lInterval;
                i++;
            }
            IGeoDataset pOutGeodataset = pReclassOp.ReclassByRemap(inRaster as IGeoDataset, pNumberRemap as IRemap, false);
            inRaster = null;
            System.GC.Collect();
            return pOutGeodataset as IRaster;
        }

/// <summary>
        /// 栅格转成等值面掩模过程
        /// </summary>
        /// <param name="pInRaster">输入栅格</param>
        /// <param name="pMaskPolygon">掩模</param>
        /// <returns></returns>
        public void RasterByMaskPolygon(IRaster pInRaster, string pMaskPolygonPath,string outRasterPath)
        {           
            //string sRasterToPolygonLayerName = "MyRasterToPolygon";
            //IFeatureClass pFeatureClass = null;
            if (pInRaster != null)
            {
                IRasterLayer pRasterLayer = new RasterLayerClass();

                pRasterLayer.CreateFromRaster(pInRaster);
                string inRasterPath = pRasterLayer.FilePath;
                ReturnRasterByMask(inRasterPath, pMaskPolygonPath, outRasterPath);
            }
            //return pFeatureClass;
        }

/// <summary>
        /// 执行掩膜操作
        /// </summary>
        /// <param name="inRasterPath">栅格输入路径</param>
        /// <param name="maskPath">多边形</param>
        /// <param name="outRasterPath">输出栅格路径</param>
        private void ReturnRasterByMask(string inRasterPath,string maskPath,string outRasterPath)
        {
            //调用Toolbox中的Reclassify
            Geoprocessor pGeoPro = new Geoprocessor();
            pGeoPro.OverwriteOutput = true;

            ExtractByMask pExtractByMask = new ExtractByMask();
            pExtractByMask.in_raster = inRasterPath;
            pExtractByMask.in_mask_data = maskPath;
            pExtractByMask.out_raster = outRasterPath;
            pGeoPro.Execute(pExtractByMask, null);
            ReturnMessages(pGeoPro);
        }

/// <summary>
        /// 栅格转换为 Polygon FeatureClass
        /// </summary>
        /// <param name="pInRasterDataset"></param>
        /// <returns></returns>
        //public IFeatureClass RasterToPolygon(IRasterDataset pInRasterDataset, IFeatureWorkspace pOutputFeatureWorkspace,string sOutFeatureClassName)
        private IFeatureClass RasterToPolygon(IRasterDataset pInRasterDataset, string sOutFeatureClassName)
        {

            IConversionOp pConversionOp = new RasterConversionOpClass();
            IGeoDataset pInGeoDataset = pInRasterDataset as IGeoDataset;

            bool b = false;
            IGeoDataset pOutGeoDataset = null;
            if (m_pTempFeatWorkspace == null)
            {
                b = DeleteDatasetIfExist(m_pMemoryWS, sOutFeatureClassName);   //如果已经存在，则删除
                pOutGeoDataset = pConversionOp.RasterDataToPolygonFeatureData(pInGeoDataset, m_pMemoryWS, sOutFeatureClassName, false);
            }
            else
            {
                b = DeleteDatasetIfExist(m_pTempFeatWorkspace as IWorkspace, sOutFeatureClassName);   //如果已经存在，则删除
                pOutGeoDataset = pConversionOp.RasterDataToPolygonFeatureData(pInGeoDataset, m_pTempFeatWorkspace as IWorkspace, sOutFeatureClassName, false);
            }

            IDataset pDataset = pOutGeoDataset as IDataset;
            if (pDataset.CanRename())
            {
                pDataset.Rename(sOutFeatureClassName);
            }
            System.GC.Collect();
            return pOutGeoDataset as IFeatureClass;
        }