From: Heng Li <lh3@live.co.uk>
Date: Fri, 26 Nov 2010 05:07:31 +0000 (+0000)
Subject:  * fixed a memory violation
X-Git-Url: https://git.donarmstrong.com/?p=samtools.git;a=commitdiff_plain;h=1a302a67a7f11225d1f523fa431fcf8ef5a6261b

 * fixed a memory violation
 * added splitchr to vcfutils.pl
---

diff --git a/bam2bcf_indel.c b/bam2bcf_indel.c
index 819c1e7..16241d0 100644
--- a/bam2bcf_indel.c
+++ b/bam2bcf_indel.c
@@ -112,7 +112,7 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
 					  const void *rghash)
 {
 	extern void ks_introsort_uint32_t(int, uint32_t*);
-	int i, s, j, k, t, n_types, *types, max_rd_len, left, right, max_ins, *score1, *score2;
+	int i, s, j, k, t, n_types, *types, max_rd_len, left, right, max_ins, *score1, *score2, max_ref2;
 	int N, K, l_run, ref_type, n_alt;
 	char *inscns = 0, *ref2, *query;
 	khash_t(rg) *hash = (khash_t(rg)*)rghash;
@@ -236,7 +236,8 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
 		free(inscns_aux);
 	}
 	// compute the likelihood given each type of indel for each read
-	ref2  = calloc(right - left + max_ins + 2, 1);
+	max_ref2 = right - left + 2 + 2 * (max_ins > -types[0]? max_ins : -types[0]);
+	ref2  = calloc(max_ref2, 1);
 	query = calloc(right - left + max_rd_len + max_ins + 2, 1);
 	score1 = calloc(N * n_types, sizeof(int));
 	score2 = calloc(N * n_types, sizeof(int));
@@ -264,6 +265,7 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
 		}
 		for (; j < right && ref[j]; ++j)
 			ref2[k++] = bam_nt16_nt4_table[bam_nt16_table[(int)ref[j]]];
+		for (; k < max_ref2; ++k) ref2[k] = 4;
 		if (j < right) right = j;
 		// align each read to ref2
 		for (s = K = 0; s < n; ++s) {
diff --git a/bcftools/vcfutils.pl b/bcftools/vcfutils.pl
index 6cc168f..cd86b0f 100755
--- a/bcftools/vcfutils.pl
+++ b/bcftools/vcfutils.pl
@@ -14,11 +14,27 @@ sub main {
   my $command = shift(@ARGV);
   my %func = (subsam=>\&subsam, listsam=>\&listsam, fillac=>\&fillac, qstats=>\&qstats, varFilter=>\&varFilter,
 			  hapmap2vcf=>\&hapmap2vcf, ucscsnp2vcf=>\&ucscsnp2vcf, filter4vcf=>\&varFilter, ldstats=>\&ldstats,
-			  gapstats=>\&gapstats);
+			  gapstats=>\&gapstats, splitchr=>\&splitchr);
   die("Unknown command \"$command\".\n") if (!defined($func{$command}));
   &{$func{$command}};
 }
 
+sub splitchr {
+  my %opts = (l=>5000000);
+  getopts('l:', \%opts);
+  my $l = $opts{l};
+  die(qq/Usage: vcfutils.pl splitchr [-l $opts{l}] <in.fa.fai>\n/) if (@ARGV == 0 && -t STDIN);
+  while (<>) {
+	my @t = split;
+	my $last = 0;
+	for (my $i = 0; $i < $t[1];) {
+	  my $e = ($t[1] - $i) / $l < 1.1? $t[1] : $i + $l;
+	  print "$t[0]:".($i+1)."-$e\n";
+	  $i = $e;
+	}
+  }
+}
+
 sub subsam {
   die(qq/Usage: vcfutils.pl subsam <in.vcf> [samples]\n/) if (@ARGV == 0);
   my ($fh, %h);