New version containing the files for tracking
[u/mrichter/AliRoot.git] / ITS / AliITS.cxx
CommitLineData
4c039060 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$Log$
55b5a0b3 18Revision 1.25 2000/11/12 22:38:05 barbera
19Added header file for the SPD Bari model
20
dfe2a93f 21Revision 1.24 2000/10/09 22:18:12 barbera
22Bug fixes from MAriana to le AliITStest.C run correctly
23
65d4384f 24Revision 1.23 2000/10/05 20:47:42 nilsen
25fixed dependencies of include files. Tryed but failed to get a root automaticly
26generates streamer function to work. Modified SetDefaults.
27
143d1056 28Revision 1.9.2.15 2000/10/04 16:56:40 nilsen
29Needed to include stdlib.h
30
31=======
32Revision 1.22 2000/10/04 19:45:52 barbera
33Corrected by F. Carminati for v3.04
34
114da3a4 35Revision 1.21 2000/10/02 21:28:08 fca
36Removal of useless dependecies via forward declarations
37
94de3818 38Revision 1.20 2000/10/02 16:31:39 barbera
39General code clean-up
40
d5da1ecf 41Revision 1.9.2.14 2000/10/02 15:43:51 barbera
42General code clean-up (e.g., printf -> cout)
43
44Revision 1.19 2000/09/22 12:13:25 nilsen
45Patches and updates for fixes to this and other routines.
46
9c74c52b 47Revision 1.18 2000/07/12 05:32:20 fca
48Correcting several syntax problem with static members
49
b669392e 50Revision 1.17 2000/07/10 16:07:18 fca
51Release version of ITS code
52
3bd79107 53Revision 1.9.2.3 2000/02/02 13:42:09 barbera
54fixed AliITS.cxx for new AliRun structure. Added ITS hits list to list of hits which will have their track numbers updated
55
56Revision 1.9.2.2 2000/01/23 03:03:13 nilsen
57//fixed FillModule. Removed fi(fabs(xl)<dx....
58
59Revision 1.9.2.1 2000/01/12 19:03:32 nilsen
60This is the version of the files after the merging done in December 1999.
61See the ReadMe110100.txt file for details
c9a71be1 62
1cedd08a 63Revision 1.9 1999/11/14 14:33:25 fca
64Correct problems with distructors and pointers, thanks to I.Hrivnacova
65
6c854497 66Revision 1.8 1999/09/29 09:24:19 fca
67Introduction of the Copyright and cvs Log
68
4c039060 69*/
70
fe4da5cc 71///////////////////////////////////////////////////////////////////////////////
58005f18 72//
73// An overview of the basic philosophy of the ITS code development
74// and analysis is show in the figure below.
fe4da5cc 75//Begin_Html
76/*
a92b2b7d 77<img src="picts/ITS/ITS_Analysis_schema.gif">
fe4da5cc 78</pre>
79<br clear=left>
80<font size=+2 color=red>
58005f18 81<p>Roberto Barbera is in charge of the ITS Offline code (1999).
fe4da5cc 82<a href="mailto:roberto.barbera@ct.infn.it">Roberto Barbera</a>.
83</font>
84<pre>
85*/
86//End_Html
58005f18 87//
88// AliITS. Inner Traking System base class.
89// This class contains the base procedures for the Inner Tracking System
90//
91//Begin_Html
92/*
a92b2b7d 93<img src="picts/ITS/AliITS_Class_Diagram.gif">
58005f18 94</pre>
95<br clear=left>
96<font size=+2 color=red>
97<p>This show the class diagram of the different elements that are part of
98the AliITS class.
99</font>
100<pre>
101*/
102//End_Html
103//
104// Version: 0
105// Written by Rene Brun, Federico Carminati, and Roberto Barbera
106//
107// Version: 1
108// Modified and documented by Bjorn S. Nilsen
109// July 11 1999
110//
3bd79107 111// Version: 2
112// Modified and documented by A. Bologna
113// October 18 1999
114//
58005f18 115// AliITS is the general base class for the ITS. Also see AliDetector for
116// futher information.
117//
fe4da5cc 118///////////////////////////////////////////////////////////////////////////////
143d1056 119#include <stdlib.h>
fe4da5cc 120#include <TMath.h>
121#include <TRandom.h>
9c74c52b 122#include <TBranch.h>
fe4da5cc 123#include <TVector.h>
3bd79107 124#include <TObjArray.h>
3bd79107 125#include <TROOT.h>
126#include <TObjectTable.h>
9c74c52b 127#include <TFile.h>
128#include <TTree.h>
d5da1ecf 129#include <TString.h>
3bd79107 130
131
bb8ac1e0 132
fe4da5cc 133#include "AliRun.h"
3bd79107 134#include "AliITS.h"
135#include "AliITSMap.h"
e8189707 136#include "AliITSDetType.h"
3bd79107 137#include "AliITSClusterFinder.h"
138#include "AliITSsimulation.h"
143d1056 139#include "AliITSresponse.h"
3bd79107 140#include "AliITSsegmentationSPD.h"
141#include "AliITSresponseSPD.h"
dfe2a93f 142#include "AliITSresponseSPDbari.h"
3bd79107 143#include "AliITSsegmentationSDD.h"
144#include "AliITSresponseSDD.h"
145#include "AliITSsegmentationSSD.h"
146#include "AliITSresponseSSD.h"
9c74c52b 147#include "AliITShit.h"
148#include "AliITSgeom.h"
149#include "AliITSdigit.h"
150#include "AliITSmodule.h"
151#include "AliITSRecPoint.h"
152#include "AliITSRawCluster.h"
94de3818 153#include "AliMC.h"
114da3a4 154#include "stdlib.h"
e8189707 155
55b5a0b3 156#include "AliITStrack.h"
157#include "AliITSiotrack.h"
158#include "AliITStracking.h"
159#include "../TPC/AliTPC.h"
160#include "../TPC/AliTPCParam.h"
161
b669392e 162const Int_t AliITS::fgkNTYPES=3;
c9a71be1 163
3bd79107 164ClassImp(AliITS)
165
fe4da5cc 166//_____________________________________________________________________________
3bd79107 167AliITS::AliITS() : AliDetector() {
fe4da5cc 168 //
169 // Default initialiser for ITS
58005f18 170 // The default constructor of the AliITS class. In addition to
171 // creating the AliITS class it zeros the variables fIshunt (a member
172 // of AliDetector class), fEuclidOut, and fIdN, and zeros the pointers
3bd79107 173 // fITSpoints, fIdSens, and fIdName. The AliDetector default constructor
174 // is also called.
fe4da5cc 175 //
3bd79107 176
9c74c52b 177
58005f18 178 fIshunt = 0;
179 fEuclidOut = 0;
e8189707 180
181 //fNDetTypes = fgkNTYPES;
58005f18 182 fIdN = 0;
183 fIdName = 0;
184 fIdSens = 0;
6c854497 185 fITSmodules = 0;
3bd79107 186 //
187 fDetTypes = 0;
3bd79107 188 //
189 fDtype = 0;
190 fNdtype = 0;
191 fCtype = 0;
192 fNctype = 0;
193 fRecPoints = 0;
194 fNRecPoints = 0;
e8189707 195 fTreeC = 0;
3bd79107 196 //
197 fITSgeom=0;
fe4da5cc 198}
3bd79107 199
fe4da5cc 200//_____________________________________________________________________________
58005f18 201AliITS::AliITS(const char *name, const char *title):AliDetector(name,title){
fe4da5cc 202 //
203 // Default initialiser for ITS
58005f18 204 // The constructor of the AliITS class. In addition to creating the
205 // AliITS class, it allocates memory for the TClonesArrays fHits and
206 // fDigits, and for the TObjArray fITSpoints. It also zeros the variables
207 // fIshunt (a member of AliDetector class), fEuclidOut, and fIdN, and zeros
208 // the pointers fIdSens and fIdName. To help in displaying hits via the ROOT
209 // macro display.C AliITS also sets the marker color to red. The variables
210 // passes with this constructor, const char *name and *title, are used by
211 // the constructor of AliDetector class. See AliDetector class for a
212 // description of these parameters and its constructor functions.
fe4da5cc 213 //
58005f18 214
9c74c52b 215
58005f18 216 fHits = new TClonesArray("AliITShit", 1560);
1cedd08a 217 gAlice->AddHitList(fHits);
3bd79107 218
e8189707 219 //fNDetTypes = fgkNTYPES;
220
221 fNdtype = new Int_t[fgkNTYPES];
222 fDtype = new TObjArray(fgkNTYPES);
223
224 fNctype = new Int_t[fgkNTYPES];
225 fCtype = new TObjArray(fgkNTYPES);
3bd79107 226
3bd79107 227
228 fRecPoints = 0;
229 fNRecPoints = 0;
230
e8189707 231 fTreeC = 0;
3bd79107 232
233 fITSmodules = 0;
58005f18 234
235 fIshunt = 0;
236 fEuclidOut = 0;
237 fIdN = 0;
238 fIdName = 0;
239 fIdSens = 0;
3bd79107 240
e8189707 241 fDetTypes = new TObjArray(fgkNTYPES);
3bd79107 242
243 Int_t i;
e8189707 244 for(i=0;i<fgkNTYPES;i++) {
3bd79107 245 (*fDetTypes)[i]=new AliITSDetType();
e8189707 246 fNdtype[i]=0;
247 fNctype[i]=0;
3bd79107 248 }
249 //
58005f18 250
fe4da5cc 251 SetMarkerColor(kRed);
58005f18 252
3bd79107 253 fITSgeom=0;
254}
255//___________________________________________________________________________
256AliITS::AliITS(AliITS &source){
e8189707 257 // copy constructor
3bd79107 258 if(this==&source) return;
d5da1ecf 259 Error("AliITS::Copy constructor",
260 "You are not allowed to make a copy of the AliITS");
3bd79107 261 exit(1);
fe4da5cc 262}
3bd79107 263//____________________________________________________________________________
264AliITS& AliITS::operator=(AliITS &source){
e8189707 265 // assignment operator
3bd79107 266 if(this==&source) return *this;
d5da1ecf 267 Error("AliITS::operator=",
268 "You are not allowed to make a copy of the AliITS");
3bd79107 269 exit(1);
e8189707 270 return *this; //fake return
3bd79107 271}
272//____________________________________________________________________________
273void AliITS::ClearModules(){
274 //clear the modules TObjArray
3bd79107 275
9c74c52b 276 if(fITSmodules) fITSmodules->Delete();
fe4da5cc 277
3bd79107 278}
fe4da5cc 279//_____________________________________________________________________________
58005f18 280AliITS::~AliITS(){
fe4da5cc 281 //
282 // Default distructor for ITS
58005f18 283 // The default destructor of the AliITS class. In addition to deleting
284 // the AliITS class it deletes the memory pointed to by the fHits, fDigits,
285 // fIdSens, fIdName, and fITSpoints.
fe4da5cc 286 //
3bd79107 287
288
fe4da5cc 289 delete fHits;
290 delete fDigits;
3bd79107 291 delete fRecPoints;
143d1056 292// delete fIdName; // TObjArray of TObjStrings
293 if(fIdName!=0) delete[] fIdName; // Array of TStrings
58005f18 294 if(fIdSens!=0) delete[] fIdSens;
3bd79107 295 if(fITSmodules!=0) {
296 this->ClearModules();
297 delete fITSmodules;
298 }// end if fITSmodules!=0
fe4da5cc 299
fe4da5cc 300 //
3bd79107 301 if(fDtype) {
9c74c52b 302 fDtype->Delete();
303 delete fDtype;
3bd79107 304 }
9c74c52b 305 delete [] fNdtype;
306 if (fCtype) {
307 fCtype->Delete();
308 delete fCtype;
3bd79107 309 }
9c74c52b 310 delete [] fNctype;
fe4da5cc 311 //
3bd79107 312
313 if (fDetTypes) {
314 fDetTypes->Delete();
315 delete fDetTypes;
316 }
317
318 if (fTreeC) delete fTreeC;
e8189707 319
9c74c52b 320 if (fITSgeom) delete fITSgeom;
321
3bd79107 322}
323
324//___________________________________________
325AliITSDetType* AliITS::DetType(Int_t id)
326{
327 //return pointer to id detector type
328 return ((AliITSDetType*) (*fDetTypes)[id]);
329
330}
331//___________________________________________
e8189707 332void AliITS::SetClasses(Int_t id, const char *digit, const char *cluster)
3bd79107 333{
334 //set the digit and cluster classes to be used for the id detector type
335 ((AliITSDetType*) (*fDetTypes)[id])->ClassNames(digit,cluster);
336
337}
338//___________________________________________
339void AliITS::SetResponseModel(Int_t id, AliITSresponse *response)
340{
341 //set the response model for the id detector type
342
343 ((AliITSDetType*) (*fDetTypes)[id])->ResponseModel(response);
344
fe4da5cc 345}
346
3bd79107 347//___________________________________________
348void AliITS::SetSegmentationModel(Int_t id, AliITSsegmentation *seg)
349{
350 //set the segmentation model for the id detector type
351
352 ((AliITSDetType*) (*fDetTypes)[id])->SegmentationModel(seg);
353
354}
355
356//___________________________________________
357void AliITS::SetSimulationModel(Int_t id, AliITSsimulation *sim)
358{
359 //set the simulation model for the id detector type
360
361 ((AliITSDetType*) (*fDetTypes)[id])->SimulationModel(sim);
362
363}
364//___________________________________________
365void AliITS::SetReconstructionModel(Int_t id, AliITSClusterFinder *reconst)
366{
367 //set the cluster finder model for the id detector type
368
369 ((AliITSDetType*) (*fDetTypes)[id])->ReconstructionModel(reconst);
58005f18 370
58005f18 371}
372
fe4da5cc 373//_____________________________________________________________________________
58005f18 374void AliITS::AddHit(Int_t track, Int_t *vol, Float_t *hits){
fe4da5cc 375 //
376 // Add an ITS hit
58005f18 377 // The function to add information to the AliITShit class. See the
378 // AliITShit class for a full description. This function allocates the
379 // necessary new space for the hit information and passes the variable
380 // track, and the pointers *vol and *hits to the AliITShit constructor
381 // function.
fe4da5cc 382 //
383 TClonesArray &lhits = *fHits;
384 new(lhits[fNhits++]) AliITShit(fIshunt,track,vol,hits);
385}
3bd79107 386//_____________________________________________________________________________
387void AliITS::AddRealDigit(Int_t id, Int_t *digits)
388{
389 // add a real digit - as coming from data
390
391 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
392 new(ldigits[fNdtype[id]++]) AliITSdigit(digits);
393
394}
395//_____________________________________________________________________________
e8189707 396void AliITS::AddSimDigit(Int_t id, AliITSdigit *d)
3bd79107 397{
398
399 // add a simulated digit
400
3bd79107 401 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
402
403 switch(id)
404 {
405 case 0:
406 new(ldigits[fNdtype[id]++]) AliITSdigitSPD(*((AliITSdigitSPD*)d));
407 break;
408 case 1:
409 new(ldigits[fNdtype[id]++]) AliITSdigitSDD(*((AliITSdigitSDD*)d));
410 break;
411 case 2:
412 new(ldigits[fNdtype[id]++]) AliITSdigitSSD(*((AliITSdigitSSD*)d));
413 break;
414 }
415
416}
417
418//_____________________________________________________________________________
e8189707 419void AliITS::AddSimDigit(Int_t id,Float_t phys,Int_t *digits,Int_t *tracks,Int_t *hits,Float_t *charges){
3bd79107 420
421 // add a simulated digit to the list
422
423 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
424 switch(id)
425 {
426 case 0:
e8189707 427 new(ldigits[fNdtype[id]++]) AliITSdigitSPD(digits,tracks,hits);
3bd79107 428 break;
429 case 1:
e8189707 430 new(ldigits[fNdtype[id]++]) AliITSdigitSDD(phys,digits,tracks,hits,charges);
3bd79107 431 break;
432 case 2:
e8189707 433 new(ldigits[fNdtype[id]++]) AliITSdigitSSD(digits,tracks,hits);
3bd79107 434 break;
435 }
436
437}
438
439//_____________________________________________________________________________
440void AliITS::AddCluster(Int_t id, AliITSRawCluster *c)
441{
442
443 // add a cluster to the list
444
3bd79107 445 TClonesArray &lcl = *((TClonesArray*)(*fCtype)[id]);
446
447 switch(id)
448 {
449 case 0:
450 new(lcl[fNctype[id]++]) AliITSRawClusterSPD(*((AliITSRawClusterSPD*)c));
451 break;
452 case 1:
453 new(lcl[fNctype[id]++]) AliITSRawClusterSDD(*((AliITSRawClusterSDD*)c));
454 break;
455 case 2:
456 new(lcl[fNctype[id]++]) AliITSRawClusterSSD(*((AliITSRawClusterSSD*)c));
457 break;
458 }
459
460}
461
462
463//_____________________________________________________________________________
464void AliITS::AddRecPoint(const AliITSRecPoint &r)
465{
466 //
467 // Add a reconstructed space point to the list
468 //
469 TClonesArray &lrecp = *fRecPoints;
470 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
471}
472
473
474//____________________________________________
475void AliITS::ResetDigits()
476{
477 //
9c74c52b 478 // Reset number of digits and the digits array for the ITS detector
3bd79107 479 //
480
481 if (!fDtype) return;
482
483 Int_t i;
e8189707 484 for (i=0;i<fgkNTYPES;i++ ) {
3bd79107 485 if ((*fDtype)[i]) ((TClonesArray*)(*fDtype)[i])->Clear();
486 if (fNdtype) fNdtype[i]=0;
487 }
488}
489
490//____________________________________________
491void AliITS::ResetDigits(Int_t i)
492{
493 //
494 // Reset number of digits and the digits array for this branch
495 //
e8189707 496 if ((*fDtype)[i]) ((TClonesArray*)(*fDtype)[i])->Clear();
497 if (fNdtype) fNdtype[i]=0;
3bd79107 498}
e8189707 499
500
3bd79107 501//____________________________________________
502void AliITS::ResetClusters()
503{
504 //
505 // Reset number of clusters and the clusters array for ITS
506 //
507
508 Int_t i;
e8189707 509 for (i=0;i<fgkNTYPES;i++ ) {
3bd79107 510 if ((*fCtype)[i]) ((TClonesArray*)(*fCtype)[i])->Clear();
511 if (fNctype) fNctype[i]=0;
512 }
513
514}
e8189707 515
3bd79107 516//____________________________________________
517void AliITS::ResetClusters(Int_t i)
518{
519 //
520 // Reset number of clusters and the clusters array for this branch
521 //
522 if ((*fCtype)[i]) ((TClonesArray*)(*fCtype)[i])->Clear();
523 if (fNctype) fNctype[i]=0;
524
525}
526
527
528//____________________________________________
529void AliITS::ResetRecPoints()
530{
531 //
532 // Reset number of rec points and the rec points array
533 //
3bd79107 534 if (fRecPoints) fRecPoints->Clear();
9c74c52b 535 fNRecPoints = 0;
3bd79107 536
537}
538
fe4da5cc 539//_____________________________________________________________________________
58005f18 540Int_t AliITS::DistancetoPrimitive(Int_t , Int_t ){
fe4da5cc 541 //
542 // Distance from mouse to ITS on the screen. Dummy routine
58005f18 543 // A dummy routine used by the ROOT macro display.C to allow for the
544 // use of the mouse (pointing device) in the macro. In general this should
545 // never be called. If it is it returns the number 9999 for any value of
546 // x and y.
fe4da5cc 547 //
548 return 9999;
549}
550
551//_____________________________________________________________________________
58005f18 552void AliITS::Init(){
fe4da5cc 553 //
554 // Initialise ITS after it has been built
58005f18 555 // This routine initializes the AliITS class. It is intended to be called
556 // from the Init function in AliITSv?. Besides displaying a banner
3bd79107 557 // indicating that it has been called it initializes the array fIdSens
558 // and sets the default segmentation, response, digit and raw cluster classes
58005f18 559 // Therefore it should be called after a call to CreateGeometry.
fe4da5cc 560 //
561 Int_t i;
143d1056 562
d5da1ecf 563 cout << endl;
564 for(i=0;i<30;i++) cout << "*";cout << " ITS_INIT ";
565 for(i=0;i<30;i++) cout << "*";cout << endl;
143d1056 566//
567 SetDefaults();
568// TObjArray of TObjStrings
569// for(i=0;i<fIdN;i++) fIdSens[i] = gMC->VolId((fIdName->At(i))->GetName());
570// Array of TStrings
58005f18 571 for(i=0;i<fIdN;i++) fIdSens[i] = gMC->VolId(fIdName[i]);
143d1056 572//
d5da1ecf 573 for(i=0;i<70;i++) cout << "*";
574 cout << endl;
fe4da5cc 575}
576
3bd79107 577//_____________________________________________________________________________
578void AliITS::SetDefaults()
579{
580 // sets the default segmentation, response, digit and raw cluster classes
581
65d4384f 582 printf("SetDefaults\n");
3bd79107 583
584 AliITSDetType *iDetType;
585
65d4384f 586
3bd79107 587 //SPD
588
65d4384f 589 iDetType=DetType(0);
143d1056 590 if (!iDetType->GetSegmentationModel()) {
65d4384f 591 AliITSsegmentationSPD *seg0=new AliITSsegmentationSPD(fITSgeom);
592 SetSegmentationModel(0,seg0);
593 }
594 if (!iDetType->GetResponseModel()) {
595 SetResponseModel(0,new AliITSresponseSPD());
596 }
3bd79107 597 // set digit and raw cluster classes to be used
65d4384f 598
599 const char *kData0=(iDetType->GetResponseModel())->DataType();
e8189707 600 if (strstr(kData0,"real")) {
3bd79107 601 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSPD");
602 } else iDetType->ClassNames("AliITSdigitSPD","AliITSRawClusterSPD");
603
65d4384f 604 // SDD //
605 iDetType=DetType(1);
606 if (!iDetType->GetResponseModel()) {
607 SetResponseModel(1,new AliITSresponseSDD());
608 }
609 AliITSresponse *resp1=iDetType->GetResponseModel();
143d1056 610 if (!iDetType->GetSegmentationModel()) {
65d4384f 611 AliITSsegmentationSDD *seg1=new AliITSsegmentationSDD(fITSgeom,resp1);
143d1056 612 SetSegmentationModel(1,seg1);
143d1056 613 }
65d4384f 614 const char *kData1=(iDetType->GetResponseModel())->DataType();
615 const char *kopt=iDetType->GetResponseModel()->ZeroSuppOption();
e8189707 616 if ((!strstr(kopt,"2D")) && (!strstr(kopt,"1D")) || strstr(kData1,"real") ) {
3bd79107 617 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSDD");
618 } else iDetType->ClassNames("AliITSdigitSDD","AliITSRawClusterSDD");
619
620 // SSD
65d4384f 621 iDetType=DetType(2);
143d1056 622 if (!iDetType->GetSegmentationModel()) {
65d4384f 623 AliITSsegmentationSSD *seg2=new AliITSsegmentationSSD(fITSgeom);
143d1056 624 SetSegmentationModel(2,seg2);
65d4384f 625 }
626 if (!iDetType->GetResponseModel()) {
627 SetResponseModel(2,new AliITSresponseSSD());
628 }
629 const char *kData2=(iDetType->GetResponseModel())->DataType();
e8189707 630 if (strstr(kData2,"real")) {
3bd79107 631 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSSD");
632 } else iDetType->ClassNames("AliITSdigitSSD","AliITSRawClusterSSD");
633
e8189707 634 if (fgkNTYPES>3) {
3bd79107 635 Warning("SetDefaults","Only the three basic detector types are initialised!");
65d4384f 636 }
3bd79107 637
638}
65d4384f 639
640
641
e8189707 642//_____________________________________________________________________________
643void AliITS::SetDefaultSimulation()
644{
645 // to be written
646
647}
648//_____________________________________________________________________________
649void AliITS::SetDefaultClusterFinders()
650{
651 // to be written
652
653}
3bd79107 654//_____________________________________________________________________________
655
656void AliITS::MakeTreeC(Option_t *option)
657{
658 // create a separate tree to store the clusters
659
d5da1ecf 660 cout << "AliITS::MakeTreeC" << endl;
9c74c52b 661
3bd79107 662 char *optC = strstr(option,"C");
663 if (optC && !fTreeC) fTreeC = new TTree("TC","Clusters in ITS");
9c74c52b 664 else return;
3bd79107 665
666 Int_t buffersize = 4000;
667 char branchname[30];
668
669 char *det[3] = {"SPD","SDD","SSD"};
670
671 // one branch for Clusters per type of detector
672 Int_t i;
e8189707 673 for (i=0; i<fgkNTYPES ;i++) {
674 if (fgkNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
3bd79107 675 else sprintf(branchname,"%sClusters%d",GetName(),i+1);
676 if (fCtype && fTreeC) {
677 TreeC()->Branch(branchname,&((*fCtype)[i]), buffersize);
d5da1ecf 678 cout << "Making Branch " << branchname;
679 cout << " for Clusters of detector type " << i+1 << endl;
3bd79107 680 }
681 }
682
683}
684
685//_____________________________________________________________________________
686void AliITS::GetTreeC(Int_t event)
687{
688
d5da1ecf 689 cout << "AliITS::GetTreeC" << endl;
9c74c52b 690
3bd79107 691 // get the clusters tree for this event and set the branch address
692 char treeName[20];
693 char branchname[30];
694
695 char *det[3] = {"SPD","SDD","SSD"};
696
697 ResetClusters();
698 if (fTreeC) {
699 delete fTreeC;
700 }
701
702 sprintf(treeName,"TreeC%d",event);
703 fTreeC = (TTree*)gDirectory->Get(treeName);
704
705
706 TBranch *branch;
707 if (fTreeC) {
708 Int_t i;
e8189707 709 for (i=0; i<fgkNTYPES; i++) {
710 if (fgkNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
711 else sprintf(branchname,"%sClusters%d",GetName(),i+1);
712 if (fCtype) {
3bd79107 713 branch = fTreeC->GetBranch(branchname);
714 if (branch) branch->SetAddress(&((*fCtype)[i]));
715 }
716 }
717 } else {
d5da1ecf 718 Error("AliITS::GetTreeC",
719 "cannot find Clusters Tree for event:%d\n",event);
3bd79107 720 }
721
722}
fe4da5cc 723//_____________________________________________________________________________
58005f18 724void AliITS::MakeBranch(Option_t* option){
fe4da5cc 725 //
3bd79107 726 // Creates Tree branches for the ITS.
fe4da5cc 727 //
3bd79107 728
729
fe4da5cc 730 Int_t buffersize = 4000;
3bd79107 731 char branchname[30];
fe4da5cc 732 sprintf(branchname,"%s",GetName());
733
734 AliDetector::MakeBranch(option);
735
fe4da5cc 736
3bd79107 737// one branch for digits per type of detector
738
739 char *det[3] = {"SPD","SDD","SSD"};
740
e8189707 741 char digclass[40];
742 char clclass[40];
3bd79107 743
744 Int_t i;
e8189707 745 for (i=0; i<fgkNTYPES ;i++) {
3bd79107 746 AliITSDetType *iDetType=DetType(i);
e8189707 747 iDetType->GetClassNames(digclass,clclass);
748 //printf("i, digclass, recclass %d %s %s\n",i,digclass,clclass);
3bd79107 749 // digits
e8189707 750 (*fDtype)[i] = new TClonesArray(digclass,10000);
3bd79107 751 // clusters
e8189707 752 (*fCtype)[i] = new TClonesArray(clclass,10000);
3bd79107 753 }
754
755
e8189707 756 for (i=0; i<fgkNTYPES ;i++) {
757 if (fgkNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
3bd79107 758 else sprintf(branchname,"%sDigits%d",GetName(),i+1);
759
760 if (fDtype && gAlice->TreeD()) {
761 gAlice->TreeD()->Branch(branchname,&((*fDtype)[i]), buffersize);
d5da1ecf 762 cout << "Making Branch " << branchname;
763 cout << " for digits of type "<< i+1 << endl;
3bd79107 764 }
765 }
766
767 // only one branch for rec points for all detector types
768 sprintf(branchname,"%sRecPoints",GetName());
769
e8189707 770 fRecPoints=new TClonesArray("AliITSRecPoint",10000);
3bd79107 771
772 if (fRecPoints && gAlice->TreeR()) {
773 gAlice->TreeR()->Branch(branchname,&fRecPoints, buffersize);
d5da1ecf 774 cout << "Making Branch " << branchname;
775 cout << " for reconstructed space points" << endl;
3bd79107 776 }
777
778
779}
780
781//___________________________________________
782void AliITS::SetTreeAddress()
783{
784
785 // Set branch address for the Trees.
786
787 char branchname[30];
788 AliDetector::SetTreeAddress();
789
790 char *det[3] = {"SPD","SDD","SSD"};
791
792 TBranch *branch;
793 TTree *treeD = gAlice->TreeD();
794 TTree *treeR = gAlice->TreeR();
795
796 Int_t i;
797 if (treeD) {
e8189707 798 for (i=0; i<fgkNTYPES; i++) {
799 if (fgkNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
3bd79107 800 else sprintf(branchname,"%sDigits%d",GetName(),i+1);
801 if (fDtype) {
802 branch = treeD->GetBranch(branchname);
803 if (branch) branch->SetAddress(&((*fDtype)[i]));
804 }
805 }
806 }
807
808
809 if (treeR) {
810 sprintf(branchname,"%sRecPoints",GetName());
811 if (fRecPoints) {
812 branch = treeR->GetBranch(branchname);
813 if (branch) branch->SetAddress(&fRecPoints);
814 }
815 }
816
817
818}
819
820//____________________________________________________________________________
821void AliITS::InitModules(Int_t size,Int_t &nmodules){
822
823 //initialize the modules array
824
e8189707 825 if(fITSmodules){
9c74c52b 826 fITSmodules->Delete();
e8189707 827 delete fITSmodules;
828 }
829
3bd79107 830 Int_t nl,indexMAX,index;
3bd79107 831
832 if(size<=0){ // default to using data stored in AliITSgeom
833 if(fITSgeom==0) {
d5da1ecf 834 Error("AliITS::InitModules",
835 "in AliITS::InitModule fITSgeom not defined\n");
3bd79107 836 return;
837 } // end if fITSgeom==0
838 nl = fITSgeom->GetNlayers();
839 indexMAX = fITSgeom->GetModuleIndex(nl,fITSgeom->GetNladders(nl),
840 fITSgeom->GetNdetectors(nl))+1;
841 nmodules = indexMAX;
842 fITSmodules = new TObjArray(indexMAX);
3bd79107 843 for(index=0;index<indexMAX;index++){
3bd79107 844 fITSmodules->AddAt( new AliITSmodule(index),index);
3bd79107 845 } // end for index
846 }else{
847 fITSmodules = new TObjArray(size);
9c74c52b 848 for(index=0;index<size;index++) {
849 fITSmodules->AddAt( new AliITSmodule(index),index);
850 }
851
3bd79107 852 nmodules = size;
853 } // end i size<=0
854}
855
856//____________________________________________________________________________
e8189707 857void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t nmodules,Option_t *option,Text_t *filename){
3bd79107 858
859 // fill the modules with the sorted by module hits; add hits from background
860 // if option=Add
861
9c74c52b 862
3bd79107 863 static TTree *trH1; //Tree with background hits
864 static TClonesArray *fHits2; //List of hits for one track only
865
866 static Bool_t first=kTRUE;
867 static TFile *file;
e8189707 868 char *addBgr = strstr(option,"Add");
3bd79107 869
870
e8189707 871 if (addBgr ) {
3bd79107 872 if(first) {
873 cout<<"filename "<<filename<<endl;
874 file=new TFile(filename);
875 cout<<"I have opened "<<filename<<" file "<<endl;
876 fHits2 = new TClonesArray("AliITShit",1000 );
877 }
878 first=kFALSE;
879 file->cd();
880 file->ls();
881 // Get Hits Tree header from file
882 if(fHits2) fHits2->Clear();
883 if(trH1) delete trH1;
884 trH1=0;
885
886 char treeName[20];
887 sprintf(treeName,"TreeH%d",bgrev);
888 trH1 = (TTree*)gDirectory->Get(treeName);
889 //printf("TrH1 %p of treename %s for event %d \n",trH1,treeName,bgrev);
890
891 if (!trH1) {
d5da1ecf 892 Error("AliITS::FillModules",
893 "cannot find Hits Tree for event:%d\n",bgrev);
3bd79107 894 }
895 // Set branch addresses
896 TBranch *branch;
897 char branchname[20];
898 sprintf(branchname,"%s",GetName());
899 if (trH1 && fHits2) {
900 branch = trH1->GetBranch(branchname);
901 if (branch) branch->SetAddress(&fHits2);
902 }
903
904 // test
905 //Int_t ntracks1 =(Int_t)TrH1->GetEntries();
906 //printf("background - ntracks1 - %d\n",ntracks1);
907 }
908
3bd79107 909 Int_t npart = gAlice->GetEvent(evnt);
910 if(npart<=0) return;
911 TClonesArray *itsHits = this->Hits();
912 Int_t lay,lad,det,index;
913 AliITShit *itsHit=0;
914 AliITSmodule *mod=0;
915
916 TTree *iTH = gAlice->TreeH();
917 Int_t ntracks =(Int_t) iTH->GetEntries();
918
919 Int_t t,h;
920 for(t=0; t<ntracks; t++){
921 gAlice->ResetHits();
922 iTH->GetEvent(t);
923 Int_t nhits = itsHits->GetEntriesFast();
9c74c52b 924 //printf("nhits %d\n",nhits);
3bd79107 925 if (!nhits) continue;
3bd79107 926 for(h=0; h<nhits; h++){
927 itsHit = (AliITShit *)itsHits->UncheckedAt(h);
928 itsHit->GetDetectorID(lay,lad,det);
9c74c52b 929 // temporarily index=det-1 !!!
930 if(fITSgeom) index = fITSgeom->GetModuleIndex(lay,lad,det);
931 else index=det-1;
932 //
3bd79107 933 mod = this->GetModule(index);
9c74c52b 934 mod->AddHit(itsHit,t,h);
3bd79107 935 } // end loop over hits
936 } // end loop over tracks
937
938 // open the file with background
939
e8189707 940 if (addBgr ) {
941 Int_t track,i;
3bd79107 942 ntracks =(Int_t)trH1->GetEntries();
943 //printf("background - ntracks1 %d\n",ntracks);
944 //printf("background - Start loop over tracks \n");
945 // Loop over tracks
946
e8189707 947 for (track=0; track<ntracks; track++) {
3bd79107 948
949 if (fHits2) fHits2->Clear();
950 trH1->GetEvent(track);
951 // Loop over hits
952 for(i=0;i<fHits2->GetEntriesFast();++i) {
953
954 itsHit=(AliITShit*) (*fHits2)[i];
955 itsHit->GetDetectorID(lay,lad,det);
9c74c52b 956 // temporarily index=det-1 !!!
957 if(fITSgeom) index = fITSgeom->GetModuleIndex(lay,lad,det);
958 else index=det-1;
959 //
3bd79107 960 mod = this->GetModule(index);
9c74c52b 961 mod->AddHit(itsHit,track,i);
3bd79107 962 } // end loop over hits
963 } // end loop over tracks
964
965 TTree *fAli=gAlice->TreeK();
966 TFile *fileAli=0;
967
968 if (fAli) fileAli =fAli->GetCurrentFile();
9c74c52b 969 fileAli->cd();
3bd79107 970
971 } // end if add
972
3bd79107 973 //gObjectTable->Print();
974
975}
976
977
978//____________________________________________________________________________
e8189707 979void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *option, Option_t *opt,Text_t *filename)
3bd79107 980{
981 // keep galice.root for signal and name differently the file for
982 // background when add! otherwise the track info for signal will be lost !
983
9c74c52b 984 // the condition below will disappear when the geom class will be
985 // initialised for all versions - for the moment it is only for v5 !
986 // 7 is the SDD beam test version
987 Int_t ver = this->IsVersion();
988 if(ver!=5 && ver!=7) return;
989
3bd79107 990 char *all = strstr(opt,"All");
991 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
3bd79107 992
993 Int_t nmodules;
994 InitModules(size,nmodules);
e8189707 995 FillModules(evNumber,bgrev,nmodules,option,filename);
3bd79107 996
9c74c52b 997 //TBranch *branch;
3bd79107 998 AliITSsimulation* sim;
9c74c52b 999 //TObjArray *branches=gAlice->TreeD()->GetListOfBranches();
3bd79107 1000 AliITSgeom *geom = GetITSgeom();
1001
1002 Int_t id,module;
9c74c52b 1003 Int_t first,last;
e8189707 1004 for (id=0;id<fgkNTYPES;id++) {
3bd79107 1005 if (!all && !det[id]) continue;
9c74c52b 1006 //branch = (TBranch*)branches->UncheckedAt(id);
3bd79107 1007 AliITSDetType *iDetType=DetType(id);
1008 sim = (AliITSsimulation*)iDetType->GetSimulationModel();
1009 if (!sim) {
1010 Error("HitsToDigits","The simulation class was not instantiated!");
1011 exit(1);
e8189707 1012 // or SetDefaultSimulation();
3bd79107 1013 }
9c74c52b 1014 if(geom) {
1015 first = geom->GetStartDet(id);
1016 last = geom->GetLastDet(id);
1017 } else first=last=0;
d5da1ecf 1018 cout << "det type " << id << " first, last "<< first << last << endl;
3bd79107 1019 for(module=first;module<=last;module++) {
1020 AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
1021 sim->DigitiseModule(mod,module,evNumber);
1022 // fills all branches - wasted disk space
1023 gAlice->TreeD()->Fill();
1024 ResetDigits();
1025 // try and fill only the branch
1026 //branch->Fill();
1027 //ResetDigits(id);
1028 } // loop over modules
1029 } // loop over detector types
1030
3bd79107 1031 ClearModules();
1032
e8189707 1033 Int_t nentries=(Int_t)gAlice->TreeD()->GetEntries();
d5da1ecf 1034 cout << "nentries in TreeD" << nentries << endl;
3bd79107 1035
1036 char hname[30];
1037 sprintf(hname,"TreeD%d",evNumber);
1038 gAlice->TreeD()->Write(hname);
1039 // reset tree
1040 gAlice->TreeD()->Reset();
1041
1042}
1043
1044
1045//____________________________________________________________________________
1046void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
1047{
e8189707 1048 // cluster finding and reconstruction of space points
3bd79107 1049
9c74c52b 1050 // the condition below will disappear when the geom class will be
1051 // initialised for all versions - for the moment it is only for v5 !
1052 // 7 is the SDD beam test version
1053 Int_t ver = this->IsVersion();
1054 if(ver!=5 && ver!=7) return;
1055
3bd79107 1056 char *all = strstr(opt,"All");
1057 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
1058
1059 static Bool_t first=kTRUE;
1060 if (first) {
1061 MakeTreeC("C");
1062 first=kFALSE;
1063 }
3bd79107 1064
9c74c52b 1065 TTree *iTC=TreeC();
1066
1067 //TBranch *branch;
3bd79107 1068 AliITSClusterFinder* rec;
1069
9c74c52b 1070 //TObjArray *branches=gAlice->TreeR()->GetListOfBranches();
3bd79107 1071 AliITSgeom *geom = GetITSgeom();
1072
1073 Int_t id,module;
e8189707 1074 for (id=0;id<fgkNTYPES;id++) {
3bd79107 1075 if (!all && !det[id]) continue;
9c74c52b 1076 //branch = (TBranch*)branches->UncheckedAt(id);
3bd79107 1077 AliITSDetType *iDetType=DetType(id);
1078 rec = (AliITSClusterFinder*)iDetType->GetReconstructionModel();
1079 if (!rec) {
1080 Error("DigitsToRecPoints","The cluster finder class was not instantiated!");
1081 exit(1);
e8189707 1082 // or SetDefaultClusterFinders();
3bd79107 1083 }
1084 TClonesArray *itsDigits = this->DigitsAddress(id);
1085
9c74c52b 1086 Int_t first,last;
1087 if(geom) {
1088 first = geom->GetStartDet(id);
1089 last = geom->GetLastDet(id);
1090 } else first=last=0;
1091 //printf("first last %d %d\n",first,last);
3bd79107 1092 for(module=first;module<=last;module++) {
3bd79107 1093 this->ResetDigits();
1094 if (all) gAlice->TreeD()->GetEvent(lastentry+module);
1095 else gAlice->TreeD()->GetEvent(lastentry+(module-first));
1096 Int_t ndigits = itsDigits->GetEntriesFast();
1097 if (ndigits) rec->FindRawClusters();
1098 gAlice->TreeR()->Fill();
1099 ResetRecPoints();
1100 iTC->Fill();
1101 ResetClusters();
1102 // try and fill only the branch
1103 //branch->Fill();
1104 //ResetRecPoints(id);
1105 } // loop over modules
1106 } // loop over detector types
1107
1108
e8189707 1109 Int_t nentries=(Int_t)gAlice->TreeR()->GetEntries();
1110 Int_t ncentries=(Int_t)iTC->GetEntries();
d5da1ecf 1111 cout << " nentries ncentries " << nentries << ncentries << endl;
3bd79107 1112
1113 char hname[30];
1114 sprintf(hname,"TreeR%d",evNumber);
1115 gAlice->TreeR()->Write(hname);
1116 // reset tree
1117 gAlice->TreeR()->Reset();
1118
1119 sprintf(hname,"TreeC%d",evNumber);
1120 iTC->Write(hname);
1121 iTC->Reset();
1122}
1123
1124
1125//____________________________________________________________________________
e8189707 1126void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t size,
3bd79107 1127Option_t *option,Option_t *opt,Text_t *filename)
1128{
1129 // keep galice.root for signal and name differently the file for
1130 // background when add! otherwise the track info for signal will be lost !
1131
9c74c52b 1132
1133 // the condition below will disappear when the geom class will be
1134 // initialised for all versions - for the moment it is only for v5 !
1135 Int_t ver = this->IsVersion();
1136 if(ver!=5) return;
1137
3bd79107 1138 char *all = strstr(opt,"All");
1139 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
1140
1141 Int_t nmodules;
e8189707 1142 InitModules(size,nmodules);
1143 FillModules(evNumber,bgrev,nmodules,option,filename);
3bd79107 1144
1145
e8189707 1146 AliITSsimulation* sim;
3bd79107 1147 AliITSgeom *geom = GetITSgeom();
1148
e8189707 1149 TRandom *random=new TRandom[9];
1150 random[0].SetSeed(111);
1151 random[1].SetSeed(222);
1152 random[2].SetSeed(333);
1153 random[3].SetSeed(444);
1154 random[4].SetSeed(555);
1155 random[5].SetSeed(666);
1156 random[6].SetSeed(777);
1157 random[7].SetSeed(888);
1158 random[8].SetSeed(999);
1159
1160
3bd79107 1161 Int_t id,module;
e8189707 1162 for (id=0;id<fgkNTYPES;id++) {
3bd79107 1163 if (!all && !det[id]) continue;
e8189707 1164 AliITSDetType *iDetType=DetType(id);
1165 sim = (AliITSsimulation*)iDetType->GetSimulationModel();
1166 if (!sim) {
1167 Error("HitsToFastPoints","The simulation class was not instantiated!");
1168 exit(1);
1169 // or SetDefaultSimulation();
1170 }
3bd79107 1171 Int_t first = geom->GetStartDet(id);
1172 Int_t last = geom->GetLastDet(id);
1173 for(module=first;module<=last;module++) {
1174 AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
e8189707 1175 sim->CreateFastRecPoints(mod,module,random);
3bd79107 1176 gAlice->TreeR()->Fill();
1177 ResetRecPoints();
1178 } // loop over modules
1179 } // loop over detector types
1180
1181
1182 ClearModules();
1183
1184 //Int_t nentries=(Int_t)gAlice->TreeR()->GetEntries();
1185
1186 char hname[30];
1187 sprintf(hname,"TreeR%d",evNumber);
1188 gAlice->TreeR()->Write(hname);
1189 // reset tree
1190 gAlice->TreeR()->Reset();
1191
e8189707 1192 delete [] random;
1193
fe4da5cc 1194}
1195
d5da1ecf 1196
e8189707 1197//____________________________________________________________________________
1198void AliITS::Streamer(TBuffer &R__b){
58005f18 1199 // Stream an object of class AliITS.
3bd79107 1200
e8189707 1201 Int_t i;
3bd79107 1202
58005f18 1203 if (R__b.IsReading()) {
235841d0 1204 Version_t R__v = R__b.ReadVersion();
1205 if (R__v) {
1206 AliDetector::Streamer(R__b);
e8189707 1207 R__b >> fIdN;
1208 R__b.ReadArray(fIdSens);
143d1056 1209 if(fIdName!=0) delete[] fIdName; // Array of TStrings
1210 fIdName = new TString[fIdN];
e8189707 1211 for(i=0;i<fIdN;i++) fIdName[i].Streamer(R__b);
235841d0 1212 R__b >> fITSgeom;
1213 R__b >> fITSmodules;
1214 R__b >> fEuclidOut;
235841d0 1215 R__b >> fMajorVersion;
1216 R__b >> fMinorVersion;
1217 R__b >> fDetTypes;
235841d0 1218 R__b >> fDtype;
1219 delete []fNdtype;
e8189707 1220 fNdtype = new Int_t[fgkNTYPES];
1221 R__b.ReadFastArray(fNdtype,fgkNTYPES);
235841d0 1222 R__b >> fCtype;
1223 delete []fNctype;
e8189707 1224 fNctype = new Int_t[fgkNTYPES];
1225 R__b.ReadFastArray(fNctype,fgkNTYPES);
235841d0 1226 R__b >> fRecPoints;
1227 R__b >> fNRecPoints;
235841d0 1228 R__b >> fTreeC;
1229 } // end if R__v
1230 } else { // writing
58005f18 1231 R__b.WriteVersion(AliITS::IsA());
1232 AliDetector::Streamer(R__b);
e8189707 1233 R__b << fIdN;
1234 R__b.WriteArray(fIdSens,fIdN);
1235 for(i=0;i<fIdN;i++) fIdName[i].Streamer(R__b);
58005f18 1236 R__b << fITSgeom;
2b2439bd 1237 R__b << fITSmodules;
58005f18 1238 R__b << fEuclidOut;
58005f18 1239 R__b << fMajorVersion;
1240 R__b << fMinorVersion;
3bd79107 1241 R__b << fDetTypes;
2b2439bd 1242 R__b << fDtype;
e8189707 1243 R__b.WriteFastArray(fNdtype,fgkNTYPES);
2b2439bd 1244 R__b << fCtype;
e8189707 1245 R__b.WriteFastArray(fNctype,fgkNTYPES);
3bd79107 1246 R__b << fRecPoints;
1247 R__b << fNRecPoints;
2b2439bd 1248 R__b << fTreeC;
235841d0 1249 } // end if
e8189707 1250
fe4da5cc 1251}
55b5a0b3 1252
1253//________________________________________________________________
1254
1255AliITStrack AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArray *fastpoints, Int_t
1256**vettid, Bool_t flagvert ) {
1257
1258//Origin A. Badala' and G.S. Pappalardo: e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
1259
1260
1261 TList *list= new TList();
1262
1263 AliITStrack tr(track);
1264
1265 list->AddLast(&tr);
1266
1267 Double_t Pt=(tr).GetPt();
1268 //cout << "\n Pt = " << Pt <<"\n";
1269
1270
1271 AliITStracking obj(list, reference, this, fastpoints,TMath::Abs(Pt),vettid, flagvert);
1272 list->Delete();
1273 delete list;
1274
1275 Int_t itot=-1;
1276 TVector VecTotLabref(18);
1277 for(Int_t lay=5; lay>=0; lay--) {
1278 TVector VecLabref(3);
1279 VecLabref=(*reference).GetLabTrack(lay);
1280 for(Int_t k=0; k<3; k++) {itot++; VecTotLabref(itot)=VecLabref(k);}
1281 }
1282 Long_t labref;
1283 Int_t freq;
1284 (*reference).Search(VecTotLabref, labref, freq);
1285
1286 if(freq < 5) labref=-labref;
1287 (*reference).SetLabel(labref);
1288
1289 return *reference;
1290
1291}
1292
1293
1294
1295//________________________________________________________________
1296
1297
1298
1299void AliITS::DoTracking(Int_t evNumber, Int_t min_t, Int_t max_t, TFile *file, Bool_t flagvert) {
1300
1301// ex macro for tracking ITS
1302
1303 printf("begin DoTracking - file %p\n",file);
1304
1305 const char *pname="75x40_100x60";
1306
1307 struct GoodTrack {
1308 Int_t lab,code;
1309 Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
1310 Bool_t flag;
1311 };
1312
1313
1314 gAlice->GetEvent(0);
1315
1316 AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC");
1317 AliTPCParam *digp = (AliTPCParam*)file->Get(pname);
1318 if (digp!=0) TPC->SetParam(digp);
1319
1320 GoodTrack gt[7000];
1321 Int_t ngood=0;
1322 ifstream in("good_tracks");
1323
1324 cerr<<"Reading good tracks...\n";
1325 while (in>>gt[ngood].lab>>gt[ngood].code
1326 >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
1327 >>gt[ngood].x >>gt[ngood].y >>gt[ngood].z
1328 >>gt[ngood].pxg >>gt[ngood].pyg >>gt[ngood].pzg
1329 >>gt[ngood].ptg >>gt[ngood].flag) {
1330 ngood++;
1331 cerr<<ngood<<'\r';
1332 if (ngood==7000) {
1333 cerr<<"Too many good tracks !\n";
1334 break;
1335 }
1336 }
1337 if (!in.eof()) cerr<<"Read error (good_tracks) !\n";
1338
1339
1340// Load tracks
1341 TFile *tf=TFile::Open("tpctracks.root");
1342 if (!tf->IsOpen()) {cerr<<"Can't open tpctracks.root !\n"; return ;}
1343 TObjArray tracks(200000);
1344 TTree *tracktree=(TTree*)tf->Get("TreeT");
1345 TBranch *tbranch=tracktree->GetBranch("tracks");
1346 Int_t nentr=(Int_t)tracktree->GetEntries();
1347 for (Int_t k=0; k<nentr; k++) {
1348 AliTPCtrack *iotrack=new AliTPCtrack;
1349 tbranch->SetAddress(&iotrack);
1350 tracktree->GetEvent(k);
1351 tracks.AddLast(iotrack);
1352 }
1353 tf->Close();
1354
1355
1356 Int_t nt = tracks.GetEntriesFast();
1357 cerr<<"Number of found tracks "<<nt<<endl;
1358
1359 TVector DataOut(9);
1360 Int_t kkk=0;
1361
1362 Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
1363
1364 ////////////////////////////// good tracks definition in TPC ////////////////////////////////
1365
1366 ofstream out1 ("AliITSTrag.out");
1367 for (Int_t i=0; i<ngood; i++) out1 << gt[i].ptg << "\n";
1368 out1.close();
1369
1370
1371 TVector vec(5);
1372 TTree *TR=gAlice->TreeR();
1373 Int_t nent=(Int_t)TR->GetEntries();
1374 TClonesArray *recPoints = RecPoints();
1375 Int_t numbpoints;
1376 Int_t totalpoints=0;
1377 Int_t *np = new Int_t[nent];
1378 Int_t **vettid = new Int_t* [nent];
1379 for (Int_t mod=0; mod<nent; mod++) {
1380 vettid[mod]=0;
1381 this->ResetRecPoints();
1382 gAlice->TreeR()->GetEvent(mod+1); //first entry in TreeR is empty
1383 numbpoints = recPoints->GetEntries();
1384 totalpoints+=numbpoints;
1385 np[mod] = numbpoints;
1386 //cout<<" mod = "<<mod<<" numbpoints = "<<numbpoints<<"\n"; getchar();
1387 vettid[mod] = new Int_t[numbpoints];
1388 for (Int_t i=0;i<numbpoints; i++) *(vettid[mod]+i)=0;
1389 }
1390
1391 AliTPCtrack *track;
1392
1393
1394 if(min_t < 0) {min_t = 0; max_t = nt-1;}
1395
1396/*
1397 ///////////////////////////////// Definition of vertex end its error ////////////////////////////
1398 ////////////////////////// In the future it will be given by a method ///////////////////////////
1399 Double_t Vx=0.;
1400 Double_t Vy=0.;
1401 Double_t Vz=0.;
1402
1403 Float_t sigmavx=0.0050; // 50 microns
1404 Float_t sigmavy=0.0050; // 50 microns
1405 Float_t sigmavz=0.010; // 100 microns
1406
1407 //Vx+=gRandom->Gaus(0,sigmavx); Vy+=gRandom->Gaus(0,sigmavy); Vz+=gRandom->Gaus(0,sigmavz);
1408 TVector vertex(3), ervertex(3)
1409 vertex(0)=Vx; vertex(1)=Vy; vertex(2)=Vz;
1410 ervertex(0)=sigmavx; ervertex(1)=sigmavy; ervertex(2)=sigmavz;
1411 /////////////////////////////////////////////////////////////////////////////////////////////////
1412*/
1413
1414 //TDirectory *savedir=gDirectory;
1415
1416 TTree tracktree1("TreeT","Tree with ITS tracks");
1417 AliITSiotrack *iotrack=0;
1418 tracktree1.Branch("ITStracks","AliITSiotrack",&iotrack,32000,0);
1419
1420 ofstream out ("AliITSTra.out");
1421
1422 for (int j=min_t; j<=max_t; j++) {
1423 track=(AliTPCtrack*)tracks.UncheckedAt(j);
1424 Int_t flaglab=0;
1425 if (!track) continue;
1426 ////// elimination of not good tracks ////////////
1427 Int_t ilab=TMath::Abs(track->GetLabel());
1428
1429 for (Int_t i=0;i<ngood;i++) {
1430 //cout<<" ilab, gt[i].lab = "<<ilab<<" "<<gt[i].lab<<"\n"; getchar();
1431 if (ilab==gt[i].lab) {
1432 flaglab=1;
1433 ptg=gt[i].ptg;
1434 pxg=gt[i].pxg;
1435 pyg=gt[i].pyg;
1436 pzg=gt[i].pzg;
1437 break;
1438 }
1439 }
1440 //cout<<" j flaglab = " <<j<<" "<<flaglab<<"\n"; getchar();
1441 if (!flaglab) continue;
1442 //cout<<" j = " <<j<<"\n"; getchar();
1443 /*
1444 ////// old propagation to the end of TPC //////////////
1445 Double_t xk=76.;
1446 track->PropagateTo(xk);
1447 xk-=0.11;
1448 track->PropagateTo(xk,42.7,2.27); //C
1449 xk-=2.6;
1450 track->PropagateTo(xk,36.2,1.98e-3); //C02
1451 xk-=0.051;
1452 track->PropagateTo(xk,42.7,2.27); //C
1453 ///////////////////////////////////////////////////
1454 */
1455
1456 ////// new propagation to the end of TPC //////////////
1457 Double_t xk=77.415;
1458 track->PropagateTo(xk, 28.94, 1.204e-3); //Ne
1459 xk -=0.01;
1460 track->PropagateTo(xk, 44.77, 1.71); //Tedlar
1461 xk -=0.04;
1462 track->PropagateTo(xk, 44.86, 1.45); //Kevlar
1463 xk -=2.0;
1464 track->PropagateTo(xk, 41.28, 0.029); //Nomex
1465 xk-=16;
1466 track->PropagateTo(xk,36.2,1.98e-3); //C02
1467 xk -=0.01;
1468 track->PropagateTo(xk, 24.01, 2.7); //Al
1469 xk -=0.01;
1470 track->PropagateTo(xk, 44.77, 1.71); //Tedlar
1471 xk -=0.04;
1472 track->PropagateTo(xk, 44.86, 1.45); //Kevlar
1473 xk -=0.5;
1474 track->PropagateTo(xk, 41.28, 0.029); //Nomex
1475
1476 ///////////////////////////////////////////////////////////////
1477
1478 ///////////////////////////////////////////////////////////////
1479 AliITStrack trackITS(*track);
1480 AliITStrack result(*track);
1481 AliITStrack primarytrack(*track);
1482
1483///////////////////////////////////////////////////////////////////////////////////////////////
1484 TVector Vgeant(3);
1485 Vgeant=result.GetVertex();
1486
1487 // Definition of Dv and Zv for vertex constraint
1488 // Double_t sigmaDv=0.0050; Double_t sigmaZv=0.010;
1489 Double_t sigmaDv=0.0015; Double_t sigmaZv=0.0015;
1490 Double_t uniform= gRandom->Uniform();
1491 Double_t signdv;
1492 if(uniform<=0.5) signdv=-1.;
1493 else
1494 signdv=1.;
1495
1496 Double_t Vr=TMath::Sqrt(Vgeant(0)*Vgeant(0)+ Vgeant(1)*Vgeant(1));
1497 Double_t Dv=gRandom->Gaus(signdv*Vr,(Float_t)sigmaDv);
1498 Double_t Zv=gRandom->Gaus(Vgeant(2),(Float_t)sigmaZv);
1499
1500 //cout<<" Dv e Zv = "<<Dv<<" "<<Zv<<"\n";
1501 trackITS.SetDv(Dv); trackITS.SetZv(Zv);
1502 trackITS.SetsigmaDv(sigmaDv); trackITS.SetsigmaZv(sigmaZv);
1503 result.SetDv(Dv); result.SetZv(Zv);
1504 result.SetsigmaDv(sigmaDv); result.SetsigmaZv(sigmaZv);
1505 primarytrack.SetDv(Dv); primarytrack.SetZv(Zv);
1506 primarytrack.SetsigmaDv(sigmaDv); primarytrack.SetsigmaZv(sigmaZv);
1507
1508/////////////////////////////////////////////////////////////////////////////////////////////////
1509
1510 primarytrack.PrimaryTrack();
1511 TVector d2=primarytrack.Getd2();
1512 TVector tgl2=primarytrack.Gettgl2();
1513 TVector dtgl=primarytrack.Getdtgl();
1514 trackITS.Setd2(d2); trackITS.Settgl2(tgl2); trackITS.Setdtgl(dtgl);
1515 result.Setd2(d2); result.Settgl2(tgl2); result.Setdtgl(dtgl);
1516 /*
1517 trackITS.SetVertex(vertex); trackITS.SetErrorVertex(ervertex);
1518 result.SetVertex(vertex); result.SetErrorVertex(ervertex);
1519 */
1520 Tracking(trackITS,&result,recPoints,vettid, flagvert);
1521
1522 cout<<" progressive track number = "<<j<<"\r";
1523 // cout<<" progressive track number = "<<j<<"\n";
1524 Long_t labITS=result.GetLabel();
1525 // cout << " ITS track label = " << labITS << "\n";
1526 int lab=track->GetLabel();
1527 //cout << " TPC track label = " << lab <<"\n";
1528 //result.GetClusters(); getchar(); //to print the cluster matrix
1529
1530//propagation to vertex
1531
1532 Double_t rbeam=3.;
1533 result.Propagation(rbeam);
1534 TMatrix *cov;
1535 cov=&result.GetCMatrix();
1536 Double_t pt=TMath::Abs(result.GetPt());
1537 Double_t Dr=result.GetD();
1538 Double_t Z=result.GetZ();
1539 Double_t tgl=result.GetTgl();
1540 Double_t C=result.GetC();
1541 Double_t Cy=C/2.;
1542 Double_t Dz=Z-(tgl/Cy)*TMath::ASin(result.arga(rbeam));
1543 Dz-=Vgeant(2);
1544 // cout<<" Dr e dz alla fine = "<<Dr<<" "<<Dz<<"\n"; getchar();
1545 Double_t phi=result.Getphi();
1546 Double_t phivertex = phi - TMath::ASin(result.argA(rbeam));
1547 Double_t duepi=2.*TMath::Pi();
1548 if(phivertex>duepi) phivertex-=duepi;
1549 if(phivertex<0.) phivertex+=duepi;
1550 Double_t Dtot=TMath::Sqrt(Dr*Dr+Dz*Dz);
1551//////////////////////////////////////////////////////////////////////////////////////////
1552 Int_t NumofCluster, idmodule,idpoint;
1553 NumofCluster=result.GetNumClust();
1554 if(NumofCluster >=5) {
1555
1556
1557 AliITSiotrack outtrack;
1558
1559 iotrack=&outtrack;
1560
1561 iotrack->SetStatePhi(phi);
1562 iotrack->SetStateZ(Z);
1563 iotrack->SetStateD(Dr);
1564 iotrack->SetStateTgl(tgl);
1565 iotrack->SetStateC(C);
1566 Double_t radius=result.Getrtrack();
1567 iotrack->SetRadius(radius);
1568 Int_t charge;
1569 if(C>0.) charge=-1; else charge=1;
1570 iotrack->SetCharge(charge);
1571
1572
1573
1574 iotrack->SetCovMatrix(cov);
1575
1576 Double_t px=pt*TMath::Cos(phi);
1577 Double_t py=pt*TMath::Sin(phi);
1578 Double_t pz=pt*tgl;
1579
1580 Double_t xtrack=Dr*TMath::Sin(phi);
1581 Double_t ytrack=Dr*TMath::Cos(phi);
1582 Double_t ztrack=Dz+Vgeant(2);
1583
1584
1585 iotrack->SetPx(px);
1586 iotrack->SetPy(py);
1587 iotrack->SetPz(pz);
1588 iotrack->SetX(xtrack);
1589 iotrack->SetY(ytrack);
1590 iotrack->SetZ(ztrack);
1591 iotrack->SetLabel(labITS);
1592
1593 for( Int_t il=0;il<6; il++){
1594 iotrack->SetIdPoint(il,result.GetIdPoint(il));
1595 iotrack->SetIdModule(il,result.GetIdModule(il));
1596 }
1597 tracktree1.Fill();
1598
1599 //cout<<" labITS = "<<labITS<<"\n";
1600 //cout<<" phi z Dr tgl C = "<<phi<<" "<<Z<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n"; getchar();
1601
1602 DataOut(kkk) = ptg; kkk++; DataOut(kkk)=labITS; kkk++; DataOut(kkk)=lab; kkk++;
1603
1604 for (Int_t il=0;il<6;il++) {
1605 idpoint=result.GetIdPoint(il);
1606 idmodule=result.GetIdModule(il);
1607 *(vettid[idmodule]+idpoint)=1;
1608 iotrack->SetIdPoint(il,idpoint);
1609 iotrack->SetIdModule(il,idmodule);
1610 }
1611
1612 Double_t difpt= (pt-ptg)/ptg*100.;
1613 DataOut(kkk)=difpt; kkk++;
1614 Double_t lambdag=TMath::ATan(pzg/ptg);
1615 Double_t lam=TMath::ATan(tgl);
1616 Double_t diflam = (lam - lambdag)*1000.;
1617 DataOut(kkk) = diflam; kkk++;
1618 Double_t phig=TMath::ATan(pyg/pxg);
1619 Double_t phi=phivertex;
1620 Double_t difphi = (phi - phig)*1000.;
1621 DataOut(kkk)=difphi; kkk++;
1622 DataOut(kkk)=Dtot*1.e4; kkk++;
1623 DataOut(kkk)=Dr*1.e4; kkk++;
1624 DataOut(kkk)=Dz*1.e4; kkk++;
1625 for (int r=0; r<9; r++) { out<<DataOut(r)<<" ";}
1626 out<<"\n";
1627 kkk=0;
1628
1629
1630 } // end if on NumofCluster
1631 //gObjectTable->Print(); // stampa memoria
1632 } // end for (int j=min_t; j<=max_t; j++)
1633
1634 out.close();
1635
1636
1637 static Bool_t first=kTRUE;
1638 static TFile *tfile;
1639
1640 if(first) {
1641 tfile=new TFile("itstracks.root","RECREATE");
1642 //cout<<"I have opened itstracks.root file "<<endl;
1643 }
1644 first=kFALSE;
1645 tfile->cd();
1646 tfile->ls();
1647
1648 char hname[30];
1649 sprintf(hname,"TreeT%d",evNumber);
1650
1651 tracktree1.Write(hname);
1652
1653
1654
1655 TTree *fAli=gAlice->TreeK();
1656 TFile *fileAli=0;
1657
1658 if (fAli) fileAli =fAli->GetCurrentFile();
1659 fileAli->cd();
1660
1661 ////////////////////////////////////////////////////////////////////////////////////////////////
1662
1663 printf("delete vectors\n");
1664 if(np) delete [] np;
1665 if(vettid) delete [] vettid;
1666
1667}