Skip to content

Instantly share code, notes, and snippets.

@sickel
Last active January 24, 2018 13:46
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save sickel/b051a7938d0147b3a83031ac05b32dd1 to your computer and use it in GitHub Desktop.
Save sickel/b051a7938d0147b3a83031ac05b32dd1 to your computer and use it in GitHub Desktop.
An R script to make a plot of the number of nuclear tests per country and year
if(!exists("lang")){lang="EN"}
directory="https://gist.githubusercontent.com/sickel/8accd0d2ca32c9f7624e10f4837095f5/raw/f70630c17010bc33407caa3d08dc8d6283bc7f23/"
# or download https://gist.githubusercontent.com/sickel/8accd0d2ca32c9f7624e10f4837095f5/raw/f70630c17010bc33407caa3d08dc8d6283bc7f23/NuclearTests.csv
# and set the directory to where the file is stored.
# Data source: http://www.ldeo.columbia.edu/~richards/my_papers/WW_nuclear_tests_IASPEI_HB.pdf and https://en.wikipedia.org/wiki/List_of_nuclear_weapons_tests_of_North_Korea
# In theory, it should be possible to source this file directly into R from the "Raw" url, for some reason this fails,
# so copy the script and past it into an editor.
# These are the parts that need to be translated if that is needed
countries=c("USA","USSR","UK","France","China","India","Pakistan","North Korea","Others") # List of countries that have been doing tests
cntrhead="Country" # Head of the legend
ylab="Number of tests"     # Label of y axis
atmos="Atmospheric"        # Text in upper left corner
uground="Underground"      # Text in lower left corner
# Add on translated fields and wrap them in if(lang==...) and set the lang variable before running the script
if(lang=='NO'){
 countries=c("USA","Sovjetunionen","Storbritannia","Frankrike","Kina","India","Pakistan","Nord-Korea","Andre")
 cntrhead="Nasjon"
 ylab="Antall tester"
 atmos="Atmosfærisk"
 uground="Underjordisk"
}
datafile="NuclearTests.csv"
# Should not be changed:
cntrcol="Country"
yrcol="Year"
library(ggplot2)
nucTests <- read.csv(file=paste(directory,datafile,sep=""), header=TRUE, sep=",")
nucTests[is.na(nucTests)]=0
atm=nucTests[c(1,2)]
atm[3]=countries[1]
names(atm)[3]=cntrcol
for (i in c(4,6,8,10,12,14,16,18)){
temp=nucTests[c(1,i)]
 temp[3]=countries[i/2]
 names(temp)[2]=names(atm)[2]
 names(temp)[3]=cntrcol
 atm=rbind(atm,temp)
}
names(atm)[1]=yrcol
names(atm)[2]="N"
atm$Country=factor(atm$Country,levels=countries)
ugd=nucTests[c(1,3)]
ugd[3]=countries[1]
names(ugd)[3]=cntrcol
for (i in c(5,7,9,11,13,15,17,19)){
 temp=nucTests[c(1,i)]
 temp[3]=countries[(i-1)/2]
 names(temp)[2]=names(ugd)[2]
 names(temp)[3]=cntrcol
 ugd=rbind(ugd,temp)
}
names(ugd)[1]=yrcol
names(ugd)[2]="N"
ugd$Country=factor(ugd$Country,levels=countries)
ugd[2]=-1*ugd[2]
p=ggplot()
width=0.8
p=p+geom_bar(data=atm, aes(x=Year,y=N, fill=Country),stat="identity", width=width)
p=p+geom_bar(data=ugd, aes(x=Year,y=N, fill=Country),stat="identity", width=width)
p=p+xlab("")+labs(fill=cntrhead)+ylab(ylab)
p=p+theme(axis.text.x = element_text(colour="grey20",size=17,angle=0,hjust=.5,vjust=.5,face="plain"),
          axis.text.y = element_text(colour="grey20",size=17,angle=0,hjust=1,vjust=.5,face="plain"), 
          axis.title.y = element_text(colour="grey20",size=17,angle=90,hjust=.5,vjust=.5,face="plain"))
p=p+annotate("text", x = 1950, y = 125, label = atmos,size=5)
p=p+annotate("text", x = 1950, y = -100, label = uground,size=5)
p=p+theme(legend.text=element_text(size=15))
p=p+theme(legend.title=element_text(size=17))
p=p+scale_fill_brewer(palette="Paired")
# Use RColorBrewer::display.brewer.all() to show all default palettes
# abs to get positive numbers for the underground detonations
p=p+scale_y_continuous(label=abs)
p=p+scale_x_continuous(breaks=seq(1940,2020,10))
if(names(dev.cur())=="null device"){
# If no plotting window is open, open one with the correct proportions
width=300
height=width/1.618
x11(width=width,height=height)
}
cat("p to plot")
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment