waves.R (1877B)
1 #!/usr/bin/Rscript 2 3 arma.bits.qweibull <- function (x, param, sh) { 4 sh <- 2.0 5 sc = param$mean / gamma(1.0 + 1.0 / sh) 6 qweibull(x, shape=sh, scale=sc) 7 } 8 arma.bits.qnorm <- function (x, param) { 9 qnorm(x, mean=param$mean, sd=param$sd) 10 } 11 12 arma.QFUNCTIONS <- list( 13 elevation = arma.bits.qnorm, 14 heights_x = function (x, param) { arma.bits.qweibull(x, param, 2.0) }, 15 heights_y = function (x, param) { arma.bits.qweibull(x, param, 2.0) }, 16 lengths_x = function (x, param) { arma.bits.qweibull(x, param, 2.3) }, 17 lengths_y = function (x, param) { arma.bits.qweibull(x, param, 2.3) }, 18 periods = function (x, param) { arma.bits.qweibull(x, param, 3.0) } 19 ) 20 21 arma.ALL_WAVE_CHARACTERISTICS <- list( 22 "elevation", 23 "heights_x", 24 "heights_y", 25 "lengths_x", 26 "lengths_y", 27 "periods" 28 ) 29 30 arma.load_wave_parameters <- function ( 31 dir = "output", 32 wave_params = arma.ALL_WAVE_CHARACTERISTICS, 33 qfuncs = arma.QFUNCTIONS 34 ) { 35 sapply(wave_params, function (filename) { 36 param <- list() 37 param$filename = filename 38 param$data <- read.table(file.path(dir, filename))[[1]] 39 param$mean <- mean(param$data) 40 param$sd <- sd(param$data) 41 param$qfunc <- function (x) { qfuncs[[filename]](x, param) } 42 param$min <- param$qfunc(.01) 43 param$max <- param$qfunc(.99) 44 elem <- list() 45 elem[[filename]] <- param 46 elem 47 }) 48 } 49 50 arma.qqplot <- function (param, nsamples=100, ttl, ...) { 51 qdata <- param$qfunc(ppoints(nsamples)) 52 qqplot( 53 qdata, 54 param$data, 55 asp=1, 56 xlim=c(param$min, param$max), 57 ylim=c(param$min, param$max), 58 bty='n', 59 ... 60 ) 61 box(col='black') 62 if (class(ttl) == 'character') { 63 title(ttl, line=-2) 64 } else { 65 title(ttl$title, line=-1, adj=ttl$adjust) 66 } 67 qqline(param$data, distribution=param$qfunc) 68 } 69 70 #wave_params <- arma.load_wave_parameters() 71 #par(pty="s", mfrow=c(3, 3)) 72 #for (name in names(wave_params)) { 73 # print(name) 74 # arma.qqplot(wave_params[[name]]) 75 #}