2 // This class is the task to check the
3 // Propagation and Update method used in the
4 // 1. AliExternalTrackParam
7 // Pure Monte-Carlo data used, not influence of detectors
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
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
34 #include <TTreeStream.h>
35 #include <AliAnalysisManager.h>
36 #include <AliESDInputHandler.h>
38 #include "AliMCEvent.h"
39 #include "AliMCEventHandler.h"
42 #include "AliTrackComparison.h"
43 #include "AliMCTrackingTestTask.h"
44 #include "AliGenInfoMaker.h"
46 #include "AliTrackPointArray.h"
47 #include "AliESDCaloCluster.h"
50 #include "AliMCInfo.h"
51 #include "AliComparisonObject.h"
52 #include "AliESDRecInfo.h"
53 #include "AliTPCParamSR.h"
54 #include "AliTracker.h"
55 #include "AliTPCseed.h"
57 #include "AliCDBManager.h"
58 #include "AliGRPManager.h"
59 #include "AliGeomManager.h"
66 ClassImp(AliMCTrackingTestTask)
68 #define USE_STREAMER 0
71 //________________________________________________________________________
72 AliMCTrackingTestTask::AliMCTrackingTestTask() :
74 fMCinfo(0), //! MC event handler
86 // Default constructor (should not be used)
90 AliMCTrackingTestTask::AliMCTrackingTestTask(const AliMCTrackingTestTask& /*info*/) :
92 fMCinfo(0), //! MC event handler
104 // Default constructor
110 //________________________________________________________________________
111 AliMCTrackingTestTask::AliMCTrackingTestTask(const char *name) :
112 AliAnalysisTask(name, "AliMCTrackingTestTask"),
113 fMCinfo(0), //! MC event handler
125 // Normal constructor
127 // Input slot #0 works with a TChain
128 DefineInput(0, TChain::Class());
129 DefineOutput(0, TList::Class());
131 // create the list for comparison objects
132 fCompList = new TList;
137 AliMCTrackingTestTask::~AliMCTrackingTestTask(){
141 if (fDebugLevel>0) printf("AliMCTrackingTestTask::~AliMCTrackingTestTask\n");
142 if (fDebugStreamer) delete fDebugStreamer;
144 if (fOutList) delete fOutList; fOutList = 0;
145 if (fCompList) delete fCompList; fCompList = 0;
149 //_____________________________________________________________________________
150 Bool_t AliMCTrackingTestTask::AddComparisonObject(AliTrackComparison *cObj)
152 // add comparison object to the list
154 Printf("ERROR: Could not add comparison object");
158 // add object to the list
159 fCompList->AddLast(cObj);
163 //________________________________________________________________________
164 void AliMCTrackingTestTask::ConnectInputData(Option_t *)
167 // Connect the input data
170 cout << "AnalysisTaskTPCCluster::ConnectInputData()" << endl;
172 TTree* tree=dynamic_cast<TTree*>(GetInputData(0));
174 //Printf("ERROR: Could not read chain from input slot 0");
177 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
179 //Printf("ERROR: Could not get ESDInputHandler");
182 fESD = esdH->GetEvent();
183 //Printf("*** CONNECTED NEW EVENT ****");
186 AliMCEventHandler* mcinfo = (AliMCEventHandler*) (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
187 mcinfo->SetReadTR(kTRUE);
189 fMCinfo = mcinfo->MCEvent();
192 //________________________________________________________________________
193 void AliMCTrackingTestTask::CreateOutputObjects()
196 // Connect the output objects
199 cout << "AnalysisTaskTPCCluster::CreateOutputObjects()" << endl;
202 fOutList = new TList;
203 fOutList->SetOwner(kTRUE);
204 fPitList = fOutList->MakeIterator();
206 // create output list
208 // add comparison objects to the output
209 AliTrackComparison *cObj=0;
211 TIterator *pitCompList = fCompList->MakeIterator();
212 pitCompList->Reset();
213 while(( cObj = (AliTrackComparison*)pitCompList->Next()) != NULL) {
217 Printf("UserCreateOutputObjects(): Number of output comparison objects: %d \n", count);
219 PostData(0, fOutList);
223 //________________________________________________________________________
224 void AliMCTrackingTestTask::Exec(Option_t *) {
226 // Execute analysis for current event
230 printf("New event!\n");
234 cout << "AliMCTrackingTestTask::Exec()" << endl;
238 // If MC has been connected
240 cout << "Not MC info\n" << endl;
247 PostData(0, fOutList);
251 //________________________________________________________________________
252 void AliMCTrackingTestTask::Terminate(Option_t *) {
257 printf("AliMCTrackingTestTask: Terminate() \n");
259 if (fDebugLevel>0) printf("AliMCtrackingTestTask::Terminate\n");
260 if (fDebugStreamer) delete fDebugStreamer;
263 // AliTrackComparison *cObj=0;
264 // TIterator *pitCompList = fCompList->MakeIterator();
265 // pitCompList->Reset();
266 // while(( cObj = (AliTrackComparison*)pitCompList->Next()) != NULL) {
267 // for(Int_t i=0; i<6; i++)
268 // cObj->MakeDistortionMap(438,i);
275 //________________________________________________________________________
276 TTreeSRedirector *AliMCTrackingTestTask::GetDebugStreamer(){
278 // Get Debug streamer
279 // In case debug streamer not yet initialized and StreamLevel>0 create new one
281 if (fStreamLevel==0) return 0;
282 if (fDebugStreamer) return fDebugStreamer;
285 dsName+="Debug.root";
287 printf(" get debug streamer \n");
289 dsName.ReplaceAll(" ","");
290 fDebugStreamer = new TTreeSRedirector(dsName.Data());
291 return fDebugStreamer;
294 //________________________________________________________________________
295 AliExternalTrackParam * AliMCTrackingTestTask::MakeTrack(const AliTrackReference* ref, TParticle*part)
298 // Make track out of the track ref
299 // part - TParticle used to determine chargr
300 // the covariance matrix - equal 0 - starting from ideal MC position
301 Double_t xyz[3]={ref->X(),ref->Y(),ref->Z()};
302 Double_t pxyz[3]={ref->Px(),ref->Py(),ref->Pz()};
304 for (Int_t i=0; i<21;i++) cv[i]=0;
305 if (!part->GetPDG()) return 0;
306 AliExternalTrackParam * param = new AliExternalTrackParam(xyz,pxyz,cv,TMath::Nint(part->GetPDG()->Charge()/3.));
310 //________________________________________________________________________
311 Bool_t AliMCTrackingTestTask::PropagateToPoint(AliExternalTrackParam *param, Double_t *xyz, Double_t mass, Float_t step){
313 // Propagate track to point xyz using
314 // AliTracker::PropagateToBxByBz functionality
316 // param - track parameters
317 // xyz - position to propagate
318 // mass - particle mass
319 // step - step to be used
320 Double_t radius=TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
321 Double_t alpha = TMath::ATan2(xyz[1],xyz[0]);
323 AliTracker::PropagateTrackToBxByBz(param, radius+step, mass, step, kTRUE,0.99,-1);
324 param->Rotate(alpha);
325 Bool_t isOK = param->PropagateTo(radius,AliTracker::GetBz());
330 //________________________________________________________________________
331 void AliMCTrackingTestTask::ProcessMCInfo(){
337 printf("ProcessMCInfo\n");
340 const AliESDVertex *pVertex = fESD->GetPrimaryVertex();
341 if(!pVertex) AliError("No primary vertex found!\n");
343 pVertex->GetXYZ(vPos);
345 TParticle * particle= new TParticle();
346 TClonesArray * trefs = new TClonesArray("AliTrackReference");
347 const Double_t kPcut=0.1;
352 Int_t npart = fMCinfo->GetNumberOfTracks();
353 if (npart==0) return;
354 Double_t vertex[4]={0,0,0,0};
355 fMCinfo->GetParticleAndTR(0, particle, trefs);
357 vertex[0]=particle->Vx();
358 vertex[1]=particle->Vy();
359 vertex[2]=particle->Vz();
360 vertex[3]=particle->R();
366 for (Int_t ipart=0;ipart<npart;ipart++){
367 Int_t status = fMCinfo->GetParticleAndTR(ipart, particle, trefs);
368 if (status<0 || !particle || !trefs) continue;
369 Int_t nref = trefs->GetEntries();
370 if (nref<5) continue;
371 FitTrackRefs(particle,trefs);
373 AliTrackReference * tpcIn=0;
374 AliTrackReference * tpcOut=0;
375 AliTrackReference * trdIn=0;
376 AliTrackReference * trdOut=0;
377 AliTrackReference * itsIn=0;
378 AliTrackReference * itsOut=0;
380 AliTrackReference * tofIn=0;
381 AliTrackReference * tofOut=0;
382 AliTrackReference * hmpidIn=0;
383 AliTrackReference * hmpidOut=0;
384 AliTrackReference * emcalIn=0;
385 AliTrackReference * emcalOut=0;
389 for (Int_t iref=0;iref<nref;iref++){
390 AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
393 Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
395 if (dir<0) break; // oposite direction - looping track - return back
396 if (ref->P()<kPcut) continue;
397 if (ref->R()<rmax) break;
398 //if (ref->R()<rmin) break;
400 if (ref->DetectorId()==AliTrackReference::kTPC){
404 if (ref->R()>tpcIn->R()) tpcOut = ref;
408 if (ref->DetectorId()==AliTrackReference::kITS){
412 if (ref->R()>itsIn->R()) itsOut = ref;
416 if (ref->DetectorId()==AliTrackReference::kTRD){
420 if (ref->R()>trdIn->R()) trdOut = ref;
424 if (ref->DetectorId()==AliTrackReference::kTOF){
428 if (ref->R()>tofIn->R()) tofOut = ref;
433 if (ref->DetectorId()==AliTrackReference::kHMPID){
437 if (ref->R()>hmpidIn->R()) hmpidOut = ref;
443 if (ref->DetectorId()==AliTrackReference::kEMCAL){
447 if (ref->R()>emcalIn->R()) emcalOut = ref;
451 if (ref->R()<rmin) rmin=ref->R();
452 if (ref->R()>rmax) rmax=ref->R();
455 // -----------------------------------------------
459 if ( TMath::Abs(particle->GetPdgCode()) == 11 ) {
462 if (IsFindable(ipart,70)) {
464 // is from gamma conversion
465 Int_t motherId = particle->GetFirstMother();
467 if (motherId < npart) {
468 TParticle* mother = (fMCinfo->Stack())->Particle(motherId);
469 if (mother && TMath::Abs(mother->GetPdgCode()) == 22) {
470 Int_t nDaughters = mother->GetNDaughters();
472 for (Int_t idx=0; idx<nDaughters; ++idx) {
473 Int_t daughterId = mother->GetDaughter(idx);
474 if ( daughterId == ipart || daughterId >= npart )
477 TParticle* daughter = (fMCinfo->Stack())->Particle(daughterId);
478 if (daughter && TMath::Abs(daughter->GetPdgCode()) == 11) {
479 //Bool_t findable = IsFindable(daughterId,70);
481 TTreeSRedirector *pcstream = GetDebugStreamer();
483 (*pcstream)<<"MCgamma"<<
484 "triggerP="<<particle<< // trigger electron
485 "motherP="<<mother<< // mother gamma
486 "daughterP="<<daughter<< // daughter electron
487 "isFindable="<<findable<< // 2nd is findable
499 // -----------------------------------------------
500 if (tpcIn && tpcOut) {
501 ProcessRefTracker(tpcIn,tpcOut,particle,1);
502 ProcessRefTracker(tpcIn,tpcOut,particle,3);
504 if (itsIn && itsOut) ProcessRefTracker(itsIn, itsOut, particle, 0);
505 if (trdIn && trdOut) ProcessRefTracker(trdIn, trdOut, particle, 2);
507 if (tpcOut && trdIn) ProcessRefTracker(tpcOut,trdIn, particle, 4);
508 if (tpcOut && tofIn) ProcessRefTracker(tpcOut,tofIn, particle, 5);
509 if (tpcOut && hmpidIn) ProcessRefTracker(tpcOut,hmpidIn,particle, 6);
510 if (tpcOut && emcalIn) ProcessRefTracker(tpcOut,emcalIn,particle, 7);
512 if (tpcOut && trdOut) ProcessRefTracker(tpcOut,trdOut, particle, 8);
513 if (trdIn && tofIn) ProcessRefTracker(trdIn, tofIn, particle, 9);
514 if (tofIn && tofOut) ProcessRefTracker(tofIn, tofOut, particle,10);
517 // -----------------------------------------------
518 //Test for new base tracking class
520 AliTrackComparison *cObj=0;
521 AliTrackPoint *point=new AliTrackPoint();
522 AliESDCaloCluster *cluster=0;
523 AliExternalTrackParam *tpcOuter=0;
528 for(Int_t iCluster=0; iCluster<fESD->GetNumberOfCaloClusters(); iCluster++)
530 cluster = fESD->GetCaloCluster(iCluster);
531 if(!cluster->IsEMCAL()) continue;
532 if(cluster->GetLabel()==ipart && cluster->E()>eMax)
541 if(cluster) cluster=0;
542 cluster = fESD->GetCaloCluster(clsIndex);
543 Float_t clusterPos[3];
544 cluster->GetPosition(clusterPos);
545 point->SetXYZ(clusterPos[0],clusterPos[1],clusterPos[2],0);
549 tpcOuter = MakeTrack(tpcOut,particle);
552 Double_t mass = particle->GetMass();
553 Int_t charge = TMath::Nint(particle->GetPDG()->Charge()/3.);
555 while(( cObj = (AliTrackComparison *)fPitList->Next()) != NULL) {
556 TString objName(cObj->GetName());
557 if(!objName.CompareTo("TPCOutToEMCalInElecCls"))
559 if(TMath::Abs(particle->GetPdgCode())==11 && tpcOuter && point && cluster && emcalIn)
561 printf("\n\nTPCOutToEMCalInElecCls: ");
562 cout<<cObj->AddTracks(tpcOuter,point,mass,cluster->E(),vPos)<<endl;
566 if(!objName.CompareTo("TPCOutToEMCalInElec"))
568 if(TMath::Abs(particle->GetPdgCode())==11 && tpcOut && emcalIn)
570 printf("TPCOutToEMCalInElec: ");
571 cout<<cObj->AddTracks(tpcOut,emcalIn,mass,charge)<<endl;
575 if(!objName.CompareTo("TPCOutToEMCalInPion"))
577 if(TMath::Abs(particle->GetPdgCode())==211 && tpcOut && emcalIn)
578 cObj->AddTracks(tpcOut,emcalIn,mass,charge);
581 if(!objName.CompareTo("TPCOutToTOFIn"))
584 cObj->AddTracks(tpcOut,tofIn,mass,charge);
587 if(!objName.CompareTo("TPCOutToHMPIDIn"))
589 if(tpcOut && hmpidIn)
590 cObj->AddTracks(tpcOut,hmpidIn,mass,charge);
593 //End of the test for new base tracking class
594 // -----------------------------------------------
605 void AliMCTrackingTestTask::ProcessRefTracker(AliTrackReference* refIn, AliTrackReference* refOut, TParticle*part,Int_t type){
607 // Test propagation from In to out
611 printf("ProcessRefTracker\n");
614 AliExternalTrackParam *param = 0;
615 AliExternalTrackParam *paramMC = 0;
616 AliExternalTrackParam *paramDebug = 0;
618 // Double_t xyzIn[3]={refIn->X(),refIn->Y(), refIn->Z()};
619 Double_t xyzOut[3]={refOut->X(),refOut->Y(), refOut->Z()};
620 Double_t mass = part->GetMass();
623 param=MakeTrack(refIn,part);
624 paramMC=MakeTrack(refOut,part);
625 paramDebug=MakeTrack(refIn,part);
628 if (type!=3) PropagateToPoint(param,xyzOut, mass, step);
634 seed.Set(param->GetX(),param->GetAlpha(),param->GetParameter(),param->GetCovariance());
635 Float_t alpha= TMath::ATan2(refIn->Y(),refIn->X());
636 seed.Rotate(alpha-seed.GetAlpha());
638 for (Float_t xlayer= seed.GetX(); xlayer<refOut->R(); xlayer+=step){
639 seed.PropagateTo(xlayer);
641 seed.PropagateTo(refOut->R());
642 param->Set(seed.GetX(),seed.GetAlpha(),seed.GetParameter(),seed.GetCovariance());
647 TTreeSRedirector *pcstream = GetDebugStreamer();
651 isOK&=param->Rotate(paramMC->GetAlpha());
652 isOK&=param->PropagateTo(paramMC->GetX(),AliTracker::GetBz());
653 param->GetXYZ(gpos.GetMatrixArray());
654 param->GetPxPyPz(gmom.GetMatrixArray());
658 "type="<<type<< // detector matching type
659 "step="<<step<< // propagation step length
660 "refIn.="<<refIn<< // starting track refernce
661 "refOut.="<<refOut<< // outer track reference
662 "p.="<<part<< // particle desription (TParticle)
663 "par.="<<param<< // AliExternalTrackParam create at starting point propagated to outer track ref radius
664 "parMC.="<<paramMC<< // AliExternalTrackParam created at the outer point
665 "gpos.="<<&gpos<< // global position
666 "gmom.="<<&gmom<< // global momenta
676 Bool_t AliMCTrackingTestTask::IsFindable(Int_t label, Float_t minTrackLength ) {
678 // Find findable tracks
681 AliMCParticle *mcParticle = (AliMCParticle*) fMCinfo->GetTrack(label);
682 if(!mcParticle) return kFALSE;
685 Float_t tpcTrackLength = mcParticle->GetTPCTrackLength(AliTracker::GetBz(),0.05,counter,3.0);
686 //printf("tpcTrackLength %f \n", tpcTrackLength);
688 return (tpcTrackLength>minTrackLength);
692 void AliMCTrackingTestTask::FitTrackRefs(TParticle * part, TClonesArray * trefs){
699 printf("FitTrackRefs\n");
702 const Int_t kMinRefs=6;
703 Int_t nrefs = trefs->GetEntries();
704 if (nrefs<kMinRefs) return; // we should have enough references
708 for (Int_t iref=0; iref<nrefs; iref++){
709 AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
711 Float_t dir = ref->X()*ref->Px()+ref->Y()*ref->Py();
713 if (ref->DetectorId()!=AliTrackReference::kTPC) continue;
714 if (iref0<0) iref0 = iref;
717 if (iref1-iref0<kMinRefs) return;
719 for (Int_t icov=0; icov<15; icov++) covar[icov]=0;
726 AliTrackReference * refIn = (AliTrackReference*)trefs->At(iref0);
727 //AliTrackReference * refOut = (AliTrackReference*)trefs->At(iref1);
728 AliExternalTrackParam *paramPropagate= MakeTrack(refIn,part);
729 AliExternalTrackParam *paramUpdate = MakeTrack(refIn,part);
730 paramUpdate->AddCovariance(covar);
731 //Double_t mass = part->GetMass();
732 //Double_t charge = part->GetPDG()->Charge()/3.;
734 Float_t alphaIn= TMath::ATan2(refIn->Y(),refIn->X());
735 Float_t radiusIn= refIn->R();
736 Float_t alphaOut= TMath::ATan2(refOut->Y(),refOut->X());
737 Float_t radiusOut= refOut->R();
741 AliMagF * field = (AliMagF*) TGeoGlobalMagField::Instance()->GetField();
742 for (Int_t iref = iref0; iref<=iref1; iref++){
743 AliTrackReference * ref = (AliTrackReference*)trefs->At(iref);
744 Float_t alphaC= TMath::ATan2(ref->Y(),ref->X());
745 Double_t pos[3] = {ref->X(), ref->Y(), ref->Z()};
747 field->Field(pos,mag);
748 isOKP&=paramPropagate->Rotate(alphaC);
749 isOKU&=paramUpdate->Rotate(alphaC);
750 for (Float_t xref= paramPropagate->GetX(); xref<ref->R(); xref++){
751 isOKP&=paramPropagate->PropagateTo(xref, mag[2]);
752 isOKU&=paramUpdate->PropagateTo(xref, mag[2]);
754 isOKP&=paramPropagate->PropagateTo(ref->R(), mag[2]);
755 isOKU&=paramUpdate->PropagateTo(ref->R(), mag[2]);
756 Double_t clpos[2] = {0, ref->Z()};
757 Double_t clcov[3] = { 0.005,0,0.005};
758 isOKU&= paramUpdate->Update(clpos, clcov);
761 TTreeSRedirector *pcstream = GetDebugStreamer();
767 paramUpdate->GetXYZ(gposU.GetMatrixArray());
768 paramUpdate->GetPxPyPz(gmomU.GetMatrixArray());
769 paramPropagate->GetXYZ(gposP.GetMatrixArray());
770 paramPropagate->GetPxPyPz(gmomP.GetMatrixArray());
772 (*pcstream)<<"MCupdate"<<
780 "pP.="<<paramPropagate<<
781 "pU.="<<paramUpdate<<
789 delete paramPropagate;
796 void AliMCTrackingTestTask::FinishTaskOutput()
799 // According description in AliAnalisysTask this method is call
800 // on the slaves before sending data
803 gSystem->Exec("pwd");
804 RegisterDebugOutput();
809 void AliMCTrackingTestTask::RegisterDebugOutput(){
814 // store - copy debug output to the destination position
815 // currently ONLY for local copy
818 dsName+="Debug.root";
819 dsName.ReplaceAll(" ","");
820 TString dsName2=fDebugOutputPath.Data();
821 gSystem->MakeDirectory(dsName2.Data());
822 dsName2+=gSystem->HostName();
823 gSystem->MakeDirectory(dsName2.Data());
825 dsName2+=gSystem->BaseName(gSystem->pwd());
827 gSystem->MakeDirectory(dsName2.Data());
829 AliInfo(Form("copy %s\t%s\n",dsName.Data(),dsName2.Data()));
830 printf("copy %s\t%s\n",dsName.Data(),dsName2.Data());
831 TFile::Cp(dsName.Data(),dsName2.Data());
836 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
837 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
838 AliXRDPROOFtoolkit tool;
839 TChain * chain = tool.MakeChain("mctracking.txt","MC",0,100);
843 chain->SetAlias("pdg","(p.fPdgCode)");
844 chain->SetAlias("dPRec","(refOut.P()-par.P())/refIn.P()");
845 chain->SetAlias("dPMC","(refOut.P()-refIn.P())/refIn->P()");
846 chain->SetAlias("dPtRec","(refOut.Pt()-par.Pt())/refIn.Pt()");
847 chain->SetAlias("dPtMC","(refOut.Pt()-refIn.Pt())/refIn->Pt()");
851 chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==0&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
852 htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
853 htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
854 gPad->SaveAs("picLoss/dPcorr_ITS_step1.gif");
855 gPad->SaveAs("picLoss/dPcorr_ITS_step1.eps");
857 chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==1&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
858 htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
859 htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
860 gPad->SaveAs("picLoss/dPcorr_TPC_step1.gif");
861 gPad->SaveAs("picLoss/dPcorr_TPC_step1.eps");
864 chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==3&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
865 htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
866 htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
867 gPad->SaveAs("picLoss/dPcorr_TPCseed_step1.gif");
868 gPad->SaveAs("picLoss/dPcorr_TPCseed_step1.eps");
872 chain->Draw("-sqrt(dPRec):-sqrt(dPMC)","abs(pdg)!=11&&type==2&&p.Pt()<0.5&&abs(p.R())<1&&abs(refOut.fZ)<220");
873 htemp->SetYTitle("#sqrt{#DeltaP_{rec}/P}");
874 htemp->SetXTitle("#sqrt{#DeltaP_{mc}/P}");
875 gPad->SaveAs("picLoss/dPcorr_TRD_step1.gif");
876 gPad->SaveAs("picLoss/dPcorr_TRD_step1.eps");
881 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");
882 his->SetXTitle("(P_{trec}-P_{tmc})/P_{tmc}");
883 gPad->SaveAs("picLoss/dPtcorr_TPCseed_step1_1D.eps");
884 gPad->SaveAs("picLoss/dPtcorr_TPCseed_step1_1D.gif");
886 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");
887 his->SetXTitle("(P_{rec}-P_{mc})/P_{mc}");
888 gPad->SaveAs("picLoss/dPcorr_TPCseed_step1_1D.eps");
889 gPad->SaveAs("picLoss/dPcorr_TPCseed_step1_1D.gif");