Commit 050ff095 authored by GODARD Vincent's avatar GODARD Vincent
Browse files

reorganize folders

parent da3c984d
......@@ -8,27 +8,28 @@ source("functions/utils.R")
# Define UI for app that draws a histogram ----
ui <- fluidPage(
img(src='amu.png', align = "right", width=100),
# App title ----
titlePanel("Cosmos profile"),
# Sidebar layout with input and output definitions ----
sidebarLayout(
# Sidebar panel for inputs ----
sidebarPanel(
selectInput(inputId = "site",
label = "Site",
choices = c("Manosque upper terrace (Siame et al. 2004, France)" = "siame2004local",
"Test" = "siame2004local2")
),
# Input: Slider pour les parametres
sliderInput(inputId = "age",
label = "Age (a):",
......@@ -36,72 +37,72 @@ ui <- fluidPage(
max = 100000,
step = 1000,
value = 10000),
sliderInput(inputId = "ero",
label = "Denudation rate (m/Ma):",
min = 0,
max = 200,
step = 10,
value = 50),
value = 50),
sliderInput(inputId = "rho",
label = "Density (g/cm3):",
min = 1.5,
max = 3,
step = 0.05,
value = 2.6)
value = 2.6)
),
# Main panel for displaying outputs ----
mainPanel(
# Output: Histogram ----
plotOutput(outputId = "distPlot"),
plotOutput(outputId = "distPlot2")
)
)
)
# Define server logic required to draw a histogram ----
server <- function(input, output) {
output$distPlot <- renderPlot({
#unité longeurs: cm
#unité poids: g
#unité temps: an
#paramètres d'atténuation
Lspal=160. #longueur d'attenuation spalation (g/cm2)
Lslow=1500. #longueur d'attenuation stopping muons (g/cm2)
Lfast=4320. #longueur d'attenuation fast muons (g/cm2)
Lambda=c(Lspal,Lslow,Lfast)
#param production et decay 10Be, 26Al et 14C
prm=data.frame(matrix(NA, nrow=4, ncol=3))
# 10Be
prm[1,1]=4.49 # production spalation slhl (at/g)
prm[2,1]=0.012 # production stopping muons slhl (at/g)
prm[3,1]=0.039 # production fast muons slhl (at/g)
prm[2,1]=0.012 # production stopping muons slhl (at/g)
prm[3,1]=0.039 # production fast muons slhl (at/g)
prm[4,1]=log(2)/1.36e6 #decay constant (1/an)
# 26Al
prm[1,2]=prm[1,1]*6.61 # production spalation slhl (at/g)
prm[2,2]=0.84 # production stopping muons slhl (at/g)
prm[3,2]=0.081 # production fast muons slhl (at/g)
prm[2,2]=0.84 # production stopping muons slhl (at/g)
prm[3,2]=0.081 # production fast muons slhl (at/g)
prm[4,2]=log(2)/0.717e6 #decay constant (1/an)
# 14C
prm[1,3]=12.3 # production spalation slhl (at/g)
prm[2,3]=3.31 # production spalation slhl (at/g)
prm[3,3]=0 # production spalation slhl (at/g)
prm[4,3]=log(2)/5730 #decay constant (1/an)
data=read.table("data/sample_data.dat",header=TRUE)
#récupération des variables
data = data[data$study==input$site,]
alt = data$altitude[1]
......@@ -111,17 +112,17 @@ server <- function(input, output) {
age = input$age
zmax = max(data$depth)*rho*1.5 # cm -> g/cm2
C0 = 0
# scaling factors
Pression = stone2000pressure(alt)
Sn = stone2000neutrons(Pression,lat)
Sm = stone2000muons(Pression,lat)
# depth profil (g/cm2)
nz=200 # number of increment
dz=zmax/(nz-1)
z=seq(0,zmax,length.out=nz)
# plot data
plot(data$concentration,data$depth,
ylim=rev(c(0,max(data$depth)*1.5)),xlim=c(0,max(data$concentration)*1.5),
......@@ -129,16 +130,16 @@ server <- function(input, output) {
grid()
arrows(data$concentration-data$error,data$depth,data$concentration+data$error,data$depth,length=0)
points(data$concentration,data$depth,pch=21,bg="lightblue",cex=2)
# compute profile
Be10 = concentration_euler(z,ero,age,C0,prm[,1],c(Sn,Sm),Lambda)
# plot profile
lines(Be10,z/rho,type="l",col="red",lwd=4)
})
}
# Create Shiny app ----
......
library(shiny)
data = read.table("data/sample_data.csv",header=T,sep=",")
ch1 = data[data$label==1,]$study
names(ch1) <- data[data$label==1,]$title
# Define UI for app that draws a histogram ----
ui <- fluidPage(
img(src='amu.png', align = "right", width=100),
# App title ----
titlePanel("Cosmos profile"),
# Sidebar layout with input and output definitions ----
sidebarLayout(
sidebarPanel(
htmlOutput("study_selector"), #add selectinput boxes from objects created in server
htmlOutput("samples_selector"),
htmlOutput("age_selector"),
htmlOutput("ero_selector")
),
# Main panel for displaying outputs ----
mainPanel(
# Output:
plotOutput(outputId = "distPlot")
)
)
)
# Define server logic required to draw a histogram ----
server <- function(input, output) {
output$study_selector = renderUI({ #creates study select box object called in ui
selectInput(inputId = "site",
label = "Site",
choices = ch1)
})
output$samples_selector = renderUI({ #creates sample check box object called in ui
sp = unique(data[data$study == input$site,]$label)
names(sp) <- unique(data[data$study == input$site,]$label)
checkboxGroupInput(inputId = "sel_sp",
label = "Selected samples",
choices = sp,
selected = sp,
inline = TRUE)
})
output$age_selector = renderUI({#creates age slider called in ui
amin = data[data$study == input$site & data$label==1,]$amin
amax = data[data$study == input$site & data$label==1,]$amax
sliderInput(inputId = "age",
label = "Age (ka):",
min = amin,
max = amax,
step = 1,
value = amin)
})
output$ero_selector = renderUI({#creates erosion slider called in ui
emin = data[data$study == input$site & data$label==1,]$emin
emax = data[data$study == input$site & data$label==1,]$emax
sliderInput(inputId = "ero",
label = "Denudation rate (m/Ma):",
min = emin,
max = emax,
step = 1,
value = emin)
})
output$distPlot <- renderPlot({
site <- input$site
print(site)
data_sel = data[data$study == site,]
#
#print(data_sel)
plot(NA,xlim=range(data_sel$C),ylim=range(data_sel$depth))
#print(input$age)
#plot(data_sel$C,data_sel$depth)
#plot(NA,xlim=range(data_sel$C),ylim=range(data_sel$depth))
#points(data_sel$C,data_sel$depth)
})
}
# Create Shiny app ----
shinyApp(ui = ui, server = server)
# TODO starting values
# TODO choice of colors
library("shiny")
library("TCNtools")
library("akima")
data = read.table("data/sample_data.csv",header=T,sep=",")
ch1 = data[data$label==1,]$study
names(ch1) <- data[data$label==1,]$title
start = "siame2004local"
# Define UI for app that draws a histogram ----
ui <- fluidPage(
img(src='amu.png', align = "right", width=100),
# App title ----
titlePanel("Cosmos profile"),
# Sidebar layout with input and output definitions ----
sidebarLayout(
sidebarPanel(
selectInput(inputId = "study",
label = "Site",
choices = ch1,
selected = ch1[1]),
checkboxGroupInput(inputId = "sel_sp",
label = "Selected samples",
choices = unique(data[data$study == start,]$label),
selected = unique(data[data$study == start,]$label),
inline = TRUE),
# Input: Slider age
sliderInput(inputId = "age",
label = "Age (ka):",
min = data[data$study==start & data$label==1,]$amin,
max = data[data$study==start & data$label==1,]$amax,
step = 1,
value = data[data$study==start & data$label==1,]$amin),
# Input: Slider erosion
sliderInput(inputId = "ero",
label = "Denudation rate (m/Ma):",
min = data[data$study==start & data$label==1,]$emin,
max = data[data$study==start & data$label==1,]$emax,
step = 1,
value = data[data$study==start & data$label==1,]$emin),
# Input: Slider density
sliderInput(inputId = "rho",
label = "Density (g/cm3):",
min = 1.2,
max = 2.7,
step = 0.1,
value = 2),
# Input: Slider inheritance
sliderInput(inputId = "C0",
label = "Inheritance (x1000 at/g):",
min = 0,
max = 1000,
step = 10,
value = 0)
),
# Main panel for displaying outputs ----
mainPanel(
# Output:
plotOutput(outputId = "distPlot")
)
)
)
# Define server logic required to draw a histogram ----
server <- function(input, output,session) {
# Attenuation lengths
Lambda = c(160,1500,4320) # g/cm2
names(Lambda) <- c("Lspal","Lstop","Lfast") # we just give names to the element of the vector
# Production and decay parameters
prm = matrix(c( 4.01 , 0.012 , 0.039 , log(2)/1.36e6,
27.93 , 0.84 , 0.081 , log(2)/0.717e6,
12.22 , 3.31 , 0 , log(2)/5730),
nrow = 4,ncol=3 )
colnames(prm) <- c("Be10","Al26","C14") # we just give names to the columns of the matrix
rownames(prm) <- c("Pspal","Pstop","Pfast","lambda") # we just give names to the rows of the matrix
observe({
# updating inputs
updateCheckboxGroupInput(session, "sel_sp",choices = unique(data[data$study == input$study ,]$label),
selected = unique(data[data$study == input$study,]$label),inline=TRUE)
#
updateSliderInput(session, "age", min = data[data$study==input$study & data$label==1,]$amin,
max = data[data$study==input$study & data$label==1,]$amax,
value = data[data$study==input$study & data$label==1,]$amin)
#
updateSliderInput(session, "ero", min = data[data$study==input$study & data$label==1,]$emin,
max = data[data$study==input$study & data$label==1,]$emax,
value = data[data$study==input$study & data$label==1,]$emin)
})
output$distPlot <- renderPlot({
# compute scaling -> optimize
P = atm_pressure(alt=data[data$study==input$study & data$label==1,]$altitude,model="stone2000")
S = as.numeric(scaling_st(P,data[data$study==input$study & data$label==1,]$latitude))
layout(matrix(c(1,2), nrow = 1, ncol = 2, byrow = TRUE), widths =c(1,1))
# plot data
print(as.numeric(input$sel_sp))
data_st = data[data$study == input$study,]
data_sel = data_st[data_st$label %in% as.numeric(input$sel_sp),]
plot(NA,xlim=range(data_st$C),ylim=rev(range(data_st$depth,0)))
arrows(data_st$C-data_st$C_e,data_st$depth,data_st$C+data_st$C_e,data_st$depth,col="grey",length = 0)
points(data_st$C,data_st$depth,pch=21,col="grey",bg="grey",cex=2)
points(data_sel$C,data_sel$depth,pch=21,col="red",bg="red",cex=2)
arrows(data_sel$C-data_sel$C_e,data_sel$depth,data_sel$C+data_sel$C_e,data_sel$depth,col="red",length = 0)
text(data_st$C,data_st$depth,data_st$label,col="white")
#
z = seq(0,max(data_st$depth),length.out = 20)*input$rho
C = solv_conc_eul(z,input$ero*100/1e6*input$rho,input$age*1000,0,prm[,"Be10"],S,Lambda) + input$C0*1000
lines(C,z/input$rho,lty=2,lwd=3)
# explo
n = 50
v_ero = seq(data[data$study==input$study & data$label==1,]$emin,data[data$study==input$study & data$label==1,]$emax,length.out = n) / 1e6 *100 *input$rho
v_age = seq(data[data$study==input$study & data$label==1,]$amin,data[data$study==input$study & data$label==1,]$amax,length.out = n) *1000
res = expand.grid(v_ero,v_age)
colnames(res) <- c("ero","age")
for (i in 1:nrow(res)){
data_sel$Cmod = solv_conc_eul(data_sel$depth*input$rho,res$ero[i],res$age[i],0,prm[,"Be10"],S,Lambda) + input$C0*1000
res$chi2[i] = sum((data_sel$Cmod-data_sel$C)^2/data_sel$C_e^2)
}
imax = which.min(res$chi2)
C = solv_conc_eul(z,res$ero[imax],res$age[imax],0,prm[,"Be10"],S,Lambda) + input$C0*1000
lines(C,z/input$rho,lty=1,lwd=3,col="red")
# plot surface
spline<-interp(res$ero/max(res$ero),res$age/max(res$age),log10(res$chi2),duplicate="mean",nx=200,ny=200)
contour(spline$x*max(res$ero)*1e6/100/input$rho,spline$y*max(res$age)/1000,spline$z)
points(input$ero,input$age,pch=21,bg="blue")
points(res$ero[imax]*1e6/100/input$rho,res$age[imax]/1000,pch=21,bg="red")
})
}
# Create Shiny app ----
shinyApp(ui = ui, server = server)
#install.packages( "maps", dependencies = TRUE) #run this to install R package maps
################################- warning this will update existing packages if already installed
#*save the following code in a file named app.R *
library(shiny)
library(maps)
##Section 1 ____________________________________________________
#load your data or create a data table as follows:
countyData = read.table(
text = "State County
Delaware Kent
Delaware 'New Castle'
Delaware Sussex
'Rhode Island' Bristol
'Rhode Island' Kent
'Rhode Island' Newport
'Rhode Island' Providence
'Rhode Island' Washington",
header = TRUE, stringsAsFactors = FALSE)
##Section 2 ____________________________________________________
#set up the user interface
ui = shinyUI(
fluidPage( #allows layout to fill browser window
titlePanel("Reactive select input boxes"),
#adds a title to page and browser tab
#-use "title = 'tab name'" to name browser tab
sidebarPanel( #designates location of following items
htmlOutput("state_selector"),#add selectinput boxs
htmlOutput("county_selector")# from objects created in server
),
mainPanel(
plotOutput("plot1") #put plot item in main area
)
) )
##Section 3 ____________________________________________________
#server controls what is displayed by the user interface
server = shinyServer(function(input, output) {
#creates logic behind ui outputs ** pay attention to letter case in names
output$state_selector = renderUI({ #creates State select box object called in ui
selectInput(inputId = "state", #name of input
label = "State:", #label displayed in ui
choices = as.character(unique(countyData$State)),
# calls unique values from the State column in the previously created table
selected = "Delaware") #default choice (not required)
})
output$county_selector = renderUI({#creates County select box object called in ui
data_available = countyData[countyData$State == input$state, "County"]
#creates a reactive list of available counties based on the State selection made
selectInput(inputId = "county", #name of input
label = "County:", #label displayed in ui
choices = unique(data_available), #calls list of available counties
selected = unique(data_available)[1])
})
output$plot1 = renderPlot({ #creates a the plot to go in the mainPanel
map('county', region = input$state)
#uses the map function based on the state selected
map('county', region =paste(input$state,input$county, sep=','),
add = T, fill = T, col = 'red')
#adds plot of the selected county filled in red
})
})#close the shinyServer
##Section 4____________________________________________________
shinyApp(ui = ui, server = server) #need this if combining ui and server into one file.
study,sample,label,depth,C,C_e,latitude,altitude,title,amin,amax,emin,emax
siame2004local,PCAR97-9,1,75,0.37e5,0.11e5,43.8482,330,"Manosque upper terrace, France (Siame et al. 2004, France)",0,300,30,100
siame2004local,PCAR97-9,2,130,0.24e5,0.09e5,NA,NA,NA,NA,NA,NA,NA
siame2004local,PCAR97-9,3,175,0.16e5,0.04e5,NA,NA,NA,NA,NA,NA,NA
siame2004local,PCAR97-9,4,300,0.06e5,0.05e5,NA,NA,NA,NA,NA,NA,NA
hidy2010geologically,GC-04-LF-404.30s,1,27.5,568744,17062,36.853,985,"Lees Ferry M4 terrace, Arizona (Hidy et al., 2010)",40,120,0,10
hidy2010geologically,GC-04-LF-404.60s,2,57.5,406713,12201,NA,NA,NA,NA,NA,NA,NA
hidy2010geologically,GC-04-LF-404.100s,3,97.5,292243,8767,NA,NA,NA,NA,NA,NA,NA
hidy2010geologically,GC-04-LF-404.140s,4,137.5,203072,6092,NA,NA,NA,NA,NA,NA,NA
hidy2010geologically,GC-04-LF-404.180s,5,177.5,157209,4716,NA,NA,NA,NA,NA,NA,NA
hidy2010geologically,GC-04-LF-404.220s,6,217.5,134198,4025,NA,NA,NA,NA,NA,NA,NA
......@@ -59,7 +59,7 @@ S = as.numeric(scaling_st(P,mean(data$latitude)))
S=c(1,1)
data$C = solv_conc_euler(data$depth*rho,ero,age,0,prm[,1],S,Lambda) + C0
data$C = solv_conc_eul(data$depth*rho,ero,age,0,prm[,1],S,Lambda) + C0
n = 50
v_ero = seq(25,90,length.out = n) / 1e6 *100 *rho
......@@ -68,12 +68,17 @@ res = expand.grid(v_ero,v_age)
colnames(res) <- c("ero","age")
for (i in 1:nrow(res)){
data$C = solv_conc_euler(data$depth*rho,res$ero[i],res$age[i],0,prm[,1],S,Lambda) + C0
data$C = solv_conc_eul(data$depth*rho,res$ero[i],res$age[i],0,prm[,1],S,Lambda) + C0
res$chi2[i] = sum((data$C-data$concentration)^2/data$error^2)
}
spline<-interp(res$ero/max(res$ero),res$age/max(res$age),log10(res$chi2),duplicate="mean",nx=200,ny=200)
contour(spline$x*max(res$ero)*1e6/100/rho,spline$y*max(res$age),spline$z)
res$p = dchisq(res$chi2,2)
spline<-interp(res$ero/max(res$ero),res$age/max(res$age),1-res$p,duplicate="mean",nx=200,ny=200)
contour(spline$x*max(res$ero)*1e6/100/rho,spline$y*max(res$age),spline$z)
......
library(RColorBrewer)
display.brewer.all(colorblindFriendly = TRUE)
display.brewer.pal(n = 8, name = 'Dark2')
devtools::install_git(url = "https://gitlab.osupytheas.fr/vgodard/tcntools")
library("TCNtools")
# Attenuation lengths
Lambda = c(160,1500,4320) # g/cm2
names(Lambda) <- c("Lspal","Lstop","Lfast") # we just give names to the element of the vector
# Production and decay parameters
prm = matrix(c( 4.01 , 0.012 , 0.039 , log(2)/1.36e6,
27.93 , 0.84 , 0.081 , log(2)/0.7