# Calculating first quartile of a raster by zonal statistics?

I want to calculate the 1st quartile (25th percentile) of a raster (beneath polygons) using zonal statistics but it does not this option.

Zonal statistics has: MIN, MAX, MEAN, MAJORITY, RANGE, and does not have 1st quartile. How can I do that?

• I guess you are using Arcgis? do you have access to other software? – dof1985 Jul 17 '16 at 15:15
• Also, see this thread from Geonet; if it works you can post the code as an answer. – dof1985 Jul 17 '16 at 15:19
• thanks for your consideration. I am using ArcGIS. I looked at that script it is used for percentile not calculation the 1th quarter. I have written some codes using arcpy and it is preferable to use ARCGIS. – Khandan Jul 17 '16 at 18:44
• 1st quarter = 25th percentile. Change line 25 from `lst_perc = [2, 5, 10, 90, 95, 98] ` to `lst_perc =  ` – user2856 Jul 17 '16 at 19:30
• thanks for your precious comment. Is it possible to calculate the average of 1th quarter of values? – Khandan Jul 17 '16 at 20:51

## 1 Answer

Dear dof1985 thanks for your comment I did it by the following code :

``````import arcpy

from arcpy import env

from arcpy.sa import *

import arcpy, os

import arcpy.da

env.workspace = "E:/spring/cases/main/9_jan_2016/New Folder/output"
output1=r"E:\main\9_jan_2016\New Folder"
output=r"E:\main\9_jan_2016\New Folder\output\test"
output2=r"E:\main\9_jan_2016\New Folder\output\reza1.gdb"
output3=r"E:\main\9_jan_2016\New Folder\output"

fc=output3+os.sep+"1800.shp"
inraster=output1+os.sep+"1800TIR1BTC.tif"
fc2=output2+os.sep+"percent07"
g=0
arcpy.CreateTable_management (output2,"percent07",       output3+os.sep+str(1800)+".shp")
arcpy.CheckOutExtension("Spatial")
inraster=Int(inraster)

for row in arcpy.da.SearchCursor(fc, ["FID","SHAPE@"]):
print(row)
average=[]
outfc = os.path.join(output, "fc" + str(g))
arcpy.CopyFeatures_management(row, outfc)
# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")

# Execute ExtractByPolygon

TabulateArea(inraster,"Value",output+os.sep+ "fc" + str(g)+".shp","Id",output+os.sep+"reza"+str(g)+".dbf")
A=int(arcpy.GetCount_management(output+os.sep+"reza"+str(g)+".dbf").getOutput(0))
print("A"+str(A))
counti=A*0.25
print("counti"+str(counti))
b=0
i=0
if (A!=0):

arcpy.Sort_management(output+os.sep+"reza"+str(g)+".dbf", output+os.sep+"rezasorted"+str(g)+".dbf",[["Value","ASCENDING"]])
for row in arcpy.da.SearchCursor(output+os.sep+"rezasorted"+str(g)+".dbf", ["Value"]):
if  i<=counti:
#cnt_sum = sum(extPolygonOut.values())
print("row"+str(row))
b=b+row
i=i+1
else :
print("loop  finished")
break

a=b/i

average.append(a)
with arcpy.da.InsertCursor(fc2,["TIR1"])as cursor :
for row in average:

cursor.insertRow(average)

arcpy.Delete_management(output+os.sep+ "fc" + str(g)+".shp")
arcpy.Delete_management(output+os.sep+"reza"+str(g)+".dbf")
arcpy.Delete_management(output+os.sep+"rezasorted"+str(g)+".dbf")
g=g+1
print(average)
#del(average)
else :
arcpy.Delete_management(output+os.sep+ "fc" + str(g)+".shp")
arcpy.Delete_management(output+os.sep+"reza"+str(g)+".dbf")
arcpy.Delete_management(output+os.sep+"rezasorted"+str(g)+".dbf")
g=g+1
``````