diff --git a/libraries/base/common/python/try_round_weight.py b/libraries/base/common/python/try_round_weight.py
index 7b8a3ec60a04224192f8610714d77899e66a34e3..b14463638e288fc3d3b855875cb76d49af1956c6 100644
--- a/libraries/base/common/python/try_round_weight.py
+++ b/libraries/base/common/python/try_round_weight.py
@@ -28,16 +28,27 @@
# . quantized subbands --> sigma_qq
# . unquantized subbands --> sigma_sq
# Preliminary conclusion:
-# . for small input noise with sigma < 2 the output sigma gets disturbed
+# . For small input noise with sigma < 2 the output sigma gets disturbed
# due to the weighting if the weighting is applied after the subband
# quantisation
-# . increasing -N improves the results, for LOFAR subbands N = 195312
+# . Increasing -N improves the results, for LOFAR subbands N = 195312
# . it may be preferred to apply the subband weights to the unquantized
# WPFB output.
+# . Choosing sufficient intermediate resolution (-r) before applying
+# weights:
+# - Rounding noise changes the sigma of the noise anyway, as shown by
+# sigmas_ratio_sq_input_T. Therefore choose resolution such that
+# jumps in sigma due to rounding weighted noise, as shown by
+# sigmas_ratio_qq_sq_T, is much smaller than sigmas_ratio_sq_input_T
+# ==> -r 2.
+# - For input sigma >= 1
+# . Choose sigmas_ratio_qq_sq_T < 10% ==> -r 3
+# . Choose sigmas_ratio_qq_sq_T < 1% ==> -r 6
+# . Choose sigmas_ratio_qq_sq_T < 0.1% ==> -r 9
# Note:
-# . For values exactly halfway between rounded decimal values, NumPy rounds to
-# the nearest even value. Thus 1.5 and 2.5 round to 2.0, -0.5 and 0.5 round
-# to 0.0, etc.
+# . For values exactly halfway between rounded decimal values, NumPy of
+# Python3 rounds to the nearest even value. Thus 1.5 and 2.5 round to 2.0,
+# -0.5 and 0.5 round to 0.0, etc. (Python2 rounds half away from zero).
# Usage:
# > python3 try_round_weight.py -N 195312
@@ -74,6 +85,7 @@ _parser = argparse.ArgumentParser(
# Use -r = 6 to see effect of having more resolution before rounding
> python try_round_weight.py --w_lo 0.7 --w_hi 0.8 --w_step 0.0001 --s_lo 1 --s_hi 10 --s_step 1 -N 195312 -S 0 -r 6
+ > python try_round_weight.py --w_lo 0.3 --w_hi 1.1 --w_step 0.001 --s_lo 1 --s_hi 10 --s_step 1 -N 195312 -S 0 -r 6
\n""")),
formatter_class=argparse.RawTextHelpFormatter)
_parser.add_argument('-S', default=0, type=int, help='Random number seed')
@@ -118,7 +130,8 @@ resolution = args.r
resolution_factor = 2**resolution
# Determine weighted rounded noise sigma / weighted noise sigma for range of weights and input noise sigmas
-sigmas_ratio = np.nan * np.zeros((N_weights, N_sigmas)) # w rows, s cols
+sigmas_ratio_qq_sq = np.nan * np.zeros((N_weights, N_sigmas)) # w rows, s cols
+sigmas_ratio_sq_input = np.nan * np.zeros((N_weights, N_sigmas)) # w rows, s cols
sigmas_qq = np.zeros((N_weights, N_sigmas))
sigmas_sq = np.zeros((N_weights, N_sigmas))
for s, sigma in enumerate(sigmas):
@@ -143,21 +156,25 @@ for s, sigma in enumerate(sigmas):
sigmas_qq[w][s] = s_qq
sigmas_sq[w][s] = s_sq
if s_sq != 0:
- sigmas_ratio[w][s] = s_qq / s_sq # weighted rounded noise sigma / weighted noise sigma
+ sigmas_ratio_qq_sq[w][s] = s_qq / s_sq # weighted rounded noise sigma / weighted noise sigma
+ sigmas_ratio_sq_input[w][s] = s_sq / sigma # weighted noise sigma / input noise sigma
# Transpose [w][s] to have index ranges [s][w]
-sigmas_ratio_T = sigmas_ratio.transpose()
+sigmas_ratio_qq_sq_T = sigmas_ratio_qq_sq.transpose()
+sigmas_ratio_sq_input_T = sigmas_ratio_sq_input.transpose()
sigmas_qq_T = sigmas_qq.transpose()
sigmas_sq_T = sigmas_sq.transpose()
# Plot results
figNr = 0
+s_colors = plt.cm.jet(np.linspace(0, 1, N_sigmas))
+w_colors = plt.cm.jet(np.linspace(0, 1, N_weights))
figNr += 1
plt.figure(figNr)
for s, sigma in enumerate(sigmas):
# Plot sigma_qq of twice quantized noise as function of weight for
# different input sigmas
- plt.plot(weights, sigmas_qq_T[s], label='s = %4.2f' % sigma)
+ plt.plot(weights, sigmas_qq_T[s], color=s_colors[s], label='s = %4.2f' % sigma)
plt.title("Sigma of weighted quantized noise")
plt.xlabel("Weight")
plt.ylabel("Sigma_qq")
@@ -171,14 +188,28 @@ for s, sigma in enumerate(sigmas):
# different input sigmas.
# Normalize the sigma_qq by the weight, so that it can be compared with
# the input sigma that is shown by the horizontal sigma reference lines.
- plt.plot(weights, sigmas_qq_T[s] / weights, label='s = %4.2f' % sigma)
- plt.plot(weights, sigmas[s]*np.ones(N_weights)) # add sigma reference lines
+ plt.plot(weights, sigmas_qq_T[s] / weights, color=s_colors[s], label='s = %4.2f' % sigma)
+ plt.plot(weights, sigmas[s]*np.ones(N_weights), '--', color=s_colors[s]) # add sigma reference lines
plt.title("Sigma of weighted quantized noise, normalized for weight")
plt.xlabel("Weight")
plt.ylabel("Sigma_qq")
plt.legend(loc='upper right')
plt.grid()
+figNr += 1
+plt.figure(figNr)
+for s, sigma in enumerate(sigmas):
+ # Plot ratio of sigma_sq / sigma as function of weight for different
+ # input sigma. The ratio deviation from 1 tells how much quantized
+ # weighted noise deviates from the input noise. This shows that rounding
+ # noise cause a change in sigma even when weight is 1.
+ plt.plot(weights, sigmas_ratio_sq_input_T[s] / weights, color=s_colors[s], label='s = %4.2f' % sigma)
+plt.title("Relative sigma difference of weighted quantized data / input data")
+plt.xlabel("Weight")
+plt.ylabel("Relative sigma difference")
+plt.legend(loc='upper right')
+plt.grid()
+
figNr += 1
plt.figure(figNr)
for s, sigma in enumerate(sigmas):
@@ -186,8 +217,8 @@ for s, sigma in enumerate(sigmas):
# input sigma. The ratio deviation from 1 tells how much the twice
# quantized noise deviates from the noise that is only quantized after
# the weighting.
- plt.plot(weights, sigmas_ratio_T[s], label='s = %4.2f' % sigma)
-plt.title("Relative sigma difference of weighting after / before quantisation")
+ plt.plot(weights, sigmas_ratio_qq_sq_T[s], color=s_colors[s], label='s = %4.2f' % sigma)
+plt.title("Relative sigma difference of weighting before / after quantisation")
plt.xlabel("Weight")
plt.ylabel("Relative sigma difference")
plt.legend(loc='upper right')
@@ -198,8 +229,8 @@ plt.figure(figNr)
for w, weight in enumerate(weights):
# Plot ratio of sigma_qq / sigma_sq as function of input sigma for
# different weights
- plt.plot(sigmas, sigmas_ratio[w], label='w = %4.2f' % weight)
-plt.title("Relative sigma difference of weighting after / before quantisation")
+ plt.plot(sigmas, sigmas_ratio_qq_sq[w], color=w_colors[w], label='w = %4.2f' % weight)
+plt.title("Relative sigma difference of weighting before / after quantisation")
plt.xlabel("Sigma")
plt.ylabel("Relative sigma difference (s_qq / s_sq)")
plt.legend(loc='upper right')
@@ -207,7 +238,7 @@ plt.grid()
figNr += 1
plt.figure(figNr)
-plt.imshow(sigmas_ratio, origin='lower', interpolation='none', aspect='auto', extent=[sigma_lo, sigma_hi, weight_lo, weight_hi])
+plt.imshow(sigmas_ratio_qq_sq, origin='lower', interpolation='none', aspect='auto', extent=[sigma_lo, sigma_hi, weight_lo, weight_hi])
plt.colorbar()
plt.title("Relative sigma difference of weighting after / before quantisation")
plt.xlabel("Sigma")