]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCCalibKr.cxx
Updated flags for low flux case (A. Dainese)
[u/mrichter/AliRoot.git] / TPC / AliTPCCalibKr.cxx
index 07ab8473fa6fc3016953df7b46f053411ed247a0..650b33026cea9c4f03960c6e012fe0bc6724bafd 100644 (file)
 // The AliTPCCalibKr class description (TPC Kr calibration).
 //
 //
-// The AliTPCCalibKr fills the array of TH3F histograms (TPC_max_padraw,TPC_max_pad,TPC_ADC_cluster),
-// its data memebers.   
+// The AliTPCCalibKr keeps the array of TH3F histograms (TPC_max_padraw,TPC_max_pad,TPC_ADC_cluster),
+// its data memebers and is filled by AliTPCCalibKrTask under conditions (Accept()).   
 // 
-// As the input it requires the tree with reconstructed Kr clusters (AliTPCclusterKr objects). 
-// The AliTPCCalibKr objects containing an array of TH3F histograms are stored (by default) in the 
-// ouptut (outHistFile.root) file.
-//
 // The ouput TH3F histograms are later used to determine the calibration parameters of TPC chambers.
 // These calculations are done by using AliTPCCalibKr::Analyse() function. The ouput calibration 
 // parameters (details in AliTPCCalibKr::Analyse()) are stored in the calibKr.root file for each TPC pad.
 // In addition the debugCalibKr.root file with debug information is created.
 //
+
+/*
 // Usage example:
 //
-// 1. Create outHistFile.root histogram file:
-//
-// -- Load libXrdClient.so if data on Xrd cluster e.g. (GSI)
-// gSystem->Load("/usr/local/grid/XRootd/GSI/lib64/libXrdClient.so");
-//
-// -- Load toolkit
-// gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
-// gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+");
-// AliXRDPROOFtoolkit tool;
-//
-// -- Make chain of files
-// TChain * chain = tool.MakeChain("KrClusters.txt","Kr","",1000,0);
-//
-// -- Run AliTPCCalibKr task (Only TPC C side)
-// AliTPCCalibKr *task = new AliTPCCalibKr;
-// task->SetInputChain(chain);
-// task->SetASide(kFALSE);
-//
-// task->Process();
-// 
-// 2. Analyse output histograms:
-//
-// TFile f("outHistFile.root");
-// AliTPCCalibKr.Analyse();
-//
-// 3. See calibration parameters e.g.:
-//
-// TFile f("calibKr.root");
-// spectrMean->GetCalROC(70)->GetValue(40,40);
-// fitMean->GetCalROC(70)->GetValue(40,40);
-//
-// 4. See debug information e.g.:
-//
-// TFile f("debugCalibKr.root");
-// .ls;
-//
+
+// 1. Analyse output histograms:
+TFile f("outHistFile.root");
+AliTPCCalibKr *obj = (AliTPCCalibKr*) cOutput.FindObject("AliTPCCalibKr")
+obj->Analyse();
+
+// 2. See calibration parameters e.g.:
+TFile f("calibKr.root");
+spectrMean->GetCalROC(70)->GetValue(40,40);
+fitMean->GetCalROC(70)->GetValue(40,40);
+
+// 3. See debug information e.g.:
+TFile f("debugCalibKr.root");
+.ls;
+
 // -- Print calibKr TTree content 
-// calibKr->Print();
-//
+calibKr->Print();
+
 // -- Draw calibKr TTree variables
-// calibKr.Draw("fitMean");
-//
+calibKr.Draw("fitMean");
+
+*/
+
+
 //
 // Author: Jacek Otwinowski (J.Otwinowski@gsi.de) and Stafan Geartner (S.Gaertner@gsi.de)
 //-----------------------------------------------------------------------------
@@ -112,34 +93,58 @@ ClassImp(AliTPCCalibKr)
 
 AliTPCCalibKr::AliTPCCalibKr() : 
   TObject(),
-  
-  bOutputHisto(kTRUE),
-  bASide(kTRUE),
-  bCSide(kTRUE),
-  fClusters(0),
-  fClustKr(0),
-  fTree(0),
-  fHistoKrArray(72)
+  fASide(kTRUE),
+  fCSide(kTRUE),
+  fHistoKrArray(72),
+  fADCOverClustSizeMin(0.0),
+  fADCOverClustSizeMax(1.0e9),
+  fMaxADCOverClustADCMin(0.0),
+  fMaxADCOverClustADCMax(1.0e9),
+  fTimeMin(0.0),
+  fTimeMax(1.0e9),
+  fClustSizeMin(0.0),
+  fClustSizeMax(1.0e9)
 {
   //
   // default constructor
   //
+
+  // TObjArray with histograms
+  fHistoKrArray.SetOwner(kTRUE); // is owner of histograms
+  fHistoKrArray.Clear();
+  // init cuts
+  SetADCOverClustSizeRange(7,200);
+  SetMaxADCOverClustADCRange(0.01,0.4);
+  SetTimeRange(200,600);
+  SetClustSizeRange(6,200);
+  // init histograms
+  Init();
 }
 
 //_____________________________________________________________________
 AliTPCCalibKr::AliTPCCalibKr(const AliTPCCalibKr& pad) : 
   TObject(pad),
   
-  bOutputHisto(pad.bOutputHisto),
-  bASide(pad.bASide),
-  bCSide(pad.bCSide),
-  fClusters(pad.fClusters),
-  fClustKr(pad.fClustKr),
-  fTree(pad.fTree),
-  fHistoKrArray(72)
+  fASide(pad.fASide),
+  fCSide(pad.fCSide),
+  fHistoKrArray(72),
+  fADCOverClustSizeMin(pad.fADCOverClustSizeMin),
+  fADCOverClustSizeMax(pad.fADCOverClustSizeMax),
+  fMaxADCOverClustADCMin(pad.fMaxADCOverClustADCMin),
+  fMaxADCOverClustADCMax(pad.fMaxADCOverClustADCMax),
+  fTimeMin(pad.fTimeMin),
+  fTimeMax(pad.fTimeMax),
+  fClustSizeMin(pad.fClustSizeMin),
+  fClustSizeMax(pad.fClustSizeMax)
 {
   // copy constructor
  
+  // TObjArray with histograms
+  fHistoKrArray.SetOwner(kTRUE); // is owner of histograms
+  fHistoKrArray.Clear();
+
   for (Int_t iSec = 0; iSec < 72; ++iSec) 
   {
     TH3F *hOld = pad.GetHistoKr(iSec);
@@ -156,9 +161,10 @@ AliTPCCalibKr::~AliTPCCalibKr()
   //
   // destructor
   //
-  if(fClustKr)  delete fClustKr; fClustKr = 0;
-  if(fClusters) delete fClusters; fClusters = 0;
-  if(fTree)     delete fTree; fTree = 0;
+
+  // for (Int_t iSec = 0; iSec < 72; ++iSec) {
+  //     if (fHistoKrArray.At(iSec)) delete fHistoKrArray.RemoveAt(iSec);
+  // }
   fHistoKrArray.Delete();
 }
 
@@ -177,105 +183,38 @@ AliTPCCalibKr& AliTPCCalibKr::operator = (const  AliTPCCalibKr &source)
 void AliTPCCalibKr::Init()
 {
   // 
-  // init input tree and output histograms 
+  // init output histograms 
   //
 
-  // set input tree
-  if(!fTree) { 
-   Printf("ERROR: Could not read chain from input");
-  }
-  else {
-   fTree->SetBranchStatus("*",1); 
-  }
-
-  // set branch address
-  fClusters = new TClonesArray("AliTPCclusterKr");
-
-  if(!fTree->GetBranch("fClusters")) {
-    Printf("ERROR: Could not get fClusters branch from input");
-  } else {
-   fTree->GetBranch("fClusters")->SetAddress(&fClusters);
-  }
-  
-  // create output TObjArray
-  fHistoKrArray.Clear();
-
   // add histograms to the TObjArray
   for(Int_t i=0; i<72; ++i) {
     
        // C - side
-       if( IsCSide(i) == kTRUE && bCSide == kTRUE) {
+       if(IsCSide(i) == kTRUE && fCSide == kTRUE) {
       TH3F *hist = CreateHisto(i);
       if(hist) fHistoKrArray.AddAt(hist,i);
        }
     
        // A - side
-       if(IsCSide(i) == kFALSE && bASide == kTRUE) {
+       if(IsCSide(i) == kFALSE && fASide == kTRUE) {
       TH3F *hist = CreateHisto(i);
       if(hist) fHistoKrArray.AddAt(hist,i);
        }
-
-  }
-}
-
-//_____________________________________________________________________
-Bool_t AliTPCCalibKr::ReadEntry(Int_t evt)
-{
-  // 
-  // read entry from the tree
-  //
-  Long64_t centry = fTree->LoadTree(evt);
-  if(centry < 0) return kFALSE;
-
-  if(!fTree->GetBranch("fClusters")) 
-  {
-    Printf("ERROR: Could not get fClusters branch from input");
-       return kFALSE;
-  } else {
-   fTree->GetBranch("fClusters")->SetAddress(&fClusters);
   }
-
-  fTree->GetEntry(evt);
-
-return kTRUE;
 }
  
 //_____________________________________________________________________
-Bool_t AliTPCCalibKr::Process()
+Bool_t AliTPCCalibKr::Process(AliTPCclusterKr *cluster)
 {
   //
   // process events 
   // call event by event
   //
 
-  // init tree
-  Init();
-
-  // get events
-  if(!fTree) return kFALSE;
-  Int_t nEvents = fTree->GetEntries();
-
-  // fill histograms 
-  for(Int_t i=0; i<nEvents; ++i)
-  {
-    if(ReadEntry(i) == kFALSE) return kFALSE;
-
-    if(!(i%10000)) cout << "evt: " << i << endl; 
-
-    // get TClonesArray entries
-    fClustKr = 0;
-    Int_t entries = fClusters->GetEntries();
-    for(Int_t j=0; j < entries; ++j)
-       {
-         fClustKr = (AliTPCclusterKr*)fClusters->At(j);
-
-      if(fClustKr) Update(fClustKr);
-         else return kFALSE;
-       }
-  }
-
-  // write output 
-  return Terminate();
+  if(cluster) Update(cluster);
+  else return kFALSE;
+ return kTRUE;
 }
 
 //_____________________________________________________________________
@@ -291,9 +230,9 @@ TH3F* AliTPCCalibKr::CreateHisto(Int_t chamber)
 
     if( IsIROC(chamber) == kTRUE ) 
        {
-          h = new TH3F(name,name,63,0,63,100,0,100,150,100,3000);
+          h = new TH3F(name,name,63,0,63,108,0,108,200,100,2500);
        } else {
-          h = new TH3F(name,name,96,0,96,100,0,100,150,100,3000);
+          h = new TH3F(name,name,96,0,96,140,0,140,200,100,1700);
        }
     h->SetXTitle("padrow");
     h->SetYTitle("pad");
@@ -340,6 +279,7 @@ Bool_t AliTPCCalibKr::Update(AliTPCclusterKr  *cl)
   return kTRUE;
 }
 
+//_____________________________________________________________________
 Bool_t AliTPCCalibKr::Accept(AliTPCclusterKr  *cl){
   //
   // cuts
@@ -349,25 +289,53 @@ Bool_t AliTPCCalibKr::Accept(AliTPCclusterKr  *cl){
     TCut cutR1("cutR1","fADCcluster/fSize>7");          // cosmic tracks and noise removal
     TCut cutR2("cutR2","fMax.fAdc/fADCcluster<0.4");    // digital noise removal
     TCut cutR3("cutR3","fMax.fAdc/fADCcluster>0.01");   // noise removal
+    TCut cutR4("cutR4","fMax.fTime>200");   // noise removal
+    TCut cutR5("cutR5","fMax.fTime<600");   // noise removal
     TCut cutS1("cutS1","fSize<200");    // adjust it according v seetings - cosmic tracks
-    TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutS1;
+    TCut cutAll = cutR0+cutR1+cutR2+cutR3+cutR4+cutR5+cutS1;
   */
+  /*
   //R0
-  if (cl->GetADCcluster()/ cl->GetSize() >200)        return kFALSE;
+  if ((float)cl->GetADCcluster()/ cl->GetSize() >200)        return kFALSE;
   // R1
-  if (cl->GetADCcluster()/ cl->GetSize() <7)          return kFALSE;
+  if ((float)cl->GetADCcluster()/ cl->GetSize() <7)          return kFALSE;
   //R2
-  if (cl->GetMax().GetAdc()/ cl->GetADCcluster() >0.4)  return kFALSE;
+  if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() >0.4)  return kFALSE;
   //R3
-  if (cl->GetMax().GetAdc()/ cl->GetADCcluster() <0.01) return kFALSE;
+  if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() <0.01) return kFALSE;
+  //R4
+  if (cl->GetMax().GetTime() < 200) return kFALSE;
+  //R5
+  if (cl->GetMax().GetTime() > 600) return kFALSE;
   //S1
   if (cl->GetSize()>200) return kFALSE;
   if (cl->GetSize()<6)  return kFALSE;
-  return kTRUE;
 
-}
+  SetADCOverClustSizeRange(7,200);
+  SetMaxADCOverClustADCRange(0.01,0.4);
+  SetTimeRange(200,600);
+  SetClustSizeRange(6,200);
+  */
+
+  //R0
+  if ((float)cl->GetADCcluster()/ cl->GetSize() >fADCOverClustSizeMax)        return kFALSE;
+  // R1
+  if ((float)cl->GetADCcluster()/ cl->GetSize() <fADCOverClustSizeMin)          return kFALSE;
+  //R2
+  if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() >fMaxADCOverClustADCMax)  return kFALSE;
+  //R3
+  if ((float)cl->GetMax().GetAdc()/ cl->GetADCcluster() <fMaxADCOverClustADCMin) return kFALSE;
+  //R4
+  if (cl->GetMax().GetTime() > fTimeMax) return kFALSE;
+  //R5
+  if (cl->GetMax().GetTime() < fTimeMin) return kFALSE;
+  //S1
+  if (cl->GetSize()>fClustSizeMax) return kFALSE;
+  if (cl->GetSize()<fClustSizeMin)  return kFALSE;
 
+  return kTRUE;
 
+}
 
 //_____________________________________________________________________
 TH3F* AliTPCCalibKr::GetHistoKr(Int_t chamber) const
@@ -375,38 +343,6 @@ TH3F* AliTPCCalibKr::GetHistoKr(Int_t chamber) const
   // get histograms from fHistoKrArray
   return (TH3F*) fHistoKrArray.At(chamber);
 }
-
-//_____________________________________________________________________
-Bool_t AliTPCCalibKr::Terminate() 
-{
-  //
-  // store AliTPCCalibKr in the output file 
-  //
-  if(bOutputHisto) {
-    TFile *outFile = new TFile("outHistFile.root","RECREATE"); 
-   
-    if(outFile) 
-       {
-         outFile->cd();
-
-         for(int i=0; i<72; ++i) {
-            if( IsCSide(i) == kTRUE && bCSide == kTRUE)
-              printf("C side chamber: %d, 3D histo entries: %10.f \n",i,((TH3F*)fHistoKrArray.At(i))->GetEntries());
-
-            if( IsCSide(i) == kFALSE && bASide == kTRUE)
-              printf("A side chamber: %d, 3D histo entries: %10.f \n",i,((TH3F*)fHistoKrArray.At(i))->GetEntries());
-         }
-         this->Write();
-         outFile->Close();
-
-         return kTRUE;
-       }
-       else 
-         return kFALSE;
-  }
-
-return kFALSE;
-}
  
 //_____________________________________________________________________
 void AliTPCCalibKr::Analyse() 
@@ -432,8 +368,9 @@ void AliTPCCalibKr::Analyse()
   Double_t windowFrac = 0.12;
   // the 3d histogram will be projected on the pads given by the following window size
   // set the numbers to 0 if you want to do a pad-by-pad calibration
-  UInt_t rowRadius = 5;
-  UInt_t padRadius = 10;
+  UInt_t rowRadius = 2;
+  UInt_t padRadius = 4;
+  
   // the step size by which pad and row are incremented is given by the following numbers
   // set them to 1 if you want the finest granularity
   UInt_t rowStep = 1;    // formerly: 2*rowRadius
@@ -481,9 +418,16 @@ void AliTPCCalibKr::Analyse()
     
         // find maximum in spectrum to define a range (given by windowFrac) for which a Gauss is fitted
         Double_t maxEntries = projH->GetBinCenter(projH->GetMaximumBin());
+               Int_t minBin = projH->FindBin((1.-windowFrac) * maxEntries);
+               Int_t maxBin = projH->FindBin((1.+windowFrac) * maxEntries);
+        Double_t integCharge =  projH->Integral(minBin,maxBin); 
+
         Int_t fitResult = projH->Fit("gaus", "Q0", "", (1.-windowFrac) * maxEntries, (1.+windowFrac) * maxEntries);
+
         if (fitResult != 0) {
-          Error("Analyse", "Error while fitting spectrum for chamber %i, rows %i - %i, pads %i - %i.", chamber, rowLow, rowUp, padLow, padUp);
+          Error("Analyse", "Error while fitting spectrum for chamber %i, rows %i - %i, pads %i - %i, integrated charge %f.", chamber, rowLow, rowUp, padLow, padUp, integCharge);
+          //printf("[ch%i r%i, p%i] entries = %f, maxEntries = %f,  %f\n", chamber, iRow, iPad, entries, maxEntries);
+    
           delete projH;
           continue;
         }
@@ -498,7 +442,7 @@ void AliTPCCalibKr::Analyse()
         delete gausFit;
         delete projH;
         if (fitMean <= 0) continue;
-        printf("[ch%i r%i, p%i] entries = %f, maxEntries = %f, fitMean = %f, fitRMS = %f\n", chamber, iRow, iPad, entries, maxEntries, fitMean, fitRMS);
+        //printf("[ch%i r%i, p%i] entries = %f, maxEntries = %f, fitMean = %f, fitRMS = %f\n", chamber, iRow, iPad, entries, maxEntries, fitMean, fitRMS);
     
         // write the calibration parameters for each pad that the 3d histogram was projected onto
         // (with considering the step size) to the CalPads
@@ -506,9 +450,9 @@ void AliTPCCalibKr::Analyse()
         // rowStep (padStep) even: fill s/2 rows (pads) in ascending direction, s/2-1 in descending direction
         for (Int_t r = iRow - (rowStep/2 - (rowStep+1)%2); r <= (Int_t)(iRow + rowStep/2); r++) {
           if (r < 0 || r >= (Int_t)nRows) continue;
-          UInt_t nPads = roc.GetNPads(r);
+          UInt_t nPadsR = roc.GetNPads(r);
           for (Int_t p = iPad - (padStep/2 - (padStep+1)%2); p <= (Int_t)(iPad + padStep/2); p++) {
-            if (p < 0 || p >= (Int_t)nPads) continue;
+            if (p < 0 || p >= (Int_t)nPadsR) continue;
             spectrMeanCalPad->GetCalROC(chamber)->SetValue(r, p, histMean);
             spectrRMSCalPad->GetCalROC(chamber)->SetValue(r, p, histRMS);
             fitMeanCalPad->GetCalROC(chamber)->SetValue(r, p, fitMean);
@@ -582,3 +526,41 @@ TH1D* AliTPCCalibKr::ProjectHisto(TH3F* histo3D, const char* name, Int_t xMin, I
   projH->SetEntries((Long64_t)entries);
   return projH;
 }
+
+//_____________________________________________________________________
+Long64_t AliTPCCalibKr::Merge(TCollection* list) {
+// merge function 
+//
+cout << "Merge " << endl;
+
+if (!list)
+return 0;
+
+if (list->IsEmpty())
+return 1;
+
+TIterator* iter = list->MakeIterator();
+TObject* obj = 0;
+
+    iter->Reset();
+    Int_t count=0;
+       while((obj = iter->Next()) != 0)
+       {
+         AliTPCCalibKr* entry = dynamic_cast<AliTPCCalibKr*>(obj);
+         if (entry == 0) continue; 
+
+               for(int i=0; i<72; ++i) { 
+                 if(IsCSide(i) == kTRUE && fCSide == kTRUE) { 
+                 ((TH3F*)fHistoKrArray.At(i))->Add( ((TH3F*)entry->fHistoKrArray.At(i)) );  
+                 } 
+
+                 if(IsCSide(i) == kFALSE && fASide == kTRUE) { 
+                   ((TH3F*)fHistoKrArray.At(i))->Add( ((TH3F*)entry->fHistoKrArray.At(i)) );  
+                 } 
+               } 
+
+         count++;
+       }
+
+return count;
+}