#-------------------------------------------------------
# example32.py
# pi0 -> gamma gamma decay calculations
# all values are in MeV
#-------------------------------------------------------
from ROOT import TLorentzVector
from numpy import sqrt

# input data
p0 = 500 # MeV
m  = 135 # MeV
E  = sqrt(p0*p0 + m*m) # MeV
betaz = p0 / E

print("pi0 momentum : ",p0)

# CM frame
E1, E2 = m/2, m/2
p1z = E1
p2z = -p1z
gamma1 = TLorentzVector()
gamma2 = TLorentzVector()
gamma1.SetXYZM(0,0,p1z,0)
gamma2.SetXYZM(0,0,p2z,0)
print("CM photon energies : ",gamma1.E(),gamma2.E())

# lab frame
# Boost() function modifies 4-vector components
gamma1.Boost(0,0,betaz)
gamma2.Boost(0,0,betaz)
print("Boosted photon energies : ",gamma1.E(),gamma2.E())