#!/bin/bash

#module load bioinfo-tools
#module load samtools/1.5
#module load bwa/0.7.17

bamlist=$1
bedfile=$2
refgenome=$3
bedfilter_dir=$4
snplist=$5
TMPDIR=$6
for inbam in `cat $bamlist`; do
# get just the bam reads that cover the vars of interest
samtools view -L $bedfile -h -b -u -q 30 $inbam > ${TMPDIR}/bedfilter.bam
# index this new bam
samtools index ${TMPDIR}/bedfilter.bam

# add md tag to the bam, pass to the modify script which presumably uses this info, which outputs a modified bam
samtools calmd ${TMPDIR}/bedfilter.bam $refgenome | ./modify_read_alternative.py $snplist | samtools view -b -u - > ${TMPDIR}/modified.bam

# align the reads from this new modified bam to the reference
bwa aln -l 16500 -n 0.01 -o 2 -t 15 $refgenome -b ${TMPDIR}/modified.bam \
| bwa samse $refgenome - ${TMPDIR}/modified.bam  \
| samtools view -F 4 -h - \
| samtools sort -o ${TMPDIR}/modified.remap.bam -

fbase=$(basename $inbam .bam)

# index the new remapped modified bam
samtools index ${TMPDIR}/modified.remap.bam

# add md tag to new remapped modified bam
samtools calmd -b ${TMPDIR}/modified.remap.bam $refgenome > ${TMPDIR}/modified.remap.MD.bam

# index MD version
samtools index ${TMPDIR}/modified.remap.MD.bam

# feed both the original bam (with only vars of interest) and the remapped modified bam into script that must compare them and outputs the final bam file.
samtools calmd ${TMPDIR}/bedfilter.bam $refgenome | samtools view -L $bedfile -h -q 30 - | ./filter_sam_startpos_dict.py ${TMPDIR}/modified.remap.MD.bam | samtools view -bh - > ${TMPDIR}/remapmodread.filter.bam

#index final bam
samtools index ${TMPDIR}/remapmodread.filter.bam

mv ${TMPDIR}/modified.remap.bam ${bedfilter_dir}/$fbase.modified.remap.bam
mv ${TMPDIR}/modified.remap.bam.bai ${bedfilter_dir}/$fbase.modified.remap.bam.bai

mv ${TMPDIR}/remapmodread.filter.bam ${bedfilter_dir}/$fbase.remapmodread.filter.bam
mv ${TMPDIR}/remapmodread.filter.bam.bai ${bedfilter_dir}/$fbase.remapmodread.filter.bam.bai

done

# mv ${TMPDIR}/modified.remap.posfilter.bam ${bedfilter_dir}/$fbase.modified.remap.posfilter.bam
# mv ${TMPDIR}/modified.remap.posfilter.bam.bai ${bedfilter_dir}/$fbase.modified.remap.posfilter.bam.bai

