]>
Commit | Line | Data |
---|---|---|
e1ef49ad | 1 | #include "AliITSUMatLUT.h" |
2 | #include "AliLog.h" | |
3 | #include "AliTrackerBase.h" | |
4 | #include <TString.h> | |
5 | #include <TH1.h> | |
6 | #include <TMath.h> | |
7 | #include <TRandom.h> | |
8 | ||
9 | ||
10 | //___________________________________________________________ | |
11 | AliITSUMatLUT::AliITSUMatLUT() | |
12 | :fRMin(0) | |
13 | ,fRMax(0) | |
14 | ,fDRInv(-1) | |
15 | ,fDR(0) | |
16 | ,fNBins(0) | |
17 | { | |
18 | // def c-tor | |
19 | for (int i=kNParTypes;i--;) fData[i]=0; | |
20 | } | |
21 | ||
22 | //___________________________________________________________ | |
23 | AliITSUMatLUT::~AliITSUMatLUT() | |
24 | { | |
25 | // d-tor | |
26 | for (int i=kNParTypes;i--;) delete fData[i]; | |
27 | } | |
28 | ||
29 | //___________________________________________________________ | |
30 | AliITSUMatLUT::AliITSUMatLUT(Double_t rmin,Double_t rmax,Int_t nbin) | |
31 | :fRMin(rmin) | |
32 | ,fRMax(rmax) | |
33 | ,fDRInv(0) | |
34 | ,fDR(0) | |
35 | ,fNBins(nbin) | |
36 | { | |
37 | // | |
38 | if (rmin<0 || rmax-rmin<1e-4 || nbin<1) AliFatal(Form("Illegal parameters Rmin:%f Rmax:%f Nbins:%d",rmin,rmax,nbin)); | |
39 | fDRInv = fNBins/(fRMax-fRMin); | |
40 | fDR = (fRMax-fRMin)/fNBins; | |
41 | for (int i=kNParTypes;i--;) { | |
42 | fData[i] = new Double_t[fNBins]; | |
43 | memset(fData[i],0,fNBins*sizeof(Double_t)); | |
44 | } | |
45 | // | |
46 | } | |
47 | ||
48 | //___________________________________________________________ | |
49 | AliITSUMatLUT::AliITSUMatLUT(const AliITSUMatLUT& src) | |
50 | :TObject(src) | |
51 | ,fRMin(src.fRMin) | |
52 | ,fRMax(src.fRMax) | |
53 | ,fDRInv(src.fDRInv) | |
54 | ,fDR(src.fDR) | |
55 | ,fNBins(src.fNBins) | |
56 | { | |
57 | // | |
58 | for (int i=kNParTypes;i--;) { | |
59 | if (src.fData[i]) { | |
60 | fData[i] = new Double_t[fNBins]; | |
61 | memcpy(fData[i],src.fData[i],fNBins*sizeof(Double_t)); | |
62 | } | |
63 | } | |
64 | // | |
65 | } | |
66 | ||
67 | //___________________________________________________________ | |
68 | AliITSUMatLUT & AliITSUMatLUT::operator=(const AliITSUMatLUT& src) | |
69 | { | |
70 | // copy | |
71 | if (this == &src) return *this; | |
72 | this->~AliITSUMatLUT(); | |
73 | new(this) AliITSUMatLUT(src); | |
74 | return *this; | |
75 | // | |
76 | } | |
77 | ||
78 | //___________________________________________________________ | |
79 | void AliITSUMatLUT::FillData(Int_t ntest, Double_t zmin,Double_t zmax) | |
80 | { | |
81 | // filla material data | |
82 | double start[3],stop[3],parStep[7]; | |
83 | if (fNBins<1) AliFatal("Limits are not set"); | |
84 | if (ntest<1 || zmin>zmax) AliFatal(Form("Wrong parameters Ntest:%d Zmin:%f Zmax:%f",ntest,zmin,zmax)); | |
85 | double dr = (fRMax-fRMin)/fNBins; | |
86 | AliInfo(Form("Building material table for %.3f<R<%.3f %.3f<Z<%.3f in %d bins using %d tracks",fRMin,fRMax,zmin,zmax,fNBins,ntest)); | |
eb18cdc6 | 87 | const double kAngEps = 1e-4; // tiny slope to avoid tracks strictly normal to Z axis |
88 | double *tmpAcc = new double[fNBins*kNParTypes]; | |
e1ef49ad | 89 | for (int itst=ntest;itst--;) { |
90 | double parInt[kNParTypes]={0}; | |
91 | double r = fRMin; | |
92 | double phi = gRandom->Rndm()*TMath::Pi()*2; | |
93 | double cs = TMath::Cos(phi); | |
94 | double sn = TMath::Sin(phi); | |
eb18cdc6 | 95 | double angz = 2*(gRandom->Rndm()-0.5)*kAngEps; |
e1ef49ad | 96 | stop[0] = r*cs; |
97 | stop[1] = r*sn; | |
98 | stop[2] = zmin + gRandom->Rndm()*(zmax-zmin); | |
eb18cdc6 | 99 | Bool_t fail = kFALSE; |
e1ef49ad | 100 | for (int ir=0;ir<fNBins;ir++) { |
101 | r += dr; | |
102 | for (int i=3;i--;) start[i] = stop[i]; | |
103 | stop[0] = r*cs; | |
104 | stop[1] = r*sn; | |
eb18cdc6 | 105 | stop[2] += dr*angz; // to avoid tracks normal to axis |
e1ef49ad | 106 | AliTrackerBase::MeanMaterialBudget(start,stop, parStep); |
eb18cdc6 | 107 | if (parStep[1]>999) {fail = kTRUE; printf("fail\n"); break;} |
e1ef49ad | 108 | // |
109 | parInt[kParX2X0] += parStep[1]; | |
110 | parInt[kParRhoL] += parStep[0]*parStep[4]; | |
111 | // | |
eb18cdc6 | 112 | for (int ip=kNParTypes;ip--;) tmpAcc[ir*kNParTypes+ip] = parInt[ip]; |
e1ef49ad | 113 | } |
eb18cdc6 | 114 | if (fail) {itst++; continue;} // propagation failed |
115 | for (int ir=0;ir<fNBins;ir++) for (int ip=kNParTypes;ip--;) fData[ip][ir] += tmpAcc[ir*kNParTypes+ip]; | |
e1ef49ad | 116 | } |
117 | // | |
118 | for (int ip=kNParTypes;ip--;) for (int ir=fNBins;ir--;) fData[ip][ir] /= ntest; | |
eb18cdc6 | 119 | delete[] tmpAcc; |
120 | // Print(); | |
e1ef49ad | 121 | // |
122 | } | |
123 | ||
124 | //___________________________________________________________ | |
125 | void AliITSUMatLUT::Print(Option_t*) const | |
126 | { | |
127 | // print data | |
128 | printf("Average material budget in %d bins for %.4f<R<%.4f\n",fNBins,fRMin,fRMax); | |
129 | printf(" # : rMin : rMax \t X2X0 ( incr )\t RhoL ( incr )\n"); | |
130 | for (int i=0;i<fNBins;i++) { | |
131 | double r = fRMin + fDR*i; | |
132 | printf("%4d:%7.3f:%7.3f\t%.6f(%.6f)\t%.6f(%.6f)\n",i,r,r+fDR, | |
133 | fData[kParX2X0][i],fData[kParX2X0][i]-(i==0 ? 0:fData[kParX2X0][i-1]), | |
134 | fData[kParRhoL][i],fData[kParRhoL][i]-(i==0 ? 0:fData[kParRhoL][i-1])); | |
135 | } | |
136 | ||
137 | } | |
138 | ||
139 | //___________________________________________________________ | |
140 | TH1* AliITSUMatLUT::GetHisto(const Option_t* option, const Char_t *name) const | |
141 | { | |
142 | // extract data to histo | |
143 | if (fNBins<1) return 0; | |
144 | TString nms = name; | |
145 | if (nms.IsNull()) nms = "matLUT"; | |
146 | TH1F* h = new TH1F(nms.Data(),nms.Data(),fNBins,fRMin,fRMax); | |
147 | TString opts = option; | |
148 | opts.ToLower(); | |
149 | Bool_t diff = opts.Contains("d"); | |
150 | Int_t par = opts.Contains("rhol") ? kParRhoL : kParX2X0; | |
151 | double valPrev = 0; | |
152 | for (int i=0;i<fNBins;i++) { | |
153 | double val = fData[par][i] - valPrev; | |
154 | if (diff) valPrev = fData[par][i]; | |
155 | h->SetBinContent(i+1, val); | |
156 | } | |
157 | return h; | |
158 | } | |
159 | ||
160 | //___________________________________________________________ | |
161 | Double_t AliITSUMatLUT::GetMatBudget(const Double_t *pnt0, const Double_t *pnt1, Double_t *ret) const | |
162 | { | |
163 | // query the mat.budget between 2 points | |
164 | double r0 = TMath::Sqrt(pnt0[0]*pnt0[0] + pnt0[1]*pnt0[1]); | |
165 | double r1 = TMath::Sqrt(pnt1[0]*pnt1[0] + pnt1[1]*pnt1[1]); | |
166 | if (r1<r0) {double t=r1;r1=r0;r0=t;} | |
167 | double dr = r1-r0; | |
168 | if (dr<1e-4) { | |
169 | for (int i=kNParTypes;i--;) ret[i] = 0; | |
170 | return 0; | |
171 | } | |
172 | // | |
173 | double di = pnt1[0]-pnt0[0]; | |
174 | double dst = di*di; | |
175 | di = pnt1[1]-pnt0[1]; | |
176 | dst += di*di; | |
177 | di = pnt1[2]-pnt0[2]; | |
178 | dst += di*di; | |
179 | dst = TMath::Sqrt(dst); | |
180 | double dstf = dst/dr; | |
181 | // | |
182 | double par0[kNParTypes],par1[kNParTypes]; | |
183 | GetData(r0,par0); | |
184 | GetData(r1,par1); | |
185 | for (int i=kNParTypes;i--;) ret[i] = dstf*(par1[i]-par0[i]); | |
186 | return dst; | |
187 | // | |
188 | } | |
189 | ||
190 | //___________________________________________________________ | |
191 | void AliITSUMatLUT::GetData(Double_t r, Double_t* params) const | |
192 | { | |
193 | // get integrated parameters from rMin to r | |
194 | double binF = (r-fRMin)*fDRInv; | |
195 | if (binF<0) binF=0; | |
196 | else if (binF>=fNBins) binF = fNBins-1; | |
197 | int bin = int(binF); | |
198 | double frac = 1.-(binF-bin); | |
199 | for (int i=kNParTypes;i--;) { | |
200 | double prev = bin ? fData[i][bin-1] : 0; | |
201 | params[i] = fData[i][bin] - frac*(fData[i][bin]-prev); | |
202 | } | |
203 | // | |
204 | } |