library(gplots);

inputFile="GSM_T10_RINGED_WHTR_C0_DEGREES.txt";

a=read.csv(file=inputFile,header=TRUE,sep="\t")
datasetName=strsplit(inputFile,"[.]")[[1]][1]
datasetName=strsplit(datasetName ,"[/]")
datasetName=datasetName[[1]][length(datasetName[[1]])]
# start plotting with the city with highest frequency 
# to make plotting nicer.
#
max_freq=max(a$FREQ) 
max_degree=max(a$DEGREE)
max_city=a[ a$FREQ ==max_freq,3];

# make RING categorical variable
 a$RING<-factor(a$RING); 

city_list=levels(a$RING);
plot_colors <- c(
 rgb(r=0.0,g=0.0,b=0.9)
,rgb(r=0.0,g=0.9,b=0.0)
,rgb(r=0.0,g=0.9,b=0.9)
,rgb(r=0.9,g=0.0,b=0.0)
,rgb(r=0.9,g=0.0,b=0.9)
,rgb(r=0.9,g=0.9,b=0.0)
,rgb(r=0.9,g=0.9,b=0.9)
,rgb(r=0.6,g=0.3,b=0.0)
,rgb(r=0.3,g=0.6,b=0.9)
,"yellow","blue", "red", "forestgreen","magenta"

)


plot_colors <- colorpanel(length(city_list),low="green",mid="blue",high="red");

plot_colors_index=1;

values=matrix(nrow=length(city_list),ncol=max_degree)
  
isFirts=TRUE;
eCCDFA=data.frame()
for (city in head(city_list,n=1000)){
  city_total_freq= sum(a[ a$RING ==city,2]);
  x=a[ a$RING ==city,1];
  y=a[ a$RING ==city,2];
     y=a[ a$RING ==city,2]/city_total_freq;
     max_freq=1; 
  ds<-data.frame(x=x,y=y)
#, ylim=c(0,max_freq)
	if(isFirts) 
		plot(ds, xlab = "degree", ylab = "pdf of frequency",log="xy", type="l",xlim=c(1,max_degree)
#, ylim=c(1,max_freq)
,lty = plot_colors_index, lwd=2,col = plot_colors[plot_colors_index]) 
	 else  
	lines(ds,lty = 1, lwd=2,col = plot_colors[plot_colors_index%%length(plot_colors)]) 
 

	##new
	eCCDF=data.frame(FREQ=y)
	eCCDF$OP=paste("RING",city,sep="-")
	eCCDFA=rbind(eCCDFA,eCCDF)



       isFirts = FALSE;
  	city_total_freq= sum(a[ a$RING ==city,2]);
	print(paste(city,"processed.","total freq was",city_total_freq));

    values[plot_colors_index,1:length(x)]=y;

	plot_colors_index=plot_colors_index+1; 
if(plot_colors_index >2000) break;
}

legend("bottomleft", c("size 1-333","size 334-666","size 667-1000"), cex=0.8, col=c(plot_colors[plot_colors_index*1/6],plot_colors[plot_colors_index*3/6],plot_colors[plot_colors_index*5/6]), lty=1:plot_colors_index, lwd=2, bty="n"); 




eCCDFA$OP <- as.factor(eCCDFA$OP)
k=kruskal.test(FREQ~OP,data=eCCDFA)
k


## plot statistics on plot
rp = vector('expression',3)
#rp[1] = substitute(expression(MYVALUE), 
#		list(MYVALUE = k$method))[2]
rp[1] = substitute(expression(italic("chi-squared") == MYOTHERVALUE2), 
		list(MYOTHERVALUE2 = format(k$statistic, digits = 4)))[2]

rp[2] = substitute(expression(italic(df) == MYOTHERVALUE), 
		list(MYOTHERVALUE = format(k$parameter, digits = 4)))[2]
rp[3] = substitute(expression(italic(p.value) == MYOTHERVALUE), 
		list(MYOTHERVALUE = format(k$p.value, digits = 4)))[2]

legend('topright', legend = rp, bty = 'n')

 

for (city in city_list){
 # city_total_freq= sum(a[ a$RING ==city,2]);
  print(city)
  print(city_total_freq)
}
#library(rgl)
# open3d()
#rgl.surface(x=1:max_degree, z=1:length(city_list), y=values, color="red", alpha=0.75, back="lines",log="xy")
dev.copy(postscript,paste('results/',datasetName,'.eps',sep=''));
dev.off()
dev.copy(png,paste('results/',datasetName,'.png',sep=''));
dev.off()
dev.copy(pdf,paste('results/',datasetName,'.pdf',sep=''));
dev.off()

values;