PWGCF/FLOW/Documentation/examples/plotGlauberCenVars.C

   1 ///==========================================================================
   2 ///
   3 ///    macro to plot centrality bin values
   4 ///==========================================================================
   5 ///
   6 #include <cstdlib>
   7
   8 const Int_t nbins=11;
   9 double centPercent[nbins]={5.,10.,20.,30.,40.,50.,60.,70.,80.,90.,100};
  10 TArrayI* binUp = new TArrayI(nbins);
  11 TArrayI* Multbin = new TArrayI(nbins);
  12
  13 void plotGlauberCenVars(const Char_t* file="glau_pbpb_ntuple.root")
  14 {
  15  TFile f(file);
  16  TNtuple* ntuple = f.Get("nt_Pb_Pb");
  17
  18  TCanvas* canvas = new TCanvas("Multiplicity","Multiplicity");
  19  canvas->SetLogy();
  20  ntuple->Draw("dNdEtaNBD>>htemp(500)");
  21  TH1F* htemp = (TH1F*)gPad->GetPrimitive("htemp");
  22  TH1F* hist = (TH1F*)htemp->Clone();
  23  hist->SetTitle("");
  24  hist->SetStats(0);
  25  TAxis *axis = hist->GetXaxis();
  26  axis->SetTitle("Multiplicity");
  27  //axis->CenterTitle(kTRUE);
  28  axis->SetTitleSize(0.05);
  29  axis->SetLabelSize(0.04);
  30  axis->SetTitleOffset(1.2);
  31  axis = hist->GetYaxis();
  32  axis->SetTitle("");
  33  axis->SetTitleOffset(1.3);
  34  //axis->CenterTitle(kTRUE);
  35  axis->SetTitleSize(0.05);
  36  axis->SetLabelSize(0.04);
  37  hist->DrawCopy();
  38
  39  //---------------------------------------------------
  40  getCentrality(hist);
  41  //---------------------------------------------------
  42
  43  int low=0;
  44  int high=hist->GetNbinsX();
  45  //printf("size: %i\n",binUp->GetSize());
  46  for (int i=0; i<binUp->GetSize(); i++)
  47  {
  48    low=binUp->At(i);
  49    //printf("%i, %i\n",low,high);
  50    hist->GetXaxis()->SetRange(low+1,high);
  51    hist->SetFillColor((i%2==0)?0:kGray);
  52    hist->DrawCopy("same");
  53    high=low;
  54  }
  55
  56  TCanvas* canvas2 = new TCanvas("N Wounded","N Wounded");
  57  canvas2->SetLogy();
  58  ntuple->Draw("Npart>>htemp2(200)");
  59  TH1F* htemp2 = (TH1F*)gPad->GetPrimitive("htemp2");
  60  TH1F* hist2 = (TH1F*)htemp2->Clone();
  61  hist2->SetStats(0);
  62  hist2->SetTitle("");
  63  axis = hist2->GetXaxis();
  64  axis->SetTitle("N Wounded");
  65  axis->SetTitleSize(0.05);
  66  axis->SetLabelSize(0.04);
  67  axis->SetTitleOffset(1.2);
  68  axis = hist2->GetYaxis();
  69  axis->SetTitle("");
  70  axis->SetTitleOffset(1.3);
  71  axis->SetTitleSize(0.05);
  72  axis->SetLabelSize(0.04);
  73  hist2->DrawCopy();
  74
  75
  76  int low=0;
  77  int high=5000;
  78  //printf("size: %i\n",Multbin->GetSize());
  79  TH1F* htemp10 = new TH1F();
  80  for (int i=0; i<Multbin->GetSize(); i++)
  81  {
  82    low=Multbin->At(i);
  83    //printf("%i, %i\n",low,high);
  84    char cuts[100];
  85    char histtitle1[100];
  86    char histtitle[100];
  87    sprintf(cuts,"dNdEtaNBD>%i && dNdEtaNBD<%i",low,high);
  88    //printf("dNdEtaNBD>%i && dNdEtaNBD<%i\n",low,high);
  89    sprintf(histtitle,"Npart>>htemp%i(200)",20+i);
  90    sprintf(histtitle1,"htemp%i",20+i);
  91    ntuple->Draw(histtitle,cuts,"same");
  92    htemp10 = (TH1F*)gPad->GetPrimitive(histtitle1);
  93    htemp10->SetLineColor(i+1);
  94    htemp10->DrawCopy("same");
  95    high=low;
  96  }
  97
  98  TCanvas* canvas3 = new TCanvas("N Binary","N Binary");
  99  canvas3->SetLogy();
 100  ntuple->Draw("Ncoll>>htemp3(500)");
 101  TH1F* htemp3 = (TH1F*)gPad->GetPrimitive("htemp3");
 102  TH1F* hist3 = (TH1F*)htemp3->Clone();
 103  hist3->SetStats(0);
 104  hist3->SetTitle("");
 105  axis = hist3->GetXaxis();
 106  axis->SetTitle("N Binary");
 107  axis->SetTitleSize(0.05);
 108  axis->SetLabelSize(0.04);
 109  axis->SetTitleOffset(1.2);
 110  axis = hist3->GetYaxis();
 111  axis->SetTitle("");
 112  axis->SetTitleOffset(1.3);
 113  axis->SetTitleSize(0.05);
 114  axis->SetLabelSize(0.04);
 115  hist3->DrawCopy();
 116
 117  int low=0;
 118  int high=5000;
 119  for (int i=0; i<Multbin->GetSize(); i++)
 120  {
 121    low=Multbin->At(i);
 122    char cuts[100];
 123    char histtitle[100];
 124    char histtitle1[100];
 125    sprintf(cuts,"dNdEtaNBD>%i && dNdEtaNBD<%i",low,high);
 126    sprintf(histtitle,"Ncoll>>htemp%i(500)",40+i);
 127    sprintf(histtitle1,"htemp%i",40+i);
 128    ntuple->Draw(histtitle,cuts,"same");
 129    TH1F* htemp11 = (TH1F*)gPad->GetPrimitive(histtitle1);
 130    htemp11->SetLineColor(i+1);
 131    htemp11->DrawCopy("same");
 132    high=low;
 133  }
 134
 135  TCanvas* canvas4 = new TCanvas("Impact Parameter","Impact Parameter");
 136  ntuple->Draw("B>>htemp4(100)");
 137  TH1F* htemp4 = (TH1F*)gPad->GetPrimitive("htemp4");
 138  TH1F* hist4 = (TH1F*)htemp4->Clone();
 139  hist4->SetStats(0);
 140  hist4->SetTitle("");
 141  axis = hist4->GetXaxis();
 142  axis->SetTitle("b");
 143  axis->SetTitleSize(0.05);
 144  axis->SetLabelSize(0.04);
 145  axis->SetTitleOffset(1.2);
 146  axis = hist4->GetYaxis();
 147  axis->SetTitle("");
 148  axis->SetTitleOffset(1.3);
 149  axis->SetTitleSize(0.05);
 150  axis->SetLabelSize(0.04);
 151  hist4->DrawCopy();
 152
 153  int low=0;
 154  int high=5000;
 155  for (int i=0; i<Multbin->GetSize(); i++)
 156  {
 157    low=Multbin->At(i);
 158    char cuts[100];
 159    char histtitle[100];
 160    char histtitle1[100];
 161    sprintf(cuts,"dNdEtaNBD>%i && dNdEtaNBD<%i",low,high);
 162    sprintf(histtitle,"B>>htemp%i(100)",50+i);
 163    sprintf(histtitle1,"htemp%i",50+i);
 164    ntuple->Draw(histtitle,cuts,"same");
 165    TH1F* htemp12 = (TH1F*)gPad->GetPrimitive(histtitle1);
 166    htemp12->SetLineColor(i+1);
 167    htemp12->DrawCopy("same");
 168    high=low;
 169  }
 170
 171  TCanvas* canvas5 = new TCanvas("Eccentricity","Eccentricity");
 172  canvas5->SetLogy();
 173  ntuple->Draw("VarE>>htemp5(100)");
 174  TH1F* htemp5 = (TH1F*)gPad->GetPrimitive("htemp5");
 175  TH1F* hist5 = (TH1F*)htemp5->Clone();
 176  hist5->SetStats(0);
 177  hist5->SetTitle("");
 178  axis = hist5->GetXaxis();
 179  axis->SetTitle("standard eccentricity");
 180  axis->SetTitleSize(0.05);
 181  axis->SetLabelSize(0.04);
 182  axis->SetTitleOffset(1.2);
 183  axis = hist5->GetYaxis();
 184  axis->SetTitle("");
 185  axis->SetTitleOffset(1.3);
 186  axis->SetTitleSize(0.05);
 187  axis->SetLabelSize(0.04);
 188  hist5->DrawCopy();
 189
 190  int low=0;
 191  int high=5000;
 192  for (int i=0; i<Multbin->GetSize(); i++)
 193  {
 194    low=Multbin->At(i);
 195    char cuts[100];
 196    char histtitle[100];
 197    char histtitle1[100];
 198    sprintf(cuts,"dNdEtaNBD>%i && dNdEtaNBD<%i",low,high);
 199    sprintf(histtitle,"VarE>>htemp%i(100)",60+i);
 200    sprintf(histtitle1,"htemp%i",60+i);
 201    ntuple->Draw(histtitle,cuts,"same");
 202    TH1F* htemp13 = (TH1F*)gPad->GetPrimitive(histtitle1);
 203    htemp13->SetLineColor(i+1);
 204    htemp13->DrawCopy("same");
 205    high=low;
 206  }
 207
 208  TCanvas* canvas6 = new TCanvas("EccentricityPart","EccentricityPart");
 209  canvas6->SetLogy();
 210  ntuple->Draw("VarEPart>>htemp6(100)");
 211  TH1F* htemp6 = (TH1F*)gPad->GetPrimitive("htemp6");
 212  TH1F* hist6 = (TH1F*)htemp6->Clone();
 213  hist6->SetStats(0);
 214  hist6->SetTitle("");
 215  axis = hist6->GetXaxis();
 216  axis->SetTitle("participant eccentricity");
 217  axis->SetTitleSize(0.05);
 218  axis->SetLabelSize(0.04);
 219  axis->SetTitleOffset(1.2);
 220  axis = hist6->GetYaxis();
 221  axis->SetTitle("");
 222  axis->SetTitleOffset(1.3);
 223  axis->SetTitleSize(0.05);
 224  axis->SetLabelSize(0.04);
 225  hist6->DrawCopy();
 226
 227  int low=0;
 228  int high=5000;
 229  for (int i=0; i<Multbin->GetSize(); i++)
 230  {
 231    low=Multbin->At(i);
 232    char cuts[100];
 233    char histtitle[100];
 234    char histtitle1[100];
 235    sprintf(cuts,"dNdEtaNBD>%i && dNdEtaNBD<%i",low,high);
 236    sprintf(histtitle,"VarEPart>>htemp%i(100)",70+i);
 237    sprintf(histtitle1,"htemp%i",70+i);
 238    ntuple->Draw(histtitle,cuts,"same");
 239    TH1F* htemp14 = (TH1F*)gPad->GetPrimitive(histtitle1);
 240    htemp14->SetLineColor(i+1);
 241    htemp14->DrawCopy("same");
 242    high=low;
 243  }
 244
 245 }
 246
 247
 248 void getCentrality(TH1 *histNch, Float_t ff=1.0)
 249 // histNch - histo of multiplicity distribution (better with binsize=1)
 250 // ff fraction of accepted events. All losses are assumed to occur in most
 251 // peripheral bin
 252 {
 253
 254  //double sum= histNch->GetEntries() - histNch->GetBinContent(1);
 255  double sum= histNch->Integral();
 256  int nbinsx=histNch->GetNbinsX();
 257  double frac=0.;
 258  int ic=0;
 259  for (int ib=nbinsx;ib>0;ib--){
 260    frac += histNch->GetBinContent(ib)/sum*100.*ff;
 261    if(frac > centPercent[ic]){
 262      binUp->SetAt(ib,ic);
 263      Multbin->SetAt(histNch->GetBinCenter(ib),ic);
 264      cout<<" centrality="<<centPercent[ic]<<"   mult <="<< histNch->GetBinCenter(ib) <<endl;
 265      ic++;
 266    }
 267    if(ic==nbins) break;
 268  }
 269  printf(" \n float multCent[%i] = {",nbins);
 270  // cout <<" \n float multCent[nbins] = {";
 271
 272  for (int ic=nbins-1; ic>-1; ic--){
 273    cout<< histNch->GetBinCenter(binUp->At(ic));
 274    if (ic!=0) cout<<", ";
 275  }
 276  cout<<"};\n"<<endl;
 277 }