hood
-------------------------------------------------------------
new field survey data
--------------------------------------------------------------
NOTE THAT -d is not needed with the new USPED
s.surf.rst survey.site elev=sel.100 ten=100 smo=0.2 npmin=300 dmin=2 slo=ssl.100 asp=sas.100 pcu=spc.100 tcu=stc.100 |fi=27
s.surf.rst survey.site elev=sel.150 ten=150 smo=0.05 npmin=300 dmin=2 | fi=41 prilis velky bump okolo gully 9skus sm=0.2
s.surf.rst -d survey.site elev=sel.80a ten=80 smo=0.3 npmin=300 dmin=2 slo=sdx.80 asp=sdy.80 | fi=22
s.surf.rst survey.site elev=sel.80a ten=80 smo=0.3 npmin=300 dmin=2 slo=ssl.80 asp=sas.80 pcu=spc.80 tcu=stc.80| fi=22
s.surf.rst survey.site elev=sel.55 ten=55 smo=0.4 npmin=300 dmin=2 slo=ssl.55 asp=sas.55 |fi=15
s.surf.rst survey.site elev=sel.45 ten=45 smo=1. npmin=350 dmin=4 slo=ssl.45 asp=sas.45 pcu=spc.45 tcu=stc.45 |fi=15

without gully
s.surf.rst -d survey.nogul elev=sel.80ng ten=80 smo=0.1 npmin=300 dmin=2 slo=sdx.80ng asp=sdy.80ng segmax=700| fi=27
s.surf.rst  survey.nogul elev=sel.80ng ten=80 smo=0.1 npmin=300 dmin=2 slo=ssl.80ng asp=sas.80ng pcu=spc.80ng tcu=stc.80ng segmax=700| fi=27
s.surf.rst -d gullyall.sites elev=sel.180g ten=180 smo=0.1 npmin=300 dmin=1 slo=sdx.180g asp=sdy.180g segmax=700 mask=gull.mask| fi=27
s.surf.rst  gullyall.sites elev=sel.100g ten=100 smo=0.5 npmin=300 dmin=2 segmax=700 mask=gull.mask| fi=18

sel.gf=if(gull.mask0,sel.180g,sel.80ng)
==========================================================================================
cover
s.surf.rst ifsar.img6K elev=ifsar.imgsm ten=60 smo=5 dmin=5
===============================================================================
s.surf.rst hoodel.site elev=elcs ten=150 smo=0.1 npmin=350 slo=elcs.sl asp=elcs.
as pcu=elcs.pc tcu=elcs.tc
s.surf.rst hoodel.site elev=elcs3 ten=1500 smo=0.05 npmin=300
s.surf.rst hoodn.site elev=elcs3n ten=1500 smo=0.05 npmin=300
s.surf.rst hoodel.site elev=elcs4 ten=300 smo=0.1 npmin=350 dmin=4
s.surf.rst hoodel.site elev=elcs5 ten=200 smo=1 npmin=350 dmin=5
s.surf.rst hoodn.site elev=elcs6n ten=200 smo=0.5 npmin=350 dmin=10
s.surf.rst hoodn.site elev=elcs7n ten=100 smo=0.5 npmin=350 dmin=10
s.surf.rst hoodn.site elev=elcs8n ten=200 smo=0.5 npmin=350 dmin=3| dnorm = 212, fi=42.4
s.surf.rst hoodn.site elev=elcs9n ten=l60 smo=1 npmin=400 dmin=15| dnorm = 212, fi=21.4
s.surf.rst hoodn.site elev=elcs1n ten=1000 smo=0.01 npmin=350 dmin=3| dnorm = 21 2, fi=212
s.surf.rst hoodn.site elev=elcs60 ten=60 smo=2.5 npmin=350 dmin=3 | dnorm=212,fi=12
           slo=slcs.60 asp=ascs.60 pcu=pccs.60 tcu=tccs.60
s.surf.rst hoodn.site elev=elcs70 ten=70 smo=0.1 npmin=200 dmin=1 dnorm=138,fi=9
.....................................................................

s.surf.rst elcsifsar1.8m elev=elcs.ifsm ten=60 smo=1 npmin=300 dmin=10 |dnorm = 318 
s.surf.rst elcsifsar1.8m elev=elcs.ifsm ten=45 smo=1 npmin=400 dmin=10 slo=slcs.ifs asp=ascs.ifs pcu=pccs.ifs tcu=tccs.ifs |dnorm = 366 
s.surf.rst elcsifsar1.8m elev=elcs.ifsm2m ten=50 smo=0.9 npmin=380 dmin=10 |dnorm = 342 
s.surf.rst elcsifsar1.8m elev=elcs.ifsm2m slo=slcs.ifsm2 asp=ascs.ifsm2 pcu=pccs.ifsm2 tcu=tccs.ifsm2 ten=60 smo=1.0 npmin=380 dmin=10 |dnorm = 357 
s.surf.rst elcsifsaf.5K elev=elcs.if5k slo=slcs.if5k asp=ascs.if5k pcu=pccs.if5k tcu=tccs.if5k ten=60 smo=0.5 npmin=350 dmin=10 |dnorm = 250 

two subwatersheds
s.surf.rst elcsw.ifs30K elev=elcsw.ifs30K slo=slcsw.ifs30K asp=ascsw.ifs30K pcu=pccsw.ifs30K tcu=tccsw.ifs30K ten=60 smo=0.3 npmin=350 dmin=10 | dnorm = 258.2
s.surf.rst elcsifsaf.5K elev=elcs.if5k40 slo=slcs.if5k40 asp=ascs.if5k40 pcu=pccs.if5k40
 tcu=tccs.if5k40 ten=40 smo=0.8 npmin=350 dmin=10 |dnorm = 250

erosion

lsf06=1.6*exp(dsd.9n*1,0.6)*exp(sin(slcs9n),1.3)
lsr06=1.6*exp(dsd.ifsma*2.,0.6)*exp(sin(slcs.ifs),1.3)
lsr06fp=1.6*exp(dsd.if5k*2./22.13,0.6)*exp(sin(slcs.if5k/0.0896),1.3)
w2.ls06 = exp(w2.dsd*2./22.13,0.6) * exp(sin(w2.sl)/0.0896,1.3)
huslefp=lsr06fp*kf*125*0.01
lsf16=exp(dsd.9n,1.6)*exp(sin(slcs9n),1.3)
lsr16=exp(dsd.ifsma*2.,1.6)*exp(sin(slcs.ifs),1.3)
lsr16sm2=exp(dsd.ifsm2*2.,1.6)*exp(sin(slcs.ifsm2),1.3)
lsr16fp=exp(dsd.if5k*2.,1.6)*exp(sin(slcs.if5k),1.3)
husle16fp=kf*125*0.01*exp(dsd.if5k*2.,1.6)*exp(sin(slcs.if5k),1.3)
ls1.40=dsd.if5k40*2*sin(slcs.if5k40)

r.slope.aspect ls1.40 dx=ls1.40.dx dy=ls1.40.dy
r.slope.aspect lsr16 dx=lsr16.dx dy=lsr16.dy
r.slope.aspect lsr16sm2 dx=lsr16sm2.dx dy=lsr16sm2.dy
r.slope.aspect lsf16 dx=lsf16.dx dy=lsf16.dy
r.slope.aspect lsr16fp dx=lsr16fp.dx dy=lsr16fp.dy

erdiv.40=ls1.40.dx*cos(ascs.if5k40)+ls1.40.dy*sin(ascs.if5k40)-dsd.if5k40*2*tccs.if5k40
erdiv.tr=lsr16.dx*cos(ascs.ifs)+lsr16.dy*sin(ascs.ifs)-dsd.ifsma*2*tccs.ifs
erdivsm2.tr=lsr16sm2.dx*cos(ascs.ifsm2)+lsr16sm2.dy*sin(ascs.ifsm2)-dsd.ifsm2*2*tccs.ifsm2
erdiv.tf=lsf16.dx*cos(ascs9n)+lsf16.dy*sin(ascs9n)-dsd.9n*tccs9n
erdiv.trfp=lsr16fp.dx*cos(ascs.if5k)+lsr16fp.dy*sin(ascs.if5k)-dsd.if5k*tccs.if5k

do they have slope field data?

how is imagery for C-factor - green = 0.01 short grass
                              brown = 0.5
                              CERL estimate = 0.04

R = 160 Slidel and Topsey K=0.32 Brackett K=0.17

lsr06fp=1.6*exp(dsd.if5k*2./22.13,0.6)*exp(sin(slcs.if5k/0.0896),1.3)
huslefp=kf*160*0.04*1.6*exp(dsd.if5k*2./22.13,0.6)*exp(sin(slcs.if5k/0.0896),1.3)
qs16x=kf*160*0.04*exp(dsd.if5k*2.,1.6)*exp(sin(slcs.if5k),1.3)*cos(ascs.if5k)
qs16y=kf*160*0.04*exp(dsd.if5k*2.,1.6)*exp(sin(slcs.if5k),1.3)*sin(ascs.if5k)

r.slope.aspect qs16x dx=qs16x.dx
r.slope.aspect qs16y dy=qs16y.dy

usped=qs16x.dx+qs16y.dy

R=160    C=0.04     K=0.32  - 0.17

-------------------------------------------------------------------------------
new survey
--------------------------------------------------------------------------------
R = 160 Slidel and Topsey K=0.32 Brackett K=0.17
R=160    C=0.04     K=0.32  - 0.17
imagery: green = 0.04 short grass(60%cover), brown = 0.1 bare(30%cover)
CERL C=0.04 Cf with trail =0.2/0.04 (20/60% cover)

r.flow sel.80 dsout=sdsd.80m
r.flow -u sel.80 lgout=slgu.80m

USLE LG
suslelg=kf*-160*cf.trail1*1.6*exp(slgu.80m/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelgi=kf*-160*cf.ifs*1.6*exp(slgu.80/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelg=0.32*160*0.04*1.6*exp(slgu.80/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)

USLE mod
shusle=kf*-160*cf.trail1*1.6*exp(sdsd.80m*2./22.13,0.6)*exp(sin(ssl.80/0.0896),1.6)
shusles2=kf*-160*cf.trail1*1.6*exp(sdsd.80*2./22.13,0.6)*exp(sin(ssl.80/0.0896),2)
shusles3=kf*-160*cf.trail1*1.6*exp(sdsd.80*2./22.13,0.5)*exp(sin(ssl.80/0.0896),2.5)
shuslesci=kf*-160*cf.ifs*1.6*exp(sdsd.80*2./22.13,0.6)*exp(sin(ssl.80/0.0896),1.3)

USPED p=0.6, 1,0, 1.6

p=0.5 sin**2.8 kf,cf.ifs
sqs05x=kf*160*cf.trail1*exp(sdsd.80m*2.,0.5)*exp(sin(ssl.80),2.8)*cos(sas.80)
sqs05y=kf*160*cf.trail1*exp(sdsd.80m*2.,0.5)*exp(sin(ssl.80),2.8)*sin(sas.80)
r.slope.aspect sqs05x dx=sqs05x.dx
r.slope.aspect sqs05y dy=sqs05y.dy
susped05cs28=(sqs05x.dx+sqs05y.dy)*500.
r.neighbors susped05cs28 out=susped05cs28.avg meth=average siz=3

p=0  kf,cf.trail1
sqssx=kf*160*cf.trail1*1.*exp(sin(ssl.80),1.5)*cos(sas.80)
sqssy=kf*160*cf.trail1*1.*exp(sin(ssl.80),1.5)*sin(sas.80)
r.slope.aspect sqssx dx=sqssx.dx
r.slope.aspect sqssy dy=sqssy.dy
suspeds=(sqssx.dx+sqssy.dy)
r.neighbors susped05c out=susped05c.avg meth=average siz=3


p=0.6
sqs06x=0.32*160*0.04*exp(sdsd.80*2.,0.6)*exp(sin(ssl.80),1.3)*cos(sas.80)
sqs06y=0.32*160*0.04*exp(sdsd.80*2.,0.6)*exp(sin(ssl.80),1.3)*sin(sas.80)
r.slope.aspect sqs06x dx=sqs06x.dx
r.slope.aspect sqs06y dy=sqs06y.dy
susped06=sqs06x.dx+sqs06y.dy

p=0.6 kf,cf
sqs06x=kf*160*cf.trail1*exp(sdsd.80m*2.,0.6)*exp(sin(ssl.80),1.3)*cos(sas.80)
sqs06y=kf*160*cf.trail1*exp(sdsd.80m*2.,0.6)*exp(sin(ssl.80),1.3)*sin(sas.80)
r.slope.aspect sqs06x dx=sqs06x.dx
r.slope.aspect sqs06y dy=sqs06y.dy
susped06c.200=(sqs06x.dx+sqs06y.dy)*200.
r.neighbors susped06c.200 out=susped06c.avg meth=average siz=3

p=0.6 sin=1.6 kf,cf
sqs06x=kf*160*cf.trail1*exp(sdsd.80m*2.,0.6)*exp(sin(ssl.80),1.6)*cos(sas.80)
sqs06y=kf*160*cf.trail1*exp(sdsd.80m*2.,0.6)*exp(sin(ssl.80),1.6)*sin(sas.80)
r.slope.aspect sqs06x dx=sqs06x.dx
r.slope.aspect sqs06y dy=sqs06y.dy
susped06cs16=(sqs06x.dx+sqs06y.dy)*200.
r.neighbors susped06cs16 out=susped06cs16.avg meth=average siz=3


p=0.6 kf,cfi
sqs06x=kf*160*cf.ifs*exp(sdsd.80*2.,0.6)*exp(sin(ssl.80),1.3)*cos(sas.80)
sqs06y=kf*160*cf.ifs*exp(sdsd.80*2.,0.6)*exp(sin(ssl.80),1.3)*sin(sas.80)
r.slope.aspect sqs06x dx=sqs06x.dx
r.slope.aspect sqs06y dy=sqs06y.dy
susped06ci=(sqs06x.dx+sqs06y.dy)*200.
r.neighbors susped06ci out=susped06ci.avg meth=average siz=3

p=0.6 kf,cf fi=55
sqs06x.55=kf*160*cf.trail1*exp(sdsd.55*2.,0.6)*exp(sin(ssl.55),1.3)*cos(sas.55)
sqs06y.55=kf*160*cf.trail1*exp(sdsd.55*2.,0.6)*exp(sin(ssl.55),1.3)*sin(sas.55)
r.slope.aspect sqs06x.55 dx=sqs06x55.dx
r.slope.aspect sqs06y.55 dy=sqs06y55.dy
susped06c55=(sqs06x55.dx+sqs06y55.dy)*200
r.neighbors susped06c55 out=susped06c55.avg meth=average siz=3

p=1
sqs1x=0.32*160*0.04*sdsd.80*2.*sin(ssl.80)*cos(sas.80)
sqs1y=0.32*160*0.04*sdsd.80*2.*sin(ssl.80)*sin(sas.80)
r.slope.aspect sqs1x dx=sqs1x.dx
r.slope.aspect sqs1y dy=sqs1y.dy
susped1=(sqs1x.dx+sqs1y.dy)*10.
r.colors susped1 rast=susped
r.neighbors susped1.10 out=susped1.10avg meth=average siz=3
r.neighbors susped1.50 out=susped1.50avg meth=average siz=3

p=1, Kf, Cf
sqs1x=kf*160*cf.trail1*sdsd.80m*2.*sin(ssl.80)*cos(sas.80)
sqs1y=kf*160*cf.trail1*sdsd.80m*2.*sin(ssl.80)*sin(sas.80)
r.slope.aspect sqs1x dx=sqs1x.dx
r.slope.aspect sqs1y dy=sqs1y.dy
susped1c=(sqs1x.dx+sqs1y.dy)*10.
susped1c.20=(sqs1x.dx+sqs1y.dy)*20.
r.colors susped1c rast=susped
r.neighbors susped1c out=susped1.cavg meth=average siz=3
r.neighbors susped1c.20 out=susped1.20cavg meth=average siz=3

NO GULLY p=1
sqs1xng=0.32*160*0.04*sdsd.80ng*2.*sin(ssl.80ng)*cos(sas.80ng)
sqs1yng=0.32*160*0.04*sdsd.80ng*2.*sin(ssl.80ng)*sin(sas.80ng)
r.slope.aspect sqs1xng dx=sqs1xng.dx
r.slope.aspect sqs1yng dy=sqs1yng.dy
susped1ng=(sqs1xng.dx+sqs1yng.dy)*10.
r.colors susped1ng rast=susped
r.neighbors susped1ng5 out=susped1ng5avg3 meth=average siz=3
r.neighbors susped1ng5 out=susped1ng5avg5 meth=average siz=5

NO GULLY p=1 CF,KF
sqs1xngc=kf*160*cf.trail1*sdsd.80ng*2.*sin(ssl.80ng)*cos(sas.80ng)
sqs1yngc=kf*160*cf.trail1*sdsd.80ng*2.*sin(ssl.80ng)*sin(sas.80ng)
r.slope.aspect sqs1xngc dx=sqs1xngc.dx
r.slope.aspect sqs1yngc dy=sqs1yngc.dy
susped1ngc.20=(sqs1xngc.dx+sqs1yngc.dy)*20.
r.colors susped1ngc rast=susped
r.colors susped1ngc.20 rast=susped
r.neighbors susped1ngc out=susped1ngca3 meth=average siz=3
r.neighbors susped1ngc.20 out=susped1ngca.20 meth=average siz=3

p=1.6
sqs16x=0.32*160*0.04*exp(sdsd.80*2.,1.6)*exp(sin(ssl.80),1.3)*cos(sas.80)
sqs16y=0.32*160*0.04*exp(sdsd.80*2.,1.6)*exp(sin(ssl.80),1.3)*sin(sas.80)
r.slope.aspect sqs16x dx=sqs16x.dx
r.slope.aspect sqs16y dy=sqs16y.dy
susped=sqs16x.dx+sqs16y.dy
r.neighbors susped out=susped.avg5 meth=average size=5

p=1.6 kf,cf
sqs16x=kf*160*cf.trail1*exp(sdsd.80m*2.,1.6)*exp(sin(ssl.80),1.3)*cos(sas.80)
sqs16y=kf*160*cf.trail1*exp(sdsd.80m*2.,1.6)*exp(sin(ssl.80),1.3)*sin(sas.80)
r.slope.aspect sqs16x dx=sqs16x.dx
r.slope.aspect sqs16y dy=sqs16y.dy
susped16c=sqs16x.dx+sqs16y.dy
r.neighbors susped16c out=susped16c.avg meth=average size=3

p=1.6 s=1.6 kf,cf
sqs16xs=kf*160*cf.trail1*exp(sdsd.80m*2.,1.6)*exp(sin(ssl.80),1.6)*cos(sas.80)
sqs16ys=kf*160*cf.trail1*exp(sdsd.80m*2.,1.6)*exp(sin(ssl.80),1.6)*sin(sas.80)
r.slope.aspect sqs16xs dx=sqs16xs.dx
r.slope.aspect sqs16ys dy=sqs16ys.dy
susped16cs=sqs16xs.dx+sqs16ys.dy
r.neighbors susped16cs out=susped16cs.avg meth=average size=3

p=1.6 kf,cf fi=55
sqs16x=kf*160*cf.trail1*exp(sdsd.55*2.,1.6)*exp(sin(ssl.55),1.3)*cos(sas.55)
sqs16y=kf*160*cf.trail1*exp(sdsd.55*2.,1.6)*exp(sin(ssl.55),1.3)*sin(sas.55)
r.slope.aspect sqs16x dx=sqs16x.dx
r.slope.aspect sqs16y dy=sqs16y.dy
susped16c.55=sqs16x.dx+sqs16y.dy
r.neighbors susped16c.55 out=susped1655.avg meth=average size=3

p=0.6 fi=45 too much smoothed try with p=1
sqs06x45=0.32*160*0.04*exp(sdsd.45*2.,0.6)*exp(sin(ssl.45),1.3)*cos(sas.45)
sqs06y45=0.32*160*0.04*exp(sdsd.45*2.,0.6)*exp(sin(ssl.45),1.3)*sin(sas.45)
r.slope.aspect sqs06x45 dx=sqs06x.dx45
r.slope.aspect sqs06y45 dy=sqs06y.dy45
susped0645.500=500*(sqs06x.dx45+sqs06y.dy45)

p=1 fi=55 kf,cf
sqs1x.55=kf*160*cf.trail1*sdsd.55*2.*sin(ssl.55)*cos(sas.55)
sqs1y.55=kf*160*cf.trail1*sdsd.55*2.*sin(ssl.55)*sin(sas.55)
r.slope.aspect sqs1x.55 dx=sqs1x55.dx
r.slope.aspect sqs1y.55 dy=sqs1y55.dy
susped1c=(sqs1x.dx+sqs1y.dy)*10.
susped1c55.20=(sqs1x55.dx+sqs1y55.dy)*20.
r.colors susped1c55.20 rast=susped
r.neighbors susped1c out=susped1.cavg meth=average siz=3
r.neighbors susped1c55.20 out=susped155.20cavg meth=average siz=3

averages from set of events

suspedav=(susped06c+susped06c+susped1c.20+susped16c)/4
suspedav=(susped06c.200+susped06c.200+susped1c.20+susped16c)/4
suspedav3=(susped06c+susped1c.20+susped16c)/3
r.neighbors suspedav out=suspedav.avg meth=average siz=3
r.neighbors suspedav3 out=suspedav3.avg meth=average siz=3
suspedav55=(susped06c55+susped06c55+susped1c55.20+susped16c.55)/4

SIMWE
urob ako vznika gully - aky velky dazd musi byt, zacni z vyhladeneho - zasypany gully

OLD USLE

r.average base=kf cover=slgu.80 out=slgav
slg1=if(soilbnd.matt,slgu.80)
slg2=if(soilbnd.matt,0,slgu.80)
r.sum slg1 = 2,816,212 (zle 3,501370) (20039cells)
average lg is  140  (zle 117)
r.sum slg2 = 828495 (zle 2,156000 (11922cells)
average lg is 69.
oldslgav=if(soilbnd.matt,140,69)
suslelgold=kf*-160*cf.ifs*1.6*exp(oldslgav/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelgold=kf*160*0.04*1.6*exp(oldslgav/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelgold=kf*-160*cf.trail1*1.6*exp(oldslgav/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelgco=kf*-160*cf.trail1*1.6*exp(100/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelgcos2=kf*-160*cf.trail1*1.6*exp(140/22.13,0.6)*100*sin(ssl.80)*sin(ssl.80)
suslelgco90=kf*-160*cf.trail1*1.6*exp(90/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)
suslelghom=0.24*-160*0.12*1.6*exp(140/22.13,0.6)*(65.4*sin(ssl.80)*sin(ssl.80)+4.56*sin(ssl.80)+0.0654)

----------------------------------------------------------------------------------------------------
Vysledky

with kf,cf
suslelgold suslelg shusle susped1.20cavg susped1ngca.20 susped1c55.avg susped16c.avg susped1655.avg
susped06c.avg susped06c55.avg suspedav.avg suspedav55
Opravene
suslelgold suslelg shusle
susped06c.avg,susped06cs16.avg,susped1.20cavg
susped16c.avg,susped16cs.10avg,suspedav.avg
...........................
without kf,cf
susped06.500 (ten=80, p=0.6, avg=3)
susped1.50avg (ten=80, p=1, avg=3)
susped.avg3
susped.avg3.5

RUN WITH MASK!!! 
r.what suslelgneg,susped06.500,susped06.500av,susped.avg3,susped.avg5 null=0< xycs.dat > smodels1.rep
r.what suslelgold,suslelg,shuslen,susped06.500av,susped1.50avg,susped.avg3,susped.avg3.5 null=0< xycsid.dat > smodels2.rep
with kf,cf
r.what suslelgold,suslelg,shusle,shusles2 null=0< xycsid.datma > sresult5.rep
r.what susped1.20cavg,susped1ngca.20,susped155.20cavg null=0 < xycsid.datma > sresult2.rep
r.what susped16c.avg,susped1655.avg,susped06c.avg,susped06c55.avg null=0 < xycsid.datma > sresult3.rep
r.what suspedav.avg,suspedav55,suspedav3.avg null=0< xycsid.datma > sresultav.rep
r.what susped0645.avg,susped.avg5,susped06.500av  null=0< xycsid.datma > sresulter.rep sresultopo.rep
r.what susped05c.avg,susped06ci.avg,susped06cs16.avg  null=0< xycsid.datma > sresulter.rep

New rwhat with mask and fixed flow (51 bez gully a rd30, 52 bez gully)
r.what suslelgold,suslelg,shusle< xycsid.datma > sresultusl.rep
r.what suslelgcos2,suslelgco90,suslelghom< xycsid.datma > sresultuslhom.rep
r.what susped06c.avg,susped06cs16.avg,susped1.20cavg< xycsid.datma > sresultusp06.rep
r.what susped16c.avg,susped16cs.10avg,suspedav.avg< xycsid.datma > sresultusp16.rep
r.what suslelgco,suspeds,susped05c.avg< xycsid.datma > sresultwco.rep
r.what susped05cs28,sdifu.80,sdifu.45< xycsid.datma > sresultdif.rep
r.what ssl.80neg,sdifu.80,sdifu.45< xycsid.datma > sresultsl.rep

r.what input=huslefp1mi,usped16n1,usped1c1.10,slcs.if5k,dsd.if5k<xycsid.dat | awk -F"|" '{printf "%s %s %s %.2f %.2f %.2f %.2f %d\n", $1, $2, $3, $4, $5, $6, $7, $8}' >models3.rep2

r.what input=huslefp1mi,usped16n1,csedcra.cvavg,slcs.if5k,dsd.if5k<xycsid.dat | awk -F"|" '{printf "%s %s %s %.2f %.2f %.2f %.2f %d\n", $1, $2, $3, $4, $5, $6, $7, $8}' >models3.rep5
write sitecalc to compute the statistics (with and without the gully)
ADD / MASK gully
and don't forget reverse y x
put into dig_ascii dir

v.in.ascii gull.ascii out=gully
v.support
v.to.rast gully out=gully.mask

r.support change NULL to 0

r.mapcalc
selfault=if(gylly.mask,elgul,elnogul)

rms=(cs60-suspedav.avg)*(cs60-suspedav.avg)
r.sum rms
-------------------------------------------------------

Later for publication:

Rerun with SIMWE with and without a gully and test prevention measures

Questions to Steve:

Is there any way to check in field by looking at the soil profile,
trunks/roots of bushes (if there are any) what the erosion is or
might have been on the ridges (is there realy 15t/ay as it comes out
from the Cs?) and whether there is depositied material where we have
the deposition from Cs?
The Cs measurements (these and many others, including the exeter field
and the Figure see....) have pattern of very local movement (that is why
tillage is so often used as explanation, but we don't have any tillage here)
and the pattern is close to the one which we are geting from simulations
of small and/or very short events (p=0.6, or SIMWE for short time -
see also Exeter, but I will do it for this field too). Then the question
is what is happening with Cs during the large events?
Issue of different sampling density - 100 points versus 1400
(from our experince from precision farming it makes a big difference
in correlation) New methodology is needed for comparing
spatiual distributions (surfaces, biv. functions) sampled/modeled
with different level of detail(sampling density, resolution,scale)

The concavities were lower and bigger and they are filling up

we are trying to use simulations not just to predict, but also 
to understand and explainwhat is going on.
The simulations indicate that Cs pattern is closer to pattern 
which result fro "short range" processes that to "long range" ones.
In our case the pattern is closer to the pattern produced by many
short/small rainfall events (before the water reaches steady state)
when erosion/depostion is more evenly distributed and maximums are
on the ridges, than the big events when erosion on ridges can
be susbstantial but it is still low relative to erosion in the 
center of the valleys and lower part of hillslopes.
The effect of small rainfall is similar to the effect of tillage
which is being used as explanation of Cs distribution in tilled
areas.
Because we know (do we?) that most of the soil is being moved
during the large events it seems that the linkage between Cs
and soil erosion by water is rather complex and there is something
else going on besides Cs moving with soil.

COMBINED/averaged many small (200), few middle (20) and one big event
are starting to get closer - WEPP and other continuous simulations
models do it so it seems that it is really necessary to get
the right pattern but setting the weights correctly (based on
temporal changes in rainfall and c) is really tricky

most of the continuous simulations are lumped and they did not
look at what impact it has on spatial pattern

In the WRR paper the p=0.6 gave the closest results to observed 
colluvial deposits. Distributed field measurements indicate
that the long term effects produce different patterns than
the typical rill erosion model ...
adding interril/detachment by raindrops which is idependent from
water flow (I tried it once it was very small, I will look
at it again) or diffusive term reflecting combined impact of various 
processes (zvetravanie, wind, raindrops, ...) should be bigger
or dynamics in terrain (however subtle) in long term is important
(e.g. as in hedges) - see MIT models???

Cs moves without water much more than with water
(weathering, wind, raindrop splash, ...)