fix for bug https://savannah.cern.ch/bugs/?98881

Decouple the definition of the number of dEdx slices from their physics
meaning. In the TRD tracklet the definition represents allocated memory.
The memory is filled in AliTRDtrackV1::UpdateESDtrack according to
physics definitions from AliTRDPIDResponse and conventions.

The fix has to be suplemented by changes in AliESDtrack.cxx and
AliAODPid.h as requested in https://savannah.cern.ch/bugs/?98881#comment5

diff --git a/TRD/AliTRDseedV1.cxx b/TRD/AliTRDseedV1.cxx
index 3208df689749eeb4e53153e8e8ab62586d566833..e06f1c2cf610d0cbb733cda16df8ecc4d2b92b1a 100644 (file)
--- a/TRD/AliTRDseedV1.cxx
+++ b/TRD/AliTRDseedV1.cxx
@@ -102,7 +102,7 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
   fZref[0] = 0.; fZref[1] = 0.; 
   fYfit[0] = 0.; fYfit[1] = 0.; 
   fZfit[0] = 0.; fZfit[1] = 0.; 
-  memset(fdEdx, 0, kNslices*sizeof(Float_t)); 
+  memset(fdEdx, 0, kNdEdxSlices*sizeof(Float_t));
   for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec]  = -1.;
   fLabels[0]=-1; fLabels[1]=-1; // most freq MC labels
   fLabels[2]=0;  // number of different labels for tracklet
@@ -224,7 +224,7 @@ void AliTRDseedV1::Copy(TObject &ref) const
   target.fZref[0] = fZref[0]; target.fZref[1] = fZref[1]; 
   target.fYfit[0] = fYfit[0]; target.fYfit[1] = fYfit[1]; 
   target.fZfit[0] = fZfit[0]; target.fZfit[1] = fZfit[1]; 
-  memcpy(target.fdEdx, fdEdx, kNslices*sizeof(Float_t)); 
+  memcpy(target.fdEdx, fdEdx, kNdEdxSlices*sizeof(Float_t));
   memcpy(target.fProb, fProb, AliPID::kSPECIES*sizeof(Float_t)); 
   memcpy(target.fLabels, fLabels, 3*sizeof(Int_t)); 
   memcpy(target.fRefCov, fRefCov, 7*sizeof(Double_t)); 
@@ -308,7 +308,7 @@ void AliTRDseedV1::Reset(Option_t *opt)
   fZref[0] = 0.; fZref[1] = 0.; 
   fYfit[0] = 0.; fYfit[1] = 0.; 
   fZfit[0] = 0.; fZfit[1] = 0.; 
-  memset(fdEdx, 0, kNslices*sizeof(Float_t)); 
+  memset(fdEdx, 0, kNdEdxSlices*sizeof(Float_t));
   for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) fProb[ispec]  = -1.;
   fLabels[0]=-1; fLabels[1]=-1; // most freq MC labels
   fLabels[2]=0;  // number of different labels for tracklet
@@ -417,7 +417,7 @@ void AliTRDseedV1::CookdEdx(Int_t nslices)
 // 3. cluster size
 //
 
-  memset(fdEdx, 0, kNslices*sizeof(Float_t));
+  memset(fdEdx, 0, kNdEdxSlices*sizeof(Float_t));
   const Double_t kDriftLength = (.5 * AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
 
   AliTRDcluster *c(NULL);
@@ -679,11 +679,11 @@ Float_t AliTRDseedV1::GetMomentum(Float_t *err) const
 //
 
   Double_t p = fPt*TMath::Sqrt(1.+fZref[1]*fZref[1]);
-  Double_t p2 = p*p;
-  Double_t tgl2 = fZref[1]*fZref[1];
-  Double_t pt2 = fPt*fPt;
   if(err){
-    Double_t s2 = 
+    Double_t p2 = p*p;
+    Double_t tgl2 = fZref[1]*fZref[1];
+    Double_t pt2 = fPt*fPt;
+    Double_t s2 =
       p2*tgl2*pt2*pt2*fRefCov[4]
      -2.*p2*fZref[1]*fPt*pt2*fRefCov[5]
      +p2*pt2*fRefCov[6];
@@ -746,6 +746,8 @@ Bool_t AliTRDseedV1::CookPID()
 // - Neural Network [default] - option "nn"  
 // - 2D Likelihood - option "!nn"  
 
+  AliWarning(Form("Obsolete function. Use AliTRDPIDResponse::GetResponse() instead."));
+
   AliTRDcalibDB *calibration = AliTRDcalibDB::Instance();
   if (!calibration) {
     AliError("No access to calibration data");
@@ -1143,7 +1145,7 @@ Bool_t  AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *const chamber, Bool_
   memset(s2y, 0, kNrows*kNcls*sizeof(Double_t));
   memset(blst, 0, kNrows*kNcls*sizeof(Bool_t));   //this is 8 times faster to memset than "memset(clst, 0, kNrows*kNcls*sizeof(AliTRDcluster*))"
 
-  Double_t roady(0.), s2Mean(0.), sMean(0.); Int_t ns2Mean(0);
+  Double_t roady(0.), s2Mean(0.); Int_t ns2Mean(0);
 
   // Do cluster projection and pick up cluster candidates
   AliTRDcluster *c(NULL);
@@ -1216,7 +1218,7 @@ Bool_t  AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *const chamber, Bool_
     for(Int_t ir(kNrows);ir--;) clst[ir].Clear();
     return kFALSE;
   }
-  s2Mean /= ns2Mean; sMean = TMath::Sqrt(s2Mean);
+  s2Mean /= ns2Mean; //sMean = TMath::Sqrt(s2Mean);
   //Double_t sRef(TMath::Sqrt(s2Mean+s2yTrk)); // reference error parameterization
 
   // organize row candidates
@@ -1939,9 +1941,10 @@ Bool_t AliTRDseedV1::Fit(UChar_t opt)
       kDZDX=kFALSE;
     }
 
-    Float_t w = 1.;
-    if(c->GetNPads()>4) w = .5;
-    if(c->GetNPads()>5) w = .2;
+//   TODO use this information to adjust cluster error parameterization
+//     Float_t w = 1.;
+//     if(c->GetNPads()>4) w = .5;
+//     if(c->GetNPads()>5) w = .2;
 
     // cluster charge
     qc[n]   = TMath::Abs(c->GetQ());

diff --git a/TRD/AliTRDseedV1.h b/TRD/AliTRDseedV1.h
index 3cdeabc0763d60b036af214c689280df415dfca9..c48da16d046a86c536c11f23645261d93925ebb0 100644 (file)
--- a/TRD/AliTRDseedV1.h
+++ b/TRD/AliTRDseedV1.h
@@ -54,19 +54,18 @@ public:
    ,kMask      = 0x3f   // bit mask
    ,kNtb       = 31     // max clusters/pad row
    ,kNclusters = 2*kNtb // max number of clusters/tracklet
-   ,kNslices   = 10     // max dEdx slices
+   ,kNdEdxSlices= 8     // dEdx slices allocated in reconstruction
   };
 
   // bits from 0-13 are reserved by ROOT (see TObject.h)
   enum ETRDtrackletStatus {
     kOwner      = BIT(14) // owner of its clusters
    ,kRowCross   = BIT(15) // pad row cross tracklet
-   ,kPID        = BIT(16) // PID contributor
+   ,kChmbGood   = BIT(16) // status of the detector from calibration view point
    ,kCalib      = BIT(17) // calibrated tracklet
    ,kKink       = BIT(18) // kink prolongation tracklet
    ,kStandAlone = BIT(19) // tracklet build during stand alone track finding
    ,kPrimary    = BIT(20) // tracklet from a primary track candidate
-   ,kChmbGood   = BIT(21) // status of the detector from calibration view point
   };
 
   enum ETRDtrackletError { // up to 8 bits
@@ -101,7 +100,6 @@ public:
   Bool_t    IsOwner() const          { return TestBit(kOwner);}
   Bool_t    IsKink() const           { return TestBit(kKink);}
   Bool_t    IsPrimary() const        { return TestBit(kPrimary);}
-  Bool_t    HasPID() const           { return TestBit(kPID);}
   Bool_t    HasError(ETRDtrackletError err) const
                                      { return TESTBIT(fErrorMsg, err);}
   Bool_t    IsOK() const             { return GetN() > 4 && GetNUsed() < 4;}
@@ -189,7 +187,6 @@ public:
   void      SetLabels(Int_t *lbls)   { memcpy(fLabels, lbls, 3*sizeof(Int_t)); }
   void      SetKink(Bool_t k = kTRUE){ SetBit(kKink, k);}
   void      SetPrimary(Bool_t k = kTRUE){ SetBit(kPrimary, k);}  
-  void      SetPID(Bool_t k = kTRUE) { SetBit(kPID, k);}
   void      SetStandAlone(Bool_t st) { SetBit(kStandAlone, st); }
   void      SetPt(Double_t pt)       { fPt = pt;}
   void      SetOwner();
@@ -247,7 +244,7 @@ private:
   Float_t          fS2Z;                    // estimated resolution in the z direction 
   Float_t          fC[2];                   // Curvature for standalone [0] rieman [1] vertex constrained 
   Float_t          fChi2;                   // Global chi2  
-  Float_t          fdEdx[kNslices];         // dE/dx measurements for tracklet
+  Float_t          fdEdx[kNdEdxSlices];     // dE/dx measurements for tracklet
   Float_t          fProb[AliPID::kSPECIES]; // PID probabilities
   Int_t            fLabels[3];              // most frequent MC labels and total number of different labels
   Double_t         fRefCov[7];              // covariance matrix of the track in the yz plane + the rest of the diagonal elements

diff --git a/TRD/AliTRDtrackV1.cxx b/TRD/AliTRDtrackV1.cxx
index ef5e8f43c0a2840a06b27045f5a2d6f0b887afc5..cc423344b5640d8e31d6c3941814c9a31c1b3d5e 100644 (file)
--- a/TRD/AliTRDtrackV1.cxx
+++ b/TRD/AliTRDtrackV1.cxx
@@ -53,6 +53,7 @@ AliTRDtrackV1::AliTRDtrackV1() : AliKalmanTrack()
   ,fTruncatedMean(0)
   ,fNchamberdEdx(0)
   ,fNclusterdEdx(0)
+  ,fNdEdxSlices(0)
   ,fkReconstructor(NULL)
   ,fBackupTrack(NULL)
   ,fTrackLow(NULL)
@@ -72,6 +73,7 @@ AliTRDtrackV1::AliTRDtrackV1() : AliKalmanTrack()
     fTrackletIndex[ip] = -1;
     fTracklet[ip]      = NULL;
   }
+  SetLabel(-123456789); // reset label
 }
 
 //_______________________________________________________________
@@ -82,6 +84,7 @@ AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref) : AliKalmanTrack(ref)
   ,fTruncatedMean(ref.fTruncatedMean)
   ,fNchamberdEdx(ref.fNchamberdEdx)
   ,fNclusterdEdx(ref.fNclusterdEdx)
+  ,fNdEdxSlices(ref.fNdEdxSlices)
   ,fkReconstructor(ref.fkReconstructor)
   ,fBackupTrack(NULL)
   ,fTrackLow(NULL)
@@ -115,6 +118,7 @@ AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t) : AliKalmanTrack()
   ,fTruncatedMean(0)
   ,fNchamberdEdx(0)                                                 
   ,fNclusterdEdx(0)
+  ,fNdEdxSlices(0)
   ,fkReconstructor(NULL)
   ,fBackupTrack(NULL)
   ,fTrackLow(NULL)
@@ -171,6 +175,7 @@ AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 * const trklts, const Double_t p[5], c
   ,fTruncatedMean(0)
   ,fNchamberdEdx(0)                                                 
   ,fNclusterdEdx(0)
+  ,fNdEdxSlices(0)
   ,fkReconstructor(NULL)
   ,fBackupTrack(NULL)
   ,fTrackLow(NULL)
@@ -228,6 +233,7 @@ AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 * const trklts, const Double_t p[5], c
   
   Float_t pid = 1./AliPID::kSPECIES;
   for(int is =0; is<AliPID::kSPECIES; is++) fPID[is] = pid;
+  SetLabel(-123456789); // reset label
 }
 
 //_______________________________________________________________
@@ -378,26 +384,24 @@ Bool_t AliTRDtrackV1::CookPID()
     AliError("PID Response not available");
     return kFALSE;
   }
-  Int_t nslices = pidResponse->GetNumberOfSlices();
-  Double_t dEdx[kNplane * (Int_t)AliTRDPIDResponse::kNslicesNN];
-  Float_t trackletP[kNplane];
-  memset(dEdx, 0, sizeof(Double_t) * kNplane * (Int_t)AliTRDPIDResponse::kNslicesNN);
-  memset(trackletP, 0, sizeof(Float_t)*kNplane);
+  Double_t dEdx[kNplane * (Int_t)AliTRDseedV1::kNdEdxSlices] = {0.};
+  Float_t trackletP[kNplane] = {0.};
+
+  fNdEdxSlices = pidResponse->GetNumberOfSlices();
   for(Int_t iseed = 0; iseed < kNplane; iseed++){
     if(!fTracklet[iseed]) continue;
     trackletP[iseed] = fTracklet[iseed]->GetMomentum();
-    fTracklet[iseed]->SetPID();
-    //if(pidResponse->GetPIDmethod() == AliTRDPIDResponse::kLQ1D){
+//    if(pidResponse->GetPIDmethod() == AliTRDPIDResponse::kLQ1D){
       dEdx[iseed] = fTracklet[iseed]->GetdQdl();
-   /* } else {
-      fTracklet[iseed]->CookdEdx(nslices);
+/*    } else {
+      fTracklet[iseed]->CookdEdx(fNdEdxSlices);
       const Float_t *trackletdEdx = fTracklet[iseed]->GetdEdx();
-      for(Int_t islice = 0; islice < nslices; islice++){
-        dEdx[iseed*nslices + islice] = trackletdEdx[islice];
+      for(Int_t islice = 0; islice < fNdEdxSlices; islice++){
+        dEdx[iseed*fNdEdxSlices + islice] = trackletdEdx[islice];
       }
-  } */
+    }*/
   }
-  pidResponse->GetResponse(nslices, dEdx, trackletP, fPID);
+  pidResponse->GetResponse(fNdEdxSlices, dEdx, trackletP, fPID);
 
   static Int_t nprint = 0;
   if(!nprint){
@@ -902,19 +906,28 @@ void AliTRDtrackV1::UpdateESDtrack(AliESDtrack *track)
   // pack the two numbers together and store them in the ESD
   track->SetTRDntracklets(nPID | (nTrk<<3));
   // allocate space to store raw PID signals dEdx & momentum
-  track->SetNumberOfTRDslices((AliTRDPIDResponse::kNslicesNN+3)*AliTRDgeometry::kNlayer);
+  // independent of the method used to calculate PID (see below AliTRDPIDResponse)
+  track->SetNumberOfTRDslices((AliTRDseedV1::kNdEdxSlices+2)*AliTRDgeometry::kNlayer);
   // store raw signals
   Float_t p, sp; Double_t spd;
   for (Int_t ip = 0; ip < kNplane; ip++) {
     if(fTrackletIndex[ip]<0 || !fTracklet[ip]) continue;
-    if(!fTracklet[ip]->HasPID()) continue;
+    // Fill TRD dEdx info into ESD track
+    // a. Set Summed dEdx into the first slice
+    track->SetTRDslice(fTracklet[ip]->GetdQdl(), ip, 0); 
+    // b. Set NN dEdx slices
     fTracklet[ip]->CookdEdx(AliTRDPIDResponse::kNslicesNN);
     const Float_t *dedx = fTracklet[ip]->GetdEdx();
-    for (Int_t js = 0; js < AliTRDPIDResponse::kNslicesNN; js++, dedx++){
-      track->SetTRDslice(*dedx, ip, js+1);
+    for (Int_t js(0), ks(1); js < AliTRDPIDResponse::kNslicesNN; js++, ks++, dedx++){
+      if(ks>=AliTRDseedV1::kNdEdxSlices){
+        AliError(Form("Exceed allocated space for dEdx slices."));
+        break;
+      }
+      track->SetTRDslice(*dedx, ip, ks);
     }
-    p = fTracklet[ip]->GetMomentum(&sp);
-    spd = sp; track->SetTRDmomentum(p, ip, &spd);
+    // fill TRD momentum info into ESD track
+    p = fTracklet[ip]->GetMomentum(&sp); spd = sp;
+    track->SetTRDmomentum(p, ip, &spd);
     // store global quality per tracklet instead of momentum error
     // 26.09.11 A.Bercuci
     // first implementation store no. of time bins filled in tracklet (5bits  see "y" bits) and
@@ -933,7 +946,6 @@ void AliTRDtrackV1::UpdateESDtrack(AliESDtrack *track)
     Int_t nCross(fTracklet[ip]->IsRowCross()?fTracklet[ip]->GetTBcross():0); if(nCross>3) nCross = 3;
     Char_t trackletQ = Char_t(fTracklet[ip]->GetTBoccupancy() | (nCross<<5) | (fTracklet[ip]->IsChmbGood()<<7));
     track->SetTRDTimBin(trackletQ, ip);
-    track->SetTRDslice(fTracklet[ip]->GetdQdl(), ip, 0); // Set Summed dEdx into the first slice
   }
   // store PID probabilities
   track->SetTRDpid(fPID);

diff --git a/TRD/AliTRDtrackV1.h b/TRD/AliTRDtrackV1.h
index 3b2ddd55ca8e8334eddd6849da438fa24e6ce627..7989a0bf19e5531ad1e0d88827991fdbcf8619f0 100644 (file)
--- a/TRD/AliTRDtrackV1.h
+++ b/TRD/AliTRDtrackV1.h
@@ -34,8 +34,6 @@ public:
    ,kNplane    = AliTRDgeometry::kNlayer
    ,kNcham     = AliTRDgeometry::kNstack
    ,kNsect     = AliTRDgeometry::kNsector
-   ,kNslice    =   3
-   ,kNMLPslice =   8 
    ,kMAXCLUSTERSPERTRACK = 210
   };
   
@@ -94,6 +92,7 @@ public:
   Double_t       GetPIDsignal() const   { return 0.;}
   Double_t       GetPID(Int_t is) const { return (is >=0 && is < AliPID::kSPECIES) ? fPID[is] : -1.;}
   UChar_t        GetNumberOfTrackletsPID() const;
+  Int_t          GetNumberOfPhysicsSlices() const { return fNdEdxSlices;  };
   Double_t       GetPredictedChi2(const AliTRDseedV1 *tracklet, Double_t *cov) const;
   Double_t       GetPredictedChi2(const AliCluster* /*c*/) const                   { return 0.0; }
   Int_t          GetProlongation(Double_t xk, Double_t &y, Double_t &z) const;
@@ -106,14 +105,13 @@ public:
   AliExternalTrackParam*
                  GetTrackOut() const  { return fTrackHigh;} 
   const Int_t* GetTrackletIndexes() const { return &fTrackletIndex[0];}
-  
   Bool_t         IsEqual(const TObject *inTrack) const;
-  Bool_t         IsKink() const    { return TestBit(kKink);}
-  Bool_t         IsOwner() const   { return TestBit(kOwner);};
-  Bool_t         IsPrimary() const   { return TestBit(kPrimary);};
-  Bool_t         IsStopped() const { return TestBit(kStopped);};
+  Bool_t         IsKink() const           { return TestBit(kKink);}
+  Bool_t         IsOwner() const          { return TestBit(kOwner);};
+  Bool_t         IsPrimary() const        { return TestBit(kPrimary);};
+  Bool_t         IsStopped() const        { return TestBit(kStopped);};
   Bool_t         IsElectron() const;
-  Bool_t         IsTPCseeded() const { return !TestBit(kSeeder);};
+  Bool_t         IsTPCseeded() const      { return !TestBit(kSeeder);};
   inline static Bool_t IsTrackError(ETRDtrackError error, UInt_t status);
   inline static Bool_t IsLayerError(ETRDlayerError error, Int_t layer, UInt_t status);
 
@@ -155,6 +153,7 @@ private:
   Double32_t   fTruncatedMean;         //  Truncated mean
   Int_t        fNchamberdEdx;          //  number of chambers used in calculating truncated mean
   Int_t        fNclusterdEdx;          //  number of clusters used in calculating truncated mean
+  Int_t        fNdEdxSlices;           //  number of "physics slices" fill via AliTRDPIDResponse
 
   const AliTRDReconstructor *fkReconstructor;//! reconstructor link 
   AliTRDtrackV1 *fBackupTrack;         //! Backup track
@@ -162,7 +161,7 @@ private:
   AliExternalTrackParam *fTrackLow;    // parameters of the track which enter TRD from below (TPC) 
   AliExternalTrackParam *fTrackHigh;   // parameters of the track which enter TRD from above (HMPID, PHOS) 
 
-  ClassDef(AliTRDtrackV1, 7)          // TRD track - tracklet based
+  ClassDef(AliTRDtrackV1, 8)           // TRD track - tracklet based
 };
 
 //____________________________________________________