def lat_long_rasters(ned_raster, lat_raster, lon_raster): import arcgisscripting gp = arcgisscripting.create(9.3) gp.CheckOutExtension("Spatial") gp.WorkSpace = os.path.dirname(ned_raster) gp.SingleOutputMapAlgebra_sa( "Con(IsNull({0}), {1}, $$YMap)".format(ned_raster, ned_raster), lat_raster) gp.SingleOutputMapAlgebra_sa( "Con(IsNull({0}), {1}, $$XMap)".format(ned_raster, ned_raster), lon_raster) del gp
Wednesday, August 31, 2011
Latitude & Longitude rasters in ArcGIS 10
The ability to calculate latitude and longitude rasters using the GRID commands $$YMap and $$XMap was removed in ArcGIS 10 (as well as all the GRID commands including $$COLMAP, $$ROWMAP, $$NROWS, $$NCOLS, $$CELLSIZE). To get around this you can still load the ArcGIS 9.3 geoprocessor in ArcGIS 10 and then call the functions using the Single Output Map Algebra tool. Latitude & longitude rasters can also be calculated with NumPy or the FlowAccumulation function but these aren't great options for really large rasters. Make sure to use a DEM that is in a geographic coordinates system and then project the results.
Labels:
arcgisscripting,
latitude,
longitude
Friday, August 19, 2011
Nan in NumPy array to NoData in ArcGIS raster
For smaller rasters, NumPy is a much faster alternative to some of the ArcGIS geoprocessing functions. One issue with doing this is that every raster must have exactly the same dimensions (rows, columns, cellsize, extent, etc.) since NumPy has no spatial component. Later on I will post about some of the techniques we use to make sure they are all identical.
There is one issue with converting NumPy arrays back to rasters though that I want to outline. The arcpy function NumPyArrayToRaster does not correctly convert NaN values. The values will display as NoData in the identity window, but they do not act as NoData in subsequent calculations and the 'Display NoData as' option in the Symbology tab has no impact.
To get around this, convert all NaN values to some unused value, and then set that value as the value_to_nodata parameter. For example:
This is Bug NIM063028 and is not fixed in any current or future version.
There is one issue with converting NumPy arrays back to rasters though that I want to outline. The arcpy function NumPyArrayToRaster does not correctly convert NaN values. The values will display as NoData in the identity window, but they do not act as NoData in subsequent calculations and the 'Display NoData as' option in the Symbology tab has no impact.
To get around this, convert all NaN values to some unused value, and then set that value as the value_to_nodata parameter. For example:
test_array[np.isnan(test_array)] = -9999
test_obj = arcpy.NumPyArrayToRaster(test_array, pnt, cellsize, cellsize, -9999)
test_obj.save(test_raster)
This is Bug NIM063028 and is not fixed in any current or future version.
Labels:
ArcGIS 10,
NumPy,
NumPyArrayToRaster
Pyramids documentation
There is a small issue with the documentation for the pyramids environment parameter. The explanation section for the pyramid_option has a few errors. The option is "PYRAMIDS", not "PYRAMID", and the default is to build full pyramids (-1 level), not NONE. This is Bug NIM062605 and is reported as fixed in 10.1. The example code at the bottom of the help page shows the correct usage. To force no pyramids I prefer to use: env.pyramid = "PYRAMIDS 0".
Also, the CreateRandomRaster function doesn't honor the pyramids setting. This is Bug NIM062606 and is not currently reported as fixed.
Also, the CreateRandomRaster function doesn't honor the pyramids setting. This is Bug NIM062606 and is not currently reported as fixed.
Labels:
ArcGIS 10,
CreateRandomRaster,
environment parameter,
pyramids
Wednesday, August 17, 2011
Offset of TIF rasters in ArcGIS 10, pre-SP2
In Arc10, pre-SP2, our Landsat TIF rasters were being shifted half a pixel away from where they should be. For older Landsat 5TM imagery the thermal band could be up 60m offset in both the x and y. ESRI fixed the problem in SP2 (Bug NIM060759) but in case you are loading TIF rasters with a pre-SP2 version of Arc10, it is easy to shift the rasters back to the correct location. Since we were not using the panchromatic rasters from Landsat 7, all of our half cell sizes were integers (30m -> 15m, 60m -> 30m, 120m -> 60m), but that may be an issue for other applications.
It is also fairly simple to check which ArcGIS version and service pack you have:
half_cell_width = int(0.5*arcpy.Describe(landsat_tif).meanCellWidth+0.5)
half_cell_height = int(0.5*arcpy.Describe(landsat_tif).meanCellHeight+0.5)
arcpy.Shift_management(landsat_tif, shift_raster, -half_cell_width, half_cell_height)
It is also fairly simple to check which ArcGIS version and service pack you have:
version_str = arcpy.GetInstallInfo("desktop")["SPBuild"]
version_str = "10.0.2.3200"
sp_str = version_str.split(".")[2]
sp_str = "2"
MakeFeatureLayer lock files
When calling the MakeFeatureLayer function in a python script, a lock file is created that may persist even after the script finishes. To remove the lock file, just call arcpy.Delete_management(layer) before the script terminates.
This forum post provided the answer to this problem.
This forum post provided the answer to this problem.
Labels:
ArcGIS 10,
lock file,
MakeFeatureLayer
Subscribe to:
Posts (Atom)