def run_kraken_taxon():

    def download():
        kd = KrakenDownload()
        kd.download('toydb')
    database = sequana_config_path + os.sep + "kraken_toydb"
    if os.path.exists(database) is False:
        download()

    file1 = sequana_data("Hm2_GTGAAA_L005_R1_001.fastq.gz", "data")
    file2 = sequana_data("Hm2_GTGAAA_L005_R2_001.fastq.gz", "data")

    kt = KrakenAnalysis([file1, file2], database=database)
    kt.run()

    kt = KrakenAnalysis(file2, database=database)
    kt.run()

    p = tempfile.TemporaryDirectory()

    kt = KrakenHierarchical([file1, file2], [database, database],
            output_directory=p.name, force=True)
    kt.run()

    kt = KrakenHierarchical(file1, [database, database],
output_directory=p.name, force=True)
    kt.run()

    p.cleanup()

def test_snpeff():
    # a custom refrence
    fh_log = TempFile()

    mydata = snpeff.SnpEff(reference=sequana_data("JB409847.gbk"), log=fh_log.name)
    with TempFile() as fh:
        mydata.launch_snpeff(sequana_data("JB409847.vcf"), fh.name)
    fh_log.delete()

    # cleanup
    try:
        os.remove("snpEff.config")
    except:
        pass

    try:
        os.remove("snpEff_genes.txt")
    except:
        pass

    try:
        os.remove("snpEff_summary.html")
    except:
        pass

    try:
        snpeff.SnpEff(reference="dummy")
        assert False
    except SystemExit:
        assert True
    except:
        assert False

def test_pacbio_input_bam(tmpdir):
    # we need a summary and a bunch of images
    filename = sequana_data("summary_pacbio_qc1.json")

    # mock the PNG files found in the summary
    import json
    summary = json.load(open(filename))
    pngname = sequana_data("no_data.jpg")
    summary["images"]["gc_vs_length"] = pngname
    summary["images"]["hist_gc_content"] = pngname
    summary["images"]["hist_read_length"] = pngname
    summary["images"]["hist_snr"] = pngname
    summary["images"]["hist_zmw"] = pngname

    summary_file = TempFile()
    with open(summary_file.name, "w") as ff:
        json.dump(summary, ff)

    # Now that we have this new summary file, let us use it
    # we also need an output handler
    ff = TempFile()

    from sequana.utils import config
    config.output_dir = "/tmp"
    #here, ff.name is of the form /tmp/djhfjh4dz so we need to remove the /tmp
    pacbio_input_bam.PacbioInputBAMModule(summary_file.name, ff.name.split("/")[1])

    # cleanup
    summary_file.delete()
    ff.delete()

    def __init__(self, wk=None):
        super(VariantCallingPipeline, self).__init__(wk=wk)
        # Define the data
        data = sequana_data("Hm2_GTGAAA_L005_R1_001.fastq.gz")
        input_directory = os.path.dirname(data)
        self.input_pattern = input_directory + "/Hm*gz"
        self.pipeline = "variant_calling"

        # Define the project and config file
        subprocess.check_call([
            "sequana", "--pipeline", self.pipeline,
            "--input-pattern", '%s' % self.input_pattern,
            "--working-directory", self.wk, "--force"
            ])


        cmd = ["sequana", "--pipeline", self.pipeline,
             "--input-pattern", '%s'% self.input_pattern,
             "--working-directory", self.wk, "--force"]

        if "TRAVIS_PYTHON_VERSION" in os.environ:
             cmd += ["--snakemake-jobs", "1"]

        subprocess.check_call(cmd)


        # Add reference in the config
        cfg = SequanaConfig(self.wk + "/config.yaml")
        # We added a TTTT in position 5881
        cfg._yaml_code['bwa_mem_ref']['reference'] = sequana_data("measles.fa")
        cfg.save(self.wk + '/config.yaml')

def test_sequana_data():

    try:
        sequana_data()
        assert False
    except ValueError:
        assert True

def test_add_locus_with_modification():

    # Alter the original GBK to alter the locus name
    data = open(sequana_data("JB409847.gbk"), "r").read()
    newdata = data.replace("JB409847", "DUMMY_JB409847")

    fh = TempFile(suffix="gbk")
    with open(fh.name, 'w') as fout:
        fout.write(newdata)

    # Now we read this new GBK file that has a different locus name as
    # compared to the fasta
    mydata = snpeff.SnpEff(reference=fh.name)

    # Here is the corresponding FASTA
    fasta = sequana_data("JB409847.fasta")

    with TempFile(suffix="fasta") as fh2:
        mydata.add_locus_in_fasta(fasta, fh2.name)

        # In theory, in the newly created fasta file, we should find back the
        # DUMMY tag
        # cleanup
        try:
            os.remove("snpEff.config")
        except:
            pass

        data = open(fh2.name, "r").read()
        assert "DUMMY" in data
    fh.delete()

def test_is_sam_bam():

    datatest = sequana_data("test_measles.sam", "testing")
    assert is_sam(datatest) is True

    datatest = sequana_data("test_measles.bam", "testing")
    assert is_bam(datatest) is True

def test_input():
    filename = sequana_data('virus.bed', 'data')
    df = summary.main([prog, '--file', filename])
    len(df)

    filename = sequana_data('test.fastq', "testing")
    df = summary.main([prog, '--file', filename])

def test_ChromosomeCovMultiChunk():
    filename = sequana_data('JB409847.bed')
    # using chunksize of 7000, we test odd number
    bed = bedtools.GenomeCov(filename, sequana_data('JB409847.gbk'),chunksize=7000)
    chrom = bed.chr_list[0]
    res = chrom.run(501, k=2, circular=True)
    res.get_summary()
    res.get_rois()

def test_pcrfree():
    design = sequana_data("test_index_mapper.csv")

    try:
        ad = FindAdaptersFromDesign(design, "error")
        assert False
    except Exception:
        assert True

    # Other input from PCRFree
    ad = FindAdaptersFromDesign(design, "PCRFree")

    # Test the index1/2_seq with 3 cases
    # index1 present only,
    # no index at all (None)
    # index1 and index2 present
    design1 = sequana_data("test_expdesign_hiseq.csv")
    ad1 = FindAdaptersFromDesign(design1, "PCRFree")
    ad1.check()
    res1 = ad1.get_adapters_from_sample("553-iH2-1")
    res2 = ad1.get_adapters_from_sample("539-st2")
    res3 = ad1.get_adapters_from_sample("107-st2")

    assert res1['index1']['fwd'].identifier == "NextFlex_PCR_Free_adapter8|name:8|seq:TTAGGC"
    assert res1['index1']['fwd'].name == "8"
    assert res1['index1']['rev'].name == "8"

    assert list(res2.keys()) == ["universal"]

    assert res3['index1']['fwd'].name == "9"
    assert res3['index1']['rev'].name == "9"
    assert res3['index2']['fwd'].name == "10"
    assert res3['index2']['rev'].name == "10"

    # double indexing
    # This is a double indexing for PCRFree, which has not been tested
    # since it requires 16S adapters not yet in sequana
    """design2 = sequana_data("test_expdesign_miseq_illumina2.csv")
    ad2 = FindAdaptersFromDesign(design2, "PCRFree")
    assert ad2.get_adapters_from_sample('M2')['index1']['fwd'].identifier == \
            'NextFlex_PCR_Free_adapter2|name:2|seq:TGACCA'
    assert ad2.get_adapters_from_sample('M2')['index2']['fwd'].identifier == \
            'NextFlex_PCR_Free_adapter13|name:13|seq:AGTCAA'
    """


    design = sequana_data("test_expdesign_miseq_illumina.csv")
    ad = FindAdaptersFromDesign(design, "PCRFree")
    res = ad.get_adapters_from_sample("CR81-L1236-P1")
    assert res['index1']['fwd'].identifier == 'NextFlex_PCR_Free_adapter1|name:1|seq:CGATGT'



    design1 = sequana_data("test_expdesign_miseq_illumina_1.csv")
    ad = FindAdaptersFromDesign(design1, "PCRFree")
    ad.check() # all sample names must be found
    res = ad.get_adapters_from_sample("CM-2685")['index1']['fwd']
    assert res.name == "3"

def test_atropos_paired(tmpdir):
    # This is for the new version of cutadapt with version 1.1
    directory = tmpdir.mkdir('test_module')
    config.output_dir = str(directory)
    config.sample_name = 'JB409847'
    c = CutadaptModule(sequana_data('test_atropos_pe.txt'), "TEST", "test.html")
    assert c.jinja['mode'] == 'Paired-end'
    c = CutadaptModule(sequana_data('test_atropos_se.txt'), "TEST", "test.html")
    assert c.jinja['mode'] == 'Singled-end'

def test_sequana_data_star():
    # all files in a specific directory (a list)
    f1 = sequana_data("*", "images")
    assert isinstance(f1, list)
    assert 'Institut_Pasteur.png' in f1

    # all files (return a dict)
    f1 = sequana_data("*")
    assert isinstance(f1, dict)

def test_krona_merger():

    k1 = KronaMerger(sequana_data("test_krona_k1.tsv"))
    k2 = KronaMerger(sequana_data("test_krona_k2.tsv"))
    k1 += k2

    with TempFile(suffix='.tsv') as fh:
        df = k1.to_tsv(fh.name)
    assert all(df['count'] == [14043,591,184,132])
    assert k1['Bacteria\tProteobacteria\tspecies1\n'] == 14043

def test_bcf_filter():
    vcf_output_expected = sequana_data('JB409847.filter.vcf')
    bcf = BCF_freebayes(sequana_data('JB409847.bcf'))
    filter_dict = {'freebayes_score': 200, 'frequency': 0.85, 'min_depth': 10,
                   'forward_depth': 3, 'reverse_depth': 3, 'strand_ratio': 0.3}
    filter_bcf = bcf.filter_bcf(filter_dict)
    with TempFile(suffix='.vcf') as fp:
        filter_bcf.to_vcf(fp.name)
        compare_file = filecmp.cmp(fp.name, vcf_output_expected)
        assert compare_file

def test_add_locus_no_modification():
    mydata = snpeff.SnpEff(reference=sequana_data("JB409847.gbk"))
    with TempFile() as fh:
        fastafile = sequana_data("JB409847.fasta")
        mydata.add_locus_in_fasta(fastafile, fh.name)
        # cleanup
        try:
            os.remove("snpEff.config")
        except:
            pass

def test_sequana_data():

    try:
        sequana_data()
        assert False
    except ValueError:
        assert True
    except:
        assert False 

    sequana_data("Hm2_GTGAAA_L005_R1_001.fastq.gz", "data")

def test_coverage_module(tmpdir):

    bed = bedtools.GenomeCov(sequana_data("JB409847.bed"))
    fasta = sequana_data("JB409847.fasta")
    bed.compute_gc_content(fasta)
    c = bed.chr_list[0]
    c.run(4001)

    directory = tmpdir.mkdir('test_coverage_module')
    config.output_dir = str(directory)
    config.sample_name = "JB409847"
    CoverageModule(bed)

def test_adapter_reader():
    from sequana.adapters import AdapterReader as AR
    data = sequana_data("adapters_with_duplicates.fa", "testing")
    try:
        ar = AR(data)
        ar.sanity_check()
    except ValueError:
        pass

    data1 = sequana_data("adapters_Nextera_fwd.fa")
    data2 = sequana_data("adapters_Nextera_rev.fa")
    data3 = sequana_data("adapters_Nextera_revcomp.fa")

    # try different constructors
    ar1 = AR(data1)
    ar1.__repr__()
    ar1.index_sequences
    ar_same = AR(ar1._data)    # from a list of dictionaries
    assert ar1 == ar_same
    ar_same = AR(ar1)           # from a AR instance
    assert ar1 == ar_same
    assert ar1[0]['identifier'] == 'Universal_Adapter|name:universal'
    ar1.index_names
    assert ar1.get_adapter_by_sequence("XXX") is None
    try:
        ar1.get_adapter_by_identifier("XXX")
        assert False
    except ValueError:
        assert True

    # __eq__
    assert len(ar1) == 56

    # accessors
    ar1.sequences, ar1.identifiers, ar1.comments

    ar1.get_adapter_by_sequence("ACGT")
    assert ar1.get_adapter_by_index_name("dummy") is None
    assert ar1.get_adapter_by_identifier("Nextera_index_N517")

    ar2 = AR(data2)
    ar2.reverse()
    # fails due to S516 ????????
    assert ar1 == ar2

    ar3 = AR(data3)
    ar3.reverse_complement()
    assert ar1 == ar3

    # test to_fasta method
    with TempFile() as fh:
        ar1.to_fasta(fh.name)

def test_canvasjs_linegraph():
    bed = bedtools.GenomeCov(sequana_data("JB409847.bed"))
    fasta = sequana_data("JB409847.fasta")
    bed.compute_gc_content(fasta)


    c = bed.chr_list[0]
    c.run(4001)

    df = bed[0].df
    csv = df.to_csv(columns=['pos', 'cov', 'gc'], index=False,
                    float_format='%.3g')
    # create CanvasJS stuff
    cjs = CanvasJSLineGraph(csv, 'cov', 'pos', ['cov', 'gc'])
    # set options
    cjs.set_options({'zoomEnabled': 'true',
                     'zoomType': 'x',
                     'exportEnabled': 'true'})
    # set title
    cjs.set_title("Genome Coverage")
    # set legend
    cjs.set_legend({'verticalAlign': 'bottom',
                    'horizontalAlign': 'center',
                    'cursor':'pointer'},
                    hide_on_click=True)
    # set axis
    cjs.set_axis_x({'title': "Position (bp)",
                    'labelAngle': 30,
                    'minimum': 0,
                    'maximum': len(df)})
    cjs.set_axis_y({'title': "Coverage (Count)"})
    cjs.set_axis_y2({'title': "GC content (ratio)",
                     'minimum':0,
                     'maximum': 1,
                     'lineColor': '#FFC425',
                     'titleFontColor': '#FFC425',
                     'labelFontColor': '#FFC425'})
    # set datas
    cjs.set_data(index=0, data_dict={'type': 'line',
                                     'name': "Coverage",
                                     'showInLegend': 'true',
                                     'color': '#5BC0DE',
                                     'lineColor': '#5BC0DE'})
    cjs.set_data(index=1, data_dict={'type': 'line',
                                     'axisYType': 'secondary',
                                     'name': "GC content",
                                     'showInLegend': 'true',
                                     'color': '#FFC425',
                                     'lineColor': '#FFC425'})
    # create canvasJS
    cjs.create_canvasjs()

def test_design_constructor():
    filename = sequana_data("test_expdesign_miseq_illumina.csv")
    tt = ExpDesignMiSeq(filename)

    # using existing design
    tt = ExpDesignAdapter(tt)
    # or a filename
    tt = ExpDesignAdapter(filename)

    # constructor for hiseq
    tt = ExpDesignAdapter(sequana_data("test_expdesign_hiseq.csv"))

    tt = ExpDesignAdapter(sequana_data("test_expdesign_generic.csv"))
    tt
    print(tt)