# Calculating median of all pixels of all layers of raster stack in ArcGIS for Desktop?

I have a stack of 15 layers and I need to calculate the median of all the pixels of all the layers.

How can i do so in ArcGIS?

I tried zonal statistics but it only calculates the median among pixels in one layer.

I also tried Cell Statistics but then it calculates the median of every single pixel trough time.

I need to consider all the values of all the pixels of all the layers together and calculate the median of that distribution.

• Please decide which GIS software you are trying to use and then tell us what you have tried do far in that software.
– PolyGeo
Nov 28, 2016 at 19:39
• OK. Can be both Erdas or ArgGIS. In Erdas I tried with Stck Median but it gives you the median of each pixel. In ArcGIS I tried Zonal statistics but it doesnt work with stacks, only single layers. I need to compute the median of all the pixel values of all the different layers of the stack, not per layer, not per pixel Nov 28, 2016 at 20:44
• It would be best to edit the question with the information in your comment. Nov 28, 2016 at 21:28
• It is unclear what your intended output looks like. You say that calculating a pixel based median for all the layers won't work nor will calculating the median for each layer. Are you asking how to add all the pixel values in all the data layers into one big list and calculate the median value of that list?
– Aaron
Nov 28, 2016 at 23:14
• Exacly Aaron, thats what I need: join all the pixels of all the different layers into one big distribution and calculate the median of it (Eventually I will also need the 20th and 5th percentiles...but the median is ok to start with...) Nov 28, 2016 at 23:32

The following code should work. With a few lines more you could also calculate the percentiles you need.

``````import arcpy
import numpy as np

arrays = []
arcpy.env.workspace = r"C:\..."  # The directory containing your rasters

for raster in arcpy.ListRasters():
arr = arcpy.RasterToNumPyArray(in_raster=raster)
arrays.append(arr)

stack = np.stack(arrays)   # 3d numpy array
median = np.median(stack)  # returns one single value
``````