Comprobación de si dos clases de características tienen la misma referencia espacial utilizando ArcPy

Question

Utilizando ArcPy, ¿cómo se comprueba si dos clases de características tienen la misma referencia espacial?

Comprobar simplemente si los dos son iguales no funciona:

>>> import arcpy
>>> fc1 = r"C:\Users\e1b8\Desktop\E1B8\GIS_Stackexchange\data.gdb\test"
>>> sr1 = arcpy.Describe (fc1).spatialReference 
>>> sr2 = arcpy.Describe (fc1).spatialReference
>>> sr1 == sr2
False

factoryCode no funciona, porque las proyecciones personalizadas no las tienen.

>>> fc2 = r"C:\Users\e1b8\Desktop\E1B8\GIS_Stackexchange\data.gdb\customproj"
>>> sr2 = arcpy.Describe (fc2).spatialReference
>>> sr2.factoryCode
0

Hay name pero los nombres pueden ser los mismos, pero tener diferentes unidades:

>>> sr1 = arcpy.Describe (fc1).spatialReference
>>> sr2 = arcpy.Describe (fc2).spatialReference
>>> sr1.name
u'NAD_1983_UTM_Zone_10N'
>>> sr2.name
u'NAD_1983_UTM_Zone_10N'
>>> sr1.linearUnitCode
9003
>>> sr2.linearUnitCode
9001

Así que se complica un poco. Lo mejor que se me ha ocurrido es:

>>> def CompareSRs (inFc1, inFc2):
    sr1 = arcpy.Describe (inFc1).spatialReference
    sr2 = arcpy.Describe (inFc2).spatialReference
    if not sr1.name != sr2.name:
        return False
    srType = sr1.type
    if srType != sr2.type:
        return False
    if srType == "Geographic":
        return sr1.angularUnitCode == sr2.angularUnitCode
    return sr1.linearUnitCode == sr2.linearUnitCode

Y todavía no estoy seguro de que el código anterior sea hermético.

¿Hay una forma mejor?

Answer 1

A juzgar por los comentarios, es posible que ya lo tengas :)

Podría comparar las descripciones de texto conocido (WKT) de las referencias espaciales.

sr1 = arcpy.Describe(dataset1).spatialReference
sr2 = arcpy.Describe(dataset2).spatialReference
sr1String = sr1.exportToString()
sr2String = sr2.exportToString()

matching = False

if sr1String == sr2String:
    # Exact string match
    matching = True
else:
    # difference
    pass

Answer 2

Versión de Arc: 10.3

Por si alguien sigue encontrando esto en 2019 , tuve problemas similares y quería estar lo más seguro posible de si las proyecciones coincidían. Como en la pregunta/respuesta anterior, puedes obtener la referencia espacial utilizando arcpy.Describe(dataset).spatialReference . En una biblioteca de funciones mía, integro esto en un flujo de trabajo, configurado para manejar la comparación de 2 conjuntos de datos.

Los atributos individuales de un objeto de referencia espacial de geoprocesamiento están disponibles aquí .

Las siguientes funciones deberían ayudar - siéntase libre de usar/modificar por supuesto. Vale la pena comprobar lo que se omite - algunos atributos de los sistemas de referencia espacial serán inofensivos si no coinciden, pero eso depende de usted :)

import arcpy

def check_crs(dataset):
    """Return a coordinate reference system string

    Get coordinate reference system of dataset
    """
    crs = arcpy.Describe(dataset).spatialReference
    return(crs)

def assert_crs_attribs(dataset1, dataset2, strict=False): 
    """Returns Nothing

    Asserts equality of all attributes of the provided geoprocessing spatial reference objects.
    These are generated using arcpy.Describe(your_dataset).spatialReference.
    Attributes of spatial reference object: https://pro.arcgis.com/en/pro-app/arcpy/classes/spatialreference.htm

    dataset1 - a spatial dataset with projection info e.g. shp
    dataset2 - a spatial dataset with projection info e.g. shp
    strict - boolean - if True will compare every element (default: False)
    """
    crs1=check_crs(dataset1)
    crs2=check_crs(dataset2)

    try:
        # Consider these
        assert(crs1.name==crs2.name) # The name of the spatial reference.
        assert(crs1.PCSCode==crs2.PCSCode) # The projected coordinate system code.1 
        assert(crs1.PCSName==crs2.PCSName) # The projected coordinate system name.1 
        assert(crs1.azimuth==crs2.azimuth) # The azimuth of a projected coordinate system.1 
        assert(crs1.centralMeridian==crs2.centralMeridian) # The central meridian of a projected coordinate system.1    
        assert(crs1.centralMeridianInDegrees==crs2.centralMeridianInDegrees) # The central meridian (Lambda0) of a projected coordinate system in degrees.1 
        assert(crs1.centralParallel==crs2.centralParallel) # The central parallel of a projected coordinate system.1
        assert(crs1.falseEasting==crs2.falseEasting) # The false easting of a projected coordinate system.1 
        assert(crs1.falseNorthing==crs2.falseNorthing) # The false northing of a projected coordinate system.1  
        assert(crs1.MFalseOriginAndUnits==crs2.MFalseOriginAndUnits) # The measure false origin and units.
        assert(crs1.MResolution==crs2.MResolution) # The measure resolution.
        assert(crs1.MTolerance==crs2.MTolerance) # The measure tolerance.
        assert(crs1.XYTolerance==crs2.XYTolerance) # The xy tolerance.
        assert(crs1.ZDomain==crs2.ZDomain) # The extent of the z domain.
        assert(crs1.ZFalseOriginAndUnits==crs2.ZFalseOriginAndUnits) # The z false origin and units.
        assert(crs1.factoryCode==crs2.factoryCode) # The factory code or well-known ID (WKID) of the spatial reference.
        assert(crs1.isHighPrecision==crs2.isHighPrecision) # Indicates whether the spatial reference has high precision set.
        assert(crs1.latitudeOf1st==crs2.latitudeOf1st) # The latitude of the first point of a projected coordinate system.1
        assert(crs1.latitudeOf2nd==crs2.latitudeOf2nd) # The latitude of the second point of a projected coordinate system.1    
        assert(crs1.latitudeOfOrigin==crs2.latitudeOfOrigin) # The latitude of origin of a projected coordinate system.1    
        assert(crs1.linearUnitCode==crs2.linearUnitCode) # The linear unit code.    
        assert(crs1.linearUnitName==crs2.linearUnitName) # The linear unit name.1
        assert(crs1.longitude==crs2.longitude) # The longitude value of this prime meridian.1
        assert(crs1.longitudeOf1st==crs2.longitudeOf1st) #The longitude of the first point of a projected coordinate system.1
        assert(crs1.longitudeOf2nd==crs2.longitudeOf2nd) # The longitude of the second point of a projected coordinate system.1
        assert(crs1.longitudeOfOrigin==crs2.longitudeOfOrigin) # The longitude of origin of a projected coordinate system.1
        assert(crs1.metersPerUnit==crs2.metersPerUnit) # The meters per linear unit.1
        assert(crs1.projectionCode==crs2.projectionCode) # The projection code.1
        assert(crs1.projectionName==crs2.projectionName) # The projection name.1
        assert(crs1.scaleFactor==crs2.scaleFactor) # The scale factor of a projected coordinate system.1
        assert(crs1.standardParallel1==crs2.standardParallel1) # The first parallel of a projected coordinate system.1
        assert(crs1.standardParallel2==crs2.standardParallel2) # The second parallel of a projected coordinate system.1
        assert(crs1.angularUnitCode==crs2.angularUnitCode) # The angular unit code.2
        assert(crs1.angularUnitName==crs2.angularUnitName) # The angular unit name.2
        assert(crs1.datumCode==crs2.datumCode) # The datum code.2
        assert(crs1.datumName==crs2.datumName) # The datum name.2
        assert(crs1.flattening==crs2.flattening) # The flattening ratio of this spheroid.2
        assert(crs1.longitude==crs2.longitude) # The longitude value of this prime meridian.2
        assert(crs1.primeMeridianCode==crs2.primeMeridianCode) # The prime meridian code.2

        ## Prob can be ignored
        if strict:
            assert(crs1.ZResolution==crs2.ZResolution) # The z resolution property.
            assert(crs1.ZTolerance==crs2.ZTolerance) # The z-tolerance property.
            assert(crs1.hasMPrecision==crs2.hasMPrecision) # Indicates whether m-value precision information has been defined.
            assert(crs1.hasXYPrecision==crs2.hasXYPrecision) # Indicates whether xy precision information has been defined.
            assert(crs1.hasZPrecision==crs2.hasZPrecision) # Indicates whether z-value precision information has been defined.
            assert(crs1.XYResolution==crs2.XYResolution) # The xy resolution.
            assert(crs1.domain==crs2.domain) # The extent of the xy domain.
            assert(crs1.MDomain==crs2.MDomain) # The extent of the measure domain.
            assert(crs1.remarks==crs2.remarks) # The comment string of the spatial reference.
            assert(crs1.type==crs2.type) # The type of the spatial reference. Geographic: A geographic coordinate system. Projected: A projected coordinate system.
            assert(crs1.usage==crs2.usage) # The usage notes.   
            assert(crs1.classification==crs2.classification) # The classification of a map projection.1 
            assert(crs1.GCSCode==crs2.GCSCode) # The geographic coordinate system code.2
            assert(crs1.GCSName==crs2.GCSName) # The geographic coordinate system name.2
            assert(crs1.primeMeridianName==crs2.primeMeridianName) # The prime meridian name.2
            assert(crs1.radiansPerUnit==crs2.radiansPerUnit) # The radians per angular unit.2
            assert(crs1.semiMajorAxis==crs2.semiMajorAxis) # The semi-major axis length of this spheroid.2
            assert(crs1.semiMinorAxis==crs2.semiMinorAxis) # The semi-minor axis length of this spheroid.2
            assert(crs1.spheroidCode==crs2.spheroidCode) # The spheroid code.2
            assert(crs1.spheroidName==crs2.spheroidName) # The spheroid name.2
        return(True)
    except:
        output_message="CRS differs between datasets."#\ncrs1: %s\ncrs2 : %s" %(crs1.exportToString(), crs2.exportToString())
        print(output_message)
        return(False)
    ## Differs to the falseEasting and falseNorthingUnits are odd on occasion but false eastings and northings make sense
    # crs.falseOriginAndUnits # The false origin and units.

    ## Not required
    #crs.GCS # A projected coordinate system returns a SpatialReference object for the geographic coordinate system it is based on. A geographic crs.coordinate system returns the same SpatialReference.
    #crs.SpatialReference
    #crs.VCS # If the coordinate system has a vertical coordinate system, it returns a VCS object for the vertical coordinate system it is based on.
    #crs.abbreviation # The abbreviated name of the spatial reference.
    #crs.alias # The alias of the spatial reference.

Teniendo en cuenta lo anterior, podrías utilizarlos como:

dataset1="your_vector_1.shp"
dataset2="your_vector_2.shp"

assert_crs_attribs(dataset1, dataset2)

Teniendo en cuenta tu caso de uso, es de esperar que los asserts no fallen.

Ahora integro estas funciones en muchos procesos, como cuando tengo una serie de conjuntos de datos espaciales que estoy uniendo y quiero eliminar cualquier duda de que las cosas estén mal alineadas.