Mercurial > repos > rlegendre > ving

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/ving.R	Tue Jul 29 03:15:34 2014 -0400
@@ -0,0 +1,1740 @@
+#! /usr/bin/Rscript
+
+# Il est interdit de mettre ARGUMENT comme argument long
+
+
+
+# spec contains at least 4 columns, as many as 5 columns.
+#
+# column 1: long name of flag
+#
+# column 2: short name of flag
+#
+# column 3: argument flag.  0=without argument, 1=with argument
+#
+# column 4: mode of argument.  one of "logical", "integer", "double",
+# "complex", "character", "" for arguments without value (return logical)
+#
+# column 5 (optional): description of ARGUMENT.
+#
+
+
+RetreiveAndDestroy=function(opt,root,stem,regexp,SearchNames,Out,isValue,NbFound,StockNames){
+	Bool=lapply(paste(root,SearchNames,stem,sep=""),grepl,opt)
+	names(Bool)=StockNames
+	Pos=lapply(Bool,which)
+	names(Pos)=StockNames
+	disable=c()
+	for (i in StockNames){
+		nbmatch=length(Pos[[i]])
+		if(nbmatch>0){
+			NbFound[[i]]=NbFound[[i]]+nbmatch
+			disable=c(disable,-1*Pos[[i]])
+			if(is.null(Out[[i]])){
+				if(isValue[[i]]!=0){
+					if(regexp=="next"){
+						Out[[i]]=opt[Pos[[i]]+1]
+						disable=c(disable,-1*(Pos[[i]]+1))
+					}else{
+						Out[[i]]=sub(regexp,"\\1",opt[Pos[[i]]])
+					}
+				}else{
+					Out[[i]]=TRUE
+				}
+			}else{
+				if(isValue[[i]]!=0){
+					if(regexp=="next"){
+						Out[[i]]=c(Out[[i]],opt[Pos[[i]]+1])
+						disable=c(disable,-1*(Pos[[i]]+1))
+					}else{
+						Out[[i]]=c(Out[[i]],sub(regexp,"\\1",opt[Pos[[i]]]))
+					}
+				}else{
+					Out[[i]]=c(Out[[i]],TRUE)
+				}
+			}
+		}
+	}
+	if(length(disable)>0){
+		opt=opt[disable]
+	}
+	Out[["ARGUMENT"]]=list()
+	Out[["ARGUMENT"]][["opt"]]=opt
+	Out[["ARGUMENT"]][["NbFound"]]=NbFound
+	return(Out)
+}
+
+getopt = function (spec=NULL,opt=commandArgs()) {
+
+	FindArgs=which(opt=="--args")
+	if(length(FindArgs)!=1){
+		stop(length(FindArgs)," --args found where 1 expected.",call.=F)
+	}
+	ExecName=sub("--file=","",opt[FindArgs-1])
+
+	if(FindArgs<length(opt)){
+		opt=opt[(FindArgs+1):length(opt)]
+	}else{
+		opt=""
+	}
+
+
+	min.columns=5
+	colNames=c("LongName","ShortName","Flag","Mod","Default")
+	max.columns=6
+	DimSpec=dim(spec)
+	if(DimSpec[2]>min.columns){
+		colNames=c(colNames,"Description")
+	}
+
+	if(is.null(spec) | !is.matrix(spec) | (DimSpec[2]<min.columns | DimSpec[2]>max.columns)){
+		stop('argument "spec" is required and must be a matrix with 4|5 columns.',call.=F)
+	}
+	colnames(spec)=colNames
+
+	spec=as.data.frame(spec,stringsAsFactors=F)
+	#spec validation
+	if(length(unique(c(spec$ShortName,"ARGUMENT","args")))!=DimSpec[1]+2 | length(unique(spec$LongName))!=DimSpec[1]){
+		stop('Long|Short names for flags must be unique (Long name : "ARGUMENT" and "args" forbidden).',
+				"\n","List of duplicated :",
+				"\n","Short: ",paste(spec$ShortName[duplicated(c(spec$ShortName,"ARGUMENT","args"))],collapse=" "),
+				"\n","Long:  ",paste(spec$ShortName[duplicated(spec$LongName)],collapse=" "),call.=F)
+	}
+	if(length(which(nchar(spec$ShortName)>1))!=0){
+		stop('Short names flags can\'t be longer than 1 character.')
+	}
+
+
+	#initialize
+	Out=list()
+	Short2Long=list()
+	NbFound=list()
+	isValue=list()
+	for (i in 1:DimSpec[1]){
+		Short2Long[[spec$ShortName[i]]]=spec$LongName[i]
+		NbFound[[spec$LongName[i]]]=0
+		isValue[[spec$LongName[i]]]=spec$Flag[i]
+	}
+
+	#Map, retreive and suppress ARGUMENTs and arguments
+	#Value ARGUMENT --example=value
+	Out=RetreiveAndDestroy(opt,"^--","=.+$",".+=(.+)$",spec$LongName,Out,isValue,NbFound,spec$LongName)
+	opt=Out[["ARGUMENT"]][["opt"]]
+	NbFound=Out[["ARGUMENT"]][["NbFound"]]
+	Out[["ARGUMENT"]]=NULL
+	#boolean ARGUMENT --example
+	Out=RetreiveAndDestroy(opt,"^--","$","$",spec$LongName,Out,isValue,NbFound,spec$LongName)
+	opt=Out[["ARGUMENT"]][["opt"]]
+	NbFound=Out[["ARGUMENT"]][["NbFound"]]
+	Out[["ARGUMENT"]]=NULL
+	#short name ARGUMENT -t value OR boolean -t
+	Out=RetreiveAndDestroy(opt,"^-","$","next",spec$ShortName,Out,isValue,NbFound,spec$LongName)
+	opt=Out[["ARGUMENT"]][["opt"]]
+	NbFound=Out[["ARGUMENT"]][["NbFound"]]
+	Out[["ARGUMENT"]]=NULL
+	#Warn about non mapped ARGUMENTs
+	if(length(opt)>0){
+		PosUnkArg=which(grepl("^-",opt))
+		if(length(PosUnkArg)){
+			message("Error, argument unreconized :","\n",paste(opt[PosUnkArg],collapse="\n"),"\n\n")
+		}
+		if(length(PosUnkArg)>0){
+			opt=opt[PosUnkArg*-1]
+		}
+	}
+	#Arguments
+	Out[["ARGUMENT"]]=opt
+
+	#Validation of ARGUMENTs
+	for(i in 1:DimSpec[1]){
+		if(spec$Flag[i]=="0"){#verify boolean arguments
+			NbValue=length(Out[[spec$LongName[i]]])
+			if(NbValue>1){
+				message("Warning : ",spec$LongName[i]," found ",NbValue," times")
+			}
+		}
+		if(length(Out[[spec$LongName[i]]])==0){
+			Out[[spec$LongName[i]]]=spec$Default[i]
+		}
+		library("methods")
+		Out[[spec$LongName[i]]]=as(Out[[spec$LongName[i]]],spec$Mod[i])
+	}
+
+	return(Out)
+}
+
+# column 3: argument flag.  0=no argument, 1=required argument, 2=optional argument
+
+
+
+### =========================================================================
+### Ving's methods
+### -------------------------------------------------------------------------
+###
+
+##Converts the flag numbers into binary
+integer.base.b = function(x, b=2){
+  xi <- as.integer(x)
+  if(any(is.na(xi) | ((x-xi)!=0)))
+    print(list(ERROR="x not integer", x=x))
+  N <- length(x)
+  xMax <- max(c(x,1))
+  ndigits <- 11
+  Base.b <- array(NA, dim=c(N, ndigits))
+  for(i in 1:ndigits){#i <- 1
+    Base.b[, ndigits-i+1] <- (x %% b)
+    x <- (x %/% b)
+  }
+  Base.b
+}
+
+##Checks if color used is an acceptable color
+is.acceptable.color = function(character) {
+  tmp=try(col2rgb(character),TRUE)
+  return(class(tmp)!="try-error")
+}
+
+##Combines two data frame
+listOfDataFrame2DataFrame = function(x,vertical=TRUE) {
+  if(length(x)==1) {
+    x[[1]]
+  }
+  else {
+    if(vertical) {
+      ncol = dim(x[[1]])[2]
+      vect=unlist(lapply(x,t))
+      retour = as.data.frame(t(array(vect,dim=c(ncol,length(vect)/ncol))),stringsAsFactors=FALSE)
+      names(retour)=names(x[[1]])
+    }
+    else {
+      nline=dim(x[[1]])[1]
+      vect=unlist(x)
+      retour=as.data.frame(array(vect,dim=c(nline,length(vect)/nline)))
+      names(retour) = unlist(lapply(x,names))
+    }
+    retour
+  }
+}
+
+# utilise la fonction scanBam d'une façon plus partique pour moi
+extractFromBam = function(file,which,what) {
+  return(scanBam(file, param=ScanBamParam(which=which,what=what))[[1]][[1]])
+}
+
+##Returns a list from cigar expression
+interpreteCIGAR = function(cigar) {
+  cigar_un = strsplit(unique(cigar),split="")
+  n_cigar_un = length(cigar_un)
+  taille_cigar = list()
+  analise_cigar = function(cigar_) {
+    cigar_sortie = list()
+    acc = ""
+    for(j in 1:length(cigar_)) {
+      if(sum(cigar_[j]==as.character(0:9))==1) {
+	acc=paste(acc,cigar_[j],sep="")
+      }
+      else {
+	cigar_sortie[[length(cigar_sortie)+1]]=as.integer(acc)
+	cigar_sortie[[length(cigar_sortie)+1]]=cigar_[j]
+	acc=""
+      }
+    }
+    return(cigar_sortie)
+  }
+  cigar_interprete = lapply(cigar_un,analise_cigar)
+  names(cigar_interprete) = unique(cigar)
+
+  return(cigar_interprete)
+}
+
+# prend un CIGAR splités et retourne la taille occupé par le read sur la séquence génomique (introns compris)
+calcule_longueur_cigar = function(cigar_) {
+  lon = 0
+  N=length(cigar_)
+  for(j in seq(2,N,2)) {
+    if(cigar_[[j]]!="I") {
+      lon=lon+cigar_[[j-1]]
+    }
+  }
+  return(lon)
+}
+
+# prend un CIGAR splités et retourne les positions
+calcule_junction_cigar = function(cigar_) {
+  retour=list()
+  lon = 0
+  N=length(cigar_)
+  for(j in seq(2,N,2)) {
+    if(cigar_[[j]]!="I") {
+      lon=lon+cigar_[[j-1]]
+    }
+    if(cigar_[[j]]=="N") {
+      retour[[length(retour)+1]]=c(lon-cigar_[[j-1]]+1,lon)
+    }
+  }
+  return(retour)
+}
+
+##Returns a list of numbers of single read with their coordinates
+compresse_coordonnees = function(debut,fin) {
+  if(length(debut)==0) {
+    return(list(numeric(),numeric(),numeric()))
+  }
+  else {
+    tmp=sort(paste(debut,fin,sep="_"))
+    tmp_rle=rle(tmp)
+    poids=tmp_rle$lengths
+    values_split=strsplit(tmp_rle$values,split="_")
+    doit = function(j) {
+      return(as.integer(values_split[[j]][1]))
+    }
+    debut_uni=sapply(1:length(poids),doit)
+    doit = function(j) {
+      return(as.integer(values_split[[j]][2]))
+    }
+    fin_uni=sapply(1:length(poids),doit)
+    ordre_debut = order(debut_uni)
+    return(list(debut_uni[ordre_debut],fin_uni[ordre_debut],poids[ordre_debut]))
+  }
+}
+
+RDataFileName = function(file) {
+  return(paste(file,".RData",sep=""))
+}
+
+##Function converts and extracts the infos from bamfile
+readBam_ <-function(file_,insert_max_=2000,stranded_=TRUE,ncore_=1,libraryType_=c("standard","inverse"),fileNameRData_=NA,normalized_=NULL,chrName_=NULL,from_=1,to_=NULL) {
+  suppressPackageStartupMessages(require("Rsamtools"))
+  suppressPackageStartupMessages(require("GenomicRanges"))
+##Declaration of variables
+  flagstat = numeric(11)
+  names(flagstat)=c("total","duplicates","mapped","paired","read1","read2","properly paired","itself and mate mapped","singletons","mate mapped on a different chr","QC-failed")
+  genome_info=scanBamHeader(file_)[[1]]$targets
+  noms_chromosomes = names(genome_info)
+  longueur_chromosomes = as.integer(genome_info)
+  nombre_chromosomes = length(noms_chromosomes)
+  brin_F = list()
+  brin_R = list()
+  brin_F_junction = list()
+  brin_R_junction = list()
+  pas<-c(1,2,6,7)
+  i_zone = 0
+  if(is.null(chrName_)) {
+    chrName__ = noms_chromosomes
+  }
+  else {
+    chrName__ = chrName_
+  }
+##Fragments identification
+  for(i in (1:nombre_chromosomes)) {
+    i_zone = i_zone +1
+    nom_chromo = noms_chromosomes[i]
+    lon_chromo = longueur_chromosomes[i]
+
+    if(!(nom_chromo %in% chrName__)) {
+      brin_F[[i]] = list(numeric(),numeric(),numeric())
+      brin_R[[i]] = list(numeric(),numeric(),numeric())
+      brin_F_junction[[i]] = list(numeric(),numeric(),numeric())
+      brin_R_junction[[i]] = list(numeric(),numeric(),numeric())
+    }
+    else {
+      if(is.null(to_)) {
+	to_i=lon_chromo
+      }
+      else {
+	to_i=to_[min(i_zone,length(to_))]
+      }
+      from_i=from_[min(i_zone,length(from_))]
+
+      commande = paste("RangesList(`",nom_chromo,"` = IRanges(",from_i,",",to_i,"))",sep="")
+      expr=try(parse(text=commande),TRUE)
+##Function used from GenomicRanges package
+      which = eval(expr)
+      what<-c("flag","mpos","cigar","mrnm","isize")
+      param<- ScanBamParam(what=what, which=which)
+      G_alignement<- readGAlignmentsFromBam(file=file_, param=param, use.names=TRUE)
+##Case of no reads on the chromosome
+      start=start(G_alignement)
+      if(length(start) == 0 ){
+	brin_F[[i]] = list(numeric(),numeric(),numeric())
+	brin_R[[i]] = list(numeric(),numeric(),numeric())
+	brin_F_junction[[i]] = list(numeric(),numeric(),numeric())
+	brin_R_junction[[i]] = list(numeric(),numeric(),numeric())
+      }
+      else {
+	end=end(G_alignement)
+	strand=as.character(strand(G_alignement))
+	flag=mcols(G_alignement)$flag
+	mpos=mcols(G_alignement)$mpos
+	cigar=mcols(G_alignement)$cigar
+	mrnm=mcols(G_alignement)$mrnm
+	isize=mcols(G_alignement)$isize
+
+	first_read=integer.base.b(flag)[,5]==1
+	strand[strand =="+" & !first_read ]="*"
+	strand[strand =="-" & !first_read ] ="+"
+	strand[strand =="*" & !first_read ] ="-"
+
+##Case of pairend reads
+	is_on_same_chr = mrnm==nom_chromo
+	is_on_same_chr[is.na(is_on_same_chr)] = FALSE
+	is_paired = is_on_same_chr & abs(isize) <= insert_max_
+	is_paired[first_read & strand=="+" & (isize<0 | isize>insert_max_)] = FALSE
+	is_paired[!first_read & strand=="+" & (isize>0 | isize < -insert_max_)] = FALSE
+	is_paired[first_read & strand=="-" & (isize>0 | isize < -insert_max_)] = FALSE
+	is_paired[!first_read & strand=="-" & (isize<0 | isize>insert_max_)] = FALSE
+	is_paired[is.na(is_paired)] = FALSE
+
+	debut_fragment_paired_plus<-mpos[!first_read & strand =="+" & is_paired]
+	fin_fragment_paired_plus<-end[!first_read & strand=="+" & is_paired]
+	debut_fragment_paired_moins<-mpos[first_read & strand=="-" & is_paired]
+	fin_fragment_paired_moins<-end[first_read & strand =="-" & is_paired]
+
+##Case of single reads
+	debut_fragment_singleton_plus = start[!is_paired & strand =="+"]
+	fin_fragment_singleton_plus = end[!is_paired & strand =="+"]
+	debut_fragment_singleton_moins=start[!is_paired & strand =="-"]
+	fin_fragment_singleton_moins= end[!is_paired & strand =="-"]
+
+##Fragments
+	debut_frag_plus = c(debut_fragment_paired_plus,debut_fragment_singleton_plus)
+	fin_frag_plus  = c(fin_fragment_paired_plus,fin_fragment_singleton_plus)
+	debut_frag_moins=c(debut_fragment_paired_moins,debut_fragment_singleton_moins)
+	fin_frag_moins = c(fin_fragment_paired_moins,fin_fragment_singleton_moins)
+	brin_F[[i]] = compresse_coordonnees(debut_frag_plus,fin_frag_plus)
+	brin_R[[i]] = compresse_coordonnees(debut_frag_moins,fin_frag_moins)
+
+##CIGAR's interpreter
+	cigar_interprete = interpreteCIGAR(cigar)
+	longueur_cigar = lapply(cigar_interprete,calcule_longueur_cigar)
+	junction_cigar = lapply(cigar_interprete,calcule_junction_cigar)
+	debut_junction = numeric()
+	fin_junction = numeric()
+	brin_junction = numeric()
+	i_junction = 0
+##junction read
+	for(j in 1:length(cigar)) {
+	  junctions_ = junction_cigar[[cigar[j]]]
+	  if(length(junctions_)){
+	    for(k in 1:length(junctions_)) {
+	      i_junction = i_junction + 1
+	      debut_junction[i_junction] = start[j] + junctions_[[k]][1] - 1
+	      fin_junction[i_junction] = start[j] + junctions_[[k]][2] - 1
+	      brin_junction[i_junction] = strand[j]
+	    }
+	  }
+	}
+	if(i_junction==0) {
+	  brin_F_junction[[i]] = list(numeric(),numeric(),numeric())
+	  brin_R_junction[[i]] = list(numeric(),numeric(),numeric())
+	}
+	else {
+	  brin_F_junction[[i]] = compresse_coordonnees(debut_junction[brin_junction=="+"],fin_junction[brin_junction=="+"])
+	  brin_R_junction[[i]] = compresse_coordonnees(debut_junction[brin_junction=="-"],fin_junction[brin_junction=="-"])
+	}
+##Flagstat interpreter
+	flag_bits = integer.base.b(flag)#remplie les  données stat pour un flag donné
+
+	##flagstat
+	##total
+	flagstat[1] = flagstat[1] + sum(flag_bits[,2]==0)
+	##duplicates
+	flagstat[2] = flagstat[2] + sum((flag_bits[,1]==1)&(flag_bits[,2]==0))
+	##mapped
+	flagstat[3] = flagstat[3] + sum((flag_bits[,9]==0)&(flag_bits[,2]==0))
+	##paired
+	flagstat[4] = flagstat[4] + sum((flag_bits[,11]==1)&(flag_bits[,2]==0))
+	##read1
+	flagstat[5] = flagstat[5] + sum((flag_bits[,5]==1)&(flag_bits[,2]==0))
+	##read2
+	flagstat[6] = flagstat[6] + sum((flag_bits[,4]==1)&(flag_bits[,2]==0))
+	##iself and mate mapped
+	flagstat[8] = flagstat[8] + sum((flag_bits[,11]==1)&(flag_bits[,9]==0)&(flag_bits[,8]==0)&(flag_bits[,2]==0))
+	##singletons
+	flagstat[9] = flagstat[9] + sum((flag_bits[,8]==1)&(flag_bits[,2]==0))
+	##QC-failed
+	flagstat[11] = flagstat[11] + sum(flag_bits[,2]==1)
+	##flagstat
+	##mate on a different chr
+	flagstat[10] = flagstat[10] + sum((!is_on_same_chr)&(flag_bits[,11]==1)&(flag_bits[,9]==0)&(flag_bits[,8]==0)&(flag_bits[,2]==0))
+	##flagstat
+	##properly paired
+	flagstat[7] = flagstat[7] + sum(is_paired)
+
+      }
+    }
+  }
+##Data storing
+  names(brin_F) = noms_chromosomes
+  names(brin_R) = noms_chromosomes
+  names(brin_F_junction) = noms_chromosomes
+  names(brin_R_junction) = noms_chromosomes
+
+  bamHandler = list()
+  if(libraryType_[1]=="inverse") {
+    bamHandler[[1]] = brin_R
+    bamHandler[[2]] = brin_F
+  }
+  else {
+    bamHandler[[1]] = brin_F
+    bamHandler[[2]] = brin_R
+  }
+  bamHandler[[3]] = longueur_chromosomes
+  bamHandler[[4]] = flagstat
+  bamHandler[[5]] = stranded_
+
+  if(libraryType_[1]=="inverse") {
+    bamHandler[[6]] = brin_R_junction
+    bamHandler[[7]] = brin_F_junction
+  }
+  else {
+    bamHandler[[6]] = brin_F_junction
+    bamHandler[[7]] = brin_R_junction
+  }
+ bamHandler[[8]] = FALSE
+if(!is.null(normalized_)) {
+    for( i in pas) {
+	  for(j in 1:nombre_chromosomes){
+	    bamHandler[[i]][[j]][[4]]=normalized_*bamHandler[[i]][[j]][[3]]
+	  }
+    }
+    bamHandler[[8]] = TRUE
+}
+
+  names(bamHandler) = c("F","R","chrLength","flagstat","stranded","junctions_F","junctions_R","norm")
+  if((is.null(chrName_))&(from_==1)&(is.null(to_))) {
+    if(is.null(fileNameRData_)|is.na(fileNameRData_)) {
+	save(bamHandler,file=RDataFileName(file_))
+    }
+    else {
+	save(bamHandler,file=fileNameRData_)
+    }
+  }
+  return(bamHandler)
+}
+
+chercheBamHandlerDansFichierRData = function(fichierRData) {
+  tmp=try(load(fichierRData),TRUE)
+  if(class(tmp)=="try-error") {
+    return(tmp)
+  }
+  else {
+    return(try(bamHandler,TRUE))
+  }
+}
+
+##Public Function
+readBam = function(file,insert_max=2000,stranded=TRUE,reload=FALSE,ncore=1,libraryType=c("standard","inverse"),normalized=NULL,chrName=NULL,from=1,to=NULL) {
+  if(length(file)==1) {
+    nom_fichier_RData = RDataFileName(file)
+    if(!file.exists(nom_fichier_RData)|reload) {
+      return(try(readBam_(file_=file,insert_max_=insert_max,stranded_=stranded,ncore_=ncore,libraryType_=libraryType,normalized_=normalized,chrName_=chrName,from_=from,to_=to),TRUE))
+    }
+    else {
+      load(nom_fichier_RData)
+      return(bamHandler)
+    }
+  }
+  else {
+    require("multicore")
+    chargeBamsDansFichiersRData = function(oneFile) {
+      nom_fichier_RData = RDataFileName(oneFile)
+      if(!file.exists(nom_fichier_RData)|reload) {
+	bamHandler=try(readBam_(oneFile,insert_max_=insert_max,stranded_=stranded,libraryType_=libraryType,normalized_=normalized,chrName_=chrName,from_=from,to_=to),TRUE)
+      }
+      return(nom_fichier_RData)
+    }
+    fichiersRData = mclapply(file,chargeBamsDansFichiersRData,mc.cores=ncore)
+    gc()
+    return(lapply(fichiersRData,chercheBamHandlerDansFichierRData))
+  }
+}
+
+##Extracts the signal from bamHandler objects
+extractSignal = function(bamHandlerList,chrName,from=1, to=NULL,normalized_=FALSE) {
+  forward=list()
+  reverse=list()
+  chr=which(names(bamHandlerList[[1]]$F)==chrName)
+  if(is.null(to)){
+    to=bamHandlerList[[1]]$chrLength[chr]
+  }
+  if(normalized_) {
+    end=4
+  }else{
+    end=3
+  }
+  for(i in 1:length(bamHandlerList)){
+      forward_=numeric(to-from+1)
+      which_read = which((bamHandlerList[[i]]$F[[chrName]][[2]]>=from)&(bamHandlerList[[i]]$F[[chrName]][[1]]<=to))
+      n_reads=length(which_read)
+      if(n_reads>0) {
+      for(k in which_read) {
+	debut_read = max(1,bamHandlerList[[i]]$F[[chrName]][[1]][k]-from+1)
+	fin_read = min(bamHandlerList[[i]]$F[[chrName]][[2]][k]-from+1,to-from+1)
+	forward_[debut_read:fin_read]=forward_[debut_read:fin_read]+bamHandlerList[[i]]$F[[chrName]][[end]][k]
+      }
+    }
+  which_junctions = which((bamHandlerList[[i]]$junctions_F[[chrName]][[2]]>=from)&(bamHandlerList[[i]]$junctions_F[[chrName]][[1]]<=to))
+  n_junctions=length(which_junctions)
+  if(n_junctions>0) {
+      for(k in which_junctions) {
+	debut_junction = max(1,bamHandlerList[[i]]$junctions_F[[chrName]][[1]][k]-from+1)
+	fin_junction = min(bamHandlerList[[i]]$junctions_F[[chrName]][[2]][k]-from+1,to-from+1)
+	forward_[debut_junction:fin_junction]=forward_[debut_junction:fin_junction]-bamHandlerList[[i]]$junctions_F[[chrName]][[end]][k]
+      }
+    }
+
+    reverse_=numeric(to-from+1)
+    which_read = which((bamHandlerList[[i]]$R[[chrName]][[2]]>=from)&(bamHandlerList[[i]]$R[[chrName]][[1]]<=to))
+    n_reads= length(which_read)
+  if(n_reads>0) {
+    for(k in which_read) {
+      debut_read = max(1,bamHandlerList[[i]]$R[[chrName]][[1]][k]-from+1)
+      fin_read = min(bamHandlerList[[i]]$R[[chrName]][[2]][k]-from+1,to-from+1)
+      reverse_[debut_read:fin_read]=reverse_[debut_read:fin_read]+bamHandlerList[[i]]$R[[chrName]][[end]][k]
+      }
+    }
+  which_junctions = which((bamHandlerList[[i]]$junctions_R[[chrName]][[2]]>=from)&(bamHandlerList[[i]]$junctions_R[[chrName]][[1]]<=to))
+  n_junctions=length(which_junctions)
+  if(n_junctions>0) {
+    for(k in which_junctions) {
+	debut_junction = max(1,bamHandlerList[[i]]$junctions_R[[chrName]][[1]][k]-from+1)
+	fin_junction = min(bamHandlerList[[i]]$junctions_R[[chrName]][[2]][k]-from+1,to-from+1)
+	reverse_[debut_junction:fin_junction]=reverse_[debut_junction:fin_junction]-bamHandlerList[[i]]$junctions_R[[chrName]][[end]][k]
+      }
+    }
+    forward_[forward_<0]=0
+    reverse_[reverse_<0]=0
+    if(bamHandlerList[[i]]$stranded) {
+      forward[[i]]=forward_
+      reverse[[i]]=reverse_
+    }
+    else {
+      forward[[i]]=forward_+reverse_
+      reverse[[i]]=numeric(to-from+1)
+    }
+    }
+  chr_ = list()
+  chr_$F = forward
+  chr_$R = reverse
+  return(chr_)
+}
+
+##Intern function for readGff function
+my.read.lines2=function(fname) {
+
+  s = file.info( fname )$size
+  buf = readChar( fname, s, useBytes=T)
+  strsplit( buf,"\n",fixed=T,useBytes=T)[[1]]
+}
+
+##Extracts the annotation infos from Gff file
+readGff = function(file_in, from=1, to=Inf, chr=NULL, infoName=c("ID","Name","Parent","gene","Alias","orf_classification","Ontology_term","Note","GO")) {
+  tmp=try(my.read.lines2(file_in))
+  if(!is.null(chr)) {
+    tmp1=grep(chr, tmp, value=TRUE,useBytes=T)
+  }
+  else {
+    tmp1=tmp
+  }
+  N=length(tmp1)
+  Chr = array()
+  Start= array()
+  Stop= array()
+  Strand=array()
+  Type=array()
+  info=list()
+  for(i in 1:length(infoName)) info[[i]] = array()
+  names(info)=infoName
+  j<-1
+  for (i in 1:N) {
+    if(substr(tmp1[i],1,1)!="#") {
+      line_split=unlist(strsplit(tmp1[i],"\t",fixed=T,useBytes=T))
+      if((as.integer(line_split[4])<=to) & (as.integer(line_split[5])>=from)) {
+	Chr[j] = line_split[1]
+	Start[j] = as.integer(line_split[4])
+	Stop[j] = as.integer(line_split[5])
+	Strand[j]=line_split[7]
+	Type[j] = line_split[3]
+	ninth=unlist(strsplit(line_split[9],";",fixed=T,useBytes=T))
+	element_ninth_empty=rep(TRUE,length(infoName))
+	for(element_ninth in ninth) {
+	  element_ninth_split = unlist(strsplit(element_ninth,"=",fixed=T,useBytes=T))
+	  if(length(element_ninth_split)==2) {
+	    if(element_ninth_split[1] %in% infoName) {
+	      info[[element_ninth_split[1]]][j]=element_ninth_split[2]
+	      element_ninth_empty[infoName==element_ninth_split[1]]=FALSE
+	    }
+	  }
+	}
+	for(infoName_ in infoName[element_ninth_empty]) {
+	  info[[infoName_]][j]="."
+	}
+	j<-j+1
+      }
+    }
+  }
+  retour = data.frame(Chr,Type,Start,Stop,Strand,info,stringsAsFactors = FALSE)
+  return(retour)
+}
+
+##Returns the classic visualisation
+plotRNAseq = function(forward,reverse,debut_vue=1,fin_vue=length(forward),chr=NULL,annot=NULL,style=NULL,top=NULL,bottom=NULL,x="",y="",titre="",repeated=FALSE,name_flags="",decal=0,ataxises=NULL,classic_plus_color="navyblue",classic_minus_color="mediumvioletred",stranded=TRUE) {
+  if(repeated) {
+    forward_=numeric(fin_vue)
+    forward_[debut_vue:fin_vue]=forward
+    reverse_=numeric(fin_vue)
+    reverse_[debut_vue:fin_vue]=reverse
+  }
+  else {
+    forward_=forward
+    reverse_=reverse
+  }
+  if(is.null(top)) {
+    top=max(forward_[debut_vue:fin_vue])
+  }
+  if(is.null(bottom)) {
+    bottom=max(reverse_[debut_vue:fin_vue])
+  }
+  if(is.null(ataxises)) {
+    plot(c(debut_vue,fin_vue)+decal,c(-bottom,top),ylim=c(-bottom,top),xlab=x,ylab=y,main=titre,col="white",xaxs="i",cex.main=2,yaxt="n",cex.lab=1.8)
+  }
+  else {
+    plot(c(debut_vue,fin_vue)+decal,c(-bottom,top),ylim=c(-bottom,top),xlab=x,ylab=y,main=titre,col="white",xaxs="i",xaxt="n",cex.main=2,yaxt="n",cex.lab=1.8)
+    ataxisesLabels=as.character(ataxises)
+    ataxisesLabels[((1:length(ataxises))%%2)==0]=""
+    lim=c(-bottom,top)
+    ataxises_y=pretty(lim,n=4)
+    ataxisesLabels_y=as.character(abs(ataxises_y))
+    axis(1,at=ataxises,labels=FALSE,cex.axis=2)
+    axis(1,at=ataxises,labels=ataxisesLabels,cex.axis=2,line=0.4,lwd=0)
+    axis(2, at=ataxises_y,labels=FALSE,cex.axis=2)
+    axis(2, at=ataxises_y,labels=ataxisesLabels_y,cex.axis=2,line=-0.4,lwd=0)
+  }
+  polygon(c(debut_vue,debut_vue:fin_vue,fin_vue)+decal,c(0,forward_[debut_vue:fin_vue],0),col=classic_plus_color,border=NA)
+  if(stranded){
+    text(fin_vue+(fin_vue-debut_vue)*0.01,top/2,"+",xpd=NA,cex=3)
+    text(fin_vue+(fin_vue-debut_vue)*0.01,-bottom/2,"-",xpd=NA,cex=3)
+    text(fin_vue+(fin_vue-debut_vue)*0.025,(top-bottom)/2,"Strand",xpd=NA,cex=2,srt=-90)
+    polygon(c(debut_vue,debut_vue:fin_vue,fin_vue)+decal,c(0,-reverse_[debut_vue:fin_vue],0),col=classic_minus_color,border=NA)
+    abline(h=0,lwd=2)
+    abline(h=0,col="white")
+  }
+  #if(!is.null(chr)&!is.null(annot)&!is.null(style)) {
+   # N=dim(annot)[1]
+    #if(N>0) {
+     # for(i in 1:N) {
+	#if(annot$Strand[i]=="+") {
+	 # arrows(annot$Start[i],top,annot$Stop[i],top,col=style$col[style$Type==annot$Type[i]],length=0.10)
+	#}
+	#if(annot$Strand[i]=="-") {
+	 # arrows(annot$Stop[i],-bottom,annot$Start[i],-bottom,col=style$col[style$Type==annot$Type[i]],length=0.10)
+	#}
+	#if(annot$Strand[i]==".") {
+	 # segments(annot$Start[i],top,annot$Stop[i],top,col=style$col[style$Type==annot$Type[i]])
+	  #segments(annot$Stop[i],-bottom,annot$Start[i],-bottom,col=style$col[style$Type==annot$Type[i]])
+	#}
+      #}
+    #}
+  #}
+  if(name_flags!="") {
+      flags=try(get(name_flags),TRUE)
+      if(class(flags)!="try-error") {
+	f_=flags[(flags$Chr==chr)&(flags$Stop>=debut_vue)&(flags$Start<=fin_vue),]
+	N=dim(f_)[1]
+	points(f_$Start,rep(0,N),col=2,pch=19,cex=2)
+      }
+   }
+   if(decal<=0) {
+    lines(c(0,0),c(top,-bottom),col=1)
+   }
+}
+
+##Intern function for heatmap visualisation
+paletteFromColors = function(colMin="blue",colMax="red",n=300,method=c("hsv","rgb")) {
+  colMinRGB = col2rgb(colMin)[,1]/255
+  colMaxRGB = col2rgb(colMax)[,1]/255
+  seqList = list()
+  if(method[1]=="rgb") {
+    for(i in 1:3) {
+      seqList[[i]]=seq(colMinRGB[i],colMaxRGB[i],length.out=n)
+    }
+    return(rgb(seqList[[1]],seqList[[2]],seqList[[3]]))
+  }
+  else {
+    colMinHSV = rgb2hsv(colMinRGB)
+    colMaxHSV = rgb2hsv(colMaxRGB)
+    for(i in 1:3) {
+      seqList[[i]]=seq(colMinHSV[i],colMaxHSV[i],length.out=n)
+    }
+    return(hsv(seqList[[1]],s=seqList[[2]],v=seqList[[3]]))
+  }
+}
+
+##Retuns the heatmap visualisation
+myHeatMap = function(data,debut_vue,fin_vue,ataxises=NULL,lim=NULL,heatmap_max_color="#000055",heatmap_min_color="#FFFFAA",heatmap_palette_method="hsv",textOnLeftSide="") {
+  palette=paletteFromColors(heatmap_min_color,heatmap_max_color,method=heatmap_palette_method)
+  if(is.null(lim)) {
+    image(debut_vue:fin_vue,1:dim(data)[1],t(data), col = palette,xlab="",ylab="",xaxt="n",yaxt="n")
+  }
+  else {
+    image(debut_vue:fin_vue,1:dim(data)[1],t(data), col = palette,xlab="",ylab="",xaxt="n",yaxt="n",zlim=lim)
+  }
+  box()
+  if(is.null(ataxises)) {
+    axis(1)
+  }
+  else {
+    if(sum(!is.na(ataxises))!=0) {
+      axis(1,at=ataxises,labels=FALSE)
+      axis(3,at=ataxises,labels=FALSE)
+    }
+  }
+  if(sum(!is.na(ataxises))!=0) {
+    axis(2,at=1:dim(data)[1],labels=rownames(data),las=2,cex.axis=1)
+  }
+  text(fin_vue+(fin_vue-debut_vue)*0.015,(dim(data)[1]+1)/2,textOnLeftSide,xpd=NA,cex=1,srt=-90)
+}
+
+##Returns the title of visualisation
+plotTitle=function(debut_vue=1,fin_vue=length(listForward[[1]]),chr=NULL,style=NULL) {
+    plot(c(0,1),c(0,1),col="white",ylab="",xlab="",fg="white",col.axis="white",xaxs="i",yaxs="i")
+    text(0,0.5,paste(chr,":",debut_vue,"-",fin_vue,sep=""),cex=2.2,adj=0)
+}
+
+##The main function of visualisation
+plotVisu = function(file,typeVisu="classic",listForward,listReverse,which=1:length(listForward),stranded=TRUE,
+debut_vue=1,fin_vue=length(listForward[[1]]),chr=NULL,annot=NULL,style=NULL,tops=NULL,bottoms=NULL,marks=NULL,strandMarks=NULL,
+titres="",repeated=FALSE,name_flags="",decal=0,log=TRUE,classic_plus_color="navyblue",classic_minus_color="mediumvioletred",
+heatmap_max_color="#FFFFAA",heatmap_min_color="#000055",heatmap_palette_method="hsv",
+lines_samples_colors=c(1,3,4,2)[((0:length(listForward))%%4)+1],lines_samples_type_line=((0:length(listForward))%%4)+1,
+smoothLength=trunc((fin_vue-debut_vue)/1200),annotation_color_by_strand=FALSE,annotation_placed_by_strand=FALSE) {
+  n_element_vue = length(which)
+  i_data = which[1]
+  forward_matrice=NULL
+  reverse_matrice=NULL
+  all_stranded = data[[i_data]]$stranded
+  if(!repeated) {
+    for(i in 1:n_element_vue) {
+      i_data = which[i]
+      forward_matrice = rbind(forward_matrice,forward[[i_data]][debut_vue:fin_vue])
+      reverse_matrice = rbind(reverse_matrice,reverse[[i_data]][debut_vue:fin_vue])
+      all_stranded = all_stranded & data[[i_data]]$stranded
+    }
+  }
+  else{
+    for(i in 1:n_element_vue) {
+      i_data = which[i]
+      forward_matrice = rbind(forward_matrice,forward[[i_data]])
+      reverse_matrice = rbind(reverse_matrice,reverse[[i_data]])
+      all_stranded = all_stranded & data[[i_data]]$stranded
+    }
+  }
+  rownames(forward_matrice)=titres[which]
+  rownames(reverse_matrice)=titres[which]
+  for( i in 1:n_element_vue){
+    if(smoothLength>1) {
+      lo=smooth(forward_matrice[i,],L=smoothLength)
+      forward_matrice[i,]=lo
+      los=smooth(reverse_matrice[i,],L=smoothLength)
+      reverse_matrice[i,]=los
+    }
+  }
+  if(is.null(annot)) {
+    heights_ann=NULL
+  }
+  else {
+    heights_ann=sizePlotAnnotation(annot=annot,chr=chr,debut=debut_vue,fin=fin_vue,style=style,annotation_placed_by_strand=annotation_placed_by_strand)*50
+  }
+  ataxises=pretty((debut_vue:fin_vue)+decal,n=14)
+  if(log){
+    label_scal="log2 tag densities"
+  }else{
+    label_scal="tag densities"
+  }
+  if(debut_vue==1&decal==0) {
+    ataxises[1]=1
+  }##Classic visualisation
+  if(typeVisu=="classic") {
+    height_panels = c(40,rep(200,n_element_vue),heights_ann)
+    png(file,1200,sum(height_panels))
+    prev=par(no.readonly = TRUE)
+    n_panels=length(height_panels)
+    layout(matrix(1:n_panels,n_panels,1),heights=height_panels)
+    par(mar=c(0, 5, 0, 4)+0.1)
+    plotTitle(debut_vue=debut_vue,fin_vue=fin_vue,chr=chr)
+    par(mar=c(2.5, 5, 2.5, 4)+0.1)
+    for(element in 1:n_element_vue) {
+      i_data = which[element]
+      plotRNAseq(listForward[[i_data]],listReverse[[i_data]],debut_vue,fin_vue,chr,annot,style,top=tops[min(element,length(tops))],bottom=bottoms[min(element,length(bottoms))],y=label_scal,titre=titres[min(i_data,length(titres))],name_flags=name_flags,repeated=repeated,decal=decal,ataxises=ataxises,classic_plus_color=classic_plus_color,classic_minus_color=classic_minus_color,stranded=stranded)
+      if(!is.null(marks)) {
+	if(is.null(tops)) {
+	  top=max(forward[[i_data]])
+	}
+	else {
+	  top=tops[min(element,length(tops))]
+	}
+      if(is.null(bottoms)) {
+	bottom=max(reverse[[i_data]])
+      }
+      else {
+	bottom=bottoms[min(element,length(bottoms))]
+      }
+      if(is.null(strandMarks)) {
+	segments(marks,top,y1=-bottom,col=2)
+      }
+      else {
+	if(strandMarks==0|!data[[i_data]]$stranded) {
+	   segments(marks,0,y1=top,col=2)
+	}
+	else {
+	  if(strandMarks==1) {
+	    segments(marks,0,y1=-bottom,col=2)
+	  }
+	  else {
+	    segments(marks,top,y1=-bottom,col=2)
+	  }
+	}
+      }
+    }
+   }
+   par(mar=c(0,5,1,4)+0.1)
+   if(!is.null(annot)) {
+    plot_annotation(annot,chr,debut=debut_vue,fin=fin_vue,style=style,textSize=1.5,annotation_color_by_strand=annotation_color_by_strand,annotation_placed_by_strand=annotation_placed_by_strand)
+   }
+  }
+  else {##Heatmap visualisation
+  if(typeVisu=="heatmap") {
+    if(stranded) {
+      height_panels=c(40,50*n_element_vue+28,heights_ann,50*n_element_vue+30,100)
+    }
+    else {
+      height_panels=c(40,50*n_element_vue+28,heights_ann,100)
+    }
+    png(file,1200,sum(height_panels))
+    prev=par(no.readonly = TRUE)
+    n_panels=length(height_panels)
+    layout(matrix(1:n_panels,n_panels,1),heights=height_panels)
+    par(mar=c(0, 8, 0, 2)+0.1)
+    plotTitle(debut_vue=debut_vue,fin_vue=fin_vue,chr=chr)
+    limIntensity = c(0,max(forward_matrice,reverse_matrice,na.rm=TRUE))
+    ataxisesLabels=as.character(ataxises)
+    ataxisesLabels[((1:length(ataxises))%%2)==0]=""
+    rownames(forward_matrice)=titres[which]
+    rownames(reverse_matrice)=titres[which]
+    par(mar=c(1, 8, 1.8, 2)+0.1,cex=1.1)
+    if(stranded) {
+      tmp="Plus strand"
+    }
+    else {
+      tmp="Both strands"
+    }
+    myHeatMap(forward_matrice[(dim(forward_matrice)[1]):1,],debut_vue+decal,fin_vue+decal,ataxises,lim=limIntensity,heatmap_max_color=heatmap_max_color,heatmap_min_color=heatmap_min_color,heatmap_palette_method=heatmap_palette_method,textOnLeftSide=tmp)
+    axis(3,at=ataxises,labels=FALSE,cex.axis=1.4)
+    axis(3,at=ataxises,labels=ataxisesLabels,cex.axis=1.4,line=-0.4,lwd=0)
+    if(!is.null(marks)) {
+      if(is.null(strandMarks)|!all_stranded) {
+	segments(marks,-2,y1=dim(forward_matrice)[1]+2,col=2,lwd=2,xpd=TRUE)
+      }
+    else {
+      if(strandMarks!=1) {
+	segments(marks,-2,y1=dim(forward_matrice)[1]+2,col=2,lwd=2,xpd=TRUE)
+      }
+    }
+    }
+    if(name_flags!="") {
+      flags=try(get(name_flags),TRUE)
+      if(class(flags)!="try-error") {
+	f_=flags[(flags$Chr==chr)&(flags$Stop>=debut_vue)&(flags$Start<=fin_vue),]
+	N=dim(f_)[1]
+	points(f_$Start,rep(0,N),col=2,pch=19,cex=2)
+      }
+    }
+    if(!is.null(annot)) {
+      par(mar=c(0,8,0,2)+0.1)
+      plot_annotation(annot,chr,debut=debut_vue,fin=fin_vue,style=style,textSize=0.9,annotation_color_by_strand=annotation_color_by_strand,annotation_placed_by_strand=annotation_placed_by_strand)
+    }
+    if(stranded) {
+      par(mar=c(1, 8, 0, 2)+0.1)
+      myHeatMap(reverse_matrice[(dim(reverse_matrice)[1]):1,],debut_vue+decal,fin_vue+decal,ataxises,lim=limIntensity,heatmap_max_color=heatmap_max_color,heatmap_min_color=heatmap_min_color,heatmap_palette_method=heatmap_palette_method,textOnLeftSide="Minus strand")
+      if(!is.null(marks)) {
+	if(is.null(strandMarks)|!all_stranded) {
+	  segments(marks,-1,y1=dim(forward_matrice)[1]+1,col=2,lwd=2,xpd=TRUE)
+	}
+	else {
+	  if(strandMarks!=0) {
+	    segments(marks,-1,y1=dim(forward_matrice)[1]+1,col=2,lwd=2,xpd=TRUE)
+	  }
+	}
+      }
+      axis(1,at=ataxises,labels=ataxisesLabels,cex.axis=1.4)
+    }
+    palette=matrix(seq(limIntensity[1],limIntensity[2],length.out=2000),1,2000)
+    rownames(palette)="0"
+    par(mar=c(2,18,3, 12)+0.1)
+    myHeatMap(palette,1,2000,ataxises=NA,heatmap_max_color=heatmap_max_color,heatmap_min_color=heatmap_min_color,heatmap_palette_method=heatmap_palette_method)
+    labelAxisHeatmapLegend=pretty(c(limIntensity[1],limIntensity[2]),n=7)
+    atAxisHeatmapLegend=1+((labelAxisHeatmapLegend-limIntensity[1])/(limIntensity[2]-limIntensity[1]))*1999
+    axis(1,at=atAxisHeatmapLegend,labels=FALSE,cex.axis=1.4)
+    axis(1,at=atAxisHeatmapLegend,labels=labelAxisHeatmapLegend,cex.axis=1.4,line=-0.4,lwd=0)
+    text(1000,2,label_scal,xpd=NA,font=2,cex=1.4)
+
+   }##Lines visualisation
+    else if(typeVisu=="lines"){
+      legendSize = (floor(n_element_vue/2)+n_element_vue%%2)*40
+      height_panels=c(40,legendSize,400,heights_ann)
+      n_panels=length(height_panels)
+      png(file,1200,sum(height_panels))
+      prev=par(no.readonly = TRUE)
+      par(mar=c(0, 5, 0,4)+0.1,cex=1.1)
+      layout(matrix(c(1:n_panels),n_panels,4),heights=height_panels)
+      par(mar=c(0, 5, 0,4)+0.1)
+      plotTitle(debut_vue=debut_vue,fin_vue=fin_vue,chr=chr)
+      lines_legend(n_element_vue=n_element_vue,which,titres,lines_samples_colors=lines_samples_colors,lines_samples_type_line=lines_samples_type_line)
+      par(mar=c(3, 5, 0,4)+0.1)
+      plotlines(forward_matrice,reverse_matrice,which=which,debut_vue, fin_vue , chr, annot, style, tops, bottoms,marks,strandMarks,titres, repeated,name_flags,decal,ataxises=ataxises,n_element_vue=n_element_vue,y=label_scal,lines_samples_colors=lines_samples_colors,lines_samples_type_line=lines_samples_type_line,stranded=stranded)
+      par(mar=c(0,5,0,4)+0.1)
+      if(!is.null(annot)) {
+	plot_annotation(annot,chr,debut=debut_vue,fin=fin_vue,style=style,textSize=1.5,annotation_color_by_strand=annotation_color_by_strand,annotation_placed_by_strand=annotation_placed_by_strand)
+      }
+
+      }
+    }
+  invisible(dev.off())
+}
+
+##
+smooth = function(X,L=10) {
+  x_smooth=filter(X,rep(1,L)/L)
+  x_smooth[is.na(x_smooth)]=0
+  return(x_smooth)
+}
+
+##Returns lines visualisation
+plotlines = function(forward_matrice,reverse_matrice,which=1:length(forward), debut_vue = 1,fin_vue =length(forward), chr=NULL, annot=NULL, style=NULL, tops=NULL,bottoms=NULL,marks=NULL,strandMarks=NULL, titres="", repeated=FALSE,name_flags="",decal=0,ataxises=NULL,n_element_vue=length(which),y="",lines_samples_colors=c(1,3,4,2)[((0:length(forward))%%4)+1],lines_samples_type_line=((0:length(forward))%%4)+1,stranded=TRUE){
+  limIntensity = c(-max(bottoms),max(tops))
+plot(c(debut_vue,fin_vue)+decal,limIntensity,col="white",ylab=y,main="",xaxs="i",xaxt="n",yaxt="n",cex.lab=1.2,xlab="",cex.lab=1.8)
+   lty=lines_samples_type_line
+   col=lines_samples_colors
+  ataxises_y= pretty(limIntensity,n=8)
+  for(i in ataxises) abline(v=i,lty=2,col="#808080")
+  for(i in ataxises_y) abline(h=i,lty=2,col="#808080")
+  for( i in 1:n_element_vue){
+    lo=forward_matrice[i,]
+    lines((debut_vue:fin_vue)+decal,lo,type="l",lty =lty[i],col=col[i],lwd=2,xaxt="n")
+    if(stranded) {
+      los=-reverse_matrice[i,]
+      lines((debut_vue:fin_vue)+decal,los,type="l",lty= lty[i],col=col[i],lwd=2,xaxt="n")
+    }
+  }
+  ataxisesLabels=as.character(ataxises)
+  ataxisesLabels[((1:length(ataxises))%%2)==0]=""
+  ataxisesLabels_y=as.character(abs(ataxises_y))
+  ataxisesLabels_y[((1:length(ataxises_y))%%2)==0]=""
+  axis(1,at=ataxises,labels=FALSE,cex.axis=2.6)
+  axis(1,at=ataxises,labels=ataxisesLabels,cex.axis=2.6,line=0.6,lwd=0)
+  axis(2, at=ataxises_y,labels=FALSE,cex.axis=2)
+  axis(2, at=ataxises_y,labels=ataxisesLabels_y,cex.axis=2,line=-0.4,lwd=0)
+  if(stranded) {
+    abline(h=0,lwd=6)
+    abline(h=0,col="white",lwd=4)
+    text(fin_vue+(fin_vue-debut_vue)*0.01,tops/2,"+",xpd=NA,cex=3)
+    text(fin_vue+(fin_vue-debut_vue)*0.01,-bottoms/2,"-",xpd=NA,cex=3)
+    text(fin_vue+(fin_vue-debut_vue)*0.025,(tops-bottoms)/2,"Strand",xpd=NA,cex=2,srt=-90)
+  }
+}
+
+##Returns lines legend
+lines_legend=function(n_element_vue,which,titres,lines_samples_colors=c(1,3,4,2)[((0:n_element_vue)%%4)+1],lines_samples_type_line=((0:n_element_vue)%%4)+1){
+  lty=lines_samples_type_line
+  col=lines_samples_colors
+  n_y = floor(n_element_vue/2)+n_element_vue%%2
+  plot(c(0,4),c(0,-(n_y+1)),col="white", ylab="",xlab="",main="",fg="white",col.axis="white",yaxs="i")
+  i_style=0
+  for(i in 1:n_y) {
+    i_style=i_style+1
+    lines(c(1.5,2),-c(i,i),col=col[i_style],lty=lty[i_style],lwd=4)
+    text(1.48,-i,titres[which[i_style]],cex=2.6,adj=1)
+  }
+  if(n_element_vue>1) {
+    for(i in (n_y+1):n_element_vue) {
+      i_style=i_style+1
+      lines(c(3.5,4),-c(i,i)+n_y,col=col[i_style],lty=lty[i_style],lwd=4)
+      text(3.48,-i+n_y,titres[which[i_style]],cex=2.6,adj=1)
+    }
+  }
+}
+
+##Returns the shape of plain arrow for the annotation
+plain_arrow = function(left,right,y,thickness=1,pickSize=(right-left)*0.1,pickSide=c("right","left","both","none"),col="blue") {
+  if(pickSide[1]=="right") {
+    pick_point = min(right - pickSize,right)
+    polygon(c(left,pick_point,right,pick_point,left),c(y-thickness/2,y-thickness/2,y,y+thickness/2,y+thickness/2),col=col)
+  }
+  if(pickSide[1]=="left") {
+    pick_point = max(left + pickSize,left)
+    polygon(c(right,pick_point,left,pick_point,right),c(y-thickness/2,y-thickness/2,y,y+thickness/2,y+thickness/2),col=col)
+  }
+  if(pickSide[1]=="none") {
+    polygon(c(right,left,left,right),c(y-thickness/2,y-thickness/2,y+thickness/2,y+thickness/2),col=col)
+  }
+  if(pickSide[1]=="both") {
+    pick_point_1 = max(left + pickSize,left)
+    pick_point_2 = min(right - pickSize,right)
+    polygon(c(left,pick_point_1,pick_point_2,right,pick_point_2,pick_point_1),c(y,y-thickness/2,y-thickness/2,y,y+thickness/2,y+thickness/2),col=col)
+  }
+}
+
+##Returns the size of the  panel of the annotation
+sizePlotAnnotation=function(annot,chr,debut,fin,style=NULL,annotation_placed_by_strand=FALSE) {
+	left=c()
+	right=c()
+	annot_chr=annot[annot$Chr==chr,]
+	unique_ID=unique(annot_chr$ID)
+	for(j in 1:length(unique_ID)){
+	  left[j]=min(annot_chr$Start[annot_chr$ID==unique_ID[j]])
+	  right[j]=max(annot_chr$Stop[annot_chr$ID==unique_ID[j]])
+	}
+	y_plot=parking(left,right)
+      return(max(y_plot)-min(y_plot)+1)
+}
+
+##Function to organise the annotation shapes to display
+parking = function(left,right,plus=rep(TRUE,length(left)),biggestOnesInside=TRUE) {
+  y=rep(0,length(left))
+  if(sum(plus)>0) {
+    left_plus=left[plus]
+    right_plus=right[plus]
+    y_plus=rep(0,sum(plus))
+    lengths_plus=right_plus-left_plus
+    for(i in order(lengths_plus,decreasing=TRUE)) {
+      otherleft_plus=left_plus
+      otherleft_plus[i]=NA
+      otherright_plus=right_plus
+      otherright_plus[i]=NA
+      placed=FALSE
+      y_plus[i]=0.5
+      while(placed==FALSE) {
+	placed=sum((right_plus[i]>otherleft_plus[y_plus==y_plus[i]])&(left_plus[i]<otherright_plus[y_plus==y_plus[i]]),na.rm=TRUE)==0
+	if(placed==FALSE) {
+	  y_plus[i]=y_plus[i]+1
+	}
+      }
+    }
+    if(biggestOnesInside) {
+      y[plus]=y_plus
+    }
+    else {
+      y[plus]=sapply(y_plus,function(i) sort(unique(y_plus))[i==sort(unique(y_plus),decreasing=TRUE)])
+    }
+  }
+  if(sum(!plus)>0) {
+    left_minus=left[!plus]
+    right_minus=right[!plus]
+    y_minus=rep(0,sum(!plus))
+    lengths_minus=right_minus-left_minus
+    for(i in order(lengths_minus,decreasing=TRUE)) {
+      otherleft_minus=left_minus
+      otherleft_minus[i]=NA
+      otherright_minus=right_minus
+      otherright_minus[i]=NA
+      placed=FALSE
+      y_minus[i]=-0.5
+      while(placed==FALSE) {
+	placed=sum((right_minus[i]>otherleft_minus[y_minus==y_minus[i]])&(left_minus[i]<otherright_minus[y_minus==y_minus[i]]),na.rm=TRUE)==0
+	if(placed==FALSE) {
+	  y_minus[i]=y_minus[i]-1
+	}
+      }
+    }
+    if(biggestOnesInside) {
+      y[!plus]=y_minus
+    }
+    else {
+      y[!plus]=sapply(y_minus,function(i) sort(unique(y_minus))[i==sort(unique(y_minus),decreasing=TRUE)])
+    }
+  }
+  return(y)
+}
+
+##Function to
+plot_annotation = function(annot,chr,debut,fin,style=NULL,textSize=par("cex"),annotation_color_by_strand=FALSE,annotation_placed_by_strand=FALSE) {
+  left=c()
+  right=c()
+  annot_chr=annot[annot$Chr==chr,]
+  N = dim(annot_chr)[1]
+  strand=c("-","+")
+  if(N>0) {
+    unique_ID=unique(annot_chr$ID)
+    par_prev=par(no.readonly=TRUE)
+    for(j in 1:length(unique_ID)){
+      left[j]=min(annot_chr$Start[annot_chr$ID==unique_ID[j]])
+      right[j]=max(annot_chr$Stop[annot_chr$ID==unique_ID[j]])
+    }
+    if(annotation_placed_by_strand) {
+      y_plot=parking(left,right,annot_chr$Strand=="+",FALSE)
+    }
+    else {
+      y_plot=parking(left,right,biggestOnesInside=FALSE)
+    }
+    plot(c(debut,fin),c(min(y_plot)-0.5,max(y_plot)+0.5),col="white",ylab="",xlab="",fg="white",col.axis="white",xaxs="i",yaxs="i")
+    for(j in 1:length(unique_ID)){
+      annot_ID=annot_chr[annot_chr$ID==unique_ID[j],]
+      x_text<-Inf
+      for(i_annot_ID in 1:dim(annot_ID)[1]) {
+	if(annot_ID$Strand[i_annot_ID]=="-") {
+	  orientation="left"
+	}
+	else {
+	  orientation="right"
+	}
+	iDraw = FALSE
+	if(annot_ID$Strand[i_annot_ID]==".") {
+	  tmp="+"
+	}
+	else {
+	  tmp=annot_ID$Strand[i_annot_ID]
+	}
+	style_demande = style[style$Type==annot_ID$Type[i_annot_ID]&style$Strand==tmp,]
+	x_text<-min(x_text,annot_ID$Start[i_annot_ID])
+	if(style_demande$shape=="box") {
+	  plain_arrow(annot_ID$Start[i_annot_ID],annot_ID$Stop[i_annot_ID],y_plot[j],thickness=0.5, pickSide=orientation,col=style_demande$col,pickSize=(fin-debut)*0.02)
+	  iDraw=TRUE
+	}
+	if(style_demande$shape =="line") {
+	  segments(annot_ID$Start[i_annot_ID],y_plot[j],annot_ID$Stop[i_annot_ID],y_plot[j],col=style_demande$col)
+	  iDraw=TRUE
+	}
+	if(style_demande$shape == "rectangle") {
+	  plain_arrow(annot_ID$Start[i_annot_ID],annot_ID$Stop[i_annot_ID],y_plot[j],thickness=0.5, pickSide="none",col=style_demande$col)
+	  iDraw=TRUE
+	}
+	if(style_demande$shape == "arrow") {
+	  arrowHeads = pretty(debut:fin,n=50)
+	  x<-c(annot_ID$Start[i_annot_ID],arrowHeads[arrowHeads>annot_ID$Start[i_annot_ID]&arrowHeads<annot_ID$Stop[i_annot_ID]],annot_ID$Stop[i_annot_ID])
+	  if(annot_ID$Strand[i_annot_ID]=="-") {
+	    arrows(x[2:length(x)],y_plot[j],x[1:length(x)-1],col=style_demande$col,length=0.08)
+	  }
+	  else {
+	    arrows(x[1:length(x)-1],y_plot[j],x[2:length(x)],col=style_demande$col,length=0.08)
+	  }
+	  iDraw=TRUE
+	}
+	if(iDraw == FALSE) {
+	  segments(annot_ID$Start[i_annot_ID],y_plot[j],annot_ID$Stop[i_annot_ID],y_plot[j],col="black")
+	}
+      }
+      text(max(x_text,debut),y_plot[j]+0.4,adj=0,labels=unique_ID[j],cex=textSize)
+    }
+    if(annotation_placed_by_strand) {
+      abline(h=0.125,lty=2,col="black")
+    }
+  }
+}
+
+# produit la visualisation de la légende
+plotLegend= function(style,width_plot=400,length_arrows=100,type="horizontal") {
+  if(type=="horizontal") {
+    style_unique = unique(style$col)
+    n_style = length(style_unique)
+    plot(c(0,width_plot),c(1,2*n_style+1),col="white",ylab="",xlab="",xaxs="i",fg="white",col.axis="white",main="Legend")
+    for(i in 1:n_style) {
+      arrows(ceiling((width_plot-length_arrows)/2),2*(n_style-i+1),ceiling((width_plot+length_arrows)/2),2*(n_style-i+1),col=style_unique[i],length=0.10)
+      text(ceiling(width_plot/2),2*(n_style-i+1)+0.5,paste(style$Type[style$col==style_unique[i]],collapse=", "),col=style_unique[i])
+    }
+  }
+  else {
+    style_unique = unique(style$col)
+    n_style = length(style_unique)
+    length_arrows = width_plot/(n_style*1.25)
+    inbetween = length_arrows/4
+    plot(c(0,width_plot),c(0,1),col="white",ylab="",xlab="",yaxs="i",xaxt="n",yaxt="n",bty="n",xaxs="i",col.axis="white",main="")
+    text(width_plot/2,0.8,"Legend",cex=1.4,font=2)
+    for(i in 1:n_style) {
+      fleche_start = (length_arrows+inbetween)*(i-1) + inbetween/2
+      fleche_stop = (length_arrows+inbetween)*(i-1) + length_arrows + inbetween/2
+#      arrows(fleche_start,0.2,fleche_stop,0.2,col=style_unique[i],length=0.10)
+      rect(fleche_start,0.2,fleche_stop,0.6,border=style_unique[i])
+      text((fleche_start+fleche_stop)/2,0.4,paste(style$Type[style$col==style_unique[i]],collapse=", "),col=style_unique[i],cex=1.2)
+    }
+  }
+}
+
+html2rgb = function(codeHTML) {
+  chiffres=c(0:9,"A","B","C","D","E","F")
+  codeHTMLsplit=unlist(strsplit(codeHTML,split=""))
+  red = (((which(chiffres==codeHTMLsplit[2])-1)*16)+(which(chiffres==codeHTMLsplit[3])-1))/255
+  green = (((which(chiffres==codeHTMLsplit[4])-1)*16)+(which(chiffres==codeHTMLsplit[5])-1))/255
+  blue = (((which(chiffres==codeHTMLsplit[6])-1)*16)+(which(chiffres==codeHTMLsplit[7])-1))/255
+  return(c(red,green,blue))
+}
+
+rgb2hsv = function(vectRGB) {
+  Cmax=max(vectRGB)
+  Cmin=min(vectRGB)
+  delta=Cmax-Cmin
+  if(delta==0) {
+    H=0
+  }
+  else {
+    if(Cmax==vectRGB[1]) H=(1/6)*(((vectRGB[2]-vectRGB[3])/delta)%%6)
+    if(Cmax==vectRGB[2]) H=(1/6)*(((vectRGB[3]-vectRGB[1])/delta)+2)
+    if(Cmax==vectRGB[3]) H=(1/6)*(((vectRGB[1]-vectRGB[2])/delta)+4)
+  }
+  if(delta==0) {
+    S=0
+  }
+  else {
+    S=delta/Cmax
+  }
+  V=Cmax
+
+  return(c(H,S,V))
+}
+
+### =========================================================================
+### Ving's Main
+### -------------------------------------------------------------------------
+###
+{
+if(version$major != "3" & version$minor != "0.2") {
+  warning("You need 3.0.2 R Version to run the program")
+}
+tmp = unlist(strsplit(commandArgs()[4],split="="))[2]
+path_split = unlist(strsplit(tmp,split="/"))
+path = paste(path_split[1:(length(path_split)-1)],collapse="/")
+arg<-commandArgs(TRUE)
+if(length(arg)==0) {
+	cat("Ving for visualisation of RNA seq data:
+Usage:
+    ./visualiseur.R [options] <input> [<input2> <input3> ...]
+  options:
+    -o/--output                     <string>     [default:./ouput.png]
+    -v/--typeVisu                   <string>     [default: classic   ] (classic,lines,heatmap)
+    -t/--description-data           <string>     [default: filename  ]
+    -n/--normalization-coefficient  <string>     [default: none      ]
+    -i/--inverseStrand
+    -u/--unstranded
+    -c/--chromosome-name            <string>     [default: first one ]
+    -d/--start                      <integer>    [default: 1         ]
+    -f/--end                        <integer>    [default: chr length]
+    -a/--annotation                 <filename>
+    -r/--typeTranscript             <string>     [default: all       ]
+    -C/--annotation-colors          <string>     [default: fabulous  ]
+    -s/--annotation-shapes          <string>     [default: box       ] (box,line,rectangle,arrow)
+    -l/--scale-log
+    -y/--symetric-scale
+    --classic-plus-color            <string>     [default: purple    ]
+    --classic-minus-color           <string>     [default: pink      ]
+    --heatmap-max-color             <string>     [default: green     ]
+    --heatmap-min-color             <string>     [default: yellow    ]
+    --heatmap-palette-method        <string>     [default: hsv       ] (hsv,rgb)
+    --lines-samples-colors          <string>     [default: black     ]
+    --lines-samples-type-line       <string>     [default: plain     ]
+    --annotation-color-by-strand    <logical>    [default: FALSE     ]
+    --annotation-placed-by-strand   <logical>    [default: FALSE     ]
+    \n\n")
+	q("no")
+}else {
+  tmp=suppressPackageStartupMessages(require("Rsamtools"))
+  if(!tmp) {
+    stop("Package Rsamtools required !!")
+  }
+  tmp=suppressPackageStartupMessages(require("GenomicRanges"))
+  if(!tmp) {
+    stop("Package GenomicRanges required !!")
+  }
+  optArgs=getopt(
+      rbind(
+
+	  c('output','o', 1, 'character',"output.png"),
+	  c('typeVisu', 'v', 1, 'character', "classic"),
+	  c('description-data','t',1,'character',NA),
+	  c('chromosome-name', 'c', 1, 'character', NA),
+	  c('start', 'd', 1, 'numeric',1),
+	  c('end', 'f', 1, 'numeric',NA),
+	  c('annotation','a',1,'character',NA),
+	  c('typeTranscript','r',1,'character',NA),
+	  c('annotation-colors','C',1,'character',NA),
+	  c('annotation-shapes','s','1','character',NA),
+	  c('normalization-coefficient','n','1','character',NA),
+	  c('classic-plus-color',1,1,'character',"navyblue"),
+	  c('classic-minus-color',2,1,'character',"mediumvioletred"),
+	  c('heatmap-max-color',3,1,'character',"#000055"),
+	  c('heatmap-min-color',4,1,'character',"#FFFFAA"),
+	  c('heatmap-palette-method',7,1,'character',"hsv"),
+	  c('lines-samples-colors',5,1,'character',NA),
+	  c('lines-samples-type-line',6,1,'character',NA),
+	  c('scale-log', 'l',0,'logical',FALSE),
+	  c('inverseStrand','i',0,'logical', FALSE),
+	  c('unstranded','u',0,'logical', FALSE),
+	  c('symetric-scale','y',0,'logical', FALSE),
+	  c('annotation-color-by-strand',8,0,'logical',FALSE),
+	  c('annotation-placed-by-strand',9,0,'logical',FALSE)
+      )
+  )
+}
+###################
+## ARGUMENTS
+#################################################################################
+files=optArgs$ARGUMENT
+##Case file doesn't exist
+  for( i in 1:length(files)) {
+    if(!file.exists(files[i])) {
+      stop(paste(files[i],"do not exists!","\n"))
+    }
+  }
+imagefile=optArgs$output
+typeVisu=optArgs$typeVisu
+description_data=optArgs$`description-data`
+chrName = optArgs$`chromosome-name`[1]
+start=optArgs$start[1]
+end=optArgs$end[1]
+annotation=optArgs$annotation
+typeTranscript=optArgs$typeTranscript
+colors=optArgs$`annotation-colors`
+shape_data=optArgs$`annotation-shapes`
+weight=optArgs$`normalization-coefficient`
+classic_plus_color=optArgs$`classic-plus-color`
+classic_minus_color=optArgs$`classic-minus-color`
+heatmap_max_color=optArgs$`heatmap-max-color`
+heatmap_min_color=optArgs$`heatmap-min-color`
+heatmap_palette_method=optArgs$`heatmap-palette-method`
+lines_samples_colors=optArgs$`lines-samples-colors`
+lines_samples_type_line=optArgs$`lines-samples-type-line`
+log=optArgs$`scale-log`
+inverseStrand=optArgs$inverseStrand
+unstranded=optArgs$unstranded
+symetric= optArgs$`symetric-scale`
+annotation_color_by_strand=optArgs$`annotation-color-by-strand`
+annotation_placed_by_strand=optArgs$`annotation-placed-by-strand`
+###################
+## MAIN
+###################################################################################
+
+
+  genome_info=scanBamHeader(files[1])[[1]]$targets
+  noms_chromosomes = names(genome_info)
+  longueur_chromosomes = as.integer(genome_info)
+  nombre_chromosomes = length(noms_chromosomes)
+
+##Case no chromosome specified
+  if(sum(is.na(chrName))) {
+    chrName = noms_chromosomes[1]
+    cat("No chromosome specified, processing chromosome :",chrName,"\n")
+  }
+
+
+##Case false chromosome name
+  if(!(chrName %in% noms_chromosomes)) {
+    stop(paste("\"",chrName,"\" is not a proper chromosome name"))
+  }
+
+  if(is.na(end)) {
+    end = longueur_chromosomes[chrName==noms_chromosomes]
+  }
+
+  if(start > end) {
+    stop("The start is bigger than the end!")
+  }
+
+##Case asked coordinates outside the chromosome
+  if(start<0|end>longueur_chromosomes[chrName==noms_chromosomes]) {
+    stop("You are outside the chromosome")
+  }
+
+  if(sum(is.na(weight))>0) {
+    normalized_data=NULL
+    isnormalized=FALSE
+  }
+   else {
+    isnormalized=TRUE
+    normalized_data=unlist(strsplit(weight,split=","))
+    if(length(files)!=length(normalized_data)) {
+      stop("Different number of files and  weights ")
+    }
+    else {
+      normalized_data=as.numeric(normalized_data)
+      }
+  }
+
+  if(inverseStrand) {
+    libraryType="inverse"
+  }
+  else {
+    libraryType="standard"
+  }
+##Read the bam file and extract the infos
+  doit=function(i,libraryType) {
+
+   try(readBam_(files[i], libraryType=libraryType, chrName_=chrName, from_=start, to_=end,normalized_=normalized_data[i]))
+  }
+  data=lapply(1:length(files),doit,libraryType=libraryType)
+  ctrl=unlist(lapply(data,class))=="try-error"
+
+  if(sum(ctrl)>0) {
+    for(i in which(ctrl)) {
+      cat("Problem with file :",files[i],"\n")
+    }
+    stop("At least a file has problem")
+  }
+
+##Read the GFF file and extract the infos
+  if(sum(is.na(annotation))==0) {
+    gff=try(readGff(annotation[1],from=start,to=end,chr=chrName),TRUE)
+    ctrl=class(gff)=="try-error"
+    if(sum(ctrl)>0) {
+	stop(paste("Problem with  Gff file :",annotation,"\n"))
+    }
+    if(length(annotation)>1) {
+      for(i in 2:length(annotation)) {
+	gff1=try(readGff(annotation[i],from=start,to=end,chr=chrName),TRUE)
+	ctrl=class(gff1)=="try-error"
+	if(sum(ctrl)>0) {
+	  stop(paste("Problem with  Gff file :",gff1,"\n"))
+	}
+	gff=rbind(gff,gff1)
+      }
+    }
+
+  }
+  else {
+    gff=NA
+  }
+
+
+  if(sum(is.na(description_data))>0) {
+    description_data = files
+  }
+##Case of different number of files and description
+  if(length(description_data)!=length(files)) {
+    stop("Different number of files and description")
+  }
+
+  if(sum(is.na(typeTranscript)>0)) {
+    if(sum(is.na(annotation))==0) {
+      typeTranscritSplit=unique(gff$Type)
+    }else {
+      typeTranscritSplit=NA
+    }
+  }else {
+    typeTranscritSplit=unlist(strsplit(typeTranscript,split=","))
+  }
+
+  if((sum(is.na(gff))>0)|(sum(is.na(typeTranscritSplit))>0)) {
+    annot=NULL
+  }else {
+    annot = gff[gff$Type %in% typeTranscritSplit,]
+  }
+
+##check the  colors
+
+  if(sum(is.na(colors)>0)) {
+    nTypeTranscrit = length(typeTranscritSplit)
+    if(annotation_color_by_strand) {
+      colorsSplit=c(classic_plus_color,classic_minus_color)
+    }
+    else {
+      colorsSplit=rainbow(nTypeTranscrit+1)[1:nTypeTranscrit]
+    }
+  }
+  else {
+    colorsSplit=unlist(strsplit(colors,split=","))
+  }
+  for(i in 1:length(colorsSplit)) {
+    tmp=unlist(strsplit(colorsSplit[i],split=""))
+    if(length(tmp)==6|length(tmp)==8) {
+      if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+	colorsSplit[i] = paste("#",colorsSplit[i],sep="")
+      }
+    }
+  }
+
+  tmp=unlist(strsplit(classic_plus_color,split=""))
+  if(length(tmp)==6|length(tmp)==8) {
+    if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+      classic_plus_color = paste("#",classic_plus_color,sep="")
+    }
+  }
+  tmp=unlist(strsplit(classic_minus_color,split=""))
+  if(length(tmp)==6|length(tmp)==8) {
+    if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+      classic_minus_color = paste("#",classic_minus_color,sep="")
+    }
+  }
+  tmp=unlist(strsplit(heatmap_max_color,split=""))
+  if(length(tmp)==6|length(tmp)==8) {
+    if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+    heatmap_max_color = paste("#",heatmap_max_color,sep="")
+    }
+  }
+  tmp=unlist(strsplit(heatmap_min_color,split=""))
+  if(length(tmp)==6|length(tmp)==8) {
+    if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+    heatmap_min_color = paste("#",heatmap_min_color,sep="")
+    }
+  }
+  if(sum(is.na(lines_samples_colors))>0) {
+    lines_samples_colors_split=c(1,3,4,2)[((0:(length(files)-1))%%4)+1]
+  }
+  else{
+    lines_samples_colors_split=unlist(strsplit(lines_samples_colors,split=","))
+    for(i in 1:length(lines_samples_colors_split)) {
+      tmp=unlist(strsplit(lines_samples_colors_split[i],split=""))
+      if(length(tmp)==6|length(tmp)==8) {
+	if(sum(tmp %in% c(0:9,"A","B","C","D","E","F"))==length(tmp)) {
+	  lines_samples_colors_split[i] = paste("#",lines_samples_colors_split[i],sep="")
+	}
+      }
+    }
+  }
+  colorToCheck=c(colorsSplit,classic_plus_color,classic_minus_color,heatmap_max_color,heatmap_min_color,lines_samples_colors_split)
+  ctrl=sapply(colorToCheck,is.acceptable.color)
+  if(sum(!ctrl)>0) {
+    for(i in which(!ctrl)) {
+      cat("\"",colorToCheck[i],"\" is not a proper color name.\n",sep="")
+    }
+    stop("At least one color has problem")
+  }
+  if(annotation_color_by_strand){
+    if(length(colorsSplit)>2) {
+      stop("You have to specify two and only two colors!")
+    }
+  }else {
+    if(length(typeTranscritSplit)!=length(colorsSplit)) {
+      stop("Please specify the same number of transcript types and colors")
+    }
+  }
+##check the line types
+  if(sum(is.na(lines_samples_type_line))>0) {
+    lines_samples_type_line_split=((0:(length(files)-1))%/%4)+1
+  }
+  else {
+    lines_samples_type_line_split=unlist(strsplit(lines_samples_type_line,split=","))
+  }
+  if(typeVisu=="lines") {
+    ctrl=sapply(lines_samples_type_line_split,function(x) {
+	tmp=suppressWarnings(as.numeric(x))
+	if(!is.na(tmp)) {
+	  return((tmp==floor(tmp))&tmp>=1&tmp<=5)
+	}
+	else {
+	  return(FALSE)
+	}
+      }
+    )
+    if(sum(!ctrl)>0) {
+      for(i in which(!ctrl)) {
+	cat("\"",lines_samples_type_line_split[i],"\" is not a proper line style.\n",sep="")
+      }
+      stop("At least one line style has problem")
+    }
+    lines_samples_type_line_split=as.integer(lines_samples_type_line_split)
+  }
+##check the shapes
+type_shape=rep(1,length(typeTranscritSplit))
+
+if(sum(is.na(shape_data)>0)) {
+    for( i in 1:length(typeTranscritSplit)) {
+	    type_shape[i]="box"
+    }
+}else{
+    shape=unlist(strsplit(shape_data,split=","))
+    shape=as.array(shape)
+    if(length(typeTranscritSplit)!=length(shape)) {
+      stop("Please specify the same number of transcript types and shapes")
+      }else {
+      for( i in 1:length(typeTranscritSplit)){
+	type_shape[i]= shape[[i]]
+      }
+  }
+}
+
+
+##Style for the annotation
+label=rep(1,length(typeTranscritSplit))
+style=data.frame(Type=c(typeTranscritSplit,typeTranscritSplit),Strand=c(rep("+",length(typeTranscritSplit)),rep("-",length(typeTranscritSplit))),col=NA,shape=NA,label,stringsAsFactors = FALSE)
+for(i in 1:length(typeTranscritSplit)) {
+  style$shape[style$Type==typeTranscritSplit[i]]=type_shape[i]
+}
+if(annotation_color_by_strand) {
+  style$col[style$Strand=="+"]=colorsSplit[1]
+  style$col[style$Strand=="-"]=colorsSplit[2]
+}
+else {
+  for(i in 1:length(typeTranscritSplit)) {
+    style$col[style$Type==typeTranscritSplit[i]]=colorsSplit[i]
+  }
+}
+
+##Main for visualisation
+databychr=extractSignal(data,chrName,from=start,to=end,normalized_=isnormalized)
+
+{
+  reverse=list()
+  forward=list()
+  if(log) {
+  for(i in 1:length(databychr$F)){
+    forward_ = numeric()
+    tmp=log2(1+databychr$F[[i]])
+    forward_[1:length(tmp)]=tmp
+    forward[[i]]=forward_
+    reverse_ = numeric()
+    tmp=log2(1+databychr$R[[i]])
+    reverse_[1:length(tmp)]=tmp
+    reverse[[i]]=reverse_
+  }
+  }
+  else {
+      forward=databychr$F
+      reverse=databychr$R
+  }
+}
+
+smoothLength=max(trunc((end-start)/1200),1)
+if(unstranded) {
+  for(i in 1:length(databychr$F)){
+    lo=smooth(forward[[i]]+reverse[[i]],L=smoothLength)
+    forward[[i]]=lo
+    los=rep(0,length(lo))
+    reverse[[i]]=los
+  }
+}
+else {
+  for(i in 1:length(databychr$F)){
+    lo=smooth(forward[[i]],L=smoothLength)
+    forward[[i]]=lo
+    los=smooth(reverse[[i]],L=smoothLength)
+    reverse[[i]]=los
+  }
+}
+
+group_maximum = rep(1,length(databychr$F))
+max_forward = numeric(length(databychr$F))
+max_reverse = numeric(length(databychr$F))
+for(i_data in 1:length(databychr$F)) {
+  max_forward[i_data]=max(forward[[i_data]],na.rm=TRUE)
+  max_reverse[i_data]=max(reverse[[i_data]],na.rm=TRUE)
+}
+for(i_max in unique(group_maximum)) {
+  max_forward[group_maximum==i_max]=max(max_forward[group_maximum==i_max],na.rm=TRUE)
+  max_reverse[group_maximum==i_max]=max(max_reverse[group_maximum==i_max],na.rm=TRUE)
+}
+
+if(symetric){
+  for(i_data in 1:length(databychr$F)) {
+    max_forward[i_data]=max(max_forward[i_data],max_reverse[i_data],na.rm=TRUE)
+  }
+  max_reverse=max_forward
+}
+
+plotVisu(imagefile,typeVisu=typeVisu,listForward=forward,listReverse=reverse,
+debut_vue=start,fin_vue=end,chr=chrName,annot=annot,repeated=TRUE,
+titres=description_data,name_flags="",style=style,log=log,stranded=!unstranded,
+tops=max_forward,bottoms=max_reverse,
+classic_plus_color=classic_plus_color,classic_minus_color=classic_minus_color,
+heatmap_max_color=heatmap_max_color,heatmap_min_color=heatmap_min_color,heatmap_palette_method=heatmap_palette_method,
+lines_samples_colors=lines_samples_colors_split,lines_samples_type_line=lines_samples_type_line_split,
+smoothLength=1,annotation_color_by_strand=annotation_color_by_strand,annotation_placed_by_strand=annotation_placed_by_strand)
+
+}