Skip to content
Snippets Groups Projects
Commit 8d032924 authored by Thomas Faraut's avatar Thomas Faraut
Browse files

merging works now for a single file

parent 0de3d4d1
No related branches found
No related tags found
No related merge requests found
......@@ -106,20 +106,6 @@ def get_contigs(reference):
return contigs
# def chromNum(chrom):
# m = re.search(r"[cC]hr(?P<num>[\dXY]+)", chrom)
# if m is not None:
# if m.group('num') == 'X':
# num = 100
# elif m.group('num') == 'Y':
# num = 101
# else:
# num = int(m.group('num'))
# return num
# else:
# return chrom
def addContigInfosTovcf(vcffile, outvcffile, reference):
"""
Add contig infos to vcf
......@@ -150,8 +136,6 @@ def addContigInfosTovcf(vcffile, outvcffile, reference):
for ctg in contigs:
fout.write('##contig=<ID=%s,length=%d>\n'
% (ctg.name, ctg.length))
# for ctg in sorted(contigs,
# key=lambda x: chromNum(x.name)):
fout.write(line)
# bgzip file
......
......@@ -62,22 +62,25 @@ def genotype_multiple_samples(bamlist, vcf_in, vcf_out, cores=1):
pool = multiprocessing.Pool(processes=cores)
launch_ind = partial(genotype_single_sample,
vcf_in=vcf_in, out_dir=tmp_dir)
vcf_files = pool.map(launch_ind, bamfiles)
merge_cmd = "bcftools merge -m id "
if vcf_out != sys.stdout:
merge_cmd += "-O z " + " -o " + vcf_out + " "
merge_cmd += " ".join(vcf_files)
# Now merge result files
if len(bamfiles) == 1:
if vcf_out != sys.stdout:
cmd = "cp %s %s" % (vcf_files[0], vcf_out)
else:
cmd = "bcftools view %s" % vcf_files[0]
else:
cmd = "bcftools merge -m id "
if vcf_out != sys.stdout:
cmd += "-O z " + " -o " + vcf_out + " "
cmd += " ".join(vcf_files)
exit_code = call(merge_cmd, shell=True)
exit_code = call(cmd, shell=True)
if exit_code == 0:
if vcf_out != sys.stdout:
tabix_index(vcf_out, force=True, preset="vcf")
shutil.rmtree(tmp_dir)
else:
print("Failed: bcftools merge exits with status %d" % exit_code)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment