#cat << EOT > polyAT.fa
#>polyA/1
#AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA
#>polyT/1
#TTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTT
#EOT
#trimmomatic SE -phred33 raw/3_TTAGGC_L003_R1_001.fastq filtered/3_TTAGGC_L003_R1_001.fastq HEADCROP:15 ILLUMINACLIP:/usr/local/share/trimmomatic/adapters/TruSeq3-SE.fa:2:30:10 ILLUMINACLIP:polyAT.fa:2:30:10 TRAILING:30 MINLEN:75
bowtie2-build dmel-all-chromosome-r6.04.fasta dmel.bowtie.index
ln -s dmel-all-chromosome-r6.04.fasta dmel.bowtie.index.fa
tophat2 -p 8 --GTF dmel-all-r6.04.gtf --transcriptome-index dmel.transcriptome.tophat dmel.bowtie.index

for f in `ls 0_raw/*.fastq.gz`; do
  file_id="${f%%.*}"
  cutadapt --cut=15 -a file:/usr/local/share/trimmomatic/adapters/TruSeq3-SE.fa -a "A{20}" -a "T{20}" --minimum-length=50 --quality-cutoff=30 --trim-n --max-n=0 -o 1_filtered/${f} 0_raw/${f}
  tophat2 -p 8 --no-novel-juncs --GTF dmel-all-r6.04.gtf --transcriptome-index dmel.transcriptome.tophat -o 2_mapped/${file_id}.tophat dmel.bowtie.index 1_filtered/${f}
  samtools view -h 2_mapped/${file_id}.tophat/accepted_hits.bam > 2_mapped/${file_id}.tophat/accepted_hits.sam
  cufflinks -p 8 --GTF dmel-all-r6.04.gtf -o 3_assembled/${file_id}.cufflinks 2_mapped/${file_id}.tophat/accepted_hits.sam
  echo 3_assembled/${file_id}.cufflinks/transcripts.gtf >> 4_merged/assembly.gtfs
done

cuffmerge -p 8 -g dmel-all-r6.04.gtf -s dmel-all-chromosome-r6.04.fasta -o 4_merged/assembly.cufflinks 4_merged/assembly.gtfs