StRoot: StFlowAnalysisMaker/plotLYZ.C Source File

00001 
00002 //
00003 // $Id: plotLYZ.C,v 1.4 2007/02/06 19:00:53 posk Exp $
00004 //
00005 // Plot histograms from flow.LeeYang.Zeros.root
00006 //
00007 // by Markus Oldenberg and Art Poskanzer
00008 //
00010 #include "TFile.h"
00011 #include "TCanvas.h"
00012 #include "TH1F.h"
00013 #include "TLine.h"
00014 #include "TMath.h"
00015 #include <iostream.h>
00016 #include <iomanip.h>
00017 
00018 void plotLYZ(TString fileName = "flow.hist.root", Bool_t firstPass = kFALSE, char* ext = "") {
00019 
00020   Bool_t mevSim, FTPC = kFALSE;
00021   //Bool_t mevSim = kTRUE;
00022   Bool_t FTPC = kTRUE;
00023 
00024   gROOT->SetStyle("Pub");                              // set style
00025   gStyle->SetFrameLineWidth(2);
00026   gStyle->SetLineWidth(2);
00027   gStyle->SetLineColor(kBlue);
00028   gStyle->SetHistLineWidth(1);
00029   gStyle->SetFuncWidth(2);
00030 //   gStyle->SetMarkerSize(2);
00031 //   gStyle->SetMarkerStyle(kFullCircle);
00032 //   gStyle->SetMarkerColor(kRed);
00033   gROOT->ForceStyle();
00034 
00035   TFile *file = new TFile(fileName, "READ");
00036 
00037   Int_t maxSel = 2;
00038   Int_t maxHar = 4; // 4
00039   Int_t maxHarPlot = 2;
00040   float ptMax = mevSim ? 2. : 6.;
00041   float v1 = 4.0; // 4. for dir9 mevSim
00042 
00043   Float_t max, min, r0;
00044 
00045   // Get Ntheta
00046   TString histNameNtheta("FlowLYZ_r0theta_Sel1_Har2");
00047   TH1D* histNtheta = new TH1D;
00048   histNtheta = (TH1D*)file->Get(histNameNtheta);
00049   if (!histNtheta) {
00050     cout << "### Can't find file " << histNameNtheta << endl;
00051     return;
00052   }
00053   int maxTheta = histNtheta->GetNbinsX();
00054   cout << "maxTheta = " << maxTheta << endl << endl;
00055   if (maxTheta > 5) maxTheta = 5;
00056   delete histNtheta;
00057 
00058   // Make the canvases
00059   TCanvas *canMult = new TCanvas("MeanMultiplicity", "MeanMultiplicity");
00060   TPaveLabel* run = new TPaveLabel(0.1,0.01,0.2,0.03,fileName.Data());  
00061   run->Draw();
00062   TDatime now;
00063   TPaveLabel* date = new TPaveLabel(0.7,0.01,0.9,0.03,now.AsString());
00064   date->Draw();
00065   TPad* graphPad = new TPad("Graphs","Graphs",0.01,0.05,0.97,0.95);
00066   graphPad->Draw();
00067 
00068   TCanvas *canGtheta = new TCanvas("Gtheta", "Gtheta", 700, 900);
00069   canGtheta->Divide(maxHarPlot*maxSel, maxTheta);
00070 
00071   TCanvas *canGthetaZoom = new TCanvas("Gtheta_Zoom", "Gtheta_Zoom", 700, 900);
00072   canGthetaZoom->Divide(maxHarPlot*maxSel, maxTheta);
00073 
00074   TCanvas *can_r0theta = new TCanvas("r0", "r0");
00075   can_r0theta->Divide(maxSel, maxHarPlot);
00076 
00077   TCanvas *can_vEta = new TCanvas("v(eta)", "v(eta)", 700, 900);
00078   can_vEta->Divide(maxSel, maxHar);
00079 
00080   TCanvas *can_vPt = new TCanvas("v(pt)", "v(pt)", 700, 900);
00081   can_vPt->Divide(maxSel, maxHar);
00082 
00083   TCanvas *can_vr0 = new TCanvas("vr0", "vr0");
00084   can_vr0->Divide(maxSel);
00085 
00086   TCanvas *can_v = new TCanvas("v", "v");
00087   can_v->Divide(maxSel);
00088 
00089   TCanvas *can_centX = new TCanvas("centX", "centX");
00090   can_centX->Divide(maxSel, maxHarPlot);
00091 
00092   TCanvas *can_centY = new TCanvas("centY", "centY");
00093   can_centY->Divide(maxSel, maxHarPlot);
00094 
00095   TCanvas *can_centQ = new TCanvas("centQ", "centQ");
00096   can_centQ->Divide(maxSel, maxHarPlot);
00097 
00098   // Make the histograms
00099   TH1D *histMult = new TH1D;
00100   TH1F **hist_v = new TH1F*[maxSel];
00101   TH1F **hist_vr0 = new TH1F*[maxSel];
00102   TH1F **hist_v = new TH1F*[maxSel];
00103   TH1F **hist_vEta = new TH1F*[maxSel*maxHar];
00104   TH1F **hist_vPt = new TH1F*[maxSel*maxHar];
00105   TH1D **hist_r0th = new TH1D*[maxSel*maxHarPlot];
00106   TH1F **histG = new TH1F*[maxSel*maxHarPlot*maxTheta];
00107   TH1F **histGZoom = new TH1F*[maxSel*maxHarPlot*maxTheta];
00108   TH1D **hist_centX = new TH1D*[maxSel*maxHarPlot];
00109   TH1D **hist_centY = new TH1D*[maxSel*maxHarPlot];
00110   TH1D **hist_centQ = new TH1D*[maxSel*maxHarPlot];
00111 
00112   TLine **r0Line = new TLine*[maxSel*maxHarPlot*maxTheta];
00113   TLine **r0LineZoom = new TLine*[maxSel*maxHarPlot*maxTheta]; 
00114   TLine *ptZeroLine = new TLine(0., 0., ptMax, 0.);
00115   TLine *pt5Line = new TLine(0., 5., 2., 5.);
00116   TLine *ptV1Line = new TLine(0., v1, 2., v1);
00117   TLine *etaZeroLine = new TLine(-1.5, 0., 1.5, 0.);
00118   TLine *etaZeroLineFTPC = new TLine(-4.5, 0., 4.5, 0.);
00119   TLine *GZeroLine = new TLine(0., 0., 0.35, 0.);
00120   TLine *eta5Line = new TLine(-4.5, 5., 4.5, 5.);
00121   TLine *etaV1LinePos = new TLine(0., v1, 4.5, v1);
00122   TLine *etaV1LineNeg = new TLine(-4.5, -v1, 0., -v1);
00123   TLine *vLine = new TLine(0.5, 0., maxHar+0.5, 0.);
00124   TLine *recentZeroLine = new TLine(0.5, 0., 3.5, 0.);
00125   TLine *recentedZeroLine = new TLine(0.5, 0., 2.5, 0.);
00126 
00127   // Multiplicity
00128   TString histName("FlowLYZ_Mult");
00129   cout << histName << endl;
00130   graphPad->cd();
00131   histMult = (TH1D*)file->Get(histName);
00132   if (!mevSim) { histMult->Fit("gaus"); }
00133   histMult->Draw();
00134   Double_t entries = histMult->GetEntries();
00135   TString* entriesChar = new TString("entries= ");
00136   *entriesChar += (int)entries;
00137   TLatex l;
00138   l.SetNDC();
00139   l.SetTextSize(0.05);
00140   l.DrawLatex(0.65,0.8,entriesChar->Data());
00141 
00142   float _v, vErr;  
00143   for (Int_t sel = 0; sel < maxSel; sel++) {
00144 
00145     if (!firstPass) {
00146       TString histName("FlowLYZ_v_Sel");
00147       histName += sel+1;
00148       cout << histName << endl;
00149       hist_v[sel] = (TH1F*)file->Get(histName);
00150       can_v->cd(sel+1); 
00151       hist_v[sel]->SetMinimum(0.);
00152       hist_v[sel]->Draw();
00153       vLine->Draw();
00154       for (int j=1; j<=maxHar; j++) {
00155         _v = hist_v[sel]->GetBinContent(j);
00156         vErr = hist_v[sel]->GetBinError(j);
00157         cout << setprecision(3) << "Sel = " << sel+1 << ": v" << j << " from pt = (" << _v <<
00158           " +/- " << vErr << ") %" << endl;
00159       }
00160     }
00161 
00162     TString histName("FlowLYZ_vr0_Sel");
00163     histName += sel+1;
00164     cout << histName << endl;
00165     hist_vr0[sel] = (TH1F*)file->Get(histName);
00166     can_vr0->cd(sel+1); 
00167     hist_vr0[sel]->SetMinimum(0.);
00168     hist_vr0[sel]->Draw();
00169     for (int j=1; j<=maxHar; j++) {
00170       _v = hist_vr0[sel]->GetBinContent(j);
00171       vErr = hist_vr0[sel]->GetBinError(j);
00172       cout << setprecision(3) << "Sel= " << sel+1 << ": v" << j << " from r0 = (" << _v <<
00173         " +/- " << vErr << ") %" << endl;
00174     }
00175 
00176     for (Int_t har = 0; har < maxHarPlot; har++) {
00177       int n = sel + har;
00178      
00179       TString histName("FlowLYZ_r0theta_Sel");
00180       histName += sel+1;
00181       histName += "_Har";
00182       histName += har+1;
00183       cout << histName << endl;
00184       hist_r0th[n] = (TH1D*)file->Get(histName);
00185       can_r0theta->cd(1+sel+har*maxSel);
00186       //hist_r0th[n]->Fit("pol0");
00187       hist_r0th[n]->SetMinimum(0.);
00188       hist_r0th[n]->Draw();
00189       for (int th=1; th<=maxTheta; th++) {
00190         _v = hist_r0th[n]->GetBinContent(th);
00191         vErr = hist_r0th[n]->GetBinError(th);
00192         if (TMath::IsNaN(vErr)) {
00193           vErr = 0.;
00194           hist_r0th[n]->SetBinError(th, 0.);
00195         }
00196         cout << setprecision(3) << "Sel=" << sel+1 << ", Har=" << har+1 <<": r0" << th << " = "
00197              << _v << " +/- " << vErr << endl;
00198       }
00199 
00200       if (!firstPass) {
00201         TString histName("FlowCentX_Sel");
00202         histName += sel+1;
00203         histName += "_Har";
00204         histName += har+1;
00205         hist_centX[n] = (TH1D*)file->Get(histName);
00206         if (hist_centX[n]) {
00207           cout << histName << endl;
00208           can_centX->cd(1+sel+har*maxSel);
00209           hist_centX[n]->Draw();
00210           recentZeroLine->Draw();
00211         }
00212         
00213         TString histName("FlowCentY_Sel");
00214         histName += sel+1;
00215         histName += "_Har";
00216         histName += har+1;
00217         hist_centY[n] = (TH1D*)file->Get(histName);
00218         if (hist_centY[n]) {
00219           cout << histName << endl;
00220           can_centY->cd(1+sel+har*maxSel);
00221           hist_centY[n]->Draw();
00222           recentZeroLine->Draw();
00223         }
00224       }
00225 
00226       TString histName("FlowQCent_Sel");
00227       histName += sel+1;
00228       histName += "_Har";
00229       histName += har+1;
00230       hist_centQ[n] = (TH1D*)file->Get(histName);
00231       if (hist_centQ[n]) {
00232         cout << histName << endl;
00233         can_centQ->cd(1+sel+har*maxSel);
00234         hist_centQ[n]->Draw();
00235         recentedZeroLine->Draw();
00236       }
00237 
00238       min = 0.00001;
00239       float expan = 1.2;
00240       for (Int_t theta = 0; theta < maxTheta; theta++) {
00241         TString histName("FlowLYZ_Gtheta");
00242         histName += theta;
00243         histName += "_Sel";
00244         histName += sel+1;
00245         histName += "_Har";
00246         histName += har+1;
00247         cout << histName << endl;
00248         histG[theta] = (TH1F*)file->Get(histName);
00249         canGtheta->cd(1+har+sel*maxHarPlot+theta*maxSel*maxHarPlot);
00250         TVirtualPad::Pad()->SetLogy();
00251         histG[theta]->SetMinimum(min);
00252         histG[theta]->DrawCopy("PH");
00253         r0 = hist_r0th[n]->GetBinContent(theta+1);
00254         r0Line[theta] = new TLine(r0, 0., r0, 1.);
00255         r0Line[theta]->Draw();
00256         
00257         min = 0.000001;
00258         max = 1.;
00259         histGZoom[theta] = (TH1F*)file->Get(histName);
00260         canGthetaZoom->cd(1+har+sel*maxHarPlot+theta*maxSel*maxHarPlot);
00261         TVirtualPad::Pad()->SetLogy();
00262         histGZoom[theta]->SetMaximum(max);
00263         histGZoom[theta]->SetMinimum(min);
00264         histGZoom[theta]->SetAxisRange(r0/expan, r0*expan, "X");
00265         histGZoom[theta]->Draw("PH");
00266         r0LineZoom[theta] = new TLine(r0, 0., r0, max);
00267         r0LineZoom[theta]->Draw();
00268 
00269       }
00270     }
00271 
00272     if (!firstPass) {
00273       for (Int_t har = 0; har < maxHar; har++) {
00274         int n = sel + har;
00275 
00276         TString histName("FlowLYZ_vEta_Sel");
00277         histName += sel+1;
00278         histName += "_Har";
00279         histName += har+1;
00280         cout << histName << endl;
00281         hist_vEta[n] = (TH1D*)file->Get(histName);
00282         can_vEta->cd(sel+1+har*maxSel);
00283         hist_vEta[n]->SetMaximum(10.);
00284         hist_vEta[n]->SetMinimum(-10.);
00285         hist_vEta[n]->Draw("E");
00286         if (FTPC) { etaZeroLineFTPC->Draw(); }
00287         else { etaZeroLine->Draw(); }
00288         if (mevSim) {
00289           if (har==1) {
00290             eta5Line->Draw();
00291           } else if (har==0) {
00292             etaV1LinePos->Draw();
00293             etaV1LineNeg->Draw();
00294           }
00295         }
00296         
00297         TString histName("FlowLYZ_vPt_Sel");
00298         histName += sel+1;
00299         histName += "_Har";
00300         histName += har+1;
00301         cout << histName << endl;
00302         hist_vPt[n] = (TH1D*)file->Get(histName);
00303         can_vPt->cd(sel+1+har*maxSel);
00304         if (mevSim) {
00305           hist_vPt[n]->SetMaximum(10.);
00306           hist_vPt[n]->SetMinimum(-10.);
00307         } else {
00308           hist_vPt[n]->SetMaximum(20.);
00309           hist_vPt[n]->SetMinimum(-20.);
00310         }
00311         hist_vPt[n]->Draw("E");
00312         ptZeroLine->Draw();
00313         if (mevSim) {
00314           if (har==1) {
00315             pt5Line->Draw();
00316           } else if (har==0) {
00317             ptV1Line->Draw();
00318           }
00319         }
00320       }
00321     }
00322   }
00323 
00324   if (strstr(ext,"ps")) {
00325     canGtheta->SaveAs("FlowLYZ_Gtheta.ps");
00326     canGthetaZoom->SaveAs("FlowLYZ_GthetaZoom.ps");
00327     can_r0theta->SaveAs("FlowLYZ_r0.ps");
00328     canMult->SaveAs("FlowLYZ_Mult.ps");
00329     can_vEta->SaveAs("FlowLYZ_vEta.ps");
00330     can_vPt->SaveAs("FlowLYZ_vPt.ps");
00331     can_vr0->SaveAs("FlowLYZ_vro.ps");
00332     can_v->SaveAs("FlowLYZ_v.ps");
00333     can_centX->SaveAs("FlowCentX.ps");
00334     can_centY->SaveAs("FlowCentY.ps");
00335     can_centQ->SaveAs("FlowCentQ.ps");
00336   } else if (strstr(ext,"gif")) {
00337     canGtheta->SaveAs("FlowLYZ_Gtheta.gif");
00338     canGthetaZoom->SaveAs("FlowLYZ_GthetaZoom.gif");
00339     can_r0theta->SaveAs("FlowLYZ_r0.gif");
00340     canMult->SaveAs("FlowLYZ_Mult.gif");
00341     can_vEta->SaveAs("FlowLYZ_vEta.gif");
00342     can_vPt->SaveAs("FlowLYZ_vPt.gif");
00343     can_vr0->SaveAs("FlowLYZ_vro.gif");
00344     can_v->SaveAs("FlowLYZ_v.gif");
00345   }
00346 
00347   return;
00348 }
00349 
00351 //
00352 // $Log: plotLYZ.C,v $
00353 // Revision 1.4  2007/02/06 19:00:53  posk
00354 // In Lee Yang Zeros method, introduced recentering of Q vector.
00355 // Reactivated eta symmetry cut.
00356 //
00357 // Revision 1.2  2006/03/22 22:02:14  posk
00358 // Updates to macros.
00359 //
00360 //