PMD/AliPMDCalibGain.cxx

   1 //////////////////////////////////////////////////////////////////////////////
   2 //
   3 // * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4 // *                                                                        *
   5 // * Author: The ALICE Off-line Project.                                    *
   6 // * Contributors are mentioned in the code where appropriate.              *
   7 // *                                                                        *
   8 // * Permission to use, copy, modify and distribute this software and its   *
   9 // * documentation strictly for non-commercial purposes is hereby granted   *
  10 // * without fee, provided that the above copyright notice appears in all   *
  11 // * copies and that both the copyright notice and this permission notice   *
  12 // * appear in the supporting documentation. The authors make no claims     *
  13 // * about the suitability of this software for any purpose. It is          *
  14 // * provided "as is" without express or implied warranty.                  *
  15 // **************************************************************************/
  16 //
  17 //////////////////////////////////////////////////////////////////////////////
  18
  19 #include "TF1.h"
  20 #include "TFile.h"
  21 #include "TObjString.h"
  22 #include "TROOT.h"
  23 #include "TClonesArray.h"
  24 #include "TH1F.h"
  25 #include "TObjArray.h"
  26 #include "TTree.h"
  27
  28 // --- Standard library ---
  29
  30 // --- AliRoot header files ---
  31 #include "AliDAQ.h"
  32 #include "AliLog.h"
  33 #include "AliRawReader.h"
  34 #include "AliPMDRawStream.h"
  35 #include "AliPMDddldata.h"
  36 #include "AliPMDCalibGain.h"
  37
  38 ClassImp(AliPMDCalibGain)
  39
  40 AliPMDCalibGain::AliPMDCalibGain():
  41   TObject(),
  42   fpw(NULL)
  43 {
  44   // Standard Constructor
  45     for(Int_t idet = 0; idet < kDet; idet++)
  46     {
  47      fDetCount[kDet] =0.;
  48         for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
  49         {
  50             fSMIso[idet][ismn]   = 0.;
  51             fSMCount[idet][ismn] = 0.;
  52             for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
  53             {
  54                 for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
  55                 {
  56                     fCellIso[idet][ismn][jrow][kcol]   = 0.;
  57                     fCellCount[idet][ismn][jrow][kcol] = 0.;
  58                     fPedMeanRMS[idet][ismn][jrow][kcol] = 0.;
  59
  60                 }
  61             }
  62         }
  63     }
  64
  65
  66
  67 }
  68 // ------------------------------------------------------------------------ //
  69 AliPMDCalibGain::AliPMDCalibGain(const AliPMDCalibGain &pmdcalibgain):
  70   TObject(pmdcalibgain),
  71   fpw(NULL)
  72 {
  73     for(Int_t idet = 0; idet < kDet; idet++)
  74     {
  75      fDetCount[idet] = pmdcalibgain.fDetCount[idet];
  76      fDetIso[idet] = pmdcalibgain.fDetIso[idet];
  77         for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
  78         {
  79             fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn] ;
  80             fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn] ;
  81             for(Int_t jrow = 0; jrow < kMaxRow; jrow++)
  82             {
  83                 for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
  84                 {
  85                   fCellIso[idet][ismn][jrow][kcol]    =
  86                     pmdcalibgain.fCellIso[idet][ismn][jrow][kcol];
  87                   fCellCount[idet][ismn][jrow][kcol]  =
  88                     pmdcalibgain.fCellCount[idet][ismn][jrow][kcol];
  89                   fPedMeanRMS[idet][ismn][jrow][kcol] =
  90                     pmdcalibgain.fPedMeanRMS[idet][ismn][jrow][kcol];
  91
  92                 }
  93             }
  94         }
  95     }
  96
  97 }
  98 // ------------------------------------------------------------------------ //
  99 AliPMDCalibGain &AliPMDCalibGain::operator=(const AliPMDCalibGain &pmdcalibgain)
 100 {
 101     if(this != &pmdcalibgain)
 102     {
 103       this->fpw = pmdcalibgain.fpw;
 104         for(Int_t idet = 0; idet < kDet; idet++)
 105         {
 106             for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
 107             {
 108                 fSMIso[idet][ismn] = pmdcalibgain.fSMIso[idet][ismn];
 109                 fSMCount[idet][ismn] = pmdcalibgain.fSMCount[idet][ismn];
 110                 for(Int_t jrow = 0; jrow < kMaxRow;jrow++)
 111                 {
 112                     for(Int_t kcol = 0; kcol < kMaxCol; kcol++)
 113                       {
 114                         fCellIso[idet][ismn][jrow][kcol]  =
 115                           pmdcalibgain.fCellIso[idet][ismn][jrow][kcol];
 116                         fCellCount[idet][ismn][jrow][kcol]  =
 117                           pmdcalibgain.fCellCount[idet][ismn][jrow][kcol];
 118                         fPedMeanRMS[idet][ismn][jrow][kcol] =
 119                           pmdcalibgain.fPedMeanRMS[idet][ismn][jrow][kcol];
 120
 121                     }
 122                 }
 123             }
 124         }
 125     }
 126     return *this;
 127 }
 128 // ------------------------------------------------------------------------ //
 129 AliPMDCalibGain::~AliPMDCalibGain()
 130 {
 131     // dtor
 132
 133 }
 134
 135 // ------------------------------------------------------------------------ //
 136
 137 Int_t AliPMDCalibGain::ExtractPedestal()
 138 {
 139   // Pedestal extraction from the PMD_PED.root file
 140   // To be called once at the beginning
 141
 142   Int_t   det, sm, row, col;
 143   Float_t mean, rms;
 144
 145   TFile *pedfile = new TFile("PMD_PED1.root");
 146
 147   if(!pedfile)
 148     {
 149       printf("ERROR --- NO PEDESTAL (PMD_PED1.root) FILE IS FOUND --- STOP GAIN DA\n");
 150       return -3;
 151     }
 152
 153
 154   TTree *ped =(TTree*)pedfile->Get("ped");
 155
 156   ped->SetBranchAddress("det",&det);
 157   ped->SetBranchAddress("sm",&sm);
 158   ped->SetBranchAddress("row",&row);
 159   ped->SetBranchAddress("col",&col);
 160   ped->SetBranchAddress("mean",&mean);
 161   ped->SetBranchAddress("rms",&rms);
 162
 163   Int_t nentries = (Int_t)ped->GetEntries();
 164
 165   for (Int_t ient = 0; ient < nentries; ient++)
 166     {
 167       ped->GetEntry(ient);
 168       fPedMeanRMS[det][sm][row][col] = mean + 3.*rms;
 169       //printf("Mean= %f, RMS= %f, PedMeanRMS=%f\n",mean,rms,fPedMeanRMS[det][sm][row][col]);
 170
 171     }
 172
 173   pedfile->Close();
 174   delete pedfile;
 175   pedfile = 0x0;
 176
 177   return 1;
 178 }
 179 // ------------------------------------------------------------------------ //
 180
 181 void AliPMDCalibGain::ReadIntermediateFile()
 182 {
 183   // Read the variables from the file
 184
 185   fpw = fopen("interfile.dat","r");
 186
 187   Float_t detcount, detiso;
 188   Float_t smcount, smiso;
 189   Float_t cellcount, celliso;
 190
 191   for (Int_t idet = 0; idet < kDet; idet++)
 192     {
 193       fscanf(fpw,"%d %f %f",&idet,&detcount,&detiso);
 194       fDetCount[idet] = detcount;
 195       fDetIso[idet]   = detiso;
 196     }
 197
 198   for (Int_t idet = 0; idet < kDet; idet++)
 199     {
 200       for (Int_t ism = 0; ism < kMaxSMN; ism++)
 201         {
 202           fscanf(fpw,"%d %d %f %f",&idet,&ism,&smcount,&smiso);
 203
 204           fSMCount[idet][ism] = smcount;
 205           fSMIso[idet][ism]   = smiso;
 206         }
 207     }
 208
 209   for (Int_t idet = 0; idet < kDet; idet++)
 210     {
 211       for (Int_t ism = 0; ism < kMaxSMN; ism++)
 212         {
 213           for (Int_t irow = 0; irow < kMaxRow; irow++)
 214             {
 215               for (Int_t icol = 0; icol < kMaxCol; icol++)
 216                 {
 217                   fscanf(fpw,"%d %d %d %d %f %f",&idet,&ism,&irow,&icol,
 218                           &cellcount,&celliso);
 219
 220                   fCellCount[idet][ism][irow][icol] = cellcount;
 221                   fCellIso[idet][ism][irow][icol]   = celliso;
 222                 }
 223             }
 224         }
 225     }
 226
 227   fclose(fpw);
 228
 229 }
 230 // ------------------------------------------------------------------------ //
 231 void AliPMDCalibGain::WriteIntermediateFile()
 232 {
 233
 234  //Following variables to be written
 235
 236   /*
 237     fDetIso[idet] ;
 238     fSMIso[idet][ismn];
 239     fCellIso[idet][ismn][irow][icol];
 240
 241     fDetCount[idet];
 242     fSMCount[idet][ismn];
 243     fCellCount[idet][ismn][irow][icol];
 244   */
 245
 246   fpw = fopen("interfile.dat","w+");
 247
 248   for (Int_t idet = 0; idet < kDet; idet++)
 249     {
 250       fprintf(fpw,"%d %f %f\n",idet,fDetCount[idet],fDetIso[idet]);
 251     }
 252
 253   for (Int_t idet = 0; idet < kDet; idet++)
 254     {
 255       for (Int_t ism = 0; ism < kMaxSMN; ism++)
 256         {
 257           fprintf(fpw,"%d %d %f %f\n",idet,ism, fSMCount[idet][ism],fSMIso[idet][ism]);
 258         }
 259     }
 260
 261   for (Int_t idet = 0; idet < kDet; idet++)
 262     {
 263       for (Int_t ism = 0; ism < kMaxSMN; ism++)
 264         {
 265           for (Int_t irow = 0; irow < kMaxRow; irow++)
 266             {
 267               for (Int_t icol = 0; icol < kMaxCol; icol++)
 268                 {
 269                   fprintf(fpw,"%d %d %d %d %f %f\n",idet,ism,irow,icol,
 270                           fCellCount[idet][ism][irow][icol],
 271                           fCellIso[idet][ism][irow][icol]);
 272                 }
 273             }
 274         }
 275     }
 276
 277   fclose(fpw);
 278
 279 }
 280
 281 // ------------------------------------------------------------------------ //
 282
 283 Bool_t AliPMDCalibGain::ProcessEvent(AliRawReader *rawReader, TObjArray *pmdddlcont)
 284 {
 285   // Calculates the ADC of isolated cell
 286
 287   const Int_t kDDL           = AliDAQ::NumberOfDdls("PMD");
 288   const Int_t kCellNeighbour = 6;
 289   Int_t neibx[6] = {1,0,-1,-1,0,1};
 290   Int_t neiby[6] = {0,1,1,0,-1,-1};
 291
 292   Int_t id1,jd1;  //neighbour row/col
 293   Int_t isocount; //number of neighbours with 0 signal
 294
 295   Float_t d1[kDet][kMaxSMN][kMaxRow][kMaxCol];
 296
 297   for(Int_t idet = 0; idet < kDet; idet++)
 298     {
 299       for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
 300         {
 301           for(Int_t irow = 0; irow < kMaxRow; irow++)
 302             {
 303               for(Int_t icol = 0; icol < kMaxCol; icol++)
 304                 {
 305                   d1[idet][ismn][irow][icol] = 0.;
 306                 }
 307             }
 308         }
 309     }
 310
 311   AliPMDRawStream rawStream(rawReader);
 312
 313   Int_t iddl = -1;
 314
 315   Int_t numberofDDLs = 0;
 316
 317     while ((iddl = rawStream.DdlData(pmdddlcont)) >=0) {
 318       numberofDDLs++;
 319
 320       Int_t ientries = pmdddlcont->GetEntries();
 321
 322       for (Int_t ient = 0; ient < ientries; ient++)
 323       {
 324           AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont->UncheckedAt(ient);
 325
 326           Int_t idet = pmdddl->GetDetector();
 327           Int_t ismn = pmdddl->GetSMN();
 328           //Int_t mcm = pmdddl->GetMCM();
 329           //Int_t ichno = pmdddl->GetChannel();
 330           Int_t irow = pmdddl->GetRow();
 331           Int_t icol = pmdddl->GetColumn();
 332           Int_t isig = pmdddl->GetSignal();
 333
 334
 335           // Pedestal subtraction
 336
 337           if (isig>0)
 338             {
 339               d1[idet][ismn][irow][icol] =
 340                 (Float_t) isig - fPedMeanRMS[idet][ismn][irow][icol];
 341 //printf("Signal_ped_subtracted=%f, pedestal=%f\n",d1[idet][ismn][irow][icol]),fPedMeanRMS[idet][ismn][irow][icol];
 342             }
 343       }
 344       pmdddlcont->Delete();
 345   }
 346
 347   for(Int_t idet=0; idet < kDet; idet++)
 348   {
 349       for(Int_t ismn = 0; ismn < kMaxSMN; ismn++)
 350       {
 351           for(Int_t irow = 0; irow < kMaxRow; irow++)
 352           {
 353               for(Int_t icol = 0; icol < kMaxCol; icol++)
 354               {
 355                   if(d1[idet][ismn][irow][icol] > 0)
 356                   {
 357                       isocount = 0;
 358                       for(Int_t ii = 0; ii < kCellNeighbour; ii++)
 359                       {
 360                           id1 = irow + neibx[ii];
 361                           jd1 = icol + neiby[ii];
 362                           if(d1[idet][ismn][id1][jd1] == 0)
 363                           {
 364                               isocount++;
 365                               if(isocount == kCellNeighbour)
 366                               {
 367                                   fDetIso[idet] += d1[idet][ismn][irow][icol];
 368                                   fSMIso[idet][ismn] += d1[idet][ismn][irow][icol];
 369                                   fCellIso[idet][ismn][irow][icol] += d1[idet][ismn][irow][icol];
 370                                   fDetCount[idet]++;
 371                                   fSMCount[idet][ismn]++;
 372                                   fCellCount[idet][ismn][irow][icol]++;
 373
 374                               }
 375                           }
 376                       }  // neigh cell cond.
 377                   }     // d>0 cond.
 378               }
 379           }
 380       }
 381   }
 382
 383   if (numberofDDLs < kDDL)
 384       return kFALSE;
 385   return kTRUE;
 386
 387 }
 388 // ------------------------------------------------------------------------ //
 389 void AliPMDCalibGain::Analyse(TTree *gaintree)
 390 {
 391     // Calculates the mean
 392     Int_t   det, sm, row, col;
 393     Float_t gain;
 394     Float_t modmean  = 0.;
 395     Float_t cellmean = 0.;
 396     Float_t detmean =0.;
 397
 398     gaintree->Branch("det",&det,"det/I");
 399     gaintree->Branch("sm",&sm,"sm/I");
 400     gaintree->Branch("row",&row,"row/I");
 401     gaintree->Branch("col",&col,"col/I");
 402     gaintree->Branch("gain",&gain,"gain/F");
 403
 404     for(Int_t idet = 0; idet < kDet; idet++)
 405       {
 406         if (fDetCount[idet]>0 )
 407           detmean=fDetIso[idet]/fDetCount[idet];
 408         for(Int_t ism = 0; ism < kMaxSMN; ism++)
 409           {
 410             if (fSMCount[idet][ism] > 0)
 411               modmean = fSMIso[idet][ism]/fSMCount[idet][ism];
 412             for(Int_t irow = 0; irow < kMaxRow; irow++)
 413               {
 414                 for(Int_t icol = 0; icol < kMaxCol; icol++)
 415                   {
 416
 417                     if (fCellCount[idet][ism][irow][icol] > 0.)
 418                       {
 419                         cellmean = fCellIso[idet][ism][irow][icol]/fCellCount[idet][ism][irow][icol];
 420                       }
 421                     det      = idet;
 422                     sm       = ism;
 423                     row      = irow;
 424                     col      = icol;
 425                     if (cellmean > 0.0 && fCellCount[idet][ism][irow][icol]>0.)
 426                       {
 427                         gain = cellmean/detmean;
 428                       }
 429                     else
 430                       {
 431                         gain = -1.;
 432                       }
 433                     //if(fCellCount[idet][ism][irow][icol]>0.) printf("CellCount =%f, gain= %f\n",fCellCount[idet][ism][irow][icol],gain);
 434                     gaintree->Fill();
 435                   }
 436               }
 437           }
 438       }
 439
 440 }
 441 // ------------------------------------------------------------------------ //