Histogram watermarks

[janniklasrose]

For a random process with bin count/value of $c$, the estimated noise is $\pm\sqrt{c}$ [REF]. At 500 walkers per bin, this should lead to fluctuations of around 22.

N_per_bin = NP/(binedges.size-1)
print("N_per_bin = %d -- fluctuation [sqrt(N)] = %g" % (N_per_bin, math.sqrt(N_per_bin)))

histograms:

## walker positions
bins_timeaveraged = oneD.get_stats(dat1['binning_global'], 1)  # averaged over time
bins_realisationaveraged = oneD.get_stats(dat1['binning_global'], 0)  # average over realisations

# show
fig = plt.figure(figsize=(20, 5))

plt.subplot(121)
rep = 0  # only one realisation
binedges = dat1['binedges_global']
plt.hist(binedges[:-1], binedges, weights=bins_timeaveraged['min'][rep], histtype='step')  # watermark (lo)
plt.hist(binedges[:-1], binedges, weights=bins_timeaveraged['max'][rep], histtype='step')  # watermark (hi)
plt.hist(binedges[:-1], binedges, weights=bins_timeaveraged['avg'][rep], histtype='step')
plt.xlabel('$x$ (µm)')
plt.ylabel('count')
plt.legend(('lo','hi','avg',))
plt.title('time-averaged histogram (and watermarks) for one realisation')

plt.subplot(122)
for n in range(0, NT, 10):  # 10 lines
    plt.hist(binedges[:-1], binedges, weights=bins_realisationaveraged['avg'][n], histtype='step')
plt.xlabel('$x$ (µm)')
plt.ylabel('count')
plt.title('realisation-averaged histogram for 10 time points')
plt.show()