Problem with my code to run cellstatistics

BenjaminJackson · ‎11-22-2018

Update: I have managed to solve this issue, with help from Dan, colleagues and friends. I'm going to put the final working code at the bottom of this post, so skip to the end if you want to see the "answer".

I have a large radar dataset which contains the rainfall depth across the UK every 5 minutes. I want to determine the total rainfall depth for each day. As there are hundreds of files, I decided to use arcpy to solve this problem, which to be frank, I have little experience in. So firstly, here is an example the format for the radar filenames:

metoffice-c-band-rain-radar_uk_201805312155_5km-composite.gif

So if I take the two characters that indicate day (here being day 31) and use that to group all filenames on day 31 into a list, I can then call a list for each day as an input to cell statistics. I wrote some code (which you can see below), but python crashes everytime I run the code. If I print (Day1) it shows a list of all filenames that indicate Day 1, as desired, so I think this may be a problem with the Arcpy element. I previously used a for loop, but I decided to manaully type in every day so that it was more explicit. Anyway, here is my code:-

import os, sys, arcpy
from collections import defaultdict
from arcpy import env
from arcpy.sa import *
arcpy.env.cartographicCoordinateSystem = arcpy.SpatialReference("British National Grid")
arcpy.env.outputCoordinateSystem = arcpy.SpatialReference("British National Grid")
arcpy.env.extent = r"C:\Users\bwj202\OneDrive - University of Exeter\GIS\Leakage\RuralTrunkMains.shp"

# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")
env.workspace = os.getcwd()

dictdata = defaultdict(list)
for filename in os.listdir(os.getcwd()):
     if filename.endswith('.gif'):
          day = filename[37:39]
          dictdata[day].append(filename)

Day1 = dictdata["01"]
Day2 = dictdata["02"]
Day3 = dictdata["03"]
Day4 = dictdata["04"]
Day5 = dictdata["05"]
Day6 = dictdata["06"]
Day7 = dictdata["07"]
Day8 = dictdata["08"]
Day9 = dictdata["09"]
Day10 = dictdata["10"]
Day11 = dictdata["11"]
Day12 = dictdata["12"]
Day13 = dictdata["13"]
Day14 = dictdata["14"]
Day15 = dictdata["15"]
Day16 = dictdata["16"]
Day17 = dictdata["17"]
Day18 = dictdata["18"]
Day19 = dictdata["19"]
Day20 = dictdata["20"]
Day21 = dictdata["21"]
Day22 = dictdata["22"]
Day23 = dictdata["23"]
Day24 = dictdata["24"]
Day25 = dictdata["25"]
Day26 = dictdata["26"]
Day27 = dictdata["27"]
Day28 = dictdata["28"]
Day29 = dictdata["29"]
Day30 = dictdata["30"]
Day31 = dictdata["31"]

print(Day1)
#print(Day2)
#print(dictdata)

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day1, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day1.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day2, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day2.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day3, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day3.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day4, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day4.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day5, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day5.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day6, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day6.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day7, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day7.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day8, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day8.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day9, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day9.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day10, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day10.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day11, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day11.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day12, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day12.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day13, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day13.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day14, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day14.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day15, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day15.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day16, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day16.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day17, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day17.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day18, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day18.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day19, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day19.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day20, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day20.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day21, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day21.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day22, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day22.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day23, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day23.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day24, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day24.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day25, "SUM", "NODATA")

# Save the output 
outCellStatistics1.save("Day25.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day26, "SUM", "NODATA")

outCellStatistics1.save("Day26.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day27, "SUM", "NODATA")

outCellStatistics1.save("Day27.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day28, "SUM", "NODATA")

outCellStatistics1.save("Day28.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day29, "SUM", "NODATA")

outCellStatistics1.save("Day29.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day30, "SUM", "NODATA")

outCellStatistics1.save("Day30.tif")

# Execute CellStatistics
outCellStatistics1 = CellStatistics(Day31, "SUM", "NODATA")

outCellStatistics1.save("Day31.tif")

#exit()‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

My error message is as follows:

outCellStatistics.save(save)

RuntimeError: ERROR 999998: Unexpected Error.

OK so the main problem was that there were memory issues (which was odd as I was able to run cellstatistics in ArcMap) but I also needed the full path in my output. Some of the other changes are implemented purely for the code to be "cleaner". Here is the code:

import os, sys, arcpy
from collections import defaultdict
from arcpy import env
from arcpy.sa import *
arcpy.env.cartographicCoordinateSystem = arcpy.SpatialReference("British National Grid")
arcpy.env.outputCoordinateSystem = arcpy.SpatialReference("British National Grid")
arcpy.env.extent = r"C:\Users\bwj202\OneDrive - University of Exeter\GIS\Leakage\RuralTrunkMains.shp"

# Check out the ArcGIS Spatial Analyst extension license
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")
Data_folder = os.getcwd() + "/Data"
env.workspace = os.getcwd()
print (Data_folder)

dictdata = defaultdict(list)
for filename in os.listdir(Data_folder):
     if filename.endswith('.gif'):
          day = filename[37:39]
          dictdata[day].append(Data_folder + "/" + filename)

# print (dictdata)

Export_folder = os.getcwd() + "/Exports"

if os.path.exists(Export_folder):
    pass
else:
    arcpy.CreateFolder_management(os.getcwd(), "Exports")



print ("Start loop")

for key, value in dictdata.iteritems():
    label = str(key)
    # print value
    outCellStatistics = CellStatistics(value, "SUM",
                                        "DATA")  # Changed NODATA to DATA - I imageine you don't want to write off a cell value if at some point there is no rain?

    # # Save the output
    outCellStatistics.save(Export_folder + "/Day" + label + ".tif")


print ("finished")‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

DanPatterson_Retired · ‎11-22-2018

From the help... yes

Cell Statistics—Help | ArcGIS Desktop

asters_or_constants

[in_raster_or_constant,...]

A list of input rasters for which a statistic will be calculated for each cell within the Analysis window.

A number can be used as an input; however, the cell size and extent must first be set in the environment.

Raster Layer; Constant

BenjaminJackson · ‎11-22-2018

OK, thanks for the reply. I guess better wording would be whether a list in a python formatted list would work, but I guess the answer is yes in this case. I'll edit the question title.

Regardless my code returns an error, which my actual issue, if lists are acceptable then there must be something else wrong with the code.

DanPatterson_Retired · ‎11-22-2018

Cell statistics is designed to drill down through a stack of rasters.

What is your error?

related to....

outCellStatistics1 = CellStatistics([Day1], "SUM", "NODATA")

I think you are double-listing since Day1 should return a list from the dictionary and you have potentially fed it a list of lists

outCellStatistics1 = CellStatistics(Day1, "SUM", "NODATA") # no [ ]

BenjaminJackson · ‎11-22-2018

I made a few modifications based on feedback from yourself and others (changing the question after I finish these edits). I removes the [ ] from the outcellstatistics.

I receive the same error as before which is:

outCellStatistics.save(save)

RuntimeError: ERROR 999998: Unexpected Error.

I also tried changed the dictionary to a list, where:

dictdata = defaultdict(list) is changed to: dictdata = [].

dictdata[day].append(filename) is changed to: dictdata[day] = filename

Error message: dictdata = filename

TypeError: list indices must be intergers, not str

--

Incidentally, if I call Day using print(Day1), I receive a list of rasters that relate to the first day of the month.

DanPatterson_Retired · ‎11-22-2018

Not sure I follow, but if you print Day 1 like you did on line 1 below and it returns a list then

If you put Day1 into a list ….. [Day1] …. as you have done on line 6

#print(Day1)
#print(Day2)
#print(dictdata)

# Execute CellStatistics
outCellStatistics1 = CellStatistics([Day1], "SUM", "NODATA")‍‍‍‍‍‍

You now have a list of lists

 a = [1, 2, 3]

[a] # put a into a list... even though it is a list already

[[1, 2, 3]]  # you now have a list containing a list... a list of lists with one list in it‍‍‍‍‍

Line numbers might help if you have further information.

/blogs/dan_patterson/2016/08/14/script-formatting

BenjaminJackson · ‎11-23-2018

Ah, that is how you format code, useful, thank you.

Well in your first block of code, line 6, I removed the "[ ]" around Day1, I'm guessing that means it is no longer a list of lists. The truncated output for print(Day1) is:

['metoffice-c-band-rain-radar_uk_201805010000_5km-composite.gif', 'metoffice-c-b
and-rain-radar_uk_201805010005_5km-composite.gif', 'metoffice-c-band-rain-radar_
uk_201805010010_5km-composite.gif', 'metoffice-c-band-rain-radar_uk_201805010015
_5km-composite.gif', 'metoffice-c-band-rain-radar_uk_201805012355_5km-composite.gif']‍‍‍‍

DanPatterson_Retired · ‎11-23-2018

Exactly, that is what you need just the plain list, not a list of lists like you had... btw, line 6 above was from your code, which is why I suspected something was going wrong since cell statistics needs a simple list of rasters.

Then my next question was to ensure that your *.gifs were single band (I suspect they are) because they aren't supported otherwise

Raster file formats—ArcGIS Pro | ArcGIS Desktop

BenjaminJackson · ‎11-23-2018

Hi Dan,

Thanks for the reply, and I updated my question to reflect the change in code. Yes the *.gifs are single band, and I've also tried running cellstatistics in Arcmap for a couple of these tifs and they came out fine, so the problem is with the code...in theory.