#!/bin/bash

#SBATCH --nodes=1
#SBATCH --ntasks-per-node=32
#SBATCH --mem=32GB
#SBATCH --time=24:00:00
#SBATCH --output=slurm-Bowtie2-%j.log


module load Bowtie2/2.5.4-GCC-13.2.0
module load SAMtools/1.19.2-GCC-13.2.0
module load Biopython/1.84-foss-2024a

# Global parameters needed for this script
# setting global parameters
cpu=16
SOFTWAREDIR=/scratch/public/genomes/AnneTools




# =============================================================================================================================================================
# ==                                                          Functions                                                                                      ==
# =============================================================================================================================================================
function my_bowtie_SE {  # Map the reads (from the fastq file) to the reference genome ( .fna ) this will generate a SAM file (Sequence Alignment Map)
	name=$1
	forward=$2
	printf "\n\n========= Bowtie2 SE: "$name" =================================\n" >> $bamdir/02.$name.bowtie2.log
	(bowtie2 --local -U $fastqdir/$forward -p $cpu -S $bamdir/$name.sam -x $genome) 2>>$bamdir/02.$name.bowtie2.log 
	SAM_2_BAM $name
}

function my_bowtie_SE_2_fastq {  # Map the reads (from the fastq file) to the reference genome ( .fna ) this will generate a SAM file (Sequence Alignment Map)
	name=$1
	printf "\n\n========= Bowtie2 SE: "$name" =================================\n" >> $bamdir/02.$name.bowtie2.log
	(bowtie2 --local -U $fastqdir/$2,$fastqdir/$3 -p $cpu -S $bamdir/$name.sam -x $genome) 2>>$bamdir/02.$name.bowtie2.log 
	SAM_2_BAM $name
}


function SAM_2_BAM {  # convert SAM to sorted BAM (binairy SAM) and remove unmapped reads with option -F 4
	name=$1
	printf "\n========== Convert "$name" ======================================\n"	
	samtools view -S -b -F 4 -o $bamdir/$name.bam $bamdir/$name.sam
	samtools sort -@ $cpu $bamdir/$name.bam > $bamdir/$name.sorted.bam
	samtools index $bamdir/$name.sorted.bam
	rm $bamdir/$name.sam  # save space by removing the SAM file 
	rm $bamdir/$name.bam  # save space by removing the non sorted BAM file                                                                          # 3) save space by removing the SAM file 
}

function hits_summery {  # counts the multiple hits
	name=$1
	samtools view -F 4 $bamdir/$name | awk '{print $1}' | sort | uniq -c | \
	awk '{count[$1]++; sum += $1; n++} END {for (i in count) print i, count[i]; print "\nAverage hits per mapped read:", sum/n}' | sort -n > $bamdir/$name.hits_summary.txt
}


fastqdir=/scratch/hb-molgen/NGSdata/20260828_036
workdir=/scratch/p127804/PROJECTS/Douwe/Sofia_2026_036
bamdir=$workdir/bam

mkdir -p $bamdir


genome=/scratch/p127804/PROJECTS/Douwe/genomes/Bifidobacterium_adolescentis_LMG
#bowtie2-build $genome.fna $genome


my_bowtie_SE SS036_01 SS036_01_S37_L001_R1_001.fastq.gz
my_bowtie_SE SS036_02 SS036_02_S38_L001_R1_001.fastq.gz
my_bowtie_SE SS036_03 SS036_03_S39_L001_R1_001.fastq.gz
my_bowtie_SE SS036_04 SS036_04_S40_L001_R1_001.fastq.gz
my_bowtie_SE SS036_05 SS036_05_S41_L001_R1_001.fastq.gz
my_bowtie_SE SS036_06 SS036_06_S42_L001_R1_001.fastq.gz
my_bowtie_SE SS036_07 SS036_07_S43_L001_R1_001.fastq.gz
my_bowtie_SE SS036_08 SS036_08_S44_L001_R1_001.fastq.gz




############################################################################################
# # # Make a report summary

cat $bamdir/02.*.bowtie2.log > $bamdir/02.bowtie2_all.log
python3 $SOFTWAREDIR/Bowtie2LogTable.py -query $bamdir/02.bowtie2_all.log > $bamdir/02.bowtie2_all.log.table