Line data Source code
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 :
16 : #include "AliCheb2DStackF.h"
17 : #include "AliCheb2DStackS.h"
18 : #include "AliTPCChebCorr.h"
19 : #include "AliLog.h"
20 : #include <TMath.h>
21 : #include <TH1F.h>
22 : #include <TAxis.h>
23 : #include <TGraph.h>
24 : #include <TBits.h>
25 :
26 24 : ClassImp(AliTPCChebCorr)
27 :
28 : const char* AliTPCChebCorr::fgkFieldTypeName[4] = {"Any","B>0"," B<0","B=0"};
29 :
30 24 : const float AliTPCChebCorr::fgkY2XHSpan = TMath::Tan(TMath::Pi()/18);
31 :
32 : const float AliTPCChebCorr::fgkPadRowX[AliTPCChebCorr::kNRows] = {
33 : 85.225, 85.975, 86.725, 87.475, 88.225, 88.975, 89.725, 90.475, 91.225, 91.975, 92.725, 93.475, 94.225, 94.975, 95.725,
34 : 96.475, 97.225, 97.975, 98.725, 99.475,100.225,100.975,101.725,102.475,103.225,103.975,104.725,105.475,106.225,106.975,
35 : 107.725,108.475,109.225,109.975,110.725,111.475,112.225,112.975,113.725,114.475,115.225,115.975,116.725,117.475,118.225,
36 : 118.975,119.725,120.475,121.225,121.975,122.725,123.475,124.225,124.975,125.725,126.475,127.225,127.975,128.725,129.475,
37 : 130.225,130.975,131.725,135.100,136.100,137.100,138.100,139.100,140.100,141.100,142.100,143.100,144.100,145.100,146.100,
38 : 147.100,148.100,149.100,150.100,151.100,152.100,153.100,154.100,155.100,156.100,157.100,158.100,159.100,160.100,161.100,
39 : 162.100,163.100,164.100,165.100,166.100,167.100,168.100,169.100,170.100,171.100,172.100,173.100,174.100,175.100,176.100,
40 : 177.100,178.100,179.100,180.100,181.100,182.100,183.100,184.100,185.100,186.100,187.100,188.100,189.100,190.100,191.100,
41 : 192.100,193.100,194.100,195.100,196.100,197.100,198.100,199.350,200.850,202.350,203.850,205.350,206.850,208.350,209.850,
42 : 211.350,212.850,214.350,215.850,217.350,218.850,220.350,221.850,223.350,224.850,226.350,227.850,229.350,230.850,232.350,
43 : 233.850,235.350,236.850,238.350,239.850,241.350,242.850,244.350,245.850
44 : };
45 :
46 : //____________________________________________________________________
47 : AliTPCChebCorr::AliTPCChebCorr()
48 0 : : TNamed()
49 0 : ,fOnFlyInitDone(kFALSE)
50 0 : ,fFieldType(kFieldAny)
51 0 : ,fRun(-1)
52 0 : ,fNRows(0)
53 0 : ,fNStacksSect(0)
54 0 : ,fNStacksZSect(0)
55 0 : ,fNStacksZ(0)
56 0 : ,fNStacks(0)
57 0 : ,fZMaxAbs(-1)
58 0 : ,fTimeStampStart(0)
59 0 : ,fTimeStampEnd(0xffffffff)
60 0 : ,fZScaleI(0)
61 0 : ,fY2XScaleI(0)
62 0 : ,fDeadZone(0)
63 0 : ,fRowXI(0)
64 0 : ,fParams(0)
65 0 : ,fTracksRate(0)
66 0 : {
67 : // def. c-tor
68 0 : }
69 :
70 : //____________________________________________________________________
71 : AliTPCChebCorr::AliTPCChebCorr(const char* name, const char* title,
72 : int nps, int nzs, float zmaxAbs, float deadZone, const float *xrow)
73 0 : : TNamed(name,title)
74 0 : ,fOnFlyInitDone(kFALSE)
75 0 : ,fFieldType(kFieldAny)
76 0 : ,fRun(-1)
77 0 : ,fNRows(kNRows)
78 0 : ,fNStacksSect(0)
79 0 : ,fNStacksZSect(0)
80 0 : ,fNStacksZ(0)
81 0 : ,fNStacks(0)
82 0 : ,fZMaxAbs(-1)
83 0 : ,fTimeStampStart(0)
84 0 : ,fTimeStampEnd(0xffffffff)
85 0 : ,fZScaleI(0)
86 0 : ,fY2XScaleI(0)
87 0 : ,fDeadZone(deadZone)
88 0 : ,fRowXI(0)
89 0 : ,fParams(0)
90 0 : ,fTracksRate(0)
91 0 : {
92 : // c-tor, optional parameter xrow provides the X of every slice if deadZone is requested
93 : // (if xrow is absent, default will be used)
94 : //
95 0 : if (deadZone>0) {
96 0 : fRowXI = new Float_t[kNRows];
97 0 : if (xrow) for (int i=kNRows;i--;) fRowXI[i] = 1./xrow[i];
98 0 : else for (int i=kNRows;i--;) fRowXI[i] = 1.0f/fgkPadRowX[i];
99 : }
100 : //
101 0 : SetBinning(nps,nzs,zmaxAbs);
102 : //
103 0 : }
104 :
105 : //____________________________________________________________________
106 : AliTPCChebCorr::~AliTPCChebCorr()
107 0 : {
108 : // d-tor
109 0 : if (fParams) for (int i=fNStacks;i--;) delete fParams[i];
110 0 : delete[] fRowXI;
111 0 : delete[] fParams;
112 0 : delete fTracksRate;
113 0 : }
114 :
115 : //____________________________________________________________________
116 : void AliTPCChebCorr::Parameterize(stFun_t fun,int dimOut,const int np[2],const float* prec)
117 : {
118 : // build parameterizations for 2->dimout, on the same grid of np[0]xnp[1] points for
119 : // every output dimension, optionally with prec[i] precision for i-th dimension
120 : //
121 0 : if (TestBit(kParamDone)) {
122 0 : AliError("Parameterization is already done");
123 0 : return;
124 : }
125 : //
126 0 : if (fZMaxAbs<0) AliFatal("First the binning and Z limits should be set");
127 : //
128 0 : float y2xMax = GetMaxY2X(); // half-sector span
129 0 : float bmn[2],bmx[2];
130 0 : Bool_t useS = GetUseShortPrec(); // float or short representation for param
131 0 : fParams = new AliCheb2DStack*[fNStacks];
132 : //
133 0 : for (int iz=0;iz<fNStacksZ;iz++) {
134 0 : bmn[1] = iz/fZScaleI - fZMaxAbs; // boundaries in Z
135 0 : bmx[1] = iz<fNStacksZ-1 ? (iz+1)/fZScaleI - fZMaxAbs : fZMaxAbs;
136 0 : for (int isc=0;isc<kNSectors;isc++) {
137 : int isc72 = isc;
138 0 : if (iz<fNStacksZSect) isc72 += 18; // change side
139 0 : fun(isc72,0,0);
140 0 : for (int isl=0;isl<fNStacksSect;isl++) {
141 0 : float dead[2] = {0.,0.};
142 0 : if (fDeadZone>0 && isl==0) { // assign dead zone
143 0 : dead[0] = fDeadZone;
144 0 : }
145 0 : if (fDeadZone>0 && isl==fNStacksSect-1) { // assign dead zone
146 0 : dead[1] = fDeadZone;
147 0 : }
148 : //
149 0 : bmn[0] = -y2xMax+isl/fY2XScaleI; // boundaries in phi
150 0 : bmx[0] = isl<fNStacksSect-1 ? -y2xMax+(isl+1)/fY2XScaleI : y2xMax;
151 0 : int id = GetParID(iz,isc,isl);
152 0 : AliInfoF("Doing param #%03d Iz:%d Sect:%02d Slice:%d | %+.1f<Z<%+.1f %+.3f<y2x<%+.3f",
153 : id,iz,isc,isl,bmn[1],bmx[1],bmn[0],bmx[0]);
154 0 : if (useS) fParams[id] = new AliCheb2DStackS(fun,fNRows,dimOut,bmn,bmx,np,dead,fRowXI,prec);
155 0 : else fParams[id] = new AliCheb2DStackF(fun,fNRows,dimOut,bmn,bmx,np,dead,fRowXI,prec);
156 0 : }
157 : }
158 : }
159 : //
160 0 : fOnFlyInitDone = kTRUE;
161 0 : SetBit(kParamDone);
162 : //
163 0 : }
164 :
165 : //____________________________________________________________________
166 : void AliTPCChebCorr::Parameterize(stFun_t fun,int dimOut,const int np[][2],const float* prec)
167 : {
168 : // build parameterizations for 2->dimout, on the same grid of np[0]xnp[1] points for
169 : // every output dimension, optionally with prec[i] precision for i-th dimension
170 : //
171 0 : if (TestBit(kParamDone)) {
172 0 : AliError("Parameterization is already done");
173 0 : return;
174 : }
175 : //
176 0 : float bmn[2],bmx[2];
177 0 : float y2xMax = GetMaxY2X(); // half-sector span
178 : //
179 0 : fParams = new AliCheb2DStack*[fNStacks];
180 0 : Bool_t useS = GetUseShortPrec(); // float or short representation for param
181 : //
182 0 : for (int iz=0;iz<fNStacksZ;iz++) {
183 0 : bmn[1] = iz/fZScaleI-fZMaxAbs; // boundaries in Z
184 0 : bmx[1] = iz<fNStacksZ-1 ? (iz+1)/fZScaleI-fZMaxAbs : fZMaxAbs;
185 0 : for (int isc=0;isc<kNSectors;isc++) {
186 : int isc72 = isc;
187 0 : if (iz<fNStacksZSect) isc72 += 18; // change side
188 0 : fun(isc72,0,0); // this will set sector ingo to function
189 : //
190 0 : for (int isl=0;isl<fNStacksSect;isl++) {
191 0 : float dead[2] = {0.,0.};
192 0 : if (fDeadZone>0 && isl==0) { // assign dead zone
193 0 : dead[0] = fDeadZone;
194 0 : }
195 0 : if (fDeadZone>0 && isl==fNStacksSect-1) { // assign dead zone
196 0 : dead[1] = fDeadZone;
197 0 : }
198 : //
199 0 : bmn[0] = -y2xMax+isl/fY2XScaleI; // boundaries in phi
200 0 : bmx[0] = isl<fNStacksSect-1 ? -y2xMax+(isl+1)/fY2XScaleI : y2xMax;
201 0 : int id = GetParID(iz,isc,isl);
202 0 : AliInfoF("Doing param #%03d Iz:%d Sect:%02d Slice:%d | %+.1f<Z<%+.1f %+.3f<y2x<%+.3f",
203 : id,iz,isc,isl,bmn[1],bmx[1],bmn[0],bmx[0]);
204 0 : if (useS) fParams[id] = new AliCheb2DStackS(fun,fNRows,dimOut,bmn,bmx,np,dead,fRowXI,prec);
205 0 : else fParams[id] = new AliCheb2DStackF(fun,fNRows,dimOut,bmn,bmx,np,dead,fRowXI,prec);
206 0 : }
207 : }
208 : }
209 : //
210 0 : fOnFlyInitDone = kTRUE;
211 0 : SetBit(kParamDone);
212 : //
213 0 : }
214 :
215 : //____________________________________________________________________
216 : void AliTPCChebCorr::Print(const Option_t* opt) const
217 : {
218 : // print itself
219 0 : printf("%s:%s Cheb2D[%c] Param: %d slices in %+.1f<%s<%+.1f %d per sector. DeadZone: %.1fcm\n",
220 0 : GetName(),GetTitle(),GetUseFloatPrec()?'F':'S',
221 0 : fNStacksZ,-fZMaxAbs,GetUseZ2R() ? "Z/R":"Z",fZMaxAbs,fNStacksSect,fDeadZone);
222 0 : printf("Run used: %d Time span: %ld:%ld TimeDependent flag: %s Field type: %s\n",
223 0 : GetRun(),fTimeStampStart,fTimeStampEnd,GetTimeDependent() ? "ON ":"OFF", fgkFieldTypeName[fFieldType]);
224 0 : TString opts = opt; opts.ToLower();
225 0 : if (opts.Contains("p") && TestBit(kParamDone)) {
226 0 : for (int iz=0;iz<fNStacksZ;iz++) {
227 0 : for (int isc=0;isc<kNSectors;isc++) {
228 0 : for (int isl=0;isl<fNStacksSect;isl++) {
229 0 : int id = GetParID(iz,isc,isl);
230 0 : printf("Z%d Sector%02d Slice%d: ",iz,isc,isl);
231 0 : GetParam(id)->Print(opt);
232 : }
233 : }
234 : }
235 0 : }
236 0 : }
237 :
238 : //____________________________________________________________________
239 : void AliTPCChebCorr::SetBinning(int nps,int nzs, float zmxAbs)
240 : {
241 : // set binning, limits
242 0 : fNStacksSect = nps;
243 0 : fNStacksZSect = nzs;
244 0 : fNStacksZ = nzs*2; // nzs is the number of bins per side!
245 0 : fZMaxAbs = zmxAbs;
246 0 : fNStacks = fNStacksZ*fNStacksSect*kNSectors;
247 : //
248 0 : if (zmxAbs<1e-8 || nzs<1 || nps<1) AliFatalF("Wrong settings: |Zmax|:%f Nz:%d Nphi:%d",
249 : zmxAbs,nzs,nps);
250 0 : fZScaleI = nzs/zmxAbs;
251 0 : fY2XScaleI = nps/(2*TMath::Tan(TMath::Pi()/kNSectors));
252 : //
253 0 : }
254 :
255 : //____________________________________________________________________
256 : void AliTPCChebCorr::Init()
257 : {
258 : // make necessary initializations
259 0 : if (fOnFlyInitDone) return;
260 0 : AliInfo("Doing on-the-fly initialization");
261 0 : if (fRowXI && fParams) {
262 0 : for (int i=fNStacks;i--;) {
263 0 : if (fParams[i] && !fParams[i]->GetXRowInv()) fParams[i]->SetXRowInv(fRowXI);
264 : }
265 0 : }
266 0 : fOnFlyInitDone = kTRUE;
267 0 : }
268 :
269 : //____________________________________________________________________
270 : Int_t AliTPCChebCorr::GetDimOut() const
271 : {
272 : // get number of output dimensions
273 0 : const AliCheb2DStack* par = GetParam(0);
274 0 : return par ? par->GetDimOut() : 0;
275 : }
276 :
277 : //____________________________________________________________________
278 : Double_t AliTPCChebCorr::GetLuminosityCOG(TGraph* lumi, time_t tmin, time_t tmax) const
279 : {
280 : // calculate used-track-weighted luminosity COG
281 0 : if (!fTracksRate) {
282 0 : AliError("Tracks per time stamp histo is not stored");
283 0 : return 0.0;
284 : }
285 0 : TAxis* tax = fTracksRate->GetXaxis();
286 0 : time_t tmnH = tax->GetXmin();
287 0 : time_t tmxH = tax->GetXmax();
288 0 : tmax = (tmax>tmin && tmax>tmnH && tmax<tmxH) ? tmax : tmxH;
289 0 : tmin = (tmin>tmnH && tmin<tmax) ? tmin : tmnH;
290 0 : int bmn = tax->FindBin(tmin);
291 0 : int bmx = tax->FindBin(tmax);
292 : double lumiCOG = 0.0, norm = 0.0;
293 0 : for (int ib=bmn;ib<bmx;ib++) {
294 0 : double wgh = fTracksRate->GetBinContent(ib);
295 0 : lumiCOG += wgh*lumi->Eval(tax->GetBinCenter(ib));
296 0 : norm += wgh;
297 : }
298 0 : if (norm>0) lumiCOG /= norm;
299 : return lumiCOG;
300 0 : }
301 :
302 : //__________________________________________
303 : Int_t AliTPCChebCorr::GetRun() const
304 : {
305 : // get run number used for this map
306 0 : if (fRun>=0) return fRun;
307 : // try to extract from name
308 0 : TString runstr = "";
309 0 : const char* nm = GetName();
310 0 : while (*nm && !isdigit(*nm)) nm++;
311 0 : while (*nm && isdigit(*nm)){ runstr += *nm;nm++;}
312 0 : if (!runstr.IsDigit()) return -1;
313 0 : return runstr.Atoi();
314 0 : }
315 :
316 : //__________________________________________
317 : Bool_t AliTPCChebCorr::IsRowMasked(int sector72,int row) const
318 : {
319 : // check if row was masked (sector 0-71, IROC/OROC rows convention
320 0 : if (sector72>kMaxIROCSector) row += kNRowsIROC; // we are in OROC
321 0 : float tz[2] = {0.f,((sector72/kNSectors)&0x1) ? -0.1f:0.1f}; // where to query
322 0 : const AliCheb2DStack* par = GetParam(sector72,tz[0],tz[1]);
323 0 : if (!par) return kTRUE;
324 0 : return par->Eval(row,par->GetDimOut()-1,tz)<1e-6; // blocked rows have dispersion set to 0
325 0 : }
326 :
327 : //__________________________________________
328 : Int_t AliTPCChebCorr::GetNMaskedRows(int sector72,TBits* masked) const
329 : {
330 : // count masked rows, if TBits provided, set the masked rows bits
331 : int nmasked = 0;
332 0 : int nr = (sector72/kNSectorsIROC)==0 ? kNRowsIROC : kNRows-kNRowsIROC;
333 0 : for (int ir=nr;ir--;) {
334 0 : if (IsRowMasked(sector72,ir)) {nmasked++; if (masked) masked->SetBitNumber(ir);}
335 : }
336 0 : return nmasked;
337 : }
|
