# Get list of contigs
bedtools genomecov -ibam ar3_10.2_gtdb_r89.bam -bg \
| awk '$4 > 9{split($0,a,"\\."); print a[1]}' \
| sort -u > ar3_10.2_gtdb_r89.acc
# Get
head -n10 ar3_10.2_gtdb_r89.acc \
| while read LINE; do
samtools view -b -L <(grep -w -F ${LINE} /vol/scratch/DBs/bowtie2/gtdb-r89_54k_genomes.fasta.fai \
| awk '{printf("%s\t0\t%s\n",$1,$2);}') ar3_10.2_gtdb_r89.bam \
| samtools view -h \
| python pmdtools.0.60.py --writesamfield --threshold -100 --header \
| samtools view -S \
| awk '{ print $1"\t"$NF }' \
| sed -e 's/DS:Z://g';
done > test_damage.txt