+ ACC=CP000253
+ SRR=SRR397558
+ N=250000
+ mkdir -p ref
+ GB=ref/CP000253.gb
+ REF=ref/CP000253.fa
+ GFF=ref/CP000253.gff
+ efetch -db nuccore -id CP000253 -format gb
+ cat ref/CP000253.gb
+ readseq -p -f fasta
+ cat ref/CP000253.gb
+ readseq -p -f gff
+ grep CDS
+ bwa index ref/CP000253.fa
+ samtools faidx ref/CP000253.fa
+ mkdir -p reads
+ fastq-dump -X 250000 --origfmt --split-files -O reads SRR397558
+ echo '*** Sequence statistics'
+ seqkit stat reads/SRR397558_1.fastq reads/SRR397558_2.fastq
+ echo '*** '
+ SIMR1=reads/sim1.fastq
+ SIMR2=reads/sim2.fastq
+ wgsim -e 0 -r 0 -R 0 -X 0 -N 250000 ref/CP000253.fa reads/sim1.fastq reads/sim2.fastq
+ bwa mem ref/CP000253.fa reads/sim1.fastq
+ samtools sort
+ samtools index sim1.bam
+ R1=reads/SRR397558_1.fastq
+ R2=reads/SRR397558_2.fastq
+ BAM=real1.bam
+ bwa mem ref/CP000253.fa reads/SRR397558_1.fastq
+ samtools sort
+ samtools index real1.bam
+ cat ref/CP000253.gff
+ awk '$7 == "+" { print $0 }'
+ bedtools intersect -a real1.bam -b forward_cds.gff
+ samtools view -b -F 4 -f 16 forward_overlap.bam
+ samtools view -b -F 4 -F 16 forward_overlap.bam
+ xargs -n 1 samtools index
+ ls -1 forward_antisense.bam forward_overlap.bam forward_sense.bam
+ bedtools genomecov -ibam forward_sense.bam -bg
+ bedtools genomecov -ibam forward_antisense.bam -bg