]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCdataQA.cxx
Update of the cluster finder (Marian, Adam Matyja)
[u/mrichter/AliRoot.git] / TPC / AliTPCdataQA.cxx
index c7c436c73716862429ce3b04c4e179faee313db7..2e9ff87f41c007716ab36dc27c6047b08e21401b 100644 (file)
 /* $Id$ */
 
 
+// stl includes
+#include <iostream>
+
+using namespace std;
+
 //Root includes
 #include <TH1F.h>
 #include <TH2F.h>
 //date
 #include "event.h"
 #include "AliTPCCalPad.h"
+#include "AliTPCPreprocessorOnline.h"
 
 //header file
 #include "AliTPCdataQA.h"
 
-
-
 ClassImp(AliTPCdataQA)
 
-AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
-  TObject(),
+AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/  
+  TH1F("TPCRAW","TPCRAW",100,0,100),
   fFirstTimeBin(60),
   fLastTimeBin(1000),
+  fMaxTime(1100),
   fAdcMin(1),
   fAdcMax(100),
   fOldRCUformat(kTRUE),
   fROC(AliTPCROC::Instance()),
   fMapping(NULL),
+  fPedestal(0),
+  fNoise(0),
   fMaxCharge(0),
+  fMeanCharge(0),
+  fNoThreshold(0),
   fOverThreshold0(0),
   fOverThreshold5(0),
   fOverThreshold10(0),
@@ -68,14 +77,22 @@ AliTPCdataQA::AliTPCdataQA() : /*FOLD00*/
   //
   // default constructor
   //
+
+  fSectorLast  = -1;
+  fRowLast     =  0;
+  fPadLast     =  0;
+  fTimeBinLast =  0;
+  fSignalLast  =  0;
+  fNAboveThreshold = 0;
 }
 
 
 //_____________________________________________________________________
 AliTPCdataQA::AliTPCdataQA(const AliTPCdataQA &ped) : /*FOLD00*/
-  TObject(ped),
+  TH1F(ped),
   fFirstTimeBin(ped.GetFirstTimeBin()),
   fLastTimeBin(ped.GetLastTimeBin()),
+  fMaxTime(ped.fMaxTime),
   fAdcMin(ped.GetAdcMin()),
   fAdcMax(ped.GetAdcMax()),
   fOldRCUformat(ped.fOldRCUformat),
@@ -129,10 +146,11 @@ Bool_t AliTPCdataQA::ProcessEventFast(AliTPCRawStreamFast *rawStreamFast)
          Int_t isector  = rawStreamFast->GetSector();                       //  current sector
          Int_t iRow     = rawStreamFast->GetRow();                          //  current row
          Int_t iPad     = rawStreamFast->GetPad();                          //  current pad
-         Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
-          Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
 
          while ( rawStreamFast->NextBunch() ){
+           Int_t startTbin = (Int_t)rawStreamFast->GetStartTimeBin();
+           Int_t endTbin = (Int_t)rawStreamFast->GetEndTimeBin();
+
              for (Int_t iTimeBin = startTbin; iTimeBin < endTbin; iTimeBin++){
                  Float_t signal=(Float_t)rawStreamFast->GetSignals()[iTimeBin-startTbin];
                  Update(isector,iRow,iPad,iTimeBin+1,signal);
@@ -150,8 +168,13 @@ Bool_t AliTPCdataQA::ProcessEventFast(AliRawReader *rawReader)
   //
   //  Event processing loop - AliRawReader
   //
+  fSectorLast  = -1;
   AliTPCRawStreamFast *rawStreamFast = new AliTPCRawStreamFast(rawReader, (AliAltroMapping**)fMapping);
   Bool_t res=ProcessEventFast(rawStreamFast);
+  if(res)
+    fEventCounter++; // only increment event counter if there is TPC data
+                     // otherwise Analyse (called in QA) fails
+
   delete rawStreamFast;
   return res;
 }
@@ -191,9 +214,15 @@ Bool_t AliTPCdataQA::ProcessEvent(AliRawReader *rawReader)
   //
 
   // if fMapping is NULL the rawstream will crate its own mapping
+  fSectorLast  = -1;
   AliTPCRawStream rawStream(rawReader, (AliAltroMapping**)fMapping);
   rawReader->Select("TPC");
-  return ProcessEvent(&rawStream);
+  Bool_t res =  ProcessEvent(&rawStream);
+
+  if(res)
+    fEventCounter++; // only increment event counter if there is TPC data
+                     // otherwise Analyse (called in QA) fails
+  return res;
 }
 
 
@@ -253,39 +282,222 @@ Int_t AliTPCdataQA::Update(const Int_t icsector, /*FOLD00*/
   //
   if (icTimeBin<fFirstTimeBin) return 0;
   if (icTimeBin>fLastTimeBin) return 0;
+
   if (!fMaxCharge) fMaxCharge = new AliTPCCalPad("MaxCharge","MaxCharge");
+  if (!fMeanCharge) fMeanCharge = new AliTPCCalPad("MeanCharge","MeanCharge");
+  if (!fNoThreshold) fNoThreshold = new AliTPCCalPad("NoThreshold","NoThreshold");
   if (!fOverThreshold0) fOverThreshold0 = new AliTPCCalPad("OverThreshold0","OverThreshold0");
   if (!fOverThreshold5) fOverThreshold5 = new AliTPCCalPad("OverThreshold5","OverThreshold5");
   if (!fOverThreshold10) fOverThreshold10 = new AliTPCCalPad("OverThreshold10","OverThreshold10");
   if (!fOverThreshold20) fOverThreshold20 = new AliTPCCalPad("OverThreshold20","OverThreshold20");
   if (!fOverThreshold30) fOverThreshold30 = new AliTPCCalPad("OverThreshold30","OverThreshold30");
   //
-  if (csignal>fMaxCharge->GetCalROC(icsector)->GetValue(icRow, icPad)){
-    fMaxCharge->GetCalROC(icsector)->SetValue(icRow, icPad,csignal);
+
+  Int_t signal = Int_t(csignal);
+
+  // if pedestal calibrations are loaded subtract pedestals
+  if(fPedestal) {
+
+    Int_t pedestal = Int_t(fPedestal->GetCalROC(icsector)->GetValue(icRow, icPad));
+    if(pedestal<10 || pedestal>90)
+      return 0;
+    signal -= pedestal;
+  }
+
+
+  if (signal >= 0) {
+    
+    Float_t count = fNoThreshold->GetCalROC(icsector)->GetValue(icRow, icPad);
+    fNoThreshold->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
   }
+
+  // Require at least 3 ADC channels
+  if (signal < 3)
+    return 0;
+
+  // if noise calibrations are loaded require at least 3*sigmaNoise
+  if(fNoise) {
   
+    Float_t noise = fNoise->GetCalROC(icsector)->GetValue(icRow, icPad);
+
+    if(signal<noise*3)
+      return 0;
+  }
   //
-  if (csignal>0){
-    Int_t count = fOverThreshold0->GetCalROC(icsector)->GetValue(icRow, icPad);
+  // this signal is ok - now see if the previous signal was connected
+  // this is a bit ugly since the standard decoder goes down in time bins
+  // (10, 9, 8..) while the fast HLT decoder goes up in time bins (1, 2, 3..) 
+  //
+  if(fSectorLast==icsector && fRowLast==icRow && fPadLast==icPad &&
+     fTimeBinLast==icTimeBin+1 || fTimeBinLast==icTimeBin-1)
+    fNAboveThreshold++;
+  else
+    fNAboveThreshold = 1;
+    
+  if(fNAboveThreshold==2) {
+    
+    //
+    // This is the only special case, because before we did not know if we
+    // should store the information
+    //
+    UpdateSignalHistograms(fSectorLast, fRowLast, fPadLast, fTimeBinLast,
+                        fSignalLast);
+  }
+  
+  // keep the information for the next signal
+  fSectorLast  = icsector;
+  fRowLast     = icRow;
+  fPadLast     = icPad;
+  fTimeBinLast = icTimeBin;
+  fSignalLast  = signal;
+  
+  if(fNAboveThreshold==1) // we don't know if this information should be stored
+    return 1;
+  
+  UpdateSignalHistograms(fSectorLast, fRowLast, fPadLast, fTimeBinLast,
+                      fSignalLast);
+
+  return 1;
+}
+//_____________________________________________________________________
+void AliTPCdataQA::UpdateSignalHistograms(const Int_t icsector, /*FOLD00*/
+                                       const Int_t icRow,
+                                       const Int_t icPad,
+                                       const Int_t icTimeBin,
+                                       const Float_t signal)
+{
+  //
+  // Signal filling method
+  //
+  
+  {
+    Float_t charge = fMeanCharge->GetCalROC(icsector)->GetValue(icRow, icPad);
+    fMeanCharge->GetCalROC(icsector)->SetValue(icRow, icPad, charge + signal);
+  }
+  
+  if (signal>fMaxCharge->GetCalROC(icsector)->GetValue(icRow, icPad)){
+    fMaxCharge->GetCalROC(icsector)->SetValue(icRow, icPad,signal);
+  }
+  
+  if (signal>0){
+    Float_t count = fOverThreshold0->GetCalROC(icsector)->GetValue(icRow, icPad);
     fOverThreshold0->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
   };
   //
-  if (csignal>5){
-    Int_t count = fOverThreshold5->GetCalROC(icsector)->GetValue(icRow, icPad);
+  if (signal>5){
+    Float_t count = fOverThreshold5->GetCalROC(icsector)->GetValue(icRow, icPad);
     fOverThreshold5->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
   };
-  if (csignal>10){
-    Int_t count = fOverThreshold10->GetCalROC(icsector)->GetValue(icRow, icPad);
+  if (signal>10){
+    Float_t count = fOverThreshold10->GetCalROC(icsector)->GetValue(icRow, icPad);
     fOverThreshold10->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
   };
-  if (csignal>20){
-    Int_t count = fOverThreshold20->GetCalROC(icsector)->GetValue(icRow, icPad);
+  if (signal>20){
+    Float_t count = fOverThreshold20->GetCalROC(icsector)->GetValue(icRow, icPad);
     fOverThreshold20->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
   };
-  if (csignal>30){
-    Int_t count = fOverThreshold30->GetCalROC(icsector)->GetValue(icRow, icPad);
+  if (signal>30){
+    Float_t count = fOverThreshold30->GetCalROC(icsector)->GetValue(icRow, icPad);
     fOverThreshold30->GetCalROC(icsector)->SetValue(icRow, icPad,count+1);
-  };
+  };  
+}
+
+//_____________________________________________________________________
+void AliTPCdataQA::Analyse()
+{
+  //
+  //  Calculate calibration constants
+  //
+  
+  cout << "Analyse called" << endl;
+
+  if(fEventCounter==0) {
+
+    cout << "EventCounter == 0, Cannot analyse" << endl;
+    return;
+  }
+
+  Int_t nTimeBins = fLastTimeBin - fFirstTimeBin +1;
+  cout << "EventCounter: " << fEventCounter << endl
+       << "TimeBins: " << nTimeBins << endl;
+
+  if (fMeanCharge && fNoThreshold) fMeanCharge->Divide(fNoThreshold);
+
+  Float_t normalization = 1.0 / Float_t(fEventCounter * nTimeBins);
+  if (fNoThreshold)     fNoThreshold->Multiply(normalization);  
+  if (fOverThreshold0)  fOverThreshold0->Multiply(normalization);  
+  if (fOverThreshold5)  fOverThreshold5->Multiply(normalization);  
+  if (fOverThreshold10) fOverThreshold10->Multiply(normalization);  
+  if (fOverThreshold20) fOverThreshold20->Multiply(normalization);  
+  if (fOverThreshold30) fOverThreshold30->Multiply(normalization);  
+}
+
+
+void AliTPCdataQA::MakeTree(const char *fname){
+  //
+  // Export result to the tree -located in the file
+  // This file can be analyzed using AliTPCCalibViewer
+  // 
+  AliTPCdataQA *ped = this;
+  AliTPCPreprocessorOnline preprocesor;
+  if (ped->GetMaxCharge()) preprocesor.AddComponent(ped->GetMaxCharge());  
+  if (ped->GetMeanCharge()) preprocesor.AddComponent(ped->GetMeanCharge());  
+  if (ped->GetOverThreshold0()) preprocesor.AddComponent(ped->GetOverThreshold0());
+  if (ped->GetOverThreshold5()) preprocesor.AddComponent(ped->GetOverThreshold5());
+  if (ped->GetOverThreshold10()) preprocesor.AddComponent(ped->GetOverThreshold10());
+  if (ped->GetOverThreshold20()) preprocesor.AddComponent(ped->GetOverThreshold20());
+  if (ped->GetOverThreshold30()) preprocesor.AddComponent(ped->GetOverThreshold30());
+  preprocesor.DumpToFile(fname);  
+}
+
+
 
-  return 0;
+void AliTPCdataQA::MakeArrays(){
+  //
+  //
+  //
+  AliTPCROC * roc = AliTPCROC::Instance();
+  //
+  Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
+  Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1);
+  fAllBins = new Float_t*[nRowsMax];
+  fAllSigBins = new Int_t*[nRowsMax];
+  fAllNSigBins = new Int_t[nRowsMax];
+  for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
+    //
+    Int_t maxBin = fMaxTime*(nPadsMax+6);  // add 3 virtual pads  before and 3 after
+    fAllBins[iRow] = new Float_t[maxBin];
+    memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin);
+    fAllSigBins[iRow] = new Int_t[maxBin];
+    fAllNSigBins[iRow]=0;
+  }
+}
+
+
+void AliTPCdataQA::CleanArrays(){
+  //
+  //
+  //
+  AliTPCROC * roc = AliTPCROC::Instance();
+  //
+  Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
+  Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1); 
+  for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
+    //
+    Int_t maxBin = fMaxTime*(nPadsMax+6);  // add 3 virtual pads  before and 3 after
+    memset(fAllBins[iRow],0,sizeof(Float_t)*maxBin);
+    fAllNSigBins[iRow]=0;
+  }
+}
+
+Float_t* AliTPCdataQA::GetExpandDigit(Int_t row, Int_t pad, Int_t time){
+  //
+  //
+  //
+  AliTPCROC * roc = AliTPCROC::Instance();
+  Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
+  if (row<0 || row>nRowsMax) return 0;
+  Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1); 
+  
 }