[u/mrichter/AliRoot.git] / HBTAN / AliHBTPair.cxx

#include "AliHBTPair.h"
//_________________________________________________________________________
///////////////////////////////////////////////////////////////////////////
//
// class AliHBTPair
//
// class implements pair of particles and taking care of caluclation (almost)
// all of pair properties (Qinv, InvMass,...)
// 
// more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html
//
////////////////////////////////////////////////////////////////////////////

#include "AliHBTParticle.h"
#include "AliHBTWeights.h"
#include "AliHBTTrackPoints.h"

ClassImp(AliHBTPair)

/************************************************************************/
AliHBTPair::AliHBTPair(Bool_t rev):
 fPart1(0x0),
 fPart2(0x0),
 fSwapedPair(0x0),
 fQSideCMSLC(0.0),
 fQSideCMSLCNotCalc(kTRUE),
 fQOutCMSLC(0.0),
 fQOutCMSLCNotCalc(kTRUE),
 fQLongCMSLC(0.0),
 fQLongCMSLCNotCalc(kTRUE),
 fQInv(0.0),
 fQInvNotCalc(kTRUE),
 fInvMass(0.0),
 fInvMassNotCalc(kTRUE),
 fKt(0.0),
 fKtNotCalc(kTRUE),
 fKStar(0.0),
 fKStarNotCalc(kTRUE),
 fPInv(0.0),
 fQSide(0.0),
 fOut(0.0),
 fQLong(0.0),
 fMt(0.0),
 fMtNotCalc(kTRUE),
 fInvMassSqr(0.0),
 fMassSqrNotCalc(kTRUE),
 fQInvL(0.0),
 fQInvLNotCalc(kTRUE),
 fWeight(0.0),
 fWeightNotCalc(kTRUE),
 fAvarageDistance(0.0),
 fAvarageDistanceNotCalc(kTRUE),
 fPxSum(0.0),
 fPySum(0.0),
 fPzSum(0.0),
 fESum(0.0),
 fSumsNotCalc(kTRUE),
 fPxDiff(0.0),
 fPyDiff(0.0),
 fPzDiff(0.0),
 fEDiff(0.0),
 fDiffsNotCalc(kTRUE),
 fGammaCMSLC(0.0),
 fGammaCMSLCNotCalc(kTRUE),
 fChanged(kTRUE)
 {
//value of rev defines if it is Swaped
//if you pass kTRUE swpaped pair will NOT be created
//though you wont be able to get the swaped pair from this pair

  if(!rev) fSwapedPair = new AliHBTPair(kTRUE); //if false create swaped pair
  
 }
/************************************************************************/

AliHBTPair::AliHBTPair(AliHBTParticle* part1, AliHBTParticle* part2, Bool_t rev):
 fPart1(part1),
 fPart2(part2),
 fSwapedPair(0x0),
 fQSideCMSLC(0.0),
 fQSideCMSLCNotCalc(kTRUE),
 fQOutCMSLC(0.0),
 fQOutCMSLCNotCalc(kTRUE),
 fQLongCMSLC(0.0),
 fQLongCMSLCNotCalc(kTRUE),
 fQInv(0.0),
 fQInvNotCalc(kTRUE),
 fInvMass(0.0),
 fInvMassNotCalc(kTRUE),
 fKt(0.0),
 fKtNotCalc(kTRUE),
 fKStar(0.0),
 fKStarNotCalc(kTRUE),
 fPInv(0.0),
 fQSide(0.0),
 fOut(0.0),
 fQLong(0.0),
 fMt(0.0),
 fMtNotCalc(kTRUE),
 fInvMassSqr(0.0),
 fMassSqrNotCalc(kTRUE),
 fQInvL(0.0),
 fQInvLNotCalc(kTRUE),
 fWeight(0.0),
 fWeightNotCalc(kTRUE),
 fAvarageDistance(0.0),
 fAvarageDistanceNotCalc(kTRUE),
 fPxSum(0.0),
 fPySum(0.0),
 fPzSum(0.0),
 fESum(0.0),
 fSumsNotCalc(kTRUE),
 fPxDiff(0.0),
 fPyDiff(0.0),
 fPzDiff(0.0),
 fEDiff(0.0),
 fDiffsNotCalc(kTRUE),
 fGammaCMSLC(0.0),
 fGammaCMSLCNotCalc(kTRUE),
 fChanged(kTRUE)
 {
//value of rev defines if it is Swaped
//if you pass kTRUE swpaped pair will NOT be created
//though you wont be able to get the swaped pair from this pair

  if(!rev) fSwapedPair = new AliHBTPair(part2,part1,kTRUE); //if false create swaped pair
  
 }
/************************************************************************/
AliHBTPair::AliHBTPair(const AliHBTPair& in):
 TObject(in),
 fPart1(0x0),
 fPart2(0x0),
 fSwapedPair(0x0),
 fQSideCMSLC(0.0),
 fQSideCMSLCNotCalc(kTRUE),
 fQOutCMSLC(0.0),
 fQOutCMSLCNotCalc(kTRUE),
 fQLongCMSLC(0.0),
 fQLongCMSLCNotCalc(kTRUE),
 fQInv(0.0),
 fQInvNotCalc(kTRUE),
 fInvMass(0.0),
 fInvMassNotCalc(kTRUE),
 fKt(0.0),
 fKtNotCalc(kTRUE),
 fKStar(0.0),
 fKStarNotCalc(kTRUE),
 fPInv(0.0),
 fQSide(0.0),
 fOut(0.0),
 fQLong(0.0),
 fMt(0.0),
 fMtNotCalc(kTRUE),
 fInvMassSqr(0.0),
 fMassSqrNotCalc(kTRUE),
 fQInvL(0.0),
 fQInvLNotCalc(kTRUE),
 fWeight(0.0),
 fWeightNotCalc(kTRUE),
 fAvarageDistance(0.0),
 fAvarageDistanceNotCalc(kTRUE),
 fPxSum(0.0),
 fPySum(0.0),
 fPzSum(0.0),
 fESum(0.0),
 fSumsNotCalc(kTRUE),
 fPxDiff(0.0),
 fPyDiff(0.0),
 fPzDiff(0.0),
 fEDiff(0.0),
 fDiffsNotCalc(kTRUE),
 fGammaCMSLC(0.0),
 fGammaCMSLCNotCalc(kTRUE),
 fChanged(kTRUE)
{
 //cpy constructor
 in.Copy(*this);
}
/************************************************************************/

AliHBTPair& AliHBTPair::operator=(const AliHBTPair& in)
{
 //Assigment operator
 in.Copy(*this);
 return *this;
}
/************************************************************************/

Double_t AliHBTPair::GetInvMass()
{
//Returns qinv value for a pair
  if(fInvMassNotCalc)
   {
     CalculateInvMassSqr(); //method is inline so we not waste th time for jumping into method 
     
     if(fInvMassSqr<0)  fInvMass = TMath::Sqrt(-fInvMassSqr);
     else fInvMass = TMath::Sqrt(fInvMassSqr); 
     
     fInvMassNotCalc = kFALSE;
   }
  return fInvMass;
}
/************************************************************************/

Double_t AliHBTPair::GetQSideCMSLC()
{
//return Q Side in Central Of Mass System in Longitudialy Comoving Frame
 
  if (fQSideCMSLCNotCalc)
   {
    fQSideCMSLC = (fPart1->Px()*fPart2->Py()-fPart2->Px()*fPart1->Py())/GetKt();
    fQSideCMSLCNotCalc = kFALSE;
   }
  return fQSideCMSLC;
}
/************************************************************************/

Double_t AliHBTPair::GetQOutCMSLC()
{
 //caculates Qout in Center Of Mass Longitudionally Co-Moving
 if(fQOutCMSLCNotCalc)
  {
   CalculateSums();
   CalculateDiffs();

   if (fPart1->GetMass() != fPart2->GetMass())
    {
/*    
      //STAR algorithm
      Double_t beta  = fPzSum/fESum;
      Double_t gamma = GetGammaToCMSLC();
      Double_t el = gamma * (fPart1->Energy() - beta * fPart1->Pz());
      Double_t x  = ( fPart1->Px()*fPxSum + fPart1->Py()*fPySum) / ( 2.0*GetKt() );
      beta  = 2.0*GetKt()/GetMt();
      gamma = GetMt()/GetQInv();
      fQOutCMSLC = gamma * (x - beta * el);
*/

      //beta=fPzSum/fESum;    // Longit. V == beta
      Double_t beta=fPzSum/fESum;
      Double_t gamma = GetGammaToCMSLC();
      
      Double_t cosphi=fPxSum/(2.0*GetKt());  // cos(phi)
      Double_t sinphi=fPySum/(2.0*GetKt()); // sin(phi)
      
//      ROTATE(part1Px,part1Py,SPHI,CPHI,part1Px,part1Py);//ROT8
//      ROTATE(part2Px,part2Py,SPHI,CPHI,part2Px,part2Py);//ROT8
      Double_t tmp;
      tmp = fPart1->Px()*cosphi + fPart1->Py()*sinphi;
      Double_t part1Py = fPart1->Py()*cosphi - fPart1->Px()*sinphi;
      Double_t part1Px = tmp;

      tmp = fPart2->Px()*cosphi + fPart2->Py()*sinphi;
      Double_t part2Py = fPart2->Py()*cosphi - fPart2->Px()*sinphi;
      Double_t part2Px = tmp;
      
      
//      LTR(part1Pz,E1,beta,GetGammaToCMSLC(),part1Pz,E1a);
//      LTR(part2Pz,E2,beta,GetGammaToCMSLC(),part2Pz,E2a);
      Double_t part1Pz=gamma*(fPart1->Pz()-beta*fPart1->Energy());
      Double_t part2Pz=gamma*(fPart2->Pz()-beta*fPart2->Energy());

      Double_t part1P2=part1Px*part1Px+part1Py*part1Py+part1Pz*part1Pz;
      Double_t part2P2=part2Px*part2Px+part2Py*part2Py+part2Pz*part2Pz;
      Double_t part1E=TMath::Sqrt(fPart1->GetMass()*fPart1->GetMass()+part1P2);
      Double_t part2E=TMath::Sqrt(fPart2->GetMass()*fPart2->GetMass()+part2P2);
      Double_t sumE=part1E+part2E;
      Double_t sumPx=part1Px+part2Px;
      Double_t sumPy=part1Py+part2Py;
      Double_t sumPZ=part1Pz+part2Pz;
      Double_t sumP2=sumPx*sumPx+sumPy*sumPy+sumPZ*sumPZ;

      Double_t relmass=TMath::Sqrt(sumE*sumE-sumP2);
      Double_t hf = (fPart1->GetMass()*fPart1->GetMass() - fPart2->GetMass()*fPart2->GetMass())/(relmass*relmass);
      fQOutCMSLC=(part1Px-part2Px);//== id
      fQOutCMSLC=fQOutCMSLC-sumPx*hf; //sumPx == fPxSum ale po rotacji i transf
    }
   else
    {
      Double_t k2 = fPxSum*fPxDiff+fPySum*fPyDiff;
      fQOutCMSLC = 0.5*k2/GetKt();
   // if (non-id)  fQOutCMSLC=fQOutCMSLC - sumPx*HF;
    }

    
   fQOutCMSLCNotCalc = kFALSE;
  }
 return fQOutCMSLC;
}
/************************************************************************/

Double_t AliHBTPair::GetQLongCMSLC()
{
 //return Q Long in Central Of Mass System in Longitudialy Comoving Frame
 if (fQLongCMSLCNotCalc)
  {
    CalculateSums();
    CalculateDiffs();
    Double_t beta = fPzSum/fESum;
    fQLongCMSLC = GetGammaToCMSLC() * ( fPzDiff - beta*fEDiff );
    fQLongCMSLCNotCalc = kFALSE;
  }
 return fQLongCMSLC; 
}
/************************************************************************/

Double_t AliHBTPair::GetKt()
{
 //calculates the evarage momentum of the pair
  if(fKtNotCalc)
   { 
     CalculateSums();
     fKt =  0.5*TMath::Hypot(fPxSum,fPySum);
     fKtNotCalc = kFALSE;
   }
  return fKt;
}
/************************************************************************/

Double_t AliHBTPair::GetKStar()
{
  //calculates invariant velocity difference
  if (fKStarNotCalc)
   { 
    CalculateSums();

    Double_t ptrans = fPxSum*fPxSum + fPySum*fPySum;
    Double_t mtrans = fESum*fESum - fPzSum*fPzSum;
    if (ptrans > mtrans)
     {
       Error("GetKStar","Tranverse momentum bigger than transverse mass. Not normal for on-shell particles");
       Error("GetKStar","Particle1:");
       fPart1->Print();
       Error("GetKStar","Particle2:");
       fPart2->Print();
       Error("GetKStar","");
       
       fKStar = 10e5;
       fKStarNotCalc = kFALSE;
       return fKStar;
     }
    Double_t pinv =   TMath::Sqrt(mtrans - ptrans);

    Double_t q = (fPart1->GetMass()*fPart1->GetMass() - fPart2->GetMass()*fPart2->GetMass())/pinv;
    
    CalculateQInvL();
    
    q = q*q - fQInvL;
    if ( q < 0)
     {
       Info("GetKStar","Sqrt of negative number q = %f",q);
       Error("GetKStar","Particle1:");
       fPart1->Print();
       Error("GetKStar","Particle2:");
       fPart2->Print();
       fKStar = 10e5;
       fKStarNotCalc = kFALSE;
       return fKStar;
     }
     
    q = TMath::Sqrt(q);
    fKStar = q/2.;
    fKStarNotCalc = kFALSE;
   }
  return fKStar;
}
/************************************************************************/

Double_t AliHBTPair::GetQInv()
{
//returns Qinv 
//warning for non-id particles you want to use 2*KStar
  if(fQInvNotCalc)
   {
    CalculateQInvL();
    fQInv = TMath::Sqrt(TMath::Abs(fQInvL));
    fQInvNotCalc = kFALSE;
   }
  return fQInv;
}
/************************************************************************/

Double_t AliHBTPair::GetGammaToCMSLC()
{
  //calculates gamma factor of the boost to CMSLC
  if(fGammaCMSLCNotCalc)
   {
     CalculateSums();
     Double_t beta = fPzSum/fESum;
     fGammaCMSLC = 1.0/TMath::Sqrt(1.0 - beta*beta);
     fGammaCMSLCNotCalc = kFALSE;
   }
  return fGammaCMSLC;
}
/************************************************************************/

Double_t AliHBTPair::GetMt()
{
  //Calculates transverse mass of the pair
  if (fMtNotCalc)
   {
     CalculateSums();
     fMt = TMath::Sqrt(fESum*fESum - fPzSum*fPzSum);
     fMtNotCalc = kFALSE;
   }
   return fMt;
}
/************************************************************************/

Double_t AliHBTPair::GetWeight()
{
  //returns and buffers weight for this pair
  if (fWeightNotCalc)
   {
      fWeight = AliHBTWeights::Weight(this);
      fWeightNotCalc = kFALSE;  
   }
  return fWeight; 
}
/************************************************************************/

Double_t AliHBTPair::GetAvarageDistance()
{
//returns and buffers avarage distance between two tracks calculated 
// out of track points (see AliHBTTrackPoints class)

  if (fAvarageDistanceNotCalc)
   {
     fAvarageDistance = AvDistance();
     fAvarageDistanceNotCalc = kFALSE;
   }
  return fAvarageDistance;
}
/************************************************************************/

Double_t AliHBTPair::AvDistance()
{
  //returns avarage distance between two tracks in range 
  //as defined in Track-Points of AliHBTParticle
  //returns negative value if error uccured f.g. tracks do not have track-points
  AliHBTTrackPoints* tpts1 = fPart1->GetTrackPoints();
  if ( tpts1 == 0x0)
   {//it could be simulated pair
//     Warning("GetValue","Track 1 does not have Track Points. Pair NOT Passed.");
     return -1.0;
   }

  AliHBTTrackPoints* tpts2 = fPart2->GetTrackPoints();
  if ( tpts2 == 0x0)
   {
//     Warning("GetValue","Track 2 does not have Track Points. Pair NOT Passed.");
     return -1.0;
   }

  return tpts1->AvarageDistance(*tpts2);
}
Commit	Line	Data
1b446896	1	#include "AliHBTPair.h"
58ee8590	2	//_________________________________________________________________________
	3	///////////////////////////////////////////////////////////////////////////
	4	//
	5	// class AliHBTPair
	6	//
	7	// class implements pair of particles and taking care of caluclation (almost)
	8	// all of pair properties (Qinv, InvMass,...)
d012b7d0	9	//
58ee8590	10	// more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html
	11	//
	12	////////////////////////////////////////////////////////////////////////////
	13
1b446896	14	#include "AliHBTParticle.h"
dd82cadc	15	#include "AliHBTWeights.h"
087f87e7	16	#include "AliHBTTrackPoints.h"
1b446896	17
	18	ClassImp(AliHBTPair)
	19
1b446896	20	/************************************************************************/
ec6e4013	21	AliHBTPair::AliHBTPair(Bool_t rev):
ec6e4013	22	fPart1(0x0),
7261986f	23	fPart2(0x0),
	24	fSwapedPair(0x0),
	25	fQSideCMSLC(0.0),
	26	fQSideCMSLCNotCalc(kTRUE),
	27	fQOutCMSLC(0.0),
	28	fQOutCMSLCNotCalc(kTRUE),
	29	fQLongCMSLC(0.0),
	30	fQLongCMSLCNotCalc(kTRUE),
	31	fQInv(0.0),
	32	fQInvNotCalc(kTRUE),
	33	fInvMass(0.0),
	34	fInvMassNotCalc(kTRUE),
	35	fKt(0.0),
	36	fKtNotCalc(kTRUE),
	37	fKStar(0.0),
	38	fKStarNotCalc(kTRUE),
	39	fPInv(0.0),
	40	fQSide(0.0),
	41	fOut(0.0),
	42	fQLong(0.0),
	43	fMt(0.0),
	44	fMtNotCalc(kTRUE),
	45	fInvMassSqr(0.0),
	46	fMassSqrNotCalc(kTRUE),
	47	fQInvL(0.0),
	48	fQInvLNotCalc(kTRUE),
dd82cadc	49	fWeight(0.0),
dd82cadc	50	fWeightNotCalc(kTRUE),
6d7ee019	51	fAvarageDistance(0.0),
6d7ee019	52	fAvarageDistanceNotCalc(kTRUE),
7261986f	53	fPxSum(0.0),
	54	fPySum(0.0),
	55	fPzSum(0.0),
	56	fESum(0.0),
	57	fSumsNotCalc(kTRUE),
	58	fPxDiff(0.0),
	59	fPyDiff(0.0),
	60	fPzDiff(0.0),
	61	fEDiff(0.0),
	62	fDiffsNotCalc(kTRUE),
	63	fGammaCMSLC(0.0),
	64	fGammaCMSLCNotCalc(kTRUE),
	65	fChanged(kTRUE)
1b446896	66	{
	67	//value of rev defines if it is Swaped
	68	//if you pass kTRUE swpaped pair will NOT be created
	69	//though you wont be able to get the swaped pair from this pair
	70
	71	if(!rev) fSwapedPair = new AliHBTPair(kTRUE); //if false create swaped pair
1b446896	72
	73	}
	74	/************************************************************************/
ec6e4013	75
	76	AliHBTPair::AliHBTPair(AliHBTParticle* part1, AliHBTParticle* part2, Bool_t rev):
	77	fPart1(part1),
7261986f	78	fPart2(part2),
	79	fSwapedPair(0x0),
	80	fQSideCMSLC(0.0),
	81	fQSideCMSLCNotCalc(kTRUE),
	82	fQOutCMSLC(0.0),
	83	fQOutCMSLCNotCalc(kTRUE),
	84	fQLongCMSLC(0.0),
	85	fQLongCMSLCNotCalc(kTRUE),
	86	fQInv(0.0),
	87	fQInvNotCalc(kTRUE),
	88	fInvMass(0.0),
	89	fInvMassNotCalc(kTRUE),
	90	fKt(0.0),
	91	fKtNotCalc(kTRUE),
	92	fKStar(0.0),
	93	fKStarNotCalc(kTRUE),
	94	fPInv(0.0),
	95	fQSide(0.0),
	96	fOut(0.0),
	97	fQLong(0.0),
	98	fMt(0.0),
	99	fMtNotCalc(kTRUE),
	100	fInvMassSqr(0.0),
	101	fMassSqrNotCalc(kTRUE),
	102	fQInvL(0.0),
	103	fQInvLNotCalc(kTRUE),
dd82cadc	104	fWeight(0.0),
dd82cadc	105	fWeightNotCalc(kTRUE),
6d7ee019	106	fAvarageDistance(0.0),
6d7ee019	107	fAvarageDistanceNotCalc(kTRUE),
7261986f	108	fPxSum(0.0),
	109	fPySum(0.0),
	110	fPzSum(0.0),
	111	fESum(0.0),
	112	fSumsNotCalc(kTRUE),
	113	fPxDiff(0.0),
	114	fPyDiff(0.0),
	115	fPzDiff(0.0),
	116	fEDiff(0.0),
	117	fDiffsNotCalc(kTRUE),
	118	fGammaCMSLC(0.0),
	119	fGammaCMSLCNotCalc(kTRUE),
	120	fChanged(kTRUE)
ec6e4013	121	{
	122	//value of rev defines if it is Swaped
	123	//if you pass kTRUE swpaped pair will NOT be created
	124	//though you wont be able to get the swaped pair from this pair
	125
	126	if(!rev) fSwapedPair = new AliHBTPair(part2,part1,kTRUE); //if false create swaped pair
ec6e4013	127
	128	}
	129	/************************************************************************/
d012b7d0	130	AliHBTPair::AliHBTPair(const AliHBTPair& in):
	131	TObject(in),
	132	fPart1(0x0),
	133	fPart2(0x0),
	134	fSwapedPair(0x0),
	135	fQSideCMSLC(0.0),
	136	fQSideCMSLCNotCalc(kTRUE),
	137	fQOutCMSLC(0.0),
	138	fQOutCMSLCNotCalc(kTRUE),
	139	fQLongCMSLC(0.0),
	140	fQLongCMSLCNotCalc(kTRUE),
	141	fQInv(0.0),
	142	fQInvNotCalc(kTRUE),
	143	fInvMass(0.0),
	144	fInvMassNotCalc(kTRUE),
	145	fKt(0.0),
	146	fKtNotCalc(kTRUE),
	147	fKStar(0.0),
	148	fKStarNotCalc(kTRUE),
	149	fPInv(0.0),
	150	fQSide(0.0),
	151	fOut(0.0),
	152	fQLong(0.0),
	153	fMt(0.0),
	154	fMtNotCalc(kTRUE),
	155	fInvMassSqr(0.0),
	156	fMassSqrNotCalc(kTRUE),
	157	fQInvL(0.0),
	158	fQInvLNotCalc(kTRUE),
	159	fWeight(0.0),
	160	fWeightNotCalc(kTRUE),
6d7ee019	161	fAvarageDistance(0.0),
6d7ee019	162	fAvarageDistanceNotCalc(kTRUE),
d012b7d0	163	fPxSum(0.0),
	164	fPySum(0.0),
	165	fPzSum(0.0),
	166	fESum(0.0),
	167	fSumsNotCalc(kTRUE),
	168	fPxDiff(0.0),
	169	fPyDiff(0.0),
	170	fPzDiff(0.0),
	171	fEDiff(0.0),
	172	fDiffsNotCalc(kTRUE),
	173	fGammaCMSLC(0.0),
	174	fGammaCMSLCNotCalc(kTRUE),
	175	fChanged(kTRUE)
	176	{
	177	//cpy constructor
	178	in.Copy(*this);
	179	}
	180	/************************************************************************/
	181
34914285	182	AliHBTPair& AliHBTPair::operator=(const AliHBTPair& in)
d012b7d0	183	{
	184	//Assigment operator
	185	in.Copy(*this);
	186	return *this;
	187	}
34914285	188	/************************************************************************/
ec6e4013	189
1b446896	190	Double_t AliHBTPair::GetInvMass()
1b446896	191	{
b928db6c	192	//Returns qinv value for a pair
1b446896	193	if(fInvMassNotCalc)
	194	{
	195	CalculateInvMassSqr(); //method is inline so we not waste th time for jumping into method
	196
	197	if(fInvMassSqr<0) fInvMass = TMath::Sqrt(-fInvMassSqr);
	198	else fInvMass = TMath::Sqrt(fInvMassSqr);
	199
	200	fInvMassNotCalc = kFALSE;
	201	}
	202	return fInvMass;
	203	}
	204	/************************************************************************/
6d7ee019	205
1b446896	206	Double_t AliHBTPair::GetQSideCMSLC()
1b446896	207	{
6d7ee019	208	//return Q Side in Central Of Mass System in Longitudialy Comoving Frame
1b446896	209
	210	if (fQSideCMSLCNotCalc)
	211	{
	212	fQSideCMSLC = (fPart1->Px()fPart2->Py()-fPart2->Px()fPart1->Py())/GetKt();
	213	fQSideCMSLCNotCalc = kFALSE;
	214	}
	215	return fQSideCMSLC;
	216	}
	217	/************************************************************************/
6d7ee019	218
1b446896	219	Double_t AliHBTPair::GetQOutCMSLC()
1b446896	220	{
58ee8590	221	//caculates Qout in Center Of Mass Longitudionally Co-Moving
1b446896	222	if(fQOutCMSLCNotCalc)
	223	{
	224	CalculateSums();
	225	CalculateDiffs();
951aadb9	226
	227	if (fPart1->GetMass() != fPart2->GetMass())
	228	{
	229	/*
	230	//STAR algorithm
	231	Double_t beta = fPzSum/fESum;
	232	Double_t gamma = GetGammaToCMSLC();
	233	Double_t el = gamma * (fPart1->Energy() - beta * fPart1->Pz());
	234	Double_t x = ( fPart1->Px()fPxSum + fPart1->Py()fPySum) / ( 2.0*GetKt() );
	235	beta = 2.0*GetKt()/GetMt();
	236	gamma = GetMt()/GetQInv();
	237	fQOutCMSLC = gamma * (x - beta * el);
	238	*/
	239
	240	//beta=fPzSum/fESum; // Longit. V == beta
	241	Double_t beta=fPzSum/fESum;
	242	Double_t gamma = GetGammaToCMSLC();
	243
	244	Double_t cosphi=fPxSum/(2.0*GetKt()); // cos(phi)
	245	Double_t sinphi=fPySum/(2.0*GetKt()); // sin(phi)
	246
	247	// ROTATE(part1Px,part1Py,SPHI,CPHI,part1Px,part1Py);//ROT8
	248	// ROTATE(part2Px,part2Py,SPHI,CPHI,part2Px,part2Py);//ROT8
	249	Double_t tmp;
	250	tmp = fPart1->Px()cosphi + fPart1->Py()sinphi;
	251	Double_t part1Py = fPart1->Py()cosphi - fPart1->Px()sinphi;
	252	Double_t part1Px = tmp;
	253
	254	tmp = fPart2->Px()cosphi + fPart2->Py()sinphi;
	255	Double_t part2Py = fPart2->Py()cosphi - fPart2->Px()sinphi;
	256	Double_t part2Px = tmp;
	257
	258
	259	// LTR(part1Pz,E1,beta,GetGammaToCMSLC(),part1Pz,E1a);
	260	// LTR(part2Pz,E2,beta,GetGammaToCMSLC(),part2Pz,E2a);
	261	Double_t part1Pz=gamma(fPart1->Pz()-betafPart1->Energy());
	262	Double_t part2Pz=gamma(fPart2->Pz()-betafPart2->Energy());
	263
	264	Double_t part1P2=part1Pxpart1Px+part1Pypart1Py+part1Pz*part1Pz;
	265	Double_t part2P2=part2Pxpart2Px+part2Pypart2Py+part2Pz*part2Pz;
	266	Double_t part1E=TMath::Sqrt(fPart1->GetMass()*fPart1->GetMass()+part1P2);
	267	Double_t part2E=TMath::Sqrt(fPart2->GetMass()*fPart2->GetMass()+part2P2);
	268	Double_t sumE=part1E+part2E;
	269	Double_t sumPx=part1Px+part2Px;
	270	Double_t sumPy=part1Py+part2Py;
	271	Double_t sumPZ=part1Pz+part2Pz;
	272	Double_t sumP2=sumPxsumPx+sumPysumPy+sumPZ*sumPZ;
	273
	274	Double_t relmass=TMath::Sqrt(sumE*sumE-sumP2);
	275	Double_t hf = (fPart1->GetMass()fPart1->GetMass() - fPart2->GetMass()fPart2->GetMass())/(relmass*relmass);
	276	fQOutCMSLC=(part1Px-part2Px);//== id
	277	fQOutCMSLC=fQOutCMSLC-sumPx*hf; //sumPx == fPxSum ale po rotacji i transf
	278	}
	279	else
	280	{
	281	Double_t k2 = fPxSumfPxDiff+fPySumfPyDiff;
	282	fQOutCMSLC = 0.5*k2/GetKt();
	283	// if (non-id) fQOutCMSLC=fQOutCMSLC - sumPx*HF;
	284	}
	285
	286
1b446896	287	fQOutCMSLCNotCalc = kFALSE;
	288	}
	289	return fQOutCMSLC;
	290	}
	291	/************************************************************************/
6d7ee019	292
1b446896	293	Double_t AliHBTPair::GetQLongCMSLC()
1b446896	294	{
58ee8590	295	//return Q Long in Central Of Mass System in Longitudialy Comoving Frame
1b446896	296	if (fQLongCMSLCNotCalc)
	297	{
	298	CalculateSums();
	299	CalculateDiffs();
	300	Double_t beta = fPzSum/fESum;
951aadb9	301	fQLongCMSLC = GetGammaToCMSLC() * ( fPzDiff - beta*fEDiff );
1b446896	302	fQLongCMSLCNotCalc = kFALSE;
	303	}
	304	return fQLongCMSLC;
	305	}
	306	/************************************************************************/
6d7ee019	307
1b446896	308	Double_t AliHBTPair::GetKt()
1b446896	309	{
58ee8590	310	//calculates the evarage momentum of the pair
1b446896	311	if(fKtNotCalc)
	312	{
	313	CalculateSums();
	314	fKt = 0.5*TMath::Hypot(fPxSum,fPySum);
	315	fKtNotCalc = kFALSE;
	316	}
	317	return fKt;
	318	}
	319	/************************************************************************/
	320
	321	Double_t AliHBTPair::GetKStar()
	322	{
58ee8590	323	//calculates invariant velocity difference
1b446896	324	if (fKStarNotCalc)
1b446896	325	{
1b446896	326	CalculateSums();
1b446896	327
d012b7d0	328	Double_t ptrans = fPxSumfPxSum + fPySumfPySum;
d012b7d0	329	Double_t mtrans = fESumfESum - fPzSumfPzSum;
eb890d59	330	if (ptrans > mtrans)
	331	{
	332	Error("GetKStar","Tranverse momentum bigger than transverse mass. Not normal for on-shell particles");
	333	Error("GetKStar","Particle1:");
	334	fPart1->Print();
	335	Error("GetKStar","Particle2:");
	336	fPart2->Print();
	337	Error("GetKStar","");
	338
	339	fKStar = 10e5;
	340	fKStarNotCalc = kFALSE;
	341	return fKStar;
	342	}
d012b7d0	343	Double_t pinv = TMath::Sqrt(mtrans - ptrans);
1b446896	344
d012b7d0	345	Double_t q = (fPart1->GetMass()fPart1->GetMass() - fPart2->GetMass()fPart2->GetMass())/pinv;
1b446896	346
	347	CalculateQInvL();
	348
d012b7d0	349	q = q*q - fQInvL;
d012b7d0	350	if ( q < 0)
4453ae67	351	{
eb890d59	352	Info("GetKStar","Sqrt of negative number q = %f",q);
	353	Error("GetKStar","Particle1:");
	354	fPart1->Print();
	355	Error("GetKStar","Particle2:");
	356	fPart2->Print();
	357	fKStar = 10e5;
	358	fKStarNotCalc = kFALSE;
	359	return fKStar;
4453ae67	360	}
4453ae67	361
d012b7d0	362	q = TMath::Sqrt(q);
d012b7d0	363	fKStar = q/2.;
1b446896	364	fKStarNotCalc = kFALSE;
	365	}
	366	return fKStar;
	367	}
	368	/************************************************************************/
	369
	370	Double_t AliHBTPair::GetQInv()
	371	{
951aadb9	372	//returns Qinv
951aadb9	373	//warning for non-id particles you want to use 2*KStar
1b446896	374	if(fQInvNotCalc)
	375	{
	376	CalculateQInvL();
30025bb4	377	fQInv = TMath::Sqrt(TMath::Abs(fQInvL));
1b446896	378	fQInvNotCalc = kFALSE;
1b446896	379	}
1b446896	380	return fQInv;
1b446896	381	}
30025bb4	382	/************************************************************************/
951aadb9	383
	384	Double_t AliHBTPair::GetGammaToCMSLC()
	385	{
58ee8590	386	//calculates gamma factor of the boost to CMSLC
951aadb9	387	if(fGammaCMSLCNotCalc)
	388	{
	389	CalculateSums();
	390	Double_t beta = fPzSum/fESum;
	391	fGammaCMSLC = 1.0/TMath::Sqrt(1.0 - beta*beta);
	392	fGammaCMSLCNotCalc = kFALSE;
	393	}
	394	return fGammaCMSLC;
	395	}
30025bb4	396	/************************************************************************/
30025bb4	397
951aadb9	398	Double_t AliHBTPair::GetMt()
951aadb9	399	{
58ee8590	400	//Calculates transverse mass of the pair
951aadb9	401	if (fMtNotCalc)
	402	{
	403	CalculateSums();
	404	fMt = TMath::Sqrt(fESumfESum - fPzSumfPzSum);
	405	fMtNotCalc = kFALSE;
	406	}
	407	return fMt;
	408	}
30025bb4	409	/************************************************************************/
1b446896	410
dd82cadc	411	Double_t AliHBTPair::GetWeight()
47d9a058	412	{
58ee8590	413	//returns and buffers weight for this pair
dd82cadc	414	if (fWeightNotCalc)
47d9a058	415	{
dd82cadc	416	fWeight = AliHBTWeights::Weight(this);
dd82cadc	417	fWeightNotCalc = kFALSE;
47d9a058	418	}
dd82cadc	419	return fWeight;
47d9a058	420	}
087f87e7	421	/************************************************************************/
	422
	423	Double_t AliHBTPair::GetAvarageDistance()
6d7ee019	424	{
	425	//returns and buffers avarage distance between two tracks calculated
	426	// out of track points (see AliHBTTrackPoints class)
	427
	428	if (fAvarageDistanceNotCalc)
	429	{
	430	fAvarageDistance = AvDistance();
	431	fAvarageDistanceNotCalc = kFALSE;
	432	}
	433	return fAvarageDistance;
	434	}
	435	/************************************************************************/
	436
	437	Double_t AliHBTPair::AvDistance()
087f87e7	438	{
	439	//returns avarage distance between two tracks in range
	440	//as defined in Track-Points of AliHBTParticle
	441	//returns negative value if error uccured f.g. tracks do not have track-points
	442	AliHBTTrackPoints* tpts1 = fPart1->GetTrackPoints();
	443	if ( tpts1 == 0x0)
	444	{//it could be simulated pair
	445	// Warning("GetValue","Track 1 does not have Track Points. Pair NOT Passed.");
	446	return -1.0;
	447	}
	448
	449	AliHBTTrackPoints* tpts2 = fPart2->GetTrackPoints();
	450	if ( tpts2 == 0x0)
	451	{
	452	// Warning("GetValue","Track 2 does not have Track Points. Pair NOT Passed.");
	453	return -1.0;
	454	}
	455
	456	return tpts1->AvarageDistance(*tpts2);
	457	}