# R-Code zu Blatt 2 Aufgabe 2b)

# Funktion zur Darstellung der 2-dimensionalen NV-Dichte mit Varianzen
# \sigma_1, \sigma_2 und Korrelation \rho (Mittelwert \mu=(0,0))

nv2.plot <- function(sigma1,sigma2,rho){

  library(mvtnorm)

  #Definitition von Erwartungswert (0,0) und Kovarianzmatrix
  mu <- rep(0,2)
  Sigma <- matrix(c(sigma1^2,rho*sigma1*sigma2,rho*sigma1*sigma2,sigma2^2),ncol=2,nrow=2)

  #Definition eines Grids von Werten zur Berechnung der Dichte
  x1 <- x2 <- seq(-5,5,0.25)
  grid <- expand.grid(x1,x2)

  fvals <- matrix(dmvnorm(grid,mean=mu,sigma=Sigma),ncol=length(x1),nrow=length(x2))

  #3-dimensionaler Plot der Dichte
  persp(x=x1,y=x2,z=fvals,ticktype="detailed",xlab="x_1",ylab="x_2",phi=35,theta=120)
  title(main=paste("sigma1=",as.character(sigma1), ", sigma2=",as.character(sigma2), ", rho=",as.character(rho)))

  #Plot der Höhenlinien
  contour(x=x1,y=x2,z=fvals,xlab="x_1",ylab="x_2")
  title(main=paste("sigma1=",as.character(sigma1), ", sigma2=",as.character(sigma2), ", rho=",as.character(rho)))
}

par(mfrow=c(2,1))
nv2.plot(sigma1=1,sigma2=3,rho=0.1)

par(mfrow=c(2,1))
nv2.plot(sigma1=1,sigma2=3,rho=0.5)

par(mfrow=c(2,1))
nv2.plot(sigma1=1,sigma2=3,rho=0.9)