]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDgtuParam.cxx
sorry i forget header
[u/mrichter/AliRoot.git] / TRD / AliTRDgtuParam.cxx
index e73351a77eaa999f0b7126d520facf08be000fbc..6ae8e075ef8508a6aa9c5587d1e07ce39b408e03 100644 (file)
@@ -23,6 +23,7 @@
 //                                                                        //
 ////////////////////////////////////////////////////////////////////////////
 
+#include "TROOT.h"
 #include "TMath.h"
 #include "TMatrix.h"
 #include "TDecompLU.h"
@@ -64,6 +65,7 @@ AliTRDgtuParam::AliTRDgtuParam() :
   fVertexSize(20.0),
   fCurrTrackletMask(0),
   fRefLayers(0x0),
+  fMagField(0.5),
   fGeo(0x0)
 {
   // default ctor
@@ -72,7 +74,13 @@ AliTRDgtuParam::AliTRDgtuParam() :
   fRefLayers[0] = 3;
   fRefLayers[1] = 2;
   fRefLayers[2] = 1;
-  zChannelGen(); 
+  for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
+    fAki[iLayer] = 0.;
+    fBki[iLayer] = 0.;
+    fCki[iLayer] = 0.;
+  }
+
+  GenerateZChannelMap(); 
 }
 
 AliTRDgtuParam::~AliTRDgtuParam() 
@@ -115,13 +123,15 @@ Int_t AliTRDgtuParam::GetZSubchannel(Int_t stack, Int_t layer, Int_t zchannel, I
 
 Int_t AliTRDgtuParam::GetRefLayer(Int_t refLayerIdx) const 
 {
+  // returns the reference layer indexed by refLayerIdx
+
   if (refLayerIdx >= 0 && refLayerIdx < fgkNRefLayers)
     return fRefLayers[refLayerIdx];
   else 
     return -1;
 }
 
-Int_t AliTRDgtuParam::zChannelGen() 
+Int_t AliTRDgtuParam::GenerateZChannelMap() 
 {
   // generate the z-channel map
   // assuming that the tracks come from the vertex 
@@ -132,15 +142,15 @@ Int_t AliTRDgtuParam::zChannelGen()
 
   for (Int_t iStack = 0; iStack < fGeo->Nstack(); iStack++) {
 
-    Float_t X[6] = { 0 };
-    Float_t Z[6][16] = {{ 0 }};
+    Float_t x[6] = { 0 };
+    Float_t z[6][16] = {{ 0 }};
     Float_t dZ[6][16] = {{ 0 }};
     
     for (Int_t iLayer = 0; iLayer < fGeo->Nlayer(); iLayer++) {
       AliTRDpadPlane *pp = fGeo->GetPadPlane(iLayer, iStack);
-      X[iLayer]  = fGeo->GetTime0(iLayer) - fGeo->CdrHght(); // ???
+      x[iLayer]  = fGeo->GetTime0(iLayer) - fGeo->CdrHght(); // ???
       for (Int_t iRow = 0; iRow < fGeo->GetRowMax(iLayer, iStack, iSec); iRow++) {
-       Z[iLayer][iRow]  = pp->GetRowPos(iRow); // this is the right (pos. z-direction) border of the pad
+       z[iLayer][iRow]  = pp->GetRowPos(iRow); // this is the right (pos. z-direction) border of the pad
        dZ[iLayer][iRow] = pp->GetRowSize(iRow); // length of the pad in z-direction
        for (Int_t i = 0; i < fgkNZChannels; i++) 
            fZSubChannel[iStack][i][iLayer][iRow] = 0;
@@ -149,26 +159,26 @@ Int_t AliTRDgtuParam::zChannelGen()
 
     for (Int_t fixRow = 0; fixRow < fGeo->GetRowMax(fgkFixLayer, iStack, iSec); fixRow++) {
        
-      Double_t fixZmin = Z[fgkFixLayer][fixRow] - dZ[fgkFixLayer][fixRow];  
-      Double_t fixZmax = Z[fgkFixLayer][fixRow];
-      Double_t fixX    = X[fgkFixLayer] + 1.5; // ??? 1.5 from where? 
+      Double_t fixZmin = z[fgkFixLayer][fixRow] - dZ[fgkFixLayer][fixRow];  
+      Double_t fixZmax = z[fgkFixLayer][fixRow];
+      Double_t fixX    = x[fgkFixLayer] + 1.5; // ??? 1.5 from where? 
 
       for (Int_t iLayer = 0; iLayer < fGeo->Nlayer(); iLayer++) {
        Double_t leftZ, rightZ;
        
        if (iLayer <= fgkFixLayer) {
-         leftZ  = (fixZmin + fVertexSize) * (X[iLayer] + 1.5) / fixX - fVertexSize;
-         rightZ = (fixZmax - fVertexSize) * (X[iLayer] + 1.5) / fixX + fVertexSize;
+         leftZ  = (fixZmin + fVertexSize) * (x[iLayer] + 1.5) / fixX - fVertexSize;
+         rightZ = (fixZmax - fVertexSize) * (x[iLayer] + 1.5) / fixX + fVertexSize;
        }
        else {
-         leftZ  = (fixZmin - fVertexSize) * (X[iLayer] + 1.5) / fixX + fVertexSize;
-         rightZ = (fixZmax + fVertexSize) * (X[iLayer] + 1.5) / fixX - fVertexSize;
+         leftZ  = (fixZmin - fVertexSize) * (x[iLayer] + 1.5) / fixX + fVertexSize;
+         rightZ = (fixZmax + fVertexSize) * (x[iLayer] + 1.5) / fixX - fVertexSize;
        }
        
        Double_t epsilon = 0.001;
        for (Int_t iRow = 0; iRow < fGeo->GetRowMax(iLayer, iStack, iSec); iRow++) {
-         if ( (Z[iLayer][iRow] )                    > (leftZ  + epsilon) && 
-              (Z[iLayer][iRow] - dZ[iLayer][iRow] ) < (rightZ - epsilon) ) {
+         if ( (z[iLayer][iRow] )                    > (leftZ  + epsilon) && 
+              (z[iLayer][iRow] - dZ[iLayer][iRow] ) < (rightZ - epsilon) ) {
            fZChannelMap[iStack][fixRow][iLayer][iRow] = 1;
            if (fZSubChannel[iStack][fixRow % fgkNZChannels][iLayer][iRow] != 0) {
              AliError("Collision in Z-Channel assignment occured! No reliable tracking!!!");
@@ -190,7 +200,7 @@ Bool_t AliTRDgtuParam::DisplayZChannelMap(Int_t zchannel, Int_t subchannel) cons
 {
   // display the z-channel map 
 
-  if (zchannel > fgkNZChannels) {
+  if (zchannel >= fgkNZChannels) {
     AliError("Invalid Z channel!");
     return kFALSE;
   }
@@ -224,12 +234,15 @@ Bool_t AliTRDgtuParam::DisplayZChannelMap(Int_t zchannel, Int_t subchannel) cons
   }
   graph->SetMarkerStyle(kDot);
   graph->Draw("AP");
+  gROOT->Add(graph);
   for (Int_t zch = zchmin; zch < zchmax; zch++) {
     graphz[zch]->SetMarkerStyle(kCircle);
     graphz[zch]->SetMarkerColor(zch+2);
     graphz[zch]->SetMarkerSize(0.3 + zch*0.2);
     graphz[zch]->Draw("P");
+    gROOT->Add(graphz[zch]);
   }
+  delete [] graphz;
   return kTRUE;
 }
 
@@ -237,17 +250,17 @@ Int_t AliTRDgtuParam::GetCiAlpha(Int_t layer) const
 {
   // get the constant for the calculation of alpha
 
-  Int_t Ci = (Int_t) (GetChamberThickness() / fGeo->GetTime0(layer) * GetBinWidthY() / GetBinWidthdY() * (1 << (GetBitExcessAlpha() + GetBitExcessY() + 1)) );
-  return Ci;
+  Int_t ci = (Int_t) (GetChamberThickness() / fGeo->GetTime0(layer) * GetBinWidthY() / GetBinWidthdY() * (1 << (GetBitExcessAlpha() + GetBitExcessY() + 1)) );
+  return ci;
 }
 
 Int_t AliTRDgtuParam::GetCiYProj(Int_t layer) const 
 {
   // get the constant for the calculation of y_proj
 
-  Float_t Xmid = (fGeo->GetTime0(0) + fGeo->GetTime0(fGeo->Nlayer()-1)) / 2.; 
-  Int_t Ci = (Int_t) (- (fGeo->GetTime0(layer) - Xmid) / GetChamberThickness() * GetBinWidthdY() / GetBinWidthY() * (1 << GetBitExcessYProj()) );
-  return Ci;
+  Float_t xmid = (fGeo->GetTime0(0) + fGeo->GetTime0(fGeo->Nlayer()-1)) / 2.; 
+  Int_t ci = (Int_t) (- (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
@@ -259,6 +272,9 @@ Int_t AliTRDgtuParam::GetYt(Int_t stack, Int_t layer, Int_t zrow) const
 
 Bool_t AliTRDgtuParam::GenerateRecoCoefficients(Int_t trackletMask) 
 {
+  // calculate the coefficients for the straight line fit 
+  // depending on the mask of contributing tracklets
+
   fCurrTrackletMask = trackletMask;
 
   TMatrix a(GetNLayers(), 3);
@@ -385,28 +401,23 @@ Bool_t AliTRDgtuParam::GetIntersectionPoints(Int_t k, Float_t &x1, Float_t &x2)
     }
   }
 
-  x1 = fGeo->GetTime0(l1) + 10./6 * (nHits -1);
-  x2 = fGeo->GetTime0(l2) - 10./6 * (nHits -1);
-
-  return ( (l1 >= 0) && (l2 >= 0) );
+  if ( (l1 >= 0) && (l2 >= 0) ) {
+    x1 = fGeo->GetTime0(l1) + 10./6 * (nHits -1);
+    x2 = fGeo->GetTime0(l2) - 10./6 * (nHits -1);
+    return kTRUE;
+  }
+  else
+    return kFALSE;
 }
 
-Float_t AliTRDgtuParam::GetRadius(Int_t a, Float_t b, Float_t x1, Float_t x2) 
+Float_t AliTRDgtuParam::GetPt(Int_t a, Float_t /* b */, Float_t x1, Float_t x2) const
 {
-  // get the radius for the track
-  Float_t d = (1 + b * b /2 ) * (x2 - x1);
-  Float_t c1 = x1 * x2 / 2;
-//  Float_t c2 = (x1 + x2) / (x1 * x2);
-  printf("c1: %f\n", c1);
+  // returns 0.3 * B * 1/a (1/128 GeV/c)
+  // a : offset, b : slope (not used)
+
+  Float_t c1 = x1 * x2 / 2. / 10000.; // conversion cm to m
   Float_t r = 0;
   if ( (a >> 1) != 0)
-    r = (375. / 10000.) * c1 * 256 / (a >> 1);
-  return r;
-
-  Float_t y1 = a + b*x1;
-  Float_t y2 = a + b*x2;
-  Float_t alpha = TMath::Abs( TMath::ATan(y2/x2) - TMath::ATan(y1/x1) );
-  d = TMath::Sqrt( TMath::Power(x2-x1, 2) + TMath::Power(y2-y1, 2) );
-  r = d / 2. / TMath::Sin(alpha);
+    r = (0.3 * fMagField / 2. / (fgkBinWidthY/100.)) * (((Int_t) c1) << 8) / (a >> 1); //??? why shift of a?
   return r;
 }