infile<-"benchmark.tab"

outfile <- "all-2-all-plots.pdf";

outdir=".";

singleplots<-1
onlylegend<-1


#          1        2     3          4         5              6             
names <- c("family","dir","APSI Reference","SCI Reference","APSI LocARNA-P","SCI LocARNA-P",
#          7                       8                     9                  
           "Reliability1 Reference","Reliability1 LocARNA-P","Reliability1 Lara",
#          10                       11                     12                  
#           "Reliability2 reference","Reliability2 LocARNA-P","Reliability2 Lara", ## this is actually reliability 2
           "Reliability Reference","Reliability LocARNA-P","Reliability Lara", ## but call it reliability
#          13                       14                     15
           "Reliability3 Reference","Reliability3 LocARNA-P","Reliability3 Lara",
#          16                       17                     18
           "Reliability4 Reference","Reliability4 LocARNA-P","Reliability4 Lara",
#          19                    20                  
           "Dotplot-Correlation","Compalign LocARNA-P",
#          21                22          23
           "Compalign Lara", "SCI Lara", "APSI Lara",
#          24                     25
           "Reliability LocARNA by Reference", "Reliability Lara by Reference"
           )


## show all combinations of the following columns
columns <- c(20,3,4,5,6,23,22,
             #7,8,9,
             10,11,12,
             #13,14,15,
             #16,17,18,19,
             21,
             24,25
             )



pdfwidth<-5;
pdfheight<-5.4;

if (!singleplots) {
  pdf(outfile,width=2.5*(length(columns)-1),height=2.4*(length(columns)-1))
}


t <- read.table(infile)

family <- read.table(infile,sep="/")
family <- family[[2]]

t <- cbind(family,t)

## compute RATIOS OF RELIABILITIES
relratio<-t[[11]]/t[[10]];
relratio2<-t[[12]]/t[[10]];
#relratio<-t[[14]]/t[[13]];
#relratio2<-t[[15]]/t[[13]];


t <- cbind(t,relratio,relratio2)



families<-c("5_8S_rRNA","Hammerhead_3","Intron_gpII","sno_14q_I_II","U2","5S_rRNA","HCV_SLIV","IRES_HCV","SRP_bact","U6","Cobalamin","HCV_SLVII","IRES_Picorna","SRP_euk_arch","UnaL2","Entero_5_CRE","HepC_CRE","K_chan_RES","TAR","yybP-ykoY","Entero_CRE","Histone3","Lysine","T-box","Entero_OriR","HIV_FE","Retroviral_psi","THI","gcvT","HIV_GSL3","S_box","tRNA","Hammerhead_1","HIV_PBS","SECIS","U1")
families<-sort(families)

#                    c(bottom, left, top, right)
par(mar=c(4, 4, 2, 2) + 0.1)

if (!singleplots) {
  par(mfrow=c(length(columns)-1,length(columns)-1))
}

colors<-rainbow(length(families))
symbols<-1:6
symbols<-c(symbols,symbols)
symbols<-c(symbols,symbols)
symbols<-c(symbols,symbols)


if(!onlylegend) {

for (ii in 1:length(columns)) {
  for (jj in 1:length(columns)) {

    i<-columns[ii]
    j<-columns[jj]
        
    x<-t[[i]]
    y<-t[[j]]

    correlation=c(cor(x,y,method="pearson"),cor(x,y,method="kendall"),cor(x,y,method="spearman"))
    correlation <- format(correlation,digits=2)
    correlation <- paste(correlation[1],correlation[2],correlation[3]);      

    print(c(names[i],names[j],correlation));

    if (singleplots) {
      filename<-paste(names[i],"-",names[j],".pdf",sep="");
      filename<-gsub(" ","_",filename);
      filename<-gsub("/","_",filename);
      pdf(paste(outdir,filename,sep=""),height=pdfheight,width=pdfwidth);
    }
          
    c=0;
    for (f in families) {
      c=c+1;
      if (c==1) {
        plot(c(),c(),
             #main=correlation,
             xlab=names[i],ylab=names[j],xlim=c(min(x),max(x)),ylim=c(min(y),max(y)))
      }
      points(x[t[[1]]==f],y[t[[1]]==f],col=colors[c],pch=symbols[c],cex=0.7);
    }

    if (singleplots) {
      dev.off();
    }
  }
}

dev.off()
}

if (singleplots) {
  pdf(paste(outdir,"legend-horiz.pdf",sep=""),height=2.5,width=10.5);
  plot.new()
  legend("topleft",as.vector(families),col=colors,pch=symbols,cex=0.65,ncol=9)
  dev.off()
}