Turn on suggestions
Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.
Cancel
- Home
- :
- All Communities
- :
- Products
- :
- ArcGIS Spatial Analyst
- :
- ArcGIS Spatial Analyst Questions
- :
- Adding total visible pixels for each observer loca...
Options
- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Mute
- Printer Friendly Page
Adding total visible pixels for each observer location from created viewsheds
Subscribe
2093
6
11-25-2013 07:03 AM
11-25-2013
07:03 AM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi,
I have used iteration in ModelBuilder to create multiple viewsheds for each observer location in a point feature class. However, I now need to transfer the value for visible pixels back to the observer feature class, so that I can then compare the total visible areas of different observer points. I can't join the tables as there is no common field in the raster and points files.
Does anyone know a way of automating this without manually having to insert the values into a new field in the points attributes table?
Thanks!
I have used iteration in ModelBuilder to create multiple viewsheds for each observer location in a point feature class. However, I now need to transfer the value for visible pixels back to the observer feature class, so that I can then compare the total visible areas of different observer points. I can't join the tables as there is no common field in the raster and points files.
Does anyone know a way of automating this without manually having to insert the values into a new field in the points attributes table?
Thanks!
6 Replies
11-25-2013
10:21 PM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi Anne,
I think that you can use python for this. You would probably have to convert the model to python and clean it up a bit.
I suppose you loop through the observer points and for each point you are calculating a viewshed raster. This raster contains 1 for visible and 0 for not visible. Using numpy it is easy to get the sum of the raster. In case you have NoData cells, the numpy array needs to be correct for this. The code is very simple:
In the output field of the point featureclass (observers), you can use a simple updatecursor (and a whereclause) to set the calculated value.
The script could look something like this:
If you have any questions, please let me know.
Kind regards,
Xander
I think that you can use python for this. You would probably have to convert the model to python and clean it up a bit.
I suppose you loop through the observer points and for each point you are calculating a viewshed raster. This raster contains 1 for visible and 0 for not visible. Using numpy it is easy to get the sum of the raster. In case you have NoData cells, the numpy array needs to be correct for this. The code is very simple:
import numpy, arcpy
from numpy import *
myArray = arcpy.RasterToNumPyArray("Viewshed1")
# avoid erroneous sum due to NoData cells
myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
visibleCells = myArray2.sum()In the output field of the point featureclass (observers), you can use a simple updatecursor (and a whereclause) to set the calculated value.
The script could look something like this:
import arcpy, os, numpy
from numpy import *
# settings of your data
fcpoints = r'C:\Project\_Forums\Viewshed\test.gdb\sources'
fldVisible = 'VisiblePix' # will contain the number of visible pixels
dem = r'C:\Project\_Forums\cost\fgdb\test.gdb\elevutm16n'
# ouput workspace and name for viewshed
outputWS = r'C:\Project\_Forums\Viewshed\test.gdb'
prefixVS = 'Viewshed'
arcpy.env.workspace = outputWS
arcpy.env.overwriteOutput = True
# grab a SA license
if arcpy.CheckExtension("Spatial") == "Available":
arcpy.CheckOutExtension("Spatial")
else:
exit("no SA license available")
# loop through points
fields = ("OID@",fldVisible)
with arcpy.da.UpdateCursor(fcpoints, fields) as cursor:
for row in cursor:
oid = row[0]
# define output name (of oid = 1, raster will be called Viewshed1
outVSname = outputWS + os.sep + prefixVS + str(oid)
# create featurelayer for current observer point
whereclause = "OBJECTID = {0}".format(oid)
FLsources = "sources_Layer"
arcpy.MakeFeatureLayer_management(fcpoints,FLsources,whereclause,outputWS)
# calculate viewshed
arcpy.gp.Viewshed_sa(dem,FLsources,outVSname,"1","FLAT_EARTH","0,13")
# determine number of visible cells
myArray = arcpy.RasterToNumPyArray(outVSname)
myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
visPixels = myArray2.sum()
print "There are {0} visible pixels for point with OID={1}".format(visPixels,oid)
# update value in point featureclass
row[1] = visPixels
cursor.updateRow(row)
del myArray, myArray2
del row
# release license
arcpy.CheckInExtension("Spatial")If you have any questions, please let me know.
Kind regards,
Xander
11-29-2013
06:35 AM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Thanks for your reply. I am rather new to python but have tried to apply your script to my data. However, when I run it I get the error messages:
ERROR 999999: Error executing function.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
ERROR 010004: All cells in Raster G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\combraster2 have the NODATA value. Stop execution.
No valid Observation points present.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed0.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed0.
Failed to execute (Viewshed).
Don't understand why it says there is no data in the raster and no valid obs points present. Any ideas?
My script is:
>>> import arcpy, os, numpy
... from numpy import*
...
... # settings of your data
... fcpoints = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\clip16points.shp'
... fldVisible = 'VisiblePix' # will contain the number of visible pixels
... dem = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\combraster2'
...
... # ouput workspace and name for viewshed
... outputWS = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb'
... prefixVS = 'Viewshed'
... arcpy.env.workspace = outputWS
... arcpy.env.overwriteOutput = True
...
... # grab a SA license
... if arcpy.CheckExtension("Spatial") == "Available":
... arcpy.CheckOutExtension("Spatial")
... else:
... exit("no SA license available")
...
... # loop through points
... fields = ("oid@",fldVisible)
... with arcpy.da.UpdateCursor(fcpoints, fields) as cursor:
... for row in cursor:
... oid = row[0]
...
... # define output name (of oid = 1, raster will be called Viewshed1
... outVSname = outputWS + os.sep + prefixVS + str(oid)
...
... # create featurelayer for current observer point
... whereclause = "OBJECTID = {0}".format(oid)
... FLsources = "sources_Layer"
... arcpy.MakeFeatureLayer_management(fcpoints,FLsources,whereclause,outputWS)
...
... # calculate viewshed
... arcpy.gp.Viewshed_sa(dem,FLsources,outVSname,"1","FLAT_EARTH","0.13")
...
... # determine number of visible cells
... myArray = arcpy.RasterToNumPyArray(outVSname)
... myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
... VisiblePix = myArray2.sum()
... print "There are {0} visible pixels for point with oid={1}".format(VisiblePix,oid)
...
... # update value in point featureclass
... row[1] = VisiblePix
... cursor.updateRow(row)
... del myArray, myArray2
...
... del row
...
... # release license
... arcpy.CheckInExtension("Spatial")
ERROR 999999: Error executing function.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
ERROR 010004: All cells in Raster G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\combraster2 have the NODATA value. Stop execution.
No valid Observation points present.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed0.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed0.
Failed to execute (Viewshed).
Don't understand why it says there is no data in the raster and no valid obs points present. Any ideas?
My script is:
>>> import arcpy, os, numpy
... from numpy import*
...
... # settings of your data
... fcpoints = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\clip16points.shp'
... fldVisible = 'VisiblePix' # will contain the number of visible pixels
... dem = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\combraster2'
...
... # ouput workspace and name for viewshed
... outputWS = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb'
... prefixVS = 'Viewshed'
... arcpy.env.workspace = outputWS
... arcpy.env.overwriteOutput = True
...
... # grab a SA license
... if arcpy.CheckExtension("Spatial") == "Available":
... arcpy.CheckOutExtension("Spatial")
... else:
... exit("no SA license available")
...
... # loop through points
... fields = ("oid@",fldVisible)
... with arcpy.da.UpdateCursor(fcpoints, fields) as cursor:
... for row in cursor:
... oid = row[0]
...
... # define output name (of oid = 1, raster will be called Viewshed1
... outVSname = outputWS + os.sep + prefixVS + str(oid)
...
... # create featurelayer for current observer point
... whereclause = "OBJECTID = {0}".format(oid)
... FLsources = "sources_Layer"
... arcpy.MakeFeatureLayer_management(fcpoints,FLsources,whereclause,outputWS)
...
... # calculate viewshed
... arcpy.gp.Viewshed_sa(dem,FLsources,outVSname,"1","FLAT_EARTH","0.13")
...
... # determine number of visible cells
... myArray = arcpy.RasterToNumPyArray(outVSname)
... myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
... VisiblePix = myArray2.sum()
... print "There are {0} visible pixels for point with oid={1}".format(VisiblePix,oid)
...
... # update value in point featureclass
... row[1] = VisiblePix
... cursor.updateRow(row)
... del myArray, myArray2
...
... del row
...
... # release license
... arcpy.CheckInExtension("Spatial")
12-01-2013
11:35 AM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi Anne,
The description of this error is:
All of the cells in the grid are NoData. The tool requires that some cells have values other than NoData to operate on. This may happen if there is no overlap between the extent and the input(s).
Solution
Ensure that the input grid has some cell values that are not NoData. If you have set the Extent environment, make sure it falls at least partway within your input data.
Can you verify that this is not the case?
Kind regards,
Xander
The description of this error is:
All of the cells in the grid are NoData. The tool requires that some cells have values other than NoData to operate on. This may happen if there is no overlap between the extent and the input(s).
Solution
Ensure that the input grid has some cell values that are not NoData. If you have set the Extent environment, make sure it falls at least partway within your input data.
Can you verify that this is not the case?
Kind regards,
Xander
12-01-2013
09:55 PM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi Anne,
I notice a few more things. I wrote the code based on rasters and points FC in a file geodatabase. You are using grids and a point shapefile. One of the errors is "No valid Observation points present.". This is likely the cause of the wrong whereclause that is being build. The whereclause is using the field OBJECTID (standard in fgdb fc's) and the field is not delimited.
To make this a little more gentle the code should be changed a bit:
1) use the describe methode to obtain the oid field (if it has one, if not exit the procedure):
2) the whereclause should be build differently:
This will result probably in a whereclause like: "FID" = 0
So the final code in your case would be something like:
Before you run it please check you input dem and environment settings as explained in my previous post.
It this still doesn't work, try to run the code with a local copy of your dem. The 'G' drive most likely maps to a network folder and I have seen issues with running scripts with grids on network locations that failed, when running it with a local copy actually worked.
Kind regards,
Xander
I notice a few more things. I wrote the code based on rasters and points FC in a file geodatabase. You are using grids and a point shapefile. One of the errors is "No valid Observation points present.". This is likely the cause of the wrong whereclause that is being build. The whereclause is using the field OBJECTID (standard in fgdb fc's) and the field is not delimited.
To make this a little more gentle the code should be changed a bit:
1) use the describe methode to obtain the oid field (if it has one, if not exit the procedure):
desc = arcpy.Describe(fcpoints)
if desc.hasOID:
fldOID = desc.OIDFieldName
else:
exit("featureclass does not have an Object ID field")2) the whereclause should be build differently:
fld = arcpy.AddFieldDelimiters(fcpoints,fldOID)
whereclause = "{1} = {0}".format(oid,fld)This will result probably in a whereclause like: "FID" = 0
So the final code in your case would be something like:
import arcpy, os, numpy
from numpy import *
# settings of your data
fcpoints = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\clip16points.shp'
fldVisible = 'VisiblePix' # will contain the number of visible pixels
dem = r'G:\xxx\xxxx\xxxx\xxxx\xxx\xxxx\xxx\combraster2'
# ouput workspace and name for viewshed
outputWS = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb'
prefixVS = 'Viewshed'
arcpy.env.workspace = outputWS
arcpy.env.overwriteOutput = True
# grab a SA license
if arcpy.CheckExtension("Spatial") == "Available":
arcpy.CheckOutExtension("Spatial")
else:
exit("no SA license available")
# get the oid field
desc = arcpy.Describe(fcpoints)
if desc.hasOID:
fldOID = desc.OIDFieldName
else:
exit("featureclass does not have an Object ID field")
# loop through points
fields = ("oid@",fldVisible)
with arcpy.da.UpdateCursor(fcpoints, fields) as cursor:
for row in cursor:
oid = row[0]
# define output name (of oid = 1, raster will be called Viewshed1
outVSname = outputWS + os.sep + prefixVS + str(oid)
# create a whereclause
fld = arcpy.AddFieldDelimiters(fcpoints,fldOID)
whereclause = "{1} = {0}".format(oid,fld)
# create featurelayer for current observer point
FLsources = "sources_Layer"
arcpy.MakeFeatureLayer_management(fcpoints,FLsources,whereclause,outputWS)
# calculate viewshed
arcpy.gp.Viewshed_sa(dem,FLsources,outVSname,"1","FLAT_EARTH","0.13")
# determine number of visible cells
myArray = arcpy.RasterToNumPyArray(outVSname)
myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
VisiblePix = myArray2.sum()
print "There are {0} visible pixels for point with oid={1}".format(VisiblePix,oid)
# update value in point featureclass
row[1] = VisiblePix
cursor.updateRow(row)
del myArray, myArray2
del row
# release license
arcpy.CheckInExtension("Spatial")Before you run it please check you input dem and environment settings as explained in my previous post.
It this still doesn't work, try to run the code with a local copy of your dem. The 'G' drive most likely maps to a network folder and I have seen issues with running scripts with grids on network locations that failed, when running it with a local copy actually worked.
Kind regards,
Xander
12-03-2013
02:53 AM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Hi Xander,
Thanks for your help so far and sorry for my ignorance. I tried converting my points to feature classes and used a fgdb raster, then tried your initial code again, but got the same errors. The points overlap with the raster and are in the same coordinate system, and I set Union of Inputs as the processing extent (though do I actually need to do this in Python?). The vast majority of cells have values. I have attached the TIFF file from which the fgdb raster was made.
The processing output is also attached.
[HTML]Executing: MakeFeatureLayer C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clip16points sources_Layer "fldOID = 1" C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb "OBJECTID_1 OBJECTID VISIBLE NONE;Shape Shape VISIBLE NONE;Join_Count Join_Count VISIBLE NONE;TARGET_FID TARGET_FID VISIBLE NONE;FID_ FID_ VISIBLE NONE;Date Date VISIBLE NONE;Date_2 Date_2 VISIBLE NONE;Year Year VISIBLE NONE;Field_comp Field_comp VISIBLE NONE;Unique_ID Unique_ID VISIBLE NONE;X X VISIBLE NONE;Y Y VISIBLE NONE;Waypoint_I Waypoint_I VISIBLE NONE;elevation elevation VISIBLE NONE;Observers Observers VISIBLE NONE;Photo_ID Photo_ID VISIBLE NONE;Sward_heig Sward_heig VISIBLE NONE;Sward_he_1 Sward_he_1 VISIBLE NONE;Sward_he_2 Sward_he_2 VISIBLE NONE;Crop_type Crop_type VISIBLE NONE;General_co General_co VISIBLE NONE;Notes Notes VISIBLE NONE;Turbine_me Turbine_me VISIBLE NONE;Anser_drop Anser_drop VISIBLE NONE;RbG_droppi RbG_droppi VISIBLE NONE;Unknown_dr Unknown_dr VISIBLE NONE;Total_drop Total_drop VISIBLE NONE;Height_at_ Height_at_ VISIBLE NONE;Height_at1 Height_at1 VISIBLE NONE;Height_a_1 Height_a_1 VISIBLE NONE;Height_a_2 Height_a_2 VISIBLE NONE;Height_a_3 Height_a_3 VISIBLE NONE;Height_a_4 Height_a_4 VISIBLE NONE;Height_a_5 Height_a_5 VISIBLE NONE;Height_a_6 Height_a_6 VISIBLE NONE;Height_a_7 Height_a_7 VISIBLE NONE;Height_a_8 Height_a_8 VISIBLE NONE;Average_he Average_he VISIBLE NONE;Cat_No__Pl Cat_No__Pl VISIBLE NONE;Cat_No___1 Cat_No___1 VISIBLE NONE;Cat_No___2 Cat_No___2 VISIBLE NONE;Cat_No___3 Cat_No___3 VISIBLE NONE;Cat_No___4 Cat_No___4 VISIBLE NONE;Cat_No___5 Cat_No___5 VISIBLE NONE;Cat_No___6 Cat_No___6 VISIBLE NONE;Cat_No___7 Cat_No___7 VISIBLE NONE;Cat_No___8 Cat_No___8 VISIBLE NONE;Cat_No___9 Cat_No___9 VISIBLE NONE;No__densit No__densit VISIBLE NONE;No__dens_1 No__dens_1 VISIBLE NONE;No__dens_2 No__dens_2 VISIBLE NONE;No__dens_3 No__dens_3 VISIBLE NONE;No__dens_4 No__dens_4 VISIBLE NONE;Crop_colou Crop_colou VISIBLE NONE;Crop_col_1 Crop_col_1 VISIBLE NONE;Crop_col_2 Crop_col_2 VISIBLE NONE;Snow_Thin Snow_Thin VISIBLE NONE;Snow_Mediu Snow_Mediu VISIBLE NONE;Snow_Thick Snow_Thick VISIBLE NONE;_snow_cov _snow_cov VISIBLE NONE;Cat_dist_n Cat_dist_n VISIBLE NONE;Nearest_wa Nearest_wa VISIBLE NONE;Nearest__1 Nearest__1 VISIBLE NONE;Nearest__2 Nearest__2 VISIBLE NONE;Nearest__3 Nearest__3 VISIBLE NONE;Nearest__4 Nearest__4 VISIBLE NONE;No_visible No_visible VISIBLE NONE;Dist_trees Dist_trees VISIBLE NONE;Dist_turbi Dist_turbi VISIBLE NONE;Dist_roads Dist_roads VISIBLE NONE;WF1km WF1km VISIBLE NONE;WF500m WF500m VISIBLE NONE;WF250m WF250m VISIBLE NONE;WF2km WF2km VISIBLE NONE;Dist_track Dist_track VISIBLE NONE;Dist_MnRd Dist_MnRd VISIBLE NONE;Trees_1km Trees_1km VISIBLE NONE;km_roads km_roads VISIBLE NONE;km_trees km_trees VISIBLE NONE;m_track m_track VISIBLE NONE;m_roads m_roads VISIBLE NONE;km_roads_1 km_roads_1 VISIBLE NONE;m_track_1 m_track_1 VISIBLE NONE;m_roads_1 m_roads_1 VISIBLE NONE;m_trees m_trees VISIBLE NONE;m_trees_1 m_trees_1 VISIBLE NONE;km_trees_1 km_trees_1 VISIBLE NONE;km_tracks km_tracks VISIBLE NONE;km_tracks_ km_tracks_ VISIBLE NONE;Offset_A Offset_A VISIBLE NONE;VisiblePix VisiblePix VISIBLE NONE"
Start Time: Tue Dec 03 12:49:04 2013
Succeeded at Tue Dec 03 12:49:04 2013 (Elapsed Time: 0.00 seconds)
Executing: Viewshed C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clipdem30mWGS1 sources_Layer C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1 1 FLAT_EARTH 0.13
Start Time: Tue Dec 03 12:49:08 2013
ERROR 999999: Error executing function.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
ERROR 010004: All cells in Raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clipdem30mWGS1 have the NODATA value. Stop execution.
No valid Observation points present.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1.
Failed to execute (Viewshed).
Failed at Tue Dec 03 12:49:08 2013 (Elapsed Time: 0.00 seconds)[/HTML]
Any ideas?
Anne
Thanks for your help so far and sorry for my ignorance. I tried converting my points to feature classes and used a fgdb raster, then tried your initial code again, but got the same errors. The points overlap with the raster and are in the same coordinate system, and I set Union of Inputs as the processing extent (though do I actually need to do this in Python?). The vast majority of cells have values. I have attached the TIFF file from which the fgdb raster was made.
>>> import arcpy, os, numpy
... from numpy import *
...
... # settings of your data
... fcpoints = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clip16points'
... fldVisible = 'VisiblePix' # will contain the number of visible pixels
... dem = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clipdem30mWGS1'
...
... # ouput workspace and name for viewshed
... outputWS = r'C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb'
... prefixVS = 'Viewshed'
... arcpy.env.workspace = outputWS
... arcpy.env.overwriteOutput = True
...
... # grab a SA license
... if arcpy.CheckExtension("Spatial") == "Available":
... arcpy.CheckOutExtension("Spatial")
... else:
... exit("no SA license available")
... desc = arcpy.Describe(fcpoints)
... if desc.hasOID:
... fldOID = desc.OIDFieldName
... else:
... exit("featureclass does not have an Object ID field")
...
... # loop through points
... fields = ("OID@",fldVisible)
... with arcpy.da.UpdateCursor(fcpoints, fields) as cursor:
... for row in cursor:
... oid = row[0]
...
... # define output name (of oid = 1, raster will be called Viewshed1
... outVSname = outputWS + os.sep + prefixVS + str(oid)
...
... # create featurelayer for current observer point
... whereclause = "fldOID = {0}".format(oid)
... FLsources = "sources_Layer"
... arcpy.MakeFeatureLayer_management(fcpoints,FLsources,whereclause,outputWS)
...
... # calculate viewshed
... arcpy.gp.Viewshed_sa(dem,FLsources,outVSname,"1","FLAT_EARTH","0,13")
...
... # determine number of visible cells
... myArray = arcpy.RasterToNumPyArray(outVSname)
... myArray2 = where(myArray==1,ones(myArray.shape),zeros(myArray.shape))
... visPixels = myArray2.sum()
... print "There are {0} visible pixels for point with OID={1}".format(visPixels,oid)
...
... # update value in point featureclass
... row[1] = visPixels
... cursor.updateRow(row)
... del myArray, myArray2
...
... del row
...
... # release license
... arcpy.CheckInExtension("Spatial")The processing output is also attached.
[HTML]Executing: MakeFeatureLayer C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clip16points sources_Layer "fldOID = 1" C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb "OBJECTID_1 OBJECTID VISIBLE NONE;Shape Shape VISIBLE NONE;Join_Count Join_Count VISIBLE NONE;TARGET_FID TARGET_FID VISIBLE NONE;FID_ FID_ VISIBLE NONE;Date Date VISIBLE NONE;Date_2 Date_2 VISIBLE NONE;Year Year VISIBLE NONE;Field_comp Field_comp VISIBLE NONE;Unique_ID Unique_ID VISIBLE NONE;X X VISIBLE NONE;Y Y VISIBLE NONE;Waypoint_I Waypoint_I VISIBLE NONE;elevation elevation VISIBLE NONE;Observers Observers VISIBLE NONE;Photo_ID Photo_ID VISIBLE NONE;Sward_heig Sward_heig VISIBLE NONE;Sward_he_1 Sward_he_1 VISIBLE NONE;Sward_he_2 Sward_he_2 VISIBLE NONE;Crop_type Crop_type VISIBLE NONE;General_co General_co VISIBLE NONE;Notes Notes VISIBLE NONE;Turbine_me Turbine_me VISIBLE NONE;Anser_drop Anser_drop VISIBLE NONE;RbG_droppi RbG_droppi VISIBLE NONE;Unknown_dr Unknown_dr VISIBLE NONE;Total_drop Total_drop VISIBLE NONE;Height_at_ Height_at_ VISIBLE NONE;Height_at1 Height_at1 VISIBLE NONE;Height_a_1 Height_a_1 VISIBLE NONE;Height_a_2 Height_a_2 VISIBLE NONE;Height_a_3 Height_a_3 VISIBLE NONE;Height_a_4 Height_a_4 VISIBLE NONE;Height_a_5 Height_a_5 VISIBLE NONE;Height_a_6 Height_a_6 VISIBLE NONE;Height_a_7 Height_a_7 VISIBLE NONE;Height_a_8 Height_a_8 VISIBLE NONE;Average_he Average_he VISIBLE NONE;Cat_No__Pl Cat_No__Pl VISIBLE NONE;Cat_No___1 Cat_No___1 VISIBLE NONE;Cat_No___2 Cat_No___2 VISIBLE NONE;Cat_No___3 Cat_No___3 VISIBLE NONE;Cat_No___4 Cat_No___4 VISIBLE NONE;Cat_No___5 Cat_No___5 VISIBLE NONE;Cat_No___6 Cat_No___6 VISIBLE NONE;Cat_No___7 Cat_No___7 VISIBLE NONE;Cat_No___8 Cat_No___8 VISIBLE NONE;Cat_No___9 Cat_No___9 VISIBLE NONE;No__densit No__densit VISIBLE NONE;No__dens_1 No__dens_1 VISIBLE NONE;No__dens_2 No__dens_2 VISIBLE NONE;No__dens_3 No__dens_3 VISIBLE NONE;No__dens_4 No__dens_4 VISIBLE NONE;Crop_colou Crop_colou VISIBLE NONE;Crop_col_1 Crop_col_1 VISIBLE NONE;Crop_col_2 Crop_col_2 VISIBLE NONE;Snow_Thin Snow_Thin VISIBLE NONE;Snow_Mediu Snow_Mediu VISIBLE NONE;Snow_Thick Snow_Thick VISIBLE NONE;_snow_cov _snow_cov VISIBLE NONE;Cat_dist_n Cat_dist_n VISIBLE NONE;Nearest_wa Nearest_wa VISIBLE NONE;Nearest__1 Nearest__1 VISIBLE NONE;Nearest__2 Nearest__2 VISIBLE NONE;Nearest__3 Nearest__3 VISIBLE NONE;Nearest__4 Nearest__4 VISIBLE NONE;No_visible No_visible VISIBLE NONE;Dist_trees Dist_trees VISIBLE NONE;Dist_turbi Dist_turbi VISIBLE NONE;Dist_roads Dist_roads VISIBLE NONE;WF1km WF1km VISIBLE NONE;WF500m WF500m VISIBLE NONE;WF250m WF250m VISIBLE NONE;WF2km WF2km VISIBLE NONE;Dist_track Dist_track VISIBLE NONE;Dist_MnRd Dist_MnRd VISIBLE NONE;Trees_1km Trees_1km VISIBLE NONE;km_roads km_roads VISIBLE NONE;km_trees km_trees VISIBLE NONE;m_track m_track VISIBLE NONE;m_roads m_roads VISIBLE NONE;km_roads_1 km_roads_1 VISIBLE NONE;m_track_1 m_track_1 VISIBLE NONE;m_roads_1 m_roads_1 VISIBLE NONE;m_trees m_trees VISIBLE NONE;m_trees_1 m_trees_1 VISIBLE NONE;km_trees_1 km_trees_1 VISIBLE NONE;km_tracks km_tracks VISIBLE NONE;km_tracks_ km_tracks_ VISIBLE NONE;Offset_A Offset_A VISIBLE NONE;VisiblePix VisiblePix VISIBLE NONE"
Start Time: Tue Dec 03 12:49:04 2013
Succeeded at Tue Dec 03 12:49:04 2013 (Elapsed Time: 0.00 seconds)
Executing: Viewshed C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clipdem30mWGS1 sources_Layer C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1 1 FLAT_EARTH 0.13
Start Time: Tue Dec 03 12:49:08 2013
ERROR 999999: Error executing function.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
Cannot acquire a lock.
ERROR 010004: All cells in Raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\clipdem30mWGS1 have the NODATA value. Stop execution.
No valid Observation points present.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1.
ERROR 010160: Unable to open raster C:\Users\Anne.Harrison\Documents\ArcGIS\Default2.gdb\Viewshed1.
Failed to execute (Viewshed).
Failed at Tue Dec 03 12:49:08 2013 (Elapsed Time: 0.00 seconds)[/HTML]
Any ideas?
Anne
12-04-2013
03:28 AM
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Maybe this is caused by the lock on the data. What you could do is close the mxd that has references to the data and try again. I tried it with the small data sample you posted and some random points and it works:
[ATTACH=CONFIG]29562[/ATTACH]
Kind regards,
Xander
BTW; I run this as standalone script in PyScripter with ArcGIS 10.2
[ATTACH=CONFIG]29562[/ATTACH]
Kind regards,
Xander
BTW; I run this as standalone script in PyScripter with ArcGIS 10.2
