#
# Prerequisites:
# links to fastq files
# feba code linked to or downloaded

echo Compute GC contents of genes
for org in Keio MR1 Phaeo psRCH2 SB2B; do
echo Compute GC content of genes for $org
feba/bin/RegionGC.pl g/$org/genome.fna g/$org/genes.tab > g/$org/genes.GC
done

echo Processing all the reads in fastq/

echo
echo Define mutant pools using TnSeq reads
feba/bin/RunTnSeqLocal.pl Keio Keio_ML9 ezTn5_kan1 fastq/TnSeq_Keio_ML9
feba/bin/RunTnSeqLocal.pl MR1 MR1_ML3 pKMW7 fastq/TnSeq_MR1_ML3.fastq.gz
feba/bin/RunTnSeqLocal.pl Phaeo Phaeo_ML1 pKMW7 fastq/TnSeq_Phaeo_ML1.fastq.gz
feba/bin/RunTnSeqLocal.pl psRCH2 psRCH2_ML7 pHIMAR_kan fastq/TnSeq_psRCH2_ML7.fastq.gz
feba/bin/RunTnSeqLocal.pl SB2B SB2B_ML5 pHIMAR_kan fastq/TnSeq_SB2B_ML5.fastq.gz
(cd g/Keio; ln -s Keio_ML9.pool pool)
(cd g/MR1; ln -s MR1_ML3.pool pool)
(cd g/Phaeo; ln -s Phaeo_ML1.pool pool)
(cd g/psRCH2; ln -s psRCH2_ML7.pool pool)
(cd g/SB2B; ln -s SB2B_ML5.pool pool)
feba/bin/DesignRandomPool.pl  -pool g/psRCH2/pool.n2lax -minN 2 -minFrac 0 -minRatio 0 -genes g/psRCH2/genes.tab fastq/TnSeq_psRCH2_ML7..parts/*.mapped
mkdir --parents other_results
cp g/psRCH2/pool.n2lax other_results/psRCH2.pool.n2lax
feba/bin/DesignRandomPool.pl -pool g/SB2B/pool.n2lax -minN 2 -minFrac 0 -minRatio 0 -genes g/SB2B/genes.tab fastq/TnSeq_SB2B_ML5..parts/*.mapped
cp g/SB2B/pool.n2lax other_results/SB2B.pool.n2lax

echo
echo Analyzing reads for psRCH2 minipool from plate 17
for i in S1 S2 S3 S4 S5; do
(feba/bin/MapTnSeq.pl -model feba/primers/model_pHIMAR_kan -first fastq/TnSeq_minipool/*_${i}_L001_R1_001.fastq.gz -genome g/psRCH2/genome.fna > other_results/psRCH2_minipool_$i.tab) >& other_results/psRCH2_minipool_$i.log;
feba/bin/DesignRandomPool.pl -pool other_results/psRCH2_minipool_$i.pool -genes g/psRCH2/genes.tab -minN 10 other_results/psRCH2_minipool_$i.tab;
done

echo
echo Analyzing quantitative TnSeq data
rm -Rf fastq/TnSeq_test_libraries/map.cmd* >& /dev/null
(for i in `find fastq/TnSeq_test_libraries -name '*.fastq.gz' -print | grep sb2b | sed -e s/.fastq.gz//`; do
  echo "(feba/bin/MapTnSeq.pl -minQuality 0 -model feba/primers/model_pHIMAR_kan -first $i.fastq.gz -genome g/SB2B/genome.fna > $i.tab) >& $i.log";
  done;
  for i in `find fastq/TnSeq_test_libraries -name '*.fastq.gz' -print | grep psRCH2 | sed -e s/.fastq.gz//`; do
  echo "(feba/bin/MapTnSeq.pl -minQuality 0 -model feba/primers/model_pHIMAR_kan -first $i.fastq.gz -genome g/psRCH2/genome.fna > $i.tab) >& $i.log";
  done
) > fastq/TnSeq_test_libraries/map.cmds
feba/bin/submitter.pl fastq/TnSeq_test_libraries/map.cmds
for dir in `(cd fastq/TnSeq_test_libraries; ls -d Sample*)`; do
cut -f 3,4,5 fastq/TnSeq_test_libraries/$dir/*.tab | sed -e 's/\t/:/g' | scripts/tabulate.pl -skip 0 | perl -F\\t -ane 'next if $F[0] eq "Var1"; @G = split /:/, $F[0], -1; print join("\t",@G,$F[1])."\n"' > other_results/$dir.tot;
echo Wrote other_results/$dir.tot;
done
for dir in `(cd fastq/TnSeq_test_libraries; ls -d Sample* | grep sb2b)`; do
scripts/AddGeneInfo.pl other_results/$dir.tot g/SB2B/genes.tab > other_results/$dir.tot2;
done
for dir in `(cd fastq/TnSeq_test_libraries; ls -d Sample* | grep RCH2)`; do
scripts/AddGeneInfo.pl other_results/$dir.tot g/psRCH2/genes.tab > other_results/$dir.tot2;
done


echo
echo Processing BarSeq reads

feba/bin/RunBarSeqLocal.pl g/Keio Keio_ML9_set1 fastq/Keio_ML9_set1
feba/bin/RunBarSeqLocal.pl g/Keio Keio_ML9_set2 fastq/Keio_ML9_set2

feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set1 fastq/MR1_ML3_set1
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set2 fastq/MR1_ML3_set2
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set3 fastq/MR1_ML3_set3
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set3_repB fastq/MR1_ML3_set3_repB
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set4 fastq/MR1_ML3_set4
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set5 fastq/MR1_ML3_set5
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_ML3_set6 fastq/MR1_ML3_set6
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/MR1 MR1_SB2B_combined_set fastq/MR1_SB2B_combined_set

feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/Phaeo Phaeo_Dino_combined_set fastq/Phaeo_Dino_combined_set
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/Phaeo Phaeo_ML1_set2 fastq/Phaeo_ML1_set2
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/Phaeo Phaeo_ML1_set3 fastq/Phaeo_ML1_set3
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/Phaeo Phaeo_ML1_set4 fastq/Phaeo_ML1_set4
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/Phaeo Phaeo_ML1_set5 fastq/Phaeo_ML1_set5

feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set1 fastq/psRCH2_ML7_set1.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set1_seqUC fastq/psRCH2_ML7_set1_seqUC
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set2 fastq/psRCH2_ML7_set2.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set3 fastq/psRCH2_ML7_set3.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set4 fastq/psRCH2_ML7_set4.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set5 fastq/psRCH2_ML7_set5
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set10 fastq/psRCH2_ML7_set10
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set11 fastq/psRCH2_ML7_set11
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/psRCH2 psRCH2_ML7_set12 fastq/psRCH2_ML7_set12

feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/SB2B SB2B_ML5_set1 fastq/SB2B_ML5_set1.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/SB2B SB2B_ML5_set2 fastq/SB2B_ML5_set2.fastq.gz
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/SB2B SB2B_ML5_set4 fastq/SB2B_ML5_set4
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/SB2B SB2B_ML5_set5 fastq/SB2B_ML5_set5
feba/bin/RunBarSeqLocal.pl -indexes feba/primers/BarSeqPrimersH48 g/SB2B MR1_SB2B_combined_set fastq/MR1_SB2B_combined_set

echo
echo Analyzing shotgun genomic DNA sequencing reads
# see samtools, bowtie-build, and then use blockCoverage.pl (10 kb)
bowtie-build -q g/SB2B/genome.fna g/SB2B/index
samtools faidx g/SB2B/genome.fna
export GDNA=`(cd  fastq/gDNA_SB2B; ls SB2B-*.gz) | sed -e 's/.fastq.gz//'`

for i in $GDNA; do
cat <<END >  fastq/gDNA_SB2B/$i.cmds
#
zcat fastq/gDNA_SB2B/$i.fastq.gz | bowtie --quiet -S -m 1 -t g/SB2B/index - fastq/gDNA_SB2B/$i.map
samtools view -bS -o fastq/gDNA_SB2B/$i.bam fastq/gDNA_SB2B/$i.map
samtools sort fastq/gDNA_SB2B/$i.bam fastq/gDNA_SB2B/$i.sorted
samtools index fastq/gDNA_SB2B/$i.sorted.bam
samtools depth fastq/gDNA_SB2B/$i.sorted.bam > fastq/gDNA_SB2B/$i.depth
scripts/blockCoverage.pl g/SB2B/genome.fna 10 fastq/gDNA_SB2B/$i.depth > other_results/$i.blocks10
END
chmod +x fastq/gDNA_SB2B/$i.cmds;
done
(for i in $GDNA; do echo fastq/gDNA_SB2B/$i.cmds; done) > fastq/gDNA_SB2B/all.cmds
feba/bin/submitter.pl fastq/gDNA_SB2B/all.cmds

echo
echo Computing fitness values
mkdir --parents html
feba/bin/BarSeqR.pl -org Keio
feba/bin/BarSeqR.pl -org MR1
feba/bin/BarSeqR.pl -org Phaeo
feba/bin/BarSeqR.pl -org psRCH2
feba/bin/BarSeqR.pl -org SB2B
scripts/save_SB2B_halves.R

echo End of setup_rbarseq