improved star wrapper

dge_workflow/star_align.sh

@@ -5,7 +5,7 @@
 
 usage='
 Use star to align fastq files against a genome
-Usage: star_align.sh <igenome> <fastq_files>
+Usage: star_align.sh <igenome> <fastq_files>...
 
 Options:
 -c        Cache results
@@ -21,31 +21,37 @@ eval "$(echo  "$usage" | ~/bin/docopts/docopts -h - : "$@")"
 # v0.7 style
 #eval $(echo  "$usage" | /home/brandl/bin/docopts_v0.7/docopts "hallo")
 
-#eval "$(echo "$usage" | ~/bin/docopts/docopts -h -  : /projects/bioinfo/igenomes/Mus_musculus/Ensembl/GRCm38 /projects/bioinfo/holger/projects/florio_11b_2nd_batch/lanereps_pooled/arhgap11b_1.fastq.gz)"
+#eval "$(echo "$usage" | ~/bin/docopts/docopts -h -  : /projects/bioinfo/igenomes/Mus_musculus/Ensembl/GRCm38 /projects/bioinfo/holger/projects/florio_11b_2nd_batch/lanereps_pooled/arhgap11b_1.fastq.gz fastq2)"
 
-#for fastqFile in ${fastq_files[@]} ; do
-#    echo processing $fastqFile
-#done
+for fastqFile in ${fastq_files[@]} ; do
+    echo processing $fastqFile
+done
 #
 #echo $igenome
 
 
+## extract all configuration parameters
 fastqFiles=${fastq_files[@]}
-#echo $fastqFiles
-
 export star_index="${igenome}/Sequence/StarIndex"
 export gtfFile="$igenome/Annotation/Genes/genes.gtf"
+
+#echo $fastqFiles
 #head $gtfFile
 
 if [ ! -f $gtfFile ]; then
     >&2 echo "gtf '$gtfFile' does not exis"; exit 1;
 fi
 
+
+if [ -z "$(which STAR)" ]; then
+    >&2 echo "no STAR binary in PATH"; exit 1;
+fi
+
 ## basic usage tutorial
 #http://www.homolog.us/blogs/blog/2012/11/02/star-really-kick-ass-rna-seq-aligner/
 
 ## build index if not present
-if ! -d "${star_index}" ]; then
+if ! -f "${star_index}/SA" ]; then
 
 mailme "${project}: creating STAR index for $igenome"
 mkdir ${star_index}
@@ -60,17 +66,8 @@ chmod -R -w ${star_index}
 fi
 
 
-
-if [ ! -f $gtfFile ]; then
-    >&2 echo "gtf '$gtfFile' does not exis"; exit 1;
-fi
-
-if [ -z "$(which STAR)" ]; then
-    >&2 echo "no tomcat binary in PATH"; exit 1;
-fi
-
-
 echo "running STAR using igenome '${igenome}' for the following files"
+
 ll $fastqFiles
 
 #fastqFiles=$(ls $baseDir/treps_pooled/*fastq.gz)
@@ -87,19 +84,20 @@ for fastqFile in $fastqFiles ; do
 
     ## note --outSAMstrandField intronMotif is required for cuffdiff compatiblity (xs flag)
     mysub "${project}__star__${fastqBaseName}" "
-#    tophat -p6  -G $gtfFile -g1 -o $outputdir $bowtie_gindex $fastqFile
-    STAR --genomeDir $star_index --readFilesIn $fastqFile --runThreadN 6 --readFilesCommand zcat --outFileNamePrefix $fastqBaseName --outSAMtype BAM SortedByCoordinate --outSAMstrandField intronMotif --sjdbGTFfile $gtfFile
-    " -n 5 -R span[hosts=1] -q long | joblist .tophatjobs
+    STAR --genomeDir $star_index --readFilesIn $fastqFile --runThreadN 6 --readFilesCommand zcat --outFileNamePrefix ${fastqBaseName}. --outSAMtype BAM SortedByCoordinate --outSAMstrandField intronMotif --sjdbGTFfile $gtfFile
+    mv ${fastqBaseName}.Aligned.sortedByCoord.out.bam ${fastqBaseName}.bam
+    samtools index ${fastqBaseName}.bam
+    " -n 5 -R span[hosts=1] -q medium | joblist .tophatjobs
 done
 
 wait4jobs .tophatjobs
 
 
 
-dge_bam_correlate .
+dge_bam_correlate . &
 
 ## create tophat mapping report
 ## todo adjust report to STAR
 #spin.R ${NGS_TOOLS}/dge_workflow/tophat_qc.R .
 
-mailme "$project: STAR done in $(pwd)"
\ No newline at end of file
+mailme "${project}: STAR done in $(pwd)"
\ No newline at end of file