PWGLF/STRANGENESS/LambdaK0PbPb/MultYields2.C

   1
   2 #include "AliMassFitControl.h"
   3
   4 // Two functions now.
   5 // MultYields3 takes a 3D histogram and calls the 2D function MultYields2
   6 // MultYields2 could be called directly if we only have a 2D histogram
   7
   8 void MultYields3(TH3F *PtMassMult, Int_t particleMode, Int_t MultBin, Char_t* label){
   9
  10   // Leave open possibility to choose different values depending on MultBin - LSB
  11   Float_t MultLo[1] = {0};
  12   Float_t MultHi[1] = {300};
  13
  14   //Make 2D projections from the 3D histogram
  15   PtMassMult->GetZaxis()->SetRange(MultLo[MultBin],MultHi[MultBin]);
  16   TH2F* hPtMass = (TH2F*)PtMassMult->Project3D("XY");// FIX:MF
  17   hPtMass->SetTitle("PtMass");
  18
  19   //hPtMass->Draw(); // Drawing invokes default c1 canvas making it inaccessible in MultYields2 - TO FIX
  20   MultYields2(hPtMass, particleMode,MultBin,label);
  21
  22 }
  23
  24 void MultYields2(TH2F *hPtMass, Int_t particleMode, Int_t MultBin, Char_t* label){
  25
  26   hPtMass->Draw();
  27
  28     // Do .L on command line first
  29     //gROOT->LoadMacro("macros/PtMassAna2.C");
  30 /* Modifications to produce a single uncorrected spectrum in a particular mult. bin for K0 or Lambda
  31 Old ratio code preserved in MultYieldsRatio.C
  32 Dec 2003 */
  33   Char_t* part; //for name of particle
  34   Float_t BR; //branching ratio - check them
  35   if (particleMode==0){
  36     part = "K0";
  37       BR=0.686;
  38   } else if (particleMode==1){
  39     part = "Lambda";
  40     BR=0.639;
  41   } else if (particleMode==2){
  42       part = "AntiLambda";
  43       BR=0.639;
  44   } else if (particleMode==3){
  45     part = "Lambda+antiLambda";
  46     BR=0.639;
  47   } else if (particleMode==4){
  48     part = "Xi";
  49     BR=1.0; // Should be Lam ->p pi Br
  50   } else if (particleMode ==6){
  51     part = "Omega";
  52     BR=1.0; // Should be Lam->p pi  * Om -> K Lam
  53   }
  54
  55
  56
  57
  58 TString title[1]={"MinimumBias"};
  59   //Make 2D projections from the 3D histogram
  60   //Minbias (i.e. everything)
  61
  62   TObjArray *controllerArray = new TObjArray(40,1); // 2nd arg means can count from 1!
  63
  64   //Here probably need switch-case, depending on mult bin and particle
  65   // LoPt, HiPt, polynomial order, rebinning factor
  66
  67   /////  LAMBDA and LAMBDA+ANTILAMBDA (Combination)
  68
  69   if(particleMode == 1 || particleMode == 3){ // Lambda or Lambda+Anti-Lambda
  70   //  controllerArray->AddLast(new AliMassFitControl(0.1,0.2, 1,2));
  71   //  controllerArray->AddLast(new AliMassFitControl(0.2,0.3, 1,2));
  72   //  controllerArray->AddLast(new AliMassFitControl(0.3,0.4, 1,1, 1.096,1.16))
  73   //  controllerArray->AddLast(new AliMassFitControl(0.4,0.5, 2,1, 1.095,1.15));
  74   //  controllerArray->AddLast(new AliMassFitControl(0.5,0.6, 2,1, 1.095,1.15));
  75   //  controllerArray->AddLast(new AliMassFitControl(0.6,0.7, 2,1, 1.095,1.15));
  76   //  controllerArray->AddLast(new AliMassFitControl(0.7,0.8, 2,1, 1.095,1.15));
  77       controllerArray->AddLast(new AliMassFitControl(0.8,1.0, 2,1, 1.095,1.15));
  78       controllerArray->AddLast(new AliMassFitControl(1.0,1.2, 2,1, 1.095,1.17));
  79       controllerArray->AddLast(new AliMassFitControl(1.2,1.4, 2,1));
  80       controllerArray->AddLast(new AliMassFitControl(1.4,1.6, 2,1));
  81       controllerArray->AddLast(new AliMassFitControl(1.6,1.8, 2,1));
  82       controllerArray->AddLast(new AliMassFitControl(1.8,2.0, 2,1));
  83       controllerArray->AddLast(new AliMassFitControl(2.0,2.2, 2,1, 1.095,1.17));
  84       controllerArray->AddLast(new AliMassFitControl(2.2,2.4, 2,1, 1.095,1.16));
  85       controllerArray->AddLast(new AliMassFitControl(2.4,2.6, 2,1, 1.095,1.17));
  86       controllerArray->AddLast(new AliMassFitControl(2.6,2.8, 1,2, 1.095,1.17));
  87       controllerArray->AddLast(new AliMassFitControl(2.8,3.0, 1,2, 1.095,1.17));
  88       controllerArray->AddLast(new AliMassFitControl(3.0,3.5, 1,2, 1.098,1.16));
  89       controllerArray->AddLast(new AliMassFitControl(3.5,4.0, 1,2, 1.095,1.17));
  90       controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 1,2, 1.095,1.17));
  91       controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 1,2, 1.095,1.165));
  92       controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 1,2, 1.095,1.17));
  93  //     controllerArray->AddLast(new AliMassFitControl(5.5,6.0, 1,2, 1.095,1.17));
  94  //   controllerArray->AddLast(new AliMassFitControl(6.0,6.5, 1,2, 1.095,1.17));
  95  //   controllerArray->AddLast(new AliMassFitControl(6.5,7.0, 1,2, 1.095,1.17));
  96  //   controllerArray->AddLast(new AliMassFitControl(7.0,7.5, 1,2, 1.095,1.17));
  97  //   controllerArray->AddLast(new AliMassFitControl(7.5,8.0, 1,2, 1.095,1.17));
  98  //   controllerArray->AddLast(new AliMassFitControl(8.0,9.0, 1,2, 1.095,1.17));
  99  //   controllerArray->AddLast(new AliMassFitControl(9.0,10.0, 1,2, 1.090,1.17));
 100  //   controllerArray->AddLast(new AliMassFitControl(10.0,12.0, 1,2, 1.082,1.165));
 101  }
 102   /// ANTI LAMBDA ---->
 103   else if (particleMode == 2){ // Anti-Lambdas
 104   controllerArray->AddLast(new AliMassFitControl(0.3,0.4, 2,1, 1.095,1.15));
 105   controllerArray->AddLast(new AliMassFitControl(0.4,0.5, 2,1, 1.095,1.15));
 106   controllerArray->AddLast(new AliMassFitControl(0.5,0.6, 2,1, 1.095,1.15));
 107   controllerArray->AddLast(new AliMassFitControl(0.6,0.7, 2,1, 1.095,1.15));
 108   controllerArray->AddLast(new AliMassFitControl(0.7,0.8, 2,1, 1.090,1.15));
 109   controllerArray->AddLast(new AliMassFitControl(0.8,1.0, 2,1, 1.095,1.15));
 110   controllerArray->AddLast(new AliMassFitControl(1.0,1.2, 2,1, 1.090,1.17));
 111   controllerArray->AddLast(new AliMassFitControl(1.2,1.4, 2,1, 1.090,1.16));
 112   controllerArray->AddLast(new AliMassFitControl(1.4,1.6, 2,1, 1.090,1.17));
 113   controllerArray->AddLast(new AliMassFitControl(1.6,1.8, 2,1, 1.090,1.17));
 114   controllerArray->AddLast(new AliMassFitControl(1.8,2.0, 2,1, 1.090,1.17));
 115   controllerArray->AddLast(new AliMassFitControl(2.0,2.2, 2,1, 1.095,1.17));
 116   controllerArray->AddLast(new AliMassFitControl(2.2,2.4, 2,1, 1.095,1.16));
 117   controllerArray->AddLast(new AliMassFitControl(2.4,2.6, 1,1, 1.095,1.17));
 118   controllerArray->AddLast(new AliMassFitControl(2.6,2.8, 1,2, 1.095,1.17));
 119   controllerArray->AddLast(new AliMassFitControl(2.8,3.0, 1,2, 1.095,1.17));
 120   controllerArray->AddLast(new AliMassFitControl(3.0,3.5, 1,2, 1.095,1.17));
 121   controllerArray->AddLast(new AliMassFitControl(3.5,4.0, 1,2, 1.095,1.17));
 122   controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 1,2, 1.090,1.17));
 123   controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 0,4, 1.095,1.17));
 124   controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 0,4, 1.095,1.17));
 125   controllerArray->AddLast(new AliMassFitControl(5.5,6.0, 0,6, 1.085,1.17));
 126   } // end if anti-Lambda
 127   else if (particleMode == 0){ // K0s case
 128       //controllerArray->AddLast(new AliMassFitControl(0.1,0.2, 1,2, 0.44,0.56));
 129       //controllerArray->AddLast(new AliMassFitControl(0.2,0.3, 1,2, 0.44,0.56));
 130       //controllerArray->AddLast(new AliMassFitControl(0.3,0.4, 1,1, 0.442,0.559));
 131       //controllerArray->AddLast(new AliMassFitControl(0.4,0.5, 1,1, 0.44,0.56));
 132       //controllerArray->AddLast(new AliMassFitControl(0.5,0.6, 1,1, 0.442,0.559));
 133       //controllerArray->AddLast(new AliMassFitControl(0.6,0.7, 1,1, 0.442,0.558));
 134       //controllerArray->AddLast(new AliMassFitControl(0.7,0.8, 1,1, 0.44,0.56));
 135       controllerArray->AddLast(new AliMassFitControl(0.8,1.0, 2,1, 0.44,0.56));
 136       controllerArray->AddLast(new AliMassFitControl(1.0,1.2, 2,1, 0.45,0.56));
 137       controllerArray->AddLast(new AliMassFitControl(1.2,1.4, 2,1, 0.44,0.56));
 138       controllerArray->AddLast(new AliMassFitControl(1.4,1.6, 2,1, 0.45,0.56));
 139       controllerArray->AddLast(new AliMassFitControl(1.6,1.8, 2,1, 0.44,0.56));
 140       controllerArray->AddLast(new AliMassFitControl(1.8,2.0, 1,1, 0.44,0.56));
 141       controllerArray->AddLast(new AliMassFitControl(2.0,2.2, 1,1, 0.44,0.56));
 142       controllerArray->AddLast(new AliMassFitControl(2.2,2.4, 1,1, 0.44,0.56));
 143       controllerArray->AddLast(new AliMassFitControl(2.4,2.6, 1,1, 0.44,0.56));
 144       controllerArray->AddLast(new AliMassFitControl(2.6,2.8, 1,2, 0.44,0.54));
 145       controllerArray->AddLast(new AliMassFitControl(2.8,3.0, 1,2, 0.44,0.56));
 146       controllerArray->AddLast(new AliMassFitControl(3.0,3.5, 1,2, 0.44,0.56));
 147       controllerArray->AddLast(new AliMassFitControl(3.5,4.0, 1,2, 0.46,0.56));
 148       controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 1,2, 0.46,0.56));
 149       controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 1,2, 0.418,0.52));
 150       //controllerArray->AddLast(new AliMassFitControl(5.0,5.5, 0,4, 0.418,0.53));
 151 //      controllerArray->AddLast(new AliMassFitControl(5.5,6.0, 1,4, 0.418,0.54));
 152 //      controllerArray->AddLast(new AliMassFitControl(6.0,6.5, 1,4, 0.418,0.56));
 153       //controllerArray->AddLast(new AliMassFitControl(6.5,7.0, 0,4, 0.418,0.578));
 154   }  else if (particleMode == 4) { //Xi case
 155     //controllerArray->AddLast(new AliMassFitControl(0.5,0.7, 1,1, 1.28,1.45)); //signal not visible with
 156     controllerArray->AddLast(new AliMassFitControl(0.7,0.9, 1,1, 1.285,1.45));
 157     controllerArray->AddLast(new AliMassFitControl(0.9,1.1, 1,1, 1.285,1.45));
 158     controllerArray->AddLast(new AliMassFitControl(1.1,1.3, 1,1, 1.27,1.45));
 159     controllerArray->AddLast(new AliMassFitControl(1.3,1.5, 1,1, 1.27,1.45));
 160     controllerArray->AddLast(new AliMassFitControl(1.5,1.7, 1,1, 1.27,1.45));
 161     controllerArray->AddLast(new AliMassFitControl(1.7,2.0, 1,1, 1.27,1.45));
 162     controllerArray->AddLast(new AliMassFitControl(2.0,2.5, 1,1, 1.27,1.44));
 163     controllerArray->AddLast(new AliMassFitControl(2.5,3.0, 1,1, 1.27,1.45));
 164     controllerArray->AddLast(new AliMassFitControl(3.0,3.5, 1,1, 1.27,1.44));
 165     controllerArray->AddLast(new AliMassFitControl(3.5,4.0, 1,1, 1.27,1.45));
 166     controllerArray->AddLast(new AliMassFitControl(4.0,4.5, 2,1, 1.27,1.44));
 167     controllerArray->AddLast(new AliMassFitControl(4.5,5.0, 2,1, 1.27,1.45));
 168   }
 169   controllerArray->Sort();
 170
 171   // Make the proper label to pass in for use in labelling the diagnostics
 172   // saved canvas files and histos
 173   Char_t fulllabel[80];
 174   sprintf(fulllabel,"%s%s",title[MultBin].Data(),label);
 175
 176   //Slice up projection into mass histograms to extract yield
 177 TH1F* hYield = (TH1F*)PtMassAna2(hPtMass,particleMode,controllerArray->GetEntries(),controllerArray,fulllabel);
 178
 179 // CORRECTIONS
 180   Int_t Nev = 2800;  // Obviously we need a better way to access the proper number FIX ME
 181  hYield->Scale(1.0/Nev); //Divide by the number of events
 182 // Better not to hide these other corrections deep in here
 183   //hYield->Scale(Veff[MultBin]); //Multiply by the vertex efficiency effectively increaing number of events
 184                               //(since Veff<1) therefore decreases yield
 185 //hYield->Scale(1.0/BR);  //Divide by branching ratio (again increases yield since BR<1)
 186  //hYield->Scale(1.0/(2*TMath::Pi())); // Always plot 1/2pi ...
 187
 188 Char_t yieldTitle[80];
 189 sprintf(yieldTitle,"Uncorrected %s yield: %s",part,title[MultBin].Data());
 190 hYield->SetTitle(yieldTitle);
 191 hYield->SetXTitle("p_{t} / [GeV/c]");
 192 hYield->SetYTitle("1/Nev.dN/dp_{t}");
 193
 194 // Create plots
 195
 196 Char_t fileNameBase[80];
 197 sprintf(fileNameBase,"Masses%s_%s%s",part,title[MultBin].Data(),label);
 198 Char_t fileNamePng[80];
 199 sprintf(fileNamePng,"%s.png",fileNameBase);
 200 Char_t fileNameEps[80];
 201 sprintf(fileNameEps,"%s.eps",fileNameBase);
 202 Char_t fileNamePdf[80];
 203 sprintf(fileNamePdf,"%s.pdf",fileNameBase);
 204
 205 c1->SaveAs(fileNamePng);
 206 c1->SaveAs(fileNameEps);
 207 c1->SaveAs(fileNamePdf);
 208
 209 //c1->Clear();
 210
 211 //c1->SetLogy();
 212 TCanvas *cYield = new TCanvas("Yield","Corrected Yield",600,400);
 213 cYield->cd();
 214 //cYield->SetLogy();
 215 hYield->SetStats(kFALSE);
 216 hYield->Draw();
 217  cYield->SetLogy();
 218  cYield->Update();
 219 //  hRC_MB->SetMarkerStyle(20);
 220 //  hRC_MB->SetMarkerColor(4);
 221 //  hRC_MB->Scale(NBinMB/NBin3);
 222 Char_t fnametext[80];
 223 sprintf(fnametext,"Yield%s_%s%s",part,title[MultBin].Data(),label);
 224 Char_t fnamePng[80];
 225 sprintf(fnamePng,"%s.png",fnametext);
 226 c1->SaveAs(fnamePng);
 227 Char_t fnameEps[80];
 228 sprintf(fnameEps,"%s.eps",fnametext);
 229 c1->SaveAs(fnameEps);
 230
 231 // This section for yield scaled by number of binary collisions.
 232 // Could add array of values and do scaling according to 'MultBin' index
 233 TH1F* hScYield = hYield->Clone("ScYield");
 234 Char_t scYieldTitle[80];
 235 sprintf(scYieldTitle,"<N_{bin}> Scaled %s",hYield->GetTitle());
 236 hScYield->SetTitle(scYieldTitle);
 237 //SCALING for scaled yield only  divide by mean Nbin (scaled yield is therefore smaller)
 238 //hScYield->Scale(1/NBin[MultBin]);
 239
 240 Char_t fnameRoot[80];
 241 sprintf(fnameRoot,"%s.root",fnametext);
 242 TFile *YieldFile = new TFile(fnameRoot,"RECREATE");
 243 hYield->Write();
 244 hScYield->Write();
 245 YieldFile->Close();
 246
 247  controllerArray->Delete();
 248 }