PWGCF/FEMTOSCOPY/AliFemto/AliFemtoModelWeightGeneratorBasic.cxx

   1 ////////////////////////////////////////////////////////////////////////////////
   2 ///                                                                          ///
   3 /// AliFemtoModelWeightGeneratorBasic -  basic femtoscopic weight generator  ///
   4 /// only return a simple                                                          ///
   5 /// Authors: Adam Kisiel kisiel@mps.ohio-state.edu                           ///
   6 ///                                                                          ///
   7 ////////////////////////////////////////////////////////////////////////////////
   8 #ifdef __ROOT__
   9   ClassImp(AliFemtoModelWeightGeneratorBasic, 1)
  10 #endif
  11
  12 #include "AliFemtoModelWeightGeneratorBasic.h"
  13 #include "AliFemtoModelHiddenInfo.h"
  14
  15 //________________________
  16 AliFemtoModelWeightGeneratorBasic::AliFemtoModelWeightGeneratorBasic():
  17   AliFemtoModelWeightGenerator()
  18 {
  19   /* no-op */
  20 }
  21 //________________________
  22 AliFemtoModelWeightGeneratorBasic::AliFemtoModelWeightGeneratorBasic(const AliFemtoModelWeightGeneratorBasic &aModel) :
  23   AliFemtoModelWeightGenerator(aModel)
  24 {
  25   /* no-op */
  26 }
  27
  28 //________________________
  29 AliFemtoModelWeightGeneratorBasic::~AliFemtoModelWeightGeneratorBasic()
  30 {
  31   /* no-op */
  32 }
  33
  34 AliFemtoModelWeightGeneratorBasic& AliFemtoModelWeightGeneratorBasic::operator=(const AliFemtoModelWeightGeneratorBasic &aModel)
  35 {
  36   if (this != &aModel) {
  37     AliFemtoModelWeightGenerator::operator=(aModel);
  38   }
  39
  40   return *this;
  41 }
  42
  43 //________________________
  44 Double_t AliFemtoModelWeightGeneratorBasic::GenerateWeight(AliFemtoPair *aPair)
  45 {
  46   // Generate a simple femtoscopic weight coming only from
  47   // quantum statistics - symmetrization or anti-symmetrization
  48   // of the pair wave function
  49
  50   // Get hidden information pointers
  51   //AliFemtoModelHiddenInfo *inf1 = (AliFemtoModelHiddenInfo *) aPair->Track1()->HiddenInfo();
  52   //AliFemtoModelHiddenInfo *inf2 = (AliFemtoModelHiddenInfo *) aPair->Track2()->HiddenInfo();
  53   AliFemtoTrack *inf1 = (AliFemtoTrack *) aPair->Track1()->Track();
  54   AliFemtoTrack *inf2 = (AliFemtoTrack *) aPair->Track2()->Track();
  55
  56   // Calculate pair variables
  57   Double_t tPx = inf1->GetTrueMomentum()->x()+inf2->GetTrueMomentum()->x();
  58   Double_t tPy = inf1->GetTrueMomentum()->y()+inf2->GetTrueMomentum()->y();
  59   Double_t tPz = inf1->GetTrueMomentum()->z()+inf2->GetTrueMomentum()->z();
  60   //  double tE  = inf1->GetTrueMomentum()->e +inf2->GetTrueMomentum()->.e;
  61   Double_t tM1 = inf1->GetMass();
  62   Double_t tM2 = inf2->GetMass();
  63   Double_t tE1 = sqrt(tM1*tM1 + inf1->GetTrueMomentum()->Mag2());
  64   Double_t tE2 = sqrt(tM2*tM2 + inf2->GetTrueMomentum()->Mag2());
  65   Double_t tE  = tE1 + tE2;
  66   Double_t tPt = tPx*tPx + tPy*tPy;
  67   Double_t tMt = tE*tE - tPz*tPz;//mCVK;
  68   Double_t tM  = sqrt(tMt - tPt);
  69   tMt = sqrt(tMt);
  70   tPt = sqrt(tPt);
  71   Double_t tBetat = tPt/tMt;
  72
  73   // Boost to LCMS
  74   Double_t tBeta = tPz/tE;
  75   Double_t tGamma = tE/tMt;
  76   fKStarLong = tGamma * (inf1->GetTrueMomentum()->z() - tBeta * tE1);
  77   Double_t tE1L = tGamma * (tE1  - tBeta * inf1->GetTrueMomentum()->z());
  78
  79   // Transform positions to LCMS
  80 //   Double_t tP1zl = tGamma * (inf1->GetEmissionPoint()->z() - tBeta * inf1->GetEmissionPoint()->t());
  81 //   Double_t tP1tl = tGamma * (inf1->GetEmissionPoint()->t() - tBeta * inf1->GetEmissionPoint()->z());
  82
  83 //   Double_t tP2zl = tGamma * (inf2->GetEmissionPoint()->z() - tBeta * inf2->GetEmissionPoint()->t());
  84 //   Double_t tP2tl = tGamma * (inf2->GetEmissionPoint()->t() - tBeta * inf2->GetEmissionPoint()->z());
  85
  86 //   Double_t tP1pzl = tGamma * (inf1->GetTrueMomentum()->z() - tBeta * tE1);
  87 //   Double_t tP1el  = tGamma * (tE1  - tBeta * inf1->GetTrueMomentum()->z());
  88
  89 //   Double_t tP2pzl = tGamma * (inf2->GetTrueMomentum()->z() - tBeta * tE2);
  90 //   Double_t tP2el  = tGamma * (tE2  - tBeta * inf2->GetTrueMomentum()->z());
  91
  92   // Rotate in transverse plane
  93   fKStarOut  = ( inf1->GetTrueMomentum()->x()*tPx + inf1->GetTrueMomentum()->y()*tPy)/tPt;
  94   fKStarSide = (-inf1->GetTrueMomentum()->x()*tPy + inf1->GetTrueMomentum()->y()*tPx)/tPt;
  95
  96 //   Double_t tP1pxl = fKStarOut;
  97 //   Double_t tP1pyl = fKStarSide;
  98
  99 //   Double_t tP2pxl = (inf2->GetTrueMomentum()->x()*tPx + inf2->GetTrueMomentum()->y()*tPy)/tPt;
 100 //   Double_t tP2pyl = (inf2->GetTrueMomentum()->y()*tPx - inf2->GetTrueMomentum()->x()*tPy)/tPt;;
 101
 102 //   Double_t tKO = tP1pxl - tP2pxl;
 103 //   Double_t tKS = tP1pyl - tP2pyl;
 104 //   Double_t tKL = tP1pzl - tP2pzl;
 105 //   Double_t tDE = tP1el  - tP2el;
 106
 107   // save the rotated coordinates in LCMS variables
 108 //   Double_t tP1xl = ( inf1->GetEmissionPoint()->x()*tPx + inf1->GetEmissionPoint()->y()*tPy)/tPt;
 109 //   Double_t tP1yl = (-inf1->GetEmissionPoint()->x()*tPy + inf1->GetEmissionPoint()->y()*tPx)/tPt;
 110
 111 //   Double_t tP2xl = ( inf2->GetEmissionPoint()->x()*tPx + inf2->GetEmissionPoint()->y()*tPy)/tPt;
 112 //   Double_t tP2yl = (-inf2->GetEmissionPoint()->x()*tPy + inf2->GetEmissionPoint()->y()*tPx)/tPt;
 113
 114   // Boost to pair cms
 115   fKStarOut = tMt/tM * (fKStarOut - tPt/tMt * tE1L);
 116
 117   tBetat = tPt/tMt;
 118 //   Double_t tGammat = 1.0/sqrt(1.0-tBetat*tBetat);
 119
 120 //   Double_t tP1xp = tGammat*(tP1xl - tBetat*tP1tl);
 121 //   Double_t tP1tp = tGammat*(tP1tl - tBetat*tP1xl);
 122
 123 //   Double_t tP2xp = tGammat*(tP2xl - tBetat*tP2tl);
 124 //   Double_t tP2tp = tGammat*(tP2tl - tBetat*tP2xl);
 125
 126 //   Double_t tRO = (tP1xl - tP2xl)/0.197327;
 127 //   Double_t tRS = (tP1yl - tP2yl)/0.197327;
 128 //   Double_t tRL = (tP1zl - tP2zl)/0.197327;
 129 //   Double_t tDT = (tP1tl - tP2tl)/0.197327;
 130
 131   Double_t tDX = inf1->GetEmissionPoint()->x()-inf2->GetEmissionPoint()->x();
 132   Double_t tDY = inf1->GetEmissionPoint()->y()-inf2->GetEmissionPoint()->y();
 133   Double_t tRLong = inf1->GetEmissionPoint()->z()-inf2->GetEmissionPoint()->z();
 134   Double_t tDTime = inf1->GetEmissionPoint()->t()-inf2->GetEmissionPoint()->t();
 135
 136   Double_t tROut = (tDX*tPx + tDY*tPy)/tPt;
 137   Double_t tRSide = (-tDX*tPy + tDY*tPx)/tPt;
 138
 139   fRStarSide = tRSide;
 140   Double_t tRSS = fRStarSide/0.197327;
 141
 142   fRStarLong = tGamma*(tRLong - tBeta* tDTime);
 143   Double_t tDTimePairLCMS = tGamma*(tDTime - tBeta* tRLong);
 144
 145   Double_t tRLS = fRStarLong/0.197327;
 146   tBeta = tPt/tMt;
 147   tGamma = tMt/tM;
 148
 149   fRStarOut = tGamma*(tROut - tBeta* tDTimePairLCMS);
 150   Double_t tROS = fRStarOut/0.197327;
 151 //   Double_t tDTimePairCMS = tGamma*(tDTimePairLCMS - tBeta* tROut);
 152   fRStar = ::sqrt(fRStarOut*fRStarOut + fRStarSide*fRStarSide +
 153                            fRStarLong*fRStarLong);
 154   fKStar = ::sqrt(fKStarOut*fKStarOut + fKStarSide*fKStarSide + fKStarLong*fKStarLong);
 155 //   Double_t tRSt = fRStar/0.197327;
 156
 157   if (fPairType != fgkPairTypeNone) {
 158     if ((fPairType == PionPlusPionPlus()) || (fPairType == KaonPlusKaonPlus()))
 159       return 1.0 + cos (2*(fKStarOut * tROS + fKStarSide * tRSS + fKStarLong * tRLS));
 160     else if (fPairType == ProtonProton())
 161       return 1.0 - 0.5 * cos (2*(fKStarOut * tROS + fKStarSide * tRSS + fKStarLong * tRLS));
 162     else
 163       return 1.0;
 164   }
 165   else {
 166     Int_t tPairType = GetPairTypeFromPair(aPair);
 167     if ((tPairType == PionPlusPionPlus()) || (tPairType == KaonPlusKaonPlus()))
 168       return 1.0 + cos (2*(fKStarOut * tROS + fKStarSide * tRSS + fKStarLong * tRLS));
 169     else if (tPairType == ProtonProton())
 170       return 1.0 - 0.5 * cos (2*(fKStarOut * tROS + fKStarSide * tRSS + fKStarLong * tRLS));
 171     else
 172       return 1.0;
 173
 174   }
 175 }
 176
 177 //________________________
 178 void     AliFemtoModelWeightGeneratorBasic::SetPairType(Int_t aPairType)
 179 {
 180   AliFemtoModelWeightGenerator::SetPairType(aPairType);
 181 }
 182 //________________________
 183 void     AliFemtoModelWeightGeneratorBasic::SetPairTypeFromPair(AliFemtoPair *aPair)
 184 {
 185   AliFemtoModelWeightGenerator::SetPairTypeFromPair(aPair);
 186 }
 187 //________________________
 188 Int_t    AliFemtoModelWeightGeneratorBasic::GetPairType() const
 189 {
 190   return AliFemtoModelWeightGenerator::GetPairType();
 191 }
 192 //________________________
 193 AliFemtoModelWeightGenerator* AliFemtoModelWeightGeneratorBasic::Clone() const
 194 {
 195   return GetGenerator();
 196 }
 197 //________________________
 198 AliFemtoModelWeightGenerator* AliFemtoModelWeightGeneratorBasic::GetGenerator() const
 199 {
 200   AliFemtoModelWeightGeneratorBasic *tGen = new AliFemtoModelWeightGeneratorBasic(*this);
 201   return tGen;
 202 }