Select to view content in your preferred language

calculating coordinates from a table, non spatial process.

1366
8
Jump to solution
07-18-2023 01:02 PM
clt_cabq
Frequent Contributor

I have a table to process that with coordinates for records expressed in two different ways - one set are GCS, the other are State Plane coordinates. I'd like to simply calculate new coordinates into a new field so that there is one common set of coordinates for all records. I can do this through a multistep process involving converting recordsets of each coordinate type to spatial layers, projecting one of them, merging them back together and calculating the coordinate values. These data ultimately will become spatial, but in an early stage of processing, it would be nice to have a simpler (aka fewer steps) process. I suspect i could go find a formula for this somewhere along the way and implement that in Python, but this seems like it could be more of a pain that the multistep process i outlined above.

0 Kudos
1 Solution

Accepted Solutions
by Anonymous User
Not applicable

You could step through the features and use pyproj to transform. Not sure how your data is set up or anything really, but add two fields and use an update cursor.

import arcpy
import pyproj

fc = r''

# Add the lon/lat fields
arcpy.AddFields_management(fc, [['converted_lon', 'FLOAT'], ['converted_lat', 'FLOAT']])

# Convert coordinate system and add XY Coordinates in WGS 1984
with arcpy.da.UpdateCursor(fc, ['GCS_x', 'GCS_y', 'sp_x', 'sp_y', 'converted_lon', 'converted_lat']) as uCur:
    for row in uCur:
        # Already in GCS?, but if needed use the pyproj step below
        if row[0]:
            row[4] = row[0]
            row[5] = row[1]
        elif row[2]:
            # Stateplane to WGS84
            stplne = pyproj.Proj(init='epsg:32133')
            wgs84 = pyproj.Proj(init='epsg:4326')
            row[4], row[5] = pyproj.transform(stplne, wgs84, row[2], row[3])

    uCur.updateRow(row)

 

View solution in original post

8 Replies
AlfredBaldenweck
MVP Regular Contributor

Try Convert Coordinate Notation (Data Management)—ArcGIS Pro | Documentation

The only issue is that it doesn't do it in place.

0 Kudos
clt_cabq
Frequent Contributor

Thanks, I fussed with this approach a bit and  still puts me in the position of having to create intermediate tables and layers. Appreciate the suggestion! I found another post that more or less outlines how to go about doing these via formula but it looked pretty complex and lots of room for error. 

0 Kudos
TylerT
by
Frequent Contributor

Do you want to perform the projection in arcpy (gdb's) or arcgis (SEDF)?  If arcpy, you could write to a memory gdb instead of gdb directly.  This might save some of your concern about extra tables.

Alternatively, you could use the arcgis API for python and a spatially enabled dataframe (SEDF).   Two possible suggestions:
1) One SEDF.  Use pandas apply to pass the geometry project method on each geometry object.  This should work even if you have mixed starting spatial references.  You can read the in_sr with the SpatialReference method into a side column.  You might have to create geometry objects first with similar pandas apply.  Docs at https://developers.arcgis.com/python/api-reference/arcgis.geometry.html#project.


2) Multiple SEDFs, one SEDF per spatial reference group of objects.  Use the GeoSeriesAccessor to set the sr of the entire SEDF shape column then use the project_as method to your desired sr.  Finally, concatenate the SEDFs.  Docs at https://developers.arcgis.com/python/api-reference/arcgis.features.toc.html#geoseriesaccessor.

HTH
Tyler

0 Kudos
clt_cabq
Frequent Contributor

Thanks Tyler! This is work on an enterprise GDB before the layer gets published so the SEDF probably isn't my solution here, but that's an interesting approach. I need to learn more about the API for Python so maybe this will give me an excuse to poke around at it.

ct

0 Kudos
by Anonymous User
Not applicable

You could step through the features and use pyproj to transform. Not sure how your data is set up or anything really, but add two fields and use an update cursor.

import arcpy
import pyproj

fc = r''

# Add the lon/lat fields
arcpy.AddFields_management(fc, [['converted_lon', 'FLOAT'], ['converted_lat', 'FLOAT']])

# Convert coordinate system and add XY Coordinates in WGS 1984
with arcpy.da.UpdateCursor(fc, ['GCS_x', 'GCS_y', 'sp_x', 'sp_y', 'converted_lon', 'converted_lat']) as uCur:
    for row in uCur:
        # Already in GCS?, but if needed use the pyproj step below
        if row[0]:
            row[4] = row[0]
            row[5] = row[1]
        elif row[2]:
            # Stateplane to WGS84
            stplne = pyproj.Proj(init='epsg:32133')
            wgs84 = pyproj.Proj(init='epsg:4326')
            row[4], row[5] = pyproj.transform(stplne, wgs84, row[2], row[3])

    uCur.updateRow(row)

 

clt_cabq
Frequent Contributor

ah, I should have known there is a library for this, I'll check that out. I've already written the code to create two spatial layers, perform a projection on one then bring them back together, but using this library may be useful in simplifying the process.

0 Kudos
Tomasz_Tarchalski
Occasional Contributor

Hi,

I recommend strongly pyproj. 

 

Tomasz_Tarchalski_0-1690788112843.png

 

clt_cabq
Frequent Contributor

Thanks to both Tomasz and JeffK - I've implemented the pyproj module to handle this in my table and it works perfectly well. Ultimately, I think this will help streamline some downstream processes and just makes everything more straightforward. One thing I am noticing is when i run pyproj I get some warnings that I think are saying the transform method is deprecated, so the module available to import from the ArcGIs Pro package manager may be out of date. I have to read those warnings a bit closer and chase that down.

Edit: something to note here if you use this approach - the deprecated method illustrated by JeffK works, but the method shown by Tomasz is significantly faster. Working with 10K records using an insert cursor the first method took 5-6 minutes to run, the second method took under a minute, and if you read the 'gotchas' section in the pyproj docs, that holds up given some timing reported there.

0 Kudos