Steve Lynch: Is it possible to do Cokriging script?

[YumeiHu]

Hi Steven Llynch,

I need to interpolate more than 1000 precipitation maps with Cokriging method, but I don't know if it is possible to use python scripit programming, because I can't find this tool in ArcToolbox.
I saw your old reply about this topic as following:
       "Message - Create a cokriging layer using the Wizard
        - Use the Get and Set Model Parameter tool to modify parameters if required
        - use Create Geostatistical Layer tool to modify the input dataset or variable used in the cokriging layer "
However, I don't really understand, could you give me an example in python script?
I will be very grateful for this.

Thank you very much!

Yumei Hu

[liutiejun]

Hi Steve,

    I have found the mistake in my script, I have misspelled the folder name of input layer.
    However, another question puzzles me that where the outlayer is. Does it only exist in memory as an object?

Bests,
Tiejun Liu

Run the CreateGeostatisticalLayer geoprocessing tool and then look at the syntax reported in the processing dialog.

Steve

[SteveLynch]

Yes it is in memory. To persist it use the SaveToLayerFile geoprocessing tool

Steve

[liutiejun]

Dear Steve,

    Thank you very much.

Bests,
Tiejun Liu

Yes it is in memory. To persist it use the SaveToLayerFile geoprocessing tool

Steve

[liutiejun]

Dear Steve,

    Could you please tell me how to release or remove the outlayer from memory.

Bests,
Tiejun Liu

Yes it is in memory. To persist it use the SaveToLayerFile geoprocessing tool

Steve

[SteveLynch]

Please use the help and search for terms like delete, remove.

[N_G_]

Dear Steve and all,

I have read with attention posts about cokriging. I also have to produce cokriging with more than 2500 variables (Elevation+rainfall).

I'm not sure to be in the good thematic but I try, because some of you have probably solved (I think) my problem in their script.

I'm new in python and I don't manage lists. I would like to use them to create geostatical layer. I used model builder (models work) but I must write now cokriging script without model builder (export from model builder to python script is not enough...).

So, I made list of shapefiles and list of values for the field rainfall. But I have this answer in the python shell:

Traceback (most recent call last):
  File "C:\INTERP\sCRIPT PYTHON\FC4.py", line 75, in <module>
    inputDset1 = fc + NewList
TypeError: coercing to Unicode: need string or buffer, list found

I think that my code is not good and my lists must be transformed before to use them to create geost. layers. But I don't know how and where.......

If someone have an idea, a solution? It would be great.

My code:

import sys, string, os, arcgisscripting
# Create the Geoprocessor object
gp = arcgisscripting.create()

def listFcsInGDB():
    # list all Feature Classes in a geodatabase 
    gp.workspace = "C:\\INTERP\\PJ_class.gdb"

    fcs = []
    for fds in gp.ListDatasets('','feature') + ['']:
        for fc in gp.ListFeatureClasses('','',fds):
            #yield os.path.join(fds, fc)
            fcs.append(os.path.join(fds, fc))
    return fcs

fcs = listFcsInGDB()
for fc in fcs:
    print fc

# List of value in fields
input_dataset =  fc   
Atts = 'PJ_RacPJ'    #field with rainfall values       

rows = gp.searchcursor(fc)
row = rows.next()

NewList = []
for fc in fcs:
    for row in gp.SearchCursor(fc):
        fcValue = row.getvalue(Atts)
        NewList.append(fcValue)
#print NewList

arcpy.CheckOutExtension("GeoStats")

# Load required toolboxes
arcpy.ImportToolbox("C:/Documents and Settings/ArcGIS/Toolbox.tbx")

# Local variables:
Cokriging_xml = "C:\\INTERP\\COK_.xml"
CK = "CK2"
ids1 = "C:\\INTERP\\data.gdb\\dem1000"


for fc in fcs:
    for fcValue in NewList:
        inputDset1 = fc + NewList #Variables
        inputDset2 = ids1            #Covariable
        InputDset = "inputDset1;inputDset2"

# Create cokriging layer (code from model builder)

    tempEnvironment0 = gp.autoCommit
    arcpy.env.autoCommit = "1000"
    tempEnvironment1 = gp.spatialGrid1
    arcpy.env.spatialGrid1 = "0"
    arcpy.GACreateGeostatisticalLayer_ga(Cokriging_xml,InputDset,CK)
    arcpy.env.autoCommit = tempEnvironment0
    arcpy.env.spatialGrid1 = tempEnvironment1

Thanks

Naho

[N_G_]

I found the problem, and how to answer to my question above.
I put my code below if that can help ...

To use a list in the GACreateGeostatisticalLayer:

 

  #Set environment
    arcpy.env.workspace = "C:/INTERP/PJ.gdb"
    # Check out Geostatistical Analyst extension license.
    gp.CheckOutExtension("GeoStats")
    
    #Variables and loop
    inputGA = "C:/INTERP/COK.xml"
    
    # fc = list component/ fcs = list
    for fc in fcs:
    
        #1st variable /"PJ" = targetfield
        inputDset1 = fc + " "+ "X=Shape Y=Shape "+ "F1=" +"PJ"
        
        #2nd variable
        inputDset2 = "C:/INTERP/data.gdb/dem"
        VARCOVAR = inputDset1 + ";" + inputDset2
        
        #outLayer
        out = "cok"  
        outLayer = out + fc
        
        # Process : GACreateGeostatisticalLayer
        for fc in fcs:      
            gp.GACreateGeostatisticalLayer(inputGA,VARCOVAR,outLayer)
        print gp.GetMessages()

note : in this case it's not useful to make field list and value list for the target field