topic Calc Latitude, Longitude Field in Python Questions

Calc Latitude, Longitude Field

BlakeTerhune — Sat, 11 Dec 2021 20:57:27 GMT

After getting a lot of requests to get latitude and longitude coordinates for features, I created a generic script that creates two new fields in the existing feature feature class and and calcs the latitude and longitude of the feature's "true centroid." The only problem is that you would have to run this script any time a feature's geometry was changed. I did some searching for something that already existed and didn't find anything that was simple and generic enough, so please let me know if you can do it better.

import arcpy

fc = r"C:\temp\working.gdb\temp_fc"

# Define and create longitude/latitude fields
## Longitude is X, Latitude is Y
fieldsNew = ["GCS_WGS_1984_X", "GCS_WGS_1984_Y"]
fieldsExisting = {f.name: f.type for f in arcpy.ListFields(fc)}
for fNew in fieldsNew:
    # Add field if it doens't already exist
    if fNew not in fieldsExisting.keys():
        arcpy.AddField_management(
            fc,  ## in_table
            fNew,  ## field_name
            "DOUBLE",  ## field_type
        )
    # If data type of existing field is not a Double, throw Exception
    elif fieldsExisting[fNew] != "Double":
        raise Exception(
            "Existing {} field is type {}. Must be 'Double'".format(
                fNew, fieldsExisting[fNew]
            )
        )
# Update fields with GCS_WGS_1984 lat/long
print("Updating Field values")
fieldsNew.append("SHAPE@")
GCS_WGS_1984 = arcpy.SpatialReference(4326)
with arcpy.da.UpdateCursor(fc, fieldsNew, "", GCS_WGS_1984) as u_cursor:
    for row in u_cursor:
        point = row[2].trueCentroid
        row[0] = point.X
        row[1] = point.Y
        u_cursor.updateRow(row)

python snippets

Re: Calc Latitude, Longitude Field

XanderBakker — Mon, 08 Jun 2015 23:02:55 GMT

Did you look at: Add Geometry Attributes (Data Management)?

Question: what if your input featureclass does not have a sr WGS 1984 and it needs a specific transformation? Will the output coordinates be correct?

For polylines I prefer to use the midpoint of the line (positionAlongLine (0.5, True)) over the center of gravity...

Re: Calc Latitude, Longitude Field

BlakeTerhune — Mon, 08 Jun 2015 23:39:17 GMT

Great points, Xander. I don't think I saw AddGeometryAttributes_management(). Or maybe I did but didn't see that you could specify the coordinate system. I'll do some testing with that!

I have not needed to use a a specific transformation so have not tested it. What is a scenario I could use to test this? I don't understand what you mean when you say it "does not have a sr WGS 1984."

For the lines, I agree, a midpoint would be better. I just didn't want to write extra code to check the geometry type and do different things. Point features are obviously preferred so anything else requires some compromise. That would be a good enhancement though, to do different calculations of feature to point conversion based on the geometry type.

I also thought maybe instead of hard-coding the field name, I could just make a new field using the current date/time in the name. That would let me get rid of the field name and data type check and it would also hint that the field's values were calculated in a batch process separately from editing the geometry.

Re: Calc Latitude, Longitude Field

NeilAyres — Tue, 09 Jun 2015 08:34:38 GMT

If the datum (underlying GCS of the coordinate system) of the input is not WGS84, then a transformation is required if you want WGS84 based coordinates to come out.

I have never actually used the spatial reference option when opening an Update cursor. Generally, the feature exists and it knows what coordinate system it is in.... The help files here are unclear as to how this option affects the outcome.

I have in the past used the .projectAs(SR2) method on the geometry object to re-project geometries on the fly inside python. Specifying a transform is a little trickier and you have to set this via environment setting.

Re: Calc Latitude, Longitude Field

XanderBakker — Tue, 09 Jun 2015 12:37:35 GMT

Neil Ayres is right that you can project a geometry using you code and apply the transformation to it. However, it is not necessary to do this in the environment settings as you can use the the second (optional) parameter of the arcpy.Geometry::projectAs (spatial_reference, {transformation_name}) function.

An example of requiring a transformation is projecting data in Colombia from GSC_Bogota to MAGNA_Colombia_Bogota. This would require multiple steps (projections) and when going from GCS_MAGNA to GCS_Bogota, the transformation Bogota_To_MAGNA_Region_5_MB should be applied (according to specifications here).

If the transformation is not specified another transformation could be applied resulting in an inaccurate on-the-fly projection of the geometry and yielding inaccurate output coordinates.

In case these types of problems do not affect you, there is no problem. However, specifying a different spatial reference (sr) then your input as parameter of a cursor may have unexpected side-effects.

About the field names, you could use a very specific name that a normal user would unlikely create him/herself or create a tool and let the user define the names of the output fields. In case they exist verify it during validation or simply return an error.

Re: Calc Latitude, Longitude Field

BlakeTerhune — Tue, 09 Jun 2015 14:11:28 GMT

Is there an easier way to just convert state plane coordinates to decimal degrees without reprojecting? Our data is in NAD_1983_HARN_StatePlane_Arizona_Central_FIPS_0202.

EDIT:

Looks like the accepted transformation from NAD_1983 to WGS_1984 in the contiguous 48 states of the US is NAD_1983_To_WGS_1984_5. However, I don't see how to specify a transformation when using AddGeometryAttributes_management().

Re: Calc Latitude, Longitude Field

XanderBakker — Tue, 09 Jun 2015 17:02:25 GMT

If the tool itself does not provide it, you can use the environment settings:

Re: Calc Latitude, Longitude Field

BlakeTerhune — Thu, 11 Jun 2015 20:32:30 GMT

I did a test of AddGeometryAttributes_management() for "POINT_X_Y_Z_M" and it can take 15 minutes to do an SDE feature class with about 80,000 points. I'm not even trying anything fancy with reprojecting, I'm just using the feature class's assigned coordinate system. The output is just Point_X and Point_Y fields. Is this typical performance?

Re: Calc Latitude, Longitude Field

XanderBakker — Thu, 11 Jun 2015 21:10:00 GMT

Sounds slow... Will do a test on a local dataset that I have and see how it compares to using the da update cursor.

Re: Calc Latitude, Longitude Field

BlakeTerhune — Thu, 11 Jun 2015 21:37:56 GMT

I just tried to copy the data to in_memory workspace and do the AddGeometry there thinking it would be faster, but it looks like AddGeomerty won't work with in_memory data. I get an

ERROR 000732: Input Features: Dataset does not exist or is not supported

I don't see anything explicit in the documentation about this (unlike with Project, which clearly says in_memory is not a valid workspace).

Re: Calc Latitude, Longitude Field

DanPatterson_Retired — Thu, 11 Jun 2015 21:54:10 GMT

Does this tool not modify a file in-place? That could be the issue. in memory appears to be used for outputs unless the inputs are there as there as well

Re: Calc Latitude, Longitude Field

DanPatterson_Retired — Thu, 11 Jun 2015 21:57:16 GMT

see my comment above

Re: Calc Latitude, Longitude Field

BlakeTerhune — Sat, 11 Dec 2021 20:57:30 GMT

I should clarify, My plan was to create a variable that has a path "in_memory"

MyInputFC = r"C:\MySDEConn.sde\MyInputFC"
inputFC_inmem = os.path.join("in_memory", "MyInputFC_inmem")

Then I copy rows from the source to the in_memory workspace

arcpy.CopyRows_management(MyInputFC, inputFC_inmem)

Then I try to run AddGeometry on the data in memory. This is where it errors.

arcpy.AddGeometryAttributes_management(
    inputFC_inmem,  ## Input_Features
    "POINT_X_Y_Z_M",  ## Geometry_Properties
    "",  ## Length_Unit
    "",  ## Area_Unit
    GCS_WGS_1984  ## Coordinate_System
)

Then I would copy the rows back out of in_memory workspace.

Re: Calc Latitude, Longitude Field

DanPatterson_Retired — Thu, 11 Jun 2015 23:34:46 GMT

your syntax seems correct

>>> import os

>>> out = os.path.join("in_memory","myfile")

>>> out

'in_memory/myfile'

so I don't know why it can't find it

Re: Calc Latitude, Longitude Field

BlakeTerhune — Fri, 12 Jun 2015 00:09:45 GMT

It takes a good 20 seconds for the copy rows to complete, so I'm sure it's actually creating the in_memory data. My guess is that AddGeometry just can't work its magic in an in_memory workspace.

Re: Calc Latitude, Longitude Field

BrentHoskisson — Tue, 16 Jun 2015 14:44:11 GMT

Here is what I started with

PC Software Download - SPCS83

The source code is in FORTRAN, but It is easily decipherable. It took me a few days to translate it into Javascript and it works very fast (I put it in my mousemove event). No ESRI or GIS - It stays in the program.

http://www.ngs.noaa.gov/PC_PROD/pc_prod.shtml

Has a bunch of these conversion programs available for public consumption.

Brent Hoskisson

Re: Calc Latitude, Longitude Field

BlakeTerhune — Wed, 17 Jun 2015 22:05:11 GMT

I just realized I was using copy rows instead of copy features so it was trying to add geometry attributes on a table instead of a feature class. DOH!

Using in_memory workspace to calculate X/Y is about twice as fast (compared to in a file geodatabase) on a point feature class with just over 80,000 rows.