diff tools/mytools/cdf-old-not-used/cdf.py @ 0:9071e359b9a3

Uploaded

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/tools/mytools/cdf-old-not-used/cdf.py	Fri Mar 09 19:37:19 2012 -0500
@@ -0,0 +1,55 @@
+'''
+python wrapper for cdf plot
+'''
+
+import os,sys
+
+filename = sys.argv[1]
+output = sys.argv[2]
+c1, c2 = sys.argv[3].strip().split(',')
+label1,label2 = sys.argv[5].strip().split(',')
+title = sys.argv[6].strip()
+log = sys.argv[4]
+
+rf = open('tmp.r','w')
+
+rf.write("options(warn=-1)\n")
+
+rf.write("x <- read.table('"+filename+"')\n")
+
+rf.write("mycdf <- function(x1,x2) {\n")
+rf.write("       sx1 <- sort(x1)\n")
+rf.write("     sy1 <- c(1:length(sx1))/length(sx1)\n")
+rf.write("     sx2 <- sort(x2)\n")
+rf.write("     sy2 <- c(1:length(sx2))/length(sx2)\n")
+rf.write("       mi <- min(c(x1,x2))\n")
+rf.write("       ma <- max(c(x1,x2))\n")
+rf.write("       sx1 <- c(mi,sx1,ma)\n")
+rf.write("       sx2 <- c(mi,sx2,ma)\n")
+rf.write("       sy1 <- c(0,sy1,1)\n")
+rf.write("       sy2 <- c(0,sy2,1)\n")
+rf.write("       pv <- ks.test(x1,x2)$p.value\n")
+rf.write("       list(sx1,sy1,sx2,sy2,pv)\n")
+rf.write("}\n")
+
+rf.write("list1 <- x[,"+c1+"]\n")
+rf.write("list2 <- x[,"+c2+"]\n")
+
+if log == 'log':
+    rf.write("list1 <- log2(1+list1) \n")
+    rf.write("list2 <- log2(1+list2) \n")
+    #rf.write("list1[list1 == -Inf] <- min(list1[list1 != -Inf])\n")
+rf.write("res <- mycdf(list1,list2)\n")
+
+rf.write("pdf('"+output+"')\n")
+rf.write("plot(res[[1]],res[[2]],type='l',lty=1,lwd=2,col='red',main=paste('"+title+"',format(res[[5]], digit=2,scientific = T),sep=' , p='),xlab='value',ylab='cumulative frequency')\n")
+rf.write("lines(res[[3]],res[[4]],type='l',lty=2,lwd=2,col='blue')\n")
+rf.write("label1 <- paste('"+label1+"','( n1=',length(list1),')')\n")
+rf.write("label2 <- paste('"+label2+"','( n2=',length(list2),')')\n")
+rf.write("legend('bottomright',c(label1,label2),col=c('red','blue'), lty=1:2, lwd=2, bty='n')\n")
+rf.write("boxplot(list1,list2,names=c(label1,label2))\n")
+rf.write("dev.off()\n")
+
+rf.close()
+os.system("R --vanilla < tmp.r")    
+os.system('rm tmp.r')