0

I am using satellite data and by default it has a nodata value of 0 for each band.

This is fine except for when I'm making NDVIs where 0 is a real value we need to include and I don't want it to be classified as nodata.

I want to make sure all the bands/tifs I am working with have the same nodata value. So how do I convert a nodata value of 0 to something else, like -10000.

I was thinking of using arcpy's reclassify tool on my bands changing 0 to -10000 (before I make NDVIs) and then copying the raster with the copyraster tool using the setting for nodata as -10000.

This seems inefficient. Is there a tool that does this better? I'm open to open source libraries or arcpy

Emtomp
  • 543
  • 6
  • 15

3 Answers3

2

You can use gdal and numpy for this. Here is an example that will change all the 0's to -10000's for all the bands in a given set of TIFF files (tha are in the same directory):

import glob
import os

from osgeo import gdal
import numpy as np

os.chdir(r'C:\path\to\raster\files')           # change to directory with the tiff files

filenames = glob.glob('*.tif')

for fn in filenames:
    ds = gdal.Open(fn, 1)                      # pass 1 to modify the raster
    n = ds.RasterCount                         # get number of bands
    for i in range(1, n+1):
        band = ds.GetRasterBand(i)
        arr = band.ReadAsArray()               # read band as numpy array
        arr = np.where(arr == 0, -10000, arr)  # change 0 to -10000
        band.WriteArray(arr)                   # write the new array
        band.SetNoDataValue(-10000)            # set the NoData value
        band.FlushCache()                      # save changes
    del ds

Note that the data type of your TIFF files has to support negative values in order for this to work. Otherwise, the raster won't recognize the -10000 value.

Marcelo Villa
  • 5,928
  • 2
  • 19
  • 38
  • 1
    Should it be band = ds.GetRasterband(i) instead of band = ds.GetRasterband(n)? Because n will just be the last raster band in the stack? – Emtomp Feb 26 '20 at 19:45
  • 1
    Also I think it should be a capital 'B' in band? If not I get the AttributeError: type object 'object' has no attribute 'getattr'

    but this is a huge help!

     – Emtomp Feb 26 '20 at 20:10
  • Yes, you are right. I missed those details but already edited the answer. Thanks for pointing them out. – Marcelo Villa Feb 26 '20 at 22:49
  • So this worked amazingly on my first dataset, but I'm trying to run it on a new set of tifs that are also nodata = 0 and uint16, but it won't save the changes once the 0 has been changed to -10000 and it has been set as the nodatavalue. It will then say that the no data value = n/a.

    I'm not sure why it is doing this.

     – Emtomp Mar 02 '20 at 14:18
  • 1
    uint16 is an unsigned 16-bit integer. This means that your raster can only take values between 0 and 65,535. Try changing the NoData value to a positive value between that range, or change the raster data type to support negative values. – Marcelo Villa Mar 02 '20 at 14:23
  • ah yes, of course. I can't seem to find a gdal tool that automatically does this, so I should create an empty 16bit raster and then save the bands into this new one? – Emtomp Mar 02 '20 at 14:40
  • That is one option. You can also use gdal_translate to create a copy of the raster, specifying a new data type. Here is an example: https://gis.stackexchange.com/a/241448/86131. gdal_translate can be called from a Python script, either as a command line utility or using gdal.Translate(). I suggest you try this and should you have any trouble, create a new question and link it here so I can help you. – Marcelo Villa Mar 02 '20 at 14:59
  • I know I'm not suppose to comment to just say thanks, but you are a gis wizard. thankyou – Emtomp Mar 02 '20 at 15:13
2

I assume you have Spatial Analyst since you can do reclassify raster values.

Using SetNull standalone tool, or using Raster Calculator, you can use SetNull function as follows:

SetNull("RasterName"==0,"RasterName")

You need to change RasterName with the name of the raster data. usually if the data is loaded into ArcMap, you can choose the raster layer from the list.

The expression means change the raster data with zero value to Null and keep other values unchanged.

You can use SetNull standalone tool in a python code - check the help above - if you want to apply it on many raster data.

To assign -10000 to Null pixels, you need to apply another expression in raster calculator:

Con(IsNull("Null_Raster"),-10000,"Null_Raster")

Where Null_Raster is the raster data resulted from the SetNull Expression mentioned above.

ahmadhanb
  • 40,826
  • 5
  • 51
  • 105
  • That's a good idea, but I need the actual nodata value in the metadata to say nodata = -10000. This tool Sets the value as null – Emtomp Feb 26 '20 at 18:59
  • @Emtomp No need to run a specific python code to assign -10000 to NULL data. All you need is to apply another expression in raster calculator. Check the updated answer. – ahmadhanb Feb 27 '20 at 00:55
0

Solved thanks to @Marcelo Villa

I just had to make a small adjustment, litterally just change a lowercase to a capital B in GetRasterband

def nodata_update(image_list):
for i in image_list:
    # pass 1 to modify the raster
    image = gdal.Open(i, 1)

    # get number of bands
    n = image.RasterCount

    for x in range(1, n+1):
        band = image.GetRasterBand(x)

        # read band as numpy array
        arr = band.ReadAsArray()

        # change 0 to -1000
        arr = np.where(arr == 0, -1000, arr)

        # write the new array
        band.WriteArray(arr)

        # set the NoData value
        band.SetNoDataValue(-10000)

        # save changes
        band.FlushCache()

    del image
Emtomp
  • 543
  • 6
  • 15