# Tag Info

10

You can access individual bands by joining the band name to the raster path - i.e. path/to/raster/band_name. Often using path/to/raster/Band_[band number] works (i.e os.path.join(rasterpath, 'Band_1'), but not always. I use Landsat 8 imagery quite a bit and ArcGIS names the bands 'CoastalAerosol', 'Blue', 'Green', 'Red', etc... If you don't know what the ...

9

The problem was that I had run out of space to write to disk, and the map algebra commands I was using were attempting to generate and write large temporary raster files.

8

I don't use arcpy but I know Python. This only means that you need to import the Python module with Con (and it is not arcpy alone) If I look at, Con (Spatial Analyst), I see that it is a function of the Spatial Analyst module and the examples provide the solution. You have forgotten to import the module: from arcpy.sa import * But the "pythonic" way to ...

8

What you want is a conditional calculation: return the value of r whenever r and r1 are equal and otherwise set the output to NA. The cell-by-cell arithmetic operations seem to be fastest. (They are much faster than, say, using mask or the reclassification functions.) Since they do not appear to offer an actual conditional operator, use two time-honored ...

6

My initial thought is that you must have assigned the variable pente_Rclass to a string representation of a raster instead of a raster object. This would cause your first error. See below: pente_Rclass = "myRasterName" represents a string... if you are using this in arcpy methods, it will automatically assume that this string is the name of a dataset in ...

6

You could sum the rasters and test whether the total is greater than 0. list_of_rasters = [Raster(a), Raster(b), Raster(c)] z = CellStatistics(list_of_rasters, "SUM") > 0

6

Additional arguments are passed as the 6th, and so on, arguments, after the custom function signature: ST_MapAlgebra(b1.rast, -- raster to operate on 1, -- band 'generate_random_raster(double precision[], integer[], text[])'::regprocedure, -- custom function signature pixeltype, -- can be null ...

6

This is how I would do it (it requires the Spatial Analyst extension which I think you have): import arcpy from arcpy.sa import * arcpy.CheckOutExtension = "Spatial" ndvi_raster = Raster("ndvivaluefrompixel") grey_raster = Raster("grayscalepixelvalue") output_raster = Con(ndvi_raster >= -0.000005 * grey_raster + 0.314367, 2, 1) output_raster.save("C:/...

5

My favourite way to deal with this would be to use a Con(IsNull) statement... Basically, using a Con statement allows you to define a condition; if it is true, do this, if it is false, do this - including leaving the values as they were. So in this case, I would be tempted to do something like Con(IsNull(myRaster), myRaster, 1) Which basically states: ...

5

I presume you're using ArcGIS Desktop which uses python 2x, not ArcGIS Pro which uses python 3x. Side note: this is one example of why it's important to specify the software you are using. In python 2x, when both operands are integers, division returns an integer result rounded towards minus infinity. This was changed in Python 3x so division returns a ...

4

If you try to output the iterate raster tool directly into raster calculator, you will only see the last iteration in the raster calculator layers and variables list. To get around this nuisance in model builder, use Collect Values to generate a list that you can pass off to Cell Statistics to do your calculations. Simply choose the "MEAN" overlay ...

4

Have you tried using the Image analysis window? I know it gives you less control, but still allows for some useful calculations. And it is computationally more efficient - so if your dealing with a large data-set it may speed things up for you. Little Tutorial I put together. http://www.thadwester.com/1/post/2012/05/ndvi-calculation-in-arcgis-10.html

4

arcpy.ListRasters() is returning a list of raster filenames not Raster objects. In your loop you're trying to multiply a string (the raster filename) by the float value which causes the TypeError exception. You need to multiply a Raster object by the float value. Convert your file name to a Raster object using arcpy.sa.Raster from arcpy.sa import * rasters ...

4

Arcpy use the Numpy array format (hidden to the users) as PyGGIS that uses the Python GDAL module. provider = raster.dataProvider() # path the original file filePath = str(provider.dataSourceUri()) # open the original file from osgeo import gdal raster_or = gdal.Open(filePath) # create a numpy array numpy_array = raster_or.ReadAsArray() # shape ...

4

Using only open source software, you will almost certainly need to do some programming yourself. GDAL is the de facto open source raster I/O library, so you will probably be using it or one of its many wrappers. You could use Python (e.g. rasterio + numpy/scipy) or node.js, e.g. node-gdal (though beware it is currently synchronous/blocking). As for actually ...

4

I think this is an issue with ArcGIS having open file handles on the files in your temp folder. As discovered in our comments, you could delete the Python variables and then use the Delete tool to remove the temporary directory (but not the del statement, that just deletes the variable), but you would still have to explicitly delete each result object... ...

4

You don't need to write code if you have GDAL installed (if you have QGIS installed that should be the case, but you can install GDAL seperately if you want/need to) This assumes you want a separate stddev for each image. If you want the stddev across ALL the images, you'll need to do a Raster Merge to combine them into a single raster. The gdalinfo tool ...

4

The \$\$RowMap, \$\$ColMap, \$\$XMap and \$\$YMap variables (and some others) are no longer (directly) supported as of ArcGIS 10.0. You can use it via python: arcpy.env.extent = arcpy.Extent(0, 0, 10, 10) #Change to suit # Or arcpy.env.extent = "path to raster" arcpy.env.cellSize = 1 #Change to suit arcpy.gp.SingleOutputMapAlgebra_sa("\$\$ROWMAP + 1"), "...

4

Are the red_band_raster and the nir_band_raster values integers? It seems likely that you should move the float() commands around each raster variable call in your equations, and not around the result of the addition/subtraction as you have it now.

4

You can do this just by using the Mosaic to New Raster geoprocessing tool. Just choose maximum for your Mosaic Operator and the new raster will use the high value anywhere that they overlap.

3

You will need to write your own aggregate function: CREATE OR REPLACE function sum_raster_state(raster, raster) returns raster language sql as \$f\$ SELECT CASE \$1 WHEN NULL THEN \$2 ELSE ST_MapAlgebra( \$1, \$2, '([rast1] + [rast2])', NULL, 'UNION', '[rast2]', '[rast1]', NULL) END;...

3

Use this r.mapcalc expression: r.mapcalc "output = not(if(Road))*Slope" because: if(Road)is equal to 1 if Road not zero, 0 otherwise. not(if(Road))is equal to 0 if Road not zero, 1 otherwise.

3

NoData is defined different ways depending on the raster source, including numerically. (You can look it up in raster properties if you are curious what the actual value is for your particular datasets.) I suggest instead using Reclassify instead of Raster Calculator, which can either explicitly ignore NoData values (leaving them NoData) or explicitly ...

3

One method which I used quite a while ago was to print the standard deviations in the Python Console in QGIS. I would first load all the raster layers (easier if they're kept in the same directory): import os, glob layers=[] raster_path = "C:\Users\You\Desktop\Raster_folder\\" # Path to directory containing rasters for raster in glob.glob(raster_path + "...

3

Not sure if this helps, but I have a script with many nested Con statements, and I find just using indentation helps. For example: EWR_2015 = Con( "ministerial_wetlands" > 0, " ministerial_wetlands ", Con( "non_ministerial_wetlands " > 0, "WTE_using_2001_min" - 0.25, Con( "nat_veg_raster" > 0 & "...

3

You're inputting your variables as string. Try without the quotes. outRas = Raster(in_raster1) * Raster(in_raster2) * Raster(in_raster3)

3

What @Joseph says. The Raster class needs a string containing the path to your raster. You've correctly assigned this to a variable, so you then need to pass those variables to instantiate your Rasters. This code should work: import arcpy from arcpy import env from arcpy.sa import * raster1 = ("C:/raster1") raster2 = ("C:/raster2") #Note lack of quotes in ...

3

Nearly there. Couple of small problems: You've used the bash syntax instead of cmd. i.e. \$ instead of %. Your --outfile parameter isn't specified explicitly, and looks a bit jumbled. Try: for %i in (*.tif) do C:\Python27\python.exe "C:\Program Files\QGIS Wien\bin\gdal_calc.py" -A %i --A_band=3 -B %i --B_band=4 --calc="((B-A)/(B+A))" --outfile=ndvi_%i If ...

3

Union (geoprocessing) feature make a "union" of 2 dataset (C = A united to B). Probably you are tring to JOIN two dataset. Try this: Go to JOIN settings: (right click -> properties -> Join) Add Vector Join of the second layer (click on "+" green-plus icon) and set these parameters: Join Layer: "the second layer to join" Join field: "ID" Join Target: "ID" ...

2

here is how I would do it in ArcGIS without coding: build the raster attribute table of your raster (you may need to multiply by some number (e.g.1000) and round it into integers) : this will give you the count of each pixel. in this table, multiply each value by its count and compute the cum sum, that you scale between 1 and 100. reclassify your raster ...

Only top voted, non community-wiki answers of a minimum length are eligible