import smbl
import rnftools


##

#reference_basic_1   = smbl.fasta.HUMAN_HG38
reference_basic_1   = smbl.fasta.EXAMPLE_2
reference_contam_1  = smbl.fasta.EXAMPLE_1

reference_basic_2   = smbl.fasta.EXAMPLE_1
reference_contam_2  = smbl.fasta.EXAMPLE_2

##

reads_1 = "reads_a"
reads_2 = "reads_b"

##

bwa_aln_1 = "bwa_mem_1.bam"
bwa_aln_2 = "bwa_mem_2.bam"

##


#
# READS
#

# sample 1

rnftools.mishmash.sample(reads_1,ends=2)

rnftools.mishmash.DwgSim(
	fa=reference_basic_1,
	read_length_1=100,
	read_length_2=100,
	number_of_reads=10000
)


rnftools.mishmash.DwgSim(
	fa=reference_contam_1,
	read_length_1=100,
	read_length_2=100,
	number_of_reads=2000
)


# sample 2

rnftools.mishmash.sample(reads_2,ends=2)

rnftools.mishmash.DwgSim(
	fa=reference_basic_2,
	read_length_1=100,
	read_length_2=100,
	number_of_reads=10000
)


rnftools.mishmash.DwgSim(
	fa=reference_contam_2,
	read_length_1=100,
	read_length_2=100,
	number_of_reads=2000
)


#
# ALIGNMENT
#

bwa1 = smbl.prog.BwaMem(
	fasta=reference_basic_1,
	fastq_1=reads_1+".1.fq",
	fastq_2=reads_1+".2.fq",
	bam=bwa_aln_1,
)

bwa2 = smbl.prog.BwaMem(
	fasta=reference_basic_2,
	fastq_1=reads_2+".1.fq",
	fastq_2=reads_2+".2.fq",
	bam=bwa_aln_2,
)

#
# RULES
#

# required final files

rule basic:
	input:
		#rnftools.input(),
		bwa1.bam_fn(),
		bwa2.bam_fn(),



#
# BAM FILES
#

rule bwa_index_1:
	input:
		bwa1.make_index_input()
	output:
		bwa1.make_index_output()
	run:
		bwa1.make_index()

rule bwa_index_2:
	input:
		bwa2.make_index_input()
	output:
		bwa2.make_index_output()
	run:
		bwa2.make_index()

rule aln_bwa_1:
	input:
		bwa1.map_reads_input()
	output:
		bwa1.map_reads_output()
	run:
		bwa1.map_reads()

rule aln_bwa_2:
	input:
		bwa2.map_reads_input()
	output:
		bwa2.map_reads_output()
	run:
		bwa2.map_reads()


include: rnftools.include()
