function_analyses_paper <-"contribution_of_covariates_paper_interpolation_functions_10152013.R"

plot_transect_m2<-function(list_trans,r_stack,title_plot,disp=FALSE,m_layers){

  #This function creates plot of transects for stack of raster images.

  #Arguments:

  #list_trans: list of files containing the transects lines in shapefile format

  #r_stack: raster stack containing the information to extect

  #title_plot: plot title

  #disp: display and save from X11 if TRUE or plot to png file if FALSE

  #m_layers: index for layerers containing alternate units to be drawned on a differnt scale

  #RETURN:

  #list containing transect information

  nb<-length(list_trans)

  t_col<-rainbow(nb)

  t_col<-c("red","green","black")

  lty_list<-c("dashed","solid","dotted")

  list_trans_data<-vector("list",nb)

  #For scale 1

  for (i in 1:nb){

    trans_file<-list_trans[[i]][1]

    filename<-sub(".shp","",trans_file)             #Removing the extension from file.

    transect<-readOGR(dirname(filename), basename(filename))                 #reading shapefile 

    trans_data<-extract(r_stack, transect)

    if (disp==FALSE){

      png(file=paste(list_trans[[i]][2],".png",sep=""))

    }

    #Plot layer values for specific transect

    for (k in 1:ncol(trans_data[[1]])){

      y<-trans_data[[1]][,k]

      x<-1:length(y)

      m<-match(k,m_layers)

      if (k==1 & is.na(m)){

        plot(x,y,type="l",xlab="transect distance from coastal origin (km)", ylab=" maximum temperature (degree C)",

             ,cex=1.2,col=t_col[k])

        #axis(2)

      }

      if (k==1 & !is.na(m)){

        plot(x,y,type="l",col=t_col[k],lty="dotted",axes=F) #plotting fusion profile

        #axis(4,xlab="",ylab="elevation(m)")  

        axis(4,cex=1.2)

      }

      if (k!=1 & is.na(m)){

        #par(new=TRUE)              # new plot without erasing old

        lines(x,y,type="l",xlab="",ylab="",col=t_col[k],axes=F) #plotting fusion profile

        #axis(2,xlab="",ylab="tmax (in degree C)")

      }

      if (k!=1 & !is.na(m)){

        par(new=TRUE)              # key: ask for new plot without erasing old

        plot(x,y,type="l",col=t_col[k],xlab="",ylab="",lty="dotted",axes=F) #plotting fusion profile

        #axis(4,xlab="",ylab="elevation(m)")  

        axis(4,cex=1.2)

      } 

    }

    title(title_plot[i])

    legend("topleft",legend=names(r_stack)[1:2], 

           cex=1.2, col=t_col,lty=1,bty="n")

    legend("topright",legend=names(r_stack)[3], 

           cex=1.2, col=t_col[3],lty="dotted",bty="n")

    if (disp==TRUE){

      savePlot(file=paste(list_trans[[i]][2],".png",sep=""),type="png")

    }

    if (disp==FALSE){

      dev.off()

    }

    list_trans_data[[i]]<-trans_data

  }

  names(list_trans_data)<-names(list_trans)

  return(list_trans_data)

}

source(file.path(script_path,function_analyses_paper)) #source all functions used in this script.

nb_transect<-3

rast_pred_selected<-subset(rast_pred,c(1,6)) #3 is referring to FSS, plot it first because it has the

                                             # the largest range.

rast_pred2<-stack(rast_pred_selected,subset(s_raster,"elev_s"))

layers_names<-layerNames(rast_pred2)<-c("mod1","mod6","elev")

debug(plot_transect_m2)

#png(filename=paste("Comparison_daily_monthly_mean_lst",out_prefix,".png",sep=""),width=960,height=480)

#par(mfrow=c(1,2))

dev.off()

## Transects image location in OR             

png(paste("Fig7_elevation_transect_paths_",date_selected,out_prefix,".png", sep=""),

    height=480*layout_m[1],width=480*layout_m[2])

plot(elev_s)

#k<-1  #transect to plot

list_transect2[[3]]

#trans_file<-list_transect2[[k]][[1]]

#filename<-sub(".shp","",trans_file)             #Removing the extension from file.

#transect<-readOGR(".", filename)                 #reading shapefile 

#plot(transect,add=TRUE)

#title("Transect Oregon")

#dev.off()