X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWG1%2FTRD%2FAliTRDresolution.cxx;h=8fb2acbaeb8046738f7ffdba7bcd4dd2579f6f24;hb=62028b10ee4e0b01e6d55d130c1c0a6bb3c33ca6;hp=1e66a94be0b342e59ac92d069c20ea78e0bd4cf2;hpb=afca20ef3f664e102f5b3beb0284757a8cfd84b9;p=u%2Fmrichter%2FAliRoot.git diff --git a/PWG1/TRD/AliTRDresolution.cxx b/PWG1/TRD/AliTRDresolution.cxx index 1e66a94be0b..8fb2acbaeb8 100644 --- a/PWG1/TRD/AliTRDresolution.cxx +++ b/PWG1/TRD/AliTRDresolution.cxx @@ -62,6 +62,7 @@ #include #include #include +#include #include #include #include @@ -70,7 +71,10 @@ #include #include "AliPID.h" +#include "AliLog.h" #include "AliESDtrack.h" +#include "AliMathBase.h" +#include "AliTrackPointArray.h" #include "AliTRDresolution.h" #include "AliTRDgeometry.h" @@ -81,159 +85,94 @@ #include "AliTRDReconstructor.h" #include "AliTRDrecoParam.h" #include "AliTRDpidUtil.h" +#include "AliTRDinfoGen.h" #include "info/AliTRDclusterInfo.h" ClassImp(AliTRDresolution) -UChar_t const AliTRDresolution::fgNhistos[kNviews] = { - 2, 2, 5, 5, - 2, 5, 12, 2, 14 -}; UChar_t const AliTRDresolution::fgNproj[kNviews] = { - 2, 2, 5, 5, - 4, 7, 12, 2, 14 + 2, 2, 5, 5, 5, + 2, 5, 11, 11, 11 }; Char_t const * AliTRDresolution::fgPerformanceName[kNviews] = { "Charge" ,"Cluster2Track" ,"Tracklet2Track" - ,"Tracklet2TPC" + ,"Tracklet2TRDin" + ,"Tracklet2TRDout" ,"Cluster2MC" ,"Tracklet2MC" - ,"TPC2MC" - ,"TOF/HMPID2MC" + ,"TRDin2MC" + ,"TRDout2MC" ,"TRD2MC" }; -UChar_t const AliTRDresolution::fgNcomp[kNprojs] = { - 1, 1, //2, - 1, 1, //2, - 1, 1, 1, 1, 1, //5, - 1, 1, 1, 1, 1, 1, 1, //5, -// MC - 1, 1, 1, 1, //4, - 1, 1, 1, 1, 1, //5, - 1, 1, 1, 1, 1, 1, 1, 1, 11, 11, 11, 11, //12, - 1, 1, //2, - 1, 1, 1, 1, 1, 1, 1, 1, 66, 66, 66, 66, 66, 66 //14 +Char_t const * AliTRDresolution::fgParticle[11]={ + " p bar", " K -", " #pi -", " #mu -", " e -", + " No PID", + " e +", " #mu +", " #pi +", " K +", " p", }; -Char_t const *AliTRDresolution::fgAxTitle[kNprojs][4] = { - // Charge - {"Impv", "x [cm]", "I_{mpv}", "x/x_{0}"} - ,{"dI/Impv", "x/x_{0}", "#delta I/I_{mpv}", "x[cm]"} - // Clusters to Kalman - ,{"Cluster2Track residuals", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"Cluster2Track pulls", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - // TRD tracklet to Kalman fit - ,{"Tracklet2Track Y residuals", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"Tracklet2Track Y pulls", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - ,{"Tracklet2Track Z residuals", "tg(#theta)", "#mu_{z} [#mum]", "#sigma_{z} [#mum]"} - ,{"Tracklet2Track Z pulls", "tg(#theta)", "#mu_{z}", "#sigma_{z}"} - ,{"Tracklet2Track Phi residuals", "tg(#phi)", "#mu_{#phi} [mrad]", "#sigma_{#phi} [mrad]"} - // TPC track 2 first TRD tracklet - ,{"Tracklet2Track Y residuals @ TRDin", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"Tracklet2Track Y pulls @ TRDin", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - ,{"Tracklet2Track Z residuals @ TRDin", "tg(#theta)", "#mu_{z} [#mum]", "#sigma_{z} [#mum]"} - ,{"Tracklet2Track Z pulls @ TRDin", "tg(#theta)", "#mu_{z}", "#sigma_{z}"} - ,{"Tracklet2Track Phi residuals @ TRDin", "tg(#phi)", "#mu_{#phi} [mrad]", "#sigma_{#phi} [mrad]"} - // MC cluster - ,{"MC Cluster Y resolution (p_{t}<1 GeV/c)", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"MC Cluster Y resolution (13 GeV/c)", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"MC Cluster Y pulls", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - // MC tracklet - ,{"MC Tracklet Y resolution (p_{t}<1 GeV/c)", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y} [#mum]"} - ,{"MC Tracklet Y resolution (13 GeV/c)", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y}[#mum]"} - ,{"MC Tracklet Y pulls", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - ,{"MC Tracklet Cross Z resolution", "tg(#theta)", "#mu_{z} [#mum]", "#sigma_{z} [#mum]"} - ,{"MC Tracklet Cross Z pulls", "tg(#theta)", "#mu_{z}", "#sigma_{z}"} - ,{"MC Tracklet Phi resolution", "tg(#phi)", "#mu_{#phi} [mrad]", "#sigma_{#phi} [mrad]"} - // MC track TPC - ,{"Y resolution @ TRDin", "tg(#phi)", "#mu_{y} [#mum]", "#sigma_{y}[#mum]"} - ,{"Y pulls @ TRDin", "tg(#phi)", "#mu_{y}", "#sigma_{y}"} - ,{"Z resolution @ TRDin", "tg(#theta)", "#mu_{z} [#mum]", "#sigma_{z} [#mum]"} - ,{"Z pulls @ TRDin", "tg(#theta)", "#mu_{z}", "#sigma_{z}"} - ,{"Phi resolution @ TRDin", "tg(#phi)", "#mu_{#phi} [mrad]", "#sigma_{#phi} [mrad]"} - ,{"SNP pulls @ TRDin", "tg(#phi)", "#mu_{snp}", "#sigma_{snp}"} - ,{"Theta resolution @ TRDin", "tg(#theta)", "#mu_{#theta} [mrad]", "#sigma_{#theta} [mrad]"} - ,{"TGL pulls @ TRDin", "tg(#theta)", "#mu_{tgl}", "#sigma_{tgl}"} - ,{"P_{t} resolution @ TRDin", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"} - ,{"1/P_{t} pulls @ TRDin", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"} - ,{"P resolution @ TRDin", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"} - ,{"P pulls @ TRDin", "p^{MC} [GeV/c]", "1/p^{REC}-1/p^{MC}", "MC PULL: #sigma^{TPC}(#Deltap/#sigma_{p})"} - // MC track TOF - ,{"PosZ", "tg(#theta)", "MC: #mu_{z}^{TOF} [#mum]", "MC: #sigma_{z}^{TOF} [#mum]"} - ,{"PullsZ", "tg(#theta)", "MC PULL: #mu_{z}^{TOF}", "MC PULL: #sigma_{z}^{TOF}"} - // MC track in TRD - ,{"TRD track MC Y resolution", "tg(#phi)", "#mu_{y}^{Trk} [#mum]", "#sigma_{y}^{Trk} [#mum]"} - ,{"TRD track MC Y pulls", "tg(#phi)", "#mu_{y}^{Trk}", "#sigma_{y}^{Trk}"} - ,{"TRD track MC Z resolution", "tg(#theta)", "#mu_{z}^{Trk} [#mum]", "#sigma_{z}^{Trk} [#mum]"} - ,{"TRD track MC Z pulls", "tg(#theta)", "#mu_{z}^{Trk}", "#sigma_{z}^{Trk}"} - ,{"TRD track MC Phi resolution", "tg(#phi)", "#mu_{#phi}^{Trk} [mrad]", "#sigma_{#phi}^{Trk} [mrad]"} - ,{"TRD track MC SNP pulls", "tg(#phi)", "#mu_{snp}^{Trk}", "#sigma_{snp}^{Trk}"} - ,{"TRD track MC Theta resolution", "tg(#theta)", "#mu_{#theta}^{Trk} [mrad]", "#sigma_{#theta}^{Trk} [mrad]"} - ,{"TRD track MC TGL pulls", "tg(#theta)", "#mu_{tgl}^{Trk}", "#sigma_{tgl}^{Trk}"} - ,{"P_{t} resolution TRD Layer", "p_{t} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "#sigma(#Deltap_{t}/p_{t}^{MC}) [%]"} - ,{"1/P_{t} pulls TRD Layer", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC} - 1/p_{t}^{MC}", "#sigma_{1/p_{t}}"} - ,{"P resolution TRD Layer", "p [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "#sigma(#Deltap/p^{MC}) [%]"} - ,{"[SA] P_{t} resolution TRD Layer", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{Trk}(#Deltap_{t}/p_{t}^{MC}) [%]"} - ,{"[SA] 1/P_{t} pulls TRD Layer", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{Trk}"} - ,{"[SA] P resolution TRD Layer", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{Trk}(#Deltap/p^{MC}) [%]"} + +// Configure segmentation for y resolution/residuals +Int_t const AliTRDresolution::fgkNresYsegm[3] = { + AliTRDgeometry::kNsector + ,AliTRDgeometry::kNsector*AliTRDgeometry::kNstack + ,AliTRDgeometry::kNdet }; +Char_t const *AliTRDresolution::fgkResYsegmName[3] = { + "Sector", "Stack", "Detector"}; + //________________________________________________________ AliTRDresolution::AliTRDresolution() :AliTRDrecoTask() - ,fStatus(0) + ,fSegmentLevel(0) ,fIdxPlot(0) ,fIdxFrame(0) - ,fReconstructor(NULL) - ,fGeo(NULL) + ,fPtThreshold(1.) + ,fDyRange(1.5) ,fDBPDG(NULL) ,fGraphS(NULL) ,fGraphM(NULL) ,fCl(NULL) - ,fTrklt(NULL) ,fMCcl(NULL) - ,fMCtrklt(NULL) +/* ,fTrklt(NULL) + ,fMCtrklt(NULL)*/ { // // Default constructor // SetNameTitle("TRDresolution", "TRD spatial and momentum resolution"); + SetSegmentationLevel(); } //________________________________________________________ AliTRDresolution::AliTRDresolution(char* name) :AliTRDrecoTask(name, "TRD spatial and momentum resolution") - ,fStatus(0) + ,fSegmentLevel(0) ,fIdxPlot(0) ,fIdxFrame(0) - ,fReconstructor(NULL) - ,fGeo(NULL) + ,fPtThreshold(1.) + ,fDyRange(1.5) ,fDBPDG(NULL) ,fGraphS(NULL) ,fGraphM(NULL) ,fCl(NULL) - ,fTrklt(NULL) ,fMCcl(NULL) - ,fMCtrklt(NULL) +/* ,fTrklt(NULL) + ,fMCtrklt(NULL)*/ { // // Default constructor // - fReconstructor = new AliTRDReconstructor(); - fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam()); - fGeo = new AliTRDgeometry(); - InitFunctorList(); + SetSegmentationLevel(); DefineOutput(kClToTrk, TObjArray::Class()); // cluster2track - DefineOutput(kTrkltToTrk, TObjArray::Class()); // tracklet2track DefineOutput(kClToMC, TObjArray::Class()); // cluster2mc - DefineOutput(kTrkltToMC, TObjArray::Class()); // tracklet2mc +/* DefineOutput(kTrkltToTrk, TObjArray::Class()); // tracklet2track + DefineOutput(kTrkltToMC, TObjArray::Class()); // tracklet2mc*/ } //________________________________________________________ @@ -245,13 +184,10 @@ AliTRDresolution::~AliTRDresolution() if(fGraphS){fGraphS->Delete(); delete fGraphS;} if(fGraphM){fGraphM->Delete(); delete fGraphM;} - delete fGeo; - delete fReconstructor; - if(gGeoManager) delete gGeoManager; if(fCl){fCl->Delete(); delete fCl;} - if(fTrklt){fTrklt->Delete(); delete fTrklt;} if(fMCcl){fMCcl->Delete(); delete fMCcl;} - if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;} +/* if(fTrklt){fTrklt->Delete(); delete fTrklt;} + if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;}*/ } @@ -259,17 +195,28 @@ AliTRDresolution::~AliTRDresolution() void AliTRDresolution::UserCreateOutputObjects() { // spatial resolution - OpenFile(1, "RECREATE"); - fContainer = Histos(); - fCl = new TObjArray(); + AliTRDrecoTask::UserCreateOutputObjects(); + InitExchangeContainers(); +} + +//________________________________________________________ +void AliTRDresolution::InitExchangeContainers() +{ +// Init containers for subsequent tasks (AliTRDclusterResolution) + + fCl = new TObjArray(200); fCl->SetOwner(kTRUE); - fTrklt = new TObjArray(); - fTrklt->SetOwner(kTRUE); fMCcl = new TObjArray(); fMCcl->SetOwner(kTRUE); +/* fTrklt = new TObjArray(); + fTrklt->SetOwner(kTRUE); fMCtrklt = new TObjArray(); - fMCtrklt->SetOwner(kTRUE); + fMCtrklt->SetOwner(kTRUE);*/ + PostData(kClToTrk, fCl); + PostData(kClToMC, fMCcl); +/* PostData(kTrkltToTrk, fTrklt); + PostData(kTrkltToMC, fMCtrklt);*/ } //________________________________________________________ @@ -280,14 +227,33 @@ void AliTRDresolution::UserExec(Option_t *opt) // fCl->Delete(); - fTrklt->Delete(); fMCcl->Delete(); - fMCtrklt->Delete(); AliTRDrecoTask::UserExec(opt); - PostData(kClToTrk, fCl); - PostData(kTrkltToTrk, fTrklt); - PostData(kClToMC, fMCcl); - PostData(kTrkltToMC, fMCtrklt); +} + +//________________________________________________________ +Bool_t AliTRDresolution::Pulls(Double_t dyz[2], Double_t cov[3], Double_t tilt) const +{ +// Helper function to calculate pulls in the yz plane +// using proper tilt rotation +// Uses functionality defined by AliTRDseedV1. + + Double_t t2(tilt*tilt); + + // rotate along pad + Double_t cc[3]; + cc[0] = cov[0] - 2.*tilt*cov[1] + t2*cov[2]; + cc[1] = cov[1]*(1.-t2) + tilt*(cov[0] - cov[2]); + cc[2] = t2*cov[0] + 2.*tilt*cov[1] + cov[2]; + // do sqrt + Double_t sqr[3]={0., 0., 0.}; + if(AliTRDseedV1::GetCovSqrt(cc, sqr)) return kFALSE; + Double_t invsqr[3]={0., 0., 0.}; + if(AliTRDseedV1::GetCovInv(sqr, invsqr)<1.e-40) return kFALSE; + Double_t tmp(dyz[0]); + dyz[0] = invsqr[0]*tmp + invsqr[1]*dyz[1]; + dyz[1] = invsqr[1]*tmp + invsqr[2]*dyz[1]; + return kTRUE; } //________________________________________________________ @@ -299,11 +265,11 @@ TH1* AliTRDresolution::PlotCharge(const AliTRDtrackV1 *track) if(track) fkTrack = track; if(!fkTrack){ - AliDebug(2, "No Track defined."); + AliDebug(4, "No Track defined."); return NULL; } TObjArray *arr = NULL; - if(!(arr = ((TObjArray*)fContainer->At(kCharge)))){ + if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCharge)))){ AliWarning("No output container defined."); return NULL; } @@ -315,14 +281,15 @@ TH1* AliTRDresolution::PlotCharge(const AliTRDtrackV1 *track) if(!(fTracklet = fkTrack->GetTracklet(ily))) continue; if(!fTracklet->IsOK()) continue; Float_t x0 = fTracklet->GetX0(); - Float_t dq, dl; + Float_t dqdl, dl; for(Int_t itb=AliTRDseedV1::kNtb; itb--;){ if(!(c = fTracklet->GetClusters(itb))){ if(!(c = fTracklet->GetClusters(AliTRDseedV1::kNtb+itb))) continue; } - dq = fTracklet->GetdQdl(itb, &dl); + dqdl = fTracklet->GetdQdl(itb, &dl); + if(dqdl<1.e-5) continue; dl /= 0.15; // dl/dl0, dl0 = 1.5 mm for nominal vd - (h = (TH3S*)arr->At(0))->Fill(dl, x0-c->GetX(), dq); + (h = (TH3S*)arr->At(0))->Fill(dl, x0-c->GetX(), dqdl); } // if(!HasMCdata()) continue; @@ -344,31 +311,108 @@ TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track) if(track) fkTrack = track; if(!fkTrack){ - AliDebug(2, "No Track defined."); + AliDebug(4, "No Track defined."); return NULL; } TObjArray *arr = NULL; - if(!(arr = ((TObjArray*)fContainer->At(kCluster)))){ + if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCluster)))){ AliWarning("No output container defined."); return NULL; } ULong_t status = fkESD ? fkESD->GetStatus():0; - Double_t covR[7], cov[3]; - Float_t x0, y0, z0, dy, dz, dydx, dzdx; - AliTRDseedV1 *fTracklet(NULL); + Int_t sgm[3]; + Double_t covR[7], cov[3], dy[2], dz[2]; + Float_t pt, x0, y0, z0, dydx, dzdx, tilt(0.); + const AliTRDgeometry *geo(AliTRDinfoGen::Geometry()); + AliTRDseedV1 *fTracklet(NULL); TObjArray *clInfoArr(NULL); + // do LINEAR track refit if asked by the user + // it is the user responsibility to check if B=0T + Float_t param[10]; memset(param, 0, 10*sizeof(Float_t)); + Int_t np(0), nrc(0); AliTrackPoint clusters[300]; + if(HasTrackRefit()){ + Bool_t kPrimary(kFALSE); + for(Int_t ily=0; ilyGetTracklet(ily))) continue; + if(!fTracklet->IsOK()) continue; + x0 = fTracklet->GetX0(); + tilt = fTracklet->GetTilt(); + AliTRDcluster *c = NULL; + fTracklet->ResetClusterIter(kFALSE); + while((c = fTracklet->PrevCluster())){ + Float_t xyz[3] = {c->GetX(), c->GetY(), c->GetZ()}; + clusters[np].SetCharge(tilt); + clusters[np].SetClusterType(0); + clusters[np].SetVolumeID(ily); + clusters[np].SetXYZ(xyz); + np++; + } + if(fTracklet->IsRowCross()){ + Float_t xcross(0.), zcross(0.); + if(fTracklet->GetEstimatedCrossPoint(xcross, zcross)){ + clusters[np].SetCharge(tilt); + clusters[np].SetClusterType(1); + clusters[np].SetVolumeID(ily); + clusters[np].SetXYZ(xcross, 0., zcross); + np++; + nrc++; + } + } + if(fTracklet->IsPrimary()) kPrimary = kTRUE; + } + if(kPrimary){ + clusters[np].SetCharge(tilt); + clusters[np].SetClusterType(1); + clusters[np].SetVolumeID(-1); + clusters[np].SetXYZ(0., 0., 0.); + np++; + } + if(!FitTrack(np, clusters, param)) { + AliDebug(1, "Linear track Fit failed."); + return NULL; + } + if(HasTrackSelection()){ + Bool_t kReject(kFALSE); + if(fkTrack->GetNumberOfTracklets() != AliTRDgeometry::kNlayer) kReject = kTRUE; + if(!kReject && !UseTrack(np, clusters, param)) kReject = kTRUE; + if(kReject){ + AliDebug(1, "Reject track for residuals analysis."); + return NULL; + } + } + } for(Int_t ily=0; ilyGetTracklet(ily))) continue; if(!fTracklet->IsOK()) continue; x0 = fTracklet->GetX0(); - + pt = fTracklet->GetPt(); + sgm[2] = fTracklet->GetDetector(); + sgm[0] = AliTRDgeometry::GetSector(sgm[2]); + sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]); + // retrive the track angle with the chamber - y0 = fTracklet->GetYref(0); - z0 = fTracklet->GetZref(0); - dydx = fTracklet->GetYref(1); - dzdx = fTracklet->GetZref(1); + if(HasTrackRefit()){ + Float_t par[3]; + if(!FitTracklet(ily, np, clusters, param, par)) continue; + dydx = par[2];//param[3]; + dzdx = param[4]; + y0 = par[1] + dydx * (x0 - par[0]);//param[1] + dydx * (x0 - param[0]); + z0 = param[2] + dzdx * (x0 - param[0]); + } else { + y0 = fTracklet->GetYref(0); + z0 = fTracklet->GetZref(0); + dydx = fTracklet->GetYref(1); + dzdx = fTracklet->GetZref(1); + } + /*printf("RC[%c] Primary[%c]\n" + " Fit : y0[%f] z0[%f] dydx[%f] dzdx[%f]\n" + " Ref: y0[%f] z0[%f] dydx[%f] dzdx[%f]\n", fTracklet->IsRowCross()?'y':'n', fTracklet->IsPrimary()?'y':'n', y0, z0, dydx, dzdx + ,fTracklet->GetYref(0),fTracklet->GetZref(0),fTracklet->GetYref(1),fTracklet->GetZref(1));*/ + tilt = fTracklet->GetTilt(); fTracklet->GetCovRef(covR); - Float_t tilt(fTracklet->GetTilt()); + Double_t t2(tilt*tilt) + ,corr(1./(1. + t2)) + ,cost(TMath::Sqrt(corr)); AliTRDcluster *c = NULL; fTracklet->ResetClusterIter(kFALSE); while((c = fTracklet->PrevCluster())){ @@ -378,57 +422,63 @@ TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track) Float_t dx = x0 - xc; Float_t yt = y0 - dx*dydx; Float_t zt = z0 - dx*dzdx; - // calculate residuals using tilt correction -// yc -= tilt*(zc-zt); -// dy = yt - yc; - - // calculate residuals using correct covariance matrix - cov[0] = c->GetSigmaY2(); - cov[1] = c->GetSigmaYZ(); - cov[2] = c->GetSigmaZ2(); - // do rotation - Double_t sy2(cov[0]), sz2(cov[2]); - Double_t t2 = tilt*tilt; - Double_t correction = 1./(1. + t2); - cov[0] = (sy2+t2*sz2)*correction; - cov[1] = tilt*(sz2 - sy2)*correction; - cov[2] = (t2*sy2+sz2)*correction; - // do inversion - Double_t r00=cov[0]+covR[0], r01=cov[1]+covR[1], r11=cov[2]+covR[2]; - Double_t det=r00*r11 - r01*r01; - Double_t tmp=r00; r00=r11/det; r11=tmp/det; - dy = (yc - yt)*TMath::Sqrt(r00); - dz = (zc - zt)*TMath::Sqrt(r11); - - ((TH2I*)arr->At(0))->Fill(dydx, dy/*, dz*/); - ((TH2I*)arr->At(1))->Fill(dydx, dy/TMath::Sqrt(cov[0] /*+ sx2*/ + sy2)); + dy[0] = yc-yt; dz[0]= zc-zt; + + // rotate along pad + dy[1] = cost*(dy[0] - dz[0]*tilt); + dz[1] = cost*(dz[0] + dy[0]*tilt); + if(pt>fPtThreshold && c->IsInChamber()) ((TH3S*)arr->At(0))->Fill(dydx, dy[1], sgm[fSegmentLevel]); + + // tilt rotation of covariance for clusters + Double_t sy2(c->GetSigmaY2()), sz2(c->GetSigmaZ2()); + cov[0] = (sy2+t2*sz2)*corr; + cov[1] = tilt*(sz2 - sy2)*corr; + cov[2] = (t2*sy2+sz2)*corr; + // sum with track covariance + cov[0]+=covR[0]; cov[1]+=covR[1]; cov[2]+=covR[2]; + Double_t dyz[2]= {dy[1], dz[1]}; + Pulls(dyz, cov, tilt); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]); - if(DebugLevel()>=2){ - // Get z-position with respect to anode wire - Int_t istk = fGeo->GetStack(c->GetDetector()); - AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk); - Float_t row0 = pp->GetRow0(); - Float_t d = row0 - zt + pp->GetAnodeWireOffset(); - d -= ((Int_t)(2 * d)) / 2.0; - if (d > 0.25) d = 0.5 - d; - - AliTRDclusterInfo *clInfo = new AliTRDclusterInfo; - fCl->Add(clInfo); - clInfo->SetCluster(c); - Float_t covF[] = {cov[0], cov[1], cov[2]}; - clInfo->SetGlobalPosition(yt, zt, dydx, dzdx, covF); - clInfo->SetResolution(dy); - clInfo->SetAnisochronity(d); - clInfo->SetDriftLength(dx); - clInfo->SetTilt(tilt); - (*DebugStream()) << "ClusterREC" - <<"status=" << status - <<"clInfo.=" << clInfo - << "\n"; + // Get z-position with respect to anode wire + Int_t istk = geo->GetStack(c->GetDetector()); + AliTRDpadPlane *pp = geo->GetPadPlane(ily, istk); + Float_t row0 = pp->GetRow0(); + Float_t d = row0 - zt + pp->GetAnodeWireOffset(); + d -= ((Int_t)(2 * d)) / 2.0; + if (d > 0.25) d = 0.5 - d; + + AliTRDclusterInfo *clInfo(NULL); + clInfo = new AliTRDclusterInfo; + clInfo->SetCluster(c); + Float_t covF[] = {cov[0], cov[1], cov[2]}; + clInfo->SetGlobalPosition(yt, zt, dydx, dzdx, covF); + clInfo->SetResolution(dy[1]); + clInfo->SetAnisochronity(d); + clInfo->SetDriftLength(dx); + clInfo->SetTilt(tilt); + if(fCl) fCl->Add(clInfo); + else AliDebug(1, "Cl exchange container missing. Activate by calling \"InitExchangeContainers()\""); + + if(DebugLevel()>=1){ + if(!clInfoArr){ + clInfoArr=new TObjArray(AliTRDseedV1::kNclusters); + clInfoArr->SetOwner(kFALSE); + } + clInfoArr->Add(clInfo); } } + if(DebugLevel()>=1 && clInfoArr){ + (*DebugStream()) << "cluster" + <<"status=" << status + <<"clInfo.=" << clInfoArr + << "\n"; + clInfoArr->Clear(); + } } - return (TH2I*)arr->At(0); + if(clInfoArr) delete clInfoArr; + + return (TH3S*)arr->At(0); } @@ -445,54 +495,84 @@ TH1* AliTRDresolution::PlotTracklet(const AliTRDtrackV1 *track) // reference position. if(track) fkTrack = track; if(!fkTrack){ - AliDebug(2, "No Track defined."); + AliDebug(4, "No Track defined."); return NULL; } TObjArray *arr = NULL; - if(!(arr = (TObjArray*)fContainer->At(kTrackTRD ))){ + if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrack ))){ AliWarning("No output container defined."); return NULL; } + Int_t sgm[3]; Double_t cov[3], covR[7]/*, sqr[3], inv[3]*/; - Float_t x, dx, dy, dz; - AliTRDseedV1 *fTracklet = NULL; - for(Int_t il=AliTRDgeometry::kNlayer; il--;){ + Double_t pt, phi, tht, x, dx, dy[2], dz[2]; + AliTRDseedV1 *fTracklet(NULL); + for(Int_t il(0); ilGetTracklet(il))) continue; if(!fTracklet->IsOK()) continue; - x = fTracklet->GetX(); - dx = fTracklet->GetX0() - x; - // compute dy^2 and dz^2 - dy = fTracklet->GetYref(0)-dx*fTracklet->GetYref(1) - fTracklet->GetY(); - dz = fTracklet->GetZref(0)-dx*fTracklet->GetZref(1) - fTracklet->GetZ(); + sgm[2] = fTracklet->GetDetector(); + sgm[0] = AliTRDgeometry::GetSector(sgm[2]); + sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]); + x = fTracklet->GetX(); + dx = fTracklet->GetX0() - x; + pt = fTracklet->GetPt(); + phi = fTracklet->GetYref(1); + tht = fTracklet->GetZref(1); + // compute dy and dz + dy[0]= fTracklet->GetYref(0)-dx*fTracklet->GetYref(1) - fTracklet->GetY(); + dz[0]= fTracklet->GetZref(0)-dx*fTracklet->GetZref(1) - fTracklet->GetZ(); + Double_t tilt(fTracklet->GetTilt()) + ,t2(tilt*tilt) + ,corr(1./(1. + t2)) + ,cost(TMath::Sqrt(corr)); + Bool_t rc(fTracklet->IsRowCross()); + + // calculate residuals using tilt rotation + dy[1]= cost*(dy[0] - dz[0]*tilt); + dz[1]= cost*(dz[0] + dy[0]*tilt); + ((TH3S*)arr->At(0))->Fill(phi, dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]); + ((TH3S*)arr->At(2))->Fill(tht, dz[1], rc); + // compute covariance matrix fTracklet->GetCovAt(x, cov); fTracklet->GetCovRef(covR); cov[0] += covR[0]; cov[1] += covR[1]; cov[2] += covR[2]; -/* // Correct PULL calculation by considering off - // diagonal elements in the covariance matrix - // compute square root matrix - if(AliTRDseedV1::GetCovInv(cov, inv)==0.) continue; - if(AliTRDseedV1::GetCovSqrt(inv, sqr)<0.) continue; - Double_t y = sqr[0]*dy+sqr[1]*dz; - Double_t z = sqr[1]*dy+sqr[2]*dz; - ((TH3*)h)->Fill(y, z, fTracklet->GetYref(1));*/ - - ((TH2I*)arr->At(0))->Fill(fTracklet->GetYref(1), dy); - ((TH2I*)arr->At(1))->Fill(fTracklet->GetYref(1), dy/TMath::Sqrt(cov[0])); - ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), TMath::ATan((fTracklet->GetYref(1)-fTracklet->GetYfit(1))/(1-fTracklet->GetYref(1)*fTracklet->GetYfit(1)))); - if(!fTracklet->IsRowCross()) continue; - ((TH2I*)arr->At(2))->Fill(fTracklet->GetZref(1), dz); - ((TH2I*)arr->At(3))->Fill(fTracklet->GetZref(1), dz/TMath::Sqrt(cov[2])); - } + Double_t dyz[2]= {dy[1], dz[1]}; + Pulls(dyz, cov, tilt); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]); + ((TH3S*)arr->At(3))->Fill(tht, dyz[1], rc); + Double_t dphi((phi-fTracklet->GetYfit(1))/(1-phi*fTracklet->GetYfit(1))); + Double_t dtht((tht-fTracklet->GetZfit(1))/(1-tht*fTracklet->GetZfit(1))); + ((TH2I*)arr->At(4))->Fill(phi, TMath::ATan(dphi)); + if(DebugLevel()>=1){ + UChar_t err(fTracklet->GetErrorMsg()); + (*DebugStream()) << "tracklet" + <<"pt=" << pt + <<"phi=" << phi + <<"tht=" << tht + <<"det=" << sgm[2] + <<"dy0=" << dy[0] + <<"dz0=" << dz[0] + <<"dy=" << dy[1] + <<"dz=" << dz[1] + <<"dphi="<< dphi + <<"dtht="<< dtht + <<"dyp=" << dyz[0] + <<"dzp=" << dyz[1] + <<"rc=" << rc + <<"err=" << err + << "\n"; + } + } return (TH2I*)arr->At(0); } //________________________________________________________ -TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) +TH1* AliTRDresolution::PlotTrackIn(const AliTRDtrackV1 *track) { // Store resolution/pulls of Kalman before updating with the TRD information // at the radial position of the first tracklet. The following points are used @@ -505,30 +585,45 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) if(track) fkTrack = track; if(!fkTrack){ - AliDebug(2, "No Track defined."); + AliDebug(4, "No Track defined."); + return NULL; + } + TObjArray *arr = NULL; + if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackIn))){ + AliWarning("No output container defined."); return NULL; } AliExternalTrackParam *tin = NULL; - if(!(tin = fkTrack->GetTrackLow())){ + if(!(tin = fkTrack->GetTrackIn())){ AliWarning("Track did not entered TRD fiducial volume."); return NULL; } TH1 *h = NULL; Double_t x = tin->GetX(); - AliTRDseedV1 *tracklet = NULL; + AliTRDseedV1 *fTracklet = NULL; for(Int_t ily=0; ilyGetTracklet(ily))) continue; + if(!(fTracklet = fkTrack->GetTracklet(ily))) continue; break; } - if(!tracklet || TMath::Abs(x-tracklet->GetX())>1.e-3){ - AliWarning("Tracklet did not match TRD entrance."); + if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){ + AliWarning("Tracklet did not match Track."); return NULL; } + Int_t sgm[3]; + sgm[2] = fTracklet->GetDetector(); + sgm[0] = AliTRDgeometry::GetSector(sgm[2]); + sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]); + Double_t tilt(fTracklet->GetTilt()) + ,t2(tilt*tilt) + ,corr(1./(1. + t2)) + ,cost(TMath::Sqrt(corr)); + Bool_t rc(fTracklet->IsRowCross()); + const Int_t kNPAR(5); Double_t parR[kNPAR]; memcpy(parR, tin->GetParameter(), kNPAR*sizeof(Double_t)); Double_t covR[3*kNPAR]; memcpy(covR, tin->GetCovariance(), 3*kNPAR*sizeof(Double_t)); - Double_t cov[3]; tracklet->GetCovAt(x, cov); + Double_t cov[3]; fTracklet->GetCovAt(x, cov); // define sum covariances TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR); @@ -543,30 +638,37 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) //COV.Print(); PAR.Print(); //TODO Double_t dydx = TMath::Sqrt(1.-parR[2]*parR[2])/parR[2]; - Double_t dy = parR[0] - tracklet->GetY(); - TObjArray *arr = (TObjArray*)fContainer->At(kTrackTPC); - ((TH2I*)arr->At(0))->Fill(tracklet->GetYref(1), dy); - ((TH2I*)arr->At(1))->Fill(tracklet->GetYref(1), dy/TMath::Sqrt(COV(0,0)+cov[0])); - if(tracklet->IsRowCross()){ - Double_t dz = parR[1] - tracklet->GetZ(); - ((TH2I*)arr->At(2))->Fill(tracklet->GetZref(1), dz); - ((TH2I*)arr->At(3))->Fill(tracklet->GetZref(1), dz/TMath::Sqrt(COV(1,1)+cov[2])); - } - Double_t dphi = TMath::ASin(PAR[2])-TMath::ATan(tracklet->GetYfit(1)); ((TH2I*)arr->At(4))->Fill(tracklet->GetYref(1), dphi); + Double_t dy[2]={parR[0] - fTracklet->GetY(), 0.} + ,dz[2]={parR[1] - fTracklet->GetZ(), 0.} + ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1))); + // calculate residuals using tilt rotation + dy[1] = cost*(dy[0] - dz[0]*tilt); + dz[1] = cost*(dz[0] + dy[0]*tilt); + + if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]); + ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc); + ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi); + + Double_t dyz[2] = {dy[1], dz[1]}; + Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]}; + Pulls(dyz, cc, tilt); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]); + ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc); + // register reference histo for mini-task h = (TH2I*)arr->At(0); - if(DebugLevel()>=1){ + if(DebugLevel()>=2){ (*DebugStream()) << "trackIn" << "x=" << x << "P=" << &PAR << "C=" << &COV << "\n"; - Double_t y = tracklet->GetY(); - Double_t z = tracklet->GetZ(); + Double_t y = fTracklet->GetY(); + Double_t z = fTracklet->GetZ(); (*DebugStream()) << "trackletIn" << "y=" << y << "z=" << z @@ -579,7 +681,7 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) if(!HasMCdata()) return h; UChar_t s; - Float_t dx, pt0, x0=tracklet->GetX0(), y0, z0, dydx0, dzdx0; + Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0; if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h; Int_t pdg = fkMC->GetPDG(), sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index @@ -605,13 +707,17 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) // TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T())); // fill histos - arr = (TObjArray*)fContainer->At(kMCtrackTPC); + if(!(arr = (TObjArray*)fContainer->At(kMCtrackIn))) { + AliWarning("No MC container defined."); + return h; + } + // y resolution/pulls - ((TH2I*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0]); - ((TH2I*)arr->At(1))->Fill(dydx0, (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0))); + if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1))); // z resolution/pulls - ((TH2I*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1]); - ((TH2I*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1))); + ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0); + ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0); // phi resolution/snp pulls ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2])); ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2))); @@ -632,7 +738,7 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) // if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx); // fill debug for MC - if(DebugLevel()>=1){ + if(DebugLevel()>=3){ (*DebugStream()) << "trackInMC" << "P=" << &PARMC << "\n"; @@ -640,6 +746,178 @@ TH1* AliTRDresolution::PlotTrackTPC(const AliTRDtrackV1 *track) return h; } +//________________________________________________________ +TH1* AliTRDresolution::PlotTrackOut(const AliTRDtrackV1 *track) +{ +// Store resolution/pulls of Kalman after last update with the TRD information +// at the radial position of the first tracklet. The following points are used +// for comparison +// - the (y,z,snp) of the first TRD tracklet +// - the (y, z, snp, tgl, pt) of the MC track reference +// +// Additionally the momentum resolution/pulls are calculated for usage in the +// PID calculation. + + if(track) fkTrack = track; + if(!fkTrack){ + AliDebug(4, "No Track defined."); + return NULL; + } + TObjArray *arr = NULL; + if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackOut))){ + AliWarning("No output container defined."); + return NULL; + } + AliExternalTrackParam *tout = NULL; + if(!(tout = fkTrack->GetTrackOut())){ + AliDebug(2, "Track did not exit TRD."); + return NULL; + } + TH1 *h(NULL); + + Double_t x = tout->GetX(); + AliTRDseedV1 *fTracklet(NULL); + for(Int_t ily=0; ilyGetTracklet(ily))) continue; + break; + } + if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){ + AliWarning("Tracklet did not match Track position."); + return NULL; + } + Int_t sgm[3]; + sgm[2] = fTracklet->GetDetector(); + sgm[0] = AliTRDgeometry::GetSector(sgm[2]); + sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]); + Double_t tilt(fTracklet->GetTilt()) + ,t2(tilt*tilt) + ,corr(1./(1. + t2)) + ,cost(TMath::Sqrt(corr)); + Bool_t rc(fTracklet->IsRowCross()); + + const Int_t kNPAR(5); + Double_t parR[kNPAR]; memcpy(parR, tout->GetParameter(), kNPAR*sizeof(Double_t)); + Double_t covR[3*kNPAR]; memcpy(covR, tout->GetCovariance(), 3*kNPAR*sizeof(Double_t)); + Double_t cov[3]; fTracklet->GetCovAt(x, cov); + + // define sum covariances + TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR); + Double_t *pc = &covR[0], *pp = &parR[0]; + for(Int_t ir=0; irGetY(), 0., 0.} + ,dz[3]={parR[1] - fTracklet->GetZ(), 0., 0.} + ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1))); + // calculate residuals using tilt rotation + dy[1] = cost*(dy[0] - dz[0]*tilt); + dz[1] = cost*(dz[0] + dy[0]*tilt); + + if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), 1.e2*dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]); // scale to fit general residual range !!! + ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc); + ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi); + + Double_t dyz[2] = {dy[1], dz[1]}; + Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]}; + Pulls(dyz, cc, tilt); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]); + ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc); + + // register reference histo for mini-task + h = (TH2I*)arr->At(0); + + if(DebugLevel()>=2){ + (*DebugStream()) << "trackOut" + << "x=" << x + << "P=" << &PAR + << "C=" << &COV + << "\n"; + + Double_t y = fTracklet->GetY(); + Double_t z = fTracklet->GetZ(); + (*DebugStream()) << "trackletOut" + << "y=" << y + << "z=" << z + << "Vy=" << cov[0] + << "Cyz=" << cov[1] + << "Vz=" << cov[2] + << "\n"; + } + + + if(!HasMCdata()) return h; + UChar_t s; + Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0; + if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h; + Int_t pdg = fkMC->GetPDG(), + sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index + sign(0); + if(!fDBPDG) fDBPDG=TDatabasePDG::Instance(); + TParticlePDG *ppdg(fDBPDG->GetParticle(pdg)); + if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1; + + // translate to reference radial position + dx = x0 - x; y0 -= dx*dydx0; z0 -= dx*dzdx0; + Float_t norm = 1./TMath::Sqrt(1.+dydx0*dydx0); // 1/sqrt(1+tg^2(phi)) + //Fill MC info + TVectorD PARMC(kNPAR); + PARMC[0]=y0; PARMC[1]=z0; + PARMC[2]=dydx0*norm; PARMC[3]=dzdx0*norm; + PARMC[4]=1./pt0; + +// TMatrixDSymEigen eigen(COV); +// TVectorD evals = eigen.GetEigenValues(); +// TMatrixDSym evalsm(kNPAR); +// for(Int_t ir=0; irAt(kMCtrackOut))){ + AliWarning("No MC container defined."); + return h; + } + // y resolution/pulls + if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1))); + // z resolution/pulls + ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0); + ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0); + // phi resolution/snp pulls + ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2])); + ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2))); + // theta resolution/tgl pulls + ((TH2I*)arr->At(6))->Fill(dzdx0, TMath::ATan((PARMC[3]-PAR[3])/(1-PARMC[3]*PAR[3]))); + ((TH2I*)arr->At(7))->Fill(dzdx0, (PARMC[3]-PAR[3])/TMath::Sqrt(COV(3,3))); + // pt resolution\\1/pt pulls\\p resolution/pull + ((TH3S*)arr->At(8))->Fill(pt0, PARMC[4]/PAR[4]-1., sign*sIdx); + ((TH3S*)arr->At(9))->Fill(PARMC[4], (PARMC[4]-PAR[4])/TMath::Sqrt(COV(4,4)), sign*sIdx); + + Double_t p0 = TMath::Sqrt(1.+ PARMC[3]*PARMC[3])*pt0, p; + p = TMath::Sqrt(1.+ PAR[3]*PAR[3])/PAR[4]; + ((TH3S*)arr->At(10))->Fill(p0, p/p0-1., sign*sIdx); +// Float_t sp = +// p*p*PAR[4]*PAR[4]*COV(4,4) +// +2.*PAR[3]*COV(3,4)/PAR[4] +// +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4]; +// if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx); + + // fill debug for MC + if(DebugLevel()>=3){ + (*DebugStream()) << "trackOutMC" + << "P=" << &PARMC + << "\n"; + } + return h; +} + //________________________________________________________ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) { @@ -648,19 +926,23 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) // if(!HasMCdata()){ - AliWarning("No MC defined. Results will not be available."); + AliDebug(2, "No MC defined. Results will not be available."); return NULL; } if(track) fkTrack = track; if(!fkTrack){ - AliDebug(2, "No Track defined."); + AliDebug(4, "No Track defined."); + return NULL; + } + if(!fContainer){ + AliWarning("No output container defined."); return NULL; } // retriev track characteristics Int_t pdg = fkMC->GetPDG(), sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index sign(0), - det(-1), + sgm[3], label(fkMC->GetLabel()); if(!fDBPDG) fDBPDG=TDatabasePDG::Instance(); TParticlePDG *ppdg(fDBPDG->GetParticle(pdg)); @@ -672,38 +954,44 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) Float_t pt0, x0, y0, z0, dx, dy, dz, dydx0, dzdx0; Double_t covR[7]/*, cov[3]*/; - if(DebugLevel()>=1){ - TVectorD dX(12), dY(12), dZ(12), Pt(12), dPt(12), cCOV(12*15); - fkMC->PropagateKalman(&dX, &dY, &dZ, &Pt, &dPt, &cCOV); + if(DebugLevel()>=3){ + TVectorD dX(12), dY(12), dZ(12), vPt(12), dPt(12), cCOV(12*15); + fkMC->PropagateKalman(&dX, &dY, &dZ, &vPt, &dPt, &cCOV); (*DebugStream()) << "MCkalman" << "pdg=" << pdg << "dx=" << &dX << "dy=" << &dY << "dz=" << &dZ - << "pt=" << &Pt + << "pt=" << &vPt << "dpt=" << &dPt << "cov=" << &cCOV << "\n"; } - - AliTRDReconstructor rec; + AliTRDgeometry *geo(AliTRDinfoGen::Geometry()); AliTRDseedV1 *fTracklet(NULL); TObjArray *clInfoArr(NULL); for(Int_t ily=0; ilyGetTracklet(ily)))/* || !fTracklet->IsOK())*/ continue; - det = fTracklet->GetDetector(); + sgm[2] = fTracklet->GetDetector(); + sgm[0] = AliTRDgeometry::GetSector(sgm[2]); + sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]); + Double_t tilt(fTracklet->GetTilt()) + ,t2(tilt*tilt) + ,corr(1./(1. + t2)) + ,cost(TMath::Sqrt(corr)); x0 = fTracklet->GetX0(); //radial shift with respect to the MC reference (radial position of the pad plane) x= fTracklet->GetX(); + Bool_t rc(fTracklet->IsRowCross()); if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) continue; xAnode = fTracklet->GetX0(); // MC track position at reference radial position dx = x0 - x; - if(DebugLevel()>=1){ + if(DebugLevel()>=4){ (*DebugStream()) << "MC" - << "det=" << det + << "det=" << sgm[2] << "pdg=" << pdg << "sgn=" << sign << "pt=" << pt0 @@ -730,13 +1018,13 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) pt = TMath::Abs(fTracklet->GetPt()); fTracklet->GetCovRef(covR); - arr = (TObjArray*)((TObjArray*)fContainer->At(kMCtrackTRD))->At(ily); + arr = (TObjArray*)((TObjArray*)fContainer->At(kMCtrack))->At(ily); // y resolution/pulls - ((TH2I*)arr->At(0))->Fill(dydx0, dy); - ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(covR[0])); + if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(covR[0]), dz/TMath::Sqrt(covR[2])); // z resolution/pulls - ((TH2I*)arr->At(2))->Fill(dzdx0, dz); - ((TH2I*)arr->At(3))->Fill(dzdx0, dz/TMath::Sqrt(covR[2])); + ((TH3S*)arr->At(2))->Fill(dzdx0, dz, 0); + ((TH3S*)arr->At(3))->Fill(dzdx0, dz/TMath::Sqrt(covR[2]), 0); // phi resolution/ snp pulls Double_t dtgp = (dydx - dydx0)/(1.- dydx*dydx0); ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dtgp)); @@ -756,7 +1044,7 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) ((TH3S*)((TObjArray*)arr->At(10)))->Fill(p0, p/p0-1., sign*sIdx); // Fill Debug stream for Kalman track - if(DebugLevel()>=1){ + if(DebugLevel()>=4){ (*DebugStream()) << "MCtrack" << "pt=" << pt << "x=" << x @@ -773,35 +1061,29 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) AliTRDseedV1 tt(*fTracklet); tt.SetZref(0, z0 - (x0-xAnode)*dzdx0); tt.SetZref(1, dzdx0); - tt.SetReconstructor(&rec); - tt.Fit(kTRUE, kTRUE); + tt.SetReconstructor(AliTRDinfoGen::Reconstructor()); + tt.Fit(1); x= tt.GetX();y= tt.GetY();z= tt.GetZ(); dydx = tt.GetYfit(1); dx = x0 - x; ymc = y0 - dx*dydx0; zmc = z0 - dx*dzdx0; - Bool_t rc = tt.IsRowCross(); - + dy = y-ymc; + dz = z-zmc; + Float_t dphi = (dydx - dydx0); + dphi /= (1.- dydx*dydx0); + // add tracklet residuals for y and dydx arr = (TObjArray*)fContainer->At(kMCtracklet); - if(!rc){ - dy = y-ymc; - Float_t dphi = (dydx - dydx0); - dphi /= (1.- dydx*dydx0); - - ((TH3S*)arr->At(0))->Fill(dydx0, dy, pt0); - if(tt.GetS2Y()>0.) ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(tt.GetS2Y())); - ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dphi)); - } else { - // add tracklet residuals for z - dz = z-zmc; - ((TH2I*)arr->At(2))->Fill(dzdl0, dz); - if(tt.GetS2Z()>0.) ((TH2I*)arr->At(3))->Fill(dzdl0, dz/TMath::Sqrt(tt.GetS2Z())); - } + if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]); + if(tt.GetS2Y()>0. && tt.GetS2Z()>0.) ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(tt.GetS2Y()), dz/TMath::Sqrt(tt.GetS2Z())); + ((TH3S*)arr->At(2))->Fill(dzdl0, dz, rc); + if(tt.GetS2Z()>0.) ((TH3S*)arr->At(3))->Fill(dzdl0, dz/TMath::Sqrt(tt.GetS2Z()), rc); + ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dphi)); // Fill Debug stream for tracklet - if(DebugLevel()>=1){ + if(DebugLevel()>=4){ Float_t s2y = tt.GetS2Y(); Float_t s2z = tt.GetS2Z(); (*DebugStream()) << "MCtracklet" @@ -815,26 +1097,26 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) << "\n"; } - Int_t istk = AliTRDgeometry::GetStack(det); - AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk); + AliTRDpadPlane *pp = geo->GetPadPlane(ily, AliTRDgeometry::GetStack(sgm[2])); Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset(); - Float_t tilt = fTracklet->GetTilt(); //Double_t exb = AliTRDCommonParam::Instance()->GetOmegaTau(1.5); arr = (TObjArray*)fContainer->At(kMCcluster); AliTRDcluster *c = NULL; - fTracklet->ResetClusterIter(kFALSE); - while((c = fTracklet->PrevCluster())){ + tt.ResetClusterIter(kFALSE); + while((c = tt.PrevCluster())){ Float_t q = TMath::Abs(c->GetQ()); x = c->GetX(); y = c->GetY();z = c->GetZ(); dx = x0 - x; ymc= y0 - dx*dydx0; zmc= z0 - dx*dzdx0; - dy = (y - tilt*(z-zmc)) - ymc; + dy = cost*(y - ymc - tilt*(z-zmc)); + dz = cost*(z - zmc + tilt*(y-ymc)); + // Fill Histograms - if(q>20. && q<250.){ - ((TH3S*)arr->At(0))->Fill(dydx0, dy, pt0); - ((TH2I*)arr->At(1))->Fill(dydx0, dy/TMath::Sqrt(c->GetSigmaY2())); + if(q>20. && q<250. && pt0>fPtThreshold && c->IsInChamber()){ + ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]); + ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(c->GetSigmaY2()), dz/TMath::Sqrt(c->GetSigmaZ2())); } // Fill calibration container @@ -847,16 +1129,20 @@ TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track) clInfo->SetGlobalPosition(ymc, zmc, dydx0, dzdx0); clInfo->SetResolution(dy); clInfo->SetAnisochronity(d); - clInfo->SetDriftLength(dx-.5*AliTRDgeometry::CamHght()); + clInfo->SetDriftLength(dx); clInfo->SetTilt(tilt); - fMCcl->Add(clInfo); - if(DebugLevel()>=2){ - if(!clInfoArr) clInfoArr=new TObjArray(AliTRDseedV1::kNclusters); + if(fMCcl) fMCcl->Add(clInfo); + else AliDebug(1, "MCcl exchange container missing. Activate by calling \"InitExchangeContainers()\""); + if(DebugLevel()>=5){ + if(!clInfoArr){ + clInfoArr=new TObjArray(AliTRDseedV1::kNclusters); + clInfoArr->SetOwner(kFALSE); + } clInfoArr->Add(clInfo); } } // Fill Debug Tree - if(DebugLevel()>=2 && clInfoArr){ + if(DebugLevel()>=5 && clInfoArr){ (*DebugStream()) << "MCcluster" <<"clInfo.=" << clInfoArr << "\n"; @@ -880,365 +1166,961 @@ Bool_t AliTRDresolution::GetRefFigure(Int_t ifig) AliWarning("Please provide a canvas to draw results."); return kFALSE; } - Int_t first(2), // first particle species to be drawn - nspecies(7); // last particle species to be drawn - TList *l = NULL; TVirtualPad *pad=NULL; + Int_t selection[100], n(0), selStart(0); // + Int_t ly0(0), dly(5); + //Int_t ly0(1), dly(2); // used for SA + TList *l(NULL); TVirtualPad *pad(NULL); + TGraphErrors *g(NULL);TGraphAsymmErrors *ga(NULL); switch(ifig){ - case kCharge: + case 0: // charge resolution gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); ((TVirtualPad*)l->At(0))->cd(); - ((TGraphAsymmErrors*)((TObjArray*)fGraphM->At(kCharge))->At(0))->Draw("apl"); + ga=((TGraphAsymmErrors*)((TObjArray*)fGraphM->At(kCharge))->At(0)); + if(ga->GetN()) ga->Draw("apl"); ((TVirtualPad*)l->At(1))->cd(); - ((TGraphErrors*)((TObjArray*)fGraphS->At(kCharge))->At(0))->Draw("apl"); + g = ((TGraphErrors*)((TObjArray*)fGraphS->At(kCharge))->At(0)); + if(g->GetN()) g->Draw("apl"); break; - case kCluster: + case 1: // cluster2track residuals + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.; + pad = (TVirtualPad*)l->At(0); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + selStart=0; for(n=0; nAt(1); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0] = -.3; xy[1] = -200.; xy[2] = .3; xy[3] = 1000.; + xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.; + pad = (TVirtualPad*)l->At(0); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + if(!GetGraphArray(xy, kCluster, 1, 1)) break; + return kTRUE; + case 3: // kTrack y + gPad->Divide(3, 2, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0] = -.3; xy[1] = -20.; xy[2] = .3; xy[3] = 100.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kCluster, 0)) break; - xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5; + selStart=0; for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kCluster, 1)) break; + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(2))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(3))->cd(); + selStart=fgkNresYsegm[fSegmentLevel]; for(n=0; nAt(4))->cd(); + selStart=fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; nAt(5))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0] = -.3; xy[1] = -500.; xy[2] = .3; xy[3] = 1500.; + + xy[0] = -1.; xy[1] = -150.; xy[2] = 1.; xy[3] = 1000.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kTrackTRD , 0)) break; - xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5; + selection[0]=1; + if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break; + + xy[0] = -1.; xy[1] = -1500.; xy[2] = 1.; xy[3] = 10000.; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kTrackTRD , 1)) break; + selection[0]=0; + if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break; + return kTRUE; - case 3: // kTrackTRD z + case 5: // kTrack pulls gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.; + + xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kTrackTRD , 2)) break; + if(!GetGraphArray(xy, kTrack, 1, 1)) break; + xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kTrackTRD , 3)) break; + if(!GetGraphArray(xy, kTrack, 3, 1)) break; return kTRUE; - case 4: // kTrackTRD phi + case 6: // kTrack phi xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.; - if(GetGraphPlot(&xy[0], kTrackTRD , 4)) return kTRUE; + if(GetGraph(&xy[0], kTrack , 4)) return kTRUE; break; - case 5: // kTrackTPC y + case 7: // kTrackIn y + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.; + pad = ((TVirtualPad*)l->At(0)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + selStart=0; for(n=0; nAt(1)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0] = -.3; xy[1] = -500.; xy[2] = .3; xy[3] = 1500.; + xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.; pad = ((TVirtualPad*)l->At(0)); pad->cd(); pad->SetMargin(0.1, 0.1, 0.1, 0.01); - if(!GetGraphPlot(&xy[0], kTrackTPC, 0)) break; - xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5; + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1)); pad->cd(); pad->SetMargin(0.1, 0.1, 0.1, 0.01); - if(!GetGraphPlot(&xy[0], kTrackTPC, 1)) break; + if(!GetGraphArray(xy, kTrackIn, 1, 1)) break; return kTRUE; - case 6: // kTrackTPC z + case 9: // kTrackIn z gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.; pad = ((TVirtualPad*)l->At(0)); pad->cd(); pad->SetMargin(0.1, 0.1, 0.1, 0.01); - if(!GetGraphPlot(&xy[0], kTrackTPC, 2)) break; + selection[0]=1; + if(!GetGraphArray(xy, kTrackIn, 2, 1, 1, selection)) break; xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; pad = ((TVirtualPad*)l->At(1)); pad->cd(); pad->SetMargin(0.1, 0.1, 0.1, 0.01); - if(!GetGraphPlot(&xy[0], kTrackTPC, 3)) break; + if(!GetGraphArray(xy, kTrackIn, 3, 1)) break; return kTRUE; - case 7: // kTrackTPC phi + case 10: // kTrackIn phi xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.; - if(GetGraphPlot(&xy[0], kTrackTPC, 4)) return kTRUE; + if(GetGraph(&xy[0], kTrackIn, 4)) return kTRUE; break; - case 8: // kMCcluster pt <2 GeV/c - gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.; - ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCcluster, 0)) break; - ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCcluster, 1)) break; - return kTRUE; - case 9: // kMCcluster pt > 3 GeV/c + case 11: // kTrackOut y gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.; - ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCcluster, 2)) break; - xy[0] = -.3; xy[1] = -0.5; xy[2] = .3; xy[3] = 2.5; - ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCcluster, 3)) break; + xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.; + pad = ((TVirtualPad*)l->At(0)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + selStart=0; for(n=0; nAt(1)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =250.; - ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 0)) break; - ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 1)) break; + xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.; + pad = ((TVirtualPad*)l->At(0)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + if(!GetGraphArray(xy, kTrackOut, 1, 1)) break; return kTRUE; - case 11: //kMCtracklet [y] pt > 3 GeV/c + case 13: // kTrackOut z gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =250.; - ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 2)) break; - xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5; + xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.; + pad = ((TVirtualPad*)l->At(0)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + if(!GetGraphArray(xy, kTrackOut, 2, 1)) break; + xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; + pad = ((TVirtualPad*)l->At(1)); pad->cd(); + pad->SetMargin(0.1, 0.1, 0.1, 0.01); + if(!GetGraphArray(xy, kTrackOut, 3, 1)) break; + return kTRUE; + case 14: // kTrackOut phi + xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.; + if(GetGraph(&xy[0], kTrackOut, 4)) return kTRUE; + break; + case 15: // kMCcluster + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=0; for(n=0; nAt(1))->cd(); + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1))->cd(); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + if(!GetGraphArray(xy, kMCcluster, 1, 1)) break; + return kTRUE; + case 17: //kMCtracklet [y] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=0; for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 3)) break; + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1))->cd(); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + if(!GetGraphArray(xy, kMCtracklet, 1, 1)) break; + return kTRUE; + case 19: //kMCtracklet [z] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-1.; xy[1]=-100.; xy[2]=1.; xy[3] =2500.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 4)) break; + if(!GetGraphArray(xy, kMCtracklet, 2)) break; xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtracklet, 5)) break; + if(!GetGraphArray(xy, kMCtracklet, 3)) break; return kTRUE; - case 13: //kMCtracklet [phi] + case 20: //kMCtracklet [phi] xy[0]=-.3; xy[1]=-3.; xy[2]=.3; xy[3] =25.; - if(!GetGraphPlot(&xy[0], kMCtracklet, 6)) break; + if(!GetGraph(&xy[0], kMCtracklet, 4)) break; return kTRUE; - case 14: //kMCtrackTRD [y] + case 21: //kMCtrack [y] ly [0] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 0)) break; - xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 3.5; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.); for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 1)) break; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-1.; xy[1]=-700.; xy[2]=1.; xy[3] =1500.; + xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 2)) break; - xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.); for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 3)) break; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.2; xy[1]=-0.2; xy[2]=.2; xy[3] =2.; + xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 4)) break; - xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.); for(n=0; nAt(1))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.); for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 5)) break; + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.); for(n=0; nAt(1))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.); for(n=0; nAt(1))->cd(); + selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.5); for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 5.5; + ((TVirtualPad*)l->At(0))->cd(); + selStart=0; for(n=0; n<6; n++) selection[n]=selStart+n; + if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break; + ((TVirtualPad*)l->At(1))->cd(); + selStart=6; for(n=0; n<6; n++) selection[n]=selStart+n; + if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break; + return kTRUE; + case 28: //kMCtrack [z] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-1.; xy[1]=-1500.; xy[2]=1.; xy[3] =6000.; + ((TVirtualPad*)l->At(0))->cd(); + if(!GetGraphArray(xy, kMCtrack, 2)) break; + xy[0] = -1.; xy[1] = -1.5; xy[2] = 1.; xy[3] = 5.; + ((TVirtualPad*)l->At(1))->cd(); + if(!GetGraphArray(xy, kMCtrack, 3)) break; + return kTRUE; + case 29: //kMCtrack [phi/snp] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.2; xy[1]=-0.5; xy[2]=.2; xy[3] =10.; + ((TVirtualPad*)l->At(0))->cd(); + if(!GetGraphArray(xy, kMCtrack, 4)) break; + xy[0] = -.2; xy[1] = -1.5; xy[2] = .2; xy[3] = 5.; + ((TVirtualPad*)l->At(1))->cd(); + if(!GetGraphArray(xy, kMCtrack, 5)) break; + return kTRUE; + case 30: //kMCtrack [theta/tgl] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-1.; xy[1]=-0.5; xy[2]=1.; xy[3] =5.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 6)) break; + if(!GetGraphArray(xy, kMCtrack, 6)) break; xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTRD, 7)) break; + if(!GetGraphArray(xy, kMCtrack, 7)) break; return kTRUE; - case 18: //kMCtrackTRD [pt] - xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.; + case 31: //kMCtrack [pt] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 8, first, nspecies)) break; + // pi selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 8, 55+first, nspecies, kTRUE)) break; + // mu selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); return kTRUE; - case 19: //kMCtrackTRD [1/pt] pulls - xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5; + case 32: //kMCtrack [pt] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 9, first, nspecies)) break; + // p selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 9, 55+first, nspecies, kTRUE)) break; + // e selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); return kTRUE; - case 20: //kMCtrackTRD [p] - xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.; + case 33: //kMCtrack [1/pt] pulls + xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5; + //xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 4.5; // SA gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 10, first, nspecies)) break; - pad->Modified(); pad->Update(); pad->SetLogx(); + // pi selection + n=0; + for(Int_t il(ly0); ilAt(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 10, 55+first, nspecies, kTRUE)) break; - pad->Modified(); pad->Update(); pad->SetLogx(); + // mu selection + n=0; + for(Int_t il(ly0); ilDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 11, first, nspecies)) break; - pad->SetLogx(); + // p selection + n=0; + for(Int_t il(ly0); ilAt(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 11, 55+first, nspecies, kTRUE)) break; - pad->SetLogx(); + // e selection + n=0; + for(Int_t il(ly0); ilDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 12, first, nspecies)) break; + // pi selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 12, 55+first, nspecies, kTRUE)) break; + // mu selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); return kTRUE; - case 23: //kMCtrackTRD - SA [p] - xy[0] = 0.; xy[1] = -7.5; xy[2] = 12.; xy[3] = 10.5; + case 36: //kMCtrack [p] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); pad = (TVirtualPad*)l->At(0); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 13, first, nspecies)) break; - pad->SetLogx(); + // p selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrack(&xy[0], 13, 55+first, nspecies, kTRUE)) break; - pad->SetLogx(); + // e selection + n=0; + for(Int_t il(ly0); ilModified(); pad->Update(); pad->SetLogx(); return kTRUE; - case 24: // kMCtrackTPC [y] + case 37: // kMCtrackIn [y] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); - xy[0]=-.25; xy[1]=-50.; xy[2]=.25; xy[3] =800.; + xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 0)) break; - xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 2.5; + selStart=0; for(n=0; nAt(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 1)) break; + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1))->cd(); + if(!GetGraphArray(xy, kMCtrackIn, 1, 1)) break; + return kTRUE; + case 39: // kMCtrackIn [z] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =800.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 2)) break; + if(!GetGraphArray(xy, kMCtrackIn, 2, 1)) break; xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 3)) break; + if(!GetGraphArray(xy, kMCtrackIn, 3, 1)) break; return kTRUE; - case 26: // kMCtrackTPC [phi|snp] + case 40: // kMCtrackIn [phi|snp] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 4)) break; + if(!GetGraph(&xy[0], kMCtrackIn, 4)) break; xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 5)) break; + if(!GetGraph(&xy[0], kMCtrackIn, 5)) break; return kTRUE; - case 27: // kMCtrackTPC [theta|tgl] + case 41: // kMCtrackIn [theta|tgl] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.; ((TVirtualPad*)l->At(0))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 6)) break; + if(!GetGraph(&xy[0], kMCtrackIn, 6)) break; xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 1.5; ((TVirtualPad*)l->At(1))->cd(); - if(!GetGraphPlot(&xy[0], kMCtrackTPC, 7)) break; + if(!GetGraph(&xy[0], kMCtrackIn, 7)) break; return kTRUE; - case 28: // kMCtrackTPC [pt] + case 42: // kMCtrackIn [pt] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.; - pad=(TVirtualPad*)l->At(0); pad->cd(); + //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.; // SA + pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break; + pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetLogx(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrackTPC(xy, 8, first, nspecies)) break; - pad->SetLogx(); - xy[0]=0.; xy[1]=-0.5; xy[2]=2.; xy[3] =3.; + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break; + return kTRUE; + case 43: //kMCtrackIn [1/pt] pulls + xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5; + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + pad = (TVirtualPad*)l->At(0); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break; pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrackTPC(xy, 9, first, nspecies, kTRUE)) break; + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break; + return kTRUE; + case 44: // kMCtrackIn [p] + xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.; + //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.; + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackIn, 10, 1, n, selection)) break; + pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackIn, 10, 1, n, selection)) break; + return kTRUE; + case 45: // kMCtrackOut [y] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=0; for(n=0; nAt(1))->cd(); + selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; nDivide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.; + ((TVirtualPad*)l->At(0))->cd(); + selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; nAt(1))->cd(); + if(!GetGraphArray(xy, kMCtrackOut, 1, 1)) break; + return kTRUE; + case 47: // kMCtrackOut [z] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =1500.; + ((TVirtualPad*)l->At(0))->cd(); + if(!GetGraphArray(xy, kMCtrackOut, 2, 1)) break; + xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5; + ((TVirtualPad*)l->At(1))->cd(); + if(!GetGraphArray(xy, kMCtrackOut, 3, 1)) break; + return kTRUE; + case 48: // kMCtrackOut [phi|snp] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5; + ((TVirtualPad*)l->At(0))->cd(); + if(!GetGraph(&xy[0], kMCtrackOut, 4)) break; + xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5; + ((TVirtualPad*)l->At(1))->cd(); + if(!GetGraph(&xy[0], kMCtrackOut, 5)) break; + return kTRUE; + case 49: // kMCtrackOut [theta|tgl] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.; + ((TVirtualPad*)l->At(0))->cd(); + if(!GetGraph(&xy[0], kMCtrackOut, 6)) break; + xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 15.; + ((TVirtualPad*)l->At(1))->cd(); + if(!GetGraph(&xy[0], kMCtrackOut, 7)) break; return kTRUE; - case 29: // kMCtrackTPC [p] + case 50: // kMCtrackOut [pt] gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.; - pad = ((TVirtualPad*)l->At(0));pad->cd(); + pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrackTPC(xy, 10, first, nspecies)) break; - pad->SetLogx(); - xy[0]=0.; xy[1]=-0.5; xy[2]=8.; xy[3] =2.5; - pad = ((TVirtualPad*)l->At(1)); pad->cd(); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break; + pad = (TVirtualPad*)l->At(1); pad->cd();pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break; + return kTRUE; + case 51: //kMCtrackOut [1/pt] pulls + xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5; + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + pad = (TVirtualPad*)l->At(0); pad->cd(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break; + pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetMargin(0.125, 0.015, 0.1, 0.015); - if(!GetGraphTrackTPC(xy, 11, first, nspecies, kTRUE)) break; + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break; return kTRUE; - case 30: // kMCtrackTOF [z] + case 52: // kMCtrackOut [p] + gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives(); + xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.; + pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8; + if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break; + pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx(); + pad->SetMargin(0.125, 0.015, 0.1, 0.015); + n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10; + if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break; return kTRUE; } AliWarning(Form("Reference plot [%d] missing result", ifig)); return kFALSE; } -Char_t const *fgParticle[11]={ - " p bar", " K -", " #pi -", " #mu -", " e -", - " No PID", - " e +", " #mu +", " #pi +", " K +", " p", -}; -const Color_t fgColorS[11]={ -kOrange, kOrange-3, kMagenta+1, kViolet, kRed, -kGray, -kRed, kViolet, kMagenta+1, kOrange-3, kOrange -}; -const Color_t fgColorM[11]={ -kCyan-5, kAzure-4, kBlue-7, kBlue+2, kViolet+10, -kBlack, -kViolet+10, kBlue+2, kBlue-7, kAzure-4, kCyan-5 -}; -const Marker_t fgMarker[11]={ -30, 30, 26, 25, 24, -28, -20, 21, 22, 29, 29 -}; +//________________________________________________________ +void AliTRDresolution::MakeSummary() +{ +// Build summary plots + + if(!fGraphS || !fGraphM){ + AliError("Missing results"); + return; + } + Float_t xy[4] = {0., 0., 0., 0.}; + Float_t range[2]; + TH2 *h2 = new TH2I("h2SF", "", 20, -.2, .2, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5); + h2->GetXaxis()->CenterTitle(); + h2->GetYaxis()->CenterTitle(); + h2->GetZaxis()->CenterTitle();h2->GetZaxis()->SetTitleOffset(1.4); + + Int_t ih2(0), iSumPlot(0); + TCanvas *cOut = new TCanvas(Form("TRDsummary%s_%d", GetName(), iSumPlot++), "Cluster & Tracklet Resolution", 1024, 768); + cOut->Divide(3,2, 2.e-3, 2.e-3, kYellow-7); + TVirtualPad *p(NULL); + + p=cOut->cd(1); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetTitle(Form("Cluster-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kCluster))->At(0), h2); + GetRange(h2, 1, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(2); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetTitle(Form("Cluster-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kCluster))->At(0), h2); + GetRange(h2, 0, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(3); + p->SetRightMargin(0.06);p->SetTopMargin(0.06); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + GetGraphArray(xy, kCluster, 1, 1); + + p=cOut->cd(4); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetTitle(Form("Tracklet-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kTrack))->At(0), h2); + GetRange(h2, 1, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(5); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetTitle(Form("Tracklet-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kTrack))->At(0), h2); + GetRange(h2, 0, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(6); + p->SetRightMargin(0.06);p->SetTopMargin(0.06); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + GetGraphArray(xy, kTrack, 1, 1); + + cOut->SaveAs(Form("%s.gif", cOut->GetName())); + + if(!HasMCdata() || + (!fGraphS->At(kMCcluster) || !fGraphM->At(kMCcluster) || + !fGraphS->At(kMCtracklet) || !fGraphM->At(kMCtracklet))){ + delete cOut; + return; + } + cOut->Clear(); cOut->SetName(Form("TRDsummary%s_%d", GetName(), iSumPlot++)); + cOut->Divide(3, 2, 2.e-3, 2.e-3, kBlue-10); + + p=cOut->cd(1); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetTitle(Form("Cluster-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kMCcluster))->At(0), h2); + GetRange(h2, 1, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(2); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetContour(7); + h2->SetTitle(Form("Cluster-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kMCcluster))->At(0), h2); + GetRange(h2, 0, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(3); + p->SetRightMargin(0.06);p->SetTopMargin(0.06); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + GetGraphArray(xy, kMCcluster, 1, 1); + + p=cOut->cd(4); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetContour(7); + h2->SetTitle(Form("Tracklet-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kMCtracklet))->At(0), h2); + GetRange(h2, 1, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(5); + p->SetRightMargin(0.16);p->SetTopMargin(0.06); + h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++)); + h2->SetContour(7); + h2->SetTitle(Form("Tracklet-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel])); + MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kMCtracklet))->At(0), h2); + GetRange(h2, 0, range); + h2->GetZaxis()->SetRangeUser(range[0], range[1]); + h2->Draw("colz"); + h2->SetContour(7); + + p=cOut->cd(6); + p->SetRightMargin(0.06);p->SetTopMargin(0.06); + xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5; + GetGraphArray(xy, kMCtracklet, 1, 1); + + cOut->SaveAs(Form("%s.gif", cOut->GetName())); + delete cOut; +} + +//________________________________________________________ +void AliTRDresolution::GetRange(TH2 *h2, Char_t mod, Float_t *range) +{ +// Returns the range of the bulk of data in histogram h2. Removes outliers. +// The "range" vector should be initialized with 2 elements +// Option "mod" can be any of +// - 0 : gaussian like distribution +// - 1 : tailed distribution + + Int_t nx(h2->GetNbinsX()) + , ny(h2->GetNbinsY()) + , n(nx*ny); + Double_t *data=new Double_t[n]; + for(Int_t ix(1), in(0); ix<=nx; ix++){ + for(Int_t iy(1); iy<=ny; iy++) + data[in++] = h2->GetBinContent(ix, iy); + } + Double_t mean, sigm; + AliMathBase::EvaluateUni(n, data, mean, sigm, Int_t(n*.8)); + + range[0]=mean-3.*sigm; range[1]=mean+3.*sigm; + if(mod==1) range[0]=TMath::Max(Float_t(1.e-3), range[0]); + AliDebug(2, Form("h[%s] range0[%f %f]", h2->GetName(), range[0], range[1])); + TH1S h1("h1SF0", "", 100, range[0], range[1]); + h1.FillN(n,data,0); + delete [] data; + + switch(mod){ + case 0:// gaussian distribution + { + TF1 fg("fg", "gaus", mean-3.*sigm, mean+3.*sigm); + h1.Fit(&fg, "QN"); + mean=fg.GetParameter(1); sigm=fg.GetParameter(2); + range[0] = mean-2.5*sigm;range[1] = mean+2.5*sigm; + AliDebug(2, Form(" rangeG[%f %f]", range[0], range[1])); + break; + } + case 1:// tailed distribution + { + Int_t bmax(h1.GetMaximumBin()); + Int_t jBinMin(1), jBinMax(100); + for(Int_t ibin(bmax); ibin--;){ + if(h1.GetBinContent(ibin)<1.){ + jBinMin=ibin; break; + } + } + for(Int_t ibin(bmax); ibin++;){ + if(h1.GetBinContent(ibin)<1.){ + jBinMax=ibin; break; + } + } + range[0]=h1.GetBinCenter(jBinMin); range[1]=h1.GetBinCenter(jBinMax); + AliDebug(2, Form(" rangeT[%f %f]", range[0], range[1])); + break; + } + } + + return; +} + +//________________________________________________________ +void AliTRDresolution::MakeSummaryPlot(TObjArray *a, TH2 *h2) +{ +// Core functionality for MakeSummary function. + + h2->Reset(); + Double_t x,y; + TGraphErrors *g(NULL); TAxis *ax(h2->GetXaxis()); + for(Int_t iseg(0); isegAt(iseg); + for(Int_t in(0); inGetN(); in++){ + g->GetPoint(in, x, y); + h2->SetBinContent(ax->FindBin(x), iseg+1, y); + } + } +} + + //________________________________________________________ Bool_t AliTRDresolution::PostProcess() { - //fContainer = dynamic_cast(GetOutputData(0)); +// Fit, Project, Combine, Extract values from the containers filled during execution + + /*fContainer = dynamic_cast(GetOutputData(0));*/ if (!fContainer) { AliError("ERROR: list not available"); return kFALSE; } - Int_t nc(0); + + // define general behavior parameters + const Color_t fgColorS[11]={ + kOrange, kOrange-3, kMagenta+1, kViolet, kRed, + kGray, + kRed, kViolet, kMagenta+1, kOrange-3, kOrange + }; + const Color_t fgColorM[11]={ + kCyan-5, kAzure-4, kBlue-7, kBlue+2, kViolet+10, + kBlack, + kViolet+10, kBlue+2, kBlue-7, kAzure-4, kCyan-5 + }; + const Marker_t fgMarker[11]={ + 30, 30, 26, 25, 24, + 28, + 20, 21, 22, 29, 29 + }; + TGraph *gm= NULL, *gs= NULL; if(!fGraphS && !fGraphM){ TObjArray *aM(NULL), *aS(NULL); Int_t n = fContainer->GetEntriesFast(); fGraphS = new TObjArray(n); fGraphS->SetOwner(); fGraphM = new TObjArray(n); fGraphM->SetOwner(); - for(Int_t ig=0; igAddAt(aM = new TObjArray(fgNproj[ig]), ig); fGraphS->AddAt(aS = new TObjArray(fgNproj[ig]), ig); - - for(Int_t ic=0; ic1){ + + for(Int_t ic=0; ic1){ TObjArray *agS(NULL), *agM(NULL); - aS->AddAt(agS = new TObjArray(fgNcomp[nc]), ic); - aM->AddAt(agM = new TObjArray(fgNcomp[nc]), ic); - for(Int_t is=fgNcomp[nc]; is--;){ + aS->AddAt(agS = new TObjArray(fNcomp[nc]), ic); + aM->AddAt(agM = new TObjArray(fNcomp[nc]), ic); + for(Int_t is=fNcomp[nc]; is--;){ agS->AddAt(gs = new TGraphErrors(), is); - Int_t is0(is%11); + Int_t is0(is%11), il0(is/11); gs->SetMarkerStyle(fgMarker[is0]); gs->SetMarkerColor(fgColorS[is0]); gs->SetLineColor(fgColorS[is0]); - gs->SetNameTitle(Form("s_%d%02d%02d", ig, ic, is), fgParticle[is0]); + gs->SetLineStyle(il0);gs->SetLineWidth(2); + gs->SetName(Form("s_%d_%02d_%02d", ig, ic, is)); agM->AddAt(gm = new TGraphErrors(), is); gm->SetMarkerStyle(fgMarker[is0]); gm->SetMarkerColor(fgColorM[is0]); gm->SetLineColor(fgColorM[is0]); - gm->SetNameTitle(Form("m_%d%02d%02d", ig, ic, is), fgParticle[is0]); + gm->SetLineStyle(il0);gm->SetLineWidth(2); + gm->SetName(Form("m_%d_%02d_%02d", ig, ic, is)); + // this is important for labels in the legend + if(ic==0) { + gs->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel])); + gm->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel])); + } else if(ic==1) { + gs->SetTitle(Form("%s Ly[%d]", is%2 ?"z":"y", is/2)); + gm->SetTitle(Form("%s Ly[%d]", is%2?"z":"y", is/2)); + } else if(ic==2||ic==3) { + gs->SetTitle(Form("%s Ly[%d]", is%2 ?"RC":"no RC", is/2)); + gm->SetTitle(Form("%s Ly[%d]", is%2?"RC":"no RC", is/2)); + } else if(ic<=7) { + gs->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer)); + gm->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer)); + } else { + gs->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0)); + gm->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0)); + } } } else { aS->AddAt(gs = new TGraphErrors(), ic); gs->SetMarkerStyle(23); gs->SetMarkerColor(kRed); gs->SetLineColor(kRed); - gs->SetNameTitle(Form("s_%d%02d", ig, ic), ""); + gs->SetNameTitle(Form("s_%d_%02d", ig, ic), "sigma"); aM->AddAt(gm = ig ? (TGraph*)new TGraphErrors() : (TGraph*)new TGraphAsymmErrors(), ic); gm->SetLineColor(kBlack); gm->SetMarkerStyle(7); gm->SetMarkerColor(kBlack); - gm->SetNameTitle(Form("m_%d%02d", ig, ic), ""); + gm->SetNameTitle(Form("m_%d_%02d", ig, ic), "mean"); } - nc+=fgNcomp[ic]; } } } @@ -1259,23 +2141,30 @@ Bool_t AliTRDresolution::PostProcess() // Charge resolution //Process3DL(kCharge, 0, &fl); // Clusters residuals - Process2D(kCluster, 0, &fg, 1.e4); - Process2D(kCluster, 1, &fg); - fNRefFigures = 1; + Process3D(kCluster, 0, &fg, 1.e4); + Process3Dlinked(kCluster, 1, &fg); + fNRefFigures = 3; // Tracklet residual/pulls - Process2D(kTrackTRD , 0, &fg, 1.e4); - Process2D(kTrackTRD , 1, &fg); - Process2D(kTrackTRD , 2, &fg, 1.e4); - Process2D(kTrackTRD , 3, &fg); - Process2D(kTrackTRD , 4, &fg, 1.e3); - fNRefFigures = 4; - // TPC track residual/pulls - Process2D(kTrackTPC, 0, &fg, 1.e4); - Process2D(kTrackTPC, 1, &fg); - Process2D(kTrackTPC, 2, &fg, 1.e4); - Process2D(kTrackTPC, 3, &fg); - Process2D(kTrackTPC, 4, &fg, 1.e3); + Process3D(kTrack , 0, &fg, 1.e4); + Process3Dlinked(kTrack , 1, &fg); + Process3D(kTrack , 2, &fg, 1.e4); + Process3D(kTrack , 3, &fg); + Process2D(kTrack , 4, &fg, 1.e3); fNRefFigures = 7; + // TRDin residual/pulls + Process3D(kTrackIn, 0, &fg, 1.e4); + Process3Dlinked(kTrackIn, 1, &fg); + Process3D(kTrackIn, 2, &fg, 1.e4); + Process3D(kTrackIn, 3, &fg); + Process2D(kTrackIn, 4, &fg, 1.e3); + fNRefFigures = 11; + // TRDout residual/pulls + Process3D(kTrackOut, 0, &fg, 1.e3); // scale to fit - see PlotTrackOut + Process3Dlinked(kTrackOut, 1, &fg); + Process3D(kTrackOut, 2, &fg, 1.e4); + Process3D(kTrackOut, 3, &fg); + Process2D(kTrackOut, 4, &fg, 1.e3); + fNRefFigures = 15; if(!HasMCdata()) return kTRUE; @@ -1283,66 +2172,59 @@ Bool_t AliTRDresolution::PostProcess() //PROCESS MC RESIDUAL DISTRIBUTIONS // CLUSTER Y RESOLUTION/PULLS - Process3Drange(kMCcluster, 0, 0, &fg, 1.e4, 1, 1); - Process3Drange(kMCcluster, 0, 1, &fg, 1.e4, 2, 2); - Process3Drange(kMCcluster, 0, 2, &fg, 1.e4, 3, 12); - Process2D(kMCcluster, 1, &fg, 1., 3); - fNRefFigures = 10; + Process3D(kMCcluster, 0, &fg, 1.e4); + Process3Dlinked(kMCcluster, 1, &fg, 1.); + fNRefFigures = 17; // TRACKLET RESOLUTION/PULLS - Process3Drange(kMCtracklet, 0, 0, &fg, 1.e4, 1, 1); // y - Process3Drange(kMCtracklet, 0, 1, &fg, 1.e4, 2, 2); // y - Process3Drange(kMCtracklet, 0, 2, &fg, 1.e4, 3, 12); // y - Process2D(kMCtracklet, 1, &fg, 1., 3); // y pulls - Process2D(kMCtracklet, 2, &fg, 1.e4, 4); // z - Process2D(kMCtracklet, 3, &fg, 1., 5); // z pulls - Process2D(kMCtracklet, 4, &fg, 1.e3, 6); // phi - fNRefFigures = 14; + Process3D(kMCtracklet, 0, &fg, 1.e4); // y + Process3Dlinked(kMCtracklet, 1, &fg, 1.); // y pulls + Process3D(kMCtracklet, 2, &fg, 1.e4); // z + Process3D(kMCtracklet, 3, &fg, 1.); // z pulls + Process2D(kMCtracklet, 4, &fg, 1.e3); // phi + fNRefFigures = 21; // TRACK RESOLUTION/PULLS - Process2D(kMCtrackTRD, 0, &fg, 1.e4); // y - Process2D(kMCtrackTRD, 1, &fg); // y PULL - Process2D(kMCtrackTRD, 2, &fg, 1.e4); // z - Process2D(kMCtrackTRD, 3, &fg); // z PULL - Process2D(kMCtrackTRD, 4, &fg, 1.e3); // phi - Process2D(kMCtrackTRD, 5, &fg); // snp PULL - Process2D(kMCtrackTRD, 6, &fg, 1.e3); // theta - Process2D(kMCtrackTRD, 7, &fg); // tgl PULL - Process4D(kMCtrackTRD, 8, &fg, 1.e2); // pt resolution - Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 4);// pt resolution e1- @ L0 - Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 6);// pt resolution e1+ @ L0 - Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 55+4);// pt resolution e1- @ L5 - Process4D(kMCtrackTRD, 8, &fpt, 1.e2, 55+6);// pt resolution e1+ @ L5 - Process4D(kMCtrackTRD, 9, &fg); // 1/pt pulls - Process4D(kMCtrackTRD, 10, &fg, 1.e2); // p resolution - Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 4);// p resolution e1- @ L0 - Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 6);// p resolution e1+ @ L0 - Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 55+4);// p resolution e1- @ L5 - Process4D(kMCtrackTRD, 10, &fpt, 1.e2, 55+6);// p resolution e1+ @ L5 - Process4D(kMCtrackTRD, 11, &fg, 1.e2); // pt resolution stand alone - Process4D(kMCtrackTRD, 12, &fg); // 1/pt pulls stand alone - Process4D(kMCtrackTRD, 13, &fg, 1.e2); // p resolution stand alone - fNRefFigures = 24; - - // TRACK TPC RESOLUTION/PULLS - Process2D(kMCtrackTPC, 0, &fg, 1.e4);// y resolution - Process2D(kMCtrackTPC, 1, &fg); // y pulls - Process2D(kMCtrackTPC, 2, &fg, 1.e4);// z resolution - Process2D(kMCtrackTPC, 3, &fg); // z pulls - Process2D(kMCtrackTPC, 4, &fg, 1.e3);// phi resolution - Process2D(kMCtrackTPC, 5, &fg); // snp pulls - Process2D(kMCtrackTPC, 6, &fg, 1.e3);// theta resolution - Process2D(kMCtrackTPC, 7, &fg); // tgl pulls - Process3D(kMCtrackTPC, 8, &fg, 1.e2);// pt resolution - Process3D(kMCtrackTPC, 9, &fg); // 1/pt pulls - Process3D(kMCtrackTPC, 10, &fg, 1.e2);// p resolution - Process3D(kMCtrackTPC, 11, &fg); // p pulls - fNRefFigures = 30; - - // TRACK HMPID RESOLUTION/PULLS - Process2D(kMCtrackTOF, 0, &fg, 1.e4); // z towards TOF - Process2D(kMCtrackTOF, 1, &fg); // z towards TOF - fNRefFigures = 31; + Process3Darray(kMCtrack, 0, &fg, 1.e4); // y + Process3DlinkedArray(kMCtrack, 1, &fg); // y PULL + Process3Darray(kMCtrack, 2, &fg, 1.e4); // z + Process3Darray(kMCtrack, 3, &fg); // z PULL + Process2Darray(kMCtrack, 4, &fg, 1.e3); // phi + Process2Darray(kMCtrack, 5, &fg); // snp PULL + Process2Darray(kMCtrack, 6, &fg, 1.e3); // theta + Process2Darray(kMCtrack, 7, &fg); // tgl PULL + Process3Darray(kMCtrack, 8, &fg, 1.e2); // pt resolution + Process3Darray(kMCtrack, 9, &fg); // 1/pt pulls + Process3Darray(kMCtrack, 10, &fg, 1.e2); // p resolution + fNRefFigures+=16; + + // TRACK TRDin RESOLUTION/PULLS + Process3D(kMCtrackIn, 0, &fg, 1.e4);// y resolution + Process3Dlinked(kMCtrackIn, 1, &fg); // y pulls + Process3D(kMCtrackIn, 2, &fg, 1.e4);// z resolution + Process3D(kMCtrackIn, 3, &fg); // z pulls + Process2D(kMCtrackIn, 4, &fg, 1.e3);// phi resolution + Process2D(kMCtrackIn, 5, &fg); // snp pulls + Process2D(kMCtrackIn, 6, &fg, 1.e3);// theta resolution + Process2D(kMCtrackIn, 7, &fg); // tgl pulls + Process3D(kMCtrackIn, 8, &fg, 1.e2);// pt resolution + Process3D(kMCtrackIn, 9, &fg); // 1/pt pulls + Process3D(kMCtrackIn, 10, &fg, 1.e2);// p resolution + fNRefFigures+=8; + + // TRACK TRDout RESOLUTION/PULLS + Process3D(kMCtrackOut, 0, &fg, 1.e4);// y resolution + Process3Dlinked(kMCtrackOut, 1, &fg); // y pulls + Process3D(kMCtrackOut, 2, &fg, 1.e4);// z resolution + Process3D(kMCtrackOut, 3, &fg); // z pulls + Process2D(kMCtrackOut, 4, &fg, 1.e3);// phi resolution + Process2D(kMCtrackOut, 5, &fg); // snp pulls + Process2D(kMCtrackOut, 6, &fg, 1.e3);// theta resolution + Process2D(kMCtrackOut, 7, &fg); // tgl pulls + Process3D(kMCtrackOut, 8, &fg, 1.e2);// pt resolution + Process3D(kMCtrackOut, 9, &fg); // 1/pt pulls + Process3D(kMCtrackOut, 10, &fg, 1.e2);// p resolution + fNRefFigures+=8; return kTRUE; } @@ -1386,47 +2268,67 @@ void AliTRDresolution::AdjustF1(TH1 *h, TF1 *f) } //________________________________________________________ -TObjArray* AliTRDresolution::BuildMonitorContainerTracklet(const char* name) +TObjArray* AliTRDresolution::BuildMonitorContainerCluster(const char* name, Bool_t expand) { -// Build performance histograms for AliExternalTrackParam.vs TRD tracklet +// Build performance histograms for AliTRDcluster.vs TRD track or MC // - y reziduals/pulls -// - z reziduals/pulls -// - phi reziduals - TObjArray *arr = new TObjArray(5); + + TObjArray *arr = new TObjArray(2); arr->SetName(name); arr->SetOwner(); TH1 *h(NULL); char hname[100], htitle[300]; // tracklet resolution/pull in y direction - sprintf(hname, "%s_Y", name); - sprintf(htitle, "Y res @ %s;tg(#phi);#Delta y [cm];entries", name); - if(!(h = (TH2I*)gROOT->FindObject(hname))){ - h = new TH2I(hname, htitle, 21, -.33, .33, 100, -.5, .5); + sprintf(hname, "%s_%s_Y", GetNameId(), name); + sprintf(htitle, "Y res for \"%s\" @ %s;tg(#phi);#Delta y [cm];%s", GetNameId(), name, fgkResYsegmName[fSegmentLevel]); + if(!(h = (TH3S*)gROOT->FindObject(hname))){ + Int_t nybins=fgkNresYsegm[fSegmentLevel]; + if(expand) nybins*=2; + h = new TH3S(hname, htitle, + 48, -.48, .48, // phi + 60, -fDyRange, fDyRange, // dy + nybins, -0.5, nybins-0.5);// segment } else h->Reset(); arr->AddAt(h, 0); - sprintf(hname, "%s_Ypull", name); - sprintf(htitle, "Y pull @ %s;tg(#phi);#Delta y / #sigma_{y};entries", name); - if(!(h = (TH2I*)gROOT->FindObject(hname))){ - h = new TH2I(hname, htitle, 21, -.33, .33, 100, -4.5, 4.5); + sprintf(hname, "%s_%s_YZpull", GetNameId(), name); + sprintf(htitle, "YZ pull for \"%s\" @ %s;%s;#Delta y / #sigma_{y};#Delta z / #sigma_{z}", GetNameId(), name, fgkResYsegmName[fSegmentLevel]); + if(!(h = (TH3S*)gROOT->FindObject(hname))){ + h = new TH3S(hname, htitle, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5, 100, -4.5, 4.5, 100, -4.5, 4.5); } else h->Reset(); arr->AddAt(h, 1); + return arr; +} + +//________________________________________________________ +TObjArray* AliTRDresolution::BuildMonitorContainerTracklet(const char* name, Bool_t expand) +{ +// Build performance histograms for AliExternalTrackParam.vs TRD tracklet +// - y reziduals/pulls +// - z reziduals/pulls +// - phi reziduals + TObjArray *arr = BuildMonitorContainerCluster(name, expand); + arr->Expand(5); + TH1 *h(NULL); char hname[100], htitle[300]; + // tracklet resolution/pull in z direction - sprintf(hname, "%s_Z", name); - sprintf(htitle, "Z res @ %s;tg(#theta);#Delta z [cm];entries", name); - if(!(h = (TH2I*)gROOT->FindObject(hname))){ - h = new TH2I(hname, htitle, 50, -1., 1., 100, -1.5, 1.5); + sprintf(hname, "%s_%s_Z", GetNameId(), name); + sprintf(htitle, "Z res for \"%s\" @ %s;tg(#theta);#Delta z [cm];row cross", GetNameId(), name); + if(!(h = (TH3S*)gROOT->FindObject(hname))){ + h = new TH3S(hname, htitle, 50, -1., 1., 100, -1.5, 1.5, 2, -0.5, 1.5); } else h->Reset(); arr->AddAt(h, 2); - sprintf(hname, "%s_Zpull", name); - sprintf(htitle, "Z pull @ %s;tg(#theta);#Delta z / #sigma_{z};entries", name); - if(!(h = (TH2I*)gROOT->FindObject(hname))){ - h = new TH2I(hname, htitle, 50, -1., 1., 100, -5.5, 5.5); + sprintf(hname, "%s_%s_Zpull", GetNameId(), name); + sprintf(htitle, "Z pull for \"%s\" @ %s;tg(#theta);#Delta z / #sigma_{z};row cross", GetNameId(), name); + if(!(h = (TH3S*)gROOT->FindObject(hname))){ + h = new TH3S(hname, htitle, 50, -1., 1., 100, -5.5, 5.5, 2, -0.5, 1.5); + h->GetZaxis()->SetBinLabel(1, "no RC"); + h->GetZaxis()->SetBinLabel(2, "RC"); } else h->Reset(); arr->AddAt(h, 3); // tracklet to track phi resolution - sprintf(hname, "%s_PHI", name); - sprintf(htitle, "#Phi res @ %s;tg(#phi);#Delta #phi [rad];entries", name); + sprintf(hname, "%s_%s_PHI", GetNameId(), name); + sprintf(htitle, "#Phi res for \"%s\" @ %s;tg(#phi);#Delta #phi [rad];entries", GetNameId(), name); if(!(h = (TH2I*)gROOT->FindObject(hname))){ h = new TH2I(hname, htitle, 21, -.33, .33, 100, -.5, .5); } else h->Reset(); @@ -1445,69 +2347,29 @@ TObjArray* AliTRDresolution::BuildMonitorContainerTrack(const char* name) // - theta resolution, tgl pulls // - pt resolution, 1/pt pulls, p resolution - - TH1 *h(NULL); char hname[100], htitle[300]; -/* TObjArray *arr = new TObjArray(11); - arr->SetName(name); arr->SetOwner();*/ - -// // y resolution -// sprintf(hname, "%s_Y", name); -// sprintf(htitle, "Y res @ %s;tg(#phi);#Delta y [cm];entries", name); -// if(!(h = (TH2I*)gROOT->FindObject(hname))){ -// h = new TH2I(hname, htitle, 48, -.48, .48, 100, -.2, .2); -// } else h->Reset(); -// arr->AddAt(h, 0); -// // y pulls -// sprintf(hname, "%s_Ypull", name); -// sprintf(htitle, "Y pull @ %s;tg(#phi);#Delta y / #sigma_{y};entries", name); -// if(!(h = (TH2I*)gROOT->FindObject(hname))){ -// h = new TH2I(hname, htitle, 48, -.48, .48, 100, -4., 4.); -// } else h->Reset(); -// arr->AddAt(h, 1); -// -// // z resolution -// sprintf(hname, "%s_Z", name); -// sprintf(htitle, "Z res @ %s;tg(#theta);#Delta z [cm];entries", name); -// if(!(h = (TH2I*)gROOT->FindObject(hname))){ -// h = new TH2I(hname, htitle, 100, -1., 1., 100, -1., 1.); -// } else h->Reset(); -// arr->AddAt(h, 2); -// // z pulls -// sprintf(hname, "%s_Zpull", name); -// sprintf(htitle, "Z pull @ %s;tg(#theta);#Delta z / #sigma_{z};entries", name); -// if(!(h = (TH2I*)gROOT->FindObject(hname))){ -// h = new TH2I(hname, htitle, 100, -1., 1., 100, -4.5, 4.5); -// } else h->Reset(); -// arr->AddAt(h, 3); -// -// // phi resolution -// sprintf(hname, "%s_PHI", name); -// sprintf(htitle, "#Phi res @ %s;tg(#phi);#Delta #phi [rad];entries", name); -// if(!(h = (TH2I*)gROOT->FindObject(hname))){ -// h = new TH2I(hname, htitle, 60, -.3, .3, 100, -5e-3, 5e-3); -// } else h->Reset(); -// arr->AddAt(h, 4); - TObjArray *arr = BuildMonitorContainerTracklet(name); arr->Expand(11); + TH1 *h(NULL); char hname[100], htitle[300]; + TAxis *ax(NULL); + // snp pulls - sprintf(hname, "%s_SNPpull", name); - sprintf(htitle, "SNP pull @ %s;tg(#phi);#Delta snp / #sigma_{snp};entries", name); + sprintf(hname, "%s_%s_SNPpull", GetNameId(), name); + sprintf(htitle, "SNP pull for \"%s\" @ %s;tg(#phi);#Delta snp / #sigma_{snp};entries", GetNameId(), name); if(!(h = (TH2I*)gROOT->FindObject(hname))){ h = new TH2I(hname, htitle, 60, -.3, .3, 100, -4.5, 4.5); } else h->Reset(); arr->AddAt(h, 5); // theta resolution - sprintf(hname, "%s_THT", name); - sprintf(htitle, "#Theta res @ %s;tg(#theta);#Delta #theta [rad];entries", name); + sprintf(hname, "%s_%s_THT", GetNameId(), name); + sprintf(htitle, "#Theta res for \"%s\" @ %s;tg(#theta);#Delta #theta [rad];entries", GetNameId(), name); if(!(h = (TH2I*)gROOT->FindObject(hname))){ h = new TH2I(hname, htitle, 100, -1., 1., 100, -5e-3, 5e-3); } else h->Reset(); arr->AddAt(h, 6); // tgl pulls - sprintf(hname, "%s_TGLpull", name); - sprintf(htitle, "TGL pull @ %s;tg(#theta);#Delta tgl / #sigma_{tgl};entries", name); + sprintf(hname, "%s_%s_TGLpull", GetNameId(), name); + sprintf(htitle, "TGL pull for \"%s\" @ %s;tg(#theta);#Delta tgl / #sigma_{tgl};entries", GetNameId(), name); if(!(h = (TH2I*)gROOT->FindObject(hname))){ h = new TH2I(hname, htitle, 100, -1., 1., 100, -4.5, 4.5); } else h->Reset(); @@ -1516,34 +2378,40 @@ TObjArray* AliTRDresolution::BuildMonitorContainerTrack(const char* name) const Int_t kNpt(14); const Int_t kNdpt(150); const Int_t kNspc = 2*AliPID::kSPECIES+1; - Float_t Pt=0.1, DPt=-.1, Spc=-5.5; + Float_t lPt=0.1, lDPt=-.1, lSpc=-5.5; Float_t binsPt[kNpt+1], binsSpc[kNspc+1], binsDPt[kNdpt+1]; - for(Int_t i=0;iFindObject(hname))){ h = new TH3S(hname, htitle, kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc); + ax = h->GetZaxis(); + for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]); } else h->Reset(); arr->AddAt(h, 8); // 1/Pt pulls - sprintf(hname, "%s_1Pt", name); - sprintf(htitle, "1/P_{t} pull @ %s;1/p_{t}^{MC} [c/GeV];#Delta(1/p_{t})/#sigma(1/p_{t});SPECIES", name); + sprintf(hname, "%s_%s_1Pt", GetNameId(), name); + sprintf(htitle, "1/P_{t} pull for \"%s\" @ %s;1/p_{t}^{MC} [c/GeV];#Delta(1/p_{t})/#sigma(1/p_{t});SPECIES", GetNameId(), name); if(!(h = (TH3S*)gROOT->FindObject(hname))){ h = new TH3S(hname, htitle, kNpt, 0., 2., 100, -4., 4., kNspc, -5.5, 5.5); + ax = h->GetZaxis(); + for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]); } else h->Reset(); arr->AddAt(h, 9); // P resolution - sprintf(hname, "%s_P", name); - sprintf(htitle, "P res @ %s;p [GeV/c];#Delta p/p^{MC};SPECIES", name); + sprintf(hname, "%s_%s_P", GetNameId(), name); + sprintf(htitle, "P res for \"%s\" @ %s;p [GeV/c];#Delta p/p^{MC};SPECIES", GetNameId(), name); if(!(h = (TH3S*)gROOT->FindObject(hname))){ h = new TH3S(hname, htitle, kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc); + ax = h->GetZaxis(); + for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]); } else h->Reset(); arr->AddAt(h, 10); @@ -1567,14 +2435,14 @@ TObjArray* AliTRDresolution::Histos() // binnings for plots containing momentum or pt const Int_t kNpt(14), kNphi(48), kNdy(60); - Float_t Phi=-.48, Dy=-.3, Pt=0.1; + Float_t lPhi=-.48, lDy=-.3, lPt=0.1; Float_t binsPhi[kNphi+1], binsDy[kNdy+1], binsPt[kNpt+1]; - for(Int_t i=0; iAddAt(arr = new TObjArray(fgNhistos[kCharge]), kCharge); + fContainer->AddAt(arr = new TObjArray(2), kCharge); arr->SetName("Charge"); if(!(h = (TH3S*)gROOT->FindObject("hCharge"))){ h = new TH3S("hCharge", "Charge Resolution", 20, 1., 2., 24, 0., 3.6, 100, 0., 500.); @@ -1585,108 +2453,111 @@ TObjArray* AliTRDresolution::Histos() arr->AddAt(h, 0); // cluster to track residuals/pulls - fContainer->AddAt(arr = new TObjArray(fgNhistos[kCluster]), kCluster); - arr->SetName("Cl"); - if(!(h = (TH2I*)gROOT->FindObject("hCl"))){ - h = new TH2I("hCl", "Cluster Residuals", 21, -.33, .33, 100, -.5, .5); - h->GetXaxis()->SetTitle("tg(#phi)"); - h->GetYaxis()->SetTitle("#Delta y [cm]"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 0); - if(!(h = (TH2I*)gROOT->FindObject("hClpull"))){ - h = new TH2I("hClpull", "Cluster Pulls", 21, -.33, .33, 100, -4.5, 4.5); - h->GetXaxis()->SetTitle("tg(#phi)"); - h->GetYaxis()->SetTitle("#Delta y/#sigma_{y}"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 1); - + fContainer->AddAt(BuildMonitorContainerCluster("Cl"), kCluster); // tracklet to TRD track - fContainer->AddAt(BuildMonitorContainerTracklet("TrkTRD"), kTrackTRD); + fContainer->AddAt(BuildMonitorContainerTracklet("Trk", kTRUE), kTrack); // tracklet to TRDin - fContainer->AddAt(BuildMonitorContainerTracklet("TrkTRDin"), kTrackTPC); + fContainer->AddAt(BuildMonitorContainerTracklet("TrkIN", kTRUE), kTrackIn); + // tracklet to TRDout + fContainer->AddAt(BuildMonitorContainerTracklet("TrkOUT"), kTrackOut); // Resolution histos if(!HasMCdata()) return fContainer; - // cluster y resolution [0] - fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCcluster]), kMCcluster); - arr->SetName("McCl"); - if(!(h = (TH3S*)gROOT->FindObject("hMcCl"))){ - h = new TH3S("hMcCl", - "Cluster Resolution;tg(#phi);#Delta y [cm];p_{t} [GeV/c]", - kNphi, binsPhi, kNdy, binsDy, kNpt, binsPt); - } else h->Reset(); - arr->AddAt(h, 0); - if(!(h = (TH2I*)gROOT->FindObject("hMcClPull"))){ - h = new TH2I("hMcClPull", "Cluster Pulls", 48, -.48, .48, 100, -4.5, 4.5); - h->GetXaxis()->SetTitle("tg(#phi)"); - h->GetYaxis()->SetTitle("#Deltay/#sigma_{y}"); - h->GetZaxis()->SetTitle("p_{t} [GeV/c]"); - } else h->Reset(); - arr->AddAt(h, 1); - - - // TRACKLET RESOLUTION - fContainer->AddAt(arr = new TObjArray(fgNhistos[kMCtracklet]), kMCtracklet); - arr->SetName("McTrklt"); - // tracklet y resolution - if(!(h = (TH3S*)gROOT->FindObject("hMcTrkltY"))){ - h = new TH3S("hMcTrkltY", - "Tracklet Y Resolution;tg(#phi);#Delta y [cm];p_{t} [GeV/c]", - kNphi, binsPhi, kNdy, binsDy, kNpt, binsPt); - } else h->Reset(); - arr->AddAt(h, 0); - // tracklet y pulls - if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltYPull"))){ - h = new TH2I("hMcTrkltYPull", "Tracklet Pulls (Y)", 48, -.48, .48, 100, -4.5, 4.5); - h->GetXaxis()->SetTitle("tg(#phi)"); - h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 1); - // tracklet z resolution - if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltZ"))){ - h = new TH2I("hMcTrkltZ", "Tracklet Resolution (Z)", 100, -1., 1., 100, -1., 1.); - h->GetXaxis()->SetTitle("tg(#theta)"); - h->GetYaxis()->SetTitle("#Delta z [cm]"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 2); - // tracklet z pulls - if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltZPull"))){ - h = new TH2I("hMcTrkltZPull", "Tracklet Pulls (Z)", 100, -1., 1., 100, -3.5, 3.5); - h->GetXaxis()->SetTitle("tg(#theta)"); - h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 3); - // tracklet phi resolution - if(!(h = (TH2I*)gROOT->FindObject("hMcTrkltPhi"))){ - h = new TH2I("hMcTrkltPhi", "Tracklet Resolution (#Phi)", 48, -.48, .48, 100, -.15, .15); - h->GetXaxis()->SetTitle("tg(#phi)"); - h->GetYaxis()->SetTitle("#Delta #phi [rad]"); - h->GetZaxis()->SetTitle("entries"); - } else h->Reset(); - arr->AddAt(h, 4); + // cluster resolution + fContainer->AddAt(BuildMonitorContainerCluster("MCcl"), kMCcluster); + // tracklet resolution + fContainer->AddAt(BuildMonitorContainerTracklet("MCtracklet"), kMCtracklet); - // KALMAN TRACK RESOLUTION - fContainer->AddAt(arr = new TObjArray(6/*fgNhistos[kMCtrackTRD]*/), kMCtrackTRD); - arr->SetName("McTrkTRD"); - for(Int_t il(0); ilAddAt(BuildMonitorContainerTrack(Form("McTrkTRD_Ly%d", il)), il); + // track resolution + fContainer->AddAt(arr = new TObjArray(AliTRDgeometry::kNlayer), kMCtrack); + arr->SetName("MCtrk"); + for(Int_t il(0); ilAddAt(BuildMonitorContainerTrack(Form("MCtrk_Ly%d", il)), il); // TRDin TRACK RESOLUTION - fContainer->AddAt(BuildMonitorContainerTrack("McTrkTRDin"), kMCtrackTPC); + fContainer->AddAt(BuildMonitorContainerTrack("MCtrkIN"), kMCtrackIn); // TRDout TRACK RESOLUTION - fContainer->AddAt(BuildMonitorContainerTrack("McTrkTRDout"), kMCtrackTOF); + fContainer->AddAt(BuildMonitorContainerTrack("MCtrkOUT"), kMCtrackOut); return fContainer; } +//________________________________________________________ +Bool_t AliTRDresolution::Load(const Char_t *file, const Char_t *dir) +{ +// Custom load function. Used to guess the segmentation level of the data. + + if(!AliTRDrecoTask::Load(file, dir)) return kFALSE; + + // look for cluster residual plot - always available + TH3S* h3((TH3S*)((TObjArray*)fContainer->At(kClToTrk))->At(0)); + Int_t segmentation(h3->GetNbinsZ()/2); + if(segmentation==fgkNresYsegm[0]){ // default segmentation. Nothing to do + return kTRUE; + } else if(segmentation==fgkNresYsegm[1]){ // stack segmentation. + SetSegmentationLevel(1); + } else if(segmentation==fgkNresYsegm[2]){ // detector segmentation. + SetSegmentationLevel(2); + } else { + AliError(Form("Unknown segmentation [%d].", h3->GetNbinsZ())); + return kFALSE; + } + + AliDebug(2, Form("Segmentation set to level \"%s\"", fgkResYsegmName[fSegmentLevel])); + return kTRUE; +} + +//________________________________________________________ +Bool_t AliTRDresolution::Process(TH2* const h2, TGraphErrors **g, Int_t stat) +{ +// Generic function to process sigma/mean for 2D plot dy(x) + + if(!h2) { + if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : NULL pointer input container.\n"); + return kFALSE; + } + if(!Int_t(h2->GetEntries())){ + if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : Empty h[%s - %s].\n", h2->GetName(), h2->GetTitle()); + return kFALSE; + } + if(!g || !g[0]|| !g[1]) { + if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>0) printf("D-AliTRDresolution::Process() : NULL pointer output container.\n"); + return kFALSE; + } + // prepare + TAxis *ax(h2->GetXaxis()), *ay(h2->GetYaxis()); + TF1 f("f", "gaus", ay->GetXmin(), ay->GetXmax()); + Int_t n(0); + if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n); + if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n); + TH1D *h(NULL); + if((h=(TH1D*)gROOT->FindObject("py"))) delete h; + + // do actual loop + for(Int_t ix = 1, np=0; ix<=ax->GetNbins(); ix++){ + Double_t x = ax->GetBinCenter(ix); + Double_t ex= ax->GetBinWidth(ix)*0.288; // w/sqrt(12) + h = h2->ProjectionY("py", ix, ix); + if((n=(Int_t)h->GetEntries())1) printf("I-AliTRDresolution::Process() : Low statistics @ x[%f] stat[%d]=%d [%d].\n", x, ix, n, stat); + continue; + } + f.SetParameter(1, 0.); + f.SetParameter(2, 3.e-2); + h->Fit(&f, "QN"); + g[0]->SetPoint(np, x, f.GetParameter(1)); + g[0]->SetPointError(np, ex, f.GetParError(1)); + g[1]->SetPoint(np, x, f.GetParameter(2)); + g[1]->SetPointError(np, ex, f.GetParError(2)); + np++; + } + return kTRUE; +} + + //________________________________________________________ Bool_t AliTRDresolution::Process(TH2 * const h2, TF1 *f, Float_t k, TGraphErrors **g) { @@ -1698,21 +2569,30 @@ Bool_t AliTRDresolution::Process(TH2 * const h2, TF1 *f, Float_t k, TGraphErrors Int_t n = 0; if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n); if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n); + if(Int_t(h2->GetEntries())){ + AliDebug(4, Form("%s: g[%s %s]", pn, g[0]->GetName(), g[0]->GetTitle())); + } else { + AliDebug(2, Form("%s: g[%s %s]: Missing entries.", pn, g[0]->GetName(), g[0]->GetTitle())); + fIdxPlot++; + return kTRUE; + } - for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){ - Double_t x = h2->GetXaxis()->GetBinCenter(ibin); - TH1D *h = h2->ProjectionY(pn, ibin, ibin); - if(h->GetEntries()<100) continue; - //AdjustF1(h, f); - + const Int_t kINTEGRAL=1; + for(Int_t ibin = 0; ibin < Int_t(h2->GetNbinsX()/kINTEGRAL); ibin++){ + Int_t abin(ibin*kINTEGRAL+1),bbin(abin+kINTEGRAL-1),mbin(abin+Int_t(kINTEGRAL/2)); + Double_t x = h2->GetXaxis()->GetBinCenter(mbin); + TH1D *h = h2->ProjectionY(pn, abin, bbin); + if((n=(Int_t)h->GetEntries())<150){ + AliDebug(4, Form(" x[%f] range[%d %d] stat[%d] low statistics !", x, abin, bbin, n)); + continue; + } h->Fit(f, "QN"); - Int_t ip = g[0]->GetN(); + AliDebug(4, Form(" x_%d[%f] range[%d %d] stat[%d] M[%f] Sgm[%f]", ip, x, abin, bbin, n, f->GetParameter(1), f->GetParameter(2))); g[0]->SetPoint(ip, x, k*f->GetParameter(1)); g[0]->SetPointError(ip, 0., k*f->GetParError(1)); g[1]->SetPoint(ip, x, k*f->GetParameter(2)); g[1]->SetPointError(ip, 0., k*f->GetParError(2)); - /* g[0]->SetPoint(ip, x, k*h->GetMean()); g[0]->SetPointError(ip, 0., k*h->GetMeanError()); @@ -1733,8 +2613,20 @@ Bool_t AliTRDresolution::Process2D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F if(!fContainer || !fGraphS || !fGraphM) return kFALSE; // retrive containers - TH2I *h2 = idx<0 ? (TH2I*)(fContainer->At(plot)) : (TH2I*)((TObjArray*)(fContainer->At(plot)))->At(idx); - if(!h2) return kFALSE; + TH2I *h2(NULL); + if(idx<0){ + if(!(h2= (TH2I*)(fContainer->At(plot)))) return kFALSE; + } else{ + TObjArray *a0(NULL); + if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE; + if(!(h2=(TH2I*)a0->At(idx))) return kFALSE; + } + if(Int_t(h2->GetEntries())){ + AliDebug(2, Form("p[%d] idx[%d] : h[%s] %s", plot, idx, h2->GetName(), h2->GetTitle())); + } else { + AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx)); + return kFALSE; + } TGraphErrors *g[2]; if(gidx<0) gidx=idx; @@ -1755,8 +2647,20 @@ Bool_t AliTRDresolution::Process3D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F if(!fContainer || !fGraphS || !fGraphM) return kFALSE; // retrive containers - TH3S *h3 = idx<0 ? (TH3S*)(fContainer->At(plot)) : (TH3S*)((TObjArray*)(fContainer->At(plot)))->At(idx); - if(!h3) return kFALSE; + TH3S *h3(NULL); + if(idx<0){ + if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE; + } else{ + TObjArray *a0(NULL); + if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE; + if(!(h3=(TH3S*)a0->At(idx))) return kFALSE; + } + if(Int_t(h3->GetEntries())){ + AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle())); + } else { + AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx)); + return kFALSE; + } TObjArray *gm, *gs; if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE; @@ -1764,18 +2668,19 @@ Bool_t AliTRDresolution::Process3D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F TGraphErrors *g[2]; TAxis *az = h3->GetZaxis(); - for(Int_t iz=1; iz<=az->GetNbins(); iz++){ - if(!(g[0] = (TGraphErrors*)gm->At(iz-1))) return kFALSE; - if(!(g[1] = (TGraphErrors*)gs->At(iz-1))) return kFALSE; - az->SetRange(iz, iz); + for(Int_t iz(0); izGetEntriesFast(); iz++){ + if(!(g[0] = (TGraphErrors*)gm->At(iz))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(iz))) return kFALSE; + az->SetRange(iz+1, iz+1); if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE; } return kTRUE; } + //________________________________________________________ -Bool_t AliTRDresolution::Process3Drange(ETRDresolutionPlot plot, Int_t hidx, Int_t gidx, TF1 *f, Float_t k, Int_t zbin0, Int_t zbin1) +Bool_t AliTRDresolution::Process3Dlinked(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k) { // // Do the processing @@ -1784,19 +2689,38 @@ Bool_t AliTRDresolution::Process3Drange(ETRDresolutionPlot plot, Int_t hidx, Int if(!fContainer || !fGraphS || !fGraphM) return kFALSE; // retrive containers - TH3S *h3 = hidx<0 ? (TH3S*)(fContainer->At(plot)) : (TH3S*)((TObjArray*)(fContainer->At(plot)))->At(hidx); - if(!h3) return kFALSE; + TH3S *h3(NULL); + if(idx<0){ + if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE; + } else{ + TObjArray *a0(NULL); + if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE; + if(!(h3=(TH3S*)a0->At(idx))) return kFALSE; + } + if(Int_t(h3->GetEntries())){ + AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle())); + } else { + AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx)); + return kFALSE; + } + TObjArray *gm, *gs; + if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE; + if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE; TGraphErrors *g[2]; - if(!(g[0] = (TGraphErrors*)((TObjArray*)(fGraphM->At(plot)))->At(gidx))) return kFALSE; - if(!(g[1] = (TGraphErrors*)((TObjArray*)(fGraphS->At(plot)))->At(gidx))) return kFALSE; - TAxis *az = h3->GetZaxis(); - az->SetRange(zbin0, zbin1); + if(!(g[0] = (TGraphErrors*)gm->At(0))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(0))) return kFALSE; if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE; + + if(!(g[0] = (TGraphErrors*)gm->At(1))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(1))) return kFALSE; + if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE; + return kTRUE; } + //________________________________________________________ Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k) { @@ -1809,6 +2733,7 @@ Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f, // retrive containers TH3S *h3 = (TH3S*)((TObjArray*)fContainer->At(plot))->At(idx); if(!h3) return kFALSE; + AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle())); TGraphAsymmErrors *gm; TGraphErrors *gs; @@ -1841,7 +2766,45 @@ Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f, } //________________________________________________________ -Bool_t AliTRDresolution::Process4D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k, Int_t n) +Bool_t AliTRDresolution::Process2Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k) +{ + // + // Do the processing + // + + if(!fContainer || !fGraphS || !fGraphM) return kFALSE; + + // retrive containers + TObjArray *arr = (TObjArray*)(fContainer->At(plot)); + if(!arr) return kFALSE; + AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName())); + + TObjArray *gm, *gs; + if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE; + if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE; + + TGraphErrors *g[2]; TH2I *h2(NULL); TObjArray *a0(NULL); + for(Int_t ia(0); iaGetEntriesFast(); ia++){ + if(!(a0 = (TObjArray*)arr->At(ia))) continue; + + if(!(h2 = (TH2I*)a0->At(idx))) return kFALSE; + if(Int_t(h2->GetEntries())){ + AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h2->GetName(), h2->GetTitle())); + } else { + AliDebug(2, Form(" idx[%d] : Missing entries.", ia)); + continue; + } + + if(!(g[0] = (TGraphErrors*)gm->At(ia))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(ia))) return kFALSE; + if(!Process(h2, f, k, g)) return kFALSE; + } + + return kTRUE; +} + +//________________________________________________________ +Bool_t AliTRDresolution::Process3Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k) { // // Do the processing @@ -1851,100 +2814,153 @@ Bool_t AliTRDresolution::Process4D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, F //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx); // retrive containers - TObjArray *arr = (TObjArray*)((TObjArray*)(fContainer->At(plot)))->At(idx); + TObjArray *arr = (TObjArray*)(fContainer->At(plot)); if(!arr) return kFALSE; + AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName())); TObjArray *gm, *gs; if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE; if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE; - TGraphErrors *g[2]; + TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL); + Int_t in(0); + for(Int_t ia(0); iaGetEntriesFast(); ia++){ + if(!(a0 = (TObjArray*)arr->At(ia))) continue; - TH3S *h3(NULL); - for(Int_t ix=0, in=0; ixGetEntriesFast(); ix++){ - if(!(h3 = (TH3S*)arr->At(ix))) return kFALSE; + if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE; + if(Int_t(h3->GetEntries())){ + AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle())); + } else { + AliDebug(2, Form(" idx[%d] : Missing entries.", ia)); + continue; + } TAxis *az = h3->GetZaxis(); - //printf(" process ix[%d] bins[%d] in[%d]\n", ix, az->GetNbins(), in); for(Int_t iz=1; iz<=az->GetNbins(); iz++, in++){ - if(n>=0 && n!=in) continue; + if(in >= gm->GetEntriesFast()) break; if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE; if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE; - //printf(" g0[%s] g1[%s]\n", g[0]->GetName(), g[1]->GetName()); az->SetRange(iz, iz); if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE; } } + AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast())); + + return kTRUE; +} + +//________________________________________________________ +Bool_t AliTRDresolution::Process3DlinkedArray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k) +{ + // + // Do the processing + // + + if(!fContainer || !fGraphS || !fGraphM) return kFALSE; + //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx); + + // retrive containers + TObjArray *arr = (TObjArray*)(fContainer->At(plot)); + if(!arr) return kFALSE; + AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName())); + + TObjArray *gm, *gs; + if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE; + if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE; + + TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL); + Int_t in(0); + for(Int_t ia(0); iaGetEntriesFast(); ia++){ + if(!(a0 = (TObjArray*)arr->At(ia))) continue; + if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE; + if(Int_t(h3->GetEntries())){ + AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle())); + } else { + AliDebug(2, Form(" idx[%d] : Missing entries.", ia)); + continue; + } + if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE; + if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE; + in++; + + if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE; + if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE; + if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE; + in++; + } + AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast())); return kTRUE; } //________________________________________________________ -Bool_t AliTRDresolution::GetGraphPlot(Float_t *bb, ETRDresolutionPlot ip, Int_t idx) +Bool_t AliTRDresolution::GetGraph(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, const Char_t *explain) { // // Get the graphs // if(!fGraphS || !fGraphM) return kFALSE; + // axis titles look up + Int_t nref = 0; + for(Int_t jp=0; jp<(Int_t)ip; jp++) nref+=fgNproj[jp]; + UChar_t jdx = idx<0?0:idx; + for(Int_t jc=0; jcSetBorderSize(0); + leg->SetFillStyle(0); + } + // build frame + TH1S *h1(NULL); + h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]); + h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]); + h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2); + // axis range + TAxis *ax = h1->GetXaxis(); + ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2); + ax = h1->GetYaxis(); + ax->SetRangeUser(bb[1], bb[3]); + ax->CenterTitle(); ax->SetTitleOffset(1.4); + h1->Draw(); + // bounding box + TBox *b = new TBox(-.15, bb[1], .15, bb[3]); + b->SetFillStyle(3002);b->SetLineColor(0); + b->SetFillColor(ip<=Int_t(kMCcluster)?kGreen:kBlue); + b->Draw(); - //printf("plotting task[%d] gidx[%d]\n", ip, idx); TGraphErrors *gm = idx<0 ? (TGraphErrors*)fGraphM->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphM->At(ip)))->At(idx); if(!gm) return kFALSE; TGraphErrors *gs = idx<0 ? (TGraphErrors*)fGraphS->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphS->At(ip)))->At(idx); if(!gs) return kFALSE; - //printf("gs[%s] gm[%s]\n", gs->GetName(), gm->GetName()); - gs->Draw("apl"); gm->Draw("pl"); - //return kTRUE; - // titles look up - Int_t nref = 0; - for(Int_t jp=0; jp<(Int_t)ip; jp++) nref+=fgNproj[jp]; - UChar_t jdx = idx<0?0:idx; - for(Int_t jc=0; jcGetN())) { + nPlots++; + gm->Draw("pl"); if(leg) leg->AddEntry(gm, gm->GetTitle(), "pl"); PutTrendValue(Form("%s_%s", fgPerformanceName[ip], at[0]), gm->GetMean(2)); PutTrendValue(Form("%s_%sRMS", fgPerformanceName[ip], at[0]), gm->GetRMS(2)); } if((n=gs->GetN())){ + nPlots++; + gs->Draw("pl"); if(leg) leg->AddEntry(gs, gs->GetTitle(), "pl"); gs->Sort(&TGraph::CompareY); PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[ip], at[0]), gs->GetY()[0]); PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[ip], at[0]), gs->GetY()[n-1]); gs->Sort(&TGraph::CompareX); } - - // axis range - TAxis *ax(NULL); TH1 *hf(NULL); - hf = gs->GetHistogram(); - hf->SetTitle(at[0]); - ax = hf->GetXaxis(); - ax->SetRangeUser(bb[0], bb[2]); - ax->SetTitle(at[1]);ax->CenterTitle(); - - ax = hf->GetYaxis(); - ax->SetRangeUser(bb[1], bb[3]); - ax->SetTitleOffset(1.1); - ax->SetTitle(at[2]);ax->CenterTitle(); - - TGaxis *gax = NULL; - gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U"); - gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed); - //gax->SetVertical(); - gax->CenterTitle(); gax->SetTitleOffset(.7); - gax->SetTitle(at[3]); gax->Draw(); - - // bounding box - TBox *b = new TBox(-.15, bb[1], .15, bb[3]); - b->SetFillStyle(3002);b->SetLineColor(0); - b->SetFillColor(ip<=Int_t(kMCcluster)?kGreen:kBlue); - b->Draw(); + if(!nPlots) return kFALSE; + if(leg) leg->Draw(); return kTRUE; } //________________________________________________________ -Bool_t AliTRDresolution::GetGraphTrack(Float_t *bb, Int_t idx, Int_t first, Int_t n, Bool_t kLEG) +Bool_t AliTRDresolution::GetGraphArray(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, Int_t n, Int_t *sel, const Char_t *explain) { // // Get the graphs @@ -1954,25 +2970,25 @@ Bool_t AliTRDresolution::GetGraphTrack(Float_t *bb, Int_t idx, Int_t first, Int_ // axis titles look up Int_t nref(0); - for(Int_t jp=0; jpSetHeader("Mean"); - legM->SetBorderSize(1); - legM->SetFillColor(0); + legM->SetBorderSize(0); + legM->SetFillStyle(0); legS=new TLegend(.65, .6, .95, .9); legS->SetHeader("Sigma"); - legS->SetBorderSize(1); - legS->SetFillColor(0); + legS->SetBorderSize(0); + legS->SetFillStyle(0); } - Int_t layer(first/11); + // build frame TH1S *h1(NULL); - h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %d;%s;%s", at[0], layer, at[1], at[2]), 2, bb[0], bb[2]); + h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]); h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]); h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2); // axis range @@ -1982,111 +2998,192 @@ Bool_t AliTRDresolution::GetGraphTrack(Float_t *bb, Int_t idx, Int_t first, Int_ ax->SetRangeUser(bb[1], bb[3]); ax->CenterTitle(); ax->SetTitleOffset(1.4); h1->Draw(); -// TGaxis *gax = NULL; -// gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U"); -// gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed); -// //gax->SetVertical(); -// gax->CenterTitle(); //gax->SetTitleOffset(.5);gax->SetTitleSize(.06); -// gax->SetTitle(at[3]); gax->Draw(); - - TGraphErrors *gm = NULL, *gs = NULL; - TObjArray *a0 = NULL, *a1 = NULL; - a0 = (TObjArray*)((TObjArray*)fGraphM->At(kMCtrackTRD))->At(idx); - a1 = (TObjArray*)((TObjArray*)fGraphS->At(kMCtrackTRD))->At(idx); - Int_t nn(0), m(n/2); - for(Int_t is=first, is0=0; isAt(is))) return kFALSE; - if(!(gm = (TGraphErrors*)a0->At(is))) return kFALSE; + + TGraphErrors *gm(NULL), *gs(NULL); + TObjArray *a0(NULL), *a1(NULL); + a0 = (TObjArray*)((TObjArray*)fGraphM->At(ip))->At(idx); + a1 = (TObjArray*)((TObjArray*)fGraphS->At(ip))->At(idx); + if(!n) n=a0->GetEntriesFast(); + AliDebug(4, Form("Graph : Ref[%d] Title[%s] Limits{x[%f %f] y[%f %f]} Comp[%d] Selection[%c]", nref, at[0], bb[0], bb[2], bb[1], bb[3], n, sel ? 'y' : 'n')); + Int_t nn(0), nPlots(0); + for(Int_t is(0), is0(0); isAt(is0))) return kFALSE; + if(!(gm = (TGraphErrors*)a0->At(is0))) return kFALSE; if((nn=gs->GetN())){ - gs->Draw("pc"); if(legS) legS->AddEntry(gs, gs->GetTitle(), "pl"); + nPlots++; + gs->Draw("pc"); + if(legS){ + //printf("LegEntry %s [%s]%s\n", at[0], gs->GetName(), gs->GetTitle()); + legS->AddEntry(gs, gs->GetTitle(), "pl"); + } gs->Sort(&TGraph::CompareY); - PutTrendValue(Form("%s_%sSigMin%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[0]); - PutTrendValue(Form("%s_%sSigMax%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[nn-1]); + PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[0]); + PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[nn-1]); gs->Sort(&TGraph::CompareX); } if(gm->GetN()){ - gm->Draw("pc");if(legM) legM->AddEntry(gm, gm->GetTitle(), "pl"); - PutTrendValue(Form("%s_%s_%s", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gm->GetMean(2)); - PutTrendValue(Form("%s_%s_%sRMS", fgPerformanceName[kMCtrackTRD], at[0], AliPID::ParticleShortName(is0)), gm->GetRMS(2)); + nPlots++; + gm->Draw("pc"); + if(legM){ + //printf("LegEntry %s [%s]%s\n", at[0], gm->GetName(), gm->GetTitle()); + legM->AddEntry(gm, gm->GetTitle(), "pl"); + } + PutTrendValue(Form("%s_%s", fgPerformanceName[kMCtrack], at[0]), gm->GetMean(2)); + PutTrendValue(Form("%s_%sRMS", fgPerformanceName[kMCtrack], at[0]), gm->GetRMS(2)); } } - if(kLEG){legM->Draw();legS->Draw();} + if(!nPlots) return kFALSE; + if(kLEG){ + legM->Draw(); + legS->Draw(); + } return kTRUE; } +//____________________________________________________________________ +Bool_t AliTRDresolution::FitTrack(const Int_t np, AliTrackPoint *points, Float_t param[10]) +{ +// +// Fit track with a staight line using the "np" clusters stored in the array "points". +// The following particularities are stored in the clusters from points: +// 1. pad tilt as cluster charge +// 2. pad row cross or vertex constrain as fake cluster with cluster type 1 +// The parameters of the straight line fit are stored in the array "param" in the following order : +// param[0] - x0 reference radial position +// param[1] - y0 reference r-phi position @ x0 +// param[2] - z0 reference z position @ x0 +// param[3] - slope dy/dx +// param[4] - slope dz/dx +// +// Attention : +// Function should be used to refit tracks for B=0T +// -//________________________________________________________ -Bool_t AliTRDresolution::GetGraphTrackTPC(Float_t *bb, Int_t idx, Int_t ist, Int_t n, Bool_t kLEG) + if(np<40){ + if(AliLog::GetDebugLevel("PWG1", "AliTRDresolution")>1) printf("D-AliTRDresolution::FitTrack: Not enough clusters to fit a track [%d].", np); + return kFALSE; + } + TLinearFitter yfitter(2, "pol1"), zfitter(2, "pol1"); + + Double_t x0(0.); + for(Int_t ip(0); ip3) printf("D-AliTRDresolution::FitTrack: x0[%f] y0[%f] z0[%f] dydx[%f] dzdx[%f]\n", x0, y0, z0, dydx, dzdx); + return kTRUE; +} + +//____________________________________________________________________ +Bool_t AliTRDresolution::FitTracklet(const Int_t ly, const Int_t np, AliTrackPoint *points, const Float_t param[10], Float_t par[3]) { - // - // Get the graphs - // +// +// Fit tracklet with a staight line using the coresponding subset of clusters out of the total "np" clusters stored in the array "points". +// See function FitTrack for the data stored in the "clusters" array - if(!fGraphS || !fGraphM) return kFALSE; +// The parameters of the straight line fit are stored in the array "param" in the following order : +// par[0] - x0 reference radial position +// par[1] - y0 reference r-phi position @ x0 +// par[2] - slope dy/dx +// +// Attention : +// Function should be used to refit tracks for B=0T +// - // axis titles look up - Int_t nref = 0; - for(Int_t jp=0; jp1) printf("D-AliTRDresolution::FitTracklet: Not enough clusters to fit a tracklet [%d].", nly); + return kFALSE; + } + // set radial reference for fit + x0 /= Float_t(nly); - TLegend *legM(NULL), *legS(NULL); - if(kLEG){ - legM=new TLegend(.35, .6, .65, .9); - legM->SetHeader("Mean"); - legM->SetBorderSize(1); - legM->SetFillColor(0); - legS=new TLegend(.65, .6, .95, .9); - legS->SetHeader("Sigma"); - legS->SetBorderSize(1); - legS->SetFillColor(0); + // find tracklet core + Double_t mean(0.), sig(1.e3); + AliMathBase::EvaluateUni(nly, dy, mean, sig, 0); + + // simple cluster error parameterization + Float_t kSigCut = TMath::Sqrt(5.e-4 + param[3]*param[3]*0.018); + + // fit tracklet core + for(Int_t jly(0); jlykSigCut) continue; + Double_t dx(x[jly]-x0); + yfitter.AddPoint(&dx, y[jly], 1.); } - TH1S *h1(NULL); - h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), at[0], 2, bb[0], bb[2]); - h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]); - h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2); - // axis range - TAxis *ax = h1->GetXaxis(); - ax->SetTitle(at[1]);ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2); - ax = h1->GetYaxis(); - ax->SetRangeUser(bb[1], bb[3]); - ax->SetTitleOffset(1.4);//ax->SetTitleSize(.06); - ax->SetTitle(at[2]);ax->CenterTitle(); - h1->Draw(); -// TGaxis *gax = NULL; -// gax = new TGaxis(bb[2], bb[1], bb[2], bb[3], bb[1], bb[3], 510, "+U"); -// gax->SetLineColor(kRed);gax->SetLineWidth(2);gax->SetTextColor(kRed); -// //gax->SetVertical(); -// gax->CenterTitle(); //gax->SetTitleOffset(.5);gax->SetTitleSize(.06); -// gax->SetTitle(at[3]); gax->Draw(); - - Int_t nn(0), m(n/2); - TGraphErrors *gm = NULL, *gs = NULL; - TObjArray *a0 = NULL, *a1 = NULL; - a0 = (TObjArray*)((TObjArray*)fGraphM->At(kMCtrackTPC))->At(idx); - a1 = (TObjArray*)((TObjArray*)fGraphS->At(kMCtrackTPC))->At(idx); - for(Int_t is=ist, is0=0; isAt(is))) return kFALSE; - if(!(gm = (TGraphErrors*)a0->At(is))) return kFALSE; - if((nn=gs->GetN())){ - gs->Draw("pl");if(legS) legS->AddEntry(gs, gs->GetTitle(), "pl"); - gs->Sort(&TGraph::CompareY); - PutTrendValue(Form("%s_%sSigMin%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[0]); - PutTrendValue(Form("%s_%sSigMax%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gs->GetY()[nn-1]); - gs->Sort(&TGraph::CompareX); - } - if(gm->GetN()){ - gm->Draw("pl"); if(legM) legM->AddEntry(gm, gm->GetTitle(), "pl"); - PutTrendValue(Form("%s_%s_%s", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gm->GetMean(2)); - PutTrendValue(Form("%s_%s_%sRMS", fgPerformanceName[kMCtrackTPC], at[0], AliPID::ParticleShortName(is0)), gm->GetRMS(2)); - } + if(yfitter.Eval() != 0) return kFALSE; + par[0] = x0; + par[1] = yfitter.GetParameter(0); + par[2] = yfitter.GetParameter(1); + return kTRUE; +} + +//____________________________________________________________________ +Bool_t AliTRDresolution::UseTrack(const Int_t np, AliTrackPoint *points, Float_t param[10]) +{ +// +// Global selection mechanism of tracksbased on cluster to fit residuals +// The parameters are the same as used ni function FitTrack(). + + const Float_t kS(0.6), kM(0.2); + TH1S h("h1", "", 100, -5.*kS, 5.*kS); + Float_t dy, dz, s, m; + for(Int_t ip(0); ipDraw(); legS->Draw();} + TF1 fg("fg", "gaus", -5.*kS, 5.*kS); + fg.SetParameter(1, 0.); + fg.SetParameter(2, 2.e-2); + h.Fit(&fg, "QN"); + m=fg.GetParameter(1); s=fg.GetParameter(2); + if(s>kS || TMath::Abs(m)>kM) return kFALSE; return kTRUE; } @@ -2117,8 +3214,96 @@ void AliTRDresolution::GetLandauMpvFwhm(TF1 * const f, Float_t &mpv, Float_t &xm //________________________________________________________ -void AliTRDresolution::SetRecoParam(AliTRDrecoParam *r) +void AliTRDresolution::SetSegmentationLevel(Int_t l) { - - fReconstructor->SetRecoParam(r); +// Setting the segmentation level to "l" + fSegmentLevel = l; + + UShort_t const lNcomp[kNprojs] = { + 1, 1, //2, + fgkNresYsegm[fSegmentLevel], 2, //2, + 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, + 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, + 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, + // MC + fgkNresYsegm[fSegmentLevel], 2, //2, + fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5, + fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11 + fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11 + 6*fgkNresYsegm[fSegmentLevel], 6*2, 6*2, 6*2, 6, 6, 6, 6, 6*11, 6*11, 6*11 //11 + }; + memcpy(fNcomp, lNcomp, kNprojs*sizeof(UShort_t)); + + Char_t const *lAxTitle[kNprojs][4] = { + // Charge + {"Impv", "x [cm]", "I_{mpv}", "x/x_{0}"} + ,{"dI/Impv", "x/x_{0}", "#delta I/I_{mpv}", "x[cm]"} + // Clusters to Kalman + ,{"Cluster2Track residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"Cluster2Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + // TRD tracklet to Kalman fit + ,{"Tracklet2Track Y residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"Tracklet2Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"Tracklet2Track Z residuals", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"Tracklet2Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"} + ,{"Tracklet2Track Phi residuals", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + // TRDin 2 first TRD tracklet + ,{"Tracklet2Track Y residuals @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"Tracklet2Track YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"Tracklet2Track Z residuals @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"Tracklet2Track Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"} + ,{"Tracklet2Track Phi residuals @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + // TRDout 2 first TRD tracklet + ,{"Tracklet2Track Y residuals @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"Tracklet2Track YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"Tracklet2Track Z residuals @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"Tracklet2Track Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"} + ,{"Tracklet2Track Phi residuals @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + // MC cluster + ,{"MC Cluster Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"MC Cluster YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + // MC tracklet + ,{"MC Tracklet Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"} + ,{"MC Tracklet YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"MC Tracklet Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"MC Tracklet Z pulls", "tg(#theta)", "z", "#sigma_{z}"} + ,{"MC Tracklet Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + // MC track TRDin + ,{"MC Y resolution @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"} + ,{"MC YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"MC Z resolution @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"MC Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"} + ,{"MC #Phi resolution @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + ,{"MC SNP pulls @ TRDin", "tg(#phi)", "SNP", "#sigma_{snp}"} + ,{"MC #Theta resolution @ TRDin", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"} + ,{"MC TGL pulls @ TRDin", "tg(#theta)", "TGL", "#sigma_{tgl}"} + ,{"MC P_{t} resolution @ TRDin", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"} + ,{"MC 1/P_{t} pulls @ TRDin", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"} + ,{"MC P resolution @ TRDin", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"} + // MC track TRDout + ,{"MC Y resolution @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"} + ,{"MC YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"MC Z resolution @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"MC Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"} + ,{"MC #Phi resolution @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + ,{"MC SNP pulls @ TRDout", "tg(#phi)", "SNP", "#sigma_{snp}"} + ,{"MC #Theta resolution @ TRDout", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"} + ,{"MC TGL pulls @ TRDout", "tg(#theta)", "TGL", "#sigma_{tgl}"} + ,{"MC P_{t} resolution @ TRDout", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"} + ,{"MC 1/P_{t} pulls @ TRDout", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"} + ,{"MC P resolution @ TRDout", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"} + // MC track in TRD + ,{"MC Track Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"} + ,{"MC Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"} + ,{"MC Track Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"} + ,{"MC Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"} + ,{"MC Track #Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"} + ,{"MC Track SNP pulls", "tg(#phi)", "SNP", "#sigma_{snp}"} + ,{"MC Track #Theta resolution", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"} + ,{"MC Track TGL pulls", "tg(#theta)", "TGL", "#sigma_{tgl}"} + ,{"MC P_{t} resolution", "p_{t} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "#sigma(#Deltap_{t}/p_{t}^{MC}) [%]"} + ,{"MC 1/P_{t} pulls", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC} - 1/p_{t}^{MC}", "#sigma_{1/p_{t}}"} + ,{"MC P resolution", "p [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "#sigma(#Deltap/p^{MC}) [%]"} + }; + memcpy(fAxTitle, lAxTitle, 4*kNprojs*sizeof(Char_t*)); }