#
#    Example scripts in 'R' to plot Viking Seismology Data
#
#    R D Lorenz     7/22/2016


# Read in data into a 2-d array   (substitute your own path/filename  : this is the Z-axis amplitude in High rate mode)
x <- read.table("C:/Users/lorenrd1/Documents/seismometer/vl_seismic_revised/data/vpds_high_zdat.csv", sep=",", header=FALSE)

# choose record n
n <- 5397

# Extract Sol, hour, min, second  information for caption
title <- paste("Z-axis High Rate #", n,"     Sol ",x[n,3],"  ",x[n,4],":",x[n,5],":",x[n,6])

# Plot the 82 datapoints, which span just over 4 seconds
# For the event mode data, note that the records are 51 samples long, and at 1.01 Hz
plot(seq(1:83)/20.2,x[n,7:89], xlim=c(0,4.10),ylim=c(-128,128), main=title, xlab="Time (s)", ylab="Seismic Signal (DN)", typ="s") ; n<-n+1







# If you know by inspection that a set of records are consecutive in time, you can plot them
# together with this code.  Here we plot the 12 records from 16064 to 16078

nn <- 12
i <- 16064

# set up title and plot axes
title <- paste("Z-axis H #", i,":",i+nn,"     Sol ",x[i,3],"  ",x[i,4],":",x[i,5],":",x[i,6])
plot(c(0), xlim=c(0,nn*4.10),ylim=c(-128,128), main=title, xlab="Time (s)", ylab="Seismic Signal (DN)", typ="s")

# loop through the record set add each one to the plot
for (n in 0:nn) {
lines((82*n+seq(1:83))/20.2,x[n+i,7:89], typ="s")  }



# Here we explore the correlation between the X-amplitude of the event mode records with windspeed.
# We select out only those datapoints where the meteorology data was acquired within 4 seconds of the seismic data
# (absolute value of column 24 < 4)

x <- read.csv("C:/Users/lorenrd1/Documents/seismometer/vl_seismic_revised/data/event_wind_summary.tab", sep="", header=FALSE)
plot(x[abs(x[,24])<4,20],x[abs(x[,24])<4,13] , xlab="Wind Speed (m/s)", ylab="RMS X Amplitude (DN)", pch=20, cex=0.5, xlim=c(0,30))