r Documentation

Starting display moniter:
xhost + status
/sbin/ifconfig (to get the ip adress of work station)
export DISPLAY=anzak.hpc.unbc.ca:0
export DISPLAY=142.207.117.35:0.0
OR
xhost + status
ssh -X stratus

mat<-(matrix(1:2,1,2))
layout(mat)
layout.show(2)

mat<-(matrix(1:4,2,2))
layout(mat)
layout.show(4)

mat<-(matrix(1:6,2,3))
layout(mat)
layout.show(6)

hotd<-read.table(file="hotday3.out",header=F,na.strings="-999",col.name=c("dat","jd","yr","mo","dy","n","tmp","tmp3","ppt","spd","dir"))

hotJA<-read.table(file="hotday3JA.out",na.strings="-999",header=F,col.name=c("dat","jd","yr","mo","dy","n","tmp","tmp3","ppt","spd","dir"))

peak<-read.table(file="peak.out",header=F,na.strings="-999",col.name=c("dat","jd","yr","mo","dy","n","tmp","tmp3","ppt","spd","dir"))

peakJA<-read.table(file="peakJA.out",header=F,na.strings="-999",col.name=c("dat","jd","yr","mo","dy","n","tmp","tmp3","ppt","spd","dir"))

brja<-read.table(file="in.dat.HC2",header=F,col.name=c("dat","yr","mo","dy","tmp3","tmp","uk1","uk2","uk3","day3"),na.string="-999",skip=1)

plot(peak$mo,peak$n,xlab="Month",ylab="Cycle Length in Days",main="PE For All Months",xlim=c(4,10),ylim=c(0,50))
text(6,48, as.expression(substitute(max==12)))
text(6,46,as.expression(substitute(min==4)))
text(6,44,as.expression(substitute(mean==5.2)))
text(6,42,as.expression(substitute(cycles==92)))
dev.copy2eps(file="peak.eps")


plot(hotd$mo,hotd$n,xlab="Month",ylab="Cycle Lenght in Days",main="FW For All Months",xlim=c(4,10))
text(5,48, as.expression(substitute(max==48)))
text(5,46,as.expression(substitute(min==4)))
text(5,44,as.expression(substitute(mean==8.1)))
text(5,42,as.expression(substitute(cycles==388)))
dev.copy2eps(file="Fair_Weather.eps")
#dev.copy2eps(file="peak_and_FW.eps")

hist(peak$n,breaks=c(4,5,6,7,8,9,10,11,12,13,14,15,16,17,18),main="PE Cycle Lenghts for All Months",xlab="length of Cycle in Days",ylim=c(0,70))
text(14,60,as.expression(substitute(mean==5.2)))
#dev.copy2eps(file="CycleFreqPeakAll.eps")

hist(peakJA$n,breaks=c(4,5,6,7,8,9,10,11,12,13,14,15,16,17,18),main="PE Cycle Lenghts for Jul. & Aug.",xlab="Length of Cycle in Days",ylim=c(0,70))
text(14,60,as.expression(substitute(mean==5.3)))
#dev.copy2eps(file="CycleFreqPeakJA.eps")

hist(hotd$n,breaks=seq(from=0,to=50,by=1),main="FW Cycle Lenghts for All Months",xlab="Length of Cycle in Days",freq=T)
text(30,80,as.expression(substitute(mean==8.1)))
#dev.copy2eps(file="CycleFreqFWAll.eps")

hist(hotJA$n,breaks=seq(from=0,to=50,by=1),main="FW Cycle Lenghts for Jul. & Aug.",xlab="Length of Cycle in Days",freq=T,ylim=c(0,100))
text(30,80,as.expression(substitute(mean==9.1)))
#dev.copy2eps(file="CycleFreqFWJA.eps")
dev.copy2eps(file="histograms.eps")

mat<-(matrix(1:2,1,2))
layout(mat)
layout.show(2)

hist(peak$mo,breaks=c(3,4,5,6,7,8,9,10),ylab="Number of cycles",bg="blue",main="PE Cycles per Month",xlab="Month",ylim=c(0,130),axes=F)
axis(2)
axis(1,labels=c(3,4,5,6,7,8,9,10),adj="0")
text(6.5,30,"July")
text(7.5,30,"August")
text(6.5,25,"30%")
text(7.5,25,"25%")
#dev.copy2eps(file="MonthFreqPeak.eps")

hist(hotd$mo,main="FW Cycles per Month",ylab="Number of Cycles",xlab="Month",freq=T,breaks=c(3,4,5,6,7,8,9,10),ylim=c(0,130),axes=FALSE)
axis(2)
axis(1,at=c(3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5),pch=c(3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5),labels=c(NA,4,NA,5,NA,6,NA,7,NA,8,NA,9,NA,10,NA,NA))
text(6.5,80,"July")
text(7.5,80,"August")
text(6.5,75,"33%")
text(7.5,75,"23%")
#dev.copy2eps(file="MonthFreqFW2.eps")
#dev.copy2eps(file="MonthFreqBoth.eps")

boxplot(split(hotd$n,hotd$mo),main="FW Cycle Length by Month", xlab="Month", ylab="Length of Cycle in Days")
#dev.copy2eps(file="BPCycleByMonthFW.eps")

boxplot(split(peak$n,peak$mo),main="PE Cycle Length by Month", xlab="Month", ylab="Length of Cycle in Days")
#dev.copy2eps(file="BPCycleByMonthPeak.eps")

boxplot(split(peak$yr,peak$mo),main="Annual Distribution of PE Cycles by Month", xlab="Month")
#dev.copy2eps(file="BPAnnualPeak.eps")

boxplot(split(hotd$yr,hotd$mo),main="Annual Distribution of FW Cycles by Month", xlab="Month")
#dev.copy2eps(file="BPAnualFW.eps")
dev.copy2eps(file="BoxPlots.eps")

Climatology comparison:

all<-read.table(file="date.out",header=F,na.strings="-999",col.name=c("dat","jd","yr","mo","dy","tmp","tmp3","ppt","spd","dir"))

clm<-read.table(file="dateJA.out",na.strings="-999",header=F,col.name=c("dat","jd","yr","mo","dy","tmp3","tmp","ppt","spd","dir"))
yr,mo,dy,hr,tmp,tmp2,prs,spd,dir,sol

ALL HLY DATA FOR JULY AND AUGUST ONLY:

hly<-read.table(file="hlyJA.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

hlyp<-read.table(file="hlypeakJA.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

hlyfw<-read.table(file="hlyhotday3JA.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

hly15<-read.table(file="hlyJA15.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

hlyp15<-read.table(file="hlypeakJA15.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

hlyfw15<-read.table(file="hlyhotday3JA15.out",na.strings="-999",header=F,col.name=c("yr","mo","dy","hr","tmp","tmp2","prs","spd","dir","sol"))

plot(peakJA$jd,peakJA$n,main="Jul15-Aug.15 vs all Jul-Aug PE Cycle Friquency",ylab="Cycle Length in Days",xlab="Julian Date",sub=("Blue lines deliniate Jul.15 to Aug.15 and include 71% of cycles and 1/2 total time."))
abline(v=227,col="blue",untf=F)
abline(v=196,untf=F,col="blue")
text(188,11,"Cycles=11")
text(235,11,"Cycles=9")
text(205,11,"Cycles=51")
text(205,11.5,"Frequency=1.7")
text(188,11.5,"Friquency=0.7")
text(235,11.5,"Friquency=0.6")
#dev.copy2eps(file="peakJAvsJA15.eps")

Loop for counting Months:

i<-1;no<-0
while(i<=length(hotd$mo)){
if(hotd$mo[i]==8){
no<-no+1
};
i<-i+1}

mat<-(matrix(1:4,2,2))
layout(mat)
layout.show(4)

hist(clm$dir,breaks=35,main="Wind Direction",xlab="wind direction",xlim=c(1,36),ylim=c(0,200))
hist(peakJA$dir,add=T,border="green",breaks=35)
dev.copy2eps(file="WindDirClmVsPeak.eps")

hist(clm$spd,breaks=120,main="Extreem Ghust Speed",xlab="Wind Speed",ylim=c(0,110))
hist(peakJA$spd,breaks=120,add=T,border="green")
dev.copy2eps(file="WindSpdClmVsPeak.eps")

hist(clm$ppt,breaks=seq(from=0,to=500,by=5),freq=T,main="Precipitation",xlab="Precipitation in 0.1mm",ylim=c(0,130))
hist(peakJA$ppt,breaks=seq(from=0,to=500,by=5),freq=T,add=T,border="green")
dev.copy2eps(file="PptClmVsPeak.eps")

hist(peakJA$ppt,breaks=seq(from=0,to=10,by=1),ylim=c(0,10),freq=T,main="Peak Precipitation",xlab="Precipitation in 0.1mm")
dev.copy2eps(file="PptPeak.eps")
hist(hly$prs)
hist(hlyp$prs,add=T,border="green")

boxplot(split(hly$prs,hly$hr),main="Pressure",ylab="pressure in 0.01 kPa",xlab="Hour",sub="Means: Peak Emergance=9368,Climatology=9358")
boxplot(split(hlyp$prs,hlyp$hr),add=T,border="blue")
boxplot(split(hlyfw$prs,hlyfw$hr),add=T,border="purple")
#dev.copy2eps(file="PressureClmVsPeak.eps")

boxplot(split(hly$spd,hly$hr),main="Wind Speed",ylab="Speed in kph",xlab="Hour",sub="Means: Peak Emergance=,Climatology=")
boxplot(split(hlyp$spd,hlyp$hr),add=T,border="blue")
#dev.copy2eps(file="PressureClmVsPeak.eps")



boxplot(split(hly$tmp,hly$hr),plot=F)
boxplot(split(hlyp$prs,hly$hr),plot=F)
boxplot(split(hlyp$tmp,hly$hr),plot=F)
boxplot(split(hly$tmp,hly$hr),plot=F)

assign("medprs",c(9358, 9359, 9360, 9361, 9363.5, 9365, 9367, 9369, 9370,  9369,  9368,  9365, 9363, 9359, 9357, 9354, 9352, 9351, 9349, 9349 , 9351, 9354 , 9354, 9357)

assign("medprsP",c(9370, 9371.0, 9372, 9373.0, 9376.0, 9378, 9381.0, 9386.0, 9386.0 ,9383.0, 9383.0, 9382, 9377.0, 9373.0, 9370.0, 9366.0,  9362, 9358.0,  9352, 9351.0, 9354.0, 9357.0))

assign("medtmpP",c(144 , 128,  122,  111, 110.0 , 106 , 122 , 156, 178.0,   200 ,  217, 233.0, 247.0 ,258.0 ,267.0, 270.0,   272, 270.0 ,267.0, 247.0, 222.0 ,  194 ,  178 ,161.0))

assign("medtmp",c(116 , 109, 104,  100,   94 ,  94 , 103 , 117 , 133  , 150  , 161  , 176  , 185 ,  194 , 200 ,  200  , 200  , 198 ,  191 ,  178 ,  158, 142.0,   132 ,  122))

assign("medfw",c(9376, 9377.5, 9378.5, 9380, 9381, 9381.5, 9384, 9386, 9386.5, 9384.5, 9382.5,  9378, 9374, 9368.5,  9364,  9361,  9359,  9354,  9352, 9353.5, 9355.5 , 9359,  9364,  9365))

assign("medvfwtmp",c(128, 117.5,  110,  100 ,  94 ,  93, 112.5,  139 , 156 ,  178 ,  190, 205.5 ,  217 ,228,   237 ,  239, 241.0,239.0 ,  233 ,  221 ,  196 ,  172 ,  156 , 139))

plot(medprs,pch="-",ylim=c(9300,9410),type="l",xlab="Hours",ylab="Pressure",main="Diurnal Median Temperature and Pressure")
add=t
points(medprsP,pch="-",col="blue",type="l")
points(medfw,pch="-",col="purple",type="l")
par(new=T)
plot(medtmpP,pch="-",col="red",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(medtmp,pch="-",col="orange",type="l")
points(medfwtmp,pch="-",col="dark red",type="l")
text(2,300,"tem.Fair Weather",col="dark red")
text(2,290,"temp.all",col="orange")
text(2,280,"temp.peak",col="red")
text(2,270,"pressure.all")
text(2,260,"pressure.peak",col="blue")
text(2,250,"pressure.Fair Weahter",col="purple")

#dev.copy2eps(file="PrsTmpMedianAll.eps")

Same as above but with means:

value<-split(hly$prs,hly$hr)
hlym<-sapply(value,mean,na.rm=T)

value<-split(hlyp$prs,hlyp$hr)
hlypm<-sapply(value,mean,na.rm=T)

value<-split(hlyfw$prs,hlyfw$hr)
hlyfwm<-sapply(value,mean,na.rm=T)


value<-split(hly$tmp,hly$hr)
hlymt<-sapply(value,mean,na.rm=T)

value<-split(hlyp$tmp,hlyp$hr)
hlypmt<-sapply(value,mean,na.rm=T)

value<-split(hlyfw$tmp,hlyfw$hr)
hlyfwmt<-sapply(value,mean,na.rm=T)


LINE PLOT
plot(hlym,pch="-",ylim=c(9300,9410),type="l",xlab="Hours",ylab="Pressure",main="Diurnal Mean Temperature and Pressure",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypm,pch="-",col="blue",type="l")
points(hlyfwm,pch="-",col="purple",type="l")
par(new=T)
plot(hlymt,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")
text(2,300,"temp.all",col="orange")
text(2,290,"temp.peak",col="red")
text(2,280,"tem.Fair Weather",col="dark red")
text(2,270,"pressure.all")
text(2,260,"pressure.peak",col="blue")
text(2,250,"pressure.Fair Weahter",col="purple")
#dev.copy2eps(file="PrsTmpMeanAll.eps")
mtext(4,Temperature,srt=90,"Temperature")

Point PLOT
plot(hlym,pch=15,ylim=c(9300,9410),xlab="Hours",ylab="Pressure",main="Diurnal Mean Temperature and Pressure",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypm,pch=17,type="p")
points(hlyfwm,pch=20)
par(new=T)
plot(hlymt,pch=22,ylim=c(90,300),ylab="",xlab="",axes=F)
axis(side=4)
add=T
points(hlypmt,pch=24)
points(hlyfwmt,pch=21)
text(2,300,"tem.Fair Weather=")
text(4,300,pch=15)
text(2,290,"temp.all",col="orange")
text(2,280,"temp.peak",col="red")
text(2,270,"pressure.all")
text(2,260,"pressure.peak",col="blue")
text(2,250,"pressure.Fair Weahter",col="purple")
#dev.copy2eps(file="PrsTmpMeanAll.eps")

ADDING JULY15 TO AUG15
value<-split(hly15$prs,hly15$hr)
hlym15<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$prs,hlyp15$hr)
hlypm15<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$prs,hlyfw15$hr)
hlyfwm15<-sapply(value,mean,na.rm=T)
value<-split(hly15$tmp,hly15$hr)
hlymt15<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$tmp,hlyp15$hr)
hlypmt15<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$tmp,hlyfw15$hr)
hlyfwmt15<-sapply(value,mean,na.rm=T)


plot(hlym,pch="-",ylim=c(9300,9410),type="l",xlab="Hours",ylab="Pressure",sub="dashed lines = Jul15 to Aug15 data only",main="Diurnal Mean Temperature and Pressure",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypm,pch="-",col="blue",type="l")
points(hlyfwm,pch="-",col="purple",type="l")
par(new=T)
plot(hlymt,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")
text(2,300,"temp.all",col="orange")
text(2,290,"temp.peak",col="red")
text(2,280,"tem.Fair Weather",col="dark red")
text(2,270,"pressure.all")
text(2,260,"pressure.peak",col="blue")
text(2,250,"pressure.Fair Weahter",col="purple")
#dev.copy2eps(file="PrsTmpMeanAll.eps")
par(new=T)
plot(hlym15,pch="-",ylim=c(9300,9410),type="l",lty=2,xlab="Hours",ylab="Pressure",main="Diurnal Mean Temperature and Pressure",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypm15,pch="-",col="blue",type="l",lty=2)
points(hlyfwm15,pch="-",col="purple",type="l",lty=2)
par(new=T)
plot(hlymt15,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",lty=2,axes=F)
axis(side=4)
add=T
points(hlypmt15,pch="-",col="red",type="l",lty=2)
points(hlyfwmt15,pch="-",col="dark red",type="l",lty=2)
#dev.copy2eps(file="MeanPrsAll.eps")

SAME PLOT FOR SOL:

value<-split(hly$sol,hly$hr)
hlyms<-sapply(value,mean,na.rm=T)
value<-split(hlyp$sol,hlyp$hr)
hlypms<-sapply(value,mean,na.rm=T)
value<-split(hlyfw$sol,hlyfw$hr)
hlyfwms<-sapply(value,mean,na.rm=T)

value<-split(hly15$sol,hly15$hr)
hlym15s<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$sol,hlyp15$hr)
hlypm15s<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$sol,hlyfw15$hr)
hlyfwm15s<-sapply(value,mean,na.rm=T)

LINE PLOT
plot(hlyms,pch="-",ylim=c(0,3000),type="l",xlab="Hours",ylab="Global Solar Radiation",sub="dashed lines = Jul15 to Aug15 data only",main="Diurnal Mean Temperature and Solar Radiation",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypms,pch="-",col="blue",type="l")
points(hlyfwms,pch="-",col="purple",type="l")
points(hlypm15s,pch="*",col="blue",type="l",lty=2)
points(hlyfwm15s,pch="*",col="purple",type="l",lty=2)
points(hlym15s,pch="*",type="l",lty=2)
par(new=T)
plot(hlymt,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")
points(hlypmt15,pch="-",col="red",type="l",lty=2)
points(hlyfwmt15,pch="-",col="dark red",type="l",lty=2)
points(hlymt15,pch="*",type="l",lty=2)
text(2,300,"temp.all",col="orange")
text(2,290,"temp.peak",col="red")
text(2,280,"tem.Fair Weather",col="dark red")
text(2,270,"sol.all")
text(2,260,"sol.peak",col="blue")
text(2,250,"sol.Fair Weahter",col="purple")
#dev.copy2eps(file="MeanSolAll.eps")

SAME PLOT FOR DEW POINT TEMPERATURE:

value<-split(hly$tmp2,hly$hr)
hlymt2<-sapply(value,mean,na.rm=T)
value<-split(hlyp$tmp2,hlyp$hr)
hlypmt2<-sapply(value,mean,na.rm=T)
value<-split(hlyfw$tmp2,hlyfw$hr)
hlyfwmt2<-sapply(value,mean,na.rm=T)

value<-split(hly15$tmp2,hly15$hr)
hlym15t2<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$tmp2,hlyp15$hr)
hlypm15t2<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$tmp2,hlyfw15$hr)
hlyfwm15t2<-sapply(value,mean,na.rm=T)

LINE PLOT
plot(hlymt2,pch="-",ylim=c(50,300),type="l",xlab="Hours",ylab="Temperature",main="Diurnal Mean Temperature and Dew Point Temperature",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypmt2,pch="-",col="blue",type="l")
points(hlyfwmt2,pch="-",col="purple",type="l")
points(hlymt,pch="-",col="orange",type="l")
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")

par(new=T)
plot(dew,pch="-",ylim=c(25,200),col="green",ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(dewp,pch="-",col="dark green",type="l")
points(dewfw,pch="-",col="light green",type="l")

text(4,200,"temp.all",col="orange")
text(4,195,"temp.peak",col="red")
text(4,190,"tem.Fair Weather",col="dark red")
text(4,185,"dew.all")
text(4,180,"dew.peak",col="blue")
text(4,175,"dew.Fair Weahter",col="purple")
text(4,170,"differemce.all", col="green")
text(4,165,"differemce.peak",col="dark green")
text(4,160,"differemce.Fair Weahter",col="light green")
#dev.copy2eps(file="MeanTmpAll.eps")

dew<-(hlymt-hlymt2)
dewp<-(hlypmt-hlypmt2)
dewfw<-(hlyfwmt-hlyfwmt2)

par(new=T)
plot(dew,pch="-",ylim=c(30,150)col="green",ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(dewp,pch="-",col="dark green",type="l")
points(dewfw,pch="-",col="light green",type="l")

SAME PLOT FOR SPD:

value<-split(hly$spd,hly$hr)
hlymspd<-sapply(value,mean,na.rm=T)
value<-split(hlyp$spd,hlyp$hr)
hlypmspd<-sapply(value,mean,na.rm=T)
value<-split(hlyfw$spd,hlyfw$hr)
hlyfwmspd<-sapply(value,mean,na.rm=T)

value<-split(hly15$spd,hly15$hr)
hlym15spd<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$spd,hlyp15$hr)
hlypm15spd<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$spd,hlyfw15$hr)
hlyfwm15spd<-sapply(value,mean,na.rm=T)

LINE PLOT
plot(hlymspd,pch="-",type="l",xlab="Hours",ylab="Wind Speed km/hr",ylim=c(0,12),sub="dashed lines = Jul15 to Aug15 data only",main="Diurnal Mean Temperature and Wind Speed",axes=F)
axis(2)
axis(1,c(1,24:1))
box()
add=t
points(hlypmspd,pch="-",col="blue",type="l")
points(hlyfwmspd,pch="-",col="purple",type="l")
points(hlypm15spd,pch="*",col="blue",type="l",lty=2)
points(hlyfwm15spd,pch="*",col="purple",type="l",lty=2)
points(hlym15spd,pch="*",type="l",lty=2)
par(new=T)
plot(hlymt,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")
points(hlypmt15,pch="-",col="red",type="l",lty=2)
points(hlyfwmt15,pch="-",col="dark red",type="l",lty=2)
points(hlymt15,pch="*",type="l",lty=2)
text(4,300,"temp.all",col="orange")
text(4,290,"temp.peak",col="red")
text(4,280,"tem.Fair Weather",col="dark red")
text(4,270,"spd.all")
text(4,260,"spd.peak",col="blue")
text(4,250,"spd.Fair Weahter",col="purple")
#dev.copy2eps(file="MeanSpdAll.eps")

SAME PLOT FOR DIR:

plot(hlymdir,ylim=c(0,25),pch="-",,ylab="Wind Direction 10's of degrees",xlab="Hours",type="l",axes=F)
axis(2)
axis(1,c(1,24:1))
box()
add=T
points(hlypmdir,pch="-",col="blue",type="l")
points(hlyfwmdir,pch="-",col="purple",type="l")
par(new=T)
plot(hlymt,pch="-",col="orange",ylim=c(90,300),ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmt,pch="-",col="red",type="l")
points(hlyfwmt,pch="-",col="dark red",type="l")
text(4,300,"temp.all",col="orange")
text(4,290,"temp.peak",col="red")
text(4,280,"tem.Fair Weather",col="dark red")
text(4,270,"dir.all")
text(4,260,"dir.peak",col="blue")
text(4,250,"dir.Fair Weahter",col="purple")
dev.copy2eps(file="MeanDirAll.eps")

SAME PLOT FOR SPD and DIR:

value<-split(hly$dir,hly$hr)
hlymdir<-sapply(value,mean,na.rm=T)
value<-split(hlyp$dir,hlyp$hr)
hlypmdir<-sapply(value,mean,na.rm=T)
value<-split(hlyfw$dir,hlyfw$hr)
hlyfwmdir<-sapply(value,mean,na.rm=T)

value<-split(hly15$dir,hly15$hr)
hlym15dir<-sapply(value,mean,na.rm=T)
value<-split(hlyp15$dir,hlyp15$hr)
hlypm15dir<-sapply(value,mean,na.rm=T)
value<-split(hlyfw15$dir,hlyfw15$hr)
hlyfwm15dir<-sapply(value,mean,na.rm=T)

LINE PLOT - ALL too jumbled
plot(hlymspd,ylim=c(0,12),pch="-",type="l",xlab="Hours",ylab="Wind Speed km/hr",sub="dashed lines = Jul15 to Aug15 data only",main="Diurnal Mean Wind Speed and Direction",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypmspd,pch="-",col="blue",type="l")
points(hlyfwmspd,pch="-",col="purple",type="l")
points(hlypm15spd,pch="*",col="blue",type="l",lty=2)
points(hlyfwm15spd,pch="*",col="purple",type="l",lty=2)
points(hlym15spd,pch="*",type="l",lty=2)
par(new=T)
plot(hlymdir,ylim=c(0,60),pch="-",col="orange",ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmdir,pch="-",col="red",type="l")
points(hlyfwmdir,pch="-",col="dark red",type="l")
points(hlypm15dir,pch="-",col="red",type="l",lty=2)
points(hlyfwm15dir,pch="-",col="dark red",type="l",lty=2)
points(hlym15dir,pch="*",col="orange",type="l",lty=2)
text(4,25,"dir.all",col="orange")
text(4,24,"dir.peak",col="red")
text(4,23,"dir.Fair Weather",col="dark red")
text(4,22,"spd.all")
text(4,21,"spd.peak",col="blue")
text(4,20,"spd.Fair Weahter",col="purple")
#dev.copy2eps(file="MeanSpdAll.eps")

WITHOUT 15-15

plot(hlymspd,ylim=c(0,12),pch="-",type="l",xlab="Hours",ylab="Wind Speed km/hr",main="Diurnal Mean Wind Speed and Direction",axes=F)
axis(2)
axis(1,c(1,24:1))
add=t
points(hlypmspd,pch="-",col="blue",type="l")
points(hlyfwmspd,pch="-",col="purple",type="l")
par(new=T)
plot(hlymdir,ylim=c(0,25),pch="-",col="orange",ylab="",xlab="",type="l",axes=F)
axis(side=4)
add=T
points(hlypmdir,pch="-",col="red",type="l")
points(hlyfwmdir,pch="-",col="dark red",type="l")

text(4,25,"dir.all",col="orange")
text(4,24,"dir.peak",col="red")
text(4,23,"dir.Fair Weather",col="dark red")
text(4,22,"spd.all")
text(4,21,"spd.peak",col="blue")
text(4,20,"spd.Fair Weahter",col="purple")
#dev.copy2eps(file="MeanSpdAll.eps")


PLOTS FOR CONVERTED WIND DATA (WIND.OUT)

wind<-read.table(file="hourly/windJA.out",header=F,col.name=c("spd","dir"))
windfw<-read.table(file="hourly/windfwJA.out",header=F,col.name=c("spd","dir"))
windp<-read.table(file="hourly/windpeakJA.out",header=F,col.name=c("spd","dir"))

dayw<-read.table(file="hourly/day_windJA.out",header=F,col.name=c("spd","dir"))
daywfw<-read.table(file="hourly/day_windfwJA.out",header=F,col.name=c("spd","dir"))
daywp<-read.table(file="hourly/day_windpeakJA.out",header=F,col.name=c("spd","dir"))

boxplot(split(wind$spd,wind$dir),plot=F)

plot(wind$dir,wind$spd)
add=t
points(windp$dir,windp$spd,pch=10,col="red")

mat<-(matrix(1:3,3,1))
layout(mat)
layout.show(3)

boxplot(split(wind$spd,wind$dir),sub="Climatology",main="Wind Speed and direction for July and August",ylab="Wind Speed (m/s)",xlab="Wind direction (degrees)")
add=T
points(windpm,pch="-",col="blue",type="l")
boxplot(split(windp$spd,windp$dir),sub="Peak",ylim=c(0,15))
boxplot(split(windfw$spd,windfw$dir),sub="Fair Weather",ylim=c(0,15))
#dev.copy2eps(file="Wind_boxplots.eps")

mat<-(matrix(1:3,3,1))
layout(mat)
layout.show(3)
boxplot(split(dayw$spd,dayw$dir),sub="Climatology",main="Day-Time Winds for July and August",ylab="Wind Speed (m/s)",xlab="Wind direction (degrees)")
boxplot(split(daywp$spd,daywp$dir),sub="Peak",ylim=c(0,14))
boxplot(split(daywfw$spd,daywfw$dir),sub="Fair Weather",ylim=c(0,14))
#dev.copy2eps(file="wind_day_boxplots.eps")

MEANS PLOT AS ABOVE:

value<-split(dayw$spd,dayw$dir)
daywm<-sapply(value,mean,na.rm=T)
value<-split(daywp$spd,daywp$dir)
daywpm<-sapply(value,mean,na.rm=T)
value<-split(daywfw$spd,daywfw$dir)
daywfwm<-sapply(value,mean,na.rm=T)

value<-split(wind$spd,wind$dir)
windm<-sapply(value,mean,na.rm=T)
value<-split(windp$spd,windp$dir)
windpm<-sapply(value,mean,na.rm=T)
value<-split(windfw$spd,windfw$dir)
windfwm<-sapply(value,mean,na.rm=T)

LINE PLOT
plot(windm,pch="-",type="l",xlab="Wind Direction (Degrees)",ylab="Mean Wind Speed m/s",sub="dashed lines = Daily Winds",main="Mean Winds",ylim=c(0,4))
axis(2)
axis(1,1:17,c(0,20,50,70,90,110,140,160,180,200,230, 250,  270,290, 320,340,360))
box()

add=t
points(windpm,pch="-",col="blue",type="l")
points(windfwm,pch="-",col="purple",type="l")
points(daywm,pch="*",col="red",type="l",lty=2)
points(daywfwm,pch="*",col="orange",type="l",lty=2)
points(daywpm,pch="*",col="dark red",type="l",lty=2)
text(2,4,"climatolgy - Day",col="red")
text(2,3.85,"Peak -Day",col="dark red")
text(2,3.7,"Fair Weather - Day",col="orange")
text(2,3.55,"Climatology")
text(2,3.4,"peak",col="blue")
text(2,3.25,"Fair Weahter",col="purple")
dev.copy2eps(file="MeanWinds.eps")