AliTPCClusterParam.h -> simple Getter for Matrix
authormarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 15 Feb 2011 21:15:07 +0000 (21:15 +0000)
committermarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 15 Feb 2011 21:15:07 +0000 (21:15 +0000)
AliTPCRecoParam.h           -> simple Setters for ctgRange, Qmax/Qtot, UseMultiplicityCorrectionDeDx
AliTPCRecoParam.cxx         -> the new variables must also be added...

AliTPCPreprocessorOffline.h   -> new Functions pad region equalization in pass 0
AliTPCPreprocessorOffline.cxx -> Adding extraction of pad region equalization and put it to pass 0

TPC/AliTPCClusterParam.h
TPC/AliTPCPreprocessorOffline.cxx
TPC/AliTPCPreprocessorOffline.h
TPC/AliTPCRecoParam.cxx
TPC/AliTPCRecoParam.h

index f5f7a70..63c4566 100644 (file)
@@ -37,7 +37,8 @@ class AliTPCClusterParam : public TObject {
   void FitRMS(TTree * tree);
   void SetQnorm(Int_t ipad, Int_t itype,  const TVectorD *const norm); 
   void SetQnormCorr(Int_t ipad, Int_t itype, Int_t corrType, Float_t val); 
-  Double_t  GetQnormCorr(Int_t ipad, Int_t itype, Int_t corrType) const; 
+  Double_t  GetQnormCorr(Int_t ipad, Int_t itype, Int_t corrType) const;
+  TMatrixD *GetQnormCorrMatrix(){return fQNormCorr;};
   void ResetQnormCorr(); 
   //
   // Charge parameterization
index 3e8e7ae..b47cf5f 100644 (file)
@@ -72,6 +72,7 @@
 #include "AliRelAlignerKalman.h"
 #include "AliTPCParamSR.h"
 #include "AliTPCcalibTimeGain.h"
+#include "AliTPCcalibGainMult.h"
 #include "AliSplineFit.h"
 #include "AliTPCPreprocessorOffline.h"
 
@@ -96,6 +97,7 @@ AliTPCPreprocessorOffline::AliTPCPreprocessorOffline():
   fGainArray(new TObjArray),               // array to be stored in the OCDB
   fGainMIP(0),          // calibration component for MIP
   fGainCosmic(0),       // calibration component for cosmic
+  fGainMult(0),
   fSwitchOnValidation(kFALSE) // flag to switch on validation of OCDB parameters
 {
   //
@@ -417,16 +419,17 @@ void AliTPCPreprocessorOffline::AddHistoGraphs(  TObjArray * vdriftArray, AliTPC
       if (!graph) {
        printf("Graph =%s filtered out\n", name.Data());
        continue;
-      }      
+      }
+      //
       if (graph){
-       graph->SetMarkerStyle(i%8+20);
-       graph->SetMarkerColor(i%7);
-       graph->GetXaxis()->SetTitle("Time");
-       graph->GetYaxis()->SetTitle("v_{dcor}");
-       graph->SetName(graphName);
-       graph->SetTitle(graphName);
-       printf("Graph %d\t=\t%s\n", i, graphName.Data());
-       vdriftArray->Add(graph);
+        graph->SetMarkerStyle(i%8+20);
+        graph->SetMarkerColor(i%7);
+        graph->GetXaxis()->SetTitle("Time");
+        graph->GetYaxis()->SetTitle("v_{dcor}");
+        graph->SetName(graphName);
+        graph->SetTitle(graphName);
+        printf("Graph %d\t=\t%s\n", i, graphName.Data());
+        vdriftArray->Add(graph);
       }
     }
   }
@@ -792,7 +795,7 @@ void AliTPCPreprocessorOffline::CalibTimeGain(const Char_t* fileName, Int_t star
   //
   AnalyzeGain(startRunNumber,endRunNumber, 1000,1.43);
   AnalyzeAttachment(startRunNumber,endRunNumber);
-
+  AnalyzePadRegionGain();
   //
   // 3. Make control plots
   //
@@ -821,9 +824,11 @@ void AliTPCPreprocessorOffline::ReadGainGlobal(const Char_t* fileName){
   if (array){
     fGainMIP    = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGain");
     fGainCosmic = ( AliTPCcalibTimeGain *)array->FindObject("calibTimeGainCosmic");
+    fGainMult   = ( AliTPCcalibGainMult *)array->FindObject("calibGainMult");
   }else{
     fGainMIP    = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGain");
     fGainCosmic = ( AliTPCcalibTimeGain *)fcalib.Get("calibTimeGainCosmic");
+    fGainMult   = ( AliTPCcalibGainMult *)fcalib.Get("calibGainMult");
   }
   TH1 * hisT=0;
   Int_t firstBinA =0, lastBinA=0;
@@ -973,6 +978,48 @@ Bool_t AliTPCPreprocessorOffline::AnalyzeAttachment(Int_t startRunNumber, Int_t
 }
 
 
+Bool_t AliTPCPreprocessorOffline::AnalyzePadRegionGain(){
+  //
+  // Analyze gain for different pad regions - produce the calibration graphs 0,1,2
+  //
+  if (fGainMult) 
+  {
+    TH2D * histQmax = (TH2D*) fGainMult->GetHistPadEqual()->Projection(0,2);
+    TH2D * histQtot = (TH2D*) fGainMult->GetHistPadEqual()->Projection(1,2);
+    //
+    TObjArray arr;
+    histQmax->FitSlicesY(0,0,-1,0,"QNR",&arr);
+    Double_t xMax[3] = {0,1,2};
+    Double_t yMax[3]    = {((TH1D*)arr.At(1))->GetBinContent(1),
+                          ((TH1D*)arr.At(1))->GetBinContent(2),
+                          ((TH1D*)arr.At(1))->GetBinContent(3)};
+    Double_t yMaxErr[3] = {((TH1D*)arr.At(1))->GetBinError(1),
+                          ((TH1D*)arr.At(1))->GetBinError(2),
+                          ((TH1D*)arr.At(1))->GetBinError(3)};
+    TGraphErrors * fitPadRegionQmax = new TGraphErrors(3, xMax, yMax, 0, yMaxErr);
+    //
+    histQtot->FitSlicesY(0,0,-1,0,"QNR",&arr);
+    Double_t xTot[3] = {0,1,2};
+    Double_t yTot[3]    = {((TH1D*)arr.At(1))->GetBinContent(1),
+                          ((TH1D*)arr.At(1))->GetBinContent(2),
+                          ((TH1D*)arr.At(1))->GetBinContent(3)};
+    Double_t yTotErr[3] = {((TH1D*)arr.At(1))->GetBinError(1),
+                          ((TH1D*)arr.At(1))->GetBinError(2),
+                          ((TH1D*)arr.At(1))->GetBinError(3)};
+    TGraphErrors * fitPadRegionQtot = new TGraphErrors(3, xTot, yTot, 0, yTotErr);
+    //
+    fitPadRegionQtot->SetName("TGRAPHERRORS_MEANQTOT_PADREGIONGAIN_BEAM_ALL");// set proper names according to naming convention
+    fitPadRegionQmax->SetName("TGRAPHERRORS_MEANQMAX_PADREGIONGAIN_BEAM_ALL");// set proper names according to naming convention
+    //
+    fGainArray->AddLast(fitPadRegionQtot);
+    fGainArray->AddLast(fitPadRegionQmax);
+    return kTRUE;
+  } 
+  return kFALSE;
+
+}
+
+
 void AliTPCPreprocessorOffline::UpdateOCDBGain(Int_t startRunNumber, Int_t endRunNumber, const Char_t *storagePath){
   //
   // Update OCDB entry
@@ -1014,14 +1061,13 @@ void AliTPCPreprocessorOffline::MakeQAPlot(Float_t  FPtoMIPratio) {
       for(Int_t i=0; i < fGraphCosmic->GetN(); i++) {
        fGraphCosmic->GetY()[i] *= FPtoMIPratio;        
       }
-    
-      fGraphCosmic->Draw("lp");
-      grfFitCosmic->SetLineColor(2);
-      grfFitCosmic->Draw("lu");
-      fGainArray->AddLast(gainHistoCosmic);
-      //fGainArray->AddLast(canvasCosmic->Clone());
-      delete canvasCosmic;    
     }
+    fGraphCosmic->Draw("lp");
+    grfFitCosmic->SetLineColor(2);
+    grfFitCosmic->Draw("lu");
+    fGainArray->AddLast(gainHistoCosmic);
+    //fGainArray->AddLast(canvasCosmic->Clone());
+    delete canvasCosmic;    
   }
   if (fFitMIP) {
     TCanvas * canvasMIP = new TCanvas("gain MIP", "time dependent gain QA histogram MIP");
index 1676c21..1e51b7e 100644 (file)
@@ -14,6 +14,7 @@
 class TObjArray;
 class AliTPCcalibTime;
 class AliTPCcalibTimeGain;
+class AliTPCcalibGainMult;
 class AliTPCROC;
 class AliTPCParam;
 class TPad;
@@ -46,6 +47,7 @@ public:
   void MakeQAPlot(Float_t  FPtoMIPratio);
   Bool_t AnalyzeGain(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesGaussFit = 500, Float_t FPtoMIPratio = 1.43); 
   Bool_t AnalyzeAttachment(Int_t startRunNumber, Int_t endRunNumber, Int_t minEntriesFit = 2000);
+  Bool_t AnalyzePadRegionGain();
   Bool_t ValidateTimeGain(Double_t minGain=2.0, Double_t maxGain = 3.0);
   //
   // QA drawing part
@@ -81,6 +83,7 @@ private:
   TObjArray    * fGainArray;               // array to be stored in the OCDB
   AliTPCcalibTimeGain * fGainMIP;          // calibration component for MIP
   AliTPCcalibTimeGain * fGainCosmic;       // calibration component for cosmic
+  AliTPCcalibGainMult * fGainMult;         // calibration component for pad region gain equalization and multiplicity dependence
   
   Bool_t fSwitchOnValidation;  // flag to switch on validation of OCDB parameters
 
index a2ff5c3..978d1e4 100644 (file)
@@ -98,6 +98,8 @@ AliTPCRecoParam::AliTPCRecoParam():
   fUseDriftCorrectionGY(1),   // use drif correction global y
   fUseGainCorrectionTime(0),  // use gain correction time
   fUseExBCorrection(1),  // use ExB correction
+  fUseMultiplicityCorrectionDedx(kTRUE), // use Dedx multiplicity correction
+  fUseAlignmentTime(kTRUE),              // use time dependent alignment correction
   //
   fUseTotCharge(kTRUE),          // switch use total or max charge
   fMinFraction(0.01),           // truncated mean - lower threshold
index 57f1315..4ace0d2 100644 (file)
@@ -57,6 +57,11 @@ class AliTPCRecoParam : public AliDetectorRecoParam
   void SetMinMaxCutSigma(Float_t th)       {   fMinMaxCutSigma=th; }
   void SetMinLeftRightCutSigma(Float_t th) {   fMinLeftRightCutSigma=th;}  // minimal amplitude left right - PRF
   void SetMinUpDownCutSigma(Float_t th)    {   fMinUpDownCutSigma=th;}  // minimal amplitude up-down - TRF 
+  void  SetUseTotCharge(Bool_t flag) {fUseTotCharge = flag;}
+  void  SetCtgRange(Double_t ctgRange) {fCtgRange = ctgRange;}
+  void  SetUseMultiplicityCorrectionDedx(Bool_t flag) {fUseMultiplicityCorrectionDedx = flag;}
+  void  SetUseAlignmentTime(Bool_t flag) {fUseAlignmentTime = flag;}
+
   //
   Int_t    GetLastSeedRowSec()       const  { return fLastSeedRowSec;} 
   Int_t    GetSeedGapPrim()        const  { return fSeedGapPrim;} 
@@ -97,13 +102,13 @@ class AliTPCRecoParam : public AliDetectorRecoParam
   Int_t GetUseDriftCorrectionGY() const {return fUseDriftCorrectionGY;}
   Int_t GetUseGainCorrectionTime() const {return fUseGainCorrectionTime;}
   Int_t GetUseExBCorrection() const {return fUseExBCorrection;}
+  Bool_t GetUseMultiplicityCorrectionDedx() const {return fUseMultiplicityCorrectionDedx;}
+  Bool_t GetUseAlignmentTime() const {return fUseAlignmentTime;}
   //
   Bool_t   GetUseTotCharge() const {return fUseTotCharge;}          // switch use total or max charge
   Float_t  GetMinFraction() const {return fMinFraction;}           // truncated mean - lower threshold
   Float_t  GetMaxFraction() const {return fMaxFaction;}            // truncated mean - upper threshold
-  //
-  void     SetUseTotCharge(Bool_t useTotCharge) {fUseTotCharge=useTotCharge; }
-  //
+
   Bool_t   GetUseTOFCorrection() {return fUseTOFCorrection;}
 
   //
@@ -166,6 +171,8 @@ class AliTPCRecoParam : public AliDetectorRecoParam
   Int_t fUseDriftCorrectionGY;   // use drif correction global y
   Int_t fUseGainCorrectionTime;  // use gain correction time
   Int_t fUseExBCorrection;       // use ExB correction
+  Bool_t fUseMultiplicityCorrectionDedx; // use Dedx multiplicity correction
+  Bool_t fUseAlignmentTime;              // use time dependent alignment correction
   //
   // dEdx switches
   //