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 **************************************************************************/
17 ////////////////////////////////////////////////////////////////////////////
19 // TRD multiplicity //
21 // Task to select true single tracks
24 // Part A - task AliTRDmultiplicity (within the framework qaRec):
25 // For TRD standalone or TRD&TOF tracks I make a FitRiemanTilt
26 // The resulting parameterization is dumped into a DebugStreamer written to a file
28 // Part B – $ALICE_ROOT/TRD/qaRec/macros/TrackletsinTRD.C[h] (analysis macro):
29 // The clusters are read in
30 // The above produced file is read in
31 // I make a fit through the parameterization of the FitRiemanTilt -> in order to get a straight line (representing the track)
32 // the distance of each cluster with respect to the straight line is calculated
33 // If there is no cluster at a distance of 0.6 to 2 with respect to the track, the track is defined as a good track, i.e. clean single track
35 // Authors: Yvonne C. Pachmayer <pachmay@physi.uni-heidelberg.de> //
37 ////////////////////////////////////////////////////////////////////////////
39 #include <TObjArray.h>
41 #include <TTreeStream.h>
43 #include "AliESDtrack.h"
44 #include "AliExternalTrackParam.h"
45 #include "AliTRDseedV1.h"
46 #include "AliTRDtrackV1.h"
47 #include "AliTRDtrackerV1.h"
48 #include "AliTrackPointArray.h"
50 #include "AliTRDmultiplicity.h"
51 #include "AliTRDgeometry.h"
52 #include "info/AliTRDtrackInfo.h"
54 ClassImp(AliTRDmultiplicity)
56 //____________________________________________________________________
57 AliTRDmultiplicity::AliTRDmultiplicity()
62 // Default constructor
66 AliTRDmultiplicity::AliTRDmultiplicity(char* name)
67 :AliTRDrecoTask(name, "Barrel Tracking Multiplicity")
71 // Default constructor
75 //____________________________________________________________________
76 AliTRDmultiplicity::~AliTRDmultiplicity()
80 //____________________________________________________________________
81 void AliTRDmultiplicity::UserCreateOutputObjects()
84 // Create output objects
88 fContainer = new TObjArray();
89 for(Int_t is=0; is<AliTRDgeometry::kNsector; is++){
90 fContainer->Add(h = new TH1F(Form("h_sector_%i", is), " ", 100,-10,10));
94 //____________________________________________________________________
95 void AliTRDmultiplicity::UserExec(Option_t *)
104 AliTRDtrackInfo *track = 0x0;
105 AliExternalTrackParam *esd = 0x0;
106 AliTRDtrackV1 *trackTRD = 0x0;
107 Double_t x_anode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
112 for(Int_t itrk=0; itrk<fTracks->GetEntriesFast(); itrk++){
113 track = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
114 if(!track->HasESDtrack()) continue;
115 status = track->GetStatus();
117 AliTrackPoint points[6];
119 memset(xyz, 0, sizeof(Float_t) * 3);
120 Float_t trackletX[6];
121 Float_t trackletY[6];
122 Float_t trackletZ[6];
123 for(Int_t a=0;a<6;a++)
126 points[a].SetXYZ(xyz);
131 Int_t detTracklet=600;
136 if(((status&AliESDtrack::kTRDout)>0) && !((status&AliESDtrack::kTRDin)>0))
138 trackTRD = track->GetTrack();
142 AliTRDseedV1 *tracklet = 0x0;
143 Int_t counterTracklet=0;
144 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++)
146 tracklet = trackTRD->GetTracklet(itl);
152 detTracklet=tracklet->GetDetector();
153 //printf("%i %f %f \n",itl,tracklet->GetYfit(0),tracklet->GetZfit(0));
154 trackletX[itl]=tracklet->GetX0();
155 trackletY[itl]=tracklet->GetYfit(0);
156 trackletZ[itl]=tracklet->GetZfit(0);
160 // this cut is needed otherwise the order of tracklets in the fit function is not correct
161 if(counterTracklet==AliTRDgeometry::kNlayer) AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(trackTRD), 0x0, kTRUE, counterTracklet, points);
165 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
166 //printf("---------------------------------------- %i %i %f %f %f \n",counterTracklet,il,points[il].GetX(),points[il].GetY(),points[il].GetZ());
167 Double_t pointX=points[il].GetX();
168 Double_t pointY=points[il].GetY();
169 Double_t pointZ=points[il].GetZ();
172 (*DebugStream()) << "TrackletsinTRD"
173 << "standalone=" << standAlone
174 << "eventcounter=" << fEventCounter
176 << "dettracklet=" << detTracklet
177 << "xtracklet=" << trackletX[il]
178 << "xtrack=" << pointX
179 << "ytracklet=" << trackletY[il]
180 << "ytrack=" << pointY
181 << "ztracklet=" << trackletZ[il]
182 << "ztrack=" << pointZ
183 << "num_tracklets=" << counterTracklet
189 } // end TRD standalone
193 // TPC cluster selection for cosmic data
194 if((track->GetTPCncls())<40) continue;
195 // missing TPC propagation
196 if(!(status&AliESDtrack::kTPCout)) continue;
199 esd = track->GetESDinfo()->GetOuterParam();
200 // calculate sector - does not matter if track ends in TPC or TRD - sector number independent
205 alpha=esd->GetAlpha();
206 fSector = (Int_t)((alpha+ TMath::Pi())/(20*TMath::Pi()/180));
207 if((fSector>-1) && (fSector<18))
218 // only these sectors have a TRD detector at the moment
219 if(fSector==0||fSector==8||fSector==9||fSector==17)
221 trackTRD = track->GetTrack();
222 if(!trackTRD) continue;
223 AliTRDseedV1 *tracklet = 0x0;
224 Int_t counterTracklet=0;
227 for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
228 tracklet = trackTRD->GetTracklet(itl);
229 if(!tracklet || !tracklet->IsOK()) continue;
231 detTracklet=tracklet->GetDetector();
232 trackletX[itl]=tracklet->GetX0();
233 trackletY[itl]=tracklet->GetYfit(0);
234 trackletZ[itl]=tracklet->GetZfit(0);
237 AliTRDcluster *c = 0x0;
238 for(Int_t itime = 0; itime < AliTRDseedV1::kNtb; itime++){
239 c = tracklet->GetClusters(itime);
241 // Float_t cluster_x = c->GetX();
244 // isprop=AliTracker::PropagateTrackTo(esd,(Double_t)cluster_x,0.105,5,kFALSE);
247 // Int_t detector = c->GetDetector();
248 // ((TH1F*)fContainer->At(fSector))->Fill((c->GetY())-(esd->GetY()));
249 // if(c->GetY()-esd->GetY()>);
250 // printf("diff: %f\n",c->GetY()-(esd->GetY()));
255 } // loop over tracklets ends
257 if(counterTracklet==AliTRDgeometry::kNlayer)
259 AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(trackTRD), 0x0, kTRUE, counterTracklet, points);
266 for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
267 //printf("---------------------------------------- %i %i %f %f %f \n",counterTracklet,il,points[il].GetX(),points[il].GetY(),points[il].GetZ());
268 Double_t pointX=points[il].GetX();
269 Double_t pointY=points[il].GetY();
270 Double_t pointZ=points[il].GetZ();
272 (*DebugStream()) << "TrackletsinTRD"
273 << "standalone=" << standAlone
274 << "eventcounter=" << fEventCounter
276 << "dettracklet=" << detTracklet
277 << "xtracklet=" << trackletX[il]
278 << "xtrack=" << pointX
279 << "ytracklet=" << trackletY[il]
280 << "ytrack=" << pointY
281 << "ztracklet=" << trackletZ[il]
282 << "ztrack=" << pointZ
283 << "num_tracklets=" << counterTracklet
292 PostData(1, fContainer);