from CoaSim import *
import sys

if len(sys.argv) != 4:
    print 'Usage:', argv[0], 'beta noSeqs noSamples'
    sys.exit(2)

beta      = float(sys.argv[1])
noSeqs    = int(sys.argv[2])
noSamples = int(sys.argv[3])

class callback(object):
    def __init__(self):
        self.count = 0
    def recombinationEvent(self, *dummy):
        self.count += 1

# mapping seqs*beta -> rho
rho = {(20,0): 0.28, (50,0): 0.22, (20,100): 3.88, (50,100): 2.15}[(noSeqs,beta)]

counts = dict()
for i in xrange(noSamples):
    cb = callback()
    arg = simulate([], noSeqs, beta=beta, rho=rho,
                   callbacks=cb,
                   keepEmptyIntervals=True)
    try:
        counts[cb.count] += 1
    except KeyError:
        counts[cb.count] = 1


observedCounts = counts.keys()
observedCounts.sort()
total = sum(counts.values())
maxFrac = max([float(counts[oc])/total for oc in observedCounts])
for oc in observedCounts:
    print '%d: %6.2f%%' % (oc, 100.0*counts[oc]/total),
    print '*'*int(300*maxFrac*counts[oc]/total)