global name 'freseq' is not defined

Hi, I've been enjoying using your software and I think it will be perfect for testing my own software but I have come across a few issues.

One I haven't been able to solve is when using multiple hap files the code is as follows:

if len(hapfiles) == 1: #reuse hapfiles for haplotype hapseq=[pysam.Fastafile(hapfiles[0])]*nploidy elif len(hapfiles)==nploidy: hapseq=[pysam.Fastafile(hapfiles[i]) for i in range(0,nploidy)] for ref in freseq: assert ref.lengths == hapseq[0].lengths, \ "hapfile seqs must be of the same length cross haplos" else: raise Exception("incorrect number of hapfiles provided")

However 'freseq' does not appear anywhere else in the code - may just be a typo but I'm unsure what the intended object was.

Also, when inputting multiple hapfiles the manual says these should be separated by a comma but in the code they are separated by ':' - it would be useful if this was mentioned in the help option (-h).

Finally, when defining a distribution in the metafile an Assertion error was raised for the line: assert new_cnt == cnt

I added the code used when using 'fix_' as follows, just wanted to check if this was an okay fix? elif dist[0] == 'nori': #create cnt normal integer; mean = int(dist[1]) sd = int(dist[2]) values = [ random.gauss(mean, sd) for i in xrange(0,cnt) ] values = [ mean if v <= 0 else int(v) for v in values ] uniques = CountUnique(values) #ADDED BY ARODDY sorted_uniques = sorted(uniques.items(), key=operator.itemgetter(1)) #increasing order sorting of values nvalue = len(values) nunique = len(uniques) amp = cnt/nvalue rem = cnt-amp*nvalue for item in sorted_uniques: key = item[0] if rem>0: delta = uniques[key] if rem-uniques[key]>=0 else rem uniques[key] = uniques[key] * amp + delta rem = rem - delta else: delta = 0 uniques[key] = uniques[key] * amp new_cnt = new_cnt + uniques[key] ### END OF ADDED CODE assert new_cnt == cnt

Thanks for your help, Aideen

Comments (2)