[R-sig-Geo] Help: R-function to plot all CSGS realization variograms in single a plot

[Zia Ahmed]

Hi,
 I am trying to plot 100 variograms from 100 CSGS realizations  in a 
single plot with the observed  variogram model. I used following codes 
to plot them  together. But problem is this code is: I have to create 
variogram for all realizations,  ie  100  times I  have to repeat this  
code!
This their  any  way  to run following codes with in a loop to create a 
plot of variograms of all realizations with observed variogram model. 
Help will  be appreciated. 
Thanks
Zia


library(gstat)
data(meuse)
data(meuse.grid)
coordinates(meuse) = ~x + y
coordinates(meuse.grid) = ~x + y

# Variogram log Zn

lzn.vgm = variogram(log(zinc) ~ 1, meuse)
lzn.fit = fit.variogram(lzn.vgm, model = vgm(1, "Sph", 900, 1))

#Conditional simulation

lzn.sim = krige(log(zinc) ~ 1, meuse, meuse.grid, model = lzn.fit,
nmax = 30, nsim = 100)

# Variogram of all relizations

var.sim1<-variogram(sim1~1, lzn.sim)    # Realization 1
var.sim2<-variogram(sim2~1, lzn.sim)    # Realization 2
var.sim3<-variogram(sim3~1, lzn.sim)    # Realization 3
var.sim4<-variogram(sim4~1, lzn.sim)    # Realization 4
.
.
.
var.sim100<-variogram(sim100~1, lzn.sim)    # Realization 100

# Vriogram of observed data and all relaizations
##################################################

plot(lzn.vgm$gamma ~ lzn.vgm$dist, xlim = c(0, max(lzn.vgm$dist) * 1.05),
ylim = c(0, max(lzn.vgm$gamma) * 1.2), pch = 21, col = "black",cex=.8,
cex.axis = 0.85, cex.lab=.8, xlab =  "Distance (m)", ylab = "Semivariance")

lines(variogramLine(fit.variogram(lzn.vgm,  lzn.fit),  # Sim1
maxdist = max(lzn.vgm$dist)), col = "red", lwd=2)

lines(variogramLine(fit.variogram(var.sim1, lzn.fit),  # sim2
maxdist = max(lzn.vgm$dist)), col = "black", lty=1)

lines(variogramLine(fit.variogram(var.sim2, lzn.fit), # Sim3
maxdist = max(lzn.vgm$dist)), col = "black", lty=1)

lines(variogramLine(fit.variogram(var.sim3, lzn.fit), # Sim4
maxdist = max(lzn.vgm$dist)), col = "black", lty=1)
.
.
.
.
lines(variogramLine(fit.variogram(var.sim100, lzn.fit), # Sim100
maxdist = max(lzn.vgm$dist)), col = "black", lty=1) 


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20100318/795b65be/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: zua3.vcf
Type: text/x-vcard
Size: 281 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20100318/795b65be/attachment.vcf>