// //
////////////////////////////////////////////////////////////////////////////
+#include <limits>
+
#include "TROOT.h"
#include "TMath.h"
#include "TMatrix.h"
AliTRDgtuParam *AliTRDgtuParam::fgInstance = 0;
Bool_t AliTRDgtuParam::fgUseGTUconst = kTRUE;
+Bool_t AliTRDgtuParam::fgUseGTUmerge = kTRUE;
+Bool_t AliTRDgtuParam::fgLimitNoTracklets = kTRUE;
// ----- matching windows -----
Int_t AliTRDgtuParam::fgDeltaY = 19;
const Int_t AliTRDgtuParam::fgkBitExcessAlpha = 10;
const Int_t AliTRDgtuParam::fgkBitExcessYProj = 2;
-// ----- z-channel tables -----
+// pt higher than the one for smallest possible a != 0
+const Int_t AliTRDgtuParam::fgkPtInfinity = std::numeric_limits<Int_t>::max();
+
+// ----- geometry constants used in GTU -----
const Bool_t AliTRDgtuParam::fgZChannelMap[5][16][6][16] = {
{ /* --- Stack 0 --- */
{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1}}}
};
+const Float_t AliTRDgtuParam::fgkRadius[6] = { 300.65, 313.25, 325.85, 338.45, 351.05, 363.65 };
+const Float_t AliTRDgtuParam::fgkThickness = 3.;
+const Float_t AliTRDgtuParam::fgkRow0Pos[6][5] = {
+ {301, 177, 53, -57, -181},
+ {301, 177, 53, -57, -181},
+ {315, 184, 53, -57, -188},
+ {329, 191, 53, -57, -195},
+ {343, 198, 53, -57, -202},
+ {347, 200, 53, -57, -204}
+};
+const Float_t AliTRDgtuParam::fgkInnerPadLength[] = {7.5, 7.5, 8.0, 8.5, 9.0, 9.0};
+const Float_t AliTRDgtuParam::fgkOuterPadLength[] = {7.5, 7.5, 7.5, 7.5, 7.5, 8.5};
+const Float_t AliTRDgtuParam::fgkAcoeff[32][6] = {
+ {-3440, -3303, 3174, 3057, 0, 0},
+ {-3481, 0, -171, 0, 3140, 0},
+ {-2850, -1380, 0, 1277, 2441, 0},
+ {-3481, 0, -171, 0, 3140, 0},
+ { 0, -3568, -3431, 3303, 3185, 0},
+ {-2783, -1378, -136, 1275, 2510, 0},
+ {-1500, -2857, 1384, 0, 0, 2461},
+ { 0, -3609, 0, -171, 0, 3268},
+ {-3685, 0, 3400, -3276, 0, 3049},
+ { 0, -3609, 0, -171, 0, 3268},
+ {-1498, -2792, 1382, -132, 0, 2528},
+ {-1850, -1777, 0, 0, 1585, 1531},
+ {-3481, 0, -171, 0, 3140, 0},
+ { 0, -2953, -1431, 0, 1328, 2544},
+ {-1808, -1776, -89, 0, 1631, 1530},
+ {-2932, 0, 0, -1314, 2511, 1223},
+ { 0, -3609, 0, -171, 0, 3268},
+ {-1849, -1738, 0, -82, 1583, 1574},
+ { 0, 0, -3696, -3559, 3431, 3313},
+ {-2863, 0, -140, -1312, 2582, 1221},
+ { 0, -2886, -1429, -136, 1327, 2613},
+ {-1806, -1736, -89, -82, 1629, 1572},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1},
+ { -1, -1, -1, -1, -1, -1}
+};
+const Int_t AliTRDgtuParam::fgkMaskID[] = {
+ -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, 0,
+ -1, -1, -1, -1, -1, -1, -1, 1, -1, -1, -1, 2, -1, 3, 4, 5,
+ -1, -1, -1, -1, -1, -1, -1, 6, -1, -1, -1, 7, -1, 8, 9, 10,
+ -1, -1, -1, 11, -1, 12, 13, 14, -1, 15, 16, 17, 18, 19, 20, 21
+};
AliTRDgtuParam::AliTRDgtuParam() :
fVertexSize(20.0),
{
// get the constant for the calculation of y_proj
- Float_t xmid = (fGeo->GetTime0(0) + fGeo->GetTime0(fGeo->Nlayer()-1)) / 2.;
- Int_t ci = TMath::Nint(- (fGeo->GetTime0(layer) - xmid) / GetChamberThickness() * GetBinWidthdY() / GetBinWidthY() * (1 << GetBitExcessYProj()) );
+ Int_t ci = 0;
+
+ if (fgUseGTUconst) {
+ Float_t xmid = (fgkRadius[0] + fgkRadius[5]) / 2.;
+ ci = TMath::Nint(- (fgkRadius[layer] - xmid) * fgkBinWidthdY / (fgkBinWidthY * fgkThickness) * (1 << GetBitExcessYProj()));
+ } else {
+ Float_t xmid = (fGeo->GetTime0(0) + fGeo->GetTime0(fGeo->Nlayer()-1)) / 2.;
+ ci = TMath::Nint(- (fGeo->GetTime0(layer) - xmid) / GetChamberThickness() * GetBinWidthdY() / GetBinWidthY() * (1 << GetBitExcessYProj()) );
+ }
+
return ci;
}
Int_t AliTRDgtuParam::GetYt(Int_t stack, Int_t layer, Int_t zrow) const
{
- return (Int_t) (- ( (layer % 2 ? 1 : -1) *
- (GetGeo()->GetPadPlane(layer, stack)->GetRowPos(zrow) - GetGeo()->GetPadPlane(layer, stack)->GetRowSize(zrow) / 2) *
- TMath::Tan(- 2.0 / 180.0 * TMath::Pi()) ) / 0.016 );
+ // return yt for the calculation of y'
+
+ Int_t yt = 0;
+
+ if (fgUseGTUconst) {
+ yt = TMath::Nint (- ( (layer % 2 ? 1. : -1.) *
+ GetZrow(stack, layer, zrow) * TMath::Tan(- 2./180. * TMath::Pi()) / fgkBinWidthY ));
+ } else {
+ yt = TMath::Nint (- ( (layer % 2 ? 1. : -1.) *
+ (GetGeo()->GetPadPlane(layer, stack)->GetRowPos(zrow) - GetGeo()->GetPadPlane(layer, stack)->GetRowSize(zrow) / 2.) *
+ TMath::Tan(- 2./180. * TMath::Pi()) ) / fgkBinWidthY );
+ }
+
+ return yt;
}
Bool_t AliTRDgtuParam::GenerateRecoCoefficients(Int_t trackletMask)
return kTRUE;
}
-Float_t AliTRDgtuParam::GetAki(Int_t k, Int_t i)
+Int_t AliTRDgtuParam::GetAki(Int_t k, Int_t i)
{
// get A_ki for the calculation of the tracking parameters
- if (fCurrTrackletMask != k)
- GenerateRecoCoefficients(k);
-
- return -fAki[i];
+ if (fgUseGTUconst) {
+ Int_t maskId = fgkMaskID[k];
+ return fgkAcoeff[maskId][i];
+ } else {
+ if (fCurrTrackletMask != k)
+ GenerateRecoCoefficients(k);
+ return -(((Int_t) fAki[i]) << 9);
+ }
}
Float_t AliTRDgtuParam::GetBki(Int_t k, Int_t i)
// returns 0.3 * B * 1/a (1/128 GeV/c)
// a : offset, b : slope (not used)
+ // protect against division by zero, covers both cases
+ if ((a >> 2) == 0)
+ return fgkPtInfinity;
+
if (fgUseGTUconst) {
//----- calculation as in the GTU ----
const Int_t maskIdLut[64] = {
static Int_t GetBitExcessAlpha() { return fgkBitExcessAlpha; }
static Int_t GetBitExcessYProj() { return fgkBitExcessYProj; }
+ Float_t GetInnerPadLength(Int_t stack, Int_t layer) const {
+ return (stack == 2) ? 9. : fgkInnerPadLength[layer];
+ }
+ Float_t GetOuterPadLength(Int_t stack, Int_t layer) const {
+ return (stack == 2) ? 8. : fgkOuterPadLength[layer];
+ }
+ Float_t GetZrow(Int_t stack, Int_t layer, Int_t padrow) const {
+ Float_t zRowCorrected = fgkRow0Pos[layer][stack] - GetOuterPadLength(stack, layer) + GetInnerPadLength(stack, layer);
+ return zRowCorrected - (0.5 + padrow) * GetInnerPadLength(stack, layer);
+ }
+
AliTRDgeometry* GetGeo() const { return fGeo; }
Float_t GetVertexSize() const { return fVertexSize; }
Int_t GetCiAlpha(Int_t layer) const;
static void SetUseGTUconst(Bool_t b) { fgUseGTUconst = b; }
static Bool_t GetUseGTUconst() { return fgUseGTUconst; }
+ static void SetUseGTUmerge(Bool_t b) { fgUseGTUmerge = b; }
+ static Bool_t GetUseGTUmerge() { return fgUseGTUmerge; }
+
+ static void SetLimitNoTracklets(Bool_t b) { fgLimitNoTracklets = b; }
+ static Bool_t GetLimitNoTracklets() { return fgLimitNoTracklets; }
+
// z-channel map
Int_t GenerateZChannelMap(); // could have different modes (for beam-beam, cosmics, ...)
Bool_t DisplayZChannelMap(Int_t zchannel = -1, Int_t subch = 0) const;
// variables for pt-reconstruction (not used at the moment)
Bool_t GenerateRecoCoefficients(Int_t trackletMask);
- Float_t GetAki(Int_t k, Int_t i);
+ Int_t GetAki(Int_t k, Int_t i);
Float_t GetBki(Int_t k, Int_t i);
Float_t GetCki(Int_t k, Int_t i);
// Float_t GetD(Int_t k) const;
static const Int_t fgkBitExcessAlpha; // excess bits for alpha
static const Int_t fgkBitExcessYProj; // excess bits for projected y-position
+ static const Int_t fgkPtInfinity; // infinite pt as obtained when a == 0
+
protected:
static Int_t fgDeltaY; // accepted deviation in y_proj, default: 9
static Int_t fgDeltaAlpha; // accepted deviation in alpha, default: 11
static Bool_t fgUseGTUconst; // use constants as in the GTU for the calculations
// instead of geometry derived quantities
+ static Bool_t fgUseGTUmerge; // use merge algorithm exactly as in hardware
+ static Bool_t fgLimitNoTracklets; // limit the number of tracklets per layer
static const Bool_t fgZChannelMap[5][16][6][16]; // z-channel tables as in GTU
+ static const Float_t fgkRadius[6]; // layer radius as used in the GTU code
+ static const Float_t fgkThickness; // drift length as used in the GTU code
+ static const Float_t fgkRow0Pos[6][5]; // geometry constant from GTU implementation
+ static const Float_t fgkInnerPadLength[6]; // geometry constant from GTU implementation
+ static const Float_t fgkOuterPadLength[6]; // geometry constant from GTU implementation
+ static const Float_t fgkAcoeff[32][6]; // geometry constant from GTU implementation
+ static const Int_t fgkMaskID[64]; // geometry constant from GTU implementation
Float_t fVertexSize; // assumed vertex size (z-dir.) for the z-channel map
else
trk->SetAlpha(alpha);
- Int_t yproj = trk->GetdY() * (fGtuParam->GetCiYProj(layer)); //??? sign?
+ Int_t yproj = trk->GetdY() * (fGtuParam->GetCiYProj(layer));
yproj = ( ( ( (yproj >> fGtuParam->GetBitExcessYProj()) + trk->GetYbin() ) >> 2) + 1) >> 1;
trk->SetYProj(yproj);
trk->SetYPrime(trk->GetYbin() + fGtuParam->GetYt(fStack, layer, trk->GetZbin()));
- AliDebug(10, Form("0x%08x: idx: %3i, z: %2i, y: %5i, dy: %3i, y': %5i, y_proj: %5i, alpha: %3i, pid: %3i",
+ AliDebug(10, Form("0x%08x: idx: %3i, z: %2i, y: %5i, dy: %3i, y': %5i, y_proj: %5i, alpha: %3i, pid: %3i, c: %5i, yt: %5i",
trk->GetTrackletWord(), trk->GetIndex(), trk->GetZbin(), trk->GetYbin(), trk->GetdY(), trk->GetYPrime(),
- trk->GetYProj(), trk->GetAlpha(), trk->GetPID() ));
+ trk->GetYProj(), trk->GetAlpha(), trk->GetPID(),
+ fGtuParam->GetCiYProj(layer), fGtuParam->GetYt(fStack, layer, trk->GetZbin()) ));
}
return kTRUE;
}
minIdx = refLayerIdx;
// done = kFALSE;
}
- else if (trkInRefLayer[refLayerIdx]->GetZSubChannel() < trkStage0->GetZSubChannel() ||
- (trkInRefLayer[refLayerIdx]->GetZSubChannel() == trkStage0->GetZSubChannel() && trkInRefLayer[refLayerIdx]->GetYapprox() < trkStage0->GetYapprox()) ) {
+ else if ( (trkInRefLayer[refLayerIdx]->GetZSubChannel() < trkStage0->GetZSubChannel()) ||
+ ((trkInRefLayer[refLayerIdx]->GetZSubChannel() == trkStage0->GetZSubChannel()) &&
+ ((trkInRefLayer[refLayerIdx]->GetYapprox() >> 3) < (trkStage0->GetYapprox() >> 3)) ) ) {
minIdx = refLayerIdx;
trkStage0 = trkInRefLayer[refLayerIdx];
// done = kFALSE;
// ----- Merging in zchannels - 1st stage -----
- do {
- // done = kTRUE;
- trkStage0 = 0x0;
- for (Int_t zch = fGtuParam->GetNZChannels() - 1; zch > -1; zch--) {
- AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksRefUnique[zch]->First();
- if (trk == 0) {
- continue;
- }
- else if (trkStage0 == 0x0 ) {
- trkStage0 = trk;
- minIdx = zch;
- // done = kFALSE;
+ if (AliTRDgtuParam::GetUseGTUmerge()) {
+ Int_t notEmpty;
+ do {
+ Bool_t lowerThan[3] = { kFALSE, kFALSE, kFALSE };
+ AliTRDtrackGTU *trk[3] = { 0x0, 0x0, 0x0 };
+ for (Int_t iChannel = 0; iChannel < fGtuParam->GetNZChannels(); ++iChannel)
+ trk[iChannel] = (AliTRDtrackGTU*) tracksRefUnique[iChannel]->First();
+
+ for (Int_t iChannel = 0; iChannel < fGtuParam->GetNZChannels(); ++iChannel) {
+ AliTRDtrackGTU *trk1 = trk[iChannel];
+ AliTRDtrackGTU *trk2 = trk[(iChannel + 1) % 3];
+ if (trk1 && trk2) {
+ Int_t sortnum1 = (trk1->GetZChannel() + 3 * trk1->GetZSubChannel()) / 2 - 1;
+ Int_t sortnum2 = (trk2->GetZChannel() + 3 * trk2->GetZSubChannel()) / 2 - 1;
+ AliDebug(5, Form("comparing tracks %i - %i: %i - %i",
+ trk1->GetZChannel(), trk2->GetZChannel(),
+ sortnum1, sortnum2));
+ if ( (sortnum1 < sortnum2) ||
+ ((sortnum1 == sortnum2) &&
+ ((trk1->GetYapprox() >> 3) < (trk2->GetYapprox() >> 3)) ) ) {
+ lowerThan[iChannel] = kTRUE;
}
- else if ( ((trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1) < ((trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 -1 ) ||
- (((trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1) == ((trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 -1 ) && (trk->GetYapprox() < trkStage0->GetYapprox()) ) ) {
- minIdx = zch;
- trkStage0 = trk;
- // done = kFALSE;
+ }
+ }
+
+ notEmpty = (trk[2] ? (1 << 2) : 0) |
+ (trk[1] ? (1 << 1) : 0) |
+ (trk[0] ? (1 << 0) : 0);
+ Int_t pop = -1;
+
+ switch (notEmpty) {
+ // one track only
+ case 0x1:
+ pop = 0;
+ break;
+ case 0x2:
+ pop = 1;
+ break;
+ case 0x4:
+ pop = 2;
+ break;
+
+ // two tracks
+ case 0x3:
+ if (lowerThan[0])
+ pop = 0;
+ else
+ pop = 1;
+ break;
+ case 0x5:
+ if (lowerThan[2])
+ pop = 2;
+ else
+ pop = 0;
+ break;
+ case 0x6:
+ if (lowerThan[1])
+ pop = 1;
+ else
+ pop = 2;
+ break;
+
+ // three tracks
+ case 0x7:
+ if (lowerThan[0]) {
+ if (lowerThan[2])
+ pop = 2;
+ else
+ pop = 0;
+ } else {
+ if (lowerThan[1])
+ pop = 1;
+ else
+ pop = 2;
+ }
+ break;
+
+ // no tracks
+ default:
+ // nop
+ ;
+ }
+
+ if (pop > -1) {
+ tracksZMergedStage0->Add(trk[pop]);
+ tracksRefUnique[pop]->RemoveFirst();
+ }
+ } while (notEmpty);
+ }
+ else {
+ // there is still a difference between this implementation and
+ // the hardware algorithm, only for expert use
+
+ do {
+ // done = kTRUE;
+ trkStage0 = 0x0;
+ // compare tracks from all adjacent zchannels
+ // (including comparison of channels 2 and 0)
+ for (Int_t i = fGtuParam->GetNZChannels() - 1; i > -1; i--) {
+ Int_t zch = i % 3;
+ AliTRDtrackGTU *trk = (AliTRDtrackGTU*) tracksRefUnique[zch]->First();
+ if (trk == 0) {
+ continue;
+ }
+ else if (trkStage0 == 0x0 ) {
+ trkStage0 = trk;
+ minIdx = zch;
+ // done = kFALSE;
+ }
+ else {
+ Int_t sortnum1 = (trk->GetZChannel() + 3 * trk->GetZSubChannel()) / 2 - 1;
+ Int_t sortnum2 = (trkStage0->GetZChannel() + 3 * trkStage0->GetZSubChannel()) / 2 - 1;
+ AliDebug(5, Form("comparing tracks %i - %i: %i - %i",
+ trk->GetZChannel(), trkStage0->GetZChannel(),
+ sortnum1, sortnum2));
+ if ( (sortnum1 < sortnum2) ||
+ ((sortnum1 == sortnum2) &&
+ ((trk->GetYapprox() >> 3) < (trkStage0->GetYapprox() >> 3)) ) ) {
+ minIdx = zch;
+ trkStage0 = trk;
+ // done = kFALSE;
}
+ }
}
if (!trkStage0)
break;
tracksZMergedStage0->Add(trkStage0);
tracksRefUnique[minIdx]->RemoveFirst();
- } while (trkStage0 != 0);
+ } while (trkStage0 != 0);
+ }
Uniquifier(tracksZMergedStage0, tracksZUniqueStage0);
minIdx = i;
// done = kFALSE;
}
- else if ( (((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) < ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) ||
- ((((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) == ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) && (trk->GetYapprox() < trkStage0->GetYapprox()) ) ) {
+ else if ( (((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) < ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) ||
+ ((((trk->GetZChannel() + 3 * (trk->GetZSubChannel() - 1)) / 2) == ((trkStage0->GetZChannel() + 3 * (trkStage0->GetZSubChannel() - 1)) / 2)) &&
+ ((trk->GetYapprox() >> 3) < (trkStage0->GetYapprox() >> 3)) ) ) {
minIdx = i;
trkStage0 = trk;
// done = kFALSE;
if (AliTRDgtuParam::GetUseGTUconst()) {
// averaging as in GTU
+ AliDebug(1, "using GTU constants for PID calculation");
ULong64_t coeff;
// selection of coefficient for averaging
coeff &= 0x1ffff; // 17-bit constant
ULong64_t sum = 0;
+ Int_t i = 0;
for (Int_t iLayer = 0; iLayer < fGtuParam->GetNLayers(); iLayer++) {
if ((track->GetTrackletMask() >> iLayer) & 1) {
sum += track->GetTracklet(iLayer)->GetPID();
+ ++i;
}
}
+ Float_t av = 1./i * sum;
sum = sum & 0x7ff;
ULong64_t prod = (sum * coeff) & 0xfffffffffull;
ULong64_t prodFinal = ((prod >> 17) + ((prod >> 16) & 1)) & 0xff;
+ if (TMath::Abs((prodFinal & 0xff) - av) > 0.5)
+ AliError(Form("problem with PID averaging (hw <-> ar): %3lli <-> %4.1f", prodFinal & 0xff, av));
track->SetPID(prodFinal & 0xff);
return kTRUE;
AliError(Form("Could not get tracklet in layer %i\n", layer));
continue;
}
- AliDebug(10,Form(" layer %i trk yprime: %6i, aki: %6i", layer, trk->GetYPrime(), (Int_t) (2048 * fGtuParam->GetAki(track->GetTrackletMask(), layer))));
- a += (((Int_t) (2048 * fGtuParam->GetAki(track->GetTrackletMask(), layer))) * trk->GetYPrime() + 1) >> 8;
+ AliDebug(10,Form(" layer %i trk yprime: %6i, aki: %6i", layer, trk->GetYPrime(),
+ fGtuParam->GetAki(track->GetTrackletMask(), layer)));
+ a += (((fGtuParam->GetAki(track->GetTrackletMask(), layer) * trk->GetYPrime()) >> 7) + 1) >> 1;
b += fGtuParam->GetBki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
c += fGtuParam->GetCki(track->GetTrackletMask(), layer) * trk->GetYPrime() * fGtuParam->GetBinWidthY();
}
a += 3;
a = a >> 2;
- track->SetFitParams(a, b, c);
+ track->SetFitParams(a << 2, b, c);
fGtuParam->GetIntersectionPoints(track->GetTrackletMask(), x1, x2);
- AliDebug(5,Form(" Track parameters: a = %i, b = %f, c = %f, x1 = %f, x2 = %f, pt = %f (trkl mask: %i)", a, b, c, x1, x2, track->GetPt(), track->GetTrackletMask()));
+ AliDebug(5,Form(" Track parameters: a16 = %i, a18 = %i, b = %f, c = %f, x1 = %f, x2 = %f, pt = %f (trkl mask: %i)",
+ a, a << 2, b, c, x1, x2, track->GetPt(), track->GetTrackletMask()));
return kTRUE;
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff