New directory for the VMC tests (Ivana, Eva)
[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDmultiplicity.cxx
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
17 ////////////////////////////////////////////////////////////////////////////
18 //                                                                        //
19 //  TRD multiplicity                                                      //
20 //                                                                        //
21 // Task to select true single tracks
22 // 
23 // Program flow:
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
27 // 
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
34 //                                                                        //
35 //  Authors: Yvonne C. Pachmayer <pachmay@physi.uni-heidelberg.de>        //
36 //                                                                        //
37 ////////////////////////////////////////////////////////////////////////////
38
39 #include <TObjArray.h>
40 #include <TH1F.h>
41 #include <TTreeStream.h>
42
43 #include "AliESDtrack.h"
44 #include "AliExternalTrackParam.h"
45 #include "AliTRDseedV1.h"
46 #include "AliTRDtrackV1.h"
47 #include "AliTRDtrackerV1.h"
48 #include "AliTrackPointArray.h"
49
50 #include "AliTRDmultiplicity.h"
51 #include "AliTRDgeometry.h"
52 #include "info/AliTRDtrackInfo.h"
53
54 ClassImp(AliTRDmultiplicity)
55
56 //____________________________________________________________________
57 AliTRDmultiplicity::AliTRDmultiplicity()
58     :AliTRDrecoTask("Multiplicity", "Barrel Tracking Multiplicity")
59   ,fEventCounter(0)
60 {
61   //
62   // Default constructor
63     //
64 }
65
66 //____________________________________________________________________
67 AliTRDmultiplicity::~AliTRDmultiplicity()
68 {
69 }
70
71 //____________________________________________________________________
72 void  AliTRDmultiplicity::CreateOutputObjects()
73 {
74   //
75   // Create output objects
76   //
77
78   OpenFile(0, "RECREATE");
79
80   TH1 *h = 0x0;
81   fContainer = new TObjArray();
82   for(Int_t is=0; is<AliTRDgeometry::kNsector; is++){
83   fContainer->Add(h = new TH1F(Form("h_sector_%i", is), " ", 100,-10,10));
84   }
85
86
87 //____________________________________________________________________
88 void AliTRDmultiplicity::Exec(Option_t *)
89 {
90   //
91   // Do it
92   //
93
94
95
96   ULong_t status;
97   AliTRDtrackInfo     *track = 0x0;
98   AliExternalTrackParam *esd = 0x0;
99   AliTRDtrackV1 *trackTRD = 0x0;
100   Double_t x_anode[AliTRDgeometry::kNlayer] = {300.2, 312.8, 325.4, 338.0, 350.6, 363.2}; // Take the default X0
101   Int_t standAlone=1;
102   fEventCounter++;
103
104
105   for(Int_t itrk=0; itrk<fTracks->GetEntriesFast(); itrk++){
106     track = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
107     if(!track->HasESDtrack()) continue;
108     status = track->GetStatus();
109
110     AliTrackPoint points[6];
111     Float_t xyz[3];
112     memset(xyz, 0, sizeof(Float_t) * 3);
113     Float_t trackletX[6];
114     Float_t trackletY[6];
115     Float_t trackletZ[6];
116     for(Int_t a=0;a<6;a++)
117     {
118         xyz[0] = x_anode[a];
119         points[a].SetXYZ(xyz);
120         trackletX[a]=-1000;
121         trackletY[a]=-1000;
122         trackletZ[a]=-1000;
123     }
124     Int_t detTracklet=600;
125
126
127
128     // TRD standalone
129     if(((status&AliESDtrack::kTRDout)>0) && !((status&AliESDtrack::kTRDin)>0))
130     {
131         trackTRD = track->GetTrack();
132
133         if(trackTRD)
134         {
135             AliTRDseedV1 *tracklet = 0x0;
136             Int_t counterTracklet=0;
137             for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++)
138             {
139                 tracklet = trackTRD->GetTracklet(itl);
140                 if(tracklet)
141                 {
142                     if(tracklet->IsOK())
143                     {
144                         counterTracklet++;
145                         detTracklet=tracklet->GetDetector();
146                         //printf("%i %f %f \n",itl,tracklet->GetYfit(0),tracklet->GetZfit(0));
147                         trackletX[itl]=tracklet->GetX0();
148                         trackletY[itl]=tracklet->GetYfit(0);
149                         trackletZ[itl]=tracklet->GetZfit(0);
150                     }
151                 }
152             }
153             // this cut is needed otherwise the order of tracklets in the fit function is not correct
154             if(counterTracklet==AliTRDgeometry::kNlayer) AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(trackTRD), 0x0, kTRUE, counterTracklet, points);
155             else continue;
156
157             if(DebugLevel()>=1){
158               for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
159                 //printf("---------------------------------------- %i %i %f %f %f \n",counterTracklet,il,points[il].GetX(),points[il].GetY(),points[il].GetZ());
160                 Double_t pointX=points[il].GetX();
161                 Double_t pointY=points[il].GetY();
162                 Double_t pointZ=points[il].GetZ();
163
164
165                 (*DebugStream())   << "TrackletsinTRD"
166                     << "standalone=" << standAlone
167                     << "eventcounter=" << fEventCounter
168                     << "layer="     << il
169                     << "dettracklet=" << detTracklet
170                     << "xtracklet=" << trackletX[il]
171                     << "xtrack="    << pointX
172                     << "ytracklet=" << trackletY[il]
173                     << "ytrack="    << pointY
174                     << "ztracklet=" << trackletZ[il]
175                     << "ztrack="    << pointZ
176                     << "num_tracklets=" << counterTracklet
177                     << "\n";
178
179               }
180             }
181         }
182     } // end TRD standalone
183
184
185
186     // TPC cluster selection for cosmic data
187     if((track->GetTPCncls())<40) continue;
188     // missing TPC propagation
189     if(!(status&AliESDtrack::kTPCout)) continue;
190     standAlone=0;
191
192     esd = track->GetESDinfo()->GetOuterParam();
193     // calculate sector - does not matter if track ends in TPC or TRD - sector number independent
194     Double_t alpha;
195     alpha=-100;
196     Int_t fSector=100;
197     if(esd){
198         alpha=esd->GetAlpha();
199         fSector = (Int_t)((alpha+ TMath::Pi())/(20*TMath::Pi()/180));
200         if((fSector>-1) && (fSector<18))
201         {
202         }
203         else
204         {
205             continue;
206         }
207     }
208
209
210
211     // only these sectors have a TRD detector at the moment
212    if(fSector==0||fSector==8||fSector==9||fSector==17)
213     {
214         trackTRD = track->GetTrack();
215         if(!trackTRD) continue;
216         AliTRDseedV1 *tracklet = 0x0;
217         Int_t counterTracklet=0;
218
219
220         for(Int_t itl = 0; itl < AliTRDgeometry::kNlayer; itl++){
221             tracklet = trackTRD->GetTracklet(itl);
222             if(!tracklet || !tracklet->IsOK()) continue;
223             counterTracklet++;
224             detTracklet=tracklet->GetDetector();
225             trackletX[itl]=tracklet->GetX0();
226             trackletY[itl]=tracklet->GetYfit(0);
227             trackletZ[itl]=tracklet->GetZfit(0);
228
229             /*
230             AliTRDcluster *c = 0x0;
231             for(Int_t itime = 0; itime < AliTRDseedV1::kNtb; itime++){
232                 c = tracklet->GetClusters(itime);
233                 if(!c) continue;
234 //                Float_t cluster_x = c->GetX();
235 //                Bool_t isprop;
236 //                isprop=kFALSE;
237 //                isprop=AliTracker::PropagateTrackTo(esd,(Double_t)cluster_x,0.105,5,kFALSE);
238 //                if(isprop)
239 //                {
240 //                    Int_t detector = c->GetDetector();
241 //                    ((TH1F*)fContainer->At(fSector))->Fill((c->GetY())-(esd->GetY()));
242 //                    if(c->GetY()-esd->GetY()>);
243 //                    printf("diff: %f\n",c->GetY()-(esd->GetY()));
244 //                }
245             }
246             */
247
248         } // loop over tracklets ends
249
250         if(counterTracklet==AliTRDgeometry::kNlayer)
251         {
252             AliTRDtrackerV1::FitRiemanTilt(const_cast<AliTRDtrackV1 *>(trackTRD), 0x0, kTRUE, counterTracklet, points);
253         }
254         else continue;
255
256        
257
258         if(DebugLevel()>=1){
259           for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){
260             //printf("---------------------------------------- %i %i %f %f %f \n",counterTracklet,il,points[il].GetX(),points[il].GetY(),points[il].GetZ());
261             Double_t pointX=points[il].GetX();
262             Double_t pointY=points[il].GetY();
263             Double_t pointZ=points[il].GetZ();
264             
265             (*DebugStream())   << "TrackletsinTRD"
266                 << "standalone=" << standAlone
267                 << "eventcounter=" << fEventCounter
268                 << "layer="     << il
269                 << "dettracklet=" << detTracklet
270                 << "xtracklet=" << trackletX[il]
271                 << "xtrack="    << pointX
272                 << "ytracklet=" << trackletY[il]
273                 << "ytrack="    << pointY
274                 << "ztracklet=" << trackletZ[il]
275                 << "ztrack="    << pointZ
276                 << "num_tracklets=" << counterTracklet
277                 << "\n";
278           }
279         }
280
281     }
282   }
283
284   
285   PostData(0, fContainer);
286 }
287
288