Need Help Automation of Kriging using Model builder or python - Reg.

20442
53
09-12-2013 05:58 AM
AnushuyaRamakrishnan
New Contributor
I am using arcmap 10.2.

My project deals with the prediction of PM and Ozone concentrations for 2004-2006 (1098 days) using known concentrations of them at
99 grids throughout the county.

I am trying to predict the unknown concentrations at mother's residences using ordinary kriging.

This is what I am doing:

1. Use Geostatistical Wizard
2. Select kriging cokriging
3. Input the Data layers and specifications
4. Select Finish
5. Read the model parameters and click OK
6. Generate Kriged predictions
7. Right click on Kriged data and convert to raster
8. Use spatial analyst tool Extraction
9. Select "Extract values to points"
10. Generate Extracted values of kriged predictions
11. Select conversion tool
12. Choose convert From Raster to Ascii
13. Generate the Ascii table.
14. Alternatively use Sample tool and input multiple rasters to generate an output table.

My question:

1. I need some guidance for automating this kriging process using model builder or python?

2. Is there any way to directly copy and paste the raster values to excel?

Can any one provide me a lead into this
0 Kudos
53 Replies
AnushuyaRamakrishnan
New Contributor
Thank you Eric!!!

I will try this and get back to you

Anushuya
0 Kudos
AnushuyaRamakrishnan
New Contributor
Hi Eric,

I tried to run this in python:

# Import arcpy module
import arcpy

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

merge1.csv.events = arcpy.GeostatisticalDatasets("C:\\temp\\KrigingordinaryJan152.xml")
Harris_Lat_Long.csv.events = "C:\\temp\\PM04Jan16.shp"
merge1.csv,events.Field = "Y04PMJan16"
arcpy.GACreateGeostatisticalLayer_ga("C:\\temp\\KrigingordinaryJan16.xml", newDataset, "PM04Jan16")

It is throwing me an error:

Traceback (most recent call last):
  File "C:\Users\aramakrishnan\Documents\ArcGIS\scratch\Ordinary Kriging_09252013_01.py", line 19, in <module>
    merge1.csv.events = arcpy.GeostatisticalDatasets("C:\\temp\\KrigingordinaryJan152.xml")
  File "C:\Program Files (x86)\ArcGIS\Desktop10.2\arcpy\arcpy\_ga.py", line 300, in GeostatisticalDatasets
    return GADatasets(gapy.GALayerInfo(ga_model_source))
TypeError: bad argument type for built-in operation


Am I missing something?
>>>
0 Kudos
EricKrause
Esri Regular Contributor
You should change your code to something like:


X = arcpy.GeostatisticalDatasets("C:\\temp\\KrigingordinaryJan152.xml")  #Make sure your xml is saved in C:\temp for this to work
X.dataset1 = "C:\\temp\\PM04Jan16.shp" #Make sure your shp file is in C:\temp
X.dataset1Field = "Y04PMJan16"
arcpy.GACreateGeostatisticalLayer_ga("C:\\temp\\KrigingordinaryJan16.xml", X, "PM04Jan16")

You can replace the X with any name you want.  The first line creates a GeostatisticalDatasets object called "X".  It populates the necessary properties by reading your XML.  The second line sets the dataset1 property of X.  The third line sets the field for X.  And the last line uses the XML source and the GeostatisticalDatasets object ("X"), and it creates a geostatistical layer called "PM04Jan16".
0 Kudos
AnushuyaRamakrishnan
New Contributor
Hi Eric,

I could succesfully generated new geostatistical layer from an existing xml input.

Can you guide me to the next step of kriging.

For kriging using SA I use the following script:

# Process: Kriging
arcpy.gp.Kriging_sa(merge3_csv_Events_lyr, "longitude", fieldname, "Exponential 0.003512 0.255986 1.092512 0.001093", "3.5120478474024E-03", "VARIABLE 12", Output_variance_of_prediction_raster)

1. I need the script for kriging using ga

2. What is the input for this script? Is it the same as SA.

Pls. let me know.

Thank you

Anushuya
0 Kudos
EricKrause
Esri Regular Contributor
The creation of a new geostatistical layer is the automation: that new geostatistical layer is the result of kriging on the new dataset.  Once you have the geostatistical layer, you can export it to a raster with GA Layer to Grid, or you can predict directly to a point feature class with GA Layer to Points.

Now you just need to iterate all your datasets through the Create Geostatistical Layer tool.
0 Kudos
AnushuyaRamakrishnan
New Contributor
Thank you Eric... I will do that

Anushuya
0 Kudos
shabijunejo
New Contributor
hi, I am new to arcgis, I was working with georeferencing tool,,, i have followed all the process but in the and when I apply project raster the whole screen turns black and when I zoom to full extent it shows only a tiny part of it.... please help me if you may...

thanks
0 Kudos
AnushuyaRamakrishnan
New Contributor
Hi,

Do you input a raster or a layer or a data file like CSV.

Please go to add tools and add your specific input.

Then go to add xy tool and project it in arcmap.

Then perform specific geostatistical operation using the events layer and get the output layer

Next right click on layer and convert it to raster.

Use this raster to generate your output table.

In any module of arcmap the screen should not go blank.

Try these steps and check if it works.

Also try and reinstall arcmap if the problem persists.

Thanks

Anushuya
0 Kudos
AnushuyaRamakrishnan
New Contributor
Hi Eric,

I generated my ordinary kriging output successfully using the python script for generating geostatistical layer and using it to generate raster and convert it to points and generate the table.

However, I see that the output generated by scripting does not exactly match with the output that I generate manually by running
ordinary kriging for the same day using the same values of range, sill and nugget.

In the script I only set auto=true for lagsize, number of lags, range, sill and nugget values.

I attach the python script for your review.

Please suggest.

Thank you

Anushuya
0 Kudos
EricKrause
Esri Regular Contributor
Your Python script looks fine.  Can you post your XML file?
0 Kudos