In fmichelsendis/sendis: A Systematic Evaluated Nuclear Data Integral Survey

knitr::opts_chunk$set(echo = params$printcode)
knitr::opts_chunk$set(warning = FALSE)
library(knitr)
library(ggplot2)
library(plotly)
library(dplyr)
library(kableExtra) 
library(data.table)
library(scales)
#library(latex2exp)
library(devtools)

# Install and load the development version from GitHub:
#devtools::install_github("fmichelsendis/sendis")
library(sendis)

return_kable<-function(df, caption=""){

k<-kable(df, 
      escape = FALSE, 
      caption = caption,
      format = "latex", 
      booktabs = TRUE,
      linesep ="")%>%
  kable_styling(full_width = TRUE,
                font_size = 8, 
                latex_options = c("scale_down", 
                                  "striped",
                                  "hold_position",
                                  "repeat_header")
                )
  return(k)
}

# sens<-fread("../data-raw/fms_all_sens.csv")%>%
#   filter(REACTION=="total")%>%
#   select(-V1, -KSENS1, -KSENS2, -KSENS3)

Diagnostics Report at time : `r params$date`

Tables

library(sendis)

options(knitr.kable.NA = '')

df<-sendis%>%
  filter(.data$INST!="PSI")%>%
  filter(.data$LIBVER!="JENDL-4.0u")

base<-unique(ungroup(df)%>%
    filter(.data$INST=="NEA")%>%
    select(FULLID, MODEL))

df<-ungroup(df)%>%
    select(-NCASES, -CUMUL, -CHISQ)

  d<-merge(
    base,
    df,
    by = c("FULLID", "MODEL"),
    all = FALSE
           )

table<-d%>% 
  build_derived_calcs()%>%
  mutate(
    CHISQ = cell_spec(CHISQ, "latex", color = ifelse(CHISQ > 10, "red", "darkgray"))
  )%>%
  return_chitable()%>%
  rename_table_header()

return_kable(table, caption = "Total $\\chi^2$ values")%>%
  add_header_above(c(" ", " " =1, "ENDFB" = 2, "JEFF" = 7, "JENDL"=1))%>%
  row_spec(0, bold = TRUE) %>%
  row_spec(5, bold = TRUE) %>%
  column_spec(1, width = "1.5cm")%>%
  footnote(
           symbol = c("Footnote Symbol 1; ", "Footnote Symbol 2"),
           footnote_as_chunk = T, title_format = c("italic", "underline"))

Isotopes sensitivities

df<-sendis%>%
  filter(LIBVER=="JEFF-3.3"|LIBVER=="JEFF-2.2")

df<-merge(df, sens, all.x=TRUE, all.y =FALSE)%>%
  filter(KSENSTOT>0.3)
  #filter(ISOTOPE=="U233" | ISOTOPE=="Be0")

p<-ggplot(df, 
          aes(x=ISOTOPE, y=RESIDUAL)) + 
  geom_boxplot(aes(shape=LIBVER, colour= ISOTOPE)) +  
  theme_bw() 

p

Plot 0

df<-sendis%>%
  select(-FISS, -FORM, -SPEC, -LIB, -VER)%>%
  build_splitcolumns()%>%
  mutate(EVALCASE = paste0(CASETYPE,"-",CASE))%>%
  #df<-sendis%>% #build_splitcolumns()%>%
  filter(LIBVER=="JEFF-2.2"|LIBVER=="JEFF-3.3", INST=="NEA")%>% 
  #filter(INST=="NRG")%>%
  group_by(EVALCASE)%>%
  mutate(
    MEAN_EXP_ERR = mean(EXPERR),
    MEAN_EALF = mean(EALF)
    )%>%
  arrange(CASETYPE)

df<-ungroup(df)

num<-df%>%select(EVALCASE)%>%unique()%>%mutate(NUM = 1:n())

df<-merge(df, num)

df<-merge(df, sens, all.x=TRUE, all.y =FALSE)%>%select(-REACTION)


p<-ggplot(df, aes(x=EVALCASE, y=COVERE)) + 
  geom_boxplot(aes(colour= LIBVER), position=position_dodge(0.0)) + 
  geom_hline(yintercept =  1, lty = "dotted") +
  geom_ribbon(aes(x= NUM, ymin = 1-MEAN_EXP_ERR, ymax=1+MEAN_EXP_ERR),
              fill="gray",
              alpha=0.2) +
  labs(x = "", y=TeX("$k_{eff}\\,C/E$")) + 
  scale_y_continuous(
    sec.axis = dup_axis(),
    # sec.axis = dup_axis(
    #   trans = ~.*1e+5-1e+5,
    #   name = TeX("$\\Delta \\, (pcm)$")),
    breaks = seq(-0.94, 1.12, by = 0.01),
    position = "bottom"
    )+
  theme_bw()

p

Plot 1

df<-sendis%>%
  arrange(FULLID)%>%
  filter(LIBVER=="JEFF-2.2"|LIBVER=="JEFF-3.3", INST=="NEA") 

g<-ggplot(df, aes(x=FULLID, y=COVERE, shape= FISS, colour= LIBVER)) +
  geom_point()+
  geom_errorbar(aes(ymin=COVERE-TOTERR, ymax=COVERE+TOTERR))+ 
  theme_bw()+
  geom_hline(yintercept =  1, lty = "dotted") +
  scale_y_continuous(breaks = seq(-1.1, 1.1, by = 0.01))+
  coord_flip() 

  # theme(
  #   legend.title = element_blank(),
  #   panel.grid = element_blank() 
  #   )

  # +
  # geom_hline(yintercept = -3, lty = "dotted") +
  # ylim(-25,25)+
  #scale_y_continuous(breaks = seq(-25, 25, by = 5))+
  # scale_x_continuous(breaks = seq(-10, 10, by = 20))+
  # labs(y = "$\\Delta$ \n",
  #      x = "\nRatio dette publique / produit intérieur brut (%)",
  #      title = "Données Reinhart et Rogoff corrigées, 1946-2009\n") +


g

df<-sendis%>%filter(LIB=="JEFF", INST=="NEA") 
df$ZERO<-0

g<-ggplot(df, aes(x=ZERO, y= RESIDUAL, shape= FISS, colour= SPEC)) +
  geom_jitter(width=0.001, height=0) +
  facet_grid(.~VER)+
  geom_hline(yintercept =  3, lty = "dotted") +
  geom_hline(yintercept = -3, lty = "dotted") +
  ylim(-25,25)+
  #scale_y_continuous(breaks = seq(-25, 25, by = 5))+
  scale_x_continuous(breaks = seq(-10, 10, by = 20))+
  labs(y = "$\\Delta$ \n",
       x = "\nRatio dette publique / produit intérieur brut (%)",
       title = "Données Reinhart et Rogoff corrigées, 1946-2009\n") +
  theme_bw()+
  theme(
    legend.title = element_blank(),
    panel.grid = element_blank()
    # axis.title.x=element_blank(),
    # axis.text.x=element_blank(),
    # axis.ticks.x=element_blank()
    )

g

dfl<-fread("../data-raw/libs.csv")
dfl<-filter(dfl, LIB!="TENDL", LIB!="EAF", LIB!="FENDL", LIB!="ROSFOND")
g<-ggplot(dfl, aes(x=YEAR, y=NUMCHUNKS, colour=LIB))+
  geom_point(aes(size = 1, alpha=0.5))+
  geom_line(aes(size=0.5))+ 
   theme_bw()+ theme(legend.title = element_blank())+
   scale_x_continuous(breaks = seq(1990, 2021, by = 5))+
  labs(x="\nYear of release", y= "Number of data segments\n")

g

#ggsave(filename = "plot2.png", width = 16, units = "cm", dpi = 450)


## Suspicious Outliers

### Plots

df<-sendis%>%filter(LIB=="JEFF"| LIB=="JENDL" | LIB=="ENDFB", INST=="NEA")%>%
  filter(FISS!="MIX")%>%
  filter(FISS!="SPEC")%>%
  filter(FISS!="U233")%>%
  filter(!(VER %in% c("3.0", "3.1", "3.2", "3.1.2")))

g<-ggplot(df, aes(x= EALF)) +
  geom_point(aes(y= COVERE, colour= LIBVER, shape= FISS)) +
  geom_ribbon(aes(ymin = 1-EXPERR, ymax= 1+EXPERR), fill="#406dad", alpha=0.1) +
    labs(x = "Average energy of neutron causing fission (eV)", y = "C/E")+
    scale_x_log10(breaks = trans_breaks("log10", function(x) 10^x),
                  labels = trans_format("log10", math_format(10^.x)))+
    annotation_logticks(sides = "b") +
    theme_bw()+ theme(legend.title = element_blank())

g
#ggsave(plot = g,filename = "plot4.png",dpi = 400)

outliers<-sendis%>%
  return_discrepant(INST1 = "NEA", INST2 = "NRG", sigma_treshold = 3)%>%
  select(-FULLID, -EALF, -CODE.x, -CODE.y, -CALCERR.x, -CALCERR.y, -TOT_STAT_ERR_PCM)%>%
  rename( 'Delta (pcm)' = DELTA_PCM,
          'k(NEA)' = CALCVAL.x,
          'k(NRG)' = CALCVAL.y,
          Nsigmas = DELTA_SIGMA_STAT)%>%
  arrange(SHORTID, -Nsigmas)%>%
  unique()

kable(outliers, 
      escape = FALSE, 
      format = "latex", 
      booktabs = TRUE,
      row.names = FALSE,
      linesep ="")%>%
  kable_styling(
    full_width = TRUE,
    font_size = 8,
    latex_options = c("scale_down",
                      "striped",
                      "hold_position",
                      "repeat_header")
    )%>%
  row_spec(0, bold = TRUE) %>%
  column_spec(2, width = "2cm")%>%
  column_spec(3, width = "2cm")