################## ZOBRAZENI ANTROPOLOGICKYCH DAT VE 3D
DATA3D <- read.delim("07_data3d.txt", header = TRUE)
attach(DATA3D)
str(DATA3D)
dim(DATA3D)

DATA3D_f <- DATA3D[sex == "f",]
DATA3D_m <- DATA3D[sex == "m",]

# vidime, ze mame souradnice dvanacti bodu na lebce pro 58 jedincu
# nejprve zobrazime glabellu (G) - a idealne zvlast pro zeny a zvlast pro muze
library(rgl)
library(car)
# rgl.open()
open3d()
rgl.spheres(DATA3D_f$x.G,DATA3D_f$y.G,DATA3D_f$z.G, data = DATA3D_f, r = 2, color="palevioletred1")
rgl.spheres(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G, data = DATA3D_m, r = 2, color="lightskyblue")
rgl.bbox(color="darkblue")

# pro zobrazeni muzeme vyuzit i scatter3d (na grupovana data)
scatter3d(x = DATA3D$x.G, y = DATA3D$y.G, z = DATA3D$z.G, data=DATA3D, groups = DATA3D$sex)
scatter3d(x = DATA3D$x.G, y = DATA3D$y.G, z = DATA3D$z.G, data=DATA3D, groups = DATA3D$sex, fit="smooth")

# ale mame mracno bodu, bylo by vhodne najit prumer tohoto mracna
# zvlast pro muze a zeny!
prumF_x <- mean(DATA3D_f$x.G)
prumF_y <- mean(DATA3D_f$y.G)
prumF_z <- mean(DATA3D_f$z.G)

prumM_x <- mean(DATA3D_m$x.G)
prumM_y <- mean(DATA3D_m$y.G)
prumM_z <- mean(DATA3D_m$z.G)

# rgl.open()
open3d()
rgl.spheres(DATA3D_f$x.G,DATA3D_f$y.G,DATA3D_f$z.G, data = DATA3D_f, r = 2, color="palevioletred4")
rgl.spheres(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G, data = DATA3D_m, r = 2, color="royalblue")
rgl.spheres(prumF_x,prumF_y,prumF_z, r = 4, color="violetred4")
rgl.spheres(prumM_x,prumM_y,prumM_z, r = 4, color="mediumblue")

# chceme zobrazit vsechny body
# nachystame si souradnice do objektu
DATA3D_x <- DATA3D[,c(1,2,3,6,9,12,15,18,21,24,27,30,33,36)]
DATA3D_y <- DATA3D[,c(1,2,4,7,10,13,16,19,22,25,28,31,34,37)]
DATA3D_z <- DATA3D[,c(1,2,5,8,11,14,17,20,23,26,29,32,35,38)]

# funkce rgl.spheres ale vyzaduje zadavani souradnic jako vektory 
# nacteme knihovnu balicku reshape a pomoci funkce melt "spojime" pod sebe vsechny sloupce danych objektu
# install.packages("reshape")
library(reshape)
x_s <- melt(DATA3D_x, id.vars= c("id", "sex"),var = "point")
y_s <- melt(DATA3D_y, id.vars= c("id", "sex"),var = "point")
z_s <- melt(DATA3D_z, id.vars= c("id", "sex"),var = "point")

# zobrazime vsechny body
# rgl.open()
open3d()
rgl.spheres(x_s[sex == "f","value"],y_s[sex == "f","value"],z_s[sex == "f","value"], r = 2, color="palevioletred4")
rgl.spheres(x_s[sex == "m","value"],y_s[sex == "m","value"],z_s[sex == "m","value"], r = 2, color="royalblue")

# jedna se o pomerne neprehledne pole bodu, vypocitame souradnice prumernych bodu
# zvlast pro muze a zeny
prumerF_x <- colMeans(DATA3D[sex == "f",c(3,6,9,12,15,18,21,24,27,30,33,36)])
prumerF_y <- colMeans(DATA3D[sex == "f",c(4,7,10,13,16,19,22,25,28,31,34,37)])
prumerF_z <- colMeans(DATA3D[sex == "f",c(5,8,11,14,17,20,23,26,29,32,35,38)])

prumerM_x <- colMeans(DATA3D[sex == "m",c(3,6,9,12,15,18,21,24,27,30,33,36)])
prumerM_y <- colMeans(DATA3D[sex == "m",c(4,7,10,13,16,19,22,25,28,31,34,37)])
prumerM_z <- colMeans(DATA3D[sex == "m",c(5,8,11,14,17,20,23,26,29,32,35,38)])

# rgl.open()
open3d()
rgl.spheres(x_s[sex == "f","value"],y_s[sex == "f","value"],z_s[sex == "f","value"], r = 2, color="palevioletred4")
rgl.spheres(x_s[sex == "m","value"],y_s[sex == "m","value"],z_s[sex == "m","value"], r = 2, color="royalblue")
rgl.spheres(prumerF_x,prumerF_y,prumerF_z, r = 4, color="violetred4")
rgl.spheres(prumerM_x,prumerM_y,prumerM_z, r = 4, color="mediumblue")

# muzeme zobrazit pouze prumery:
# rgl.open()
open3d()
rgl.viewpoint(theta = 0, phi = -90)
rgl.spheres(prumerF_x,prumerF_y,prumerF_z, r = 4, color="violetred4")
rgl.spheres(prumerM_x,prumerM_y,prumerM_z, r = 4, color="mediumblue")
# chceme k bodum priradit nazvy bodu
colnames(DATA3D)
# vidime, ze se jedna o body: "G","P","FMOdx","FMOsin","JUdx","JUsin","B","L","OP","BA","POdx","POsin"
rgl.texts(prumerF_x,prumerF_y,prumerF_z, text = c("G","P","FMOdx","FMOsin","JUdx","JUsin","B","L","OP","BA","POdx","POsin"), 
          adj = 2.5, color = "gray15")
# chceme spojit bod G s bodem B
rgl.lines(prumerM_x[c(1,7)],prumerM_y[c(1,7)],prumerM_z[c(1,7)], color="mediumblue", lwd = 2)
rgl.lines(prumerM_x[c(7,8)],prumerM_y[c(7,8)],prumerM_z[c(7,8)], color="mediumblue", lwd = 2)

# Priklad: souradnice bodu glabella - pridame elipsoidy pro 95% interval spolehlivosti
rgl.open()
par3d(windowRect=c(20, 30, 800, 800))
rgl.viewpoint(theta = 20, phi = 45, zoom = 1) 
rgl.spheres(DATA3D_f$x.G,DATA3D_f$y.G,DATA3D_f$z.G, data = DATA3D_f, r = 2, color="plum2")
rgl.spheres(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G, data = DATA3D_m, r = 2, color="paleturquoise2")
rgl.spheres(prumF_x,prumF_y,prumF_z, r = 4, color="palevioletred1")
rgl.spheres(prumM_x,prumM_y,prumM_z, r = 4, color="paleturquoise4")
# vypocet elipsoidu, funkce elipse3d(): 
# parametr x = korelacni nebo kovariancni matice mezi x, y a z
# parametr centre = pozice stredu elipsoidu - v nasem pripade v bode o souradnicich: prumM_x,prumM_y,prumM_z (prumery)
# parametr level = interval spolehlivost
# korelacni elipsu3d vykreslujeme pomoci finkce wire3d() nebo napr shade3d()

ellips_F <- ellipse3d(cov(cbind(DATA3D_f$x.G,DATA3D_f$y.G,DATA3D_f$z.G)), centre=c(prumF_x,prumF_y,prumF_z), level = 0.95)
wire3d(ellips_F, col = "plum2",  lit = FALSE)
ellips_M <- ellipse3d(cov(cbind(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G)), centre=c(prumM_x,prumM_y,prumM_z), level = 0.95)
wire3d(ellips_M, col = "paleturquoise2",  lit = FALSE)
# chceme ulozit nahled jako png
rgl.snapshot("G_ellipse3d.png",fmt = "png")
# chceme ulozit nahled jako pdf (note - obcas komplikace v pripade polo-pruhlednych objektu)
# v pripade komplikovanejsich objektu trva prikaz dlouho (zejm. na pomalejsich pocitacich!)
rgl.postscript("G_ellipse3d.pdf", fmt = "pdf")

# vidime jeden bod u glabelly nespada do 95% intervalu => chceme zjistit ktery, pridame text. ID
rgl.open()
par3d(windowRect=c(20, 30, 800, 800))
rgl.viewpoint(theta = 20, phi = 45, zoom = 1) 
rgl.spheres(DATA3D_f$x.G,DATA3D_f$y.G,DATA3D_f$z.G, data = DATA3D_f, r = 2, color="plum2")
rgl.spheres(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G, data = DATA3D_m, r = 2, color="paleturquoise2")
rgl.spheres(prumF_x,prumF_y,prumF_z, r = 4, color="palevioletred1")
rgl.spheres(prumM_x,prumM_y,prumM_z, r = 4, color="paleturquoise4")
rgl.texts(DATA3D_m$x.G,DATA3D_m$y.G,DATA3D_m$z.G, text = DATA3D_m$id, adj = 2.5, color = "ivory")
# vidime, ze se jedna o pripad 24, otazkou je, zda jsou spatne superponovana data, zda se jedna o chybu mereni, extremni hodnotu,
# patologickeho jedince nebo jedince z jine populace

##### DOMACI UKOL 9 ###################################
# A) Zobrazte bod bregma (B) pro jedince ze souboru data3d.opr.txt: - obdobným způsobem jako bod glabella (G),
# tzn. muže a ženy zvlášť (odlište barvou).
# Změňte typ bodů (pomocí funkce: shapelist3d)
# Vypočítejte průměry a zobrazte je poloprůhledné
# Zobrazte elipsoidy pro 50% výskytu dat
# Uložte náhled okna jako .png soubor (zvolte vhodně úhel pohledu a velikost pole)
# Do ISu vložte .R script a obrázek .png ve tvaru: DU_09_UCO_Jmeno.R/png.


######### OBRAZKY ###########################
# install.packages("imager")
library(imager)
# moznost nahrani obrazku do objektu R
im <- load.image("08_britannica.jpg")
str(im)
plot(im)
# zobrazeni rovnou pomoci plot
plot(load.image("08_britannica.jpg"))
points(427,308, pch=3, col="tomato") #Brno
# lze upravit napr. osy apod.
plot(im, axes=FALSE)
box()
axis(3, at=c(100,200,300,400,500,600), labels = c(900,1000,1100,1200,1300,1400))

# install.packages("png") 
# install.packages("jpeg")
library(jpeg)
library(png)
obrazekA <- readPNG("08_britannica.png")
windows()
plot(1:2, type='n', main="Plotting Over an Image", xlab="x", ylab="y")  # nastavení oblasti plot
lim <- par()                                                            # get the plot information so the image will fill the plot box, and draw it
rasterImage(obrazekA, lim$usr[1], lim$usr[3], lim$usr[2], lim$usr[4])
points(1.62,1.28, pch=15, cex = 2, col="tomato")                        # Brno
grid()                                                                  # mrizka

obrazekB <- readJPEG("08_britannica.jpg")
windows()
plot(1:2, type='n', main="Plotting Over an Image", xlab="x", ylab="y")  # nastavení oblasti plot
lim <- par()                                                            # get the plot information so the image will fill the plot box, and draw it
rasterImage(obrazekB, lim$usr[1], lim$usr[3], lim$usr[2], lim$usr[4])
points(1.62,1.28, pch=15, cex = 2, col="tomato")                        # Brno
grid() 
