commit d973538c276e33c29b0e533a726e26607f25feb7
Author: Florent Chuffart <florent.chuffart@ens-lyon.fr>
Date:   Thu Mar 5 15:31:47 2015 +0100

    Using dplyr to filter reads
    Fixing bug concernong number of MNase replicate != 3
    Adding NEUTRAL_TRANSLATE_CUR parameter in configuration file

diff --git a/build/nucleominer_2.3.47.tar.gz b/build/nucleominer_2.3.47.tar.gz
new file mode 100644
index 0000000..3c2474a
Binary files /dev/null and b/build/nucleominer_2.3.47.tar.gz differ
diff --git a/doc/Chuffart_NM2_sweave/Chuffart_NM2_sweave.Rnw b/doc/Chuffart_NM2_sweave/Chuffart_NM2_sweave.Rnw
index 12749a5..06e6617 100644
--- a/doc/Chuffart_NM2_sweave/Chuffart_NM2_sweave.Rnw
+++ b/doc/Chuffart_NM2_sweave/Chuffart_NM2_sweave.Rnw
@@ -41,7 +41,7 @@ R.Version()$version.string
   config = fromJSON(paste(readLines(json_conf_file), collapse=""))
 
   # Read sample file
-  all_samples = get_content(config$CSV_SAMPLE_FILE, "cvs", sep=";", head=TRUE, stringsAsFactors=FALSE)  
+  all_samples = get_content(config$CSV_SAMPLE_FILE, "csv", sep=";", head=TRUE, stringsAsFactors=FALSE)  
   
   # Remove samples that seems to be eronous for some manipulation efficency reasons e.g. PCR
   all_samples = all_samples[!(all_samples$id %in% c(33,45,48,55)), ]
@@ -705,7 +705,7 @@ print(nb_iso)
 print(nb_iso/sum(check_iso_snep$snep_index))
 
 
-wpunr_mnase = read.table( paste(config$RESULTS_DIR, "/full/" ,combi[1],"_",combi[2],"_wp_mnase.tab",sep=""), header=TRUE)
+wpunr_mnase = read.table( paste(config$RESULTS_DIR, "/full/" ,combi[1],"_",combi[2],"_wpunr_mnase.tab",sep=""), header=TRUE)
 head(wpunr_mnase)
 
 sf = read.table(paste(config$RESULTS_DIR,"/full/size_factors.tab",sep=""), header=TRUE)        
diff --git a/doc/Chuffart_NM2_workdir/src/current/configurator.py b/doc/Chuffart_NM2_workdir/src/current/configurator.py
index fd9a07d..302714e 100755
--- a/doc/Chuffart_NM2_workdir/src/current/configurator.py
+++ b/doc/Chuffart_NM2_workdir/src/current/configurator.py
@@ -3,7 +3,7 @@ import csv
 import json
 
 CSV_SAMPLE_FILE = None
-"""Path to cvs file that contains sample information."""
+"""Path to csv file that contains sample information."""
 if __name__ == "__main__":
   CSV_SAMPLE_FILE = "data/samples.csv"
 
@@ -98,6 +98,12 @@ if __name__ == "__main__":
     ]
   }
 
+NEUTRAL_TRANSLATE_CUR = None
+"""Desactivate the genome translation facilities"""
+if __name__ == "__main__":
+  NEUTRAL_TRANSLATE_CUR = False
+
+
 READ_LENGTH = None
 """Length of Illumina reads."""
 if __name__ == "__main__":
diff --git a/doc/Chuffart_NM2_workdir/src/current/configurator.pyc b/doc/Chuffart_NM2_workdir/src/current/configurator.pyc
index ce8e035..692bb7d 100644
Binary files a/doc/Chuffart_NM2_workdir/src/current/configurator.pyc and b/doc/Chuffart_NM2_workdir/src/current/configurator.pyc differ
diff --git a/doc/Chuffart_NM2_workdir/src/current/extract_maps.R b/doc/Chuffart_NM2_workdir/src/current/extract_maps.R
index a92280f..210a1df 100644
--- a/doc/Chuffart_NM2_workdir/src/current/extract_maps.R
+++ b/doc/Chuffart_NM2_workdir/src/current/extract_maps.R
@@ -46,7 +46,7 @@ extract_maps = function(task) {
 }
 
 print("Running bot engine...")
-run_engine(tasks, extract_maps, debug = substr(Sys.info()[["nodename"]],1,7) == "stainer", rm_starter = TRUE, log_dir = "log", nb_proc=NULL)
+run_engine(tasks, extract_maps, debug = substr(Sys.info()[["nodename"]],1,7) == "stainer", rm_starter = FALSE, log_dir = "log", nb_proc=NULL)
 
 
 
@@ -107,18 +107,27 @@ build_wp_maps = function() {
 			for (strain in strains) {
 				print(paste("Collecting mpas for", strain, "..."))
 				partial_strain_maps = aggregated_intra_strain_nucs[[strain]][[1]]
-        tmp_nuc_map = flat_aggregated_intra_strain_nucs(partial_strain_maps, cur_index)
+        print(strain)
+        nb_tracks = sum(all_samples$strain == strain & all_samples$marker == "Mnase_Seq")
+        
+        tmp_nuc_map = flat_aggregated_intra_strain_nucs(partial_strain_maps, cur_index, nb_tracks)
 
         wp_indexes = which(tmp_nuc_map$wp == 1)
         tmp_wp_llr = apply(t(wp_indexes), 2, function(wp_index){
-          tmp_wp = partial_strain_maps[[wp_index]]  
-    			res = llr_score_nvecs(list(tmp_wp$nucs[[1]]$original_reads, tmp_wp$nucs[[2]]$original_reads, tmp_wp$nucs[[3]]$original_reads))
+          tmp_wp = partial_strain_maps[[wp_index]]
+          l = lapply(tmp_wp$nucs, function(tmp_wp_nucs){
+            tmp_wp_nucs$original_reads
+          })  
+          res = llr_score_nvecs(l)
           return(res)
         })
 
         tmp_dyad_shift = apply(t(wp_indexes), 2, function(wp_index){
           tmp_wp = partial_strain_maps[[wp_index]]  
-          ds = max(c(tmp_wp$nucs[[1]]$center, tmp_wp$nucs[[2]]$center, tmp_wp$nucs[[3]]$center)) - min(c(tmp_wp$nucs[[1]]$center, tmp_wp$nucs[[2]]$center, tmp_wp$nucs[[3]]$center))
+          c = sapply(tmp_wp$nucs, function(tmp_wp_nucs){
+            tmp_wp_nucs$center
+          })  
+          ds = max(c) - min(c)
           return(ds)
         })
 
diff --git a/doc/Chuffart_NM2_workdir/src/current/headers.R b/doc/Chuffart_NM2_workdir/src/current/headers.R
index 0da2a33..fd3a3ce 100644
--- a/doc/Chuffart_NM2_workdir/src/current/headers.R
+++ b/doc/Chuffart_NM2_workdir/src/current/headers.R
@@ -5,7 +5,7 @@ library(cachecache)
 json_conf_file = "src/current/nucleominer_config.json"
 config = fromJSON(paste(readLines(json_conf_file), collapse=""))
 
-all_samples = get_content(config$CSV_SAMPLE_FILE, "cvs", sep=";", head=TRUE, stringsAsFactors=FALSE)
+all_samples = get_content(config$CSV_SAMPLE_FILE, "csv", sep=";", head=TRUE, stringsAsFactors=FALSE)
 all_samples = all_samples[!(all_samples$id %in% c(33,45,48,55)), ]
 
 markers = unique(all_samples$marker)
diff --git a/doc/sphinx_doc/Makefile b/doc/sphinx_doc/Makefile
index 17a8f92..eaf3403 100644
--- a/doc/sphinx_doc/Makefile
+++ b/doc/sphinx_doc/Makefile
@@ -112,7 +112,7 @@ latexpdf:
 
 text:
 	$(SPHINXBUILD) -b text $(ALLSPHINXOPTS) $(BUILDDIR)/text
-	cat build/text/index.txt build/text/readme.txt build/text/tuto.txt ../../README 
+	cat build/text/index.txt build/text/readme.txt build/text/tuto.txt > ../../README 
 	@echo
 	@echo "Build finished. The text files are in $(BUILDDIR)/text."
 
diff --git a/doc/sphinx_doc/build/latex/NucleoMiner2.pdf b/doc/sphinx_doc/build/latex/NucleoMiner2.pdf
index a418027..53ae8b3 100644
Binary files a/doc/sphinx_doc/build/latex/NucleoMiner2.pdf and b/doc/sphinx_doc/build/latex/NucleoMiner2.pdf differ
diff --git a/doc/sphinx_doc/conf.py b/doc/sphinx_doc/conf.py
index 55057fe..c647fb5 100644
--- a/doc/sphinx_doc/conf.py
+++ b/doc/sphinx_doc/conf.py
@@ -50,9 +50,9 @@ copyright = u'2012, Florent CHUFFART, Jean-Baptiste VEYRIERAS, Gael YVERT'
 # built documents.
 #
 # The short X.Y version.
-version = '2.3.46'
+version = '2.3.47'
 # The full version, including alpha/beta/rc tags.
-release = '2.3.46'
+release = '2.3.47'
 
 # The language for content autogenerated by Sphinx. Refer to documentation
 # for a list of supported languages.
diff --git a/doc/sphinx_doc/rref.rst b/doc/sphinx_doc/rref.rst
index f5e2110..1323002 100644
--- a/doc/sphinx_doc/rref.rst
+++ b/doc/sphinx_doc/rref.rst
@@ -690,7 +690,8 @@ Usage
 ::
 
     filter_tf_inputs(inputs, chr, x_min, x_max, nuc_width = 160, 
-        only_f = FALSE, only_r = FALSE, filter_for_coverage = FALSE)
+        only_f = FALSE, only_r = FALSE, filter_for_coverage = FALSE, 
+        USE_DPLYR = TRUE)
 
 Arguments
 ~~~~~~~~~
@@ -727,6 +728,10 @@ Filter only R reads.
 
 Does it filter for plot coverage?
 
+``USE_DPLYR``
+
+Use dplyr lib to filter reads.
+
 Value
 ~~~~~
 
@@ -814,7 +819,8 @@ Usage
 
 ::
 
-    flat_aggregated_intra_strain_nucs(partial_strain_maps, cur_index)
+    flat_aggregated_intra_strain_nucs(partial_strain_maps, cur_index, 
+        nb_tracks = 3)
 
 Arguments
 ~~~~~~~~~
@@ -827,6 +833,10 @@ the output of aggregate\_intra\_strain\_nucs function
 
 the index of the roi involved
 
+``nb_tracks``
+
+the number of replicates
+
 Value
 ~~~~~
 
@@ -1251,13 +1261,13 @@ Details
 +---------------+---------------------------------------------------+
 | Author:       | Florent Chuffart                                  |
 +---------------+---------------------------------------------------+
-| Version:      | 2.3.46                                            |
+| Version:      | 2.3.47                                            |
 +---------------+---------------------------------------------------+
 | License:      | CeCILL                                            |
 +---------------+---------------------------------------------------+
 | Title:        | nm                                                |
 +---------------+---------------------------------------------------+
-| Depends:      | seqinr, plotrix, DESeq, cachecache                |
+| Depends:      | seqinr, plotrix, DESeq, cachecache, dplyr         |
 +---------------+---------------------------------------------------+
 
 Author(s)
diff --git a/src/DESCRIPTION b/src/DESCRIPTION
index 2ffc4fa..8a15c36 100644
--- a/src/DESCRIPTION
+++ b/src/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: nucleominer
 Maintainer: Florent Chuffart <florent.chuffart@ens-lyon.fr>
 Author: Florent Chuffart
-Version: 2.3.46
+Version: 2.3.47
 License: CeCILL 
 Title: nm
-Depends: seqinr, plotrix, DESeq, cachecache
+Depends: seqinr, plotrix, DESeq, cachecache, dplyr
 Description: It provides a set of useful functions allowing to perform quantitative analysis of nucleosomal epigenome.
diff --git a/src/R/nucleominer.R b/src/R/nucleominer.R
index a06e449..75e530b 100644
--- a/src/R/nucleominer.R
+++ b/src/R/nucleominer.R
@@ -153,20 +153,42 @@ x_max, ##<< Coordinate of the last bp observed.
 nuc_width = 160, ##<< Nucleosome width.
 only_f = FALSE, ##<< Filter only F reads.
 only_r = FALSE, ##<< Filter only R reads.
-filter_for_coverage = FALSE ##<< Does it filter for plot coverage?
+filter_for_coverage = FALSE, ##<< Does it filter for plot coverage?
+USE_DPLYR = TRUE ##<< Use dplyr lib to filter reads.
 ) {
-	if (only_f) {
-		inputs = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,3] == "F" & inputs[,2] <= x_max + nuc_width,]
-	} else if (only_r) {
-		inputs = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,3] == "R" & inputs[,2] <= x_max + nuc_width,]
-	} else {
-		inputs = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,2] <= x_max + nuc_width,]
-	}
+  n = names(inputs)
+
+  if (!USE_DPLYR) {  
+    if (only_f) {
+      inputs_out = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,3] == "F" & inputs[,2] <= x_max + nuc_width,]
+    } else if (only_r) {
+      inputs_out = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,3] == "R" & inputs[,2] <= x_max + nuc_width,]
+    } else {
+      inputs_out = inputs[inputs[,1]==chr & inputs[,2] >= x_min - nuc_width & inputs[,2] <= x_max + nuc_width,]
+    }
+  } else {
+    names(inputs) = c("chr", "pos", "str", "lev")
+    if (only_f) {
+      inputs_out = filter(inputs, chr == chr,  pos >= x_min - nuc_width, str == "F", pos <= x_max + nuc_width)
+    } else if (only_r) {
+      inputs_out = filter(inputs, chr == chr, pos >= x_min - nuc_width, str == "R" & pos <= x_max + nuc_width)
+    } else {
+      inputs_out = filter(inputs, chr == chr, pos >= x_min - nuc_width, pos <= x_max + nuc_width)
+    }
+      # if (!filter_for_coverage) {
+      #   inputs$corrected_inputs_coords = inputs[,2] + nuc_width/2 * sign_from_strand(inputs[,3])
+      #   inputs = filter(inputs, chr == chr, corrected_inputs_coords >= x_min, corrected_inputs_coords <= x_max)
+      #   inputs$corrected_inputs_coords = NULL
+      # }  
+  }
+
   if (!filter_for_coverage) {
-    corrected_inputs_coords = inputs[,2] + nuc_width/2 * sign_from_strand(inputs[,3])
-    inputs = inputs[inputs[,1]==chr & corrected_inputs_coords >= x_min & corrected_inputs_coords <= x_max,]
+    corrected_inputs_coords = inputs_out[,2] + nuc_width/2 * sign_from_strand(inputs_out[,3])
+    inputs_out = inputs_out[inputs_out[,1]==chr & corrected_inputs_coords >= x_min & corrected_inputs_coords <= x_max,]
   }
-	return(inputs)
+
+  names(inputs_out) = n 
+	return(inputs_out)
 ### Returns filtred inputs.
 }
 
@@ -205,7 +227,7 @@ coord_max=20000000 ##<< A too big value to be a coord for a nucleosome lower bou
 	  }
 	  return(TRUE)
 	}
-	store_cluster = function(clusters, new_cluster, nb_nucs_in_cluster,nuc_from_track,nb_tracks, min_nuc_center, max_nuc_center) {
+	store_cluster = function(clusters, new_cluster, nb_nucs_in_cluster, nuc_from_track, nb_tracks, min_nuc_center, max_nuc_center) {
 		if ( nb_nucs_in_cluster==nb_tracks & sum(nuc_from_track)==nb_tracks) {
 			new_cluster$wp = TRUE
 			center = (new_cluster$lower_bound + new_cluster$upper_bound) / 2
@@ -248,11 +270,12 @@ coord_max=20000000 ##<< A too big value to be a coord for a nucleosome lower bou
 		}
 		i = i+1
   }
+  nb_tracks = length(tf_outs)
 	# print(track_readers)
   new_cluster = NULL
   nb_nucs_in_cluster = 0
   nuc_from_track = c()
-  for (i in 1:length(tf_outs)){
+  for (i in 1:nb_tracks){
     nuc_from_track[i] = FALSE
   }
   # Start clustering
@@ -297,12 +320,12 @@ coord_max=20000000 ##<< A too big value to be a coord for a nucleosome lower bou
     } else {
 			if (!is.null(new_cluster)) {
         # store old cluster
-	      clusters = store_cluster(clusters, new_cluster, nb_nucs_in_cluster,nuc_from_track,length(tf_outs),min_nuc_center, max_nuc_center)
+	      clusters = store_cluster(clusters, new_cluster, nb_nucs_in_cluster, nuc_from_track, nb_tracks, min_nuc_center, max_nuc_center)
 			}
       # Reinit current cluster composition stuff
       nb_nucs_in_cluster = 0
       nuc_from_track = c()
-      for (i in 1:length(tf_outs)){
+      for (i in 1:nb_tracks){
         nuc_from_track[i] = FALSE
       }
       # create new cluster
@@ -331,7 +354,7 @@ coord_max=20000000 ##<< A too big value to be a coord for a nucleosome lower bou
   # store last cluster
   if (!is.null(new_cluster)) {
     # store old cluster
-    clusters = store_cluster(clusters, new_cluster, nb_nucs_in_cluster,nuc_from_track,length(tf_outs),min_nuc_center, max_nuc_center)
+    clusters = store_cluster(clusters, new_cluster, nb_nucs_in_cluster, nuc_from_track, nb_tracks, min_nuc_center, max_nuc_center)
   }
 	return(list(clusters, llr_scores))
 ### Returns a list of clusterized nucleosomes, and all computed llr scores.
@@ -364,7 +387,8 @@ coord_max=20000000 ##<< A too big value to be a coord for a nucleosome lower bou
 flat_aggregated_intra_strain_nucs = function(# to flat aggregate_intra_strain_nucs function output
 ### This function builds a dataframe of all clusters obtain from aggregate_intra_strain_nucs function.
 partial_strain_maps, ##<< the output of aggregate_intra_strain_nucs function
-cur_index ##<< the index of the roi involved
+cur_index, ##<< the index of the roi involved
+nb_tracks=3 ##<< the number of replicates
 ) {
 	if  (length(partial_strain_maps) == 0 ){
 		print(paste("Empty partial_strain_maps for roi", cur_index, "ands current strain." ))
@@ -386,11 +410,13 @@ cur_index ##<< the index of the roi involved
 			tmp_nuc_as_list[["nb_reads"]] = length(all_original_reads)
 			tmp_nuc_as_list[["nb_nucs"]] = length(tmp_nuc$nucs)
 			if (tmp_nuc$wp) {
-				tmp_nuc_as_list[["llr_1"]] = signif(tmp_nuc$nucs[[2]]$llr_score,5)
-				tmp_nuc_as_list[["llr_2"]] = signif(tmp_nuc$nucs[[3]]$llr_score,5)
+        for (i in 1:(nb_tracks-1)) {
+  				tmp_nuc_as_list[[paste("llr", i, sep="_")]] = signif(tmp_nuc$nucs[[i + 1]]$llr_score,5)          
+        }
 			} else {
-				tmp_nuc_as_list[["llr_1"]] = NA
-				tmp_nuc_as_list[["llr_2"]] = NA
+        for (i in 1:(nb_tracks-1)) {
+  				tmp_nuc_as_list[[paste("llr", i, sep="_")]] = NA
+        }
 			}
       return(tmp_nuc_as_list)
     })
@@ -1141,12 +1167,17 @@ strain2, ##<< The strain in wich you want the genome region of interest.
 config=NULL, ##<< GLOBAL config variable
 big_cur=NULL ##<< A largest region than roi use to filter c2c if it is needed.
 ) {
-	strain1 = roi$strain_ref
-	if (strain1 == strain2) {
-		roi$length = roi$end - roi$begin + sign(roi$end - roi$begin) * 1
+	strain1 = roi$strain_ref  
+  # Do something or nothing?
+  if (is.null(config$NEUTRAL_TRANSLATE_CUR)) {
+    config$NEUTRAL_TRANSLATE_CUR = FALSE
+  }
+	if (strain1 == strain2 | config$NEUTRAL_TRANSLATE_CUR) {
+    roi$strain_ref = strain2
+    roi$length = roi$end - roi$begin + sign(roi$end - roi$begin)
 		return(roi)
 	}
-
+  
 	# Extract c2c file
 	if (!is.null(big_cur)) {
   	# Dealing with big_cur
@@ -1964,7 +1995,7 @@ config=NULL ##<< GLOBAL config variable
       if (plot_wp_nucs | plot_fuzzy_nucs | plot_common_nucs ) {
     		replicates_wp_nucs[[replicate_rank]] = wp_nucs
         strain = samples[[1]]$strain
-        wp_maps[[strain]] = flat_aggregated_intra_strain_nucs(wp_nucs, "foo")
+        wp_maps[[strain]] = flat_aggregated_intra_strain_nucs(wp_nucs, "foo", nb_tracks)
         fuzzy_maps[[strain]] = get_intra_strain_fuzzy(wp_maps[[strain]], as.list(samples[[1]]$roi), samples[[1]]$strain, config=config)
 
         if (plot_fuzzy_nucs) {
