// * AUTHOR:      Piotr Golonka 
// * LAST UPDATE: 2000-01-17
// * COPYRIGHT:   (C) Faculty of Nuclear Physics & Techniques, UMM Cracow.

////////////////////////////////////////////////////////////////////////////
//
//  HEPParticle : abstract class for HEP particles.
//
////////////////////////////////////////////////////////////////////////////

#include "HEPParticle.H"
#include "HEPEvent.H"
#include "MC3Vector.H"
#include "MC4Vector.H"
#include <math.h>
#include <stdlib.h>
#include <list>
#include <iostream>
#include <assert.h>

using namespace std;

//#ifdef _USE_ROOT_
ClassImp(HEPParticle)
ClassImp(HEPParticleList)
ClassImp(HEPParticleListIterator)
//#endif

HEPParticleList::HEPParticleList()//:
//first(HEPParticleListIterator(0)),
//last(HEPParticleListIterator(0))
{
    data = new std::list<HEPParticle*>;

}

HEPParticleList::~HEPParticleList()
{
    if (data) {
	std::list<HEPParticle*> *plist = (std::list<HEPParticle*>*)data;
	delete plist;
	data=0;
    }
}

long int HEPParticleList::itr_begin()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    long int retval;
    std::list<HEPParticle*>::iterator &itr= *(std::list<HEPParticle*>::iterator*)(&retval);    
    itr=plist.begin();
    return retval;
}

long int HEPParticleList::itr_end()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    long int retval;
    std::list<HEPParticle*>::iterator &itr=*(std::list<HEPParticle*>::iterator*)(&retval);
    itr=plist.end();
    return retval;
}

int HEPParticleList::empty()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    return plist.empty();
}

int HEPParticleList::size()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    return plist.size();
}

void HEPParticleList::clear()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    plist.clear();
    firstParticle=0;
    lastParticle=0;
}

void HEPParticleList::pop_front()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;

    plist.pop_front();

    if (plist.empty()) {
	firstParticle=0;
	lastParticle=0;
    } else {
	firstParticle = (plist.front());
    }
}

void HEPParticleList::push_back(HEPParticle *p)
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;

    if(plist.empty()) firstParticle=p;

    plist.push_back(p);
    lastParticle=p;
}

void HEPParticleList::ls(char *option)
{
    HEPParticleListIterator itr(*this);
    for (HEPParticle *p = itr.first(); p!=0; p=itr.next() ) {
	p->ls();
    }
}

void HEPParticleList::remove(HEPParticle* p)
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    plist.remove(p);
}

void HEPParticleList::unique()
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    plist.unique();
}

void HEPParticleList::insert(HEPParticleListIterator *itr, int n, HEPParticle *p)
{
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) data;
    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(itr->ptr);
    plist.insert(itrptr,n,p);
}


bool HEPParticleList::contains(int id)
{
  HEPParticleListIterator itr(*this);
  for (HEPParticle *p = itr.first(); p!=0; p=itr.next() ){
    if(p->GetId()==id)
      return true;
  }
  return false;
}

HEPParticleListIterator::HEPParticleListIterator(HEPParticleList &l):
list(&l)
{
//    ptr=malloc(sizeof (std::list<HEPParticle*>::iterator));
    ptr=new (std::list<HEPParticle*>::iterator);

    first();
}

HEPParticleListIterator::~HEPParticleListIterator()
{
    if (ptr) {
	std::list<HEPParticle*>::iterator *itrptr = (std::list<HEPParticle*>::iterator *)(ptr);
	delete itrptr;
	itrptr=0;
    }
}

HEPParticle* HEPParticleListIterator::me()
{
    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(ptr);
    particle = *itrptr;
    return particle;
}

HEPParticle* HEPParticleListIterator::first()
{
    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(ptr);
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) list->data;

    if (plist.empty()) return 0;

    itrptr= plist.begin();

    particle = *itrptr;

    return particle;
}

HEPParticle* HEPParticleListIterator::last()
{
    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(ptr);
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) list->data;

    if (plist.empty()) return 0;

    itrptr = plist.end();
    itrptr--;

    particle = *itrptr;
    
    return particle;
}

HEPParticle* HEPParticleListIterator::next()
{

    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(ptr);

    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) list->data;

    std::list<HEPParticle*>::iterator last_particle = --(plist.end());

    if ( itrptr == last_particle ) {
      //	cout<<"this exceeds upper limit!"<<endl;
	return 0;
    }
    itrptr++;
    
    particle = *itrptr;
    
    return particle;
}

HEPParticle* HEPParticleListIterator::prev()
{
    std::list<HEPParticle*>::iterator &itrptr = *(std::list<HEPParticle*>::iterator *)(ptr);
    std::list<HEPParticle*> &plist = *(std::list<HEPParticle*>*) list->data;
    
    if ( itrptr == (plist.begin()) ) {
	
//	cout<<"you try to get before begin"<<endl;
	return 0;
    }
    
    itrptr--;
    
    particle = *itrptr;
    
    return particle;
}


HEPParticleList* HEPParticle::GetDaughterList(HEPParticleList *list)
{
    // generic, non-optimized version based on indices.
    //
    // if list is not provided, it is created.
    //
    if (!list) list=new HEPParticleList();
    
    int d1=GetFirstDaughter();
    int d2=GetLastDaughter();

    // only if properly ordered! ie d
    if ( (d1 <= d2) && (d1!=0) &&(d2!=0) ) {
    
	for (int i=d1; i<=d2; i++) {
    
	    HEPParticle *daughter=GetEvent()->GetParticle(i);
	    assert(daughter!=0);
	    // sanity check...
	    if (daughter->GetMother() == GetId())
		list->push_back(daughter);
	}
    } else {
	//printf("WARNING: HEPParticle::GetListOfDaughters - scanning for dauhgters\n");
	//ls();

	// find all daughter's that point to me...
	for (int i=1;i<=GetEvent()->GetNumOfParticles();i++) {
	    HEPParticle *particle=GetEvent()->GetParticle(i);
	    if (particle->GetMother() == GetId()) {
		list->push_back(particle);
		//printf("   ->");particle->ls();
	    }
	}
	//printf("-----\n");
    }
    
    return list;
}

HEPParticleList* HEPParticle::GetMotherList(HEPParticleList *list)
{
    // generic, non-optimized version based on indices.
    //
    // if list is not provided, it is created.
    //
    if (!list) list=new HEPParticleList();
    
    int m1=GetMother();
    int m2=GetMother2();
    
    if (m2<m1) m2=m1;
    
    // at least m1 needs to be valid...
    if (m1!=0) {
    
	for (int i=m1; i<=m2; i++) {
    
	    HEPParticle *mother=GetEvent()->GetParticle(i);
	    assert (mother!=0);
	    list->push_back(mother);
	}
    }
    return list;
}


  // we want to compare doubles with arbitrary precission.
  // this is done by EQUAL macro defined here.
  // In case we compare non-zero numbers - EQUAL compares numbers
  // with RELATIVE precission given by DEFAULT_PRECISION definition.
  // i.e. EQUAL(1e-20,1e-22) will return false, because relatively
  // difference is two significant places.
  // In case when we compare a number with zero - comparison is
  // absolute: i.e: EQUAL(1e-3.0) returns false, but
  // EQUAL (1e-15,0) returns true.

#define DEFAULT_PRECISION 1e-6

#define EQL(n1,n2)   ( ( fabs(n1-n2)      < DEFAULT_PRECISION) ? 1:0)
#define EQREL(n1,n2) ( ( fabs((n1-n2)/n1) < DEFAULT_PRECISION) ? 1:0)

#define EQUAL(n1,n2) ( (n1==0) ? EQL(n1,n2) : EQREL(n1,n2) )



HEPParticle::HEPParticle()
{
  // default constructor - dummy.
}


bool HEPParticle::Compare_WithoutId( HEPParticle &p)
{
  // Compares another HEPParticle object to this, taking into
  // account their identity (PDG codes), status, momenta
  // and vertices's coordinates.
  //
  // This function is useful for comparing particles regardless
  // of their poition in event record.
  //
  // i.e.
  // HEPParticle photon1 = MyAnalysis.GetPhotons(1);
  // HEPParticle photon2 = MyAnalysis.GetPhotons(2);
  //
  // if (photon1.Compare_WithoutId(photon2) )
  //     printf("Photons are the same\n");
  //
  // See also:  
  //           HEPParticle::operator==  
  //           HEPParticle::operator!=  

  bool result =  ( (GetPDGId()     == p.GetPDGId() ) &&
		   (GetStatus()    == p.GetStatus()) &&
		   EQUAL(GetPx()   ,  p.GetPx()    ) &&
		   EQUAL(GetPy()   ,  p.GetPy()    ) &&
		   EQUAL(GetPz()   ,  p.GetPz()    ) &&
		   EQUAL(GetE()    ,  p.GetE()     ) &&
		   EQUAL(GetM()    ,  p.GetM()     ) &&
		   EQUAL(GetVx()   ,  p.GetVx()    ) &&
		   EQUAL(GetVy()   ,  p.GetVy()    ) &&
		   EQUAL(GetVz()   ,  p.GetVz()    ) &&
		   EQUAL(GetTau()  ,  p.GetTau()   ) );
return result;
}


double HEPParticle::Spin(int id)
{
  // private function that returns the spin of some particles.
  // It is used internaly by HEPParticle::GetSpin
  
  if ( (id >=1 && id <=8) || (id>=11 && id<=18) ) // quarks:
    return 0.5;
  if ( (id >=25 && id<=100 ) || id==10 )
    return 0.0;

  switch(id) {
  case 9:
  case 21:
  case 22:
  case 23:
  case 24:
    return 1.0;     //bosons
    break;

  }
  return 0.0;
}


double HEPParticle::GetSpin(int pdgId)
{
  // returns a spin of a particle with given PDG code.
  //
  // This function is declared public,static, that means it may be
  // called without an instance of object in any place of a code.
  // i.e.
  //
  // int PDGcode=2212
  // double spin = HEPParticle::GetSpin(PDGcode);
  // printf("particle identified by code %i has spi %f \n",PDGcode,spin);
  //
  double spin=0.0;

  int pdg = abs(pdgId);

  if (pdg<=100)
    spin = Spin(pdg);
  else {
    // other particles...
    spin =  ( (pdg % 10) - 1.0)/2.0;

    // ...K_short and K_long are special !!
    if(spin < 0)
      spin = 0.0;      
  }

    return spin;
  
}

double HEPParticle::Charge(int id)
{  
  // private function that returns the charge of one of the
  // first 100 particles... It is used internaly by HEPParticle::GetCharge
  
 
  if ( (id>=18&& id<=23) || (id>=25 && id<37) || (id>=38 && id<=100) )
       return 0.0;
  
  switch (id) {
  case 0:
    return 0.0;
    break;
  case 1:
  case 3:
  case 5:
  case 7:
    //    return -1.0/3.0;
    return -0.3333333333;
    break;
  case 2:
  case 4:
  case 6:
  case 8:
    //    return 2.0/3.0;
    return  0.6666666667;
    break;
  case 9:
  case 10:
  case 12:
  case 14:
  case 16:
    return 0.0;
    break;
  case 11:
  case 13:
  case 15:
  case 17:
    return -1.0;
    break;
  case 24:
  case 37:
    return 1.0;
    break;
    
  } 
  printf("ERROR in HEPParticle::Charge(%i): pdg not found!\n",id); 
  return 0.0;
}


double HEPParticle::GetCharge(int pdgId)
{
  // returns an electric charge of a particle with given PDG code.
  //
  // This function is declared public,static, that means it may be
  // called without an instance of object in any place of a code.
  // i.e.
  //
  // int PDGcode=2212
  // double charge = HEPParticle::GetCharge(PDGcode);
  // printf("particle identified by code %i has charge %f \n",PDGcode,charge);
  //
 
  
  double charge=0;

  int pdg=abs(pdgId);

  if (pdg<=100) {
    charge=Charge(pdg);
  } else {
    // check particles made of quarks.
    // unpack their code:
    int quark3=(pdg/1000)%10;
    int quark2=(pdg/100)%10;
    int quark1=(pdg/10)%10;
    

    if (quark3==0) {
      //MESON
      if ( (quark2 % 2) == 0) 
	charge=Charge(quark2) - Charge(quark1);
      else
	charge=Charge(quark1) - Charge(quark2);
    } else {
      //DIQUARK or Baryon;
      charge=Charge(quark1) + Charge(quark2) + Charge(quark3);
    }

  }
  
  // determine sign...
  if (pdgId<0)
    charge=-charge;

  // check it once again
  if (charge*charge < 1e-6)
    charge=0.0;

  return charge;

}

char* HEPParticle::GetParticleName(int PDGCode)
{
  // returns a string containing a name of a particle
  // with given PDG code. If given PDGCode does not corespond
  // to any named particle, returned string contains a textual
  // representaion of PDGCode (i.e. GetParticleName(99999) will
  // return string "99999"
  //
  // This function is declared public,static, that means it may be
  // called without an instance of object in any place of a code.
  // i.e.
  //
  // int PDGcode=2212
  // char *name = HEPParticle::GetParticleName(PDGcode);
  // printf("particle identified by code %i is called %s \n",PDGcode,name);
  //
 
  
  static char txtbuf[21]="                    ";

  switch(PDGCode){
  case 1: 
    return ("d");
    break;
  case -1:
    return("d~");
    break;
  case 2: 
    return ("u");
    break;
  case -2:
    return("u~");
    break;
  case 3: 
    return ("s");
    break;
  case -3:
    return("s~");
    break;
  case 4: 
    return ("c");
    break;
  case -4:
    return("c~");
    break;
  case 5: 
    return ("b");
    break;
  case -5:
    return("b~");
    break;
  case 6: 
    return ("t");
    break;
  case -6:
    return("t~");
    break;
  case 7: 
    return ("l");
    break;
  case -7:
    return("l~");
    break;
  case 8: 
    return ("h");
    break;
  case -8:
    return("h~");
    break;

  case 11: 
    return ("e-");
    break;
  case -11:
    return("e+");
    break;
  case 12: 
    return ("nu_e");
    break;
  case -12:
    return("nu_e~");
    break;
  case 13: 
    return ("mu-");
    break;
  case -13:
    return("mu+");
    break;
  case 14: 
    return ("nu_mu");
    break;
  case -14:
    return("nu_mu~");
    break;
  case 15: 
    return ("tau-");
    break;
  case -15:
    return("tau+");
    break;
  case 16: 
    return ("nu_tau");
    break;
  case -16:
    return("nu_tau~");
    break;
  case 17: 
    return ("chi-");
    break;
  case -17:
    return("chi+");
    break;
  case 18: 
    return ("nu_chi");
    break;
  case -18:
    return("nu_chi~");
    break;

  case 9:
  case 21:
    return ("g");
    break;
  case 22:
    return("gamma");
    break;
  case 23: 
    return ("Z0");
    break;
  case 24:
    return("W+");
    break;
  case -24:
    return("W-");
    break;
  case 25: 
    return ("H0_1");
    break;
  case 35:
    return("H0_2");
    break;
  case 36: 
    return ("H0_3");
    break;
  case 37:
    return("H+");
    break;
  case -37:
    return("H~-");
    break;
  case 28:
    return("reggeon");
    break;
  case 29:
    return("pomeron");
    break;
  case 39:
    return("LQ");
    break;
  case 40:
    return("RQ");
    break;
    

  case 90:
    return("diquark");
    break;
  case 91:
    return("cluster");
    break;
   case 92:
    return("string");
    break;
   case 93:
    return("indep.");
    break;

  case 1103:
    return("dd_1");
    break;
  case 2101:
    return("ud_0");
    break;
  case 2103:
    return("ud_1");
    break;
  case 2203:
    return("uu_1");
    break;
  case 3101:
    return("sd_0");
    break;
  case 3103:
    return("sd_1");
    break;
  case 3201:
    return("su_0");
    break;
  case 3203:
    return("su_1");
    break;
  case 3303:
    return("ss_1");
    break;



  case -1103:
    return("dd_1~");
    break;
  case -2101:
    return("ud_0~");
    break;
  case -2103:
    return("ud_1~");
    break;
  case -2203:
    return("uu_1~");
    break;
  case -3101:
    return("sd_0~");
    break;
  case -3103:
    return("sd_1~");
    break;
  case -3201:
    return("su_0~");
    break;
  case -3203:
    return("su_1~");
    break;
  case -3303:
    return("ss_1~");
    break;

  case 211:
    return("pi+");
    break;
  case -211:
    return("pi-");
    break;
  case 111:
    return("pi0");
    break;
  case 311:
    return("K0");
    break;
  case -311:
    return("K~0");
    break;
  case 321:
    return("K+");
    break;
  case -321:
    return("K-");
    break;
  case 411:
    return("D+");
    break;
  case -411:
    return("D-");
    break;
  case 421:
    return("D0");
    break;
  case -421:
    return("D~0");
    break;
  case 431:
    return("D_s+");
    break;
  case -431:
    return("D_s-");
    break;
  case 511:
    return("B0");
    break;
  case -511:
    return("B~0");
    break;
  case 521:
    return("B+");
    break;
  case -521:
    return("B-");
    break;
  case 531:
    return("B_s0");
    break;
  case -531:
    return("B_s~0");
    break;
  case 541:
    return("B_c+");
    break;
  case -541:
    return("B_c-");
    break;
  case 221:
    return("eta");
    break;
  case 331:
    return("eta'");
    break;
  case 441:
    return("eta_c");
    break;
  case 551:
    return("eta_b");
    break;
  case 661:
    return("eta_t");
    break;
  case 130:
    return("K_L0");
    break;
  case 310:
    return("K_S0");
    break;
  case -130:
    return("K_L~0");
    break;
  case -310:
    return("K_S~0");
    break;
  case 213:
    return("rho+");
    break;
  case -213:
    return("rho-");
    break;
  case 313:
    return("K*0");
    break;
  case -313:
    return("K*~0");
    break;
  case 323:
    return("K*+");
    break;
  case -323:
    return("K*-");
    break;
  case 413:
    return("D*+");
    break;
  case -413:
    return("D*-");
    break;
  case 423:
    return("D*0");
    break;
  case -423:
    return("D*~0");
    break;
  case 433:
    return("D*_s+");
    break;
  case -433:
    return("D*_s-");
    break;
  case 513:
    return("B*0");
    break;
  case -513:
    return("B*~0");
    break;
  case 523:
    return("B*+");
    break;
  case -523:
    return("B*-");
    break;
  case 533:
    return("B*_s0");
    break;
  case -533:
    return("B*_s~0");
    break;
  case 543:
    return("B*_c+");
    break;
  case -543:
    return("B*_c-");
    break;
  case 113:
    return("rho0");
    break;
  case 223:
    return("omega");
    break;
  case 333:
    return("phi");
    break;
  case 443:
    return("Jpsi");
    break;
  case 553:
    return("Upsilon");
    break;
  case 663:
    return("Theta");
    break;

  case 10213:
    return("b_1+");
    break;
  case -10213:
    return("b_1-");
    break;
  case 10313:
    return("K_10");
    break;
  case -10313:
    return("K_1~0");
    break;
  case 10323:
    return("K_1+");
    break;
  case -10323:
    return("K_1-");
    break;
  case 10413:
    return("D_1+");
    break;
  case -10413:
    return("D_1-");
    break;
  case 10423:
    return("D_10");
    break;
  case -10423:
    return("D_1~0");
    break;
  case 10433:
    return("D_1s+");
    break;
  case -10433:
    return("D_1s-");
    break;
  case 10113:
    return("b_10");
    break;
  case 10223:
    return("h_10");
    break;
  case 10333:
    return("h'_10");
    break;
  case 10443:
    return("h_1c0");
    break;
  case 10211:
    return("a_0+");
    break;
  case -10211:
    return("a_0-");
    break;
  case 10311:
    return("K*_00");
    break;
  case -10311:
    return("K*_0~0");
    break;
  case 10321:
    return("K*_0+");
    break;
  case -10321:
    return("K*_0-");
    break;
  case 10411:
    return("D*_0+");
    break;
  case -10411:
    return("D*_0-");
    break;
  case 10421:
    return("D*_00");
    break;
  case -10421:
    return("D*_0~0");
    break;
  case 10431:
    return("D*_0s+");
    break;
  case -10431:
    return("D*_0s-");
    break;
  case 10111:
    return("a_00");
    break;
  case 10221:
    return("f_00");
    break;
  case 10331:
    return("f'_00");
    break;
  case 10441:
    return("chi_0c0");
    break;
  case 20213:
    return("a_1+");
    break;
  case -20213:
    return("a_1-");
    break;
  case 20313:
    return("K*_10");
    break;
  case -20313:
    return("K*_1~0");
    break;
  case 20323:
    return("K*_1+");
    break;
  case -20323:
    return("K*_1-");
    break;
  case 20413:
    return("D*_1+");
    break;
  case -20413:
    return("D*_1-");
    break;
  case 20423:
    return("D*_10");
    break;
  case -20423:
    return("D*_1~0");
    break;
  case 20433:
    return("D*_1s+");
    break;
  case -20433:
    return("D*_1s-");
    break;
  case 20113:
    return("a_10");
    break;
  case 20223:
    return("f_10");
    break;
  case 20333:
    return("f'_10");
    break;
  case 20443:
    return("chi_1c0");
    break;
  case 215:
    return("a_2+");
    break;
  case -215:
    return("a_2-");
    break;
  case 315:
    return("K*_20");
    break;
  case -315:
    return("K*_2~0");
    break;
  case 325:
    return("K*_2+");
    break;
  case -325:
    return("K*_2-");
    break;
  case 415:
    return("D*_2+");
    break;
  case -415:
    return("D*_2-");
    break;
  case 425:
    return("D*_20");
    break;
  case -425:
    return("D*_2~0");
    break;
  case 435:
    return("D*_2s+");
    break;
  case -435:
    return("D*_2s-");
    break;
  case 115:
    return("a_20");
    break;
  case 225:
    return("f_20");
    break;
  case 335:
    return("f'_20");
    break;
  case 445:
    return("chi_2c0");
    break;
  case 30443:
    return("psi'");
    break;
  case 30553:
    return("Upsilon'");
    break;

  case 2112:
    return("n0");
    break;
  case 2212:
    return("p+");
    break;
  case 3112:
    return("Sigma-");
    break;
  case 3122:
    return("Lambda0");
    break;
  case 3212:
    return("Sigma0");
    break;
  case 3222:
    return("Sigma+");
    break;
  case 3312:
    return("Xi-");
    break;
  case 3322:
    return("Xi0");
    break;
  case 4112:
    return("Sigma_c0");
    break;
  case 4122:
    return("Lambda_c+");
    break;
  case 4212:
    return("Sigma_c+");
    break;
  case 4222:
    return("Sigma_c++");
    break;
  case 4132:
    return("Xi_c0");
    break;
  case 4312:
    return("Xi'_c0");
    break;
  case 4232:
    return("Xi_c+");
    break;
  case 4322:
    return("Xi'_c+");
    break;
  case 4332:
    return("Omega_c0");
    break;
  case 5112:
    return("Sigma_b-");
    break;
  case 5222:
    return("Sigma_b+");
    break;
  case 1114:
    return("Delta-");
    break;
  case 2114:
    return("Delta0");
    break;
  case 2214:
    return("Delta+");
    break;
  case 2224:
    return("Delta++");
    break;
  case 3114:
    return("Sigma*-");
    break;
  case 3224:
    return("Sigma*+");
    break;
  case 3314:
    return("Xi*-");
    break;
  case 3324:
    return("Xi*0");
    break;
  case 3334:
    return("Omega-");
    break;
  case 4114:
    return("Sigma*_c0");
    break;
  case 4214:
    return("Sigma*_c+");
    break;
  case 4224:
    return("Sigma*_c++");
    break;
  case 4314:
    return("Xi*_c0");
    break;
  case 4324:
    return("Xi*_c+");
    break;
  case 4334:
    return("Omega*_c0");
    break;
  case 5114:
    return("Sigma*_b-");
    break;
  case 5214:
    return("Sigma*_b0");
    break;
  case 5224:
    return("Sigma*_b+");
    break;




  case -2112:
    return("n~0");
    break;
  case -2212:
    return("p~-");
    break;
  case -3112:
    return("Sigma~+");
    break;
  case -3122:
    return("Lambda~0");
    break;
  case -3212:
    return("Sigma~0");
    break;
  case -3222:
    return("Sigma~-");
    break;
  case -3312:
    return("Xi~+");
    break;
  case -3322:
    return("Xi~0");
    break;
  case -4112:
    return("Sigma_c~0");
    break;
  case -4122:
    return("Lambda_c~-");
    break;
  case -4212:
    return("Sigma_c~-");
    break;
  case -4222:
    return("Sigma_c~--");
    break;
  case -4132:
    return("Xi_c~0");
    break;
  case -4312:
    return("Xi'_c~0");
    break;
  case -4232:
    return("Xi_c~-");
    break;
  case -4322:
    return("Xi'_c~-");
    break;
  case -4332:
    return("Omega_c~0");
    break;
  case -5112:
    return("Sigma_b~+");
    break;
  case -5222:
    return("Sigma_b~-");
    break;
  case -1114:
    return("Delta~+");
    break;
  case -2114:
    return("Delta~0");
    break;
  case -2214:
    return("Delta~-");
    break;
  case -2224:
    return("Delta~--");
    break;
  case -3114:
    return("Sigma*~+");
    break;
  case -3224:
    return("Sigma*~-");
    break;
  case -3314:
    return("Xi*~+");
    break;
  case -3324:
    return("Xi*~0");
    break;
  case -3334:
    return("Omega~+");
    break;
  case -4114:
    return("Sigma*_c~0");
    break;
  case -4214:
    return("Sigma*_c~-");
    break;
  case -4224:
    return("Sigma*_c~--");
    break;
  case -4314:
    return("Xi*_c~0");
    break;
  case -4324:
    return("Xi*_c~-");
    break;
  case -4334:
    return("Omega*_c~0");
    break;
  case -5114:
    return("Sigma*_b~+");
    break;
  case -5214:
    return("Sigma*_b~0");
    break;
  case -5224:
    return("Sigma*_b~-");
    break;



  default:
    sprintf(txtbuf,"%i",PDGCode);
    return (txtbuf);
    break;
  }


}



char* HEPParticle::GetLaTeXName(int PDGCode)
{
  
  static char txtbuf[21]="                    ";

  switch(PDGCode){
  case 1: 
    return ("d");
    break;
  case -1:
    return("\\widetilde{d}");
    break;
  case 2: 
    return ("u");
    break;
  case -2:
    return("\\widetilde{u}");
    break;
  case 3: 
    return ("s");
    break;
  case -3:
    return("\\widetilde{s}");
    break;
  case 4: 
    return ("c");
    break;
  case -4:
    return("\\widetilde{c}");
    break;
  case 5: 
    return ("b");
    break;
  case -5:
    return("\\widetilde{b}");
    break;
  case 6: 
    return ("t");
    break;
  case -6:
    return("\\widetilde{t}");
    break;
  case 7: 
    return ("l");
    break;
  case -7:
    return("\\widetilde{l}");
    break;
  case 8: 
    return ("h");
    break;
  case -8:
    return("\\widetilde{h}");
    break;

  case 11: 
    return ("e^{-}");
    break;
  case -11:
    return("e^{+}");
    break;
  case 12: 
    return ("\\nu_{e}");
    break;
  case -12:
    return("\\widetilde{\\nu_{e}}");
    break;
  case 13: 
    return ("\\mu^{-}");
    break;
  case -13:
    return("\\mu^{+}");
    break;
  case 14: 
    return ("\\nu_{\\mu}");
    break;
  case -14:
    return("\\widetilde{\\nu_{\\mu}}");
    break;
  case 15: 
    return ("\\tau^{-}");
    break;
  case -15:
    return("\\tau^{+}");
    break;
  case 16: 
    return ("\\nu_{\\tau}");
    break;
  case -16:
    return("\\widetilde{\\nu_{\\tau}}");
    break;
  case 17: 
    return ("\\chi^{-}");
    break;
  case -17:
    return("\\chi^{+}");
    break;
  case 18: 
    return ("\\nu_{\\chi}");
    break;
  case -18:
    return("\\widetilde{\\nu_{\\chi}}");
    break;

  case 9:
  case 21:
    return ("g");
    break;
  case 22:
    return("\\gamma");
    break;
  case 23: 
    return ("Z^{0}");
    break;
  case 24:
    return("W^{+}");
    break;
  case -24:
    return("W^{-}");
    break;
  case 25: 
    return ("H^{0}_{1}");
    break;
  case 35:
    return("H^{0}_{2}");
    break;
  case 36: 
    return ("H^{0}_{3}");
    break;
  case 37:
    return("H^{+}");
    break;
  case -37:
    return("\\widetilde{H^{-}}");
    break;
  case 28:
    return("reggeon");
    break;
  case 29:
    return("pomeron");
    break;
  case 39:
    return("LQ");
    break;
  case 40:
    return("RQ");
    break;
    

  case 90:
    return("diquark");
    break;
  case 91:
    return("cluster");
    break;
   case 92:
    return("string");
    break;
   case 93:
    return("indep.");
    break;

  case 1103:
    return("dd_{1}");
    break;
  case 2101:
    return("ud_{0}");
    break;
  case 2103:
    return("ud_{1}");
    break;
  case 2203:
    return("uu_{1}");
    break;
  case 3101:
    return("sd_{0}");
    break;
  case 3103:
    return("sd_{1}");
    break;
  case 3201:
    return("su_{0}");
    break;
  case 3203:
    return("su_{1}");
    break;
  case 3303:
    return("ss_{1}");
    break;



  case -1103:
    return("\\widetilde{dd_{1}}");
    break;
  case -2101:
    return("\\widetilde{ud_{0}}");
    break;
  case -2103:
    return("\\widetilde{ud_{1}}");
    break;
  case -2203:
    return("\\widetilde{uu_{1}}");
    break;
  case -3101:
    return("\\widetilde{sd_{0}}");
    break;
  case -3103:
    return("\\widetilde{sd_{1}}");
    break;
  case -3201:
    return("\\widetilde{su_{0}}");
    break;
  case -3203:
    return("\\widetilde{su_{1}}");
    break;
  case -3303:
    return("\\widetilde{ss_{1}}");
    break;

  case 211:
    return("\\pi^{+}");
    break;
  case -211:
    return("\\pi^{-}");
    break;
  case 111:
    return("\\pi^{0}");
    break;
  case 311:
    return("K^{0}");
    break;
  case -311:
    return("\\widetilde{K^{0}}");
    break;
  case 321:
    return("K^{+}");
    break;
  case -321:
    return("K^{-}");
    break;
  case 411:
    return("D^{+}");
    break;
  case -411:
    return("D^{-}");
    break;
  case 421:
    return("D^{0}");
    break;
  case -421:
    return("\\widetilde{D^{0}}");
    break;
  case 431:
    return("D_{s}^{+}");
    break;
  case -431:
    return("D_{s}^{-}");
    break;
  case 511:
    return("B^{0}");
    break;
  case -511:
    return("\\widetilde{B^{0}}");
    break;
  case 521:
    return("B^{+}");
    break;
  case -521:
    return("B^{-}");
    break;
  case 531:
    return("B_{s}^{0}");
    break;
  case -531:
    return("\\widetilde{B_{s}^{0}}");
    break;
  case 541:
    return("B_{c}^{+}");
    break;
  case -541:
    return("B_{c}^{-}");
    break;
  case 221:
    return("\\eta");
    break;
  case 331:
    return("\\eta'");
    break;
  case 441:
    return("\\eta_{c}");
    break;
  case 551:
    return("\\eta_{b}");
    break;
  case 661:
    return("\\eta_{t}");
    break;
  case 130:
    return("K_{L}^{0}");
    break;
  case 310:
    return("K_{S}^{0}");
    break;
  case -130:
    return("\\widetilde{K_{L}^{0}}");
    break;
  case -310:
    return("\\widetilde{K_{S}^{0}}");
    break;
  case 213:
    return("\\rho^{+}");
    break;
  case -213:
    return("\\rho^{-}");
    break;
  case 313:
    return("K^{*0}");
    break;
  case -313:
    return("\\widetilde{K^{*0}}");
    break;
  case 323:
    return("K^{*+}");
    break;
  case -323:
    return("K^{*-}");
    break;
  case 413:
    return("D^{*+}");
    break;
  case -413:
    return("D^{*-}");
    break;
  case 423:
    return("D^{*0}");
    break;
  case -423:
    return("\\widetilde{D^{*0}}");
    break;
  case 433:
    return("{D^{*}_{s}}^{+}");
    break;
  case -433:
    return("{D^{*}_{s}}^{-}");
    break;
  case 513:
    return("B^{*0}");
    break;
  case -513:
    return("\\widetilde{B^{*0}}");
    break;
  case 523:
    return("B^{*+}");
    break;
  case -523:
    return("B^{*-}");
    break;
  case 533:
    return("{B^{*}_{s}}^{0}");
    break;
  case -533:
    return("\\widetilde{{B^{*}_{s}}^{0}}");
    break;
  case 543:
    return("{B^{*}_{c}}^{+}");
    break;
  case -543:
    return("{B^{*}_{c}}^{-}");
    break;
  case 113:
    return("\\rho^{0}");
    break;
  case 223:
    return("\\omega");
    break;
  case 333:
    return("\\phi");
    break;
  case 443:
    return("J/\\psi");
    break;
  case 553:
    return("\\Upsilon");
    break;
  case 663:
    return("\\Theta");
    break;

  case 10213:
    return("b_{1}^{+}");
    break;
  case -10213:
    return("b_{1}^{-}");
    break;
  case 10313:
    return("K_{1}^{0}");
    break;
  case -10313:
    return("\\widetilde{K_{1}^{0}}");
    break;
  case 10323:
    return("K_{1}^{+}");
    break;
  case -10323:
    return("K_{1}^{-}");
    break;
  case 10413:
    return("D_{1}^{+}");
    break;
  case -10413:
    return("D_{1}^{-}");
    break;
  case 10423:
    return("D_{1}^{0}");
    break;
  case -10423:
    return("\\widetilde{D_{1}^{0}}");
    break;
  case 10433:
    return("D_{1s}^{+}");
    break;
  case -10433:
    return("D_{1s}^{-}");
    break;
  case 10113:
    return("b_{1}^{0}");
    break;
  case 10223:
    return("h_{1}^{0}");
    break;
  case 10333:
    return("{h'}_{1}^{0}");
    break;
  case 10443:
    return("h_{1c}^{0}");
    break;
  case 10211:
    return("a_{0}^{+}");
    break;
  case -10211:
    return("a_{0}^{-}");
    break;
  case 10311:
    return("{K^{*}_{0}}^{0}");
    break;
  case -10311:
    return("\\widetilde{{K^{*}_{0}}^{0}}");
    break;
  case 10321:
    return("{K^{*}_{0}}^{+}");
    break;
  case -10321:
    return("{K^{*}_{0}}^{-}");
    break;
  case 10411:
    return("{D^{*}_{0}}^{+}");
    break;
  case -10411:
    return("{D^{*}_{0}}^{-}");
    break;
  case 10421:
    return("{D^{*}_{0}}^{0}");
    break;
  case -10421:
    return("\\widetilde{{D^{*}_{0}}^{0}}");
    break;
  case 10431:
    return("{D^{*}_{0s}}^{+}");
    break;
  case -10431:
    return("{D^{*}_{0s}}^{-}");
    break;
  case 10111:
    return("a_{0}^{0}");
    break;
  case 10221:
    return("f_{0}^{0}");
    break;
  case 10331:
    return("f'_{0}^{0}");
    break;
  case 10441:
    return("\\chi_{0c}^{0}");
    break;
  case 20213:
    return("a_{1}^{+}");
    break;
  case -20213:
    return("a_{1}^{-}");
    break;
  case 20313:
    return("{K^{*}_{1}}^{0}");
    break;
  case -20313:
    return("\\widetilde{{K^{*}_{1}}^{0}}");
    break;
  case 20323:
    return("{K^{*}_{1}}^{+}");
    break;
  case -20323:
    return("{K^{*}_{1}}^{-}");
    break;
  case 20413:
    return("{D^{*}_{1}}^{+}");
    break;
  case -20413:
    return("{D^{*}_{1}}^{-}");
    break;
  case 20423:
    return("{D^{*}_{1}}^{0}");
    break;
  case -20423:
    return("\\widetilde{{D^{*}_{1}}^{0}}");
    break;
  case 20433:
    return("{D^{*}_{1s}}^{+}");
    break;
  case -20433:
    return("{D^{*}_{1s}}^{-}");
    break;
  case 20113:
    return("a_{1}^{0}");
    break;
  case 20223:
    return("f_{1}^{0}");
    break;
  case 20333:
    return("f'_{1}^{0}");
    break;
  case 20443:
    return("\\chi_{1c}^{0}");
    break;
  case 215:
    return("a_{2}^{+}");
    break;
  case -215:
    return("a_{2}^{-}");
    break;
  case 315:
    return("{K^{*}_{2}}^{0}");
    break;
  case -315:
    return("\\widetilde{{K^{*}_{2}}^{0}}");
    break;
  case 325:
    return("{K^{*}_{2}}^{+}");
    break;
  case -325:
    return("{K^{*}_{2}}^{-}");
    break;
  case 415:
    return("{D^{*}_{2}}^{+}");
    break;
  case -415:
    return("{D^{*}_{2}}^{-}");
    break;
  case 425:
    return("{D^{*}_{2}}^{0}");
    break;
  case -425:
    return("\\widetilde{{D^{*}_{2}}^{0}}");
    break;
  case 435:
    return("{D^{*}_{2s}}^{+}");
    break;
  case -435:
    return("{D^{*}_{2s}}^{-}");
    break;
  case 115:
    return("a_{2}^{0}");
    break;
  case 225:
    return("f_{2}^{0}");
    break;
  case 335:
    return("{f'_{2}}^{0}");
    break;
  case 445:
    return("\\chi_{2c}^{0}");
    break;
  case 30443:
    return("\\psi'");
    break;
  case 30553:
    return("\\Upsilon'");
    break;

  case 2112:
    return("n^{0}");
    break;
  case 2212:
    return("p^{+}");
    break;
  case 3112:
    return("\\Sigma^{-}");
    break;
  case 3122:
    return("\\Lambda^{0}");
    break;
  case 3212:
    return("\\Sigma^{0}");
    break;
  case 3222:
    return("\\Sigma^{+}");
    break;
  case 3312:
    return("\\Xi^{-}");
    break;
  case 3322:
    return("\\Xi^{0}");
    break;
  case 4112:
    return("\\Sigma_{c}^{0}");
    break;
  case 4122:
    return("\\Lambda_{c}^{+}");
    break;
  case 4212:
    return("\\Sigma_{c}^{+}");
    break;
  case 4222:
    return("\\Sigma_{c}^{++}");
    break;
  case 4132:
    return("\\Xi_{c}^{0}");
    break;
  case 4312:
    return("\\Xi'_{c}^{0}");
    break;
  case 4232:
    return("\\Xi_{c}^{+}");
    break;
  case 4322:
    return("\\Xi'_{c}^{+}");
    break;
  case 4332:
    return("\\Omega_{c}^{0}");
    break;
  case 5112:
    return("\\Sigma_{b}^{-}");
    break;
  case 5222:
    return("\\Sigma_{b}^{+}");
    break;
  case 1114:
    return("\\Delta^{-}");
    break;
  case 2114:
    return("\\Delta^{0}");
    break;
  case 2214:
    return("\\Delta^{+}");
    break;
  case 2224:
    return("\\Delta^{++}");
    break;
  case 3114:
    return("\\Sigma^{*}^{-}");
    break;
  case 3224:
    return("\\Sigma^{*}^{+}");
    break;
  case 3314:
    return("\\Xi^{*}^{-}");
    break;
  case 3324:
    return("Xi^{*}^{0}");
    break;
  case 3334:
    return("\\Omega^{-}");
    break;
  case 4114:
    return("\\Sigma^{*}_{c}^{0}");
    break;
  case 4214:
    return("\\Sigma^{*}_{c}^{+}");
    break;
  case 4224:
    return("\\Sigma^{*}_{c}^{++}");
    break;
  case 4314:
    return("\\Xi^{*}_{c}^{0}");
    break;
  case 4324:
    return("\\Xi^{*}_{c}^{+}");
    break;
  case 4334:
    return("\\Omega^{*}_{c}^{0}");
    break;
  case 5114:
    return("\\Sigma^{*}_{b}^{-}");
    break;
  case 5214:
    return("\\Sigma^{*}_{b}^{0}");
    break;
  case 5224:
    return("\\Sigma^{*}_{b}^{+}");
    break;




  case -2112:
    return("\\widetilde{n^{0}}");
    break;
  case -2212:
    return("\\widetilde{p^{-}}");
    break;
  case -3112:
    return("\\widetilde{\\Sigma^{+}}");
    break;
  case -3122:
    return("\\widetilde{\\Lambda^{0}}");
    break;
  case -3212:
    return("\\widetilde{\\Sigma^{0}}");
    break;
  case -3222:
    return("\\widetilde{\\Sigma^{-}}");
    break;
  case -3312:
    return("\\widetilde{\\Xi^{+}}");
    break;
  case -3322:
    return("\\widetilde{\\Xi^{0}}");
    break;
  case -4112:
    return("\\widetilde{\\Sigma_{c}^{0}}");
    break;
  case -4122:
    return("\\widetilde{\\Lambda_{c}^{-}}");
    break;
  case -4212:
    return("\\widetilde{\\Sigma_{c}^{-}}");
    break;
  case -4222:
    return("\\widetilde{\\Sigma_{c}^{--}}");
    break;
  case -4132:
    return("\\widetilde{\\Xi_{c}^{0}}");
    break;
  case -4312:
    return("\\widetilde{\\Xi'_{c}^{0}}");
    break;
  case -4232:
    return("\\widetilde{\\Xi_{c}^{-}}");
    break;
  case -4322:
    return("\\widetilde{\\Xi'_{c}^{-}}");
    break;
  case -4332:
    return("\\widetilde{\\Omega_{c}^{0}}");
    break;
  case -5112:
    return("\\widetilde{\\Sigma_{b}^{+}}");
    break;
  case -5222:
    return("\\widetilde{\\Sigma_{b}^{-}}");
    break;
  case -1114:
    return("\\widetilde{\\Delta^{+}}");
    break;
  case -2114:
    return("\\widetilde{\\Delta^{0}}");
    break;
  case -2214:
    return("\\widetilde{\\Delta^{-}}");
    break;
  case -2224:
    return("\\widetilde{\\Delta^{--}}");
    break;
  case -3114:
    return("\\widetilde{\\Sigma^{*}^{+}}");
    break;
  case -3224:
    return("\\widetilde{\\Sigma^{*}^{-}}");
    break;
  case -3314:
    return("\\widetilde{\\Xi^{*}^{+}}");
    break;
  case -3324:
    return("\\widetilde{\\Xi^{*}^{0}}");
    break;
  case -3334:
    return("\\widetilde{\\Omega^{+}}");
    break;
  case -4114:
    return("\\widetilde{\\Sigma^{*}_{c}^{0}}");
    break;
  case -4214:
    return("\\widetilde{\\Sigma^{*}_}c}^{-}}");
    break;
  case -4224:
    return("\\widetilde{\\Sigma^{*}_{c}^{--}}");
    break;
  case -4314:
    return("\\widetilde{\\Xi^{*}_{c}^{0}}");
    break;
  case -4324:
    return("\\widetilde{\\Xi^{*}_{c}^{-}}");
    break;
  case -4334:
    return("\\widetilde{\\Omega^{*}_{c}^{0}}");
    break;
  case -5114:
    return("\\widetilde{\\Sigma^{*}_{b}^{+}}");
    break;
  case -5214:
    return("\\widetilde{\\Sigma^{*}_{b}^{0}}");
    break;
  case -5224:
    return("\\widetilde{\\Sigma^{*}_{b}^{-}}");
    break;



  default:
    sprintf(txtbuf,"%i",PDGCode);
    return (txtbuf);
    break;
  }


}



char*   HEPParticle::GetParticleName()
{
  // returns a string containig particle's name.

  return GetParticleName( GetPDGId() );
}

char*   HEPParticle::GetLaTeXName()
{

  return GetLaTeXName( GetPDGId() );
}


double  HEPParticle::GetCharge()
{
  // returns particle's electric charge.

  return GetCharge( GetPDGId() );
}



double  HEPParticle::GetSpin()
{
  // returns particle's spin
  return GetSpin( GetPDGId() );
}

void HEPParticle::Diff(HEPParticle &p)
{
  // Shows a differences between this particle and another.
  // For floats absolute and relative difference is calculated
  // (the latter given in percents). If compared values is equal
  // to 0.0 relative difference (x2-x1)/x1 cannot be calculated.
  // In this case displayed relative difference will be -100%
  //
  // example:
  //
  // HEPParticle p1 = MyEvent1.GetParticle(1);
  // HEPParticle p2 = MyEvent2.GetParticle(1);
  // /* print differences*/
  // p1.Diff(p2);
  
  if (GetPDGId() != p.GetPDGId() )
    printf("PDG:            %15i | %15i\n",GetPDGId(),p.GetPDGId());
  
  if (GetStatus() != p.GetStatus()   )
    printf("Status:         %15i | %15i\n",GetStatus(),p.GetStatus());
  
  if (!EQUAL(GetPx(), p.GetPx())) 
    printf("Px:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetPx(),p.GetPx(),GetPx()-p.GetPx(),
	   (GetPx()!=0.0) ? 100*fabs((GetPx()-p.GetPx())/GetPx()): -100);
  
  if (!EQUAL(GetPy(), p.GetPy())) 
    printf("Py:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetPy(),p.GetPy(),GetPy()-p.GetPy(),
	   (GetPy()!=0.0) ? 100*fabs((GetPy()-p.GetPy())/GetPy()): -100);
  
  if (!EQUAL(GetPz(), p.GetPz())) 
    printf("Pz:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetPz(),p.GetPz(),GetPz()-p.GetPz(),
	   (GetPz()!=0.0) ? 100*fabs((GetPz()-p.GetPz())/GetPz()): -100);
  
  if (!EQUAL(GetE(), p.GetE())) 
    printf("E:              %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetE(),p.GetE(),GetE()-p.GetE(),
	   (GetE()!=0.0) ? 100*fabs((GetE()-p.GetE())/GetE()): -100);
  
  if (!EQUAL(GetM() , p.GetM())) 
    printf("M:              %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetM(),p.GetM(),GetM()-p.GetM(),
	   (GetM()!=0.0) ? 100*fabs((GetM()-p.GetM())/GetM()): -100);
  
  if (!EQUAL(GetVx() , p.GetVx())) 
    printf("Vx:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetVx(),p.GetVx(),GetVx()-p.GetVx(),
	   (GetVx()!=0.0) ? 100*fabs((GetVx()-p.GetVx())/GetVx()): -100);
  
  if (!EQUAL(GetVy() , p.GetVy())) 
    printf("Vy:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetVy(),p.GetVy(),GetVy()-p.GetVy(),
	   (GetVy()!=0.0) ? 100*fabs((GetVy()-p.GetVy())/GetVy()): -100);
  
  if (!EQUAL(GetVz() ,p.GetVz())) 
    printf("Vz:             %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetVz(),p.GetVz(),GetVz()-p.GetVz(),
	   (GetVz()!=0.0) ? 100*fabs((GetVz()-p.GetVz())/GetVz()): -100);
  
  if (!EQUAL(GetTau() , p.GetTau())) 
    printf("tau:            %15.9G | %15.9G   D=%15.9G (REL=%12.7f%%)\n",
	   GetTau(),p.GetTau(),GetTau()-p.GetTau(),
	   (GetTau()!=0.0) ? 100*fabs((GetTau()-p.GetTau())/GetTau()): -100);
  
  
  

}


void HEPParticle::ls(char *option)
{
   if (option==0) {
   
      printf("%3i%7s [%i] <%3i,%3i|%3i,%3i> (% 8.3f,% 8.3f,% 8.3f)%8.3f %8.3f\n",
          GetId(), GetParticleName(),GetStatus(),
	  GetMother(), GetMother2(), 
	  GetFirstDaughter(), GetLastDaughter(),
	  GetPx(), GetPy(), GetPz(), GetE(), GetM()
	  );
   
   }
   
   
}

HEPParticle* HEPParticle::Assign (HEPParticle &p)
{
//    SetId(p.GetId());
    SetMother(p.GetMother());
    SetMother2(p.GetMother2());
    SetFirstDaughter(p.GetFirstDaughter());
    SetLastDaughter(p.GetLastDaughter());
    SetE(p.GetE());
    SetPx(p.GetPx());
    SetPy(p.GetPy());
    SetPz(p.GetPz());
    SetM(p.GetM());
    SetPDGId(p.GetPDGId());
    SetStatus(p.GetStatus());
    SetVx(p.GetVx());
    SetVy(p.GetVy());
    SetVz(p.GetVz());
    SetTau(p.GetTau());

    return this;
}




MC4Vector HEPParticle::GetP4()
{
    return MC4Vector(GetE(), GetPx(), GetPy(), GetPz());
}

MC3Vector HEPParticle::GetP3()
{
    return MC3Vector(GetPx(), GetPy(), GetPz());
}

MC3Vector HEPParticle::GetV3()
{
    return MC3Vector(GetVx(), GetVy(), GetVz());
}

void HEPParticle::SetP4(MC4Vector &v)
{
    SetE(v.GetX0());
    SetPx(v.GetX1());
    SetPy(v.GetX2());
    SetPz(v.GetX3());
}
void HEPParticle::SetP3(MC3Vector &v)
{
    SetPx(v.GetX0());
    SetPy(v.GetX1());
    SetPz(v.GetX2());
}

void HEPParticle::SetV3(MC3Vector &v)
{
    SetVx(v.GetX0());
    SetVy(v.GetX1());
    SetVz(v.GetX2());
}




#ifdef _USE_ROOT_
void HEPParticle::Streamer(TBuffer &)
{
  // streamer class for ROOT compatibility
}
void HEPParticleList::Streamer(TBuffer &)
{
  // streamer class for ROOT compatibility
}
void HEPParticleListIterator::Streamer(TBuffer &)
{
  // streamer class for ROOT compatibility
}
#endif