Using R to write SWAT input was:Re: [SWAT-user:1363] Re: Running SWAT without its interface !

[willem vervoort]

Jim,
It actually is much simpler than I thought. Attached is a little
script that rewrites the HRU files, let me know if this causes any
problems.

The pseudo code is:
set directory

read in file names from directory that end in ".hru"

cycle through loop, open each file and read the text

rewrite the lines you want

write the file back out

close file

I think this can be done more elegantly using "apply" than with a for
loop and that will make it faster. I will rewrite and post that too.

Willem

On Fri, Aug 8, 2008 at 8:36 AM, willem vervoort
<willemv...@gmail.com> wrote:
> Jim,
> That is a great idea, I am always happy to point out the advantages of
> R. I will do that, just give me some time to organise a good example.
> WIlle
>
> On Fri, Aug 8, 2008 at 3:45 AM, Jim Almendinger <din...@smm.org> wrote:
>> Willem and others --
>> I've dabbled in R for simple analyses and plotting -- and I know I'm missing
>> much of its many capabilities. If you find it to be useful for your
>> analyses, it would be instructive to post a sample script, along with a
>> short example data file -- especially when you've used it for
>> pre-processing. I assume you're reading in an input file, searching for the
>> lines with the parameters you want to change, parsing the line to find that
>> parameter, changing the parameter, and then re-writing the file in the
>> format needed by SWAT. I don't know any of these commands, and I'm always
>> happy to try and learn some new tricks.
>> As an aside, we've used small Fortran programs to post-process files, to
>> extract out just the data we want. But for pre-processing, the most I've
>> done is to use update queries in Access to change selected parameters in the
>> input tables -- which still requires the interface to re-write the text
>> input files.
>> Thanks,
>> -- Jim
>>
>> On Aug 6, 2008, at 5:51 PM, willem vervoort wrote:
>>
>> Dear Jaswant,
>> running SWAT without the interface works fine. I have used R
>> (www.r-project.org) for pre and post processing of data files for
>> running different scenarios, but you could use any other math command
>> line structure. I have only used the AVSWAT interface to set up the
>> watershed as I haven't yet bothered working out how to write the fig
>> file, but I am sure that can be done as well. You can run SWAT2005.exe
>> just from the dos command prompt, or via another command structure.
>> That way you can also link SWAT to PEST if that is what you want to
>> do.
>> The main trick is to pay very good attention to the manual in terms of
>> whether the input files are FORTRAN formatted or just FREE.
>> Willem
>>
>> Dr. James E. Almendinger, Senior Scientist
>>
>> St. Croix Watershed Research Station
>>
>> Science Museum of Minnesota
>>
>> 16910 152nd St. N
>>
>> Marine on St. Croix, MN 55047
>>
>> tel: 651-433-5953 X 19
>>
>> fax: 651-433-5924
>>
>> email: din...@smm.org
>>
>> web: www.smm.org/SCWRS/
>>
>>
>
>
>
> --
> http://blogs.usyd.edu.au/waterhydrosu/
>

--
http://blogs.usyd.edu.au/waterhydrosu/

[willem vervoort]

Jim,
here is another version and a file to rewrite the hru files using
sapply. I did not notice any change in speed directly, but have not
tested that thoroughly. Note that I read in from the "txtinout" dir
and read out to the "Pest" dir. That way I don't overwrite my original
files. I had to create the Pest dir first manually (Only the dir not
the files). I have about 488 hru's so this works quite rapidly.

I use tinn-R for scripting, but you can use others.
(http://www.sciviews.org/Tinn-R/ and
http://sourceforge.net/projects/tinn-r)

Attached is also a file to plot observed vs predicted from the output
of the autocalibration (which you can also run from DOS), you can
write a similar file to process the Pest output.
Willem

On Fri, Aug 8, 2008 at 8:46 AM, willem vervoort

--
http://blogs.usyd.edu.au/waterhydrosu/

[willem vervoort]

Oops I think I attached an older version of the read output file, just
to make sure here is the latest version
Willem

On Fri, Aug 8, 2008 at 9:10 AM, willem vervoort

--
http://blogs.usyd.edu.au/waterhydrosu/

[Naresh]

You can use

dir.create('C:\\willem\\research\\swat2005\\Pest)

to create the Pest directory automatically.

>  ObsvsPred.r
> 1KDownload- Hide quoted text -
>
> - Show quoted text -

[Rohan Sadler]

Hi All,

NB: under linux: system("mkdir C:\\willem\\research\\swat2005\\Pest")
works as well.

Here is a very rough code hack in R for writing .mgt, .hru, .sol and
.gw files as part of manual calibration specific to our project.

We had 5500 hru's over 23 subbasins, I wasn't able to specify unique
hrus through soil x land-use combinations in ArcSWAT as each hru was a
paddock whose spatial location we wanted to keep. Ideally arcswat
would have dealt with 5500 unque (dummy) land uses, but it couldn't.
We had to write out the hru files outside of ArcSWAT.

I am not going to explain it all here, the code works but is currently
tiresome to fiddle with (and from your point of view to interpret -
the project was very much time limited so the hacks were cheap and
nasty). It was built so as to rewrite landuses to see effects of
landuse change on recharge rates. The second half of the "gen" codes
is more useful (i.e. general) than the first half. Swat is run from
inside gen_sol_zana_2.R using specified parameters and saves output
files to separate directories for latter collection and compilation
(ultimately into plots).

From the discussion, it sounds like it would be useful to aim towards
putting together a SWAT-R package under R to do specified tasks, which
would need to be discussed and planned so as to have broad usability.
The key is to get an autocalibration running. It would certainly have
saved me many hours of pain. I am sorry to say but I can't help overly
in this, as I no longer work on the project where I was regularly
using SWAT.

1. hru/mgt (gen_man_2.R)
g.fun(z=3,hru.yes,fd=1,party=1,n.fert=40,cn2.mult=1,nveg=FALSE)
* runs from a database with every .hru linked to an arcgis layer (ie.
hru10.dbf).
* according to a specified distribution, assigns land uses randomly
(so we can investigate the roll of increasing perennial vegetation on
recharge rates: z variable below relates to the percentage of alfalfa
in the landscape. nveg=TRUE means all of landscape is native
vegetation
* accesses cn2 numbers; writes to different directories to be run
simultaneously by a quad processor (fd); (party is a redundant
variable).
* can manipulate nitrogen levels; overwrites land use templates for
management operations such as awht_mgt.txt.
* I was going to read EPSCO and ESCO from a table for different land
uses, and include that as arguments int he function, so the prototype
of that exercise is in the function.

2. sol/gw (gen_sol_zana_2.R)
solgw(fd=1,party=1,gw.yes,sol.yes,sdepth=15000,gw_delay=30,alpha_bf=0.5,
gwqmin=500,gw_revap=0.05,revapmin=1000,rchrg_dp=0.1)
* runs swat from inside the function.
* gw.yes=TRUE means write gw files, similarly sol.yes
* arguments are different gw parameters.
* sdepth: we specified our soils (deep laterites) up to 15m depth; and
then experimented with altering how deeply we went with our .sol
parameters (e.g. restricted all soils to 2m depth).
* read from a .csv spreadsheet of soil parameters co-opted from
arcswat soil database table.

3. simlist1.R
* a batch file to run different parameters in different working
directories, run through a SSH tunnel to our quad machine:
R CMD BATCH simlist1.R
* doing a OA-latin hypercube would entail writing a script that writes
a host of files such as simlist1.R.

I do have some working proto-code for plotting sensitivity analyses in
R from output.* files, and responses over time from different
parameter runs. But it is rough, and perhaps incommunicable (it is
targetted towards the proportion of landscape as alfalfa). I also have
code for rewriting Australian Bureau of meterology files as multiple
stations to a pcp.pcp file by finding the nearest weather station to
the centroid of each basin; and to extract summaries of the weather
data so as to write the .wgn file (apart from dewpoints and the like
which I had to find from other places).

Regards
Rohan Sadler

2008/8/8 Naresh <nares...@gmail.com>:

--
Research Associate
School of Plant Biology
M089 The University of Western Australia
35 Stirling Highway,
Crawley, 6009

Office: +61 8 6488 3491
Fax: +61 8 6488 1098
Mobile: 0433 192 600

[Jim Almendinger]

Willem, Rohan, Naresh, et al. --

Many thanks for the R scripts -- there's a lot to be learned by working through the examples.  I'm going to have to start with Willem's basic scripts, which are about at the level I can digest for now.  But Rohan's long and intricate scripts are good examples of just how far one can use R as a basic programming language to manipulate the data files.  

I think eventually a SWAT-R package would be very useful, especially since it would be open-source and freely available internationally.  I don't fully understand what constitutes such a package.  My instincts say that it is better, at least for some and in the beginning, to have the simplest, most basic array of tools that allow user flexibility to tailor the scripts to their particular problem.  Perhaps for now, just a set of example scripts demonstrating the most useful commands would be enough.  I can see that eventually, a full package that provides pre- and post-processing capabilities would be useful, especially given the statistical and graphical power of R (which is, of course, what it was designed for).  The concern is -- perhaps misplaced -- that the more complicated and user-friendly such a a program becomes, the less flexible it is to the user, and the more of a black-box it becomes.  And we don't want to re-invent the wheel that Srini and others have already put into the ArcSWAT inferface, which has continued to improve the usability of SWAT.  

It will be interesting to watch how R usage among SWAT users will (or won't) grow, and whether it will evolve into a SWAT-R package.  I have to think that R usage will grow -- I mean, cutting and pasting data and hand-parsing text files into Excel and Access is like using a rock to pound in nails.  Makes me feel like a Neanderthal.  

Thanks again for the advice and examples,

-- Jim

-- 

BEGIN-ANTISPAM-VOTING-LINKS

------------------------------------------------------

Teach CanIt if this mail (ID 5130981) is spam:

Spam:        https://canit.smm.org/canit/b.php?i=5130981&m=3540b52ec6e5&c=s

Not spam:    https://canit.smm.org/canit/b.php?i=5130981&m=3540b52ec6e5&c=n

Forget vote: https://canit.smm.org/canit/b.php?i=5130981&m=3540b52ec6e5&c=f

------------------------------------------------------

END-ANTISPAM-VOTING-LINKS

 .mgt file Subbasin:1 HRU:1 Luse:AWHT Soil: LIMERICK Slope: 0-9999 13/09/2007 12:00:00 AM ARCGIS-SWAT2003 interface MAVZ

               0    | NMGT:Management code

Initial Plant Growth Parameters

               0    | IGRO: Land cover status: 0-none growing; 1-growing

               0    | PLANT_ID: Land cover ID number (IGRO = 1)

            0.00    | LAI_INIT: Initial leaf are index (IGRO = 1)

            0.00    | BIO_INIT: Initial biomass (kg/ha) (IGRO = 1)

            0.00    | PHU_PLT: Number of heat units to bring plant to maturity (IGRO = 1)

General Management Parameters

            0.20    | BIOMIX: Biological mixing efficiency

           81.00    | CN2: Initial SCS CN II value

            1.00    | USLE_P: USLE support practice factor

            0.00    | BIO_MIN: Minimum biomass for grazing (kg/ha)

           0.000    | FILTERW: width of edge of field filter strip (m)

Urban Management Parameters

               0    | IURBAN: urban simulation code, 0-none, 1-USGS, 2-buildup/washoff

               0    | URBLU: urban land type

Irrigation Management Parameters

               0    | IRRSC: irrigation code

               0    | IRRNO: irrigation source location

           0.000    | FLOWMIN: min in-stream flow for irr diversions (m^3/s)

           0.000    | DIVMAX: max irrigation diversion from reach (+mm/-10^4m^3)

           0.000    | FLOWFR: : fraction of flow allowed to be pulled for irr

Tile Drain Management Parameters

           0.000    | DDRAIN: depth to subsurface tile drain (mm)

           0.000    | TDRAIN: time to drain soil to field capacity (hr)

           0.000    | GDRAIN: drain tile lag time (hr)

Management Operations:

               1    | NROT: number of years of rotation

Operation Schedule:

  5 20           1   98          2200.00000            

  5 21           4  177             1.00000    1.00

  5 22          11    7             0.75000             

  5 23           6    4 

 12  7           5

<gen_man_2.R><gen_sol_zana_2.R><simlist1.R>