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