1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // lu@physi.uni-heidelberg.de
23 grep " AliTRDdEdxUtils::" AliTRDdEdxUtils.cxx | grep "=" -v | grep -v "[6]" | grep -v printf |wc
24 grep "(" AliTRDdEdxUtils.h | grep ";" | grep -v grep | grep -v ClassDef | grep -v "{" | grep -v typedef | wc
28 #ifndef ALITRDDEDXUTILS_H
29 #define ALITRDDEDXUTILS_H
36 #include "THnSparse.h"
40 #include "TTreeStream.h"
57 //===================================================================================
59 //===================================================================================
60 static void GetCDFCuts(const TH1D *hh, const Int_t ncut, Double_t cuts[], const Double_t cdfs[], const Double_t thres);
61 static Double_t GetMeanRMS(const Double_t nn, const Double_t sum, const Double_t w2s, Double_t * grms=0x0, Double_t * gerr=0x0);
62 static Double_t TruncatedMean(const Int_t nx, const Double_t xdata[], const Double_t lowfrac, const Double_t highfrac, Double_t * grms=0x0, Double_t * gerr=0x0, Double_t *wws=0x0);
63 static Double_t TruncatedMean(const TH1 *hh, const Double_t lowfrac, const Double_t highfrac, Double_t * grms=0x0, Double_t * gerr=0x0);
64 static void FitSlicesY(const TH2D *hh, TH1D *&hnor, TH1D *&hmpv, TH1D *&hwid, TH1D *&hres, const Double_t thres, const Double_t lowfrac, const Double_t highfrac);
66 //===================================================================================
67 // TRD Analysis Fast Tool
68 //===================================================================================
69 static Int_t GetNtracklet(const AliESDEvent *esd);
70 static AliTRDtrackV1 * GetTRDtrackV1(const AliESDtrack * esdtrack);
71 static Bool_t IsInSameStack(const AliTRDtrackV1 *trdtrack);
72 static Bool_t GetFirstSectorStackMomentum(const AliTRDtrackV1 *trdtrack, Int_t & isec, Int_t & istk, Double_t & mom);
74 //===================================================================================
76 //===================================================================================
77 static void SetCalibFile(const TString file){fgCalibFile = file;}
78 static void DeleteCalibObj();
79 static void IniCalibObj();
81 static void SetObjPHQ(TObjArray * obj){fgObjPHQ = obj;}
82 static Bool_t GenerateDefaultPHQOCDB(const TString path="local://./");
84 static Double_t GetCalibTPCscale(const Int_t tpcncls, const Double_t tpcsig);
86 static void DeleteCalibHist();
87 static void IniCalibHist(TList *list, const Bool_t kPHQonly=kFALSE);
88 static Bool_t ReadCalibHist(const TString filename, const TString listname);
90 static TObjArray * GetObjPHQ();
91 static TObjArray * GetHistPHQ(){return fgHistPHQ;}
92 static TObjArray * GetObjPHQ(const Bool_t kinvq, const Double_t mag, const Int_t charge);
93 static THnSparseD * GetHistPHQ(const Bool_t kinvq, const Double_t mag, const Int_t charge);
95 static void FillCalibHist(const Int_t ncls, const TVectorD *arrayQ, const TVectorD *arrayX, THnSparseD * hcalib, const Double_t scale);
96 static void FillCalibHist(const AliTRDtrackV1 *trdv1, const Bool_t kinvq, const Double_t mag, const Int_t charge, const Double_t scale);
98 static void CalibOutput(const TList *l0, Int_t run);
99 static TObjArray* GetCalibObj(const THnSparseD *hh, Int_t run=-999, TList *lout=0x0, TTreeSRedirector *calibStream=0x0);
101 static THnSparseD * GetHistGain(){return fgHistGain;}
102 static THnSparseD * GetHistT0(){return fgHistT0;}
103 static THnSparseD * GetHistVd(){return fgHistVd;}
105 //===================================================================================
107 //===================================================================================
108 static Double_t ToyCook(const Bool_t kinvq, Int_t &ncluster, TVectorD *arrayQ, TVectorD *arrayX, const TObjArray *cobj=0x0);
109 static Double_t CombineddEdx(const Bool_t kinvq, Int_t &concls, TVectorD *coarrayQ, TVectorD *coarrayX, const Int_t tpcncls, const TVectorD *tpcarrayQ, const TVectorD *tpcarrayX, const Int_t trdncls, const TVectorD *trdarrayQ, const TVectorD *trdarrayX);
111 //===================================================================================
112 // dEdx Getter and Setter
113 //===================================================================================
114 static Int_t GetArrayClusterQ(const Bool_t kinvq, TVectorD *arrayQ, TVectorD *arrayX, const AliTRDtrackV1 *trdtrack, Int_t timeBin0=-1, Int_t timeBin1=1000, Int_t tstep=1);
115 static Int_t UpdateArrayX(const Int_t ncls, TVectorD* arrayX);
116 static void SetChamberQT(const AliTRDtrackV1 *trdtrack, const Int_t kcalib, THnSparseD * hgain=0x0, THnSparseD * ht0=0x0, THnSparseD * hvd=0x0);
118 static Int_t GetNchamber() {return fgNchamber;}
119 static Double_t GetChamberQ(const Int_t ich) {return fgChamberQ[ich];}
120 static Double_t GetChamberTmean(const Int_t ich) {return fgChamberTmean[ich];}
121 static Double_t GetTrackTmean() {return fgTrackTmean;}
123 //===================================================================================
124 // dEdx Parameterization
125 //===================================================================================
127 static Double_t MeandEdx(const Double_t * xx, const Double_t * par);
128 static Double_t MeanTR(const Double_t * xx, const Double_t * par);
129 static Double_t MeandEdxTR(const Double_t * xx, const Double_t * par);
131 static Double_t QMeanTPC(const Double_t bg);
132 static Double_t Q0MeanTRDpp(const Double_t bg);
133 static Double_t Q1MeanTRDpp(const Double_t bg);
134 static Double_t Q0MeanTRDPbPb(const Double_t bg);
135 static Double_t Q1MeanTRDPbPb(const Double_t bg);
137 typedef Double_t (*FFunc)(const Double_t *xx, const Double_t *par);
139 static Double_t MeandEdxLogx(const Double_t * xx, const Double_t * par){return ToLogx(MeandEdx, xx, par);}
140 static Double_t MeanTRLogx(const Double_t * xx, const Double_t * par){return ToLogx(MeanTR, xx, par);}
141 static Double_t MeandEdxTRLogx(const Double_t * xx, const Double_t * par){return ToLogx(MeandEdxTR, xx, par);}
143 //===================================================================================
144 // Detector, Data and Control Constant
145 //===================================================================================
146 static Int_t NTRDchamber(){return 18*5*6;} //540
147 static Int_t NTRDtimebin(){return NTRDchamber()*31;} //16740
148 static Int_t ToDetector(const Int_t gtb);
149 static Int_t ToTimeBin(const Int_t gtb);
150 static Int_t ToSector(const Int_t gtb);
151 static Int_t ToStack(const Int_t gtb);
152 static Int_t ToLayer(const Int_t gtb);
154 static TString GetRunType(const Int_t run);
156 static void SetPadGainOn(const Bool_t kon){ fgPadGainOn = kon; }
157 static void SetExBOn(const Bool_t kon){ fgExBOn = kon; }
158 static void SetQScale(const Double_t scale){ fgQScale = scale; }
159 static void SetQ0Frac(const Double_t q0){ fgQ0Frac = q0; }
160 static void SetQ1Frac(const Double_t q1){ fgQ1Frac = q1; }
161 static void SetTimeBinCountCut(const Double_t tbc){ fgTimeBinCountCut = tbc; }
162 static void SetCalibTPCnclsCut(const Int_t tpc){ fgCalibTPCnclsCut = tpc; }
164 static Bool_t IsPadGainOn(){return fgPadGainOn;}
165 static Bool_t IsExBOn(){return fgExBOn;}
166 static Double_t QScale(){return fgQScale;}
167 static Double_t Q0Frac(){return fgQ0Frac;}
168 static Double_t Q1Frac(){return fgQ1Frac;}
169 static Double_t TimeBinCountCut(){return fgTimeBinCountCut;}
170 static Int_t CalibTPCnclsCut(){return fgCalibTPCnclsCut;}
172 static void PrintControl();
173 //===================================================================================
174 //===================================================================================
178 //dEdx Getter and Setter
179 static Double_t GetAngularCorrection(const AliTRDseedV1 *seed);
180 static Double_t GetPadGain(const Int_t det, const Int_t icol, const Int_t irow);
181 static Double_t GetRNDClusterQ(AliTRDcluster *cl);
182 static Double_t GetClusterQ(const Bool_t kinvq, const AliTRDseedV1 * seed, const Int_t itb);
184 //dEdx Parameterization
185 static Double_t ToLogx(FFunc func, const Double_t * xx, const Double_t * par);
188 static Int_t ApplyCalib(const Int_t nc0, TVectorD *arrayQ, TVectorD *arrayX, const TObjArray *cobj);
189 static void GetPHCountMeanRMS(const TH1D *hnor, TH1D *&hmean);
190 static Int_t GetPHQIterator(const Bool_t kinvq, const Double_t mag, const Int_t charge);
191 static TString GetPHQName(const Bool_t kobj, const Int_t iter);
193 //Detector, Data and Control Constant
194 static THnSparseD *fgHistGain;//PH hist
195 static THnSparseD *fgHistT0;//PH hist
196 static THnSparseD *fgHistVd;//PH hist
197 static TObjArray * fgHistPHQ;//array containing 8 THnSparseD!
199 static TString fgCalibFile; //private path for calibration object
201 static TObjArray * fgObjGain;//chamber gain obj
202 static TObjArray * fgObjT0;//t0 obj
203 static TObjArray * fgObjVd;//Vd obj
204 static TObjArray * fgObjPHQ;//array containing 8 TObjArray!
206 static Int_t fgNchamber; //number of chamber used in cookdedx
207 static Double_t fgChamberQ[6]; //dqdl in chamber [i]
208 static Double_t fgChamberTmean[6]; //Q-weighted timebin \sum Q*T / \sum Q
210 static Double_t fgTrackTmean; //mean timebin over track
212 static Bool_t fgPadGainOn; //pad gain
213 static Bool_t fgExBOn; //exbon
214 static Double_t fgQScale; //Qscale
215 static Double_t fgQ0Frac; //q0frac
216 static Double_t fgQ1Frac; //q1frac
217 static Double_t fgTimeBinCountCut; //tbcut
218 static Int_t fgCalibTPCnclsCut; //tpccut