# Getting BAM files.
home.dir <- "../.."
bam.files <- list.files(file.path(home.dir, "bam"), pattern="bam$", full=TRUE)

# Assembling an annotation file.
anno.dir <- "~/lustre/annotation"
system(sprintf("cat %s/ERCC92.gtf %s/mm10.gtf > temp.gtf", anno.dir, anno.dir))

# Running featureCounts.
require(Rsubread)
out <- featureCounts(bam.files, annot.ext="temp.gtf", isGTFAnnotationFile=TRUE,	minMQS=10, nthreads=4)

# Saving counts to file, with gene names.
colnames(out$counts) <- sub("\\.bam$", "", basename(bam.files))
final <- data.frame(GeneID=rownames(out$counts), Length=out$annotation$Length, out$counts)
write.table(file="genic_counts.tsv", final, col.names=TRUE, row.names=FALSE, quote=FALSE, sep="\t")

# Augmenting the stats.
my.stats <- read.table(file.path(home.dir, "all_qual.tsv"), header=TRUE)
m <- match(my.stats$Sample, colnames(out$counts))
my.stats$Genic <- as.integer(out$stat[1,-1][m])
write.table(file="my_qual.tsv", my.stats, col.names=TRUE, row.names=FALSE, quote=FALSE, sep="\t")

# Saving the session information.
unlink("temp.gtf")
sessionInfo()