PMD/AliPMDRecpointRead.C

   1 // ----------------------------------------------------//
   2 //                                                     //
   3 //       This macro reads the PMD clusters which       //
   4 //       are stored in the file "PMD.RecPoints.root"   //
   5 //                                                     //
   6 // ----------------------------------------------------//
   7
   8 #include <Riostream.h>
   9 #include "TBranch.h"
  10 #include "TStopwatch.h"
  11
  12 extern AliRun *gAlice;
  13
  14 Int_t AliPMDRecpointRead(Int_t nevent = 1)
  15 {
  16   if (gAlice)
  17     {
  18       delete AliRunLoader::Instance();
  19       delete gAlice;//if everything was OK here it is already NULL
  20       gAlice = 0x0;
  21     }
  22   AliRunLoader *fRunLoader = AliRunLoader::Open("galice.root","Event","update");
  23   if (!fRunLoader)
  24     {
  25       cerr<<"Can't load RunLoader"<<endl;
  26       return 1;
  27     }
  28   AliLoader *pmdloader = fRunLoader->GetLoader("PMDLoader");
  29   Int_t nevent = fRunLoader->GetNumberOfEvents();
  30   cout << " * *********** nevent = " << nevent << endl;
  31
  32   if (pmdloader == 0x0)
  33     {
  34       cerr<<" ===> Can not find PMD or PMDLoader <===\n";
  35       delete fRunLoader;
  36       return 2;
  37     }
  38
  39   pmdloader->LoadRecPoints("READ");
  40
  41
  42   AliCDBManager *man = AliCDBManager::Instance();
  43   man->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
  44   man->SetRun(0);
  45
  46   TClonesArray *fRecpoints;
  47   AliPMDUtility *cc = new AliPMDUtility();
  48
  49   cc->ApplyAlignment();
  50
  51   TH2F *h2 = new TH2F("h2"," ",100,-100.,100.,100,-100.,100.);
  52
  53   FILE *fpw = fopen("junk_rec.dat","w");
  54
  55   for (Int_t ievt = 0; ievt < nevent; ievt++)
  56     {
  57       fRunLoader->GetEvent(ievt);
  58       TTree *treeR = pmdloader->TreeR();
  59       if (treeR == 0x0)
  60         {
  61           cout << " Can not get TreeR" << endl;
  62           return 3;
  63         }
  64
  65       AliPMDrecpoint1  *pmdrecpoint;
  66       TBranch *branch1 = treeR->GetBranch("PMDRecpoint");
  67       branch1->SetAddress(&fRecpoints);
  68       /**********************************************************************
  69        *    det   : Detector, 0: PRE & 1:CPV                                *
  70        *    smn   : Serial Module Number from 0 to 23 for both detector     *
  71        *    xpos  : x-position of the cluster                               *
  72        *    ypos  : y-position of the cluster                               *
  73        *            THESE xpos & ypos are not the true xpos and ypos        *
  74        *            for some of the unit modules. They are rotated.         *
  75        *    adc   : ADC contained in the cluster                            *
  76        *    ncell : Number of cells contained in the cluster                *
  77        *    rad   : radius of the cluster (1d fit)                          *
  78        *    xpad  : TRUE x-position of the cluster                          *
  79        *    ypad  : TRUE y-position of the cluster                          *
  80        **********************************************************************/
  81
  82       Int_t   det,smn;
  83       Float_t xpos,ypos, xpad, ypad;
  84       Float_t adc, ncell, sigx, sigy;
  85       Float_t xx, yy;
  86       Int_t   nmodules = branch1->GetEntries();
  87       cout << " nmodules = " << nmodules << endl;
  88       for (Int_t imodule = 0; imodule < nmodules; imodule++)
  89         {
  90           branch1->GetEntry(imodule);
  91           Int_t nentries = fRecpoints->GetLast();
  92           for(Int_t ient = 0; ient < nentries+1; ient++)
  93             {
  94               pmdrecpoint = (AliPMDrecpoint1*)fRecpoints->UncheckedAt(ient);
  95               det   = (Int_t) pmdrecpoint->GetDetector();
  96               smn   = (Int_t) pmdrecpoint->GetSMNumber();
  97               xpos  = pmdrecpoint->GetClusX();
  98               ypos  = pmdrecpoint->GetClusY();
  99               adc   = pmdrecpoint->GetClusADC();
 100               ncell = pmdrecpoint->GetClusCells();
 101               sigx  = pmdrecpoint->GetClusSigmaX();
 102               sigy  = pmdrecpoint->GetClusSigmaY();
 103
 104               //
 105               // User has to plug in his analysis code here
 106               //
 107
 108               fprintf(fpw,"%d %d %d %d\n",
 109                       det,smn,xpos,ypos);
 110               //
 111               // Plot the cluster centroid to see the PMD geometry
 112               // using the PMD Utility class
 113               //
 114               if (det == 1)
 115                 {
 116                   // Draw only for PRE plane
 117                   cc->RectGeomCellPos(smn,xpos,ypos,xx,yy);
 118                   h2->Fill(xx,yy);
 119                 }
 120
 121               //
 122               // End of the User code
 123               //
 124             }
 125         }
 126
 127     }
 128
 129   h2->Draw();
 130
 131   fclose(fpw);
 132   return 0;
 133 }
 134