with db

sim/gen-SZ.R

+library(ggplot2)
+library(R4CouchDB)
+
+cdb <- cdbIni()
+cdb$DBName <- "sz_mc"
+for(randSdRel in c(1e-3, 3e-3, 5e-3, 9e-3, 1e-2, 3e-2, 5e-2, 9e-2,1e-1, 3e-1, 5e-1, 9e-1)){
+for( k in 1:100){
+
+### --------------------------------------------------------
+    d <- list(
+        M              = 100,        # Anz. Messpunkte pro SZ
+        N              = 5,          # Anz. der SZ
+        dyn            = 1.6e-3,     # Dynamikbereich des SZ
+        p = list(
+            fill        = 44.6,      # mbar
+            randSdRel   = randSdRel   ),# 
+        m = list(
+            mean         =  -1.8e-3),# mbar/s
+        t = list(
+            between      = 5.0,      # s
+            randSdRel    = 0)        #
+        )
+### -------------------------------------------------------
+
+    d$p$meanSz       <- rep(0.0, d$N) # Druckmittelwert der SZ
+    d$p$typeSz       <- "const" # Druckmittelwert der SZ
+
+    d$m$sz           <- rep(d$m$mean, d$N) # Ausgangssteigung der SZ
+    idx              <- seq(1, d$M)
+
+    d$p$D            <- d$dyn * d$p$fill
+    d$p$d            <- -d$p$D / d$M
+
+    d$p$upper        <- d$p$meanSz + d$p$D / 2
+    d$p$lower        <- d$p$meanSz - d$p$D / 2
+    d$t$d            <- d$p$d / d$m$sz
+
+    d$t$randSd       <- abs(d$t$d * d$t$randSdRel)
+    d$p$randSd       <- abs(d$p$D * d$p$randSdRel)
+
+    ## neue knstliche SZs
+    ## erzeugen
+
+    d$p$pure         <- matrix(NA, ncol = d$M , nrow = d$N)
+    d$p$rand         <- matrix(NA, ncol = d$M , nrow = d$N)
+
+    d$t$pure         <- matrix(NA, ncol = d$M , nrow = d$N)
+    d$t$rand         <- matrix(NA, ncol = d$M , nrow = d$N)
+
+    d$m$rand         <- rep(NA, d$N)
+    d$c$rand         <- rep(NA, d$N)
+
+    d$p$randMean     <- rep(NA, d$N)
+    d$t$randMean     <- rep(NA, d$N)
+
+    for(i in 1:d$N){
+
+        if(i == 1){
+            d$t$pure[i,]  <- d$t$d[i] * idx
+        }else{
+            d$t$pure[i,]  <- d$t$pure[i-1, d$M] + d$t$between + d$t$d[i] * idx
+        }
+
+        d$p$pure[i,]      <- seq(d$p$upper[i],
+                             d$p$lower[i],
+                             length.out = d$M)
+
+        d$p$rand[i,]      <- d$p$pure[i,] + rnorm(d$M, 0, d$p$randSd)
+        d$t$rand[i,]      <- d$t$pure[i,] + rnorm(d$M, 0, d$t$randSd)
+
+        temp.lm           <- lm( d$p$rand[i,] ~ d$t$rand[i,])
+        d$m$rand[i]       <- coefficients(temp.lm)[2]
+        d$p$randMean[i]   <- mean(d$p$rand[i,])
+        d$t$randMean[i]   <- mean(d$t$rand[i,])
+    }
+
+    d$p$randIntercept     <- d$p$randMean - d$m$rand * d$t$randMean
+    d$t$rand0             <- (mean(d$p$randMean) - d$p$randIntercept) / d$m$rand
+
+    d$t$DeltaSd           <- sd(diff(d$t$rand0))
+    d$t$DeltaMean         <- mean(diff(d$t$rand0))
+    d$t$DeltaSdRel        <- sd(diff(d$t$rand0)) / d$t$DeltaMean
+
+
+    cdb$dataList <- d
+    dn <- cdbAddDoc(cdb)
+}
+
+}
+pt.df <- data.frame(d.p.pure  = as.vector(d$p$pure),
+                    d.t.pure  = as.vector(d$t$pure),
+                    d.p.rand = as.vector(d$p$rand),
+                    d.t.rand = as.vector(d$t$rand))
+ggplot(pt.df) +
+    geom_point(aes(x = d.t.pure, y = d.p.pure), color='black') +
+    geom_point(aes(x = d.t.rand, y = d.p.rand), color='red') +
+    geom_abline(aes(intercept = d$p$randIntercept, slope = d$m$rand), color='red') +
+    geom_vline(xintercept = d$t$rand0, linetype = "longdash") +
+    geom_point(aes(y = d$p$randMean, x = d$t$randMean), color="black", size=4) +
+    ggtitle("dp vs. dt")
+
+

sim/plot-SZ.R

+
+cdb <- cdbInit()
+cdb$DBName <- "sz_mc"
+cdb$id <- "1f817eefa918f5bb1a8620eaa01f9aac"
+
+
+d <- cdbGetDoc(cdb)$res
+
+
+pt.df <- data.frame(d.p.pure  = unlist(d$p$pure),
+                    d.t.pure  = unlist(d$t$pure),
+                    d.p.rand  = unlist(d$p$rand),
+                    d.t.rand  = unlist(d$t$rand))
+ggplot(pt.df) +
+    geom_point(aes(x = d.t.pure, y = d.p.pure), color='black') +
+    geom_point(aes(x = d.t.rand, y = d.p.rand), color='red') +
+    geom_abline(aes(intercept = unlist(d$p$randIntercept), slope = unlist(d$m$rand)), color='red') +
+    geom_vline(xintercept = unlist(d$t$rand0), linetype = "longdash") +
+    geom_point(aes(y = unlist(d$p$randMean), x = unlist(d$t$randMean)), color="black", size=4) +
+    ggtitle("dp vs. dt")

sim/theo.R

-library(ggplot2)
-library(R4CouchDB)
-
-### --------------------------------------------------------
-d <- list(
-  M              = 100,       # Anz. Messpunkte pro SZ
-  N              = 5,         # Anz. der SZ
-  dyn            = 1.6e-3,    # dynamikbereich des SZ
-  p = list(                   
-    fill        = 44.6,       # mbar
-    randSdRel   = 0   ),          #1e-8
-      # relativer hub der SZ
-  m = list(
-    mean         =  -1.8e-3), # mbar/s
-  t = list(
-    between      = 5.0,       # s
-    randSdRel    = 1.0e-3)    #
-  )
-### -------------------------------------------------------
-  
-  d$p$meanSz     <- rep(0.0, d$N) # Druckmittelwert der SZ
-  d$m$sz           <- rep(d$m$mean, d$N) # Ausgangssteigung der SZ
-  d$idx            <- seq(1, d$M)
-
-  d$p$D            <- d$dyn * d$p$fill
-  d$p$d            <- -d$p$D / d$M
-
-  d$p$upper        <- d$p$meanSz + d$p$D / 2
-  d$p$lower        <- d$p$meanSz - d$p$D / 2
-  d$t$d            <- d$p$d / d$m$sz
-
-  d$t$randSd      <- abs(d$t$d * d$t$randSdRel)
-  d$p$randSd      <- abs(d$p$fill * d$t$randSdRel)
-
-  ## neue knstliche SZs
-  ## erzeugen
-
-  d$p$pure  <- matrix(NA, ncol = d$M , nrow = d$N)
-  d$p$rand <- matrix(NA, ncol = d$M , nrow = d$N)
-    
-  d$t$pure  <- matrix(NA, ncol = d$M , nrow = d$N)
-  d$t$rand <- matrix(NA, ncol = d$M , nrow = d$N)
-
-  d$m$rand       <- rep(NA, d$N)
-  d$c$rand       <- rep(NA, d$N)
-
-  d$p$randMean   <- rep(NA, d$N)
-  d$t$randMean   <- rep(NA, d$N)
-
-  for(i in 1:d$N){
-
-    if(i == 1){
-      d$t$pure[i,]  <- d$t$d[i] * d$idx
-    }else{
-      d$t$pure[i,]  <- d$t$pure[i-1, d$M] + d$t$between + d$t$d[i] * d$idx
-    }
-
-    d$p$pure[i,]  <- seq(d$p$upper[i],
-                      d$p$lower[i],
-                      length.out = d$M)
- 
-    d$p$rand[i,] <- d$p$pure[i,] + rnorm(d$M, 0, d$p$randSd)
-    d$t$rand[i,] <- d$t$pure[i,] + rnorm(d$M, 0, d$t$randSd)
- 
-    temp.lm        <- lm( d$p$rand[i,] ~ d$t$rand[i,])
-    d$c$rand[i]      <- coefficients(temp.lm)[1]
-    d$m$rand[i]      <- coefficients(temp.lm)[2]
-    d$p$randMean[i] <- mean(d$p$rand[i,])
-    d$t$randMean[i] <- mean(d$t$rand[i,])
-  }
-
-  d$ci      <- d$p$randMean - d$m$rand * d$t$randMean
-  d$t$rand0 <- (mean(d$p$randMean) - d$ci)/ d$m$rand
-
-  d$t$DeltaSd <- sd(diff(d$t$rand0))
-  d$t.DeltaMean <- mean(diff(d$t$rand0))
-  d$t.DeltaSdRel <- sd(diff(d$t$rand0)) / d$t$DeltaMean
-
-### --- plot ---
-
-  pt.df <- data.frame(d.p.pure  = as.vector(d$p$pure),
-                      d.t.pure  = as.vector(d$t$pure),
-                      d.p.rand = as.vector(d$p$rand),
-                      d.t.rand = as.vector(d$t$rand))
-  ggplot(pt.df) +
-  geom_point(aes(x = d.t.pure, y = d.p.pure), color='black') +
-  geom_point(aes(x = d.t.rand, y = d.p.rand), color='red') +
-  geom_abline(aes(intercept = d$c$rand, slope = d$m$rand), color='red') +
-  geom_vline(xintercept = d$t$rand0, linetype = "longdash") +
-  geom_point(aes(y = d$p$randMean, x = d$t$randMean), color="black", size=4) +
-  ggtitle("dp vs. dt")