FMD/AliFMDPedestalDA.cxx

   1 /**************************************************************************
   2  * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /** @file    AliFMDPedestalDA.cxx
  17     @author  Hans Hjersing Dalsgaard <canute@nbi.dk>
  18     @date    Mon Mar 10 09:46:05 2008
  19     @brief   Derived class for the pedestal detector algorithm.
  20 */
  21 //
  22 // This class implements the virtual functions of the AliFMDBaseDA
  23 // class.  The most important of these functions, FillChannels(..) and
  24 // Analyse(...) collect and analyse the data of each channel. The
  25 // resulting pedestal and noise values are written to a comma
  26 // separated values (csv) file on the go. The csv files produced in
  27 // this way are the basic input to the AliFMDPreprocessor.
  28 //
  29
  30 #include "AliFMDPedestalDA.h"
  31 #include "AliFMDAltroMapping.h"
  32 #include <iostream>
  33 #include <fstream>
  34 #include "AliLog.h"
  35 #include "TF1.h"
  36 #include "TObject.h"
  37 #include "TMath.h"
  38 #include <TSystem.h>
  39 #include <TDatime.h>
  40 #include <TH2.h>
  41
  42 //_____________________________________________________________________
  43 ClassImp(AliFMDPedestalDA)
  44
  45 //_____________________________________________________________________
  46 AliFMDPedestalDA::AliFMDPedestalDA()
  47   : AliFMDBaseDA(),
  48     fCurrentChannel(1),
  49     fPedSummary("PedestalSummary","pedestals",51200,0,51200),
  50     fNoiseSummary("NoiseSummary","noise",51200,0,51200),
  51     fZSfileFMD1(),
  52     fZSfileFMD2(),
  53     fZSfileFMD3(),
  54     fMinTimebin(3 * 4 * 3 * 16), // 3 ddls, 4 FECs, 3 Altros, 16 channels
  55     fMaxTimebin(3 * 4 * 3 * 16), // 3 ddls, 4 FECs, 3 Altros, 16 channels
  56     fSummaryFMD1i(0),
  57     fSummaryFMD2i(0),
  58     fSummaryFMD2o(0),
  59     fSummaryFMD3i(0),
  60     fSummaryFMD3o(0)
  61 {
  62   // Default constructor
  63   Rotate("peds.csv", 3);
  64   fOutputFile.open("peds.csv");
  65   Rotate("ddl3072.csv", 10);
  66   fZSfileFMD1.open("ddl3072.csv");
  67   Rotate("ddl3073.csv", 10);
  68   fZSfileFMD2.open("ddl3073.csv");
  69   Rotate("ddl3074.csv", 10);
  70   fZSfileFMD3.open("ddl3074.csv");
  71 }
  72
  73 //_____________________________________________________________________
  74 AliFMDPedestalDA::AliFMDPedestalDA(const AliFMDPedestalDA & pedDA) :
  75   AliFMDBaseDA(pedDA),
  76   fCurrentChannel(1),
  77   fPedSummary("PedestalSummary","pedestals",51200,0,51200),
  78   fNoiseSummary("NoiseSummary","noise",51200,0,51200),
  79   fZSfileFMD1(),
  80   fZSfileFMD2(),
  81   fZSfileFMD3(),
  82   fMinTimebin(pedDA.fMinTimebin),
  83   fMaxTimebin(pedDA.fMaxTimebin),
  84   fSummaryFMD1i(pedDA.fSummaryFMD1i),
  85   fSummaryFMD2i(pedDA.fSummaryFMD2i),
  86   fSummaryFMD2o(pedDA.fSummaryFMD2o),
  87   fSummaryFMD3i(pedDA.fSummaryFMD3i),
  88   fSummaryFMD3o(pedDA.fSummaryFMD3o)
  89 {
  90   // Copy constructor
  91 }
  92
  93 //_____________________________________________________________________
  94 AliFMDPedestalDA::~AliFMDPedestalDA()
  95 {
  96   // Destructor.
  97 }
  98
  99 //_____________________________________________________________________
 100 void AliFMDPedestalDA::Init()
 101 {
 102   // Initialise
 103   SetRequiredEvents(1000);
 104   fMinTimebin.Reset(1024);
 105   fMaxTimebin.Reset(-1);
 106
 107
 108 }
 109
 110 //_____________________________________________________________________
 111 void AliFMDPedestalDA::AddChannelContainer(TObjArray* sectorArray,
 112                                            UShort_t det,
 113                                            Char_t   ring,
 114                                            UShort_t sec,
 115                                            UShort_t strip)
 116 {
 117   // Add a channel to the containers.
 118   //
 119   // Parameters:
 120   //     sectorArray  Array of sectors
 121   //     det          Detector
 122   //     ring         Ring
 123   //     sec          Sector
 124   //     strip        Strip
 125   AliFMDParameters* pars        = AliFMDParameters::Instance();
 126   UInt_t            samples     = pars->GetSampleRate(det, ring, sec, strip);
 127   TObjArray*        sampleArray = new TObjArray(samples);
 128   sampleArray->SetOwner();
 129   for (UInt_t sample = 0; sample < samples; sample++) {
 130     TH1S* hSample = new TH1S(Form("FMD%d%c[%02d,03%d]_%d",
 131                                   det,ring,sec,strip,sample),
 132                              Form("FMD%d%c[%02d,%03%d]_%d",
 133                                   det,ring,sec,strip),
 134                              1024,-.5,1023.5);
 135     hSample->SetXTitle("ADC");
 136     hSample->SetYTitle("Events");
 137     hSample->SetDirectory(0);
 138     sampleArray->AddAt(hSample, sample);
 139   }
 140   sectorArray->AddAtAndExpand(sampleArray, strip);
 141 }
 142
 143 //_____________________________________________________________________
 144 void AliFMDPedestalDA::FillChannels(AliFMDDigit* digit)
 145 {
 146   // Fill ADC values from a digit into the corresponding histogram.
 147   //
 148   // Parameters:
 149   //    digit Digit to fill ADC values for.
 150   UShort_t det   = digit->Detector();
 151   Char_t   ring  = digit->Ring();
 152   UShort_t sec   = digit->Sector();
 153   UShort_t strip = digit->Strip();
 154
 155   AliFMDParameters* pars     = AliFMDParameters::Instance();
 156   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
 157   for (UInt_t sample = 0; sample < samples; sample++) {
 158     TH1S* hSample = GetChannel(det, ring, sec, strip, sample);
 159     hSample->Fill(digit->Count(sample));
 160   }
 161
 162 }
 163
 164 //_____________________________________________________________________
 165 void AliFMDPedestalDA::MakeSummary(UShort_t det, Char_t ring)
 166 {
 167   std::cout << "Making summary for FMD" << det << ring << " ..." << std::endl;
 168   switch (det) {
 169   case 1:
 170     fSummaryFMD1i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
 171     break;
 172   case 2:
 173     switch (ring) {
 174     case 'I': case 'i':
 175       fSummaryFMD2i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
 176       break;
 177     case 'O': case 'o':
 178       fSummaryFMD2o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
 179       break;
 180     }
 181     break;
 182   case 3:
 183     switch (ring) {
 184     case 'I': case 'i':
 185       fSummaryFMD3i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
 186       break;
 187     case 'O': case 'o':
 188       fSummaryFMD3o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
 189       break;
 190     }
 191     break;
 192   }
 193 }
 194
 195 //_____________________________________________________________________
 196 void AliFMDPedestalDA::Analyse(UShort_t det,
 197                                Char_t   ring,
 198                                UShort_t sec,
 199                                UShort_t strip)
 200 {
 201   // Analyse a strip.  That is, compute the mean and spread of the ADC
 202   // spectra for all strips.  Also output on files the values.
 203   //
 204   // Parameters:
 205   //     det   Detector
 206   //     ring  Ring
 207   //     sec   Sector
 208   //     strip Strip.
 209   AliFMDParameters* pars     = AliFMDParameters::Instance();
 210   TH2* summary = 0;
 211   switch (det) {
 212   case 1: summary = fSummaryFMD1i; break;
 213   case 2:
 214     switch (ring) {
 215     case 'I':  summary = fSummaryFMD2i; break;
 216     case 'O':  summary = fSummaryFMD2o; break;
 217     }
 218     break;
 219   case 3:
 220     switch (ring) {
 221     case 'I':  summary = fSummaryFMD3i; break;
 222     case 'O':  summary = fSummaryFMD3o; break;
 223     }
 224     break;
 225   }
 226   static bool first = true;
 227   if (summary && first) {
 228     std::cout << "Filling summary " << summary->GetName() << std::endl;
 229     first = false;
 230   }
 231
 232   // Float_t           factor   = pars->GetPedestalFactor();
 233   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
 234   for (UShort_t sample = 0; sample < samples; sample++) {
 235
 236     TH1S* hChannel = GetChannel(det, ring, sec, strip,sample);
 237     if(hChannel->GetEntries() == 0) {
 238       //AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
 239       //            det,ring,sec,strip));
 240       return;
 241     }
 242
 243     AliDebug(50, Form("Fitting FMD%d%c_%d_%d with %d entries",
 244                       det,ring,sec,strip, hChannel->GetEntries()));
 245     TF1 fitFunc("fitFunc","gausn",0,300);
 246     fitFunc.SetParameters(100,100,1);
 247     hChannel->Fit("fitFunc","Q0+","",10,200);
 248
 249     Float_t mean = hChannel->GetMean();
 250     Float_t rms  = hChannel->GetRMS();
 251
 252
 253
 254     hChannel->GetXaxis()->SetRangeUser(mean-5*rms,mean+5*rms);
 255
 256     mean = hChannel->GetMean();
 257     rms  = hChannel->GetRMS();
 258
 259
 260     UShort_t ddl, board, altro, channel;
 261     UShort_t timebin;
 262
 263     pars->Detector2Hardware(det,ring,sec,strip,sample,
 264                             ddl,board,altro,channel,timebin);
 265     Int_t idx = HWIndex(ddl, board, altro, channel);
 266     if (idx >= 0) {
 267       fMinTimebin[idx] = TMath::Min(Short_t(timebin),   fMinTimebin[idx]);
 268       fMaxTimebin[idx] = TMath::Max(Short_t(timebin+1), fMaxTimebin[idx]);
 269     }
 270
 271     std::ostream* zsFile = 0;
 272     switch(det) {
 273     case 1:  zsFile = &fZSfileFMD1; break;
 274     case 2:  zsFile = &fZSfileFMD2; break;
 275     case 3:  zsFile = &fZSfileFMD3; break;
 276     default: AliWarning("Unknown sample!"); break;
 277
 278     }
 279     *zsFile << board   << ','
 280             << altro   << ','
 281             << channel << ','
 282             << timebin << ','
 283             << mean    << ','
 284             << rms     << "\n";
 285
 286     Float_t chi2ndf = 0;
 287
 288
 289     if(fitFunc.GetNDF())
 290       chi2ndf = fitFunc.GetChisquare() / fitFunc.GetNDF();
 291
 292
 293     Int_t sampleToWrite = 2;
 294     if      (samples == 2) sampleToWrite = 1;
 295     else if (samples <  2) sampleToWrite = 0;
 296
 297     if(sample != sampleToWrite) continue;
 298
 299
 300     fOutputFile << det                         << ','
 301                 << ring                        << ','
 302                 << sec                         << ','
 303                 << strip                       << ','
 304                 << mean                        << ','
 305                 << rms                         << ','
 306                 << fitFunc.GetParameter(1)     << ','
 307                 << fitFunc.GetParameter(2)     << ','
 308                 << chi2ndf                     <<"\n";
 309
 310     if (summary) {
 311       Int_t bin = summary->FindBin(sec, strip);
 312       summary->SetBinContent(bin, mean);
 313       summary->SetBinError(bin, rms);
 314     }
 315
 316     if(fSaveHistograms  ) {
 317       gDirectory->cd(GetSectorPath(det, ring, sec, kTRUE));
 318       TH1F* sumPed   = dynamic_cast<TH1F*>(gDirectory->Get("Pedestals"));
 319       TH1F* sumNoise = dynamic_cast<TH1F*>(gDirectory->Get("Noise"));
 320       Int_t nStr = (ring == 'I' ? 512 : 256);
 321       if (!sumPed) {
 322         sumPed = new TH1F("Pedestals",
 323                           Form("Summary of pedestals in FMD%d%c[%02d]",
 324                                det, ring, sec),
 325                           nStr, -.5, nStr-.5);
 326         sumPed->SetXTitle("Strip");
 327         sumPed->SetYTitle("Pedestal [ADC]");
 328         sumPed->SetDirectory(gDirectory);
 329       }
 330       if (!sumNoise) {
 331         sumNoise = new TH1F("Noise",
 332                             Form("Summary of noise in FMD%d%c[%02d]",
 333                                  det, ring, sec),
 334                             nStr, -.5, nStr-.5);
 335         sumNoise->SetXTitle("Strip");
 336         sumNoise->SetYTitle("Noise [ADC]");
 337
 338         sumNoise->SetDirectory(gDirectory);
 339       }
 340       sumPed->SetBinContent(strip+1, mean);
 341       sumPed->SetBinError(strip+1, rms);
 342       sumNoise->SetBinContent(strip+1, rms);
 343
 344       if(sumNoise->GetEntries() == nStr)
 345         sumNoise->Write(sumNoise->GetName(),TObject::kOverwrite);
 346       if(sumPed->GetEntries() == nStr)
 347         sumPed->Write(sumPed->GetName(),TObject::kOverwrite);
 348
 349       fPedSummary.SetBinContent(fCurrentChannel,mean);
 350
 351       fNoiseSummary.SetBinContent(fCurrentChannel,rms);
 352       fCurrentChannel++;
 353
 354       gDirectory->cd(GetStripPath(det, ring, sec, strip, kTRUE));
 355       hChannel->GetXaxis()->SetRange(1,1024);
 356
 357       hChannel->Write();
 358     }
 359   }
 360 }
 361
 362 //_____________________________________________________________________
 363 void AliFMDPedestalDA::Terminate(TFile* diagFile)
 364 {
 365   // Called at the end of a job.  Fills in missing time-bins and
 366   // closes output files
 367   if(fSaveHistograms) {
 368     diagFile->cd();
 369
 370     fPedSummary.Write();
 371     fNoiseSummary.Write();
 372   }
 373   AliFMDAltroMapping* map = AliFMDParameters::Instance()->GetAltroMap();
 374   for (Int_t i = 0; i < 3; i++) {
 375     std::ofstream& out = (i == 0 ? fZSfileFMD1 :
 376                           i == 1 ? fZSfileFMD2 :
 377                           fZSfileFMD3);
 378     if (out.is_open() && fSeenDetectors[i]) {
 379       FillinTimebins(out, map->Detector2DDL(i+1));
 380     }
 381     if (!fSeenDetectors[i]) {
 382       TString n(Form("ddl%d.csv",3072+map->Detector2DDL(i+1)));
 383       gSystem->Unlink(n.Data());
 384     }
 385   }
 386
 387 }
 388
 389 //_____________________________________________________________________
 390 void AliFMDPedestalDA::FillinTimebins(std::ofstream& out, UShort_t /*ddl*/)
 391 {
 392 #if 0
 393   unsigned short  boards[] = { 0x0, 0x1, 0x10, 0x11, 0xFFFF };
 394   unsigned short* board    = boards;
 395   while ((*boards) != 0xFFFF) {
 396     for (UShort_t altro = 0; altro < 3; altro++) {
 397       for (UShort_t channel = 0; channel < 16; channel++) {
 398         Int_t idx = HWIndex(ddl, *board, altro, channel);
 399         if (idx < 0) {
 400           AliWarning(Form("Invalid index for %4d/0x%02x/0x%x/0x%x: %d",
 401                           ddl, *board, altro, channel, idx));
 402           continue;
 403         }
 404         Short_t min = fMinTimebin[idx];
 405         Short_t max = fMaxTimebin[idx];
 406
 407         // Channel not seen at all.
 408         if (min > 1023 || max < 0) continue;
 409
 410         out << "# Extra timebins for 0x" << std::hex
 411             << board << ',' << altro << ',' << channel
 412             << " got time-bins " << min << " to " << max-1
 413             << std::dec << std::endl;
 414
 415         for (UShort_t t = 15; t < min; t++)
 416           // Write a phony line
 417           out << board << "," << altro << "," << channel << ","
 418               << t << "," << 1023 << "," << 0 << std::endl;
 419
 420         for (UShort_t t = max; t < 1024; t++)
 421           // Write a phony line
 422           out << board << "," << altro << "," << channel << ","
 423               << t << "," << 1023 << "," << 0 << std::endl;
 424       } // channel loop
 425     } // altro loop
 426   } // board loop
 427   // Write trailer, and close
 428 #endif
 429   out.write("# EOF\n", 6);
 430   out.close();
 431 }
 432
 433 //_____________________________________________________________________
 434 void AliFMDPedestalDA::WriteHeaderToFile()
 435 {
 436   // Write headers to output files
 437   AliFMDParameters* pars       = AliFMDParameters::Instance();
 438   fOutputFile.write(Form("# %s \n",pars->GetPedestalShuttleID()),13);
 439   TDatime now;
 440   fOutputFile << "# This file created from run # " << fRunno
 441               << " @ " << now.AsString() << std::endl;
 442   fOutputFile.write("# Detector, "
 443                     "Ring, "
 444                     "Sector, "
 445                     "Strip, "
 446                     "Pedestal, "
 447                     "Noise, "
 448                     "Mu, "
 449                     "Sigma, "
 450                     "Chi2/NDF \n", 71);
 451
 452   std::ostream* zss[] = { &fZSfileFMD1, &fZSfileFMD2, &fZSfileFMD3, 0 };
 453   for (size_t i = 0; i < 3; i++)  {
 454     *(zss[i]) << "# FMD " << (i+1) << " pedestals \n"
 455               << "# board, "
 456               << "altro, "
 457               << "channel, "
 458               << "timebin, "
 459               << "pedestal, "
 460               << "noise\n";
 461     *(zss[i]) << "# This file created from run # " << fRunno
 462               << " @ " << now.AsString() << std::endl;
 463   }
 464 }
 465
 466 //_____________________________________________________________________
 467 TH1S* AliFMDPedestalDA::GetChannel(UShort_t det,
 468                                    Char_t   ring,
 469                                    UShort_t sec,
 470                                    UShort_t strip,
 471                                    UInt_t   sample)
 472 {
 473   // Get the histogram corresponding to a strip sample.
 474   //
 475   // Parameters:
 476   //     det    Detector
 477   //     ring   Ring
 478   //     sec    Sector
 479   //     strip  Strip
 480   //     sample Sample
 481   //
 482   // Return:
 483   //     ADC spectra of a strip.
 484   UShort_t   iring       = (ring == 'O' ? 0 : 1);
 485   TObjArray* detArray    = static_cast<TObjArray*>(fDetectorArray.At(det));
 486   TObjArray* ringArray   = static_cast<TObjArray*>(detArray->At(iring));
 487   TObjArray* secArray    = static_cast<TObjArray*>(ringArray->At(sec));
 488   TObjArray* sampleArray = static_cast<TObjArray*>(secArray->At(strip));
 489   TH1S*      hSample     = static_cast<TH1S*>(sampleArray->At(sample));
 490   return hSample;
 491
 492 }
 493
 494 //_____________________________________________________________________
 495 //
 496 //EOF
 497 //