DV4 (https://physikmix.lewinb.net for more!)

1b.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Apr 30 13:17:31 2019

@author: Bormann, Lewin; Jürgens, Kira; Leung, Emily
"""

import random

rng = (155,160)

def random_sample(tau, n):
    return [random.expovariate(1/tau) for i in range(0, n)]

def count_in(sample, rng):
    return sum([1 if (rng[0] <= e and e <= rng[1]) else 0 for e in sample])

print("There are", count_in(random_sample(100, 1000), rng), "samples in {}".format(rng))

1c.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Apr 30 13:17:31 2019

@author: Bormann, Lewin; Jürgens, Kira; Leung, Emily
"""

import numpy as np
import matplotlib.pyplot as plt

import random

rng = (155,160)
super_samples = 100

def random_sample(tau, n):
    return [random.expovariate(1/tau) for i in range(0, n)]

def count_in(sample, rng):
    return sum([1 if (rng[0] <= e and e <= rng[1]) else 0 for e in sample])

sample_of_samples = np.array([count_in(random_sample(100, 1000), rng) for i in range(0, super_samples)], dtype=np.int)

fig = plt.figure()
ax = fig.add_subplot(111)
ax.hist(sample_of_samples, bins=range(sample_of_samples.min(), sample_of_samples.max()))

plt.show()

1d.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Apr 30 13:17:31 2019

@author: Bormann, Lewin; Jürgens, Kira; Leung, Emily
"""

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats

from math import exp
import random

tau = 100
rng = (155,160)
samples = 1000
super_samples = 100

def random_sample(tau, n):
    return [random.expovariate(1/tau) for i in range(0, n)]

def count_in(sample, rng):
    return sum([1 if (rng[0] <= e and e <= rng[1]) else 0 for e in sample])

def probability_in_range(tau, rng):
    return - exp(-1/tau * rng[1]) + exp(-1/tau * rng[0])

sample_of_samples = np.array([count_in(random_sample(tau, samples), rng) for i in range(0, super_samples)], dtype=np.int)

fig = plt.figure()
ax = fig.add_subplot(111)

super_sample_range = range(sample_of_samples.min(), sample_of_samples.max())
ax.hist(sample_of_samples, bins=super_sample_range)

# Theoretisch erwartete Verteilung

def probability_for_r_events(r, n, p):
    return stats.binom.pmf(r, n, p)

predicted = stats.binom.pmf(super_sample_range, samples, probability_in_range(tau, rng)) * super_samples
xaxis = np.array(super_sample_range) + 1/2
ax.plot(xaxis, predicted, 'ro-')

plt.show()

1e.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Apr 30 13:17:31 2019

@author: Bormann, Lewin; Jürgens, Kira; Leung, Emily
"""

import numpy as np
import matplotlib.pyplot as plt

import math
import random

rng = (155,160)
super_samples = 100

def random_sample(tau, n):
    return [random.expovariate(1/tau) for i in range(0, n)]

def count_in(sample, rng):
    return sum([1 if (rng[0] <= e and e <= rng[1]) else 0 for e in sample])

sample_of_samples = np.array([count_in(random_sample(100, 1000), rng) for i in range(0, super_samples)], dtype=np.int)

def sample_variance(sample):
    avg = sample.sum()/sample.size
    return sum([(x - avg)**2 for x in sample])/sample.size

def stat_error(sample):
    return math.sqrt(1/((sample.size-1)) * sample_variance(sample))

print("Average of sample:", sample_of_samples.mean())
print("Std deviation of sample:", math.sqrt(sample_variance(sample_of_samples)))
print("Statistical error of sample:", stat_error(sample_of_samples))

fig = plt.figure()
ax = fig.add_subplot(111)
ax.hist(sample_of_samples, bins=range(sample_of_samples.min(), sample_of_samples.max()))

plt.show()

3abcd.py

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Apr 30 20:46:26 2019

@author: Bormann, Lewin; Jürgens, Kira; Leung, Emily
"""

import matplotlib.pyplot as plt
import numpy as np

import math

data = np.loadtxt("DateiBlatt4.csv")
bins_per_int = 3


# 3a)
fig = plt.figure(dpi=150)
ax = fig.add_subplot(111)
ax.grid(True)

bins = np.linspace(data.min(), data.max(), int((data.max()-data.min())*bins_per_int))
ax.hist(data, bins=bins)

# Wir berechnen die Werte einmal händisch, dann mit numpy-Methoden.
# 3b)

mean = data.sum() / data.size
print("Average:", data.mean(), mean)
err = math.sqrt(1/(data.size-1) * 1/data.size * ((data - mean)**2).sum())
print("Statistical error:", math.sqrt(1/(data.size) * data.var(ddof=1)), err)

# 3c)

stddev = math.sqrt(1/(data.size-1) * ((data -  mean)**2).sum())
print("Standard deviation:", data.std(ddof=1), stddev)

# Wähle Parameter σ = stddev, µ = mean

def gauss_fn(σ, µ):
    def gauss(x):
        return 1/(math.sqrt(2*math.pi)*σ) * math.exp(-(x - µ)**2/(2*σ**2))
    return np.vectorize(gauss)

def gauss_prob(min, max, σ, µ):
    """Calculates the probability of events between min and max.

    min and max can be None to express +/- Inf."""
    upper = 1 if max is None else math.erf((max - µ)/(math.sqrt(2)*σ))
    lower = 0 if min is None else math.erf((min - µ)/(math.sqrt(2)*σ))
    return 1/2 * (upper - lower)

# Berechne Integral für jede Balkenbreite, um die Zahlen direkt zu vergleichen.
x = np.zeros(bins.size-1)
y = np.zeros(x.size)
for (i, min) in enumerate(bins[:-1]):
    mid = (bins[i]+bins[i+1])/2
    x[i] = mid
    delta = (bins[i+1]-bins[i])/2
    y[i] = gauss_prob(mid-delta, mid+delta, stddev, mean) * data.size

ax.plot(x, y, 'rx-')
ax.plot(x, gauss_fn(stddev, mean)(x) * data.size/3, 'go-')

# 3d)

print("P(X > 3.5) =", gauss_prob(3.5, None, stddev, mean))