TPC2-Muon-Teresa

# Python 3

from random import random
from math   import log, sqrt

# Constantes físicas
m    = 0.106		# massa do muon
tau0 = 2.2e-6		
c    = 3.0e8		# velocidade da luz

# Funções

def rnd_exponencial(a):
    """
    Devolve um número aleatório n distribuído de acordo com exp(-a*n).
    """
    return -1/a*log(random())
    

def prob_deteccao(L, E0, num_iterecoes):
	"""
	Devolve a probabilidade de que um muon de Energia E0 produzido à altura L
	seja detectado na terra.
	"""

	# Contador da quantidade de muões que sobrevivem
	contador=0

	# Ciclo que repete o algoritmo num_iterecoes vezes
	for i in range(num_iterecoes):

		# Gerar aleatoriamente uma energia
		Ek = rnd_exponencial(1/E0)

		# Cálculo de gamma e tau
		gamma = (Ek+m)/m
		tau   = gamma*tau0

		# Geramos um tempo de vida aleatório
		t = rnd_exponencial(1/tau) 

		# Calculamos a distância que este muon pode percorrer
		beta = sqrt( abs(Ek**2 - m**2) ) /Ek
		l    = beta*c*t

		# Determinamos se este muon sobreviveu
		if l>L:
			contador+=1


	# Calculamos a probabilidade
	p = contador/num_iterecoes

	return p


# Código principal

num_iterecoes = 1000

parametros = [ [5000, 2.5], [4000, 3.0], [3000, 1.5], [9000, 2.0], [7500, 1.0] ]

for L, E0 in parametros:
    p=prob_deteccao(L, E0, num_iterecoes)
    print("L=%d E0=%.1f --> p=%0.3f" % (L, E0, p) )