AliRoot Core  v5-06-15 (45dab64)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTPCCalibPulser.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 
164 
165 //Root includes
166 #include <TObjArray.h>
167 #include <TH1F.h>
168 #include <TH2F.h>
169 #include <TH2S.h>
170 #include <TString.h>
171 #include <TVectorF.h>
172 #include <TMath.h>
173 #include <TMap.h>
174 #include <TDirectory.h>
175 #include <TSystem.h>
176 #include <TROOT.h>
177 #include <TFile.h>
178 
179 //AliRoot includes
180 #include "AliRawReader.h"
181 #include "AliRawReaderRoot.h"
182 #include "AliRawReaderDate.h"
183 #include "AliTPCCalROC.h"
184 #include "AliTPCCalPad.h"
185 #include "AliTPCROC.h"
186 #include "AliTPCParam.h"
187 #include "AliTPCCalibPulser.h"
188 #include "AliTPCcalibDB.h"
189 #include "AliMathBase.h"
190 #include "AliLog.h"
191 #include "TTreeStream.h"
192 
193 //date
194 #include "event.h"
195 
196 
197 
198 
202 
205 fNbinsT0(200),
206 fXminT0(-2),
207 fXmaxT0(2),
208 fNbinsQ(200),
209 fXminQ(10),
210 fXmaxQ(40),
211 fNbinsRMS(100),
212 fXminRMS(0.1),
213 fXmaxRMS(5.1),
214 fPeakIntMinus(2),
215 fPeakIntPlus(2),
216 fIsZeroSuppressed(kFALSE),
217 fLastSector(-1),
218 fParam(new AliTPCParam),
219 fPedestalTPC(0x0),
220 fPadNoiseTPC(0x0),
221 fOutliers(0x0),
222 fPedestalROC(0x0),
223 fPadNoiseROC(0x0),
224 fCalRocArrayT0(72),
225 fCalRocArrayQ(72),
226 fCalRocArrayRMS(72),
227 fCalRocArrayOutliers(72),
228 fHistoQArray(72),
229 fHistoT0Array(72),
230 fHistoRMSArray(72),
231 fHMeanTimeSector(0x0),
232 fVMeanTimeSector(72),
233 fPadTimesArrayEvent(72),
234 fPadQArrayEvent(72),
235 fPadRMSArrayEvent(72),
236 fPadPedestalArrayEvent(72),
237 fCurrentChannel(-1),
238 fCurrentSector(-1),
239 fCurrentRow(-1),
240 fCurrentPad(-1),
241 fMaxPadSignal(-1),
242 fMaxTimeBin(-1),
243 fPadSignal(1024),
244 fPadPedestal(0),
245 fPadNoise(0),
246 fVTime0Offset(72),
247 fVTime0OffsetCounter(72)
248 {
249  //
250  // AliTPCSignal default constructor
251  //
252  SetNameTitle("AliTPCCalibPulser","AliTPCCalibPulser");
253  fFirstTimeBin=60;
254  fLastTimeBin=900;
255  fParam->Update();
256 }
257 //_____________________________________________________________________
259 AliTPCCalibRawBase(sig),
260 fNbinsT0(sig.fNbinsT0),
261 fXminT0(sig.fXminT0),
262 fXmaxT0(sig.fXmaxT0),
263 fNbinsQ(sig.fNbinsQ),
264 fXminQ(sig.fXminQ),
265 fXmaxQ(sig.fXmaxQ),
266 fNbinsRMS(sig.fNbinsRMS),
267 fXminRMS(sig.fXminRMS),
268 fXmaxRMS(sig.fXmaxRMS),
269 fPeakIntMinus(sig.fPeakIntMinus),
270 fPeakIntPlus(sig.fPeakIntPlus),
271 fIsZeroSuppressed(sig.fIsZeroSuppressed),
272 fLastSector(-1),
273 fParam(new AliTPCParam),
274 fPedestalTPC(0x0),
275 fPadNoiseTPC(0x0),
276 fOutliers(0x0),
277 fPedestalROC(0x0),
278 fPadNoiseROC(0x0),
279 fCalRocArrayT0(72),
280 fCalRocArrayQ(72),
281 fCalRocArrayRMS(72),
282 fCalRocArrayOutliers(72),
283 fHistoQArray(72),
284 fHistoT0Array(72),
285 fHistoRMSArray(72),
286 fHMeanTimeSector(0x0),
287 fVMeanTimeSector(72),
288 fPadTimesArrayEvent(72),
289 fPadQArrayEvent(72),
290 fPadRMSArrayEvent(72),
291 fPadPedestalArrayEvent(72),
292 fCurrentChannel(-1),
293 fCurrentSector(-1),
294 fCurrentRow(-1),
295 fCurrentPad(-1),
296 fMaxPadSignal(-1),
297 fMaxTimeBin(-1),
298 fPadSignal(1024),
299 fPadPedestal(0),
300 fPadNoise(0),
301 fVTime0Offset(72),
302 fVTime0OffsetCounter(72)
303 {
305 
306  for (Int_t iSec = 0; iSec < 72; ++iSec){
307  const AliTPCCalROC *calQ = (AliTPCCalROC*)sig.fCalRocArrayQ.UncheckedAt(iSec);
308  const AliTPCCalROC *calT0 = (AliTPCCalROC*)sig.fCalRocArrayT0.UncheckedAt(iSec);
309  const AliTPCCalROC *calRMS = (AliTPCCalROC*)sig.fCalRocArrayRMS.UncheckedAt(iSec);
310  const AliTPCCalROC *calOut = (AliTPCCalROC*)sig.fCalRocArrayOutliers.UncheckedAt(iSec);
311 
312  const TH2S *hQ = (TH2S*)sig.fHistoQArray.UncheckedAt(iSec);
313  const TH2S *hT0 = (TH2S*)sig.fHistoT0Array.UncheckedAt(iSec);
314  const TH2S *hRMS = (TH2S*)sig.fHistoRMSArray.UncheckedAt(iSec);
315 
316  if ( calQ != 0x0 ) fCalRocArrayQ.AddAt(new AliTPCCalROC(*calQ), iSec);
317  if ( calT0 != 0x0 ) fCalRocArrayT0.AddAt(new AliTPCCalROC(*calT0), iSec);
318  if ( calRMS != 0x0 ) fCalRocArrayRMS.AddAt(new AliTPCCalROC(*calRMS), iSec);
319  if ( calOut != 0x0 ) fCalRocArrayOutliers.AddAt(new AliTPCCalROC(*calOut), iSec);
320 
321  if ( hQ != 0x0 ){
322  TH2S *hNew = new TH2S(*hQ);
323  hNew->SetDirectory(0);
324  fHistoQArray.AddAt(hNew,iSec);
325  }
326  if ( hT0 != 0x0 ){
327  TH2S *hNew = new TH2S(*hT0);
328  hNew->SetDirectory(0);
329  fHistoQArray.AddAt(hNew,iSec);
330  }
331  if ( hRMS != 0x0 ){
332  TH2S *hNew = new TH2S(*hRMS);
333  hNew->SetDirectory(0);
334  fHistoQArray.AddAt(hNew,iSec);
335  }
336  fVMeanTimeSector[iSec]=sig.fVMeanTimeSector[iSec];
337  }
338 
339  if ( sig.fHMeanTimeSector ) fHMeanTimeSector=(TH2F*)sig.fHMeanTimeSector->Clone();
340  fParam->Update();
341 }
344 fNbinsT0(200),
345 fXminT0(-2),
346 fXmaxT0(2),
347 fNbinsQ(200),
348 fXminQ(10),
349 fXmaxQ(40),
350 fNbinsRMS(100),
351 fXminRMS(0.1),
352 fXmaxRMS(5.1),
353 fPeakIntMinus(2),
354 fPeakIntPlus(2),
355 fIsZeroSuppressed(kFALSE),
356 fLastSector(-1),
357 fParam(new AliTPCParam),
358 fPedestalTPC(0x0),
359 fPadNoiseTPC(0x0),
360 fOutliers(0x0),
361 fPedestalROC(0x0),
362 fPadNoiseROC(0x0),
363 fCalRocArrayT0(72),
364 fCalRocArrayQ(72),
365 fCalRocArrayRMS(72),
366 fCalRocArrayOutliers(72),
367 fHistoQArray(72),
368 fHistoT0Array(72),
369 fHistoRMSArray(72),
370 fHMeanTimeSector(0x0),
371 fVMeanTimeSector(72),
372 fPadTimesArrayEvent(72),
373 fPadQArrayEvent(72),
374 fPadRMSArrayEvent(72),
375 fPadPedestalArrayEvent(72),
376 fCurrentChannel(-1),
377 fCurrentSector(-1),
378 fCurrentRow(-1),
379 fCurrentPad(-1),
380 fMaxPadSignal(-1),
381 fMaxTimeBin(-1),
382 fPadSignal(1024),
383 fPadPedestal(0),
384 fPadNoise(0),
385 fVTime0Offset(72),
386 fVTime0OffsetCounter(72)
387 {
389 
390  SetNameTitle("AliTPCCalibPulser","AliTPCCalibPulser");
391  fFirstTimeBin=60;
392  fLastTimeBin=900;
393  if (config->GetValue("FirstTimeBin")) fFirstTimeBin = ((TObjString*)config->GetValue("FirstTimeBin"))->GetString().Atoi();
394  if (config->GetValue("LastTimeBin")) fLastTimeBin = ((TObjString*)config->GetValue("LastTimeBin"))->GetString().Atoi();
395  if (config->GetValue("NbinsT0")) fNbinsT0 = ((TObjString*)config->GetValue("NbinsT0"))->GetString().Atoi();
396  if (config->GetValue("XminT0")) fXminT0 = ((TObjString*)config->GetValue("XminT0"))->GetString().Atof();
397  if (config->GetValue("XmaxT0")) fXmaxT0 = ((TObjString*)config->GetValue("XmaxT0"))->GetString().Atof();
398  if (config->GetValue("NbinsQ")) fNbinsQ = ((TObjString*)config->GetValue("NbinsQ"))->GetString().Atoi();
399  if (config->GetValue("XminQ")) fXminQ = ((TObjString*)config->GetValue("XminQ"))->GetString().Atof();
400  if (config->GetValue("XmaxQ")) fXmaxQ = ((TObjString*)config->GetValue("XmaxQ"))->GetString().Atof();
401  if (config->GetValue("NbinsRMS")) fNbinsRMS = ((TObjString*)config->GetValue("NbinsRMS"))->GetString().Atoi();
402  if (config->GetValue("XminRMS")) fXminRMS = ((TObjString*)config->GetValue("XminRMS"))->GetString().Atof();
403  if (config->GetValue("XmaxRMS")) fXmaxRMS = ((TObjString*)config->GetValue("XmaxRMS"))->GetString().Atof();
404  if (config->GetValue("PeakIntMinus")) fPeakIntMinus = (Int_t)((TObjString*)config->GetValue("PeakIntMinus"))->GetString().Atof();
405  if (config->GetValue("PeakIntPlus")) fPeakIntPlus = (Int_t)((TObjString*)config->GetValue("PeakIntPlus"))->GetString().Atof();
406  if (config->GetValue("IsZeroSuppressed")) fIsZeroSuppressed = (Bool_t)((TObjString*)config->GetValue("IsZeroSuppressed"))->GetString().Atoi();
407 
408  fParam->Update();
409 }
410 //_____________________________________________________________________
412 {
414 
415  if (&source == this) return *this;
416  new (this) AliTPCCalibPulser(source);
417 
418  return *this;
419 }
420 //_____________________________________________________________________
422 {
424 
425  Reset();
426  delete fParam;
427 }
429 {
431 
432  fCalRocArrayT0.Delete();
433  fCalRocArrayQ.Delete();
434  fCalRocArrayRMS.Delete();
435  fCalRocArrayOutliers.Delete();
436 
437  fHistoQArray.Delete();
438  fHistoT0Array.Delete();
439  fHistoRMSArray.Delete();
440 
441  fPadTimesArrayEvent.Delete();
442  fPadQArrayEvent.Delete();
443  fPadRMSArrayEvent.Delete();
444  fPadPedestalArrayEvent.Delete();
445 
447 }
448 //_____________________________________________________________________
449 Int_t AliTPCCalibPulser::Update(const Int_t icsector,
450  const Int_t icRow,
451  const Int_t icPad,
452  const Int_t icTimeBin,
453  const Float_t csignal)
454 {
458 
459  if (icRow<0) return 0;
460  if (icPad<0) return 0;
461  if (icTimeBin<0) return 0;
462  if ( (icTimeBin>fLastTimeBin) || (icTimeBin<fFirstTimeBin) ) return 0;
463 
464  if ( icRow<0 || icPad<0 ){
465  AliWarning("Wrong Pad or Row number, skipping!");
466  return 0;
467  }
468 
469  Int_t iChannel = fROC->GetRowIndexes(icsector)[icRow]+icPad; // global pad position in sector
470 
471  //init first pad and sector in this event
472  if ( fCurrentChannel == -1 ) {
473  fCurrentChannel = iChannel;
474  fCurrentSector = icsector;
475  fCurrentRow = icRow;
476  fCurrentPad = icPad;
477  }
478 
479  //process last pad if we change to a new one
480  if ( iChannel != fCurrentChannel ){
481  ProcessPad();
483  fCurrentChannel = iChannel;
484  fCurrentSector = icsector;
485  fCurrentRow = icRow;
486  fCurrentPad = icPad;
487  }
488 
489  //fill signals for current pad
490  fPadSignal[icTimeBin]=csignal;
491  if ( csignal > fMaxPadSignal ){
492  fMaxPadSignal = csignal;
493  fMaxTimeBin = icTimeBin;
494  }
495  return 0;
496 }
497 //_____________________________________________________________________
499 {
502 
503  Bool_t noPedestal = kTRUE;;
505  //use pedestal database
506  //only load new pedestals if the sector has changed
507  if ( fCurrentSector!=fLastSector ){
510  }
511 
512  if ( fPedestalROC&&fPadNoiseROC ){
515  noPedestal = kFALSE;
516  }
517 
518  }
519 
520  //if we are not running with pedestal database, or for the current sector there is no information
521  //available, calculate the pedestal and noise on the fly
522  if ( noPedestal ) {
523  const Int_t kPedMax = 100; //maximum pedestal value
524  Float_t max = 0;
525  Int_t median = -1;
526  Int_t count0 = 0;
527  Int_t count1 = 0;
528  //
529  Float_t padSignal=0;
530  //
531  UShort_t histo[kPedMax];
532  memset(histo,0,kPedMax*sizeof(UShort_t));
533 
534  for (Int_t i=fFirstTimeBin; i<=fLastTimeBin; ++i){
535  padSignal = fPadSignal.GetMatrixArray()[i];
536  if (padSignal<=0) continue;
537  if (padSignal>max && i>10) {
538  max = padSignal;
539  }
540  if (padSignal>kPedMax-1) continue;
541  histo[Int_t(padSignal+0.5)]++;
542  count0++;
543  }
544  //
545  for (Int_t i=1; i<kPedMax; ++i){
546  if (count1<count0*0.5) median=i;
547  count1+=histo[i];
548  }
549  // truncated mean
550  //
551  // what if by chance histo[median] == 0 ?!?
552  Float_t count=histo[median] ,mean=histo[median]*median, rms=histo[median]*median*median ;
553  //
554  for (Int_t idelta=1; idelta<10; ++idelta){
555  if (median-idelta<=0) continue;
556  if (median+idelta>kPedMax) continue;
557  if (count<part*count1){
558  count+=histo[median-idelta];
559  mean +=histo[median-idelta]*(median-idelta);
560  rms +=histo[median-idelta]*(median-idelta)*(median-idelta);
561  count+=histo[median+idelta];
562  mean +=histo[median+idelta]*(median+idelta);
563  rms +=histo[median+idelta]*(median+idelta)*(median+idelta);
564  }
565  }
566  fPadPedestal = 0;
567  fPadNoise = 0;
568  if ( count > 0 ) {
569  mean/=count;
570  rms = TMath::Sqrt(TMath::Abs(rms/count-mean*mean));
571  fPadPedestal = mean;
572  fPadNoise = rms;
573  }
574  }
575  fPadPedestal*=(Float_t)(!fIsZeroSuppressed);
576 }
577 //_____________________________________________________________________
578 void AliTPCCalibPulser::FindPulserSignal(TVectorD &param, Float_t &qSum)
579 {
583 
584  Float_t ceQmax =0, ceQsum=0, ceTime=0, ceRMS=0;
585  Int_t cemaxpos = fMaxTimeBin;
586  Float_t ceSumThreshold = 10.*TMath::Max(fPadNoise,Float_t(1.)); // threshold for the signal sum
587  Float_t ceMaxThreshold = 5.*TMath::Max(fPadNoise,Float_t(1.)); // threshold for the signal max
588  const Int_t kCemin = fPeakIntMinus; // range for the analysis of the ce signal +- channels from the peak
589  const Int_t kCemax = fPeakIntPlus;
590  param[0] = ceQmax;
591  param[1] = ceTime;
592  param[2] = ceRMS;
593  qSum = ceQsum;
594 
595  if (cemaxpos>0){
596  ceQmax = fPadSignal.GetMatrixArray()[cemaxpos]-fPadPedestal;
597  if ( ceQmax<ceMaxThreshold ) return;
598  for (Int_t i=cemaxpos-kCemin; i<=cemaxpos+kCemax; ++i){
599  Float_t signal = fPadSignal.GetMatrixArray()[i]-fPadPedestal;
600  if ( (i>fFirstTimeBin) && (i<fLastTimeBin) && (signal>0) ){
601  ceTime+=signal*(i+0.5);
602  ceRMS +=signal*(i+0.5)*(i+0.5);
603  ceQsum+=signal;
604  }
605  }
606  }
607  if (ceQsum>ceSumThreshold) {
608  ceTime/=ceQsum;
609  ceRMS = TMath::Sqrt(TMath::Abs(ceRMS/ceQsum-ceTime*ceTime));
610  ceTime-=GetL1PhaseTB();
611  //only fill the Time0Offset if pad was not marked as an outlier!
612  if ( !fOutliers ){
613  //skip edge pads for calculating the mean time
615  fVTime0Offset.GetMatrixArray()[fCurrentSector]+=ceTime; // mean time for each sector
616  fVTime0OffsetCounter.GetMatrixArray()[fCurrentSector]++;
617  GetHistoTSec()->Fill(ceTime,fCurrentSector);
618  }
619  } else {
621  fVTime0Offset.GetMatrixArray()[fCurrentSector]+=ceTime; // mean time for each sector
622  fVTime0OffsetCounter.GetMatrixArray()[fCurrentSector]++;
623  }
624  }
625 
626  //Normalise Q to the 'cell-size': The wire density is the same in the IROC and OROC, therefore the
627  // the pick-up signal should scale with the pad area. In addition
628  // the signal should decrease with the wire distance (4mm in IROC, 6mm in OROC),
629  // ratio 2/3. The pad area we express in mm2 (factor 100). We normalise the signal
630  // to the OROC signal (factor 2/3 for the IROCs).
632  if ( fCurrentSector<fParam->GetNInnerSector() ) norm*=3./2.;
633 
634  ceQsum/=norm;
635  } else {
636  ceQmax=0;
637  ceTime=0;
638  ceRMS =0;
639  ceQsum=0;
640  }
641  param[0] = ceQmax;
642  param[1] = ceTime;
643  param[2] = ceRMS;
644  qSum = ceQsum;
645 }
646 //_____________________________________________________________________
648 {
650 
651  FindPedestal();
652  TVectorD param(3);
653  Float_t qSum;
654  FindPulserSignal(param, qSum);
655 
656  Double_t meanT = param[1];
657  Double_t sigmaT = param[2];
658 
659 
660  //Fill Event T0 counter
661  (*GetPadTimesEvent(fCurrentSector,kTRUE)).GetMatrixArray()[fCurrentChannel] = meanT;
662 
663  //Fill Q histogram
664 // GetHistoQ(fCurrentSector,kTRUE)->Fill( TMath::Sqrt(qSum), fCurrentChannel ); //use linear scale, needs also a change in the Analyse funciton.
665  GetHistoQ(fCurrentSector,kTRUE)->Fill( qSum, fCurrentChannel );
666 
667  //Fill RMS histogram
668  GetHistoRMS(fCurrentSector,kTRUE)->Fill( sigmaT, fCurrentChannel );
669 
670 
671  //Fill debugging info
672  if ( GetStreamLevel()>0 ){
673  TTreeSRedirector *streamer=GetDebugStreamer();
674  if ( GetStreamLevel() == 1 ){
675  if ( streamer ) {
677  (*streamer) << "PadSignals" <<
678  "Event=" <<fNevents <<
679  "Sector=" <<fCurrentSector<<
680  "Row=" <<fCurrentRow<<
681  "Pad=" <<fCurrentPad<<
682  "PadC=" <<padc<<
683  "Channel="<<fCurrentChannel<<
684  "Sum=" <<qSum<<
685  "params.="<<&param<<
686  "signal.=" <<&fPadSignal<<
687  "\n";
688  }
689  } else { //debug > 1
693  }
694  }
695  ResetPad();
696 }
697 //_____________________________________________________________________
699 {
701 
702  //check if last pad has allready been processed, if not do so
703  if ( fMaxTimeBin>-1 ) ProcessPad();
704 
705  //loop over all ROCs, fill Time0 histogram corrected for the mean Time0 of each ROC
706  for ( Int_t iSec = 0; iSec<72; ++iSec ){
707  TVectorF *vTimes = GetPadTimesEvent(iSec);
708  if ( !vTimes || fVTime0OffsetCounter[iSec]==0 ) continue;
709  Float_t time0 = fVTime0Offset[iSec]/fVTime0OffsetCounter[iSec];
710  for ( UInt_t iChannel=0; iChannel<fROC->GetNChannels(iSec); ++iChannel ){
711  Float_t time = (*vTimes).GetMatrixArray()[iChannel];
712 
713  GetHistoT0(iSec,kTRUE)->Fill( time-time0,iChannel );
714  //GetHistoT0(iSec,kTRUE)->Fill( time,iChannel );
715 
716 
717  //Debug start
718  if ( GetStreamLevel()>1 ){
719  TTreeSRedirector *streamer=GetDebugStreamer();
720  if ( streamer ) {
721  Int_t row=0;
722  Int_t pad=0;
723  Int_t padc=0;
724 
725  Float_t q = (*GetPadQEvent(iSec)).GetMatrixArray()[iChannel];
726  Float_t rms = (*GetPadRMSEvent(iSec)).GetMatrixArray()[iChannel];
727 
728  UInt_t channel=iChannel;
729  Int_t sector=iSec;
730 
731  while ( channel > (fROC->GetRowIndexes(sector)[row]+fROC->GetNPads(sector,row)-1) ) row++;
732  pad = channel-fROC->GetRowIndexes(sector)[row];
733  padc = pad-(fROC->GetNPads(sector,row)/2);
734 
735  (*streamer) << "DataPad" <<
736  "Event=" << fNevents <<
737  "Sector="<< sector <<
738  "Row=" << row<<
739  "Pad=" << pad <<
740  "PadC=" << padc <<
741  "PadSec="<< channel <<
742  "Time0=" << time0 <<
743  "Time=" << time <<
744  "RMS=" << rms <<
745  "Sum=" << q <<
746  "\n";
747  }
748  }
749  //Debug end
750  }
751  }
752  ++fNevents;
753 }
754 //_____________________________________________________________________
755 TH2S* AliTPCCalibPulser::GetHisto(Int_t sector, TObjArray *arr,
756  Int_t nbinsY, Float_t ymin, Float_t ymax,
757  const Char_t *type, Bool_t force)
758 {
761 
762  if ( !force || arr->UncheckedAt(sector) )
763  return (TH2S*)arr->UncheckedAt(sector);
764 
765  // if we are forced and histogram doesn't yes exist create it
766  // new histogram with Q calib information. One value for each pad!
767  TH2S* hist = new TH2S(Form("hCalib%s%.2d",type,sector),Form("%s calibration histogram sector %.2d",type,sector),
768  nbinsY, ymin, ymax,
769  fROC->GetNChannels(sector),0,fROC->GetNChannels(sector));
770  hist->SetDirectory(0);
771  arr->AddAt(hist,sector);
772  return hist;
773 }
774 //_____________________________________________________________________
775 TH2S* AliTPCCalibPulser::GetHistoT0(Int_t sector, Bool_t force)
776 {
779 
780  TObjArray *arr = &fHistoT0Array;
781  return GetHisto(sector, arr, fNbinsT0, fXminT0, fXmaxT0, "T0", force);
782 }
783 //_____________________________________________________________________
784 TH2S* AliTPCCalibPulser::GetHistoQ(Int_t sector, Bool_t force)
785 {
788 
789  TObjArray *arr = &fHistoQArray;
790  return GetHisto(sector, arr, fNbinsQ, fXminQ, fXmaxQ, "Q", force);
791 }
792 //_____________________________________________________________________
793 TH2S* AliTPCCalibPulser::GetHistoRMS(Int_t sector, Bool_t force)
794 {
797 
798  TObjArray *arr = &fHistoRMSArray;
799  return GetHisto(sector, arr, fNbinsRMS, fXminRMS, fXmaxRMS, "RMS", force);
800 }
801 //_____________________________________________________________________
803 {
806 
807  if ( !fHMeanTimeSector )
808  fHMeanTimeSector = new TH2F("fHMeanTimeSector","Abs mean time per sector",
810  72,0,72);
811  return fHMeanTimeSector;
812 }
813 //_____________________________________________________________________
814 TVectorF* AliTPCCalibPulser::GetPadInfoEvent(Int_t sector, TObjArray *arr, Bool_t force)
815 {
818 
819  if ( !force || arr->UncheckedAt(sector) )
820  return (TVectorF*)arr->UncheckedAt(sector);
821 
822  TVectorF *vect = new TVectorF(fROC->GetNChannels(sector));
823  arr->AddAt(vect,sector);
824  return vect;
825 }
826 //_____________________________________________________________________
827 TVectorF* AliTPCCalibPulser::GetPadTimesEvent(Int_t sector, Bool_t force)
828 {
831 
833  return GetPadInfoEvent(sector,arr,force);
834 }
835 //_____________________________________________________________________
836 TVectorF* AliTPCCalibPulser::GetPadQEvent(Int_t sector, Bool_t force)
837 {
841 
842  TObjArray *arr = &fPadQArrayEvent;
843  return GetPadInfoEvent(sector,arr,force);
844 }
845 //_____________________________________________________________________
846 TVectorF* AliTPCCalibPulser::GetPadRMSEvent(Int_t sector, Bool_t force)
847 {
851 
853  return GetPadInfoEvent(sector,arr,force);
854 }
855 //_____________________________________________________________________
856 TVectorF* AliTPCCalibPulser::GetPadPedestalEvent(Int_t sector, Bool_t force)
857 {
861 
863  return GetPadInfoEvent(sector,arr,force);
864 }
865 //_____________________________________________________________________
866 AliTPCCalROC* AliTPCCalibPulser::GetCalRoc(Int_t sector, TObjArray* arr, Bool_t force) const
867 {
870 
871  if ( !force || arr->UncheckedAt(sector) )
872  return (AliTPCCalROC*)arr->UncheckedAt(sector);
873 
874  // if we are forced and histogram doesn't yes exist create it
875 
876  // new AliTPCCalROC for T0 information. One value for each pad!
877  AliTPCCalROC *croc = new AliTPCCalROC(sector);
878  arr->AddAt(croc,sector);
879  return croc;
880 }
881 //_____________________________________________________________________
882 AliTPCCalROC* AliTPCCalibPulser::GetCalRocT0(Int_t sector, Bool_t force)
883 {
886 
887  TObjArray *arr = &fCalRocArrayT0;
888  return GetCalRoc(sector, arr, force);
889 }
890 //_____________________________________________________________________
891 AliTPCCalROC* AliTPCCalibPulser::GetCalRocQ(Int_t sector, Bool_t force)
892 {
895 
896  TObjArray *arr = &fCalRocArrayQ;
897  return GetCalRoc(sector, arr, force);
898 }
899 //_____________________________________________________________________
900 AliTPCCalROC* AliTPCCalibPulser::GetCalRocRMS(Int_t sector, Bool_t force)
901 {
904 
905  TObjArray *arr = &fCalRocArrayRMS;
906  return GetCalRoc(sector, arr, force);
907 }
908 //_____________________________________________________________________
910 {
913 
915  return GetCalRoc(sector, arr, force);
916 }
917 //_____________________________________________________________________
919 {
921 
922  fLastSector=-1;
923  fCurrentSector=-1;
924  fCurrentRow=-1;
925  fCurrentPad=-1;
926  fCurrentChannel=-1;
927 
928  ResetPad();
929 
930  fPadTimesArrayEvent.Delete();
931 
932  fPadQArrayEvent.Delete();
933  fPadRMSArrayEvent.Delete();
934  fPadPedestalArrayEvent.Delete();
935 
936  for ( Int_t i=0; i<72; ++i ){
937  fVTime0Offset[i]=0;
939  }
940 }
941 //_____________________________________________________________________
943 {
945 
946  for (Int_t i=fFirstTimeBin; i<fLastTimeBin+1; ++i)
947  fPadSignal[i] = 0;
948  fMaxTimeBin = -1;
949  fMaxPadSignal = -1;
950  fPadPedestal = -1;
951  fPadNoise = -1;
952 }
953 //_____________________________________________________________________
954 Bool_t AliTPCCalibPulser::IsEdgePad(Int_t sector, Int_t row, Int_t pad)
955 {
957 
958  Int_t edge1 = 0;
959  Int_t edge2 = fROC->GetNPads(sector,row)-1;
960  if ( pad == edge1 || pad == edge2 ) return kTRUE;
961 
962  return kFALSE;
963 }
964 //_____________________________________________________________________
966 {
968 
969  MergeBase(sig);
970  //merge histograms
971  for (Int_t iSec=0; iSec<72; ++iSec){
972  TH2S *hRefQmerge = sig->GetHistoQ(iSec);
973  TH2S *hRefT0merge = sig->GetHistoT0(iSec);
974  TH2S *hRefRMSmerge = sig->GetHistoRMS(iSec);
975 
976 
977  if ( hRefQmerge ){
978  TDirectory *dir = hRefQmerge->GetDirectory(); hRefQmerge->SetDirectory(0);
979  TH2S *hRefQ = GetHistoQ(iSec);
980  if ( hRefQ ) hRefQ->Add(hRefQmerge);
981  else {
982  TH2S *hist = new TH2S(*hRefQmerge);
983  hist->SetDirectory(0);
984  fHistoQArray.AddAt(hist, iSec);
985  }
986  hRefQmerge->SetDirectory(dir);
987  }
988  if ( hRefT0merge ){
989  TDirectory *dir = hRefT0merge->GetDirectory(); hRefT0merge->SetDirectory(0);
990  TH2S *hRefT0 = GetHistoT0(iSec);
991  if ( hRefT0 ) hRefT0->Add(hRefT0merge);
992  else {
993  TH2S *hist = new TH2S(*hRefT0merge);
994  hist->SetDirectory(0);
995  fHistoT0Array.AddAt(hist, iSec);
996  }
997  hRefT0merge->SetDirectory(dir);
998  }
999  if ( hRefRMSmerge ){
1000  TDirectory *dir = hRefRMSmerge->GetDirectory(); hRefRMSmerge->SetDirectory(0);
1001  TH2S *hRefRMS = GetHistoRMS(iSec);
1002  if ( hRefRMS ) hRefRMS->Add(hRefRMSmerge);
1003  else {
1004  TH2S *hist = new TH2S(*hRefRMSmerge);
1005  hist->SetDirectory(0);
1006  fHistoRMSArray.AddAt(hist, iSec);
1007  }
1008  hRefRMSmerge->SetDirectory(dir);
1009  }
1010 
1011  }
1012  if ( sig->fHMeanTimeSector ){
1013  TDirectory *dir = sig->fHMeanTimeSector->GetDirectory(); sig->fHMeanTimeSector->SetDirectory(0);
1015  else {
1016  fHMeanTimeSector = new TH2F(*sig->fHMeanTimeSector);
1017  fHMeanTimeSector->SetDirectory(0);
1018  }
1019  sig->fHMeanTimeSector->SetDirectory(dir);
1020  }
1021 }
1022 
1023 
1024 //_____________________________________________________________________
1025 Long64_t AliTPCCalibPulser::Merge(TCollection * const list)
1026 {
1028 
1029  Long64_t nmerged=1;
1030 
1031  TIter next(list);
1032  AliTPCCalibPulser *ce=0;
1033  TObject *o=0;
1034 
1035  while ( (o=next()) ){
1036  ce=dynamic_cast<AliTPCCalibPulser*>(o);
1037  if (ce){
1038  Merge(ce);
1039  ++nmerged;
1040  }
1041  }
1042 
1043  return nmerged;
1044 }
1045 
1046 //_____________________________________________________________________
1048 {
1050 
1051  TVectorD paramQ(3);
1052  TVectorD paramT0(3);
1053  TVectorD paramRMS(3);
1054  TMatrixD dummy(3,3);
1055  //calculate mean time for each sector and mean time for each side
1056  TH1F hMeanTsec("hMeanTsec","hMeanTsec",20*(fLastTimeBin-fFirstTimeBin),fFirstTimeBin,fLastTimeBin);
1057  fVMeanTimeSector.Zero();
1058 
1059  for (Int_t iSec=0; iSec<72; ++iSec){
1060  TH2S *hT0 = GetHistoT0(iSec);
1061  if (!hT0 ) continue;
1062  //calculate sector mean T
1063  if ( fHMeanTimeSector ){
1064  Int_t nbinsT = fHMeanTimeSector->GetNbinsX();
1065  Int_t offset = (nbinsT+2)*(iSec+1);
1066  Float_t *arrP=fHMeanTimeSector->GetArray()+offset;
1067  Int_t entries=0;
1068  for ( Int_t i=0; i<nbinsT; i++ ) entries+=(Int_t)arrP[i+1];
1069  hMeanTsec.Set(nbinsT+2,arrP);
1070  hMeanTsec.SetEntries(entries);
1071  paramT0.Zero();
1072  // truncated mean: remove lower 5% and upper 5%
1073  if ( entries>0 ) AliMathBase::TruncatedMean(&hMeanTsec,&paramT0,0.05,.95);
1074  fVMeanTimeSector[iSec]=paramT0[1];
1075  }
1076 
1077  AliTPCCalROC *rocQ = GetCalRocQ (iSec,kTRUE);
1078  AliTPCCalROC *rocT0 = GetCalRocT0 (iSec,kTRUE);
1079  AliTPCCalROC *rocRMS = GetCalRocRMS(iSec,kTRUE);
1080  AliTPCCalROC *rocOut = GetCalRocOutliers(iSec,kTRUE);
1081 
1082  TH2S *hQ = GetHistoQ(iSec);
1083  TH2S *hRMS = GetHistoRMS(iSec);
1084 
1085  Short_t *arrayhQ = hQ->GetArray();
1086  Short_t *arrayhT0 = hT0->GetArray();
1087  Short_t *arrayhRMS = hRMS->GetArray();
1088 
1089  UInt_t nChannels = fROC->GetNChannels(iSec);
1090  Float_t meanTsec = fVMeanTimeSector[iSec];
1091 
1092  //debug
1093  Int_t row=0;
1094  Int_t pad=0;
1095  Int_t padc=0;
1097 
1098  for (UInt_t iChannel=0; iChannel<nChannels; ++iChannel){
1099 
1100  Float_t cogTime0 = -1000;
1101  Float_t cogQ = -1000;
1102  Float_t cogRMS = -1000;
1103  Float_t cogOut = 0;
1104 
1105  Int_t offsetQ = (fNbinsQ+2)*(iChannel+1)+1;
1106  Int_t offsetT0 = (fNbinsT0+2)*(iChannel+1)+1;
1107  Int_t offsetRMS = (fNbinsRMS+2)*(iChannel+1)+1;
1108 /*
1109  AliMathBase::FitGaus(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ,&paramQ,&dummy);
1110  AliMathBase::FitGaus(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0,&paramT0,&dummy);
1111  AliMathBase::FitGaus(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS,&paramRMS,&dummy);
1112  cogQ = paramQ[1];
1113  cogTime0 = paramT0[1];
1114  cogRMS = paramRMS[1];
1115 */
1116  cogQ = AliMathBase::GetCOG(arrayhQ+offsetQ,fNbinsQ,fXminQ,fXmaxQ);
1117  cogTime0 = AliMathBase::GetCOG(arrayhT0+offsetT0,fNbinsT0,fXminT0,fXmaxT0);
1118  cogRMS = AliMathBase::GetCOG(arrayhRMS+offsetRMS,fNbinsRMS,fXminRMS,fXmaxRMS);
1119 
1120  /*
1121  if ( (cogQ < ??) && (cogTime0 > ??) && (cogTime0<??) && ( cogRMS>??) ){
1122  cogOut = 1;
1123  cogTime0 = 0;
1124  cogQ = 0;
1125  cogRMS = 0;
1126  }
1127 */
1128 // rocQ->SetValue(iChannel, cogQ*cogQ); // changed to linear scale again
1129  rocQ->SetValue(iChannel, cogQ);
1130  rocT0->SetValue(iChannel, cogTime0+meanTsec); //offset by mean time of the sector
1131  rocRMS->SetValue(iChannel, cogRMS);
1132  rocOut->SetValue(iChannel, cogOut);
1133 
1134  // in case a channel has no data set the value to 0
1135  if (TMath::Abs(cogTime0-fXminT0)<1e-10){
1136  rocQ->SetValue(iChannel, 0);
1137  rocT0->SetValue(iChannel, 0); //offset by mean time of the sector
1138  rocRMS->SetValue(iChannel, 0);
1139  }
1140 
1141  //debug
1142  if ( GetStreamLevel() > 2 ){
1143  TTreeSRedirector *streamer=GetDebugStreamer();
1144  if ( streamer ) {
1145  while ( iChannel > (fROC->GetRowIndexes(iSec)[row]+fROC->GetNPads(iSec,row)-1) ) row++;
1146  pad = iChannel-fROC->GetRowIndexes(iSec)[row];
1147  padc = pad-(fROC->GetNPads(iSec,row)/2);
1148 
1149  (*streamer) << "DataEnd" <<
1150  "Sector=" << iSec <<
1151  "Pad=" << pad <<
1152  "PadC=" << padc <<
1153  "Row=" << row <<
1154  "PadSec=" << iChannel <<
1155  "Q=" << cogQ <<
1156  "T0=" << cogTime0 <<
1157  "RMS=" << cogRMS <<
1158  "\n";
1159  }
1160  }
1162  }
1163 
1164 
1165  }
1166 }
1167 //_____________________________________________________________________
1168 //_________________________ Test Functions ___________________________
1169 //_____________________________________________________________________
1171 {
1178 
1179 
1180  TObjArray *histArray = new TObjArray(6);
1181  const Char_t *type[] = {"T0","Q","RMS"};
1182  Int_t fNbins[3] = {fNbinsT0,fNbinsQ,fNbinsRMS};
1183 
1184  for (Int_t itype = 0; itype<3; ++itype){
1185  for (Int_t imode=0; imode<2; ++imode){
1186  Int_t icount = itype*2+imode;
1187  histArray->AddAt(new TH1F(Form("hTestBinning%s%d",type[itype],imode),
1188  Form("Test Binning of '%s', mode - %d",type[itype],imode),
1189  72,0,72),
1190  icount);
1191  }
1192  }
1193 
1194 
1195  TH2S *hRef=0x0;
1196  Short_t *array=0x0;
1197  for (Int_t itype = 0; itype<3; ++itype){
1198  for (Int_t iSec=0; iSec<72; ++iSec){
1199  if ( itype == 0 ) hRef = GetHistoT0(iSec);
1200  if ( itype == 1 ) hRef = GetHistoQ(iSec);
1201  if ( itype == 2 ) hRef = GetHistoRMS(iSec);
1202  if ( hRef == 0x0 ) continue;
1203  array = (hRef->GetArray());
1204  UInt_t nChannels = fROC->GetNChannels(iSec);
1205 
1206  Int_t nempty=0;
1207  for (UInt_t iChannel=0; iChannel<nChannels; ++iChannel){
1208  Int_t nfilled=0;
1209  Int_t offset = (fNbins[itype]+2)*(iChannel+1)+1;
1210  Int_t c1 = 0;
1211  Int_t c2 = 0;
1212  for (Int_t iBin=0; iBin<fNbins[itype]; ++iBin){
1213  if ( array[offset+iBin]>0 ) {
1214  nfilled++;
1215  if ( c1 && c2 ) nempty++;
1216  else c1 = 1;
1217  }
1218  else if ( c1 ) c2 = 1;
1219 
1220  }
1221  ((TH1F*)histArray->At(itype*2))->Fill(nfilled);
1222  }
1223  ((TH1F*)histArray->At(itype*2+1))->Fill(iSec,nempty);
1224  }
1225  }
1226  return histArray;
1227 }
Float_t fXmaxRMS
xmax of T0 reference histogram
Float_t GetPadPitchLength(Int_t isector=0, Int_t padrow=0) const
Definition: AliTPCParam.h:299
Double_t GetL1PhaseTB() const
Int_t fNbinsT0
Number of bins for T0 reference histogram.
UInt_t GetNPads(UInt_t sector, UInt_t row) const
Definition: AliTPCROC.h:29
AliTPCCalROC * GetCalRocRMS(Int_t sector, Bool_t force=kFALSE)
Bool_t fIsZeroSuppressed
if data is zero suppressed
Int_t fCurrentChannel
! current channel processed
AliTPCCalROC * GetCalROC(Int_t sector) const
Definition: AliTPCCalPad.h:42
AliTPCCalROC * GetCalRocOutliers(Int_t sector, Bool_t force=kFALSE)
AliTPCCalROC * fPedestalROC
! Pedestal Information for current ROC
TH2F * fHMeanTimeSector
Timing distribution per sector.
Float_t fXmaxQ
xmax of T0 reference histogram
TVectorF * GetPadTimesEvent(Int_t sector, Bool_t force=kFALSE)
#define TObjArray
Int_t fNevents
Number of processed events.
Float_t fPadPedestal
! Pedestal Value of current pad
TH2S * GetHisto(Int_t sector, TObjArray *arr, Int_t nbinsY, Float_t ymin, Float_t ymax, const Char_t *type, Bool_t force)
Manager and of geomety classes for set: TPC.
Definition: AliTPCParam.h:18
TObjArray fHistoQArray
Calibration histograms for Charge distribution.
void MergeBase(const AliTPCCalibRawBase *calib)
Float_t GetValue(UInt_t row, UInt_t pad) const
Definition: AliTPCCalROC.h:33
Float_t fPadNoise
! Noise Value of current pad
Float_t fXminQ
xmin of T0 reference histogram
TObjArray fPadTimesArrayEvent
! Pad Times for the event, before mean Time0 corrections
TObjArray fPadPedestalArrayEvent
! Signal width for the event, only needed for debugging streamer
TH2S * GetHistoQ(Int_t sector, Bool_t force=kFALSE)
UInt_t GetNChannels(UInt_t sector) const
Definition: AliTPCROC.h:28
AliTPCCalPad * fPadNoiseTPC
! Pad noise Information whole TPC
Int_t fLastTimeBin
Last Time bin used for analysis.
Float_t GetPadPitchWidth(Int_t isector=0) const
Definition: AliTPCParam.h:297
Base class for the calibration algorithms using raw data as input.
AliTPCCalROC * fPadNoiseROC
! Pad noise Information for current ROC
AliTPCParam * fParam
! TPC information
TTreeSRedirector * GetDebugStreamer()
Float_t fMaxPadSignal
! maximum bin of current pad
TObjArray * array
Definition: AnalyzeLaser.C:12
ClassImp(TPCGenInfo)
Definition: AliTPCCmpNG.C:254
Implementation of the TPC pulser calibration.
AliTPCCalPad * fPedestalTPC
! Pedestal Information
TObjArray fCalRocArrayT0
Array of AliTPCCalROC class for Time0 calibration.
AliTPCROC * fROC
! ROC information
TVectorF fPadSignal
! signal of current Pad
Int_t fFirstTimeBin
First Time bin used for analysis.
Float_t fXminRMS
xmin of T0 reference histogram
AliTPCCalROC * GetCalRocQ(Int_t sector, Bool_t force=kFALSE)
Int_t fCurrentRow
! current row processed
const UInt_t * GetRowIndexes(UInt_t sector) const
Definition: AliTPCROC.h:32
AliTPCCalROC * GetCalRocT0(Int_t sector, Bool_t force=kFALSE)
void SetValue(UInt_t row, UInt_t pad, Float_t vd)
Definition: AliTPCCalROC.h:35
Int_t fLastSector
! Last sector processed
virtual Int_t Update(const Int_t isector, const Int_t iRow, const Int_t iPad, const Int_t iTimeBin, const Float_t signal)
TH2S * GetHistoT0(Int_t sector, Bool_t force=kFALSE)
Int_t fPeakIntMinus
Peak integral range for COG determination. Bins used before max bin.
Float_t fXmaxT0
xmax of T0 reference histogram
TObjArray fPadRMSArrayEvent
! Signal width for the event, only needed for debugging streamer
TPC calibration base class for one ROC.
Definition: AliTPCCalROC.h:19
void Merge(AliTPCCalibPulser *const sig)
TVectorF * GetPadRMSEvent(Int_t sector, Bool_t force=kFALSE)
virtual Bool_t Update()
AliTPCCalROC * GetCalRoc(Int_t sector, TObjArray *arr, Bool_t force) const
TVectorF * GetPadPedestalEvent(Int_t sector, Bool_t force=kFALSE)
TVectorF * GetPadQEvent(Int_t sector, Bool_t force=kFALSE)
TH2S * GetHistoRMS(Int_t sector, Bool_t force=kFALSE)
TVectorF fVTime0Offset
! Time0 Offset from preprocessing for each sector;
Int_t fNbinsRMS
Number of bins for T0 reference histogram.
TObjArray fCalRocArrayRMS
Array of AliTPCCalROC class for signal width calibration.
AliTPCCalibPulser & operator=(const AliTPCCalibPulser &source)
void FindPedestal(Float_t part=.6)
TObjArray fCalRocArrayQ
Array of AliTPCCalROC class for Charge calibration.
TVectorF fVTime0OffsetCounter
! Time0 Offset from preprocessing for each sector;
Int_t fPeakIntPlus
Peak integral range for COG determination. Bins used after max bin.
virtual void EndEvent()
TObjArray fHistoT0Array
Calibration histograms for Time0 distribution.
Int_t GetStreamLevel() const
TVectorF fVMeanTimeSector
Mean time per sector from analysis of fHMeanTimeSector.
TVectorF * GetPadInfoEvent(Int_t sector, TObjArray *arr, Bool_t force=kFALSE)
Int_t fCurrentSector
! current sector processed
virtual void ResetEvent()
TObjArray fCalRocArrayOutliers
Array of AliTPCCalROC class for signal outliers.
Bool_t IsEdgePad(Int_t sector, Int_t row, Int_t pad)
Int_t fNbinsQ
Number of bins for T0 reference histogram.
TObjArray fPadQArrayEvent
! Charge for the event, only needed for debugging streamer
void FindPulserSignal(TVectorD &param, Float_t &qSum)
Float_t fXminT0
xmin of T0 reference histogram
AliTPCCalPad * fOutliers
! Outlier information. Those will not be used for calculating the T0
TObjArray fHistoRMSArray
Calibration histograms for signal width distribution.
Int_t fCurrentPad
! current pad processed
Int_t fMaxTimeBin
! time bin with maximum value