Trouver la distance minimale entre les bords des polygones à l'aide d'ArcGIS Desktop?

9

J'ai une carte d'environ 3000 polygones dans ArcGIS 10. Je cherche à trouver la distance entre chacun d'eux. Je sais comment le faire en utilisant les coordonnées lat et longue du centroïde, mais je recherche la distance en ligne droite la plus courte entre le bord le plus proche d'un polygone et le bord le plus proche de l'autre polygone. Des idées?

Kevin
la source

Réponses:

11

C'est un bon morceau de code, mais pas aussi agréable que (en supposant que votre table est en coordonnées géographiques, sinon supprimez simplement les transformations géographiques)

CREATE TABLE mytable_distances AS
SELECT a.id, b.id, ST_Distance(a.geom::geography, b.geom::geography) as distance
FROM mytable a, mytable b;

Ai-je mentionné que les bases de données spatiales sont géniales? Ils font. Oh, ils le font.

Paul Ramsey
la source
Cela permettra de trouver la distance entre les sommets les plus proches, mais pas les bords eux-mêmes - il ne semble pas que GEOS expose cette réponse plus précise. Pourtant, assez pratique!
scw
1
Désolé scw, vous vous trompez à bien des égards. PostGIS a des calculs de distance natifs. GOES n'y participe pas. Deuxièmement, il donne absolument la distance la plus proche entre les bords, non seulement les sommets à la fois dans la distance géométrique et dans le calcul de la distance sphéroïdale de type géographique. Paul l'a écrit.
Nicklas Avén
Pour le voir visuellement pour la géométrie, vous pouvez utiliser st_shortestline qui renvoie la ligne qui donne la distance.
Nicklas Avén
1
Nik a raison, tant en géométrie qu'en géographie, la fonction distance renvoie la distance entre les bords. Par exemple, sélectionnez st_distance ('LINESTRING (0 0, 0 100)', 'LINESTRING (50 1, 51 1)')
Paul Ramsey
2
wow, les bases de données spatiales font du rock! je calcule la distance entre un ensemble de ~ 8200 polygones et le plus proche voisin dans un autre ensemble de ~ 8400 polygones. dans arcgis 10, l'outil 'générer près de la table' avec un rayon de recherche de 10000 m a pris 1 heure et 15 minutes (sur le bureau i7 quad-core 3,4 GHz). la même requête dans PostGIS n'a pris que 3,5 minutes, et c'était sur un ordinateur plus lent (un macbook pro i7 dual-core à 2,7 GHz).
pistachionut
8

La distance de A à B est identique à B à A, et la distance de A à A est nulle, donc une demi-matrice vous fera économiser du travail.

IProximityOperator renvoie la distance à partir du bord. Le code ci-dessous utilise une projection azimutale centrée sur le centroïde de chaque polygone (devrait également fonctionner avec les lignes). Si les polygones ne sont pas trop complexes (ou si vous avez beaucoup de mémoire) charger toutes les géométries en mémoire et les projeter serait plus rapide. (Ce n'est pas testé à fond).

public class Pair
{
    public int Oid1;
    public int Oid2;
    public double Dist;
    public static void TestGetDistances()
    {
        IWorkspaceFactory wsf = new ESRI.ArcGIS.DataSourcesGDB.FileGDBWorkspaceFactoryClass();

        string path = @"C:\Program Files\ArcGIS\DeveloperKit10.0\Samples\data\Usa\USA.gdb";
        var fws = wsf.OpenFromFile(path, 0) as IFeatureWorkspace;
        IFeatureClass fc = fws.OpenFeatureClass("states");
        var halfMatrix = Pair.GetPairs(fc);

    }
    /// <summary>
    /// key is oid of each feature, value is pairs for features with smaller oids.
    /// </summary>
    /// <param name="fc"></param>
    /// <returns></returns>
    public static SortedList<int, List<Pair>> GetPairs(IFeatureClass fc)
    {
        ISpatialReferenceFactory3 srf = new SpatialReferenceEnvironmentClass();
        IProjectedCoordinateSystem pcs = 
        srf.CreateProjectedCoordinateSystem((int)esriSRProjCSType.esriSRProjCS_WGS1984N_PoleAziEqui);

        var outList = new SortedList<int, List<Pair>>();
        IFeatureCursor fCur = fc.Search(null, true);
        IFeature f;
        while ((f = fCur.NextFeature()) != null)
        {
            var pairs = GetDistances(f, pcs);
            Debug.Print("{0} has {1} pairs", f.OID, pairs.Count);
            outList.Add(f.OID, pairs);
        }
        System.Runtime.InteropServices.Marshal.FinalReleaseComObject(fCur);
        return outList;
    }

    private static IPoint GetGCSCentroid(IGeometry geom)
    {
        if (geom.SpatialReference is IProjectedCoordinateSystem)
        {
            geom.Project(((IProjectedCoordinateSystem)geom.SpatialReference).GeographicCoordinateSystem);
        }
        IArea a = geom is IArea ? geom as IArea : geom.Envelope as IArea;
        return a.Centroid;
    }

    /// <summary>
    /// return a list of all other features whose OID is lesser than f1
    /// </summary>
    /// <param name="f1"></param>
    /// <param name="pcs"></param>
    /// <returns></returns>
    private static List<Pair> GetDistances(IFeature f1, IProjectedCoordinateSystem pcs)
    {
        IPoint centroid = GetGCSCentroid(f1.ShapeCopy);

        pcs.set_CentralMeridian(true, centroid.X);
        ((IProjectedCoordinateSystem2)pcs).LatitudeOfOrigin = centroid.Y;
        var g1 = f1.ShapeCopy;
        g1.Project(pcs);

        var outList = new List<Pair>();
        var fc = f1.Class as IFeatureClass;
        var proxOp = g1 as IProximityOperator;
        IFeatureCursor fCur = fc.Search(null, true);
        IFeature f2 = null;
        while ((f2 = fCur.NextFeature()) != null)
        {
            if (f2.OID < f1.OID)
            {
                var g2 = f2.ShapeCopy;
                g2.Project(pcs);
                outList.Add(new Pair()
                {
                    Oid1 = f1.OID,
                    Oid2 = f2.OID,
                    Dist = proxOp.ReturnDistance(g2)
                });
            }
        }
        System.Runtime.InteropServices.Marshal.FinalReleaseComObject(fCur);
        return outList;
    }
}
Kirk Kuykendall
la source
c'est un joli morceau de code. Je ne connaissais pas IproximityOperator, et j'ai fini par coder quelque chose comme ça moi-même (évidemment c'est plus lent)
George Silva
2

Je pense que l' outil de table proche fonctionnerait pour ce que vous voulez:

Détermine les distances entre chaque entité dans les entités en entrée et une ou plusieurs entités proches dans les entités proches, dans le rayon de recherche. Les résultats sont enregistrés dans le tableau de sortie.

Laissez simplement le rayon de recherche vide.

Jason Scheirer
la source
C'est la solution que j'essaierais en premier, mais elle nécessite un niveau de licence ArcInfo pour déverrouiller l'outil Générer une table proche (analyse).
PolyGeo