perform IDW interpolation to select multiple Z field from attribute table using ArcPy

ShouvikJha · ‎08-21-2016

I have 12 point feature datasets for different month. I want to perform IDW interpolation.

I have created the python code but its taking only one "Z" field for different month feature datasets from the attribute. Below I have mentioned the code. How to assign here multiple Z field from each month folder (JANUARY, FEBRUARY......DECEMBER).

Point feature class from different month folder (e.g. For January Merged_001, for FEBRUARY Merged_002.......for DECEMBER Merged_012 ) and want to perform IDW interpolation and save it as on same month folder as it Z field name (Max_Temper, Min_Temper, Precipitatat, Wind, Relative_H, Solar) with month (e.g for January 001_Max_Temper, 001_Min_Temper, 001_Precipitatat, 001_Wind, 001_Relative_H, 001_Solar).

After Whatever Solar raster for different month (for January 001_Solar....for DECEMBER 012_Solar) we will get, I want to do some calculation (001_Solar * 30) *0.5 and save it 001_sr ........012_sr in Month wise folder.

Below i have attached Point shapefile also for reference purpose

Here is my code.

import arcpy
from arcpy import env
from arcpy.sa import *
arcpy.env.parallelProcessingFactor = "100%"
arcpy.env.overwriteOutput = True
# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")
env.workspace = "D:\SWAT-WEATHER-DATA2\APRIL"
# Set local variables
inPointFeatures = "D:\SWAT-WEATHER-DATA2\APRIL\Merged_004.shp"
zField = "Min_Temper" #Z_Filed = Max_Temper, Min_Temper, Precipitatat, Wind, Relative_H, Solar
cellSize = 0.002298707671
power = 2
searchRadius = RadiusVariable(10, 150000)
#Mask region of interest
mask="D:\Gujarta Shape file\GUJARATSTATE.shp"

# Execute IDW
outIDW = Idw(inPointFeatures, zField, cellSize, power, searchRadius)
# Execute Mask
IDWMASk = ExtractByMask(outIDW, mask)
# Save output, except Solar Radiation raster
IDWMASk.save("004_Min_Temp.tif")
# 004_Max_Temp.tif, 004_Min_Temp.tif....Solar
#Only for Solar Radiation raster
#PAR = (IDWMASk * 30) * 0.5
#Save output, Only for Solar Radiation raster
#PAR.save("003_sr.tif")
#print done
print 'done'‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

DanPatterson_Retired · ‎08-30-2016

fix the potential errors I posted about regarding file paths etc

ShouvikJha · ‎08-30-2016

Dan_Patterson yes i have added r in file path but still the same error

DanPatterson_Retired · ‎08-30-2016

I gave an extensive list! not just the r

Move your files to comply and fix all on the list

Read Xander's last post regarding the input data containing duplicates.

Take your time and work through each step

ShouvikJha · ‎08-30-2016

I am trying to fix it.

ShouvikJha · ‎08-30-2016

Earlier i had post the issues of taking multiple Z value to perfrom IDW. From IDW interpolation its seems to be not much accurate. I am trying to execute Universal Kriging (LINEARDRIFT). I have written one code to perform kriging, but its not running, And i am confused to select some variable for Kriging. Below i have attached point shape file, which i am performing Kriging.

Error massage from below code while performing Kriging

Message File Name Line Position
Traceback
 <module> D:\Arc-GIS-Python-Script\KRIGING_PAR1.py 31
 Kriging C:\Program Files\ArcGIS\Desktop10.1\arcpy\arcpy\sa\Functions.py 2405
 swapper C:\Program Files\ArcGIS\Desktop10.1\arcpy\arcpy\sa\Utils.py 47
 wrapper C:\Program Files\ArcGIS\Desktop10.1\arcpy\arcpy\sa\Functions.py 2397
 <lambda> C:\Program Files\ArcGIS\Desktop10.1\arcpy\arcpy\geoprocessing\_base.py 484
RuntimeError: Object: Error in executing tool‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Below is the code

import arcpy
from arcpy import env
from arcpy.sa import *
arcpy.env.parallelProcessingFactor = "100%"
arcpy.env.overwriteOutput = True

# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")

env.workspace = r"D:\SWAT-WEATHER-DATA"

# Set local variables
inPointFeatures = r"D:\SWAT-WEATHER-DATA\SELECTED-EX-MERGE.shp"

zFields = ["Max_Temper", "Min_Temper", "Precipitat", "Wind", "Relative_H", "Solar"]

#Kriging Veriable
cellSize = 0.002298707671
lagSize = 0.5780481172534
majorRange = 6
partialSill = 3.304292110
nugget = 0.002701348
kRadius = RadiusFixed(20000, 1)

#Mask region of interest
mask="D:\Gujarta Shape file\GUJARATSTATE.shp"

# Execute Kriging
for zField in zFields:
    kModelUniversalObj = KrigingModelUniversal("LINEARDRIFT", lagSize, majorRange, partialSill, nugget)
    OutKriging = Kriging(inPointFeatures, zField, kModelUniversalObj, kRadius)
    # Execute Mask
    #IDWMASk = ExtractByMask(outIDW, mask)
    IDWKrig = ExtractByMask(OutKriging, mask)

    if not zField == "Solar":
    # Save output, except Solar Radiation raster
     IDWKrig.save("012_{}.tif".format(zField))

    else:
     #Only for Solar Radiation raster
     PAR = (IDWKrig * 30) * 0.5
     #Save output, Only for Solar Radiation raster
     PAR.save("012_{}.tif".format(zField))

#print done
print 'done'‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

XanderBakker · ‎08-30-2016

I have a number of observations:

You state that you are confused about the lagSize, majorRange, partialSill and nugget. They are all part of the semivariogram. Tools to model these parameters are available with the GeoStatistical Analyst extension. More on Kriging here: How Kriging works—Help | ArcGIS for Desktop and more specifically see the topic "Understanding a semivariogram—Range, sill, and nugget". For GeoStatistical Analyst and the sensitivity of the semivariogram, see this topic: Semivariogram Sensitivity—Help | ArcGIS for Desktop
If you don't have the tools to define the lagSize, majorRange, partialSill and nugget correctly, it is better not to specify them.
Your fixed radius definition should be specified in map units. If you have decimal degrees a distance of 20000 decimal degrees will be "rather" big.

When I executed the tool with the default values provided by the tool (but changing to Linear drift). The snippets of the tool executions was:

arcpy.gp.Kriging_sa("Merged_004", "Max_Temper", "D:/Xander/GeoNet/Kriging/tif/krimaxtemp01.TIF", "LinearDrift 0,023729", "2,37291946411132E-02", "FIXED 0,118646 20", "D:/Xander/GeoNet/Kriging/tif/var01.TIF")

With this result.

Personally, I don't think it makes much sense, but an important reason for this are the overlapping points. It would be better to get the average of the month into a single point for each location and use that instead.

ShouvikJha · ‎08-30-2016

Xander Bakker‌. Thank you. as you mentioned take average of the month into a single point for each parameters,

By default interpolation techniques takes the mean value of each coincident points for each parameters . I have added one snapshot of this.

ShouvikJha · ‎08-30-2016

Xander Bakker‌. As you said to skip those variable which i don't have tool to find out. If i skip those variable is there any effect on accuracy of output.

Please can you tell me, how to skip those variable by defining default values for Kriging. Kriging is tested by several researcher therefor would like to employ Kriging.

XanderBakker · ‎08-30-2016

When you skip variables you will affect accuracy. However, specifying wrong variables, will affect accuracy even more. Best option is to se Geostatistical Analyst. Just scroll through the exercises of the tutorial to get a general idea what a normal workflow with Geostatistical Analyst involves: Exercise 1: Creating a surface using default parameters—Help | ArcGIS for Desktop (scroll down to 50% of the page) and more specifically Exercise 2: Exploring your data—Help | ArcGIS for Desktop and Exercise 3: Mapping ozone concentration—Help | ArcGIS for Desktop to get the general idea.

ShouvikJha · ‎08-30-2016

Xander Bakker Following the same i have assigned the variable for Kriging. But had some confusion i am on right track or not. I have attached some step to take the variable which i used for kriging.