ex-6: wrote something about ex-6 in the readme file and put in order ex.6/plots

ex-6/plots/emd-round.py renamed into ex-6/plots/emd.py and made interactive.
This commit is contained in:
Giù Marcer 2020-05-22 15:27:16 +02:00 committed by rnhmjoj
parent ea6fccfa5a
commit 76d5f3bebd
9 changed files with 210 additions and 221 deletions

View File

@ -162,19 +162,42 @@ To plot the results do:
`ex-6/bin/main` simulates a Fraunhöfer diffraction experiment. The program
prints to `stdout` the bin counts of the intensity as a function of the
diffraction angle. To plot a histogram simply pipe the output to the
program `ex-6/plot.py`.
diffraction angle. To plot a histogram do:
$ ex-6/bin/main | ex-6/plot.py
The program convolves the original signal with a Gaussian kernel (`-s` to
change the σ), optionally adds a Poisson noise (`-m` to change the mean μ) and
performs either a naive deconvolution by a FFT (`-m fft` mode) or applying the
Richard-Lucy deconvolution algorithm (`-m rl` mode), which is expected to
perform optimally in this case.
change the kernel σ), optionally adds a Gaussian noise (`-n` to change the
noise σ) and performs either a naive deconvolution by a FFT (`-m fft` mode)
or applying the Richardson-Lucy deconvolution algorithm (`-m rl` mode).
The `-c` and `-d` options control whether the convolved or deconvolved
histogram counts should be printed to `stdout`. For more options
The `-o`, `-c` and `-d` options control whether the original, convolved or
deconvolved histogram counts should be printed to `stdout`. For more options
run the program with `-h` to see the usage screen.
`ex-6/bin/test` simulates a customizable number of experiments and prints
to `stdout` the histograms of the distribution of the EMD of the convolved
signal together with the deconvolved signal EMD histograms with both FFT and
Richardson-Lucy procedures. It also prints to `stderr` the average, standard
deviation and skewness of the histograms. To plot the results, do:
$ ex-6/bin/test | ex-6/dist-plot.py
The program accepts some parameters to control the histogram and number of
events, run it with `-h` to see their usage.
`ex-6/plots/emd.py` plots the content of the files `ex-6/plots/emd-noisy.txt`
and `ex-6/plots/emd-noiseless.txt` depending on the argument passed to it from
stdin. Do:
$ ex-6/plots/emd.py 'noisy'
to plot the content of the first file and do:
$ ex-6/plots/emd.py 'noiseless'
to plot the content of the second file.
### Exercise 7

View File

@ -30,6 +30,4 @@ hist(a[2*n:], insert(b[2*n:], 0, a[2*n]), weights=f[2*n:],
ticklabel_format(style='sci', axis='x', scilimits=(0, 0), useMathText=True)
tight_layout()
name = (sys.argv[1] if len(sys.argv) > 1 else "prova")
savefig('notes/images/' + name + '.pdf' )
show()

View File

@ -14,5 +14,4 @@ xlabel(r'$\theta$ (radians)')
ylabel(r'$I(\theta)$ (a.u.)')
tight_layout()
name = sys.argv[2] if len(sys.argv) > 2 else "prova"
savefig('notes/images/' + name + '.pdf' )
show()

View File

@ -0,0 +1,75 @@
#rounds mean stddev
#σ = 1
1 1.88e-04 1.4e-05
2 3.17e-04 1.4e-05
3 1.80e-04 1.4e-05
4 2.02e-04 1.5e-05
5 1.70e-04 1.4e-05
6 1.64e-04 1.4e-05
7 1.75e-04 1.4e-05
8 1.67e-04 1.4e-05
9 1.61e-04 1.4e-05
10 1.57e-04 1.4e-05
11 1.59e-04 1.4e-05
12 1.56e-04 1.4e-05
13 1.54e-04 1.4e-05
14 1.53e-04 1.4e-05
15 1.51e-04 1.4e-05
16 1.50e-04 1.3e-05
17 1.49e-04 1.3e-05
18 1.48e-04 1.3e-05
19 1.46e-04 1.3e-05
20 1.47e-04 1.3e-05
40 1.35e-04 1.3e-05
60 1.28e-04 1.3e-05
80 1.24e-04 1.3e-05
100 1.21e-04 1.3e-05
400 1.02e-04 1.2e-05
800 9.37e-05 1.1e-05
1600 8.61e-05 1.0e-05
3800 7.98e-05 9.2e-06
7600 7.97e-05 9.2e-06
10000 8.06e-05 9.2e-06
15000 8.16e-05 9.2e-06
#σ = 0.5
1 1.17e-04 8.0e-06
2 2.47e-05 3.0e-06
3 3.65e-05 3.8e-06
4 5.85e-05 5.4e-06
5 1.48e-05 3.3e-06
6 1.82e-05 3.0e-06
7 1.22e-05 3.8e-06
8 1.31e-05 3.6e-06
9 1.19e-05 3.9e-06
10 1.17e-05 3.9e-06
11 1.16e-05 4.0e-06
12 1.16e-05 4.0e-06
13 1.16e-05 4.0e-06
14 1.16e-05 4.0e-06
15 1.16e-05 4.0e-06
16 1.16e-05 4.0e-06
17 1.16e-05 4.0e-06
18 1.16e-05 4.0e-06
19 1.16e-05 4.0e-06
20 1.16e-05 4.0e-06
#σ = 0.1
1 4.23e-16 3.1e-16
2 4.23e-16 3.2e-16
3 4.24e-16 3.1e-16
4 4.24e-16 3.1e-16
5 4.25e-16 3.2e-16
6 4.25e-16 3.2e-16
7 4.24e-16 3.1e-16
8 4.24e-16 3.1e-16
9 4.24e-16 3.1e-16
10 4.25e-16 3.2e-16
11 4.25e-16 3.2e-16
12 4.24e-16 3.1e-16
13 4.25e-16 3.2e-16
14 4.24e-16 3.1e-16
15 4.25e-16 3.2e-16
16 4.24e-16 3.1e-16
17 4.25e-16 3.2e-16
18 4.24e-16 3.1e-16
19 4.25e-16 3.2e-16
20 4.24e-16 3.1e-16

64
ex-6/plots/emd-noisy.txt Normal file
View File

@ -0,0 +1,64 @@
#rounds mean stddev
#σ=0.005
1 2.92e-04 4.7e-05
2 2.18e-04 4.7e-05
3 2.28e-04 4.5e-05
4 2.12e-04 4.8e-05
5 2.08e-04 4.8e-05
6 2.02e-04 4.9e-05
7 2.05e-04 4.9e-05
8 2.00e-04 4.9e-05
9 1.98e-04 5.0e-05
10 1.95e-04 5.0e-05
11 1.97e-04 5.0e-05
12 1.94e-04 5.0e-05
13 1.93e-04 5.0e-05
14 1.92e-04 5.0e-05
15 1.91e-04 5.1e-05
16 1.90e-04 5.1e-05
17 1.90e-04 5.1e-05
18 1.89e-04 5.1e-05
19 1.88e-04 5.1e-05
20 1.88e-04 5.1e-05
#σ = 0.01
1 3.47e-04 1.0e-04
2 3.96e-04 1.1e-04
3 3.39e-04 1.1e-04
4 3.42e-04 1.1e-04
5 3.34e-04 1.1e-04
6 3.35e-04 1.1e-04
7 3.37e-04 1.1e-04
8 3.33e-04 1.1e-04
9 3.33e-04 1.1e-04
10 3.32e-04 1.1e-04
11 3.32e-04 1.1e-04
12 3.31e-04 1.1e-04
13 3.31e-04 1.1e-04
14 3.31e-04 1.1e-04
15 3.31e-04 1.1e-04
16 3.31e-04 1.1e-04
17 3.31e-04 1.1e-04
18 3.31e-04 1.1e-04
19 3.31e-04 1.1e-04
20 3.31e-04 1.1e-04
#σ=0.05
1 1.49e-03 5.7e-04
2 1.51e-03 5.6e-04
3 1.50e-03 5.7e-04
4 1.50e-03 5.6e-04
5 1.52e-03 5.6e-04
6 1.53e-03 5.6e-04
7 1.52e-03 5.6e-04
8 1.51e-03 5.6e-04
9 1.54e-03 5.5e-04
10 1.53e-03 5.5e-04
11 1.53e-03 5.5e-04
12 1.54e-03 5.5e-04
13 1.54e-03 5.5e-04
14 1.54e-03 5.5e-04
15 1.55e-03 5.5e-04
16 1.55e-03 5.5e-04
17 1.55e-03 5.5e-04
18 1.56e-03 5.5e-04
19 1.56e-03 5.4e-04
20 1.56e-03 5.4e-04

View File

@ -1,75 +0,0 @@
#rounds mean stddev skew
#σ = 1
1 1.88e-04 1.4e-05 0.26
2 3.17e-04 1.4e-05 0.19
3 1.80e-04 1.4e-05 0.27
4 2.02e-04 1.5e-05 0.21
5 1.70e-04 1.4e-05 0.27
6 1.64e-04 1.4e-05 0.27
7 1.75e-04 1.4e-05 0.27
8 1.67e-04 1.4e-05 0.27
9 1.61e-04 1.4e-05 0.26
10 1.57e-04 1.4e-05 0.25
11 1.59e-04 1.4e-05 0.26
12 1.56e-04 1.4e-05 0.25
13 1.54e-04 1.4e-05 0.25
14 1.53e-04 1.4e-05 0.25
15 1.51e-04 1.4e-05 0.25
16 1.50e-04 1.3e-05 0.25
17 1.49e-04 1.3e-05 0.25
18 1.48e-04 1.3e-05 0.25
19 1.46e-04 1.3e-05 0.25
20 1.47e-04 1.3e-05 0.25
40 1.35e-04 1.3e-05 0.29
60 1.28e-04 1.3e-05 0.31
80 1.24e-04 1.3e-05 0.33
100 1.21e-04 1.3e-05 0.35
400 1.02e-04 1.2e-05 0.40
800 9.37e-05 1.1e-05 0.36
1600 8.61e-05 1.0e-05 0.27
3800 7.98e-05 9.2e-06 0.26
7600 7.97e-05 9.2e-06 0.28
10000 8.06e-05 9.2e-06 0.26
15000 8.16e-05 9.2e-06 0.25
#σ = 0.5
1 1.17e-04 8.0e-06 0.25
2 2.47e-05 3.0e-06 0.32
3 3.65e-05 3.8e-06 0.35
4 5.85e-05 5.4e-06 0.34
5 1.48e-05 3.3e-06 0.51
6 1.82e-05 3.0e-06 0.44
7 1.22e-05 3.8e-06 0.42
8 1.31e-05 3.6e-06 0.47
9 1.19e-05 3.9e-06 0.38
10 1.17e-05 3.9e-06 0.36
11 1.16e-05 4.0e-06 0.35
12 1.16e-05 4.0e-06 0.34
13 1.16e-05 4.0e-06 0.34
14 1.16e-05 4.0e-06 0.34
15 1.16e-05 4.0e-06 0.34
16 1.16e-05 4.0e-06 0.34
17 1.16e-05 4.0e-06 0.34
18 1.16e-05 4.0e-06 0.34
19 1.16e-05 4.0e-06 0.34
20 1.16e-05 4.0e-06 0.34
#σ = 0.1
1 4.23e-16 3.1e-16 0.81
2 4.23e-16 3.2e-16 0.85
3 4.24e-16 3.1e-16 0.87
4 4.24e-16 3.1e-16 0.88
5 4.25e-16 3.2e-16 0.82
6 4.25e-16 3.2e-16 0.82
7 4.24e-16 3.1e-16 0.87
8 4.24e-16 3.1e-16 0.87
9 4.24e-16 3.1e-16 0.87
10 4.25e-16 3.2e-16 0.82
11 4.25e-16 3.2e-16 0.82
12 4.24e-16 3.1e-16 0.87
13 4.25e-16 3.2e-16 0.82
14 4.24e-16 3.1e-16 0.87
15 4.25e-16 3.2e-16 0.82
16 4.24e-16 3.1e-16 0.87
17 4.25e-16 3.2e-16 0.82
18 4.24e-16 3.1e-16 0.87
19 4.25e-16 3.2e-16 0.82
20 4.24e-16 3.1e-16 0.87

View File

@ -1,92 +0,0 @@
#rounds mean stddev skew
#σ = 0
1 2.43e-4 1.3e-5 0.22
2 1.67e-4 1.3e-5 0.27
3 1.52e-4 1.3e-5 0.29
4 1.43e-4 1.3e-5 0.29
5 1.36e-4 1.2e-5 0.28
6 1.31e-4 1.2e-5 0.29
7 1.26e-4 1.2e-5 0.29
8 1.22e-4 1.2e-5 0.30
9 1.19e-4 1.2e-5 0.31
10 1.16e-4 1.2e-5 0.33
11 1.13e-4 1.2e-5 0.34
12 1.11e-4 1.2e-5 0.34
13 1.09e-4 1.2e-5 0.36
14 1.07e-4 1.1e-5 0.36
15 1.05e-4 1.1e-5 0.36
16 1.03e-4 1.1e-5 0.37
17 1.01e-4 1.1e-5 0.37
18 9.97e-5 1.1e-5 0.37
19 9.82e-5 1.1e-5 0.37
20 9.68e-5 1.1e-5 0.37
40 7.74e-5 1.0e-5 0.35
60 6.57e-5 9.4e-6 0.32
80 5.78e-5 9.1e-6 0.30
100 5.20e-5 9.0e-6 0.30
400 3.61e-5 1.0e-5 0.33
800 3.60e-5 1.0e-5 0.33
1600 3.60e-5 1.0e-5 0.33
#σ=0.005
1 2.92e-04 4.7e-05 1.43
2 2.18e-04 4.7e-05 1.72
3 2.28e-04 4.5e-05 1.73
4 2.12e-04 4.8e-05 1.72
5 2.08e-04 4.8e-05 1.71
6 2.02e-04 4.9e-05 1.69
7 2.05e-04 4.9e-05 1.70
8 2.00e-04 4.9e-05 1.69
9 1.98e-04 5.0e-05 1.68
10 1.95e-04 5.0e-05 1.67
11 1.97e-04 5.0e-05 1.68
12 1.94e-04 5.0e-05 1.67
13 1.93e-04 5.0e-05 1.67
14 1.92e-04 5.0e-05 1.67
15 1.91e-04 5.1e-05 1.66
16 1.90e-04 5.1e-05 1.66
17 1.90e-04 5.1e-05 1.66
18 1.89e-04 5.1e-05 1.66
19 1.88e-04 5.1e-05 1.66
20 1.88e-04 5.1e-05 1.66
#σ = 0.01
1 3.47e-04 1.0e-04 1.58
2 3.96e-04 1.1e-04 1.53
3 3.39e-04 1.1e-04 1.56
4 3.42e-04 1.1e-04 1.57
5 3.34e-04 1.1e-04 1.55
6 3.35e-04 1.1e-04 1.55
7 3.37e-04 1.1e-04 1.55
8 3.33e-04 1.1e-04 1.55
9 3.33e-04 1.1e-04 1.55
10 3.32e-04 1.1e-04 1.55
11 3.32e-04 1.1e-04 1.54
12 3.31e-04 1.1e-04 1.54
13 3.31e-04 1.1e-04 1.54
14 3.31e-04 1.1e-04 1.54
15 3.31e-04 1.1e-04 1.55
16 3.31e-04 1.1e-04 1.55
17 3.31e-04 1.1e-04 1.55
18 3.31e-04 1.1e-04 1.55
19 3.31e-04 1.1e-04 1.55
20 3.31e-04 1.1e-04 1.55
#σ=0.05
1 1.49e-03 5.7e-04 1.38
2 1.51e-03 5.6e-04 1.41
3 1.50e-03 5.7e-04 1.39
4 1.50e-03 5.6e-04 1.40
5 1.52e-03 5.6e-04 1.42
6 1.53e-03 5.6e-04 1.43
7 1.52e-03 5.6e-04 1.43
8 1.51e-03 5.6e-04 1.42
9 1.54e-03 5.5e-04 1.45
10 1.53e-03 5.5e-04 1.44
11 1.53e-03 5.5e-04 1.44
12 1.54e-03 5.5e-04 1.45
13 1.54e-03 5.5e-04 1.45
14 1.54e-03 5.5e-04 1.46
15 1.55e-03 5.5e-04 1.46
16 1.55e-03 5.5e-04 1.46
17 1.55e-03 5.5e-04 1.47
18 1.56e-03 5.5e-04 1.47
19 1.56e-03 5.4e-04 1.47
20 1.56e-03 5.4e-04 1.48

30
ex-6/plots/emd-round.py → ex-6/plots/emd.py Normal file → Executable file
View File

@ -1,5 +1,8 @@
#!/usr/bin/env python
import numpy as np
import matplotlib.pyplot as plt
import sys
def plot(table, title='', log=False):
@ -24,23 +27,18 @@ def plot(table, title='', log=False):
twin.set_ylabel('standard deviation', color='gray')
twin.ticklabel_format(style='sci', axis='y',
scilimits=(0, 0), useMathText=True)
# plt.subplot(212)
# if log:
# plt.xscale('log')
# plt.title('skewness', loc='right')
# plt.xlabel('RL rounds')
# plt.plot(table[0], table[3], color='xkcd:gray')
plt.tight_layout()
table = np.loadtxt('ex-6/plots/emd-round-noise.txt')
file = sys.argv[1] if len(sys.argv) > 1 else 'noiseless'
table = np.loadtxt('ex-6/plots/emd-' + file + '.txt')
if file == 'noiseless':
plot(table[:31].T, title=r'noise at $\sigma_N = 0.005$')
plot(table[31:51].T, title=r'noise at $\sigma_N = 0.005$')
plot(table[51:].T, title=r'noise at $\sigma_N = 0.01$')
else:
plot(table[:20].T, title=r'noise at $\sigma_N = 0.005$')
plot(table[20:40].T, title=r'noise at $\sigma_N = 0.005$')
plot(table[40:].T, title=r'noise at $\sigma_N = 0.01$')
plt.show()
# plot(table[:27].T, title='noiseless', log=True)
plot(table[27:47].T, title=r'noise at $\sigma_N = 0.005$')
plt.savefig('notes/images/6-rounds-noise-0.005.pdf')
plot(table[47:67].T, title=r'noise at $\sigma_N = 0.01$')
plt.savefig('notes/images/6-rounds-noise-0.01.pdf')
plot(table[67:].T, title=r'noise at $\sigma_N = 0.05$')
plt.savefig('notes/images/6-rounds-noise-0.05.pdf')

View File

@ -647,31 +647,6 @@ histogram as follow.
![Example of Noisy histogram,
$\sigma_N = 0.05$.](images/6-noisy.pdf){#fig:noisy}
<div id="fig:rounds-noise">
![](images/6-noise-rounds-0.005.pdf){#fig:rnoise-0.005}
![](images/6-noise-rounds-0.01.pdf){#fig:rnoise-0.01}
![](images/6-noise-rounds-0.05.pdf){#fig:rnoise-0.05}
EMD as a function of RL rounds for different noise $\sigma_N$ values with the
kernel $\sigma = 0.8 \Delta \theta$. The average is shown in red and the
standard deviation in grey. Noisy results.
</div>
<div id="fig:emd-noisy">
![$\sigma_N = 0.005$](images/6-noise-emd-0.005.pdf){#fig:enoise-0.005}
![$\sigma_N = 0.01$](images/6-noise-emd-0.01.pdf){#fig:enoise-0.01}
![$\sigma_N = 0.05$](images/6-noise-emd-0.05.pdf){#fig:enoise-0.05}
EMD distributions for different noise $\sigma_N$ values. The plots on the left
show the results for the FFT deconvolution, the central column the results for
the RL deconvolution and the third one shows the EMD for the convolved signal.
Noisy results.
</div>
For each bin, once the convolved histogram was computed, a value $v_N$ was
randomly sampled from a Gaussian distribution with standard deviation
$\sigma_N$, and the value $v_n \cdot b$ was added to the bin itself, where $b$
@ -724,4 +699,28 @@ noise. When the kernel is not known a-priori, either of them turns out to be as
good as the FFT in the aforementioned situation: only a poor approximation of
the original signal can be derived.
<div id="fig:rounds-noise">
![](images/6-noise-rounds-0.005.pdf){#fig:rnoise-0.005}
![](images/6-noise-rounds-0.01.pdf){#fig:rnoise-0.01}
![](images/6-noise-rounds-0.05.pdf){#fig:rnoise-0.05}
EMD as a function of RL rounds for different noise $\sigma_N$ values with the
kernel $\sigma = 0.8 \Delta \theta$. The average is shown in red and the
standard deviation in grey. Noisy results.
</div>
<div id="fig:emd-noisy">
![$\sigma_N = 0.005$](images/6-noise-emd-0.005.pdf){#fig:enoise-0.005}
![$\sigma_N = 0.01$](images/6-noise-emd-0.01.pdf){#fig:enoise-0.01}
![$\sigma_N = 0.05$](images/6-noise-emd-0.05.pdf){#fig:enoise-0.05}
EMD distributions for different noise $\sigma_N$ values. The plots on the left
show the results for the FFT deconvolution, the central column the results for
the RL deconvolution and the third one shows the EMD for the convolved signal.
Noisy results.
</div>