]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/vertexingHF/AliAnalysisTaskSEImproveITS.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Protection added for empty fTree->GetTree() (H.Beck)
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskSEImproveITS.cxx
CommitLineData
2d4517ec 1/*************************************************************************
2* Copyright(c) 1998-2011, 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 <TObjArray.h>
17#include <TClonesArray.h>
18#include <TGraph.h>
19#include <TFile.h>
20#include <TList.h>
21#include <TNtuple.h>
10371f90 22
23#include "AliVertex.h"
24#include "AliVVertex.h"
25#include "AliESDVertex.h"
26#include "AliVertexerTracks.h"
27#include "AliAODEvent.h"
28#include "AliAODTrack.h"
29#include "AliAODMCParticle.h"
30#include "AliExternalTrackParam.h"
04de9699 31#include "AliAODRecoDecayHF2Prong.h"
10371f90 32#include "AliAODRecoDecayHF3Prong.h"
2d4517ec 33#include "AliAnalysisTaskSEImproveITS.h"
34
35//
36// Implementation of the "hybrid-approach" for ITS upgrade studies.
37// The tastk smears the track parameters according to estimations
38// from single-track upgrade studies. Afterwards it recalculates
39// the parameters of the reconstructed decays.
40//
41// WARNING: This will affect all tasks in a train after this one
42// (which is typically desired, though).
43//
44
45AliAnalysisTaskSEImproveITS::AliAnalysisTaskSEImproveITS()
46 :AliAnalysisTaskSE(),
47 fD0ZResPCur (0),
48 fD0ZResKCur (0),
49 fD0ZResPiCur (0),
50 fD0RPResPCur (0),
51 fD0RPResKCur (0),
52 fD0RPResPiCur(0),
53 fPt1ResPCur (0),
54 fPt1ResKCur (0),
55 fPt1ResPiCur (0),
56 fD0ZResPUpg (0),
57 fD0ZResKUpg (0),
58 fD0ZResPiUpg (0),
59 fD0RPResPUpg (0),
60 fD0RPResKUpg (0),
61 fD0RPResPiUpg(0),
62 fPt1ResPUpg (0),
63 fPt1ResKUpg (0),
64 fPt1ResPiUpg (0),
b6427aac 65 fRunInVertexing(kFALSE),
47736b6e 66 fImproveTracks(kTRUE),
2d4517ec 67 fDebugOutput (0),
68 fDebugNtuple (0),
69 fDebugVars (0),
70 fNDebug (0)
71{
72 //
73 // Default constructor.
74 //
75}
76
77AliAnalysisTaskSEImproveITS::AliAnalysisTaskSEImproveITS(const char *name,
78 const char *resfileCurURI,
79 const char *resfileUpgURI,
b6427aac 80 Bool_t isRunInVertexing,
81 Int_t ndebug)
2d4517ec 82 :AliAnalysisTaskSE(name),
83 fD0ZResPCur (0),
84 fD0ZResKCur (0),
85 fD0ZResPiCur (0),
86 fD0RPResPCur (0),
87 fD0RPResKCur (0),
88 fD0RPResPiCur(0),
89 fPt1ResPCur (0),
90 fPt1ResKCur (0),
91 fPt1ResPiCur (0),
92 fD0ZResPUpg (0),
93 fD0ZResKUpg (0),
94 fD0ZResPiUpg (0),
95 fD0RPResPUpg (0),
96 fD0RPResKUpg (0),
97 fD0RPResPiUpg(0),
98 fPt1ResPUpg (0),
99 fPt1ResKUpg (0),
100 fPt1ResPiUpg (0),
b6427aac 101 fRunInVertexing(isRunInVertexing),
47736b6e 102 fImproveTracks(kTRUE),
2d4517ec 103 fDebugOutput (0),
104 fDebugNtuple (0),
105 fDebugVars (0),
106 fNDebug (ndebug)
107{
108 //
109 // Constructor to be used to create the task.
110 // The the URIs specify the resolution files to be used.
111 // They are expected to contain TGraphs with the resolutions
112 // for the current and the upgraded ITS (see code for details).
113 // One may also specify for how many tracks debug information
114 // is written to the output.
115 //
116 TFile *resfileCur=TFile::Open(resfileCurURI);
117 fD0RPResPCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0RPResP" )));
118 fD0RPResKCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0RPResK" )));
119 fD0RPResPiCur=new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0RPResPi")));
120 fD0ZResPCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0ZResP" )));
121 fD0ZResKCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0ZResK" )));
122 fD0ZResPiCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("D0ZResPi" )));
123 fPt1ResPCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("Pt1ResP" )));
124 fPt1ResKCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("Pt1ResK" )));
125 fPt1ResPiCur =new TGraph(*static_cast<TGraph*>(resfileCur->Get("Pt1ResPi" )));
126 fD0RPResPCur ->SetName("D0RPResPCur" );
127 fD0RPResKCur ->SetName("D0RPResKCur" );
128 fD0RPResPiCur->SetName("D0RPResPiCur");
129 fD0ZResPCur ->SetName("D0ZResPCur" );
130 fD0ZResKCur ->SetName("D0ZResKCur" );
131 fD0ZResPiCur ->SetName("D0ZResPiCur" );
132 fPt1ResPCur ->SetName("Pt1ResPCur" );
133 fPt1ResKCur ->SetName("Pt1ResKCur" );
134 fPt1ResPiCur ->SetName("Pt1ResPiCur" );
135 delete resfileCur;
136 TFile *resfileUpg=TFile::Open(resfileUpgURI );
137 fD0RPResPUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0RPResP" )));
138 fD0RPResKUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0RPResK" )));
139 fD0RPResPiUpg=new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0RPResPi")));
140 fD0ZResPUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0ZResP" )));
141 fD0ZResKUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0ZResK" )));
142 fD0ZResPiUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("D0ZResPi" )));
143 fPt1ResPUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("Pt1ResP" )));
144 fPt1ResKUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("Pt1ResK" )));
145 fPt1ResPiUpg =new TGraph(*static_cast<TGraph*>(resfileUpg->Get("Pt1ResPi" )));
146 fD0RPResPUpg ->SetName("D0RPResPUpg" );
147 fD0RPResKUpg ->SetName("D0RPResKUpg" );
148 fD0RPResPiUpg->SetName("D0RPResPiUpg");
149 fD0ZResPUpg ->SetName("D0ZResPUpg" );
150 fD0ZResKUpg ->SetName("D0ZResKUpg" );
151 fD0ZResPiUpg ->SetName("D0ZResPiUpg" );
152 fPt1ResPUpg ->SetName("Pt1ResPUpg" );
153 fPt1ResKUpg ->SetName("Pt1ResKUpg" );
154 fPt1ResPiUpg ->SetName("Pt1ResPiUpg" );
155 delete resfileUpg;
156
157 DefineOutput(1,TNtuple::Class());
158}
159
160AliAnalysisTaskSEImproveITS::~AliAnalysisTaskSEImproveITS() {
161 //
162 // Destructor.
163 //
164 delete fDebugOutput;
165}
166
167void AliAnalysisTaskSEImproveITS::UserCreateOutputObjects() {
168 //
169 // Creation of user output objects.
170 //
171 fDebugOutput=new TList();
172 fDebugOutput->SetOwner();
173 fDebugNtuple=new TNtuple("fDebugNtuple","Smearing","pdg:ptmc:d0rpo:d0zo:pt1o:sd0rpo:sd0zo:spt1o:d0rpn:d0zn:pt1n:sd0rpn:sd0zn:spt1n:d0rpmc:d0zmc:pt1mc");
174 fDebugVars=new Float_t[fDebugNtuple->GetNvar()];
175
176 fDebugOutput->Add(fDebugNtuple );
177
178 fDebugOutput->Add(fD0RPResPCur );
179 fDebugOutput->Add(fD0RPResKCur );
180 fDebugOutput->Add(fD0RPResPiCur);
181 fDebugOutput->Add(fD0ZResPCur );
182 fDebugOutput->Add(fD0ZResKCur );
183 fDebugOutput->Add(fD0ZResPiCur );
184 fDebugOutput->Add(fPt1ResPCur );
185 fDebugOutput->Add(fPt1ResKCur );
186 fDebugOutput->Add(fPt1ResPiCur );
187 fDebugOutput->Add(fD0RPResPUpg );
188 fDebugOutput->Add(fD0RPResKUpg );
189 fDebugOutput->Add(fD0RPResPiUpg);
190 fDebugOutput->Add(fD0ZResPUpg );
191 fDebugOutput->Add(fD0ZResKUpg );
192 fDebugOutput->Add(fD0ZResPiUpg );
193 fDebugOutput->Add(fPt1ResPUpg );
194 fDebugOutput->Add(fPt1ResKUpg );
195 fDebugOutput->Add(fPt1ResPiUpg );
196
197 PostData(1,fDebugOutput);
198}
199
200void AliAnalysisTaskSEImproveITS::UserExec(Option_t*) {
201 //
202 // The event loop
203 //
b6427aac 204 AliAODEvent *ev=0x0;
205 if(!fRunInVertexing) {
206 ev=dynamic_cast<AliAODEvent*>(InputEvent());
207 } else {
208 if(AODEvent() && IsStandardAOD()) ev = dynamic_cast<AliAODEvent*> (AODEvent());
209 }
a828d135 210 if(!ev) return;
2d4517ec 211 Double_t bz=ev->GetMagneticField();
212
318e7f18 213
214
215
2d4517ec 216 // Smear all tracks
217 TClonesArray *mcs=static_cast<TClonesArray*>(ev->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
218 if (!mcs) return;
47736b6e 219 if (fImproveTracks) {
220 for(Int_t itrack=0;itrack<ev->GetNumberOfTracks();++itrack) {
221 SmearTrack(ev->GetTrack(itrack),mcs);
222 }
223 }
2d4517ec 224
225 // TODO: recalculated primary vertex
226 AliVVertex *primaryVertex=ev->GetPrimaryVertex();
227
228 // Recalculate all candidates
04de9699 229 // D0->Kpi
230 TClonesArray *array2Prong=static_cast<TClonesArray*>(ev->GetList()->FindObject("D0toKpi"));
231 if (array2Prong) {
232 for (Int_t icand=0;icand<array2Prong->GetEntries();++icand) {
233 AliAODRecoDecayHF2Prong *decay=static_cast<AliAODRecoDecayHF2Prong*>(array2Prong->At(icand));
234
235 // recalculate vertices
236 AliVVertex *oldSecondaryVertex=decay->GetSecondaryVtx();
237
238
239 AliExternalTrackParam et1; et1.CopyFromVTrack(static_cast<AliAODTrack*>(decay->GetDaughter(0)));
240 AliExternalTrackParam et2; et2.CopyFromVTrack(static_cast<AliAODTrack*>(decay->GetDaughter(1)));
241
242 TObjArray ta12;
243
244 ta12.Add(&et1); ta12.Add(&et2);
245 AliESDVertex *v12 =RecalculateVertex(oldSecondaryVertex,&ta12 ,bz);
246
247
248 // update secondary vertex
249 Double_t pos[3];
250 v12->GetXYZ(pos);
251
252 decay->GetSecondaryVtx()->SetPosition(pos[0],pos[1],pos[2]);
253 decay->GetSecondaryVtx()->SetChi2perNDF(v12->GetChi2toNDF());
254
255 //!!!!TODO: covariance matrix
256
257 // update d0
318e7f18 258 Double_t d0z0[2],covd0z0[3];
04de9699 259 Double_t d0[2],d0err[2];
260 et1.PropagateToDCA(primaryVertex,bz,100.,d0z0,covd0z0);
261 d0[0]=d0z0[0];
262 d0err[0] = TMath::Sqrt(covd0z0[0]);
263 et2.PropagateToDCA(primaryVertex,bz,100.,d0z0,covd0z0);
264 d0[1]=d0z0[0];
265 d0err[1] = TMath::Sqrt(covd0z0[0]);
266 decay->Setd0Prongs(2,d0);
267 decay->Setd0errProngs(2,d0err);
268 //
269
270
271 Double_t xdummy=0.,ydummy=0.;
272 Double_t dca;
273 dca=et1.GetDCA(&et2,bz,xdummy,ydummy);
274 decay->SetDCA(dca);
275
276
277 delete v12;
278
279 Double_t px[2],py[2],pz[2];
280 for (Int_t i=0;i<2;++i) {
281 const AliAODTrack *t=static_cast<AliAODTrack*>(decay->GetDaughter(i));
282 px[i]=t->Px();
283 py[i]=t->Py();
284 pz[i]=t->Pz();
285 }
286 decay->SetPxPyPzProngs(2,px,py,pz);
287 }
288 }
289
290
291
292 // Three prong
2d4517ec 293 TClonesArray *array3Prong=static_cast<TClonesArray*>(ev->GetList()->FindObject("Charm3Prong"));
294 if (array3Prong) {
295 for (Int_t icand=0;icand<array3Prong->GetEntries();++icand) {
296 AliAODRecoDecayHF3Prong *decay=static_cast<AliAODRecoDecayHF3Prong*>(array3Prong->At(icand));
297
298 // recalculate vertices
299 AliVVertex *oldSecondaryVertex=decay->GetSecondaryVtx();
4ec1ca0f 300 AliExternalTrackParam et1; et1.CopyFromVTrack(static_cast<AliAODTrack*>(decay->GetDaughter(0)));
301 AliExternalTrackParam et2; et2.CopyFromVTrack(static_cast<AliAODTrack*>(decay->GetDaughter(1)));
302 AliExternalTrackParam et3; et3.CopyFromVTrack(static_cast<AliAODTrack*>(decay->GetDaughter(2)));
2d4517ec 303 TObjArray ta123,ta12,ta23;
304 ta123.Add(&et1);ta123.Add(&et2);ta123.Add(&et3);
305 ta12. Add(&et1);ta12 .Add(&et2);
306 ta23 .Add(&et2);ta23 .Add(&et3);
307 AliESDVertex *v123=RecalculateVertex(oldSecondaryVertex,&ta123,bz);
308 AliESDVertex *v12 =RecalculateVertex(oldSecondaryVertex,&ta12 ,bz);
309 AliESDVertex *v23 =RecalculateVertex(oldSecondaryVertex,&ta23 ,bz);
310
311 // update secondary vertex
312 Double_t pos[3];
313 v123->GetXYZ(pos);
314 decay->GetSecondaryVtx()->SetPosition(pos[0],pos[1],pos[2]);
315 decay->GetSecondaryVtx()->SetChi2perNDF(v123->GetChi2toNDF());
316 //TODO: covariance matrix
317
318 // update d0 for all progs
319 Double_t d0z0[2],covd0z0[3];
320 Double_t d0[3],d0err[3];
321 et1.PropagateToDCA(primaryVertex,bz,100.,d0z0,covd0z0);
322 d0[0]=d0z0[0];
323 d0err[0] = TMath::Sqrt(covd0z0[0]);
324 et2.PropagateToDCA(primaryVertex,bz,100.,d0z0,covd0z0);
325 d0[1]=d0z0[0];
326 d0err[1] = TMath::Sqrt(covd0z0[0]);
327 et3.PropagateToDCA(primaryVertex,bz,100.,d0z0,covd0z0);
328 d0[2]=d0z0[0];
329 d0err[2] = TMath::Sqrt(covd0z0[0]);
330 decay->Setd0Prongs (3,d0 );
331 decay->Setd0errProngs(3,d0err);
332 // TODO: setter missing
333
334 // update dca for prong combinations
335 Double_t xdummy=0.,ydummy=0.;
336 Double_t dca[3];
337 dca[0]=et1.GetDCA(&et2,bz,xdummy,ydummy);
338 dca[1]=et3.GetDCA(&et2,bz,xdummy,ydummy);
339 dca[2]=et1.GetDCA(&et3,bz,xdummy,ydummy);
340 decay->SetDCAs(3,dca);
341
342 // update dist12 and dist23
343 primaryVertex->GetXYZ(pos);
344 decay->SetDist12toPrim(TMath::Sqrt((v12->GetX()-pos[0])*(v12->GetX()-pos[0])
345 +(v12->GetY()-pos[1])*(v12->GetY()-pos[1])
346 +(v12->GetZ()-pos[2])*(v12->GetZ()-pos[2])));
347 decay->SetDist23toPrim(TMath::Sqrt((v23->GetX()-pos[0])*(v23->GetX()-pos[0])
348 +(v23->GetY()-pos[1])*(v23->GetY()-pos[1])
349 +(v23->GetZ()-pos[2])*(v23->GetZ()-pos[2])));
350
351 delete v123;delete v12;delete v23;
352
353 Double_t px[3],py[3],pz[3];
354 for (Int_t i=0;i<3;++i) {
355 const AliAODTrack *t=static_cast<AliAODTrack*>(decay->GetDaughter(i));
356 px[i]=t->Px();
357 py[i]=t->Py();
358 pz[i]=t->Pz();
359 }
360 decay->SetPxPyPzProngs(3,px,py,pz);
361 }
362 }
363}
364
365void AliAnalysisTaskSEImproveITS::SmearTrack(AliAODTrack *track,const TClonesArray *mcs) {
366 // Early exit, if this track has nothing in common with the ITS
318e7f18 367
2d4517ec 368 if (!(track->HasPointOnITSLayer(0) || track->HasPointOnITSLayer(1)))
369 return;
370
47736b6e 371 // Check if the track was already "improved" (this is done with a trick using layer 7 (ie the 8th))
372 if (TESTBIT(track->GetITSClusterMap(),7)) return;
373 //
374
375
2d4517ec 376 // Get reconstructed track parameters
4ec1ca0f 377 AliExternalTrackParam et; et.CopyFromVTrack(track);
2d4517ec 378 Double_t *param=const_cast<Double_t*>(et.GetParameter());
379//TODO: Double_t *covar=const_cast<Double_t*>(et.GetCovariance());
380
381 // Get MC info
382 Int_t imc=track->GetLabel();
383 if (imc<=0) return;
384 const AliAODMCParticle *mc=static_cast<AliAODMCParticle*>(mcs->At(imc));
385 Double_t mcx[3];
386 Double_t mcp[3];
387 Double_t mccv[36]={0.};
388 Short_t mcc;
389 mc->XvYvZv(mcx);
390 mc->PxPyPz(mcp);
391 mcc=mc->Charge();
392 AliExternalTrackParam mct(mcx,mcp,mccv,mcc);
393 const Double_t *parammc=mct.GetParameter();
394//TODO: const Double_t *covermc=mct.GetCovariance();
395 AliVertex vtx(mcx,1.,1);
396
397 // Correct reference points and frames according to MC
398 // TODO: B-Field correct?
399 // TODO: failing propagation....
400 et.PropagateToDCA(&vtx,track->GetBz(),10.);
401 et.Rotate(mct.GetAlpha());
402
403 // Select appropriate smearing functions
404 Double_t ptmc=TMath::Abs(mc->Pt());
405 Double_t sd0rpn=0.;
406 Double_t sd0zn =0.;
407 Double_t spt1n =0.;
408 Double_t sd0rpo=0.;
409 Double_t sd0zo =0.;
410 Double_t spt1o =0.;
411 switch (mc->GetPdgCode()) {
412 case 2212: case -2212:
413 sd0rpo=EvalGraph(fD0RPResPCur,ptmc);
414 sd0zo =EvalGraph(fD0ZResPCur ,ptmc);
415 spt1o =EvalGraph(fPt1ResPCur ,ptmc);
416 sd0rpn=EvalGraph(fD0RPResPUpg,ptmc);
417 sd0zn =EvalGraph(fD0ZResPUpg ,ptmc);
418 spt1n =EvalGraph(fPt1ResPUpg ,ptmc);
419 break;
420 case 321: case -321:
421 sd0rpo=EvalGraph(fD0RPResKCur,ptmc);
422 sd0zo =EvalGraph(fD0ZResKCur ,ptmc);
423 spt1o =EvalGraph(fPt1ResKCur ,ptmc);
424 sd0rpn=EvalGraph(fD0RPResKUpg,ptmc);
425 sd0zn =EvalGraph(fD0ZResKUpg ,ptmc);
426 spt1n =EvalGraph(fPt1ResKUpg ,ptmc);
427 break;
428 case 211: case -211:
429 sd0rpo=EvalGraph(fD0RPResPiCur,ptmc);
430 sd0zo =EvalGraph(fD0ZResPiCur ,ptmc);
431 spt1o =EvalGraph(fPt1ResPiCur ,ptmc);
432 sd0rpn=EvalGraph(fD0RPResPiUpg,ptmc);
433 sd0zn =EvalGraph(fD0ZResPiUpg ,ptmc);
434 spt1n =EvalGraph(fPt1ResPiUpg ,ptmc);
435 break;
436 default:
437 return;
438 }
439
440 // Use the same units (i.e. cm and GeV/c)! TODO: pt!
441 sd0rpo*=1.e-4;
442 sd0zo *=1.e-4;
443 sd0rpn*=1.e-4;
444 sd0zn *=1.e-4;
445
446 // Apply the smearing
447 Double_t d0zo =param [1];
448 Double_t d0zmc =parammc[1];
449 Double_t d0rpo =param [0];
450 Double_t d0rpmc=parammc[0];
451 Double_t pt1o =param [4];
452 Double_t pt1mc =parammc[4];
453 Double_t dd0zo =d0zo-d0zmc;
454 Double_t dd0zn =dd0zo *(sd0zo >0. ? (sd0zn /sd0zo ) : 1.);
455 Double_t d0zn =d0zmc+dd0zn;
456 Double_t dd0rpo=d0rpo-d0rpmc;
457 Double_t dd0rpn=dd0rpo*(sd0rpo>0. ? (sd0rpn/sd0rpo) : 1.);
458 Double_t d0rpn =d0rpmc+dd0rpn;
459 Double_t dpt1o =pt1o-pt1mc;
460 Double_t dpt1n =dpt1o *(spt1o >0. ? (spt1n /spt1o ) : 1.);
461 Double_t pt1n =pt1mc+dpt1n;
462 param[0]=d0rpn;
463 param[1]=d0zn ;
464 param[4]=pt1n ;
465
466 // Copy the smeared parameters to the AOD track
467 Double_t x[3];
468 Double_t p[3];
469 et.GetXYZ(x);
470 et.GetPxPyPz(p);
471 track->SetPosition(x,kFALSE);
472 track->SetP(p,kTRUE);
473
47736b6e 474
475 // Mark the track as "improved" with a trick (this is done with a trick using layer 7 (ie the 8th))
476 UChar_t itsClusterMap = track->GetITSClusterMap();
477 SETBIT(itsClusterMap,7);
478 track->SetITSClusterMap(itsClusterMap);
479 //
480
2d4517ec 481 // write out debug infos
482 if (fDebugNtuple->GetEntries()<fNDebug) {
483 Int_t idbg=0;
484 fDebugVars[idbg++]=mc->GetPdgCode();
485 fDebugVars[idbg++]=ptmc ;
486 fDebugVars[idbg++]=d0rpo ;
487 fDebugVars[idbg++]=d0zo ;
488 fDebugVars[idbg++]=pt1o ;
489 fDebugVars[idbg++]=sd0rpo;
490 fDebugVars[idbg++]=sd0zo ;
491 fDebugVars[idbg++]=spt1o ;
492 fDebugVars[idbg++]=d0rpn ;
493 fDebugVars[idbg++]=d0zn ;
494 fDebugVars[idbg++]=pt1n ;
495 fDebugVars[idbg++]=sd0rpn;
496 fDebugVars[idbg++]=sd0zn ;
497 fDebugVars[idbg++]=spt1n ;
498 fDebugVars[idbg++]=d0rpmc;
499 fDebugVars[idbg++]=d0zmc ;
500 fDebugVars[idbg++]=pt1mc ;
501 fDebugNtuple->Fill(fDebugVars);
502 PostData(1,fDebugOutput);
503 }
504}
505
506AliESDVertex* AliAnalysisTaskSEImproveITS::RecalculateVertex(const AliVVertex *old,TObjArray *tracks,Double_t bField) {
507 //
508 // Helper function to recalculate a vertex.
509 //
510
511 static UShort_t ids[]={1,2,3}; //TODO: unsave...
512 AliVertexerTracks vertexer(bField);
513 vertexer.SetVtxStart(old->GetX(),old->GetY(),old->GetZ());
514 AliESDVertex *vertex=vertexer.VertexForSelectedTracks(tracks,ids);
515 return vertex;
516}
517
518Double_t AliAnalysisTaskSEImproveITS::EvalGraph(const TGraph *graph,Double_t x) const {
519 //
520 // Evaluates a TGraph without linear extrapolation. Instead the last
521 // valid point of the graph is used when out of range.
522 // The function assumes an ascending order of X.
523 //
524
525 // TODO: find a pretty solution for this:
526 Int_t n =graph->GetN();
527 Double_t xmin=graph->GetX()[0 ];
528 Double_t xmax=graph->GetX()[n-1];
529 if (x<xmin) return graph->Eval(xmin);
530 if (x>xmax) return graph->Eval(xmax);
531 return graph->Eval(x);
532}
533
534ClassImp(AliAnalysisTaskSEImproveITS);
535