#Running WASP (mappability filtering) with STAR alignmnet. Sample: sample_id | xthreads   
#Script run for each sample, base edits to commands can be made to this file and then recall "/home/asiimwe/projects/run_env/alpha_star_wasp_comparison/Run_Dirs_STAR_WASP.py"
#Replacements per run are made to: sample_id, xthreads, "runThreadN x", and the "sample snp directory - SSD"
#Directory STAR_WASP_Runs needs to be created prior to running the python file


WASP=/home/asiimwe/WASP
PYTHON=/home/asiimwe/miniconda3/bin/python3.9
export PATH=/usr/bin/samtools/:$PATH
STAR=/usr/bin/STAR
vcfFileDir=/scratch/asiimwe/STAR-WASP_FASTQs_VCFs/VCF
baseFileDir=/scratch/asiimwe/STAR-WASP_FASTQs_VCFs/FASTQ
ulimit -n 10000 

STARdir=/home/asiimwe/projects/run_env/alpha_star_wasp_comparison/STAR_WASP_RunsHet/sample_id/xthreads

STARpar="--runThreadN x --genomeDir /home/asiimwe/projects/run_env/alpha_star_wasp_comparison/genome_directory/  --genomeLoad NoSharedMemory --outSAMtype BAM SortedByCoordinate --outSAMattributes NH HI AS nM NM MD jM jI  --alignEndsType EndToEnd --outSAMunmapped Within --outFilterMultimapNmax 1"

STARreads="--readFilesCommand gunzip -c --readFilesIn $baseFileDir/sample_id/R1.fastq.gz  $baseFileDir/sample_id/R2.fastq.gz"

WASPvcf="--varVCFfile $vcfFileDir/ssd_input_snp_dir/ssd.vcf.snv1het" 

$STAR $STARpar $STARreads $WASPvcf

mv $STARdir/Aligned.sortedByCoord.out.bam $STARdir/A_sorted.bam
samtools index $STARdir/A_sorted.bam $STARdir/A_sorted.bai
$PYTHON $WASP/mapping/find_intersecting_snps.py --is_paired_end --is_sorted --snp_dir $vcfFileDir/ssd_input_snp_dir/SNPdir --output_dir ./ A_sorted.bam 
$STAR $STARpar $WASPvcf --readFilesCommand gunzip -c --readFilesIn  A_sorted.remap.fq1.gz A_sorted.remap.fq2.gz  
mv Aligned.sortedByCoord.out.bam Aligned.out.bam
samtools index Aligned.out.bam Aligned.out.bai
$PYTHON $WASP/mapping/filter_remapped_reads.py A_sorted.to.remap.bam  Aligned.out.bam A_sorted.keep.bam