### RX -o 'ps<-"x"' -i RFtest.all.R | R > xa.out & ### RX -o 'testlevel<-3;ps<-"x"' -i RFtest.all.R | R ## R --no-save < RFtest.all.R # source("RFtest.all.R") source("./RFtest.R") if (!exists("testlevel")) testlevel <- 1 maxtestlevel <- 5 switch (testlevel, {pointnumber<-20; valuerepet<-3; pointrepet<-1; NP<-1; duration<-"70 sec"; random.parameters <- TRUE; PrintLevel <- 1; lines<-100;repeatscript <- 1 }, {pointnumber<-50; valuerepet<-30; pointrepet<-4; NP<-1; duration<-"2.5 h"; random.parameters <- TRUE; PrintLevel <- 2; lines<-300;repeatscript <- 1 }, {pointnumber<-50; valuerepet<-200; pointrepet<-8; NP<-3; duration<-"?? d"; random.parameters <- TRUE; PrintLevel <- 3; lines<-500;repeatscript <- 1}, {pointnumber<-200; valuerepet<-300; pointrepet<-12; NP<-3; duration<-"?? w"; random.parameters <- TRUE; PrintLevel <- 3; lines<-1000;repeatscript <- 3 }, {pointnumber<-200; valuerepet<-300; pointrepet<-12; NP<-3; duration<-"?? w"; random.parameters <- FALSE; PrintLevel <- 3; lines<-1000;repeatscript <- 3 } # further test: random point number range [0,10] ) #****************************************************************** #* * #* This test tries to simulate 1, 2, and 3 dimensional random * #* fields of any defined covariance function with any defined * #* method. The random fields may be defined on a grid or not. * #* * #* For each combination, several random fields are generated * #* and the empirical variogram is calculated and (graphically) * #* compared with the true variogram. * #* In some dimensions or for some (randomly chosen) parameters * #* the variogram (thus the random field) does not exist. In this * #* case, a single red cross in the middle of the picture is * #* shown. * #* For some parameter values it may happen that none of the * #* simulation methods are suitable. In this case the true * #* variogram is shown in form of red crosses. * #* In case the random field can be simulated by a certain method * #* the estimated variogram is shown. The latter may nonetheless * #* deviate largely from the true variogram (depending on the * #* parameters and the precision chosen) * #* * #* Please ignore the strange titles and labels, which code al l * #* the relevant parameter values, and some summary statistic * #* of the simulation. * #* * #* As there are about 20 covariance functions, some of them * #* allowing for additional parameters, and 5 methods implemented * #* the number of combinations is quite large and the running * #* time can be pretty high! * #* Therefore different test levels exist. Low level means quick * #* testing, but the estimated variogram may be very eratic. * #* * #* The estimated time for running the test on a 500Mhz PC is * #* given by the variable `duration' above * #* * #******************************************************************" if (exists("ps") && !is.null(ps)) { ## for author's use only ps <- paste(ps,pid(),sep=".") save <- TRUE; } else { ## for users ps <- NULL; histo<-FALSE; save <- FALSE; } ### minor adjustments ... fieldsize<-3 endofbins <- 1.3 numberbins <-31 if (!is.null(ps)) {histo <- FALSE} ## Da R CMD check die Variable `ps' nicht ## definiert hat, if (histo) {X11(); X11();} ## wird X11() nicht mehr aufgerufen RFparameters(Storing=TRUE,PrintLevel=PrintLevel, TBM2.lines=6/5*lines,TBM2.linesimufactor=0,TBM2.linesimustep=0.01, TBM3D2.lines=lines,TBM3D2.linesimufactor=0,TBM3D2.linesimustep=0.01, TBM3D3.lines=lines,TBM3D3.linesimufactor=0,TBM3D3.linesimustep=0.01, spectral.lines=lines) ENVIR <- environment() randomize <- function(){ assign("quadraticgrid", runif(1)<0.5, envir=ENVIR) ##quadraticgrid <- FALSE assign("Mean", runif(1,0,10), envir=ENVIR) assign("scaling", if (runif(1)<0.5) runif(1,1,10) else runif(1,0.2,1), envir=ENVIR) assign("variance", if (runif(1)<0.5) runif(1,1,10) else runif(1,0.2,1), envir=ENVIR) assign("nugget", if (runif(1)<0.5) runif(1,1,10) else runif(1,0.2,1), envir=ENVIR) } simplemodels <- c("circular","cubic", "exponential","gauss", "gneiting","penta", "spherical","wave"); if (random.parameters){ models <- list( list(model="bessel", kappa1=runif(NP,-1,3), kappa2=NULL), list(model="cauchy", kappa1=runif(NP,-0.5,3), kappa2=NULL), list(model="damped", kappa1=runif(NP,-0.5,2.5), kappa2=NULL), list(model="FD", kappa1=runif(NP,-1,1), kappa2=NULL), list(model="fractg", kappa1=runif(NP,-0.5,2.5), kappa2=NULL), list(model="gencau", kappa1=runif(NP,-0.5,2.5), kappa2=runif(NP,-0.5,3)), list(model="gengn", kappa1=round(runif(NP,-0.5,4)), kappa2=runif(NP,-0.5,3)), list(model="lgd1", kappa1=runif(NP,-0.5,3.5), kappa2=runif(NP,-0.5,1.5)), list(model="power", kappa1=runif(NP,-0.5,3), kappa2=NULL), list(model="qexpone", kappa1=runif(NP,-0.2,1.2), kappa2=NULL), list(model="stable", kappa1=runif(NP,-0.5,2.5), kappa2=NULL), list(model="whittle", kappa1=runif(NP,-0.5,3), kappa2=NULL), list(model="fractalB",kappa1=runif(NP,-0.5,2.5), kappa2=NULL), ) largemodels <- list( list(model="cone", kappa1=runif(NP,0.5,1), kappa2=runif(NP,0.5,3), kappa3=runif(NP,0.5,3)), list(model="cone", kappa1=runif(NP,-0.5,2.5), kappa2=runif(NP,-0.5,3), kappa3=runif(NP,-0.5,3)), list(model="cauchytbm", kappa1=runif(NP,-0.5,2.5), kappa2=runif(NP,-0.5,3), kappa3=round(runif(NP,0,50)/10)), list(model="hyper", kappa1=runif(NP,-0.5,3), kappa2=runif(NP,-0.5,3), kappa3=runif(NP,-0.5,3)) ) ## very large models are not tested yet: ## nsst, nsst2 } else { ## fixed random parameters models <- list( list(model="bessel", kappa1=c(-1,-0.5,0,0.5,4),kappa2=NULL), list(model="cauchy", kappa1=c(-1,0,0.1,5), kappa2=NULL), list(model="damped", kappa1=c(-1,0,1,2), kappa2=NULL), list(model="FD", kappa1=c(-1,-0.5,0,0.5), kappa2=NULL), list(model="fractg", kappa1=c(0,1,2,3), kappa2=NULL), list(model="gencau", kappa1=c(-1,0,1,2), kappa2=c(-1,0,0.1,5)), list(model="gengn", kappa1=c(-1,0,0.1,1,8), kappa2=c(-1,0,0.1,4)), list(model="lgd1", kappa1=runif(0,1,2,3,4), kappa2=c(-1,0,1)), list(model="power", kappa1=c(-1,0.1,1,2), kappa2=NULL), list(model="qexpon", kappa1=c(-1,0,0.5,1), kappa2=NULL), list(model="stable", kappa1=c(-1,0,0.1,1,2), kappa2=NULL), list(model="whittle", kappa1=c(-1,0,0.1,1,8), kappa2=NULL), list(model="fractalB",kappa1=c(-1,0,1,2,3), kappa2=NULL), ) largemodels <- list( list(model="cauchytbm", kappa1=c(-1,0,0.1,1,8), kappa2=c(-1,0,0.1,5), kappa3=c(0.5,1,1.5,2,3)), list(model="cone", kappa1=c(-1,0,8), kappa2=c(-0.5,0,5), kappa3=c(-0.5,0,3)), list(model="hyper", kappa1=c(-0.1,0,0.1,5), kappa2=c(-5,-0.1,0,0.1,5), kappa3=c(-0.1,0,0.1,5)) ) ## very large models are not tested yet: ## nsst, nsst2 } for (i in 1:repeatscript) { randomize() print("nugget only simulates the nugget effect. Therefore the value of variance is not taken into account"); RFcontrol("nugget",pointnumber=pointnumber,valuerepet=valuerepet, pointrepet=pointrepet,scal=1,var=0,nug=nugget,mean=Mean, field=fieldsize,endof=endofbins,numberb=numberbins,histo=histo, ps=ps,quadraticgrid=quadraticgrid,save=save) for (naturalscaling in FALSE:TRUE) { RFparameters(PracticalRange=as.logical(naturalscaling)) for (model in models) {## randomize() RFcontrol(model$model,pointnumber=pointnumber,valuerepet=valuerepet, pointrepet=pointrepet,kappa1=model$kappa1,kappa2=model$kappa2, scaling=scaling,var=variance,nug=nugget,mean=Mean, field=fieldsize,endof=endofbins,numberb=numberbins, histo=histo,ps=ps,quadraticgrid=quadraticgrid,save=save) } for (model in largemodels) { randomize() RFcontrol(model$model,pointnumber=pointnumber,valuerepet=valuerepet, pointrepet=pointrepet,kappa1=model$kappa1,kappa2=model$kappa2, kappa3=model$kappa3,scaling=scaling,var=variance,nug=nugget, mean=Mean,field=fieldsize,endof=endofbins,numberb=numberbins, histo=histo,ps=ps,quadraticgrid=quadraticgrid,save=save) print("OK") } for (model in simplemodels) { randomize() RFcontrol(model,pointnumber=pointnumber,valuerepet=valuerepet, pointrepet=pointrepet,scaling=scaling,var=variance,nug=nugget, mean=Mean,field=fieldsize,endof=endofbins,numberb=numberbins, histo=histo,ps=ps,quadraticgrid=quadraticgrid,save=save) } } }