Select to view content in your preferred language

Trying to project a layer in a Python script

7222
10
02-10-2011 06:14 AM
ChristinaHerrick1
New Contributor
I am trying to write a python script to do a batch extraction of pixel values.  I have a feature class in a file geodatabase that has about 6500 records in it, and I also have a list of rasters.  I am trying to write code that will select one point at a time based on an attribute called 'cellnum', scroll through the list of rasters and see if the point overlaps with any of the rasters, and if it does, it will extract the pixel value.  However, the points are in WGS84 and my rasters are all in different UTM zones.  I need to reproject the selected point to match the raster in question (the script doesn't do on-the-fly projections).

However, when I try to reproject the selected point, I get errors.

import arcpy
from arcpy import env

env.overwriteOutput = True
env.workspace = "C:/Data/Project.gdb"
env.scratchWorkspace = env.workspace
sort = "points"
inlyr = "in_memory\\sort"
newlyr = "in_memory\\sort2"
hypdir = arcpy.ListRasters("")

arcpy.MakeFeatureLayer_management(sort, inlyr)

rows = arcpy.SearchCursor(inlyr)
for row in rows:
    center = row.getValue("cellnum")
    Expression = "\"cellnum\" = " + str(center)
    
    if center == 0:
        arcpy.SelectLayerByAttribute_management(inlyr, "NEW_SELECTION", Expression)
        
        point = arcpy.SearchCursor(inlyr)
        for i in point:
            for scene in hypdir:
                des = arcpy.Describe(scene)
                descr = arcpy.Describe(inlyr)
                sr1 = des.spatialReference
                sr2 = descr.spatialReference
                            
                if sr1 != sr2:
                    print "Matching projections..."
                    arcpy.Project_management(inlyr, newlyr, str(sr1.exportToString()))


I get the error message:

Failed to execute. Parameters are not valid.
ERROR 000944: Output feature class cannot be in the in_memory workspace.
Failed to execute (Project).

I understand that.  So I change 'newlyr' to be "C:/Data/temp_sort.lyr", and I get this error:

ERROR 999999: Error executing function.
Create output feature class failed
Failed to execute (Project).

I've tried this also, newlyr = "C:/Data/Project.gdb/temp_sort", but still get an error:

Runtime error <type 'exceptions.IOError'>: C:/Data/Project.gdb/temp_sort does not exist

I'd rather not reproject the entire point file each time if I don't have to.  This script will be run in the future with other datasets that could be much larger. 

I'm at my wit's end.  Can anyone shed some light?

Thanks a million in advance, this has been giving me headaches for over a week now.
0 Kudos
10 Replies
JasonScheirer
Regular Contributor II
0 Kudos
KevinGooss
Occasional Contributor
In my experience this just does not work.
And it doesn't appear that it was intended to work that way.
in the documentation for the SearchCursor the spatialReference parameter is defined as:

The Spatial Reference of the feature class


if that param is designed to reproject the fc on the cursor then the doc is poor.

Additionally, in the topic named "Setting a cursor's spatial reference" there is the following:

By default, the spatial reference of the geometry returned from a search cursor is the same as the feature class opened by the cursor. You can also set the spatial reference on an update or insert cursor.


which doesn't sound promising for SearchCursor.
Based on some tests i did the spatial reference is not considered when creating the SearchCursor.
I wish that it were jscheirer.
I'm trying to get UTM15N measures for features that are in Web Mercator and that would be perfect for me.
0 Kudos
JasonScheirer
Regular Contributor II
Please file a bug.
0 Kudos
KG
by
New Contributor
I am experiencing this behavior as well. I have noticed that when you drop the sort parameter the issue is not present.

I thought the forums were for the exchange of ideas to arrive at solutions?
having someone from esri reply and say 'file a bug' is not very helpful nor is it in the spirit of what these forums are used for.

If KG22 files a bug that doesn't help cfw2, myself or any other esri software or forum user.
What we need is a workaround for the issue and that typically comes from a lively discussion of the issue on these forums.

It's disconcerting to have an esri representative chime in with a solution that doesn't work. Then not take the time to offer a helpful tip or suggestion to arrive at an answer.

Am I alone on this one?
0 Kudos
JasonScheirer
Regular Contributor II
Please file an issue with Esri tech support so it can be investigated, tested and fixed for a future SP, Quick Fix, or release.
0 Kudos
ChristinaHerrick1
New Contributor
Using the spatial reference information in the search cursor seems to have worked for me.  I'm using a point feature class within a file geodatabase that is in WGS84, and I'm checking each point for overlap with a list of rasters.  I have something that looks like this:

import arcpy
from arcpy import env
from arcpy import management as dm

inspace = "C:/test"
env.workspace = inspace

sort = "pointfile"
hypdir = arcpy.ListRasters("","TIF")

dm.MakeFeatureLayer(sort,inlyr)
for scene in hypdir:
   des = arcpy.Describe(scene)
   sr1 = des.spatialReference

   rows = arcpy.SearchCursor(inlyr, "", sr1)
   for row in rows:
        ID = row.getValue("ID")
        Expression = "\"ID \" = " + str(ID)

        dm.SelectLayerByAttribute(inlyr, "NEW_SELECTION", Expression


... and so on.  Since my rasters are all in different UTM zones, the extracted values would all be 'Null' without the 'sr1' parameter in the search cursor (that's what prompted me to originally post).   🙂
0 Kudos
KimOllivier
Regular Contributor II
Web Mercator? Isn't this just a sop for a simple sphere that G**gle insists that must be used to map the world? Its not a proper projection because it has no datum. I think that is why there is no projection using a cursor option. Try pretending it is WGS84 and suck hard on the registration errors, which is evident on every road network overlay of the images in GE or GM.

I use the projection option on cursors and it works for me. I can even get the transform between two datums to work as well if I match the projection names exactly to suit the transform.

I could not find an easy way to project a shape object. I had to write out a temporary featureclass to the scratch workspace, not in_memory.
(I got the shape object from an interactive pick using a tool, and then I needed the pick coordinate in lat/long to query a server.)
0 Kudos
StevePeaslee
Occasional Contributor
My experience has been that when a feature layer has a selection set, then the coordinates returned by a cursor fail to project. No error, they just aren't projected.

-Steve
0 Kudos
MattDeitemeyer
New Contributor III
This seems like related problem, I'm updating an arcgisscripting 9.3 python script to arcpy.  The In_Memory layer errors out on the Project geoprocess.

ERROR 000944: Output feature class cannot be in the in_memory workspace.

Here's the test code, minus the bogus projection info string.

import sys, string, os, shutil, arcpy

# Create the Geoprocessor object and set variables
arcpy.env.overwriteOutput = 1

Ticket_Envelope = r"C:\Temp\Matt.shp"
Ticket_Layer = r"In_Memory\Ticket_Layer"
arcpy.Project_management(Ticket_Envelope, Ticket_Layer, "blah blah blah", "", "")

So "In_Memory" feature classes seem to work like normal, with other geoprocesses in ArcPy.  However, Projection appears to be an exception And\Or a bug.  Any insight?
0 Kudos