/* text_ks.aml
/* Script for producing three VIC layers of Ksat based on processed
/* STATSGO soils data from Penn State.
/*
/* this is run with &run text_ks from the arc prompt.
/* The required input files that must be in the active workspace are:
/*   mo_stext (this is the grid of soil texture layers clipped to basin)
/*   relate environment relksat.rel must be available
/*   info file text_soil must be available
/*
/* First, remove old versions of the output files
/*
&if [exists ksat_l1_cell -grid] &then
   kill ksat_l1_cell all
&if [exists ksat_l2_cell -grid] &then
   kill ksat_l2_cell all
&if [exists ksat_l3_cell -grid] &then
   kill ksat_l3_cell all
/*
/* Next restore the relate environment. The relate environment relsoil.rel
/* consists of 9 individual relates (txsoil1, txsoil2,..., txsoil9),
/* each of which links a column (item) in one file (the grid containing
/* the soil texture data mo_stext) with an item in a related file (the 
/* INFO file text_soil functioning as a lookup table). 
/*
/* The INFO file was created using the Arctools/Edittools then /File/Database/
/* Table/... The INFO file has columns tcode (0-16), ksat, wcr, and wp.
/* The relates were created in Arctools/EditTools then Arctools/manage/relates
/* where the values to input are: Database (INFO file name), Item (from .VAT
/* file of grid file), Column (common column in INFO file to which Item will
/* be related). Since the info file is sorted min to max, Ordered is 
/* appropriate, though linear also works, for type of relate.
/* The original database has items called Ln_TCODE (n = STATSGO layer number)
/* which correspond to the "tcode" column of the INFO file. 
/*

relate restore relsoil.rel

grid

/* Pixel level computations of Ksat for VIC layer 1 (pixels at 30 arc sec).
/* mo_stext is the full soil texture grid with data on all 11 layers.
/* since there is "no data" for layers 10 and 11 for nearly the entire 
/* database of the grid, I only include layers 1-9.
/* txsoil1, txsoil2,..., txsoil9 are the names of the individual relates
/* which are used to relate soil texture to any other attribute in INFO file

/* note that, if space is a problem, the first three commands can be combined
/* into one. They are broken into separate commands for clarity here.
 
ksat_l1_1 = mo_stext.txsoil1//ksat

ksat_l1_2 = mo_stext.txsoil2//ksat

ksat_l1_pix = mean (ksat_l1_1, ksat_l1_2)

kill ksat_l1_1 all
kill ksat_l1_2 all

/* aggregate to 1/8 degree cell size

ksat_l1_cell = aggregate (ksat_l1_pix, 15, mean, expand, data)
kill ksat_l1_pix all

/* display result for debugging
display 9999
mape ksat_l1_cell
grids ksat_l1_cell

/*************************************************************
/* Continue with merging STATSGO layers 3 through 7 into VIC layer 2.

ksat_l2_3 = mo_stext.txsoil3//ksat

ksat_l2_4 = mo_stext.txsoil4//ksat

ksat_l2_5 = mo_stext.txsoil5//ksat

ksat_l2_6 = mo_stext.txsoil6//ksat

ksat_l2_7 = mo_stext.txsoil7//ksat


ksat_l2_pix = mean (ksat_l2_3, ksat_l2_4, ksat_l2_5, ksat_l2_6, ksat_l2_7)

kill ksat_l2_3 all
kill ksat_l2_4 all
kill ksat_l2_5 all
kill ksat_l2_6 all
kill ksat_l2_7 all

/* aggregate to 1/8 degree cell size

ksat_l2_cell = aggregate (ksat_l2_pix, 15, mean, expand, data)
kill ksat_l2_pix all


/* display result for debugging
/* display 9999
/* mape ksat_l2_cell
/* grids ksat_l2_cell

/***************************************************************
/* VIC layer 3 is comprised of STATSGO layers 8 and 9

ksat_l3_8 = mo_stext.txsoil8//ksat

ksat_l3_9 = mo_stext.txsoil9//ksat

ksat_l3_pix = mean (ksat_l3_8, ksat_l3_9)

kill ksat_l3_8 all
kill ksat_l3_9 all

/* aggregate to 1/8 degree cell size

ksat_l3_cell = aggregate (ksat_l3_pix, 15, mean, expand, data)
kill ksat_l3_pix all

/* ksat_l3_cell is the VIC layer 3.

/* display result for debugging
/* display 9999
/* mape ksat_l3_cell
/* grids ksat_l3_cell