0

I am working on a single tiff with over 20k bands, each representing daily temperature, and I'll be extracting and preforming simple raster math on roughly 3500 of the bands. The problem I'm running into is that each band is taking ~2 min to process which would equate to a 4 day running time. I am sure there must be a more efficient manner for processing this.

*edit: When running for 1 band in range (1, 2) I am still getting a ~2min process time.

*edit 2: Relevant information of file.
Columns, Rows: 179, 195
Number_of_Bands: 20454
Cell_Size__X.Y: 0.0625, 0.0625
Format: TIFF
Pixel_Type: floating point
Pixel_Depth: 32 Bit

I have run the same script on a subset with ~200 bands from the file above and the processing time per-band is nearly 2 seconds each.

>>import took 0:00:20.362000
>>set up environment took 0:00:08.220000
>>a full band loop took 0:01:43.342000
>>a full band loop took 0:01:42.302000
>>a full band loop took 0:01:41.980000
>>a full band loop took 0:01:41.359000
>>a full band loop took 0:01:41.072000
>>one season loop took 0:08:32.306000

import datetime
t1 = datetime.datetime.now()

import arcpy, os, arcinfo
from arcpy.sa import *

t2 = datetime.datetime.now()
print "import took %s" %  ( t2-t1)  
arcpy.env.overwriteOutput = True

# Set the current workspace
arcpy.env.workspace = "C:\\raster\\folder"
arcpy.CheckOutExtension('Spatial')

# Get and print a list from the workspace
rasters = arcpy.ListRasters('*', 'TIF')
cutoffTemp = float(277.2)

t3 = datetime.datetime.now()
print "set up environment took %s" %  ( t3-t2)

for raster in rasters:
    t4 = datetime.datetime.now()
    print(raster)
    d = arcpy.Describe(raster)
    number_of_bands = d.bandCount #number of bands from Hist2005
    print(number_of_bands)

    arcpy.AddMessage("Processing raster {} for {} bands".format(raster, number_of_bands))
    listcomp90 = [] #create a list of output tifs, empty for each input raster

    for band_number in range(1, 6): #small range for test
        t6 = datetime.datetime.now()
        listcomp91.append(band_number)
        Hist = arcpy.Raster("{}\\Band_{}".format(raster, band_number)) #gets inidiviual bands for input
        arcpy.AddMessage('in the band number cycle and the current band number is {}'.format( band_number))
        output = arcpy.sa.Con(Hist >= cutoffTemp, sapHist) #if band values >= cuttoffTemp then temp     
        Rasname, Rasext = os.path.splitext(raster)
        output.save(os.path.join(arcpy.env.workspace, "New{}_91_comp_{}.tif".format(Rasname, band_number)))
        listcomp91.append(os.path.join(arcpy.env.workspace, "New{}_91_comp_{}.tif".format(Rasname, band_number)))
        t7 = datetime.datetime.now()
        print "a full band loop took %s" % (t7-t6)

t5 = datetime.datetime.now()
print "one season loop took %s" % (t5-t4)
  • 1
    Welcome to GIS SE! As a new user please take the tour to learn about our focused Q&A format. What happens if you extract a single band into a separate file and process just the one - How long does it take? Please edit your question to include this information. – Midavalo Apr 15 '17 at 19:49
  • Thanks for the information Midavalo! I reran it and updated the original post. I am still getting a ~2 minute processing time. – dil_lockd Apr 15 '17 at 20:20
  • 1
    Please edit the question to provide information on the size of the rasters (rows & cols, and datatype & bit-depth) – Vince Apr 15 '17 at 21:31
  • 1
    Band organization is also a significant factor -- BIL or BIP require far more I/O. – Vince Apr 15 '17 at 23:32
  • 1
    If your processing is still 2 minutes on a single band outside of a loop, then I'd suggest that the looping isn't the problem – Midavalo Apr 16 '17 at 0:08
0

It looks like the size of file or the number of bands is what was slowing the processing time down. To fix this I broke the script into two parts: first, composite the bands for each season; then pass that into a cleaned up version of the file above. Additionally, I chose not to re-composite the final output. Whole process for the two scripts on ~200 bands is about 8min.

Example of band composite script:

for raster in rasters  
d = arcpy.Describe(raster)
number_of_bands = d.bandCount #number of bands
listcomp = [] #create a list of output tifs, empty for each input raster

    for band_number in range(3013, 3227):

        listcomp.append(os.path.join(arcpy.env.workspace, "{}\\Band_{}".format(raster, band_number)))  #adds raster bands without needlessly re-saving them

#collect the output tifs into a composite
arcpy.CompositeBands_management(listcomp, os.path.join(arcpy.env.workspace,"{}_Composite.tif".format(raster)))

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.