adding next sector check (to be improved)
[u/mrichter/AliRoot.git] / TPC / Upgrade / AliToyMCEventGenerator.cxx
CommitLineData
526ddf0e 1#include <iostream>
a1a695e5 2
526ddf0e 3#include <TDatabasePDG.h>
4#include <TRandom.h>
5#include <TH2F.h>
a1a695e5 6#include <TGeoGlobalMagField.h>
32438f4e 7#include <TSpline.h>
cd8ed0ac 8#include <TObjString.h>
a1a695e5 9
b0c26b55 10#include <AliLog.h>
a1a695e5 11#include <AliTPCROC.h>
12#include <AliTrackPointArray.h>
526ddf0e 13#include <AliTrackerBase.h>
14#include <AliCDBManager.h>
15#include <AliTPCParam.h>
16#include <AliGeomManager.h>
17#include <AliTPCcalibDB.h>
526ddf0e 18#include <AliTPCclusterMI.h>
a1a695e5 19#include <AliTPCSpaceCharge3D.h>
32438f4e 20#include <AliTPCROC.h>
a1a695e5 21
22#include "AliToyMCEvent.h"
23#include "AliToyMCTrack.h"
24
25#include "AliToyMCEventGenerator.h"
26
de0014b7 27ClassImp(AliToyMCEventGenerator);
526ddf0e 28
29
de0014b7 30AliToyMCEventGenerator::AliToyMCEventGenerator()
526ddf0e 31 :TObject()
32 ,fTPCParam(0x0)
a1a695e5 33 ,fEvent(0x0)
223d9e38 34 ,fCurrentTrack(0)
35 ,fTPCCorrection(0x0)
36 ,fCorrectionFile("$ALICE_ROOT/TPC/Calib/maps/SC_NeCO2_eps5_50kHz_precal.lookup.root")
a1a695e5 37 ,fOutputFileName("toyMC.root")
38 ,fOutFile(0x0)
39 ,fOutTree(0x0)
d92b8630 40 ,fUseStepCorrection(kFALSE)
d1cf83f5 41 ,fUseMaterialBudget(kFALSE)
526ddf0e 42{
526ddf0e 43 fTPCParam = AliTPCcalibDB::Instance()->GetParameters();
44 fTPCParam->ReadGeoMatrices();
223d9e38 45 gRandom->SetSeed();
526ddf0e 46}
47//________________________________________________________________
de0014b7 48AliToyMCEventGenerator::AliToyMCEventGenerator(const AliToyMCEventGenerator &gen)
526ddf0e 49 :TObject(gen)
a1a695e5 50 ,fTPCParam(gen.fTPCParam)
51 ,fEvent(0x0)
223d9e38 52 ,fCurrentTrack(0)
53 ,fTPCCorrection(gen.fTPCCorrection)
54 ,fCorrectionFile(gen.fCorrectionFile)
a1a695e5 55 ,fOutputFileName(gen.fOutputFileName)
56 ,fOutFile(0x0)
57 ,fOutTree(0x0)
d92b8630 58 ,fUseStepCorrection(gen.fUseStepCorrection)
d1cf83f5 59 ,fUseMaterialBudget(gen.fUseMaterialBudget)
526ddf0e 60{
61 //
223d9e38 62 gRandom->SetSeed();
526ddf0e 63}
64//________________________________________________________________
de0014b7 65AliToyMCEventGenerator::~AliToyMCEventGenerator()
526ddf0e 66{
223d9e38 67 delete fTPCCorrection;
526ddf0e 68}
32438f4e 69
526ddf0e 70//________________________________________________________________
a1a695e5 71Bool_t AliToyMCEventGenerator::DistortTrack(AliToyMCTrack &trackIn, Double_t t0)
72{
73 //
74 //
75 //
76
526ddf0e 77 if(!fTPCParam) {
78 fTPCParam = AliTPCcalibDB::Instance()->GetParameters();
79 fTPCParam->ReadGeoMatrices();
526ddf0e 80 }
32438f4e 81
1e62e876 82 MakeITSClusters(trackIn/*,t0*/);
83 MakeTPCClusters(trackIn, t0);
84 MakeTRDClusters(trackIn/*,t0*/);
85
86 return 1;
87}
88//________________________________________________________________
89void AliToyMCEventGenerator::MakeITSClusters(AliToyMCTrack &trackIn/*, Double_t t0*/)
90{
91 //Upgrade ITS parameters
92 const Int_t nITSLayers = 7;
93 const Double_t ITSRadii[nITSLayers] = {2.2, 2.8, 3.6, 20.0, 22.0, 41.0, 43.0};
94 const Double_t lengthITS[nITSLayers] = {22.4, 24.2, 26.8, 78.0, 83.6, 142.4, 148.6};
95
96 const Double_t kMaxSnp = 0.85;
97 const Double_t kMaxZ0 = fTPCParam->GetZLength();
98 const Double_t kMass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
99
100 AliToyMCTrack track(trackIn);
101 Double_t xyz[3] = {0.,0.,0.};
102
103 if (!AliTrackerBase::PropagateTrackTo(&track,ITSRadii[0],kMass,5,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
104 AliError(Form("Propagation to %.2f failed\n",ITSRadii[0]));
105 return;
106 }
107
108 for(Int_t iLayer = 0; iLayer<nITSLayers; iLayer++){
109
110 if (!AliTrackerBase::PropagateTrackTo(&track,ITSRadii[iLayer],kMass,1,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
111 AliError(Form("Propagation to %.2f failed\n",ITSRadii[iLayer]));
112 continue;
113 }
114 track.GetXYZ(xyz);
115
116 if (TMath::Abs(track.GetZ())>kMaxZ0) continue;
117 if (TMath::Abs(track.GetZ())>lengthITS[iLayer]/2) continue;
118
119 AliCluster* tempCl = new AliCluster();
120
121 tempCl->SetX(xyz[0]);
122 tempCl->SetY(xyz[1]);
123 tempCl->SetZ(xyz[2]);
124
125 Double_t sigmaY = 0.0004;
126 Double_t sigmaZ = 0.0004;
127
128 tempCl->SetSigmaY2(sigmaY*sigmaY);
129 tempCl->SetSigmaZ2(sigmaZ*sigmaZ);
130
131 trackIn.AddITSPoint(*tempCl);
132
133 }
134
135}
136//________________________________________________________________
137void AliToyMCEventGenerator::MakeTRDClusters(AliToyMCTrack &trackIn/*, Double_t t0*/)
138
139{
140 //Uses current TRD parameters
141 const Int_t nTRDLayers = 6;
142 const Double_t distToMid = 3.2 + 30./2; //dist to middle of drift region (radiator + half drift region)
143 const Double_t TRDRadii[nTRDLayers] = {294.5 + distToMid, 307.1 + distToMid, 319.7 + distToMid, 332.3 + distToMid, 344.9 + distToMid, 357.5 + distToMid};
144 const Double_t lengthTRD[nTRDLayers] = {604.0, 634.0, 656.0, 686.0, 700.0, 700.0};
145
146
147 const Double_t kMaxSnp = 0.85;
148 const Double_t kMaxZ0 = fTPCParam->GetZLength();
149 const Double_t kMass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
150
151 AliToyMCTrack track(trackIn);
152 Double_t xyz[3] = {0.,0.,0.};
153
154 if (!AliTrackerBase::PropagateTrackTo(&track,TRDRadii[0],kMass,5,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
155 AliError(Form("Propagation to %.2f failed\n",TRDRadii[0]));
156 return;
157 }
158
159 for(Int_t iLayer = 0; iLayer<nTRDLayers; iLayer++){
160
161 if (!AliTrackerBase::PropagateTrackTo(&track,TRDRadii[iLayer],kMass,1,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
162 AliError(Form("Propagation to %.2f failed\n",TRDRadii[iLayer]));
163 continue;
164 }
165 track.GetXYZ(xyz);
166
167 if (TMath::Abs(track.GetZ())>kMaxZ0) continue;
168 if (TMath::Abs(track.GetZ())>lengthTRD[iLayer]/2) continue;
169
170 AliCluster* tempCl = new AliCluster();
171
172 tempCl->SetX(xyz[0]);
173 tempCl->SetY(xyz[1]);
174 tempCl->SetZ(xyz[2]);
175
176 Double_t sigmaY = 0.06;
177 Double_t sigmaZ = 0.2;
178
179 tempCl->SetSigmaY2(sigmaY*sigmaY);
180 tempCl->SetSigmaZ2(sigmaZ*sigmaZ);
181
182 trackIn.AddTRDPoint(*tempCl);
183
184 }
185
186}
187//________________________________________________________________
188void AliToyMCEventGenerator::MakeTPCClusters(AliToyMCTrack &trackIn, Double_t t0)
189{
190
32438f4e 191 // make it big enough to hold all points
192 // store real number of generated points in the unique id
193 const Int_t nMaxPoints=3000;
194 AliTrackPointArray pointArray0(nMaxPoints); //undistorted
195 AliTrackPointArray pointArray1(nMaxPoints); //distorted
1e62e876 196
32438f4e 197 //Create space point of undistorted and distorted clusters along the propagated track trajectory
198 CreateSpacePoints(trackIn,pointArray0,pointArray1);
199 //Convert the space points into clusters in the local frame
200 //for undistorted and distorted clusters using the same function
201 ConvertTrackPointsToLocalClusters(pointArray0,trackIn,t0,0);
202 ConvertTrackPointsToLocalClusters(pointArray1,trackIn,t0,1);
526ddf0e 203
526ddf0e 204
1e62e876 205}
32438f4e 206//________________________________________________________________
207void AliToyMCEventGenerator::CreateSpacePoints(AliToyMCTrack &trackIn,
208 AliTrackPointArray &arrUdist,
209 AliTrackPointArray &arrDist)
210{
211 //
212 // sample the track from the inner to the outer wall of the TPC
213 // a graph is filled in local coordinates for later
214 //
215
216 const Double_t kMaxSnp = 0.85;
217 const Double_t kMaxZ0 = fTPCParam->GetZLength();
218 const Double_t kMass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
219
220 const Double_t iFCRadius = 83.5; //radius constants found in AliTPCCorrection.cxx
221 const Double_t oFCRadius = 254.5;
222
223 AliToyMCTrack track(trackIn);
224 //!!! TODO: make this adjustable perhaps
225 const Double_t stepSize=0.1;
226 Double_t xyz[3] = {0.,0.,0.};
227 Float_t xyzf[3] = {0.,0.,0.};
228
a1a695e5 229 //!!! when does the propagation not work, how often does it happen?
d1cf83f5 230 if (!AliTrackerBase::PropagateTrackTo(&track,iFCRadius,kMass,5,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
b0c26b55 231 AliError(Form("Propagation to IFC: %.2f failed\n",iFCRadius));
32438f4e 232 return;
b0c26b55 233 }
32438f4e 234
526ddf0e 235 Int_t npoints=0;
b0c26b55 236
32438f4e 237 for (Double_t radius=iFCRadius; radius<oFCRadius; radius+=stepSize){
238
b0c26b55 239 //!!! changed from return 0 to continue -> Please check
d1cf83f5 240 if (!AliTrackerBase::PropagateTrackTo(&track,radius,kMass,1,kTRUE,kMaxSnp,0,kFALSE,fUseMaterialBudget)) {
32438f4e 241 AliError(Form("Propagation to %.2f failed\n",radius));
b0c26b55 242 continue;
243 }
526ddf0e 244 track.GetXYZ(xyz);
32438f4e 245
b0c26b55 246 //!!! Why is this smeared
05da1b4e 247// xyz[0]+=gRandom->Gaus(0,0.000005);
248// xyz[1]+=gRandom->Gaus(0,0.000005);
249// xyz[2]+=gRandom->Gaus(0,0.000005);
32438f4e 250
251 xyzf[0]=Float_t(xyz[0]);
252 xyzf[1]=Float_t(xyz[1]);
253 xyzf[2]=Float_t(xyz[2]);
254
526ddf0e 255 if (TMath::Abs(track.GetZ())>kMaxZ0) continue;
32438f4e 256 // if (TMath::Abs(track.GetX())<iFCRadius) continue;
257 // if (TMath::Abs(track.GetX())>oFCRadius) continue;
258
259
260 AliTrackPoint pUdist; // undistorted space point
261 AliTrackPoint pDist; // distorted space point
262 // Set undistorted point
263 SetPoint(xyzf,pUdist);
264 arrUdist.AddPoint(npoints, &pUdist);
265 Int_t sector=pUdist.GetVolumeID();
32438f4e 266
267 // set distorted point
526ddf0e 268 Float_t distPoint[3]={xyz[0],xyz[1],xyz[2]};
d92b8630 269 Float_t dxyz[3]={0.,0.,0.};
270 if (!fUseStepCorrection){
223d9e38 271 fTPCCorrection->DistortPoint(distPoint, sector);
d92b8630 272 } else {
223d9e38 273 fTPCCorrection->GetCorrectionIntegralDz(distPoint,sector,dxyz,5);
d92b8630 274 distPoint[0]-=dxyz[0];
275 distPoint[1]-=dxyz[1];
276 distPoint[2]-=dxyz[2];
277 }
32438f4e 278 SetPoint(distPoint, pDist);
279 arrDist.AddPoint(npoints, &pDist);
32438f4e 280
281 ++npoints;
526ddf0e 282 }
32438f4e 283
284 arrUdist.SetUniqueID(npoints);
285 arrDist.SetUniqueID(npoints);
286}
526ddf0e 287
32438f4e 288//________________________________________________________________
289void AliToyMCEventGenerator::SetPoint(Float_t xyz[3], AliTrackPoint &point)
290{
291 //
292 // make AliTrackPoint out of AliTPCclusterMI
293 //
294
295 //covariance at the local frame
296 //assume 1mm distortion in y and z
297 Int_t i[3]={0,0,0};
298 const Double_t kSigmaY=0.1;
299 const Double_t kSigmaZ=0.1;
300 Float_t cov[6]={0,0,0, kSigmaY*kSigmaY,0,kSigmaZ*kSigmaZ};
301
302 const Float_t alpha = -TMath::ATan2(xyz[1],xyz[0]);
303 const Float_t sin = TMath::Sin(alpha), cos = TMath::Cos(alpha);
304
305 Float_t newcov[6];
306 newcov[0] = cov[0]*cos*cos + 2*cov[1]*sin*cos + cov[3]*sin*sin;
307 newcov[1] = cov[1]*(cos*cos-sin*sin) + (cov[3]-cov[0])*sin*cos;
308 newcov[2] = cov[2]*cos + cov[4]*sin;
309 newcov[3] = cov[0]*sin*sin - 2*cov[1]*sin*cos + cov[3]*cos*cos;
310 newcov[4] = cov[4]*cos - cov[2]*sin;
311 newcov[5] = cov[5];
312
313 // voluem ID to add later ....
314 point.SetXYZ(xyz);
315 point.SetCov(newcov);
d92b8630 316 // abuse volume ID for the sector number
32438f4e 317 point.SetVolumeID(fTPCParam->Transform0to1(xyz,i));
b0c26b55 318
32438f4e 319 // TODO: Add sampled dE/dx (use SetCharge)
320}
b0c26b55 321
32438f4e 322//________________________________________________________________
323void AliToyMCEventGenerator::ConvertTrackPointsToLocalClusters(AliTrackPointArray &arrPoints,
324 AliToyMCTrack &tr, Double_t t0, Int_t type)
325{
326 //
327 //
328 //
526ddf0e 329
526ddf0e 330
32438f4e 331 const Int_t npoints=Int_t(arrPoints.GetUniqueID());
332 Int_t secOld=-1;
333
334 // create an array for graphs which are used for local interpolation
335 // we need a new graph if the sector changed since then the local frame will also change
336 TObjArray arrGraphsXY(72);
337 arrGraphsXY.SetOwner();
338 TObjArray arrGraphsXZ(72);
339 arrGraphsXZ.SetOwner();
340 AliTrackPoint p;
341 //create initial graph
342 TGraph *grXY=0x0;
343 TGraph *grXZ=0x0;
344 //row -> sector mapping
345 Int_t rowMap[159];
346 for (Int_t irow=0; irow<159; ++irow) rowMap[irow]=-1;
347
348 // 1. Step
349 // Make from the list of global space points graphs in the local frame
350 // one graph per sector is needed
351 for (Int_t ipoint=0; ipoint<npoints; ++ipoint){
352 arrPoints.GetPoint(p,ipoint);
353 Float_t xyz[3] = {p.GetX(),p.GetY(),p.GetZ()};
354 Int_t index[3] = {0,0,0};
355 Int_t row = fTPCParam->GetPadRow(xyz,index);
356 // rotate space point to local frame
357 // the angle is given by the VolumeID which was set in CrateSpacePoints
358 const Int_t sec=p.GetVolumeID();
359 if (row<0 || (sec<36 && row>62) || row>95 ) continue;
360 Double_t angle=((sec%18)*20.+10.)/TMath::RadToDeg();
361 AliTrackPoint pRot=p.Rotate(angle);
362 const Int_t secrow=row+(sec>35)*63;
363 if (rowMap[secrow]==-1) rowMap[secrow]=sec;
364 // check if we need a new graph (sector change)
365 if (secOld!=sec){
366 grXY=new TGraph;
367 grXZ=new TGraph;
368 arrGraphsXY.AddAt(grXY,sec);
369 arrGraphsXZ.AddAt(grXZ,sec);
370 }
b0c26b55 371
32438f4e 372 //add coordinates in local frame for later interpolation
373 grXY->SetPoint(grXY->GetN(), pRot.GetX(), pRot.GetY());
374 grXZ->SetPoint(grXZ->GetN(), pRot.GetX(), pRot.GetZ());
375 secOld=sec;
376 }
b0c26b55 377
32438f4e 378 // 2. Step
379 // create in the center of each row a space point by using the graph to interpolate
380 // the the center of the row. This is done in xy and xz
381 TSpline3 *splXY=0x0;
382 TSpline3 *splXZ=0x0;
383 AliTPCclusterMI tempCl;
384 secOld=-1;
385 for (Int_t irow=0; irow<159; ++irow ){
386 const Int_t sec = rowMap[irow];
387 if (sec==-1) continue;
388 const Int_t secrow = irow<63?irow:irow-63;
389 Double_t localX = fTPCParam->GetPadRowRadii(sec,secrow);
390 // get graph for the current row
391 if (sec!=secOld){
392 delete splXY;
393 splXY=0x0;
394 delete splXZ;
395 splXZ=0x0;
396
397 grXY=(TGraph*)arrGraphsXY.At(sec);
398 grXZ=(TGraph*)arrGraphsXZ.At(sec);
399 if (!grXY) continue;
b0c26b55 400
05da1b4e 401// if(grXY->GetN()>1 && grXZ->GetN()>1) { //causes segmentation violation if N==1
402// splXY=new TSpline3("splXY",grXY);
403// splXZ=new TSpline3("splXZ",grXZ);
404// }
405// else {
544076ea 406 //TODO: make a cluster also in the sector w only one space point?
05da1b4e 407// continue;
544076ea 408 // Double_t tempX=0., tempY = 0., tempZ = 0.;
409
410 // grXY->GetPoint(0,tempX,localY);
411 // grXZ->GetPoint(0,tempX,localZ);
05da1b4e 412// }
544076ea 413
32438f4e 414 }
415 secOld=sec;
b0c26b55 416
32438f4e 417 // check we are in an active area
05da1b4e 418// if (splXY->FindX(localX)<1 || splXZ->FindX(localX)<1) continue;
419 if ( localX<grXY->GetX()[0] || localX>grXY->GetX()[grXY->GetN()-1] || localX<grXZ->GetX()[0] || localX>grXZ->GetX()[grXZ->GetN()-1]) continue;
b0c26b55 420
32438f4e 421 //get interpolated value at the center for the pad row
422 // using splines
05da1b4e 423// const Double_t localY=splXY->Eval(localX/*,0x0,"S"*/);
424// const Double_t localZ=splXZ->Eval(localX/*,0x0,"S"*/);
425 const Double_t localY=grXY->Eval(localX/*,0x0,"S"*/);
426 const Double_t localZ=grXZ->Eval(localX/*,0x0,"S"*/);
32438f4e 427 Float_t xyz[3]={localX,localY,localZ};
526ddf0e 428
32438f4e 429 if (!SetupCluster(tempCl,xyz,sec,t0)) continue;
d1cf83f5 430 tempCl.SetLabel(tr.GetUniqueID(), 0);
32438f4e 431
432 if (type==0) tr.AddSpacePoint(tempCl);
433 else tr.AddDistortedSpacePoint(tempCl);
434// printf("SetupCluster %3d: (%.2f, %.2f, %.2f), %d, %.2f\n",irow,xyz[0],xyz[1],xyz[2],sec,t0);
435 }
526ddf0e 436
32438f4e 437 delete splXY;
438 splXY=0x0;
439 delete splXZ;
440 splXZ=0x0;
441
442}
526ddf0e 443
32438f4e 444//________________________________________________________________
445Bool_t AliToyMCEventGenerator::SetupCluster(AliTPCclusterMI &tempCl, Float_t xyz[3], Int_t sec, Double_t t0)
446{
447 //
448 //
449 //
526ddf0e 450
32438f4e 451 const Double_t kSigmaY = 0.1;
452 const Double_t kSigmaZ = 0.1;
453 const Double_t kMaxZ0 = fTPCParam->GetZLength();
454 //TODO: Get this from the OCDB at some point?
455 const Double_t kDriftVel = fTPCParam->GetDriftV();
526ddf0e 456
32438f4e 457 tempCl.SetX(xyz[0]);
458 tempCl.SetY(xyz[1]);
459 tempCl.SetZ(xyz[2]);
526ddf0e 460
32438f4e 461 tempCl.SetSigmaY2(kSigmaY*kSigmaY);
462 tempCl.SetSigmaZ2(kSigmaZ*kSigmaZ);
463
464 // transform from the local coordinates to the coordinates expressed in pad coordinates
465 Int_t index[3] = {0,sec,0};
466 fTPCParam->Transform2to3(xyz,index);
467 fTPCParam->Transform3to4(xyz,index);
468 fTPCParam->Transform4to8(xyz,index);
469
470 const Int_t row = index[2];
471 const Int_t nPads = fTPCParam->GetNPads(sec, row);
472 // pad is fractional, but it needs to be shifted from the center
473 // to the edge of the row
474 const Float_t pad = xyz[1] + nPads/2;
475
476 tempCl.SetRow(row);
477 tempCl.SetPad(pad);
478 Float_t timeBin=Float_t(t0 + (kMaxZ0-TMath::Abs(tempCl.GetZ()))/kDriftVel);
479 tempCl.SetTimeBin(timeBin); // set time as t0 + drift time from dist z
480 tempCl.SetDetector(sec);
481
482 //check if we are in the active area
483 if (pad<0 || pad>=nPads) return kFALSE;
484
485 return kTRUE;
526ddf0e 486}
32438f4e 487
a1a695e5 488//________________________________________________________________
489Bool_t AliToyMCEventGenerator::ConnectOutputFile()
490{
491 //
492 // Create the output file name and tree and connect the event
493 //
494
495 fOutFile = new TFile(fOutputFileName.Data(),"recreate");
496
497 if (!fOutFile || !fOutFile->IsOpen()){
498 delete fOutFile;
499 fOutFile=0x0;
500 return kFALSE;
501 }
502
503 fOutTree = new TTree("toyMCtree","Tree with toyMC simulation");
504 fOutTree->Branch("event","AliToyMCEvent",&fEvent);
505
506 return kTRUE;
507}
508
509//________________________________________________________________
510Bool_t AliToyMCEventGenerator::CloseOutputFile()
511{
512 //
513 // close the output file
514 //
515 if (!fOutFile) return kFALSE;
516 fOutFile->Write();
517 fOutFile->Close();
518 delete fOutFile;
519 fOutFile=0x0;
520
521 return kTRUE;
522}
523
524//________________________________________________________________
525void AliToyMCEventGenerator::FillTree()
526{
527 // fill the tree
528 if (fOutTree) fOutTree->Fill();
529}
530
cd8ed0ac 531//________________________________________________________________
223d9e38 532void AliToyMCEventGenerator::SetSpaceCharge(EEpsilon epsilon, EGasType gasType/*=kNeCO2_9010*/,
533 ECollRate collRate/*=k50kHz*/, ECorrection corrType/*=kLookup*/)
cd8ed0ac 534{
535 //
536 // Set the space charge conditions
537 //
223d9e38 538 fCorrectionFile="$ALICE_ROOT/TPC/Calib/maps/SC";
cd8ed0ac 539 switch (gasType) {
540 case kNeCO2_9010:
223d9e38 541 fCorrectionFile.Append("_NeCO2");
cd8ed0ac 542 break;
543 }
544 switch (epsilon) {
545 case kEps5:
223d9e38 546 fCorrectionFile.Append("_eps5");
cd8ed0ac 547 break;
548 case kEps10:
223d9e38 549 fCorrectionFile.Append("_eps10");
cd8ed0ac 550 break;
551 case kEps20:
223d9e38 552 fCorrectionFile.Append("_eps20");
cd8ed0ac 553 break;
554 }
555 switch (collRate) {
556 case k50kHz:
223d9e38 557 fCorrectionFile.Append("_50kHz");
558 break;
559 }
560 switch (corrType) {
561 case kLookup:
562 fCorrectionFile.Append("_precal.lookup.root");
563 break;
564 case kSpaceChargeFile:
565 fCorrectionFile.Append("_precal.root");
cd8ed0ac 566 break;
567 }
cd8ed0ac 568}
569
570//________________________________________________________________
571void AliToyMCEventGenerator::InitSpaceCharge()
572{
573 //
574 // init the space charge conditions
575 // this should be called after the tree was connected
576 //
577
223d9e38 578 AliInfo(Form("Using space charge map file: '%s'",fCorrectionFile.Data()));
579
580 TString corrName("map");
581
582 // allow for specifying an object name for the AliTPCCorrection in the file name
583 // separated by a ':'
584 TObjArray *arr=fCorrectionFile.Tokenize(":");
585 if (arr->GetEntriesFast()>1) {
586 fCorrectionFile=arr->At(0)->GetName();
587 corrName=arr->At(1)->GetName();
588 }
589 delete arr;
590
cd8ed0ac 591
223d9e38 592 TFile f(fCorrectionFile.Data());
593 fTPCCorrection=(AliTPCSpaceCharge3D*)f.Get("map");
cd8ed0ac 594
595 if (fOutTree){
596 AliInfo("Attaching space charge map file name to the tree");
223d9e38 597 fOutTree->GetUserInfo()->Add(new TObjString(fCorrectionFile.Data()));
cd8ed0ac 598 }
599}