Example of handling overlapping polygons when using ArcObjects IZonalOps

overlapping_polygon_example.cs

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

using System.Collections.Specialized;
using System.Runtime.InteropServices;

using ESRI.ArcGIS.esriSystem;
using ESRI.ArcGIS.Server;
using ESRI.ArcGIS.Geometry;
using ESRI.ArcGIS.Geodatabase;
using ESRI.ArcGIS.Carto;
using ESRI.ArcGIS.SOESupport;
using ESRI.ArcGIS.SpatialAnalyst;

using ESRI.ArcGIS.GeoAnalyst;
using ESRI.ArcGIS.DataSourcesRaster;

using Newtonsoft.Json;
using Newtonsoft.Json.Linq;
using Newtonsoft.Json.Converters;
using Newtonsoft.Json.Serialization;
using Newtonsoft.Json.Bson;
using Newtonsoft.Json.Schema;
using Newtonsoft.Json.Utilities;

namespace suitability_stats_soe
{
    /*
     * Darn Esri Zonal Tools don't handle overlapping polygons....I also wanted results as JSON.NET JArrays...
     * 
     */
    public class ZonalOps2 : IDisposable
    {
        private IZonalOp zonalOp;

        public ZonalOps2()
        {
            this.zonalOp = new RasterZonalOpClass();
        }

        public void Dispose()
        {
            Marshal.FinalReleaseComObject(this.zonalOp);
        }

        public JArray ZonalStatisticsAsTable2(IMapServer4 mapServer, IFeatureClass zones, int valueRasterIndex, string zoneField, bool ignoreNoData = true)
        {
            List<object> comObjectsToRelease = new List<object>();

            IGeoDataset valueRaster = null;
            comObjectsToRelease.Add(valueRaster);

            try
            {
                IMapServerDataAccess dataAccess = (IMapServerDataAccess)mapServer;
                valueRaster = dataAccess.GetDataSource(mapServer.DefaultMapName, valueRasterIndex) as IGeoDataset;
                return this.ZonalStatisticsAsTable2(zones, valueRaster, zoneField, ignoreNoData);
            }

            finally
            {
                this.ReleaseComObjects(comObjectsToRelease);
            }
        }

        //divides input features into non-overlapping groups, runs ZonalStatisticsAsTable for each group, and returns single JSON.NET JArray
        public JArray ZonalStatisticsAsTable2(IFeatureClass zonesFeatureClass, IGeoDataset valueRaster, string zoneField, bool ignoreNoData = true, IWorkspace outputWorkspace = null)
        {
            JArray resultStats = new JArray();

            List<object> comObjectsToRelease = new List<object>();

            IFeatureClassDescriptor featureClassDescriptor = new FeatureClassDescriptorClass() as IFeatureClassDescriptor;
            comObjectsToRelease.Add(featureClassDescriptor);

            ITable areaTable = null;
            comObjectsToRelease.Add(areaTable);

            IQueryFilter queryFilter = null;
            comObjectsToRelease.Add(queryFilter);

            try
            {
                List<List<int>> nonOverlappingGroups = this.CreateNonOverlappingGroups(zonesFeatureClass, zoneField);

                IGeoDataset zonesGeodatset = (IGeoDataset)zonesFeatureClass;

                this.SetRasterEnvironmentVaribles(zonesGeodatset.Extent, -1, zonesGeodatset, valueRaster, outputWorkspace, zonesGeodatset.SpatialReference);

                foreach (List<int> group in nonOverlappingGroups)
                {
                    queryFilter = new QueryFilterClass();
                    queryFilter.WhereClause = String.Format("\"{0}\" IN ({1})", zoneField, string.Join(",", group.Select(n => n.ToString()).ToArray()));
                    featureClassDescriptor.Create(zonesFeatureClass, queryFilter, zoneField);

                    areaTable = this.zonalOp.ZonalStatisticsAsTable((IGeoDataset)featureClassDescriptor, valueRaster, ignoreNoData);
                    JArray tableArray = this.ITableToJArray(areaTable);
                    for (int r = 0; r < tableArray.Count; r++)
                    {
                        resultStats.Add(tableArray[r]);
                    }
                }

                return resultStats;
            }

            finally
            {
                this.ReleaseComObjects(comObjectsToRelease);
            }
        }

        //divides input features into non-overlapping groups, runs TabulateArea for each group, and returns single JSON.NET JArray
        public JArray TabulateArea2(IFeatureClass zonesFeatureClass, IGeoDataset valueRaster, string uniqueIntegerIdField, IWorkspace outputWorkspace = null)
        {
            JArray resultStats = new JArray();

            List<object> comObjectsToRelease = new List<object>();

            IFeatureClassDescriptor featureClassDescriptor = new FeatureClassDescriptorClass() as IFeatureClassDescriptor;
            comObjectsToRelease.Add(featureClassDescriptor);

            ITable areaTable = null;
            comObjectsToRelease.Add(areaTable);

            IQueryFilter queryFilter = null;
            comObjectsToRelease.Add(queryFilter);

            try
            {
                List<List<int>> nonOverlappingGroups = this.CreateNonOverlappingGroups(zonesFeatureClass, uniqueIntegerIdField);

                IGeoDataset zonesGeodatset = (IGeoDataset)zonesFeatureClass;

                this.SetRasterEnvironmentVaribles(zonesGeodatset.Extent, -1, zonesGeodatset, valueRaster, outputWorkspace, zonesGeodatset.SpatialReference);

                foreach (List<int> group in nonOverlappingGroups)
                {
                    queryFilter = new QueryFilterClass();
                    queryFilter.WhereClause = String.Format("\"{0}\" IN ({1})", uniqueIntegerIdField, string.Join(",", group.Select(n => n.ToString()).ToArray()));
                    featureClassDescriptor.Create(zonesFeatureClass, queryFilter, uniqueIntegerIdField);

                    areaTable = this.zonalOp.TabulateArea((IGeoDataset)featureClassDescriptor, valueRaster);
                    JArray tableArray = this.ITableToJArray(areaTable);
                    for (int r = 0; r < tableArray.Count; r++)
                    {
                        resultStats.Add(tableArray[r]);
                    }
                }

                return resultStats;
            }

            finally
            {
                this.ReleaseComObjects(comObjectsToRelease);
            }


        }

        //converts Esri ITable to JSON.NET JArray
        public JArray ITableToJArray(ITable table)
        {
            JArray resultArray = new JArray();

            List<object> comObjectsToRelease = new List<object>();

            IFields fields = table.Fields;
            comObjectsToRelease.Add(fields);

            ICursor cursor = table.Search(null, false);
            comObjectsToRelease.Add(cursor);

            IRow row = null;
            comObjectsToRelease.Add(row);

            IField field = null;
            comObjectsToRelease.Add(field);

            try
            {
                while ((row = cursor.NextRow()) != null)
                {
                    JObject jobject = new JObject();
                    for (int f = 0; f < fields.FieldCount; f++)
                    {
                        field = fields.get_Field(f);
                        if (field.Type == esriFieldType.esriFieldTypeInteger || field.Type == esriFieldType.esriFieldTypeSmallInteger || field.Type == esriFieldType.esriFieldTypeOID)
                        {
                            jobject[field.Name] = (int)row.get_Value(f);
                        }

                        else if (field.Type == esriFieldType.esriFieldTypeDouble)
                        {
                            jobject[field.Name] = (double)row.get_Value(f);
                        }

                        else if (field.Type == esriFieldType.esriFieldTypeSingle)
                        {
                            jobject[field.Name] = (float)row.get_Value(f);
                        }

                        else if (field.Type == esriFieldType.esriFieldTypeString)
                        {
                            jobject[field.Name] = (string)row.get_Value(f);
                        }

                        else
                        {
                            throw new Exception(String.Format("Unable to handle field type {0} ==> {1}", field.Name, field.Type));
                        }
                    }

                    resultArray.Add(jobject);
                }

                return resultArray;
            }

            finally
            {
                this.ReleaseComObjects(comObjectsToRelease);
            }
        }

        //returns arrays of non-overlapping groups indicated by the ids from the uniqueId Field: 
        public List<List<int>> CreateNonOverlappingGroups(IFeatureClass zonesFeatureClass, string uniqueIntegerIdField)
        {
            List<List<int>> finalNonOverlappingGroups = new List<List<int>>();

            List<object> comObjectsToRelease = new List<object>();

            IGeometry geometryBag = new GeometryBagClass();
            comObjectsToRelease.Add(geometryBag);

            IGeoDataset zonesGeodatset = (IGeoDataset)zonesFeatureClass;
            geometryBag.SpatialReference = zonesGeodatset.SpatialReference;

            IFeatureCursor featureCursor = zonesFeatureClass.Search(null, false);
            comObjectsToRelease.Add(featureCursor);

            IGeometryCollection geometryCollection = geometryBag as IGeometryCollection;
            comObjectsToRelease.Add(geometryCollection);

            IFeature currentFeature = featureCursor.NextFeature();
            comObjectsToRelease.Add(currentFeature);

            try
            {
                if (zonesFeatureClass == null)
                {
                    return null;
                }

                int uniqueIdFieldIndex = -1;

                List<int> featureIds = new List<int>();
                while (currentFeature != null)
                {
                    if (uniqueIdFieldIndex == -1)
                    {
                        uniqueIdFieldIndex = currentFeature.Fields.FindField(uniqueIntegerIdField);
                        if (uniqueIdFieldIndex == -1)
                        {
                            throw new Exception(String.Format("GetOverlappingGroupsById :: uniqueIntegerIdField not found: {0}", uniqueIntegerIdField));
                        }

                        if (currentFeature.Fields.get_Field(uniqueIdFieldIndex).Type != esriFieldType.esriFieldTypeInteger)
                        {
                            throw new Exception(String.Format("GetOverlappingGroupsById :: Field {0} must be of type Integer", uniqueIntegerIdField));
                        }
                    }

                    object missing = Type.Missing;
                    geometryCollection.AddGeometry(currentFeature.Shape, ref missing, ref missing);
                    featureIds.Add(Convert.ToInt32(currentFeature.get_Value(uniqueIdFieldIndex)));
                    currentFeature = featureCursor.NextFeature();
                }

                IRelationalOperatorNxM overlapOperator = (IRelationalOperatorNxM)geometryBag;
                comObjectsToRelease.Add(overlapOperator);

                IRelationResult overlappingResult = overlapOperator.Overlaps((IGeometryBag)geometryBag);
                comObjectsToRelease.Add(overlappingResult);

                int count = overlappingResult.RelationElementCount;
                int left, right;
                Dictionary<int, List<int>> overlapIndex = new Dictionary<int, List<int>>();
                for (int i = 0; i < count; i++)
                {
                    overlappingResult.RelationElement(i, out left, out right);
                    if (!overlapIndex.ContainsKey(left))
                    {
                        overlapIndex[left] = new List<int>();
                    }
                    overlapIndex[left].Add(right);
                }

                List<KeyValuePair<List<int>, List<int>>> nonOverlappingGroups = new List<KeyValuePair<List<int>, List<int>>>();
                nonOverlappingGroups.Add(CreateEmptyOverlapEntry());

                for (int j = 0; j < featureIds.Count; j++)
                {
                    if (!overlapIndex.ContainsKey(j))
                    {
                        nonOverlappingGroups[0].Key.Add(featureIds[j]);
                        continue;
                    }
                    else
                    {
                        bool groupFound = false;
                        foreach (KeyValuePair<List<int>, List<int>> pair in nonOverlappingGroups)
                        {
                            if (!pair.Value.Contains(j))
                            {
                                pair.Key.Add(featureIds[j]);
                                pair.Value.AddRange(overlapIndex[j]);
                                groupFound = true;
                                break;
                            }
                        }

                        if (!groupFound)
                        {
                            KeyValuePair<List<int>, List<int>> newGroup = CreateEmptyOverlapEntry();
                            newGroup.Key.Add(featureIds[j]);
                            newGroup.Value.AddRange(overlapIndex[j]);
                            nonOverlappingGroups.Add(newGroup);
                        }
                    }
                }

                foreach (KeyValuePair<List<int>, List<int>> pair in nonOverlappingGroups)
                {
                    finalNonOverlappingGroups.Add(pair.Key);
                }

                return finalNonOverlappingGroups;
            }

            finally
            {
                this.ReleaseComObjects(comObjectsToRelease);
            }
        }

        public void SetRasterEnvironmentVaribles(IEnvelope envelope_Extent = null,
                                                double nCellSize = -1,
                                                IGeoDataset geoDataset_Mask = null,
                                                IGeoDataset SnapRaster = null,
                                                IWorkspace workspace = null,
                                                ISpatialReference spatialReference = null)
        {
            IRasterAnalysisEnvironment rasterAnalysisEnvironment = (IRasterAnalysisEnvironment)this.zonalOp;

            //extent.
            if (envelope_Extent != null)
            {
                object extentProvider = (object)envelope_Extent;
                object snapRasterData = (object)SnapRaster;
                rasterAnalysisEnvironment.SetExtent
                    (esriRasterEnvSettingEnum.esriRasterEnvValue, ref extentProvider,
                    ref snapRasterData);
            }

            //cell size.
            if (nCellSize > 0)
            {
                object cellSizeProvider = (object)nCellSize;
                rasterAnalysisEnvironment.SetCellSize
                    (esriRasterEnvSettingEnum.esriRasterEnvValue, ref cellSizeProvider);
            }

            //mask.
            if (geoDataset_Mask != null)
            {
                rasterAnalysisEnvironment.Mask = geoDataset_Mask;
            }

            //output workspace.
            if (workspace != null)
            {
                rasterAnalysisEnvironment.OutWorkspace = workspace;
            }

            //spatial reference.
            if (spatialReference != null)
            {
                rasterAnalysisEnvironment.OutSpatialReference = spatialReference;
            }
        }

        public bool WithinAreaThreshold(IArea area, int propsedCellSize, int minCellCount = 25)
        {
            return true;
        }

        //KeyValuePair< <ids in group> : <ids forbidden from group> >
        private KeyValuePair<List<int>, List<int>> CreateEmptyOverlapEntry()
        {
            List<int> groupsIds = new List<int>();
            List<int> outlawIds = new List<int>();
            return new KeyValuePair<List<int>, List<int>>(groupsIds, outlawIds);
        }

        private void ReleaseComObjects(List<object> comObjects)
        {
            foreach (object item in comObjects)
            {
                if (item != null && Marshal.IsComObject(item))
                {
                    Marshal.FinalReleaseComObject(item);
                }
            }

            comObjects.Clear();
        }
    }
}