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