PWG1/TPC/AliMCTrackingTestTask.cxx

   1 //
   2 // This class is the task to check the
   3 // Propagation and Update method used in the
   4 //               1. AliExternalTrackParam
   5 //               2. AliTracker
   6 //
   7 // Pure Monte-Carlo data used, not influence of detectors
   8
   9 // Input - TParticle + Array of track references - (Points alogn track trajectories)
  10 // Output - Trees with track references - no histograms
  11 //          MC tree -  test for material budget correction
  12 //                     see function ProcessRefTracker
  13 //          MCupdate tree - test for correctness of propagation and update
  14 //                     see function AliMCTrackingTestTask::FitTrackRefs
  15 //
  16 // Principle - Creates AliExternalTrackParam form 1 Track Refernece -
  17 //             Propagate it to other
  18 // Magnetic field and the geometry has to be created before using it
  19
  20 //
  21 //
  22 //
  23
  24 //
  25 // ROOT includes
  26 #include <TChain.h>
  27 #include <TMath.h>
  28 #include <TVectorD.h>
  29 #include <TSystem.h>
  30 #include <TFile.h>
  31
  32 // ALIROOT includes
  33 #include <TTreeStream.h>
  34 #include <AliAnalysisManager.h>
  35 #include <AliESDInputHandler.h>
  36 #include "AliStack.h"
  37 #include "AliMCEvent.h"
  38 #include "AliMCEventHandler.h"
  39
  40 #include <AliESD.h>
  41 #include "AliMCTrackingTestTask.h"
  42 #include "AliGenInfoMaker.h"
  43 #include "AliHelix.h"
  44
  45 //
  46 #include "AliMCInfo.h"
  47 #include "AliESDRecInfo.h"
  48 #include "AliTPCParamSR.h"
  49 #include "AliTracker.h"
  50 #include "AliTPCseed.h"
  51
  52 // STL includes
  53 #include <iostream>
  54
  55 using namespace std;
  56
  57 ClassImp(AliMCTrackingTestTask)
  58
  59 //________________________________________________________________________
  60 AliMCTrackingTestTask::AliMCTrackingTestTask() :
  61   AliAnalysisTask(),
  62   fMCinfo(0),     //! MC event handler
  63   fESD(0),
  64   fDebugStreamer(0),
  65   fStreamLevel(0),
  66   fDebugLevel(0),
  67   fDebugOutputPath()
  68 {
  69   //
  70   // Default constructor (should not be used)
  71   //
  72 }
  73
  74 AliMCTrackingTestTask::AliMCTrackingTestTask(const AliMCTrackingTestTask& /*info*/) :
  75   AliAnalysisTask(),
  76   fMCinfo(0),     //! MC event handler
  77   fESD(0),
  78   //
  79   fDebugStreamer(0),
  80   fStreamLevel(0),
  81   fDebugLevel(),
  82   fDebugOutputPath()
  83 {
  84   //
  85   // Default constructor
  86   //
  87 }
  88
  89
  90
  91 //________________________________________________________________________
  92 AliMCTrackingTestTask::AliMCTrackingTestTask(const char *name) :
  93   AliAnalysisTask(name, "AliMCTrackingTestTask"),
  94   fMCinfo(0),     //! MC event handler
  95   fESD(0),
  96   fDebugStreamer(0),
  97   fStreamLevel(0),
  98   fDebugLevel(0),
  99   fDebugOutputPath()
 100 {
 101   //
 102   // Normal constructor
 103   //
 104   // Input slot #0 works with a TChain
 105   DefineInput(0, TChain::Class());
 106   // Output slot #0 writes into a TList
 107   DefineOutput(0, TObjArray::Class());
 108   //
 109   //
 110 }
 111
 112 AliMCTrackingTestTask::~AliMCTrackingTestTask(){
 113   //
 114   //
 115   //
 116   if (fDebugLevel>0)  printf("AliMCTrackingTestTask::~AliMCTrackingTestTask\n");
 117   if (fDebugStreamer) delete fDebugStreamer;
 118   fDebugStreamer=0;
 119 }
 120
 121
 122 //________________________________________________________________________
 123 void AliMCTrackingTestTask::ConnectInputData(Option_t *)
 124 {
 125   //
 126   // Connect the input data
 127   //
 128   if(fDebugLevel>3)
 129     cout << "AnalysisTaskTPCCluster::ConnectInputData()" << endl;
 130
 131   TTree* tree=dynamic_cast<TTree*>(GetInputData(0));
 132   if (!tree) {
 133     //Printf("ERROR: Could not read chain from input slot 0");
 134   }
 135   else {
 136     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 137     if (!esdH) {
 138       //Printf("ERROR: Could not get ESDInputHandler");
 139     }
 140     else {
 141       fESD = esdH->GetEvent();
 142       //Printf("*** CONNECTED NEW EVENT ****");
 143     }
 144   }
 145   AliMCEventHandler* mcinfo = (AliMCEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 146   mcinfo->SetReadTR(kTRUE);
 147
 148   fMCinfo = mcinfo->MCEvent();
 149
 150
 151 }
 152
 153
 154
 155
 156
 157
 158 //________________________________________________________________________
 159 void AliMCTrackingTestTask::CreateOutputObjects()
 160 {
 161   //
 162   // Connect the output objects
 163   //
 164   if(fDebugLevel>3)
 165     cout << "AnalysisTaskTPCCluster::CreateOutputObjects()" << endl;
 166
 167 }
 168
 169
 170 //________________________________________________________________________
 171 void AliMCTrackingTestTask::Exec(Option_t *) {
 172   //
 173   // Execute analysis for current event
 174   //
 175
 176   if(fDebugLevel>3)
 177     cout << "AliMCTrackingTestTask::Exec()" << endl;
 178
 179
 180   // If MC has been connected
 181
 182   if (!fMCinfo){
 183     cout << "Not MC info\n" << endl;
 184   }else{
 185     ProcessMCInfo();
 186     //mcinfo->Print();
 187     //DumpInfo();
 188   }
 189   //
 190 }
 191
 192
 193
 194
 195 //________________________________________________________________________
 196 void AliMCTrackingTestTask::Terminate(Option_t *) {
 197     //
 198     // Terminate loop
 199     //
 200   if(fDebugLevel>3)
 201     printf("AliMCTrackingTestTask: Terminate() \n");
 202   //
 203   if (fDebugLevel>0) printf("AliMCtrackingTestTask::Terminate\n");
 204   if (fDebugStreamer) delete fDebugStreamer;
 205   fDebugStreamer = 0;
 206   return;
 207 }
 208
 209
 210
 211 TTreeSRedirector *AliMCTrackingTestTask::GetDebugStreamer(){
 212   //
 213   // Get Debug streamer
 214   // In case debug streamer not yet initialized and StreamLevel>0 create new one
 215   //
 216   if (fStreamLevel==0) return 0;
 217   if (fDebugStreamer) return fDebugStreamer;
 218   TString dsName;
 219   dsName=GetName();
 220   dsName+="Debug.root";
 221   dsName.ReplaceAll(" ","");
 222   fDebugStreamer = new TTreeSRedirector(dsName.Data());
 223   return fDebugStreamer;
 224 }
 225
 226
 227
 228
 229
 230
 231 AliExternalTrackParam * AliMCTrackingTestTask::MakeTrack(const AliTrackReference* ref, TParticle*part)
 232 {
 233   //
 234   // Make track out of the track ref
 235   // part - TParticle used to determine chargr
 236   // the covariance matrix - equal 0 - starting from ideal MC position
 237   Double_t xyz[3]={ref->X(),ref->Y(),ref->Z()};
 238   Double_t pxyz[3]={ref->Px(),ref->Py(),ref->Pz()};
 239   Double_t cv[21];
 240   for (Int_t i=0; i<21;i++) cv[i]=0;
 241   if (!part->GetPDG()) return 0;
 242   AliExternalTrackParam * param = new AliExternalTrackParam(xyz,pxyz,cv,TMath::Nint(part->GetPDG()->Charge()/3.));
 243   return param;
 244 }
 245
 246 Bool_t  AliMCTrackingTestTask::PropagateToPoint(AliExternalTrackParam *param, Double_t *xyz, Double_t mass, Float_t step){
 247   //
 248   // Propagate track to point xyz using
 249   // AliTracker::PropagateTo functionality
 250   //
 251   //  param - track parameters
 252   //  xyz   - position to propagate
 253   //  mass  - particle mass
 254   //  step  - step to be used
 255   Double_t radius=TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
 256   AliTracker::PropagateTrackTo(param, radius+step, mass, step, kTRUE,0.99);
 257   AliTracker::PropagateTrackTo(param, radius+0.5, mass, step*0.1, kTRUE,0.99);
 258   Double_t sxyz[3]={0,0,0};
 259   AliESDVertex vertex(xyz,sxyz);
 260   Bool_t isOK = param->PropagateToDCA(&vertex,AliTracker::GetBz(),10);
 261   return isOK;
 262 }
 263
 264
 265 void  AliMCTrackingTestTask::ProcessMCInfo(){
 266   //
 267   //
 268   //
 269    //
 270   TParticle * particle= new TParticle;
 271   TClonesArray * trefs = new TClonesArray("AliTrackReference");
 272   const Double_t kPcut=0.1;
 273   //
 274   //
 275   // Process tracks
 276   //
 277   Int_t npart = fMCinfo->GetNumberOfTracks();
 278   if (npart==0) return;
 279   Double_t vertex[4]={0,0,0,0};
 280   fMCinfo->GetParticleAndTR(0, particle, trefs);
 281   if (particle){
 282     vertex[0]=particle->Vx();
 283     vertex[1]=particle->Vy();
 284     vertex[2]=particle->Vz();
 285     vertex[3]=particle->R();
 286   }
 287   //
 288   //
 289
 290   for (Int_t ipart=0;ipart<npart;ipart++){
 291     Int_t status = fMCinfo->GetParticleAndTR(ipart, particle, trefs);
 292     if (status<0) continue;
 293     if (!particle) continue;
 294     if (!trefs) continue;
 295     Int_t nref = trefs->GetEntries();
 296     if (nref<5) continue;
 297     FitTrackRefs(particle,trefs);
 298     AliTrackReference * tpcIn=0;
 299     AliTrackReference * tpcOut=0;
 300     AliTrackReference * trdIn=0;
 301     AliTrackReference * trdOut=0;
 302     AliTrackReference * itsIn=0;
 303     AliTrackReference * itsOut=0;
 304     Double_t rmax=0;
 305     Double_t rmin=1000;
 306     for (Int_t iref=0;iref<nref;iref++){
 307       AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
 308       if (!ref) continue;
 309
 310       Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
 311
 312       if (dir<0) break; // oposite direction - looping track - return back
 313       if (ref->P()<kPcut) continue;
 314       if (ref->R()<rmax) break;
 315       //if (ref->R()<rmin)  break;
 316       //TPC
 317       if (ref->DetectorId()==AliTrackReference::kTPC){
 318         if (!tpcIn) {
 319           tpcIn  = ref;
 320         }else{
 321           if (ref->R()>tpcIn->R()) tpcOut = ref;
 322         }
 323       }
 324       //ITS
 325       if (ref->DetectorId()==AliTrackReference::kITS){
 326         if (!itsIn) {
 327           itsIn  = ref;
 328         }else{
 329           if (ref->R()>itsIn->R()) itsOut = ref;
 330         }
 331       }
 332       //TRD
 333       if (ref->DetectorId()==AliTrackReference::kTRD){
 334         if (!trdIn) {
 335           trdIn  = ref;
 336         }else{
 337           if (ref->R()>trdIn->R()) trdOut = ref;
 338         }
 339       }
 340       if (ref->R()<rmin) rmin=ref->R();
 341       if (ref->R()>rmax) rmax=ref->R();
 342     }
 343     if (tpcIn && tpcOut) {
 344       ProcessRefTracker(tpcIn,tpcOut,particle,1);
 345       ProcessRefTracker(tpcIn,tpcOut,particle,3);
 346     }
 347     if (itsIn && itsOut) ProcessRefTracker(itsIn,itsOut,particle,0);
 348     if (trdIn && trdOut) ProcessRefTracker(trdIn,trdOut,particle,2);
 349   }
 350 }
 351
 352
 353
 354 void AliMCTrackingTestTask::ProcessRefTracker(AliTrackReference* refIn,  AliTrackReference* refOut, TParticle*part,Int_t type){
 355   //
 356   // Test propagation from In to out
 357   //
 358   AliExternalTrackParam *param = 0;
 359   AliExternalTrackParam *paramMC = 0;
 360   Double_t xyzIn[3]={refIn->X(),refIn->Y(), refIn->Z()};
 361   Double_t mass = part->GetMass();
 362   Double_t step=1;
 363   //
 364   param=MakeTrack(refOut,part);
 365   paramMC=MakeTrack(refOut,part);
 366   if (!param) return;
 367   if (type<3) PropagateToPoint(param,xyzIn, mass, step);
 368   if (type==3) {
 369     AliTPCseed seed;
 370     seed.Set(param->GetX(),param->GetAlpha(),param->GetParameter(),param->GetCovariance());
 371     Float_t alpha= TMath::ATan2(refIn->Y(),refIn->X());
 372     if(seed.Rotate(alpha-seed.GetAlpha())==kFALSE) return;
 373     seed.SetMass(mass);
 374     for (Float_t xlayer= seed.GetX(); xlayer>refIn->R(); xlayer-=step){
 375       seed.PropagateTo(xlayer);
 376     }
 377     seed.PropagateTo(refIn->R());
 378     param->Set(seed.GetX(),seed.GetAlpha(),seed.GetParameter(),seed.GetCovariance());
 379   }
 380   TTreeSRedirector *pcstream = GetDebugStreamer();
 381   TVectorD gpos(3);
 382   TVectorD gmom(3);
 383   param->GetXYZ(gpos.GetMatrixArray());
 384   param->GetPxPyPz(gmom.GetMatrixArray());
 385   if (pcstream){
 386     (*pcstream)<<"MC"<<
 387       "type="<<type<<
 388       "step="<<step<<
 389       "refIn.="<<refIn<<
 390       "refOut.="<<refOut<<
 391       "p.="<<part<<
 392       "par.="<<param<<
 393       "parMC.="<<paramMC<<
 394       "gpos.="<<&gpos<<
 395       "gmom.="<<&gmom<<
 396       "\n";
 397   }
 398 }
 399
 400
 401 void  AliMCTrackingTestTask::FitTrackRefs(TParticle * part, TClonesArray * trefs){
 402   //
 403   //
 404   //
 405   //
 406   const Int_t kMinRefs=6;
 407   Int_t nrefs = trefs->GetEntries();
 408   if (nrefs<kMinRefs) return; // we should have enough references
 409   Int_t iref0 =-1;
 410   Int_t iref1 =-1;
 411
 412   for (Int_t iref=0; iref<nrefs; iref++){
 413     AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
 414     if (!ref) continue;
 415     Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
 416     if (dir<0) break;
 417     if (ref->DetectorId()!=AliTrackReference::kTPC) continue;
 418     if (iref0<0) iref0 = iref;
 419     iref1 = iref;
 420   }
 421   if (iref1-iref0<kMinRefs) return;
 422   Double_t covar[15];
 423   for (Int_t icov=0; icov<15; icov++) covar[icov]=0;
 424   covar[0]=1;
 425   covar[2]=1;
 426   covar[5]=1;
 427   covar[9]=1;
 428   covar[14]=1;
 429
 430   AliTrackReference * refIn = (AliTrackReference*)trefs->At(iref0);
 431   AliTrackReference * refOut = (AliTrackReference*)trefs->At(iref1);
 432   AliExternalTrackParam *paramPropagate= MakeTrack(refIn,part);
 433   AliExternalTrackParam *paramUpdate   = MakeTrack(refIn,part);
 434   paramUpdate->AddCovariance(covar);
 435   Double_t mass = part->GetMass();
 436   Double_t charge = part->GetPDG()->Charge()/3.;
 437 /*
 438   Float_t alphaIn= TMath::ATan2(refIn->Y(),refIn->X());
 439   Float_t radiusIn= refIn->R();
 440   Float_t alphaOut= TMath::ATan2(refOut->Y(),refOut->X());
 441   Float_t radiusOut= refOut->R();
 442 */
 443   Bool_t isOKP=kTRUE;
 444   Bool_t isOKU=kTRUE;
 445   AliMagF * field = (AliMagF*) TGeoGlobalMagField::Instance()->GetField();
 446   for (Int_t iref = iref0; iref<=iref1; iref++){
 447     AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
 448     Float_t alphaC= TMath::ATan2(ref->Y(),ref->X());
 449     Double_t pos[3] = {ref->X(), ref->Y(), ref->Z()};
 450     Double_t mag[3];
 451     field->Field(pos,mag);
 452     isOKP&=paramPropagate->Rotate(alphaC);
 453     isOKU&=paramUpdate->Rotate(alphaC);
 454     for (Float_t xref= paramPropagate->GetX(); xref<ref->R(); xref++){
 455       isOKP&=paramPropagate->PropagateTo(xref, mag[2]);
 456       isOKU&=paramUpdate->PropagateTo(xref, mag[2]);
 457     }
 458     isOKP&=paramPropagate->PropagateTo(ref->R(), mag[2]);
 459     isOKU&=paramUpdate->PropagateTo(ref->R(), mag[2]);
 460     Double_t clpos[2] = {0, ref->Z()};
 461     Double_t clcov[3] = { 0.005,0,0.005};
 462     isOKU&= paramUpdate->Update(clpos, clcov);
 463   }
 464   TTreeSRedirector *pcstream = GetDebugStreamer();
 465   if (pcstream){
 466     TVectorD gposU(3);
 467     TVectorD gmomU(3);
 468     TVectorD gposP(3);
 469     TVectorD gmomP(3);
 470     paramUpdate->GetXYZ(gposU.GetMatrixArray());
 471     paramUpdate->GetPxPyPz(gmomU.GetMatrixArray());
 472     paramPropagate->GetXYZ(gposP.GetMatrixArray());
 473     paramPropagate->GetPxPyPz(gmomP.GetMatrixArray());
 474
 475      (*pcstream)<<"MCupdate"<<
 476        "isOKU="<<isOKU<<
 477        "isOKP="<<isOKP<<
 478        "m="<<mass<<
 479        "q="<<charge<<
 480        "part.="<<part<<
 481        "refIn.="<<refIn<<
 482        "refOut.="<<refOut<<
 483        "pP.="<<paramPropagate<<
 484        "pU.="<<paramUpdate<<
 485        "gposU.="<<&gposU<<
 486        "gmomU.="<<&gmomU<<
 487        "gposP.="<<&gposP<<
 488        "gmomP.="<<&gmomP<<
 489        "\n";
 490    }
 491 }
 492
 493
 494
 495
 496 void AliMCTrackingTestTask::FinishTaskOutput()
 497 {
 498   //
 499   // According description in AliAnalisysTask this method is call
 500   // on the slaves before sending data
 501   //
 502   Terminate("slave");
 503   gSystem->Exec("pwd");
 504   RegisterDebugOutput();
 505
 506 }
 507
 508
 509 void AliMCTrackingTestTask::RegisterDebugOutput(){
 510   //
 511   //
 512   //
 513   //
 514   // store  - copy debug output to the destination position
 515   // currently ONLY for local copy
 516   TString dsName;
 517   dsName=GetName();
 518   dsName+="Debug.root";
 519   dsName.ReplaceAll(" ","");
 520   TString dsName2=fDebugOutputPath.Data();
 521   gSystem->MakeDirectory(dsName2.Data());
 522   dsName2+=gSystem->HostName();
 523   gSystem->MakeDirectory(dsName2.Data());
 524   dsName2+="/";
 525   dsName2+=gSystem->BaseName(gSystem->pwd());
 526   dsName2+="/";
 527   gSystem->MakeDirectory(dsName2.Data());
 528   dsName2+=dsName;
 529   AliInfo(Form("copy %s\t%s\n",dsName.Data(),dsName2.Data()));
 530   printf("copy %s\t%s\n",dsName.Data(),dsName2.Data());
 531   TFile::Cp(dsName.Data(),dsName2.Data());
 532 }
 533
 534
 535 /*
 536   gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
 537   gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
 538   AliXRDPROOFtoolkit tool;
 539   TChain * chain = tool.MakeChain("mctracking.txt","MC",0,100);
 540   chain->Lookup();
 541   //
 542   //
 543   chain->SetAlias("pdg","(p.fPdgCode)");
 544   chain->SetAlias("dPRec","(refOut.P()-par.P())/refIn.P()");
 545   chain->SetAlias("dPMC","(refOut.P()-refIn.P())/refIn->P()");
 546   chain->SetAlias("dPtRec","(refOut.Pt()-par.Pt())/refIn.Pt()");
 547   chain->SetAlias("dPtMC","(refOut.Pt()-refIn.Pt())/refIn->Pt()");
 548
 549
 550   // ITS
 551   chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==0&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 552   htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
 553   htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
 554   gPad->SaveAs("picLoss/dPcorr_ITS_step1.gif");
 555   gPad->SaveAs("picLoss/dPcorr_ITS_step1.eps");
 556   // TPC
 557   chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==1&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 558   htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
 559   htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
 560   gPad->SaveAs("picLoss/dPcorr_TPC_step1.gif");
 561   gPad->SaveAs("picLoss/dPcorr_TPC_step1.eps");
 562   //
 563    // TPC
 564   chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==3&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 565   htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
 566   htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
 567   gPad->SaveAs("picLoss/dPcorr_TPCseed_step1.gif");
 568   gPad->SaveAs("picLoss/dPcorr_TPCseed_step1.eps");
 569
 570
 571   // TRD
 572   chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==2&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 573   htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
 574   htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
 575   gPad->SaveAs("picLoss/dPcorr_TRD_step1.gif");
 576   gPad->SaveAs("picLoss/dPcorr_TRD_step1.eps");
 577
 578   //
 579   //
 580   //
 581   chain->Draw("(par.Pt()-refIn.Pt())/refIn.Pt()>>his(100,-0.02,0.02)","abs(pdg)!=11&&type==3&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 582   his->SetXTitle("(P_{trec}-P_{tmc})/P_{tmc}");
 583   gPad->SaveAs("picLoss/dPtcorr_TPCseed_step1_1D.eps");
 584   gPad->SaveAs("picLoss/dPtcorr_TPCseed_step1_1D.gif");
 585
 586   chain->Draw("(par.P()-refIn.P())/refIn.P()>>his(100,-0.02,0.02)","abs(pdg)!=11&&type==3&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
 587   his->SetXTitle("(P_{rec}-P_{mc})/P_{mc}");
 588   gPad->SaveAs("picLoss/dPcorr_TPCseed_step1_1D.eps");
 589   gPad->SaveAs("picLoss/dPcorr_TPCseed_step1_1D.gif");
 590 */