#include<iostream>
#include<cmath>

using namespace std;

class LorentzVector {
  public:
    //Constructor:
    LorentzVector(double t, double x, double y, double z) 
        {_t = t; _x = x; _y = y; _z = z;}

    //Access functions:
    double GetT() {return _t;}
    double GetX() {return _x;}
    double GetY() {return _y;}
    double GetZ() {return _z;}

    //Methods:
    void Print();
    LorentzVector operator+(LorentzVector vec);
    double Contract();
    LorentzVector LorentzTransformation(double beta);

  private:
    //Encapsulation of space-time coordinates:
    double _t, _x, _y, _z;
};

int main() {
  //Make a couple of Lorentz vectors and print them
  LorentzVector A(-5,6,-1,2);
  A.Print();
  LorentzVector B(4,-6,4,-5);
  B.Print();

  //Test the overloaded + operator
  LorentzVector C = A + B;
  C.Print();

  //Test the Lorentz boost method
  cout << "Contraction before boost = " << C.Contract() << endl;
  LorentzVector D = C.LorentzTransformation(0.7);
  D.Print();
  cout << "Contraction after boost = " << D.Contract() << endl;

  return 0;
}

//Print space-time coordinates
void LorentzVector::Print() {
  cout << "Lorentz vector is (";
  cout << "t=" << _t;
  cout << ", x=" << _x;
  cout << ", y=" << _y;
  cout << ", z=" << _z;
  cout << ")" << endl;
}

//Overload + operator
LorentzVector LorentzVector::operator+(LorentzVector vec) {
  return LorentzVector(vec.GetT()+_t, vec.GetX()+_x, vec.GetY()+_y, vec.GetZ()+_z);
}

//Contract Lorentz 4-vector
double LorentzVector::Contract() {
  return _t*_t  - _x*_x - _y*_y - _z*_z;
}
    
//Lorentz boost - return a boosted Lorentz vector
//(Alternative: boost the class itself)
LorentzVector LorentzVector::LorentzTransformation(double beta) {
  double gamma = 1. / sqrt(1. - beta*beta);
  double T = gamma * (_t - beta*_x);
  double X = gamma * (_x - beta*_t);
  return LorentzVector(T, X, _y, _z);
}