PWG1/AliComparisonRes.cxx

   1 //------------------------------------------------------------------------------\r
   2 // Implementation of AliComparisonRes class. It keeps information from \r
   3 // comparison of reconstructed and MC particle tracks. In addtion, \r
   4 // it keeps selection cuts used during comparison. The comparison \r
   5 // information is stored in the ROOT histograms. Analysis of these \r
   6 // histograms can be done by using Analyse() class function. The result of \r
   7 // the analysis (histograms) are stored in the output picture_res.root file.\r
   8 //  \r
   9 // Author: J.Otwinowski 04/02/2008 \r
  10 //------------------------------------------------------------------------------\r
  11 \r
  12 /*\r
  13   //after running analysis, read the file, and get component\r
  14   gSystem->Load("libPWG1.so");\r
  15   TFile f("Output.root");\r
  16   AliComparisonRes * comp = (AliComparisonRes*)f.Get("AliComparisonRes");\r
  17 \r
  18   // analyse comparison data (output stored in pictures_res.root)\r
  19   comp->Analyse();\r
  20   \r
  21   // paramtetrisation of the TPC track length (for information only) \r
  22   TF1 fl("fl","((min(250./(abs(x+0.000001)),250)-90))",0,2);  // TPC track length function\r
  23   TF1 fl2("fl2","[0]/((min(250./(abs(x+0.000001)),250)-90))^[1]",0,2);\r
  24   fl2.SetParameter(1,1);\r
  25   fl2.SetParameter(0,1);\r
  26 */\r
  27 \r
  28 #include <iostream>\r
  29 \r
  30 #include "TFile.h"\r
  31 #include "TCint.h"\r
  32 #include "TH3F.h"\r
  33 #include "TH2F.h"\r
  34 #include "TF1.h"\r
  35 #include "TProfile.h"\r
  36 #include "TProfile2D.h"\r
  37 #include "TGraph2D.h"\r
  38 #include "TCanvas.h"\r
  39 #include "TGraph.h"\r
  40 \r
  41 #include "AliESDEvent.h"   \r
  42 #include "AliESD.h"\r
  43 #include "AliESDfriend.h"\r
  44 #include "AliESDfriendTrack.h"\r
  45 #include "AliESDVertex.h"\r
  46 #include "AliRecInfoCuts.h" \r
  47 #include "AliMCInfoCuts.h" \r
  48 #include "AliLog.h" \r
  49 #include "AliTracker.h" \r
  50 \r
  51 #include "AliMathBase.h"\r
  52 #include "AliTreeDraw.h" \r
  53 \r
  54 #include "AliMCInfo.h" \r
  55 #include "AliESDRecInfo.h" \r
  56 #include "AliComparisonRes.h" \r
  57 \r
  58 using namespace std;\r
  59 \r
  60 ClassImp(AliComparisonRes)\r
  61 \r
  62 //_____________________________________________________________________________\r
  63 AliComparisonRes::AliComparisonRes():\r
  64   TNamed("AliComparisonRes","AliComparisonRes"),\r
  65 \r
  66   // Resolution \r
  67   fPtResolLPT(0),        // pt resolution - low pt\r
  68   fPtResolHPT(0),        // pt resolution - high pt \r
  69   fPtPullLPT(0),         // pt resolution - low pt\r
  70   fPtPullHPT(0),         // pt resolution - high pt \r
  71   //\r
  72   // Resolution constrained param\r
  73   //\r
  74   fCPhiResolTan(0),  // angular resolution -  constrained\r
  75   fCTanResolTan(0),  // angular resolution -  constrained\r
  76   fCPtResolTan(0),   // pt resolution      -  constrained\r
  77   fCPhiPullTan(0),   // angular resolution -  constrained\r
  78   fCTanPullTan(0),   // angular resolution -  constrained\r
  79   fCPtPullTan(0),    // pt resolution      -  constrained\r
  80 \r
  81   //\r
  82   // Parametrisation histograms\r
  83   //\r
  84   fC1Pt2Resol1PtTPC(0),\r
  85   fC1Pt2Resol1PtTPCITS(0),\r
  86   fCYResol1PtTPC(0),\r
  87   fCYResol1PtTPCITS(0),\r
  88   fCZResol1PtTPC(0),\r
  89   fCZResol1PtTPCITS(0),\r
  90   fCPhiResol1PtTPC(0),\r
  91   fCPhiResol1PtTPCITS(0),\r
  92   fCThetaResol1PtTPC(0),\r
  93   fCThetaResol1PtTPCITS(0),\r
  94 \r
  95   // vertex\r
  96   fVertex(0),\r
  97  \r
  98   // Cuts \r
  99   fCutsRC(0),  \r
 100   fCutsMC(0)  \r
 101 {\r
 102   InitHisto();\r
 103   InitCuts();\r
 104   \r
 105   // vertex (0,0,0)\r
 106   fVertex = new AliESDVertex();\r
 107   fVertex->SetXv(0.0);\r
 108   fVertex->SetYv(0.0);\r
 109   fVertex->SetZv(0.0);\r
 110 }\r
 111 \r
 112 //_____________________________________________________________________________\r
 113 AliComparisonRes::~AliComparisonRes(){\r
 114   \r
 115   // Resolution histograms\r
 116   if(fPtResolLPT) delete  fPtResolLPT; fPtResolLPT=0;     \r
 117   if(fPtResolHPT) delete  fPtResolHPT; fPtResolHPT=0;    \r
 118   if(fPtPullLPT)  delete  fPtPullLPT;  fPtPullLPT=0;    \r
 119   if(fPtPullHPT)  delete  fPtPullHPT;  fPtPullHPT=0;   \r
 120 \r
 121   // Resolution histograms (constrained param)\r
 122   if(fCPhiResolTan) delete fCPhiResolTan; fCPhiResolTan=0;\r
 123   if(fCTanResolTan) delete fCTanResolTan; fCTanResolTan=0;\r
 124   if(fCPtResolTan)  delete fCPtResolTan;  fCPtResolTan=0; \r
 125   if(fCPhiPullTan)  delete fCPhiPullTan;  fCPhiPullTan=0;\r
 126   if(fCTanPullTan)  delete fCTanPullTan;  fCTanPullTan=0;\r
 127   if(fCPtPullTan)   delete fCPtPullTan;   fCPtPullTan=0;\r
 128 \r
 129   // Parametrisation histograms\r
 130   if(fC1Pt2Resol1PtTPC) delete fC1Pt2Resol1PtTPC; fC1Pt2Resol1PtTPC=0;\r
 131   if(fC1Pt2Resol1PtTPCITS) delete fC1Pt2Resol1PtTPCITS; fC1Pt2Resol1PtTPCITS=0;\r
 132   if(fCYResol1PtTPC) delete fCYResol1PtTPC; fCYResol1PtTPC=0;\r
 133   if(fCYResol1PtTPCITS) delete fCYResol1PtTPCITS; fCYResol1PtTPCITS=0;\r
 134   if(fCZResol1PtTPC) delete fCZResol1PtTPC; fCZResol1PtTPC=0;\r
 135   if(fCZResol1PtTPCITS) delete fCZResol1PtTPCITS; fCZResol1PtTPCITS=0;\r
 136   if(fCPhiResol1PtTPC) delete fCPhiResol1PtTPC; fCPhiResol1PtTPC=0;\r
 137   if(fCPhiResol1PtTPCITS) delete fCPhiResol1PtTPCITS; fCPhiResol1PtTPCITS=0;\r
 138   if(fCThetaResol1PtTPC) delete fCThetaResol1PtTPC; fCThetaResol1PtTPC=0;\r
 139   if(fCThetaResol1PtTPCITS) delete fCThetaResol1PtTPCITS; fCThetaResol1PtTPCITS=0;\r
 140 \r
 141   if(fVertex) delete fVertex; fVertex=0;\r
 142 \r
 143 }\r
 144 \r
 145 //_____________________________________________________________________________\r
 146 void AliComparisonRes::InitHisto(){\r
 147 \r
 148   // Init histograms\r
 149   fCPhiResolTan = new TH2F("CPhiResolTan","CPhiResolTan",50, -2,2,200,-0.025,0.025);   \r
 150   fCPhiResolTan->SetXTitle("tan(#theta)");\r
 151   fCPhiResolTan->SetYTitle("#Delta#phi");\r
 152 \r
 153   fCTanResolTan = new TH2F("CTanResolTan","CTanResolTan",50, -2,2,200,-0.025,0.025);\r
 154   fCTanResolTan->SetXTitle("tan(#theta)");\r
 155   fCTanResolTan->SetYTitle("#Delta#theta");\r
 156 \r
 157   fCPtResolTan=new TH2F("CPtResol","CPtResol",50, -2,2,200,-0.2,0.2);    \r
 158   fCPtResolTan->SetXTitle("Tan(#theta)");\r
 159   fCPtResolTan->SetYTitle("#Deltap_{t}/p_{t}");\r
 160 \r
 161   fCPhiPullTan = new TH2F("CPhiPullTan","CPhiPullTan",50, -2,2,200,-5,5);   \r
 162   fCPhiPullTan->SetXTitle("Tan(#theta)");\r
 163   fCPhiPullTan->SetYTitle("#Delta#phi/#Sigma");\r
 164 \r
 165   fCTanPullTan = new TH2F("CTanPullTan","CTanPullTan",50, -2,2,200,-5,5);\r
 166   fCTanPullTan->SetXTitle("Tan(#theta)");\r
 167   fCTanPullTan->SetYTitle("#Delta#theta/#Sigma");\r
 168 \r
 169   fCPtPullTan=new TH2F("CPtPull","CPtPull",50, -2,2,200,-5,5);    \r
 170   fCPtPullTan->SetXTitle("Tan(#theta)");\r
 171   fCPtPullTan->SetYTitle("(1/mcp_{t}-1/p_{t})/#Sigma");\r
 172 \r
 173   fPtResolLPT = new TH2F("Pt_resol_lpt","pt resol",10, 0.1,3,200,-0.2,0.2);\r
 174   fPtResolLPT->SetXTitle("p_{t}");\r
 175   fPtResolLPT->SetYTitle("#Deltap_{t}/p_{t}");\r
 176 \r
 177   fPtResolHPT = new TH2F("Pt_resol_hpt","pt resol",10, 2,100,200,-0.3,0.3);  \r
 178   fPtResolHPT->SetXTitle("p_{t}");\r
 179   fPtResolHPT->SetYTitle("#Deltap_{t}/p_{t}");\r
 180 \r
 181   fPtPullLPT = new TH2F("Pt_pull_lpt","pt pull",10, 0.1,3,200,-6,6);\r
 182   fPtPullLPT->SetXTitle("p_{t}");\r
 183   fPtPullLPT->SetYTitle("#Deltap_{t}/#Sigma");\r
 184 \r
 185   fPtPullHPT = new TH2F("Pt_pull_hpt","pt pull",10,2,100,200,-6,6);  \r
 186   fPtPullHPT->SetXTitle("p_{t}");\r
 187   fPtPullHPT->SetYTitle("#Deltap_{t}/#Sigma");\r
 188 \r
 189   // Parametrisation histograms\r
 190   fC1Pt2Resol1PtTPC = new TH2F("fC1Pt2Resol1PtTPC","(1/mcpt-1/pt)/(1+1/mcpt)^2 vs 1/pt)",100,0,10,200,-0.020,0.020);  \r
 191   fC1Pt2Resol1PtTPC->SetXTitle("1/mcp_{t}");\r
 192   fC1Pt2Resol1PtTPC->SetYTitle("(1/mcp_{t}-1/p_{t})/(1+1/mcp_{t})^2)");\r
 193 \r
 194   fC1Pt2Resol1PtTPCITS = new TH2F("fC1Pt2Resol1PtTPCITS","(1/mcpt-1/pt)/(1+1/mcpt)^2 vs 1/pt)",100,0,10,200,-0.020,0.020);  \r
 195   fC1Pt2Resol1PtTPCITS->SetXTitle("1/mcp_{t}");\r
 196   fC1Pt2Resol1PtTPCITS->SetYTitle("(1/mcp_{t}-1/p_{t})/(1+1/mcp_{t})^2)");\r
 197 \r
 198   fCYResol1PtTPC = new TH2F("fCYResol1PtTPC","fCYResol1PtTPC",100, 0,10,200,-0.010,0.010);   \r
 199   fCYResol1PtTPC->SetXTitle("1/mcpt");\r
 200   fCYResol1PtTPC->SetYTitle("#DeltaY");\r
 201 \r
 202   fCYResol1PtTPCITS = new TH2F("fCYResol1PtTPCITS","fCYResol1PtTPCITS",100, 0,10,200,-0.010,0.010);   \r
 203   fCYResol1PtTPCITS->SetXTitle("1/mcpt");\r
 204   fCYResol1PtTPCITS->SetYTitle("#DeltaY");\r
 205 \r
 206   fCZResol1PtTPC = new TH2F("fCZResol1PtTPC","fCZResol1PtTPC",50, 0,10,200,-0.020,0.020);   \r
 207   fCZResol1PtTPC->SetXTitle("1/mcpt");\r
 208   fCZResol1PtTPC->SetYTitle("#DeltaZ");\r
 209 \r
 210   fCZResol1PtTPCITS = new TH2F("fCZResol1PtTPCITS","fCZResol1PtTPCITS",50, 0,10,200,-0.020,0.020);   \r
 211   fCZResol1PtTPCITS->SetXTitle("1/mcpt");\r
 212   fCZResol1PtTPCITS->SetYTitle("#DeltaZ");\r
 213 \r
 214   fCPhiResol1PtTPC = new TH2F("fCPhiResol1PtTPC","fCPhiResol1PtTPC",50, 0,10,200,-0.005,0.005);   \r
 215   fCPhiResol1PtTPC->SetXTitle("1/mcpt");\r
 216   fCPhiResol1PtTPC->SetYTitle("#Delta#phi");\r
 217 \r
 218   fCPhiResol1PtTPCITS = new TH2F("fCPhiResol1PtTPCITS","fCPhiResol1PtTPCITS",50, 0,10,200,-0.005,0.005);   \r
 219   fCPhiResol1PtTPCITS->SetXTitle("1/mcpt");\r
 220   fCPhiResol1PtTPCITS->SetYTitle("#Delta#phi");\r
 221 \r
 222   fCThetaResol1PtTPC = new TH2F("fCThetaResol1PtTPC","fCThetaResol1PtTPC",50, 0,10,200,-0.005,0.005);   \r
 223   fCThetaResol1PtTPC->SetXTitle("1/mcpt");\r
 224   fCThetaResol1PtTPC->SetYTitle("#Delta#theta");\r
 225 \r
 226   fCThetaResol1PtTPCITS = new TH2F("fCThetaResol1PtTPCITS","fCThetaResol1PtTPCITS",50, 0,10,200,-0.005,0.005);   \r
 227   fCThetaResol1PtTPCITS->SetXTitle("1/mcpt");\r
 228   fCThetaResol1PtTPCITS->SetYTitle("#Delta#theta");\r
 229 }\r
 230 \r
 231 //_____________________________________________________________________________\r
 232 void AliComparisonRes::InitCuts()\r
 233 {\r
 234   // Init cuts \r
 235   if(!fCutsMC) \r
 236     AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");\r
 237   if(!fCutsRC) \r
 238     AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");\r
 239 }\r
 240 \r
 241 //_____________________________________________________________________________\r
 242 void AliComparisonRes::Process(AliMCInfo* infoMC, AliESDRecInfo *infoRC)\r
 243 {\r
 244   // Fill resolution comparison information \r
 245   Float_t mcpt = infoMC->GetParticle().Pt();\r
 246 \r
 247   // distance to Prim. vertex \r
 248   const Double_t* dv = infoMC->GetVDist(); \r
 249 \r
 250   Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
 251 \r
 252   // Check selection cuts\r
 253   if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
 254   if (!isPrim) return;\r
 255   if (infoRC->GetStatus(1)==0) return;\r
 256   if (!infoRC->GetESDtrack()) return;  \r
 257   if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;\r
 258 \r
 259   // calculate and set prim. vertex\r
 260   fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
 261   fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
 262   fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
 263 \r
 264   Float_t deltaPt= (mcpt-infoRC->GetESDtrack()->Pt())/mcpt;  \r
 265   Float_t pullPt= (1/mcpt-infoRC->GetESDtrack()->OneOverPt())/TMath::Sqrt(infoRC->GetESDtrack()->GetSigma1Pt2());  \r
 266 \r
 267   // Fill histograms\r
 268   fPtResolLPT->Fill(mcpt,deltaPt);\r
 269   fPtResolHPT->Fill(mcpt,deltaPt);\r
 270   fPtPullLPT->Fill(mcpt,pullPt);\r
 271   fPtPullHPT->Fill(mcpt,pullPt);  \r
 272 }\r
 273 \r
 274 //_____________________________________________________________________________\r
 275 void AliComparisonRes::ProcessConstrained(AliMCInfo* infoMC, AliESDRecInfo *infoRC)\r
 276 {\r
 277   // Fill resolution comparison information (constarained parameters) \r
 278   //\r
 279   AliExternalTrackParam *track = 0;\r
 280   Double_t kRadius    = 3.0;      // beam pipe radius\r
 281   Double_t kMaxStep   = 5.0;      // max step\r
 282   Double_t field      = AliTracker::GetBz(); // nominal Bz field [kG]\r
 283   Double_t kMaxD      = 123456.0; // max distance\r
 284 \r
 285   Int_t clusterITS[200];\r
 286   Double_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
 287  \r
 288   Float_t mcpt = infoMC->GetParticle().Pt();\r
 289   Float_t tantheta = TMath::Tan(infoMC->GetParticle().Theta()-TMath::Pi()*0.5);\r
 290 \r
 291   // distance to Prim. vertex \r
 292   const Double_t* dv = infoMC->GetVDist(); \r
 293   Bool_t isPrim = TMath::Sqrt(dv[0]*dv[0] + dv[1]*dv[1])<fCutsMC->GetMaxR() && TMath::Abs(dv[2])<fCutsMC->GetMaxVz();\r
 294   \r
 295   // Check selection cuts\r
 296   if (fCutsMC->IsPdgParticle(TMath::Abs(infoMC->GetParticle().GetPdgCode())) == kFALSE) return; \r
 297   if (!isPrim) return;\r
 298   if (infoRC->GetStatus(1)!=3) return;\r
 299   if (!infoRC->GetESDtrack()) return;  \r
 300   if (infoRC->GetESDtrack()->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return;\r
 301   if (!infoRC->GetESDtrack()->GetConstrainedParam()) return;\r
 302 \r
 303   // calculate and set prim. vertex\r
 304   fVertex->SetXv( infoMC->GetParticle().Vx() - dv[0] );\r
 305   fVertex->SetYv( infoMC->GetParticle().Vy() - dv[1] );\r
 306   fVertex->SetZv( infoMC->GetParticle().Vz() - dv[2] );\r
 307 \r
 308   // constrained parameters resolution\r
 309   const AliExternalTrackParam * cparam = infoRC->GetESDtrack()->GetConstrainedParam();\r
 310   Float_t deltaCPt= (mcpt-cparam->Pt())/mcpt;  \r
 311   Float_t pullCPt= (1/mcpt-cparam->OneOverPt())/\r
 312     TMath::Sqrt(cparam->GetSigma1Pt2());          \r
 313   Float_t deltaPhi = TMath::ATan2(cparam->Py(),cparam->Px())-\r
 314     TMath::ATan2(infoMC->GetParticle().Py(),infoMC->GetParticle().Px());\r
 315   Float_t pullPhi = deltaPhi/TMath::Sqrt(cparam->GetSigmaSnp2()); \r
 316 \r
 317   Float_t deltaTan = TMath::ATan2(cparam->Pz(),cparam->Pt())-TMath::ATan2(infoMC->GetParticle().Pz(),infoMC->GetParticle().Pt());\r
 318   Float_t pullTan = deltaPhi/TMath::Sqrt(cparam->GetSigmaSnp2()); \r
 319 \r
 320 \r
 321   //Float_t delta1Pt2 = (1/mcpt-infoRC->GetESDtrack()->OneOverPt())/TMath::Power(1+1/mcpt,2));       \r
 322   Float_t delta1Pt2 = (1/mcpt-cparam->OneOverPt())/TMath::Power(1+1/mcpt,2);       \r
 323 \r
 324   Float_t deltaY1Pt = (infoMC->GetParticle().Vy()-cparam->GetY()) / (0.2+1/mcpt);\r
 325   Float_t deltaZ1Pt = (infoMC->GetParticle().Vz()-cparam->GetZ()) / (0.2+1/mcpt);\r
 326   Float_t deltaPhi1Pt = deltaPhi   / (0.1+1/mcpt);\r
 327   Float_t deltaTheta1Pt = deltaTan / (0.1+1/mcpt);\r
 328 \r
 329   // calculate track parameters at vertex\r
 330   if (infoRC->GetESDtrack()->GetTPCInnerParam())\r
 331   {\r
 332     if ((track = new AliExternalTrackParam(*infoRC->GetESDtrack()->GetTPCInnerParam())) != 0 )\r
 333     {\r
 334       Bool_t bStatus = AliTracker::PropagateTrackTo(track,kRadius,infoMC->GetMass(),kMaxStep,kTRUE);\r
 335       Bool_t bDCAStatus = track->PropagateToDCA(fVertex,field,kMaxD,dca,cov);\r
 336 \r
 337       // Fill parametrisation histograms (only TPC track)\r
 338       if(bStatus && bDCAStatus) \r
 339           {\r
 340         fC1Pt2Resol1PtTPC->Fill(1/mcpt,delta1Pt2);\r
 341         fCYResol1PtTPC->Fill(1/mcpt,deltaY1Pt);\r
 342         fCZResol1PtTPC->Fill(1/mcpt,deltaZ1Pt);\r
 343         fCPhiResol1PtTPC->Fill(1/mcpt,deltaPhi1Pt);\r
 344         fCThetaResol1PtTPC->Fill(1/mcpt,deltaTheta1Pt);\r
 345           }\r
 346           delete track;\r
 347     }\r
 348   }\r
 349 \r
 350  // TPC and ITS (nb. of clusters >2) in the system\r
 351   if(infoRC->GetStatus(1)==3 && infoRC->GetESDtrack()->GetITSclusters(clusterITS)>2) \r
 352   {\r
 353       fC1Pt2Resol1PtTPCITS->Fill(1/mcpt,delta1Pt2);\r
 354       fCYResol1PtTPCITS->Fill(1/mcpt,deltaY1Pt);\r
 355       fCZResol1PtTPCITS->Fill(1/mcpt,deltaZ1Pt);\r
 356       fCPhiResol1PtTPCITS->Fill(1/mcpt,deltaPhi1Pt);\r
 357       fCThetaResol1PtTPCITS->Fill(1/mcpt,deltaTheta1Pt);\r
 358   }\r
 359 \r
 360   // Fill histograms\r
 361   fCPtResolTan->Fill(tantheta,deltaCPt);\r
 362   fCPtPullTan->Fill(tantheta,pullCPt);\r
 363   fCPhiResolTan->Fill(tantheta,deltaPhi);\r
 364   fCPhiPullTan->Fill(tantheta,pullPhi);\r
 365   fCTanResolTan->Fill(tantheta,deltaTan);\r
 366   fCTanPullTan->Fill(tantheta,pullTan);\r
 367 }\r
 368  \r
 369 //_____________________________________________________________________________\r
 370 void AliComparisonRes::Exec(AliMCInfo* infoMC, AliESDRecInfo *infoRC){\r
 371   \r
 372   // Process comparison information \r
 373   Process(infoMC,infoRC);\r
 374   ProcessConstrained(infoMC,infoRC);\r
 375 }\r
 376 \r
 377 //_____________________________________________________________________________\r
 378 TH1F* AliComparisonRes::MakeResol(TH2F * his, Int_t integ, Bool_t type){\r
 379   // Create resolution histograms\r
 380  \r
 381   TH1F *hisr, *hism;\r
 382   if (!gPad) new TCanvas;\r
 383   hisr = AliTreeDraw::CreateResHistoI(his,&hism,integ);\r
 384   if (type) return hism;\r
 385   else \r
 386     return hisr;\r
 387 }\r
 388 \r
 389 //_____________________________________________________________________________\r
 390 void AliComparisonRes::Analyse(){\r
 391   // Analyse comparison information and store output histograms \r
 392   // in the "pictures_res.root" file \r
 393   \r
 394   AliComparisonRes * comp=this;\r
 395   TH1F *hiss=0;\r
 396 \r
 397   TFile *fp = new TFile("pictures_res.root","recreate");\r
 398   fp->cd();\r
 399 \r
 400   TCanvas * c = new TCanvas("Phi resol Tan","Phi resol Tan");\r
 401   c->cd();\r
 402   //\r
 403   hiss = comp->MakeResol(comp->fCPtResolTan,1,0);\r
 404   hiss->SetXTitle("Tan(#theta)");\r
 405   hiss->SetYTitle("#sigmap_{t}/p_{t}");\r
 406   hiss->Draw(); \r
 407   hiss->Write("CptResolTan");\r
 408   //\r
 409   hiss = comp->MakeResol(comp->fCPhiResolTan,1,0);\r
 410   hiss->SetXTitle("Tan(#theta)");\r
 411   hiss->SetYTitle("#sigma#phi (rad)");\r
 412   hiss->Draw();\r
 413   hiss->Write("PhiResolTan");\r
 414   //\r
 415   hiss = comp->MakeResol(comp->fCTanResolTan,1,0);\r
 416   hiss->SetXTitle("Tan(#theta)");\r
 417   hiss->SetYTitle("#sigma#theta (rad)");\r
 418   hiss->Draw();\r
 419   hiss->Write("ThetaResolTan");\r
 420   //\r
 421   hiss = comp->MakeResol(comp->fCPtPullTan,1,0);\r
 422   hiss->SetXTitle("Tan(#theta)");\r
 423   hiss->SetYTitle("1/mcp_{t}-1/p_{t}/#Sigma(1/p_{t})");\r
 424   hiss->Draw();\r
 425   hiss->Write("CptPullTan");\r
 426   //\r
 427   hiss = comp->MakeResol(comp->fC1Pt2Resol1PtTPC,1,0);\r
 428   hiss->SetXTitle("1/mcp_{t}");\r
 429   hiss->SetYTitle("1/mcp_{t}-1/p_{t}/(1+1/p_{t})^2");\r
 430   hiss->Draw();\r
 431   hiss->Write("Pt2Resol1PtTPC");\r
 432   fC1Pt2Resol1PtTPC->Write();\r
 433 \r
 434   hiss = comp->MakeResol(comp->fC1Pt2Resol1PtTPCITS,1,0);\r
 435   hiss->SetXTitle("1/mcp_{t}");\r
 436   hiss->SetYTitle("1/mcp_{t}-1/p_{t}/(1+1/p_{t})^2");\r
 437   hiss->Draw();\r
 438   hiss->Write("Pt2Resol1PtTPCITS");\r
 439   fC1Pt2Resol1PtTPCITS->Write();\r
 440   //\r
 441   hiss = comp->MakeResol(comp->fCYResol1PtTPC,1,0);\r
 442   hiss->SetXTitle("1/mcp_{t}");\r
 443   hiss->SetYTitle("(mcy-y)/(0.2+1/mcp_{t})");\r
 444   hiss->Draw();\r
 445   hiss->Write("CYResol1PtTPC");\r
 446   fCYResol1PtTPC->Write();\r
 447 \r
 448   hiss = comp->MakeResol(comp->fCYResol1PtTPCITS,1,0);\r
 449   hiss->SetXTitle("1/mcp_{t}");\r
 450   hiss->SetYTitle("(mcy-y)/(0.2+1/mcp_{t})");\r
 451   hiss->Draw();\r
 452   hiss->Write("CYResol1PtTPCITS");\r
 453   fCYResol1PtTPCITS->Write();\r
 454   //\r
 455   hiss = comp->MakeResol(comp->fCZResol1PtTPC,1,0);\r
 456   hiss->SetXTitle("1/mcp_{t}");\r
 457   hiss->SetYTitle("(mcz-z)/(0.2+1/mcp_{t})");\r
 458   hiss->Draw();\r
 459   hiss->Write("CZResol1PtTPC");\r
 460   fCZResol1PtTPC->Write();\r
 461 \r
 462   hiss = comp->MakeResol(comp->fCZResol1PtTPCITS,1,0);\r
 463   hiss->SetXTitle("1/mcp_{t}");\r
 464   hiss->SetYTitle("(mcz-z)/(0.2+1/mcp_{t})");\r
 465   hiss->Draw();\r
 466   hiss->Write("CZResol1PtTPCITS");\r
 467   fCZResol1PtTPCITS->Write();\r
 468   //\r
 469   hiss = comp->MakeResol(comp->fCPhiResol1PtTPC,1,0);\r
 470   hiss->SetXTitle("1/mcp_{t}");\r
 471   hiss->SetYTitle("(mc#phi-#phi)/(0.1+1/mcp_{t})");\r
 472   hiss->Draw();\r
 473   hiss->Write("CPhiResol1PtTPC");\r
 474   fCPhiResol1PtTPC->Write();\r
 475 \r
 476   hiss = comp->MakeResol(comp->fCPhiResol1PtTPCITS,1,0);\r
 477   hiss->SetXTitle("1/mcp_{t}");\r
 478   hiss->SetYTitle("(mc#phi-#phi)/(0.1+1/mcp_{t})");\r
 479   hiss->Draw();\r
 480   hiss->Write("CPhiResol1PtTPCITS");\r
 481   fCPhiResol1PtTPCITS->Write();\r
 482   //\r
 483   hiss = comp->MakeResol(comp->fCThetaResol1PtTPC,1,0);\r
 484   hiss->SetXTitle("1/mcp_{t}");\r
 485   hiss->SetYTitle("(mc#theta-#theta)/(0.1+1/mcp_{t})");\r
 486   hiss->Draw();\r
 487   hiss->Write("CThetaResol1PtTPC");\r
 488   fCThetaResol1PtTPC->Write();\r
 489 \r
 490   hiss = comp->MakeResol(comp->fCThetaResol1PtTPCITS,1,0);\r
 491   hiss->SetXTitle("1/mcp_{t}");\r
 492   hiss->SetYTitle("(mc#theta-#theta)/(0.1+1/mcp_{t})");\r
 493   hiss->Draw();\r
 494   hiss->Write("CThetaResol1PtTPCITS");\r
 495   fCThetaResol1PtTPCITS->Write();\r
 496 \r
 497   fp->Close();\r
 498 }\r
 499 \r
 500 //_____________________________________________________________________________\r
 501 Long64_t AliComparisonRes::Merge(TCollection* list) \r
 502 {\r
 503   // Merge list of objects (needed by PROOF)\r
 504 \r
 505   if (!list)\r
 506   return 0;\r
 507 \r
 508   if (list->IsEmpty())\r
 509   return 1;\r
 510 \r
 511   TIterator* iter = list->MakeIterator();\r
 512   TObject* obj = 0;\r
 513 \r
 514   // collection of generated histograms\r
 515   Int_t count=0;\r
 516   while((obj = iter->Next()) != 0) \r
 517   {\r
 518   AliComparisonRes* entry = dynamic_cast<AliComparisonRes*>(obj);\r
 519   if (entry == 0) continue; \r
 520 \r
 521   fPtResolLPT->Add(entry->fPtResolLPT);\r
 522   fPtResolHPT->Add(entry->fPtResolHPT);\r
 523   fPtPullLPT->Add(entry->fPtPullLPT);\r
 524   fPtPullHPT->Add(entry->fPtPullHPT);\r
 525 \r
 526   // Resolution histograms (constrained param)\r
 527   fCPhiResolTan->Add(entry->fCPhiResolTan);\r
 528   fCTanResolTan->Add(entry->fCTanResolTan);\r
 529   fCPtResolTan->Add(entry->fCPtResolTan);\r
 530   fCPhiPullTan->Add(entry->fCPhiPullTan);\r
 531   fCTanPullTan->Add(entry->fCTanPullTan);\r
 532   fCPtPullTan->Add(entry->fCPtPullTan);\r
 533 \r
 534   //  Histograms for 1/pt parameterisation\r
 535   fC1Pt2Resol1PtTPC->Add(entry->fC1Pt2Resol1PtTPC);\r
 536   fCYResol1PtTPC->Add(entry->fCYResol1PtTPC);\r
 537   fCZResol1PtTPC->Add(entry->fCZResol1PtTPC);\r
 538   fCPhiResol1PtTPC->Add(entry->fCPhiResol1PtTPC);\r
 539   fCThetaResol1PtTPC->Add(entry->fCThetaResol1PtTPC);\r
 540 \r
 541   fC1Pt2Resol1PtTPCITS->Add(entry->fC1Pt2Resol1PtTPCITS);\r
 542   fCYResol1PtTPCITS->Add(entry->fCYResol1PtTPCITS);\r
 543   fCZResol1PtTPCITS->Add(entry->fCZResol1PtTPCITS);\r
 544   fCPhiResol1PtTPCITS->Add(entry->fCPhiResol1PtTPCITS);\r
 545   fCThetaResol1PtTPCITS->Add(entry->fCThetaResol1PtTPCITS);\r
 546 \r
 547   count++;\r
 548   }\r
 549 \r
 550 return count;\r
 551 }\r