PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoCutMonitorParticlePtPDG.cxx

   1 ////////////////////////////////////////////////////////////////////////////////
   2 //                                                                            //
   3 // AliFemtoCutMonitorParticlePtPDG - the cut monitor for particles to study    //
   4 // the difference between reconstructed and true momentum                     //
   5 //                                                                            //
   6 ////////////////////////////////////////////////////////////////////////////////
   7 #include "AliFemtoCutMonitorParticlePtPDG.h"
   8 #include "AliFemtoModelHiddenInfo.h"
   9 #include <TH1D.h>
  10 #include <TH2D.h>
  11 #include <TList.h>
  12 #include <TMath.h>
  13
  14 AliFemtoCutMonitorParticlePtPDG::AliFemtoCutMonitorParticlePtPDG():
  15   fPtPDG(0),ftpcHist(0),fPtGoodPi(0),fPtFakePi(0),fPtGoodK(0),fPtFakeK(0),
  16   fPtGoodP(0),fPtFakeP(0),fPtRPi(0),fPtRK(0),fPtRP(0),
  17   fPtContP(0),
  18   fPtContPi(0),
  19   fPtContMup(0),
  20   fPtContElp(0),
  21   fMass(0.13957)
  22 {
  23   // Default constructor
  24   fPtPDG = new TH2D("PtPDG", "PDG vs Pt",  5, .0, 5.0, 100, 0.1, 2.0);
  25   ftpcHist=new TH2D("tpcHist","TPC dE/dX vs momentum",100,0.1,2.7,100,0.,6.);
  26   fPtGoodPi = new TH1D("PtGoodPi", "good pions Pt",              100, 0.1, 2.0);
  27   fPtRPi = new TH1D("PtRPi", "right pdg pions Pt",               100, 0.1, 2.0);
  28   fPtFakePi = new TH1D("PtFakePi", "fake pions Pt",              100, 0.1, 2.0);
  29   fPtGoodK = new TH1D("PtGoodK", "good kaons Pt",                100, 0.1, 2.0);
  30   fPtRK = new TH1D("PtRK", "right pdg kaons Pt",                 100, 0.1, 2.0);
  31   fPtFakeK = new TH1D("PtFakeK", "fake kaons Pt",                100, 0.1, 2.0);
  32   fPtGoodP = new TH1D("PtGoodP", "good protons Pt",              100, 0.1, 2.0);
  33   fPtFakeP = new TH1D("PtFakeP", "fake protons Pt",              100, 0.1, 2.0);
  34   fPtRP = new TH1D("PtRP", "right pdg protons Pt",               100, 0.1, 2.0);
  35
  36   fPtContP = new TH1D("PtContP", "contamination Pt",               100, 0.1, 2.0);
  37   fPtContPi = new TH1D("PtContPi", "contamination Pt",               100, 0.1, 2.0);
  38   fPtContMup = new TH1D("PtContMup", "contamination Pt",               100, 0.1, 2.0);
  39   fPtContElp = new TH1D("PtContElp", "contamination Pt",               100, 0.1, 2.0);
  40
  41 }
  42
  43 AliFemtoCutMonitorParticlePtPDG::AliFemtoCutMonitorParticlePtPDG(const char *aName, float aMass):
  44   AliFemtoCutMonitor(),
  45   fPtPDG(0),ftpcHist(0),fPtGoodPi(0),fPtFakePi(0),fPtGoodK(0),fPtFakeK(0),
  46   fPtGoodP(0),fPtFakeP(0),fPtRPi(0),fPtRK(0),fPtRP(0),
  47   fPtContP(0),
  48   fPtContPi(0),
  49   fPtContMup(0),
  50   fPtContElp(0),
  51   fMass(aMass)
  52 {
  53   // Normal constructor
  54   char name[200];
  55   snprintf(name, 200, "PtPDG%s", aName);
  56   fPtPDG = new TH2D(name, "PDG vs Pt", 10, 0.0, 5.0, 100, 0.1, 2.0);
  57   snprintf(name, 200, "tpcHist%s", aName);
  58   ftpcHist=new TH2D(name,"TPC dE/dX vs momentum",100,0.1,2.7,100,0.,6.);
  59   snprintf(name, 200, "PtGoodPi%s", aName);
  60   fPtGoodPi = new TH1D(name, "good pions Pt",                    100, 0.1, 2.0);
  61   snprintf(name, 200, "PtFakePi%s", aName);
  62   fPtFakePi = new TH1D(name, "fake pions Pt",              100, 0.1, 2.0);
  63   snprintf(name, 200, "PtRPi%s", aName);
  64   fPtRPi = new TH1D(name, "right pdg pions Pt",               100, 0.1, 2.0);
  65   snprintf(name, 200, "PtGoodK%s", aName);
  66   fPtGoodK = new TH1D(name, "good kaons Pt",                     100, 0.1, 2.0);
  67   snprintf(name, 200, "PtFakeK%s", aName);
  68   fPtFakeK = new TH1D(name, "fake kaons Pt",                100, 0.1, 2.0);
  69   snprintf(name, 200, "PtRK%s", aName);
  70   fPtRK = new TH1D(name, "right pdg kaons Pt",                 100, 0.1, 2.0);
  71    snprintf(name, 200, "PtGoodP%s", aName);
  72   fPtGoodP = new TH1D(name, "good protons Pt",                     100, 0.1, 2.0);
  73   snprintf(name, 200, "PtFakeP%s", aName);
  74   fPtFakeP = new TH1D(name, "fake protons Pt",                100, 0.1, 2.0);
  75   snprintf(name, 200, "PtRP%s", aName);
  76   fPtRP = new TH1D(name, "right pdg protons Pt",                 100, 0.1, 2.0);
  77
  78   snprintf(name, 200, "PtContP%s", aName);
  79   fPtContP = new TH1D(name, "contamination",                 100, 0.1, 2.0);
  80   snprintf(name, 200, "PtContPi%s", aName);
  81   fPtContPi = new TH1D(name, "contamination",                 100, 0.1, 2.0);
  82   snprintf(name, 200, "PtContMup%s", aName);
  83   fPtContMup = new TH1D(name, "contamination",                 100, 0.1, 2.0);
  84   snprintf(name, 200, "PtContElp%s", aName);
  85   fPtContElp = new TH1D(name, "contamination",                 100, 0.1, 2.0);
  86
  87 }
  88
  89 AliFemtoCutMonitorParticlePtPDG::AliFemtoCutMonitorParticlePtPDG(const AliFemtoCutMonitorParticlePtPDG &aCut):
  90   AliFemtoCutMonitor(),
  91   fPtPDG(0),ftpcHist(0),fPtGoodPi(0),fPtFakePi(0),fPtGoodK(0),fPtFakeK(0),
  92   fPtGoodP(0),fPtFakeP(0),fPtRPi(0),fPtRK(0),fPtRP(0),
  93   fPtContP(0),
  94   fPtContPi(0),
  95   fPtContMup(0),
  96   fPtContElp(0),
  97   fMass(0.13957)
  98 {
  99   // copy constructor
 100   if (fPtPDG) delete fPtPDG;
 101   fPtPDG = new TH2D(*aCut.fPtPDG);
 102   ftpcHist= new TH2D(*aCut.ftpcHist);
 103   fPtGoodPi= new TH1D(*aCut.fPtGoodPi);
 104   fPtFakePi= new TH1D(*aCut.fPtFakePi);
 105   fPtGoodK= new TH1D(*aCut.fPtGoodK);
 106   fPtFakeK= new TH1D(*aCut.fPtFakePi);
 107   fPtGoodP= new TH1D(*aCut.fPtGoodP);
 108   fPtFakeP= new TH1D(*aCut.fPtFakePi);
 109   fPtRPi= new TH1D(*aCut.fPtRPi);
 110   fPtRK= new TH1D(*aCut.fPtRK);
 111   fPtRP= new TH1D(*aCut.fPtRP);
 112
 113   fPtContP= new TH1D(*aCut.fPtContP);
 114   fPtContPi= new TH1D(*aCut.fPtContPi);
 115   fPtContMup= new TH1D(*aCut.fPtContMup);
 116   fPtContElp= new TH1D(*aCut.fPtContElp);
 117
 118   fMass = aCut.fMass;
 119 }
 120
 121 AliFemtoCutMonitorParticlePtPDG::~AliFemtoCutMonitorParticlePtPDG()
 122 {
 123   // Destructor
 124   delete fPtPDG;
 125   delete ftpcHist;
 126   delete fPtGoodPi;
 127   delete fPtFakePi;
 128   delete fPtGoodK;
 129   delete fPtFakeK;
 130   delete fPtGoodP;
 131   delete fPtFakeP;
 132   delete fPtRPi;
 133   delete fPtRK;
 134   delete fPtRP;
 135
 136   delete fPtContP;
 137   delete fPtContPi;
 138   delete fPtContMup;
 139   delete fPtContElp;
 140
 141 }
 142
 143 AliFemtoCutMonitorParticlePtPDG& AliFemtoCutMonitorParticlePtPDG::operator=(const AliFemtoCutMonitorParticlePtPDG& aCut)
 144 {
 145   // assignment operator
 146   if (this == &aCut)
 147     return *this;
 148
 149   if (fPtPDG) delete fPtPDG;
 150   fPtPDG = new TH2D(*aCut.fPtPDG);
 151
 152   if (ftpcHist) delete ftpcHist;
 153   ftpcHist = new TH2D(*aCut.ftpcHist);
 154
 155    if (fPtGoodPi) delete fPtGoodPi;
 156   fPtGoodPi = new TH1D(*aCut.fPtGoodPi);
 157
 158    if (fPtFakePi) delete fPtFakePi;
 159   fPtFakePi = new TH1D(*aCut.fPtFakePi);
 160
 161    if (fPtRPi) delete fPtRPi;
 162   fPtRPi = new TH1D(*aCut.fPtRPi);
 163
 164   if (fPtGoodK) delete fPtGoodK;
 165   fPtGoodK = new TH1D(*aCut.fPtGoodK);
 166
 167    if (fPtFakeK) delete fPtFakeK;
 168   fPtFakeK = new TH1D(*aCut.fPtFakeK);
 169
 170    if (fPtRK) delete fPtRK;
 171   fPtRK = new TH1D(*aCut.fPtRK);
 172
 173   if (fPtGoodP) delete fPtGoodP;
 174   fPtGoodP = new TH1D(*aCut.fPtGoodP);
 175
 176    if (fPtFakeP) delete fPtFakeP;
 177   fPtFakeP = new TH1D(*aCut.fPtFakeP);
 178
 179    if (fPtRP) delete fPtRP;
 180   fPtRP = new TH1D(*aCut.fPtRP);
 181
 182    if (fPtContP) delete fPtContP;
 183   fPtContP = new TH1D(*aCut.fPtContP);
 184
 185    if (fPtContPi) delete fPtContPi;
 186   fPtContPi = new TH1D(*aCut.fPtContPi);
 187
 188      if (fPtContMup) delete fPtContMup;
 189   fPtContMup = new TH1D(*aCut.fPtContMup);
 190
 191      if (fPtContElp) delete fPtContElp;
 192   fPtContElp = new TH1D(*aCut.fPtContElp);
 193
 194   return *this;
 195 }
 196
 197 AliFemtoString AliFemtoCutMonitorParticlePtPDG::Report(){
 198   // Prepare report from the execution
 199   string stemp = "*** AliFemtoCutMonitorParticlePtPDG report";
 200   AliFemtoString returnThis = stemp;
 201   return returnThis;
 202 }
 203
 204 void AliFemtoCutMonitorParticlePtPDG::Fill(const AliFemtoTrack* aTrack)
 205 {
 206   // Fill in the monitor histograms with the values from the current track
 207   float tEnergy = ::sqrt(aTrack->P().mag2()+fMass*fMass);
 208   float tRapidity = 0.5*::log((tEnergy+aTrack->P().z())/(tEnergy-aTrack->P().z()));
 209   float tPt = ::sqrt((aTrack->P().x())*(aTrack->P().x())+(aTrack->P().y())*(aTrack->P().y()));
 210   float tEta = -TMath::Log(TMath::Tan(aTrack->P().theta()/2.0));
 211   float tPhi = aTrack->P().phi();
 212   float tP = ::sqrt((aTrack->P().z())*(aTrack->P().z())+(aTrack->P().x())*(aTrack->P().x())+(aTrack->P().y())*(aTrack->P().y()));;
 213   float dedx = aTrack->TPCsignalN();
 214   float w[10];
 215   w[0] = aTrack->PidProbElectron();
 216   w[1] = aTrack->PidProbMuon();
 217   w[2] = aTrack->PidProbPion();
 218   w[3] = aTrack->PidProbKaon();
 219   w[4] = aTrack->PidProbProton();
 220
 221    Int_t pdg1=0;
 222   AliFemtoModelHiddenInfo *info = ( AliFemtoModelHiddenInfo *) aTrack->GetHiddenInfo();
 223   if(info)pdg1 = info->GetPDGPid();
 224
 225
 226 //most probable particle
 227   fPtGoodPi->Fill(tPt);
 228   fPtGoodP->Fill(tPt);
 229   fPtGoodK->Fill(tPt);
 230
 231 //contaminations
 232   if (abs(pdg1)!=321)fPtFakeK->Fill(tPt);
 233   if (abs(pdg1)!=211)fPtFakePi->Fill(tPt);
 234   if (abs(pdg1)!=2212)fPtFakeP->Fill(tPt);
 235
 236
 237 //contaminations for kaons
 238   if (abs(pdg1)==2212)fPtContP->Fill(tPt);
 239   if (abs(pdg1)==211)fPtContPi->Fill(tPt);
 240   if (abs(pdg1)==13)fPtContMup->Fill(tPt);
 241   if (abs(pdg1)==11)fPtContElp->Fill(tPt);
 242
 243   Float_t pdg=-1.0;
 244   if(abs(pdg1)==211){
 245    pdg=2.0;
 246    fPtRPi->Fill(tPt);
 247   }
 248
 249   if(abs(pdg1)==321){
 250    pdg=3.0;
 251    fPtRK->Fill(tPt);
 252    }
 253
 254   if(abs(pdg1)==2212){
 255    pdg=4.0;
 256    fPtRP->Fill(tPt);
 257    }
 258
 259   if(abs(pdg1)==11)pdg=0.0; //+-electron
 260   if(abs(pdg1)==13)pdg=1.0;  //+-muon
 261
 262   //cout<<"pdg from CutMonitor.."<<pdg1<<"pdg"<<pdg<<endl;
 263
 264
 265   fPtPDG->Fill(pdg, tPt);
 266   ftpcHist->Fill(tP,dedx);
 267
 268
 269 }
 270
 271 void AliFemtoCutMonitorParticlePtPDG::Write()
 272 {
 273   // Write out the relevant histograms
 274
 275   fPtPDG->Write();
 276   ftpcHist->Write();
 277   fPtGoodPi->Write();
 278   fPtFakePi->Write();
 279   fPtGoodK->Write();
 280   fPtFakeK->Write();
 281   fPtGoodP->Write();
 282   fPtFakeP->Write();
 283   fPtRPi->Write();
 284   fPtRK->Write();
 285   fPtRP->Write();
 286   fPtContP->Write();
 287   fPtContPi->Write();
 288   fPtContMup->Write();
 289   fPtContElp->Write();
 290 }
 291
 292 TList *AliFemtoCutMonitorParticlePtPDG::GetOutputList()
 293 {
 294   TList *tOutputList = new TList();
 295   tOutputList->Add(fPtPDG);
 296   tOutputList->Add(ftpcHist);
 297   tOutputList->Add(fPtGoodPi);
 298   tOutputList->Add(fPtFakePi);
 299   tOutputList->Add(fPtGoodK);
 300   tOutputList->Add(fPtFakeK);
 301   tOutputList->Add(fPtGoodP);
 302   tOutputList->Add(fPtFakeP);
 303   tOutputList->Add(fPtRPi);
 304   tOutputList->Add(fPtRK);
 305   tOutputList->Add(fPtRP);
 306
 307   tOutputList->Add(fPtContP);
 308   tOutputList->Add(fPtContPi);
 309   tOutputList->Add(fPtContMup);
 310   tOutputList->Add(fPtContElp);
 311
 312   return tOutputList;
 313 }