TRD/AliTRDdisplayDigits3D.C

   1 //_____________________________________________________________________________
   2 Int_t AliTRDdisplayDigits3D(Int_t event = 0, Int_t thresh = 0
   3                           , Bool_t sdigits = kFALSE)
   4 {
   5   //
   6   //  TRD digits display
   7   //
   8   //  Input parameter:
   9   //    <event>   : Event number
  10   //    <thresh>  : Threshold to suppress the noise
  11   //    <sdigits> : If kTRUE it will display summable digits, normal digits otherwise.
  12   //                The signal event is displayed in yellow.
  13   //
  14
  15   Char_t *inputFile = "galice.root";
  16
  17   const Int_t kNdict = 3;
  18
  19   // Define the objects
  20   AliTRDv1       *trd;
  21   AliTRDgeometry *geo;
  22
  23   Int_t           track;
  24   Int_t           idict;
  25
  26   // Connect the AliRoot file containing Geometry, Kine, Hits, and digits
  27   TFile *gafl = (TFile*) gROOT->GetListOfFiles()->FindObject(inputFile);
  28   if (!gafl) {
  29     printf("Open the ALIROOT-file %s\n",inputFile);
  30     gafl = new TFile(inputFile);
  31   }
  32   else {
  33     printf("%s is already open",inputFile);
  34   }
  35
  36   // Get AliRun object from file or create it if not on file
  37   gAlice = (AliRun *) gafl->Get("gAlice");
  38   if (!gAlice) return 1;
  39
  40   // Import the Trees for the event nEvent in the file
  41   Int_t nparticles = gAlice->GetEvent(event);
  42   if (nparticles <= 0) return 1;
  43
  44   // Get the pointer to the detector object
  45   trd = (AliTRDv1*) gAlice->GetDetector("TRD");
  46
  47   // Get the pointer to the geometry object
  48   if (trd) {
  49     geo = trd->GetGeometry();
  50   }
  51   else {
  52     printf("Cannot find the geometry\n");
  53     return 1;
  54   }
  55
  56   TCanvas *c1 = new TCanvas("digits","TRD digits display",0,0,700,730);
  57   TView   *v  = new TView(1);
  58   v->SetRange(-430,-560,-430,430,560,1710);
  59   c1->Clear();
  60   c1->SetFillColor(1);
  61   c1->SetTheta(90.0);
  62   c1->SetPhi(0.0);
  63
  64   Int_t markerColorSignal = 2;
  65   Int_t markerColorBgnd   = 7;
  66   Int_t markerColorMerged = 5;
  67   Int_t mask              = 10000000;
  68
  69   // Create the digits manager
  70   AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager();
  71   digitsManager->SetDebug(1);
  72   digitsManager->SetSDigits(sdigits);
  73
  74   // Read the digits from the file
  75   digitsManager->Open(inputFile);
  76   digitsManager->ReadDigits();
  77
  78   Int_t totalsignal = 0;
  79   Int_t totalbgnd   = 0;
  80   Int_t totalmerged = 0;
  81
  82   // Loop through all detectors
  83   for (Int_t idet = 0; idet < geo->Ndet(); idet++) {
  84
  85     printf("<AliTRDdisplayDigits3D> Loading detector %d\n",idet);
  86     AliTRDdataArrayI *digits  = digitsManager->GetDigits(idet);
  87     digits->Expand();
  88     AliTRDdataArrayI *tracks[kNdict];
  89     for (Int_t idict = 0; idict < kNdict; idict++) {
  90       tracks[idict] = digitsManager->GetDictionary(idet,idict);
  91       tracks[idict]->Expand();
  92     }
  93
  94     Int_t isec    = geo->GetSector(idet);
  95     Int_t icha    = geo->GetChamber(idet);
  96     Int_t ipla    = geo->GetPlane(idet);
  97     Int_t  rowMax = geo->GetRowMax(ipla,icha,isec);
  98     Int_t  colMax = geo->GetColMax(ipla);
  99     Int_t timeMax = geo->GetTimeMax();
 100
 101     Int_t ndigits = digits->GetOverThreshold(thresh);
 102
 103     if (ndigits > 0) {
 104
 105       TPolyMarker3D *pmSignal = new TPolyMarker3D(ndigits);
 106       TPolyMarker3D *pmBgnd   = new TPolyMarker3D(ndigits);
 107       TPolyMarker3D *pmMerged = new TPolyMarker3D(ndigits);
 108
 109       Int_t ibgnd   = 0;
 110       Int_t isignal = 0;
 111       Int_t imerged = 0;
 112
 113       for (Int_t time = 0; time < timeMax; time++) {
 114         for (Int_t  col = 0;  col <  colMax;  col++) {
 115           for (Int_t  row = 0;  row <  rowMax;  row++) {
 116
 117             Int_t type = 1;
 118
 119             Int_t amp = digits->GetDataUnchecked(row,col,time);
 120             for (idict = 0; idict < kNdict; idict++) {
 121               Int_t trk = tracks[idict]->GetDataUnchecked(row,col,time) - 1;
 122               if ((idict == 0) && (trk >= mask)) {
 123                 type = 2;
 124               }
 125               if ((type  == 1) && (trk >= mask)) {
 126                 type = 3;
 127               }
 128             }
 129
 130             if (amp > thresh) {
 131
 132               Float_t glb[3];
 133               Float_t loc[3];
 134
 135               loc[0] = row;
 136               loc[1] = col;
 137               loc[2] = time;
 138               geo->Local2Global(idet,loc,glb);
 139               Float_t x = glb[0];
 140               Float_t y = glb[1];
 141               Float_t z = glb[2];
 142
 143               if      (type == 1) {
 144                 pmSignal->SetPoint(isignal,x,y,z);
 145                 isignal++;
 146                 totalsignal++;
 147               }
 148               else if (type == 2) {
 149                 pmBgnd->SetPoint(ibgnd,x,y,z);
 150                 ibgnd++;
 151                 totalbgnd++;
 152               }
 153               else if (type == 3) {
 154                 pmMerged->SetPoint(imerged,x,y,z);
 155                 imerged++;
 156                 totalmerged++;
 157               }
 158
 159             }
 160
 161           }
 162         }
 163       }
 164
 165       digits->Compress(1,0);
 166       for (idict = 0; idict < kNdict; idict++) {
 167         tracks[idict]->Compress(1,0);
 168       }
 169
 170       pmMerged->SetMarkerSize(1);
 171       pmMerged->SetMarkerColor(markerColorMerged);
 172       pmMerged->SetMarkerStyle(1);
 173       pmMerged->Draw();
 174
 175       pmBgnd->SetMarkerSize(1);
 176       pmBgnd->SetMarkerColor(markerColorBgnd);
 177       pmBgnd->SetMarkerStyle(1);
 178       pmBgnd->Draw();
 179
 180       pmSignal->SetMarkerSize(1);
 181       pmSignal->SetMarkerColor(markerColorSignal);
 182       pmSignal->SetMarkerStyle(1);
 183       pmSignal->Draw();
 184
 185     }
 186
 187   }
 188
 189   TGeometry *geoAlice = (TGeometry *) gafl->Get("AliceGeom");
 190   TNode     *main     = (TNode *) ((geoAlice->GetListOfNodes())->First());
 191   TIter      next(main->GetListOfNodes());
 192   TNode     *module   = 0;
 193   while ((module = (TNode *) next())) {
 194     Char_t ch[100];
 195     sprintf(ch,"%s\n",module->GetTitle());
 196     if ((ch[0] == 'T') && ((ch[1] == 'R') || (ch[1] == 'P'))) {
 197       module->SetVisibility( 3);
 198     }
 199     else {
 200       module->SetVisibility(-1);
 201     }
 202   }
 203
 204   geoAlice->Draw("same");
 205
 206   c1->Modified();
 207   c1->Update();
 208
 209   gafl->Close();
 210
 211   printf("<AliTRDdisplayDigits3D> Number of digits:\n");
 212   printf("                        signal = %d, bgnd = %d, merged = %d\n"
 213         ,totalsignal,totalbgnd,totalmerged);
 214
 215   return 0;
 216
 217 }
 218