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 **************************************************************************/
18 Revision 1.12 2000/06/12 23:43:16 nilsen
19 New ITS code replacing the old structure and simulations code.
21 Revision 1.9.2.8 2000/06/12 18:05:59 barbera
22 fixed posible compilation errors on HP unix
24 Revision 1.9.2.7 2000/06/11 20:20:18 barbera
25 New AliITS base clase for the new structure.
27 Revision 1.9.2.3 2000/02/02 13:42:09 barbera
28 fixed AliITS.cxx for new AliRun structure. Added ITS hits list to list of hits which will have their track numbers updated
30 Revision 1.9.2.2 2000/01/23 03:03:13 nilsen
31 //fixed FillModule. Removed fi(fabs(xl)<dx....
33 Revision 1.9.2.1 2000/01/12 19:03:32 nilsen
34 This is the version of the files after the merging done in December 1999.
35 See the ReadMe110100.txt file for details
37 Revision 1.9 1999/11/14 14:33:25 fca
38 Correct problems with distructors and pointers, thanks to I.Hrivnacova
40 Revision 1.8 1999/09/29 09:24:19 fca
41 Introduction of the Copyright and cvs Log
45 ///////////////////////////////////////////////////////////////////////////////
47 // An overview of the basic philosophy of the ITS code development
48 // and analysis is show in the figure below.
51 <img src="picts/ITS/ITS_Analysis_schema.gif">
54 <font size=+2 color=red>
55 <p>Roberto Barbera is in charge of the ITS Offline code (1999).
56 <a href="mailto:roberto.barbera@ct.infn.it">Roberto Barbera</a>.
62 // AliITS. Inner Traking System base class.
63 // This class contains the base procedures for the Inner Tracking System
67 <img src="picts/ITS/AliITS_Class_Diagram.gif">
70 <font size=+2 color=red>
71 <p>This show the class diagram of the different elements that are part of
79 // Written by Rene Brun, Federico Carminati, and Roberto Barbera
82 // Modified and documented by Bjorn S. Nilsen
86 // Modified and documented by A. Bologna
89 // AliITS is the general base class for the ITS. Also see AliDetector for
90 // futher information.
92 ///////////////////////////////////////////////////////////////////////////////
97 #include <TObjArray.h>
98 #include <TClonesArray.h>
100 #include <TObjectTable.h>
106 #include "AliITSMap.h"
107 #include "AliITSClusterFinder.h"
108 #include "AliITSsimulation.h"
109 #include "AliITSsimulationFastPoints.h"
110 #include "AliITSsegmentationSPD.h"
111 #include "AliITSresponseSPD.h"
112 #include "AliITSsegmentationSDD.h"
113 #include "AliITSresponseSDD.h"
114 #include "AliITSsegmentationSSD.h"
115 #include "AliITSresponseSSD.h"
116 //#include "AliITStrack.h"
121 //_____________________________________________________________________________
122 AliITS::AliITS() : AliDetector() {
124 // Default initialiser for ITS
125 // The default constructor of the AliITS class. In addition to
126 // creating the AliITS class it zeros the variables fIshunt (a member
127 // of AliDetector class), fEuclidOut, and fIdN, and zeros the pointers
128 // fITSpoints, fIdSens, and fIdName. The AliDetector default constructor
152 //_____________________________________________________________________________
153 AliITS::AliITS(const char *name, const char *title):AliDetector(name,title){
155 // Default initialiser for ITS
156 // The constructor of the AliITS class. In addition to creating the
157 // AliITS class, it allocates memory for the TClonesArrays fHits and
158 // fDigits, and for the TObjArray fITSpoints. It also zeros the variables
159 // fIshunt (a member of AliDetector class), fEuclidOut, and fIdN, and zeros
160 // the pointers fIdSens and fIdName. To help in displaying hits via the ROOT
161 // macro display.C AliITS also sets the marker color to red. The variables
162 // passes with this constructor, const char *name and *title, are used by
163 // the constructor of AliDetector class. See AliDetector class for a
164 // description of these parameters and its constructor functions.
167 fHits = new TClonesArray("AliITShit", 1560);
168 gAlice->AddHitList(fHits);
189 fDetTypes = new TObjArray(fNDetTypes);
192 for(i=0;i<fNDetTypes;i++) {
193 (*fDetTypes)[i]=new AliITSDetType();
197 SetMarkerColor(kRed);
201 //___________________________________________________________________________
202 AliITS::AliITS(AliITS &source){
203 if(this==&source) return;
204 printf("Error: You are not allowed to make a copy of the AliITS\n");
207 //____________________________________________________________________________
208 AliITS& AliITS::operator=(AliITS &source){
209 if(this==&source) return *this;
210 printf("Error: You are not allowed to make a copy of the AliITS\n");
213 //____________________________________________________________________________
214 void AliITS::ClearModules(){
215 //clear the modules TObjArray
219 Int_t indSPD = fITSgeom->GetModuleIndex(2,fITSgeom->GetNladders(2),
220 fITSgeom->GetNdetectors(2));
221 Int_t indSDD = fITSgeom->GetModuleIndex(4,fITSgeom->GetNladders(4),
222 fITSgeom->GetNdetectors(4));
223 for(i=0;i<fITSmodules->GetEntriesFast();i++){
225 delete (AliITSmodule *) fITSmodules->At(i);
227 delete (AliITSmodule *) fITSmodules->At(i);
229 delete (AliITSmodule *) fITSmodules->At(i);
231 }// end if fITSmodules!=0
234 //_____________________________________________________________________________
237 // Default distructor for ITS
238 // The default destructor of the AliITS class. In addition to deleting
239 // the AliITS class it deletes the memory pointed to by the fHits, fDigits,
240 // fIdSens, fIdName, and fITSpoints.
247 if(fIdName!=0) delete[] fIdName;
248 if(fIdSens!=0) delete[] fIdSens;
250 this->ClearModules();
252 }// end if fITSmodules!=0
257 for(i=0;i<fNDetTypes;i++) {
263 for(i=0;i<fNDetTypes;i++) {
275 if (fTreeC) delete fTreeC;
279 //___________________________________________
280 AliITSDetType* AliITS::DetType(Int_t id)
282 //return pointer to id detector type
283 return ((AliITSDetType*) (*fDetTypes)[id]);
286 //___________________________________________
287 void AliITS::SetClasses(Int_t id, TString digit, TString cluster)
289 //set the digit and cluster classes to be used for the id detector type
290 ((AliITSDetType*) (*fDetTypes)[id])->ClassNames(digit,cluster);
293 //___________________________________________
294 void AliITS::SetResponseModel(Int_t id, AliITSresponse *response)
296 //set the response model for the id detector type
298 ((AliITSDetType*) (*fDetTypes)[id])->ResponseModel(response);
302 //___________________________________________
303 void AliITS::SetSegmentationModel(Int_t id, AliITSsegmentation *seg)
305 //set the segmentation model for the id detector type
307 ((AliITSDetType*) (*fDetTypes)[id])->SegmentationModel(seg);
311 //___________________________________________
312 void AliITS::SetSimulationModel(Int_t id, AliITSsimulation *sim)
314 //set the simulation model for the id detector type
316 ((AliITSDetType*) (*fDetTypes)[id])->SimulationModel(sim);
319 //___________________________________________
320 void AliITS::SetReconstructionModel(Int_t id, AliITSClusterFinder *reconst)
322 //set the cluster finder model for the id detector type
324 ((AliITSDetType*) (*fDetTypes)[id])->ReconstructionModel(reconst);
328 //_____________________________________________________________________________
329 void AliITS::AddHit(Int_t track, Int_t *vol, Float_t *hits){
332 // The function to add information to the AliITShit class. See the
333 // AliITShit class for a full description. This function allocates the
334 // necessary new space for the hit information and passes the variable
335 // track, and the pointers *vol and *hits to the AliITShit constructor
338 TClonesArray &lhits = *fHits;
339 new(lhits[fNhits++]) AliITShit(fIshunt,track,vol,hits);
341 //_____________________________________________________________________________
342 void AliITS::AddRealDigit(Int_t id, Int_t *digits)
344 // add a real digit - as coming from data
346 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
347 new(ldigits[fNdtype[id]++]) AliITSdigit(digits);
350 //_____________________________________________________________________________
351 void AliITS::AddDigit(Int_t id, AliITSdigit *d)
354 // add a simulated digit
356 // should have ctors of type AliITSdigitSDD(const AliITSdigitSDD &)
358 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
363 new(ldigits[fNdtype[id]++]) AliITSdigitSPD(*((AliITSdigitSPD*)d));
366 new(ldigits[fNdtype[id]++]) AliITSdigitSDD(*((AliITSdigitSDD*)d));
369 new(ldigits[fNdtype[id]++]) AliITSdigitSSD(*((AliITSdigitSSD*)d));
375 //_____________________________________________________________________________
376 void AliITS::AddDigit(Int_t id,Float_t phys,Int_t *digits,Int_t *tracks,Float_t *charges){
378 // add a simulated digit to the list
380 TClonesArray &ldigits = *((TClonesArray*)(*fDtype)[id]);
384 new(ldigits[fNdtype[id]++]) AliITSdigitSPD(digits,tracks);
387 new(ldigits[fNdtype[id]++]) AliITSdigitSDD(phys,digits,tracks,charges);
390 new(ldigits[fNdtype[id]++]) AliITSdigitSSD(digits,tracks);
396 //_____________________________________________________________________________
397 void AliITS::AddCluster(Int_t id, AliITSRawCluster *c)
400 // add a cluster to the list
402 // should have ctors of type AliITSRawClusterSDD(const AliITSRawClusterSDD &)
404 TClonesArray &lcl = *((TClonesArray*)(*fCtype)[id]);
409 new(lcl[fNctype[id]++]) AliITSRawClusterSPD(*((AliITSRawClusterSPD*)c));
412 new(lcl[fNctype[id]++]) AliITSRawClusterSDD(*((AliITSRawClusterSDD*)c));
415 new(lcl[fNctype[id]++]) AliITSRawClusterSSD(*((AliITSRawClusterSSD*)c));
422 //_____________________________________________________________________________
423 void AliITS::AddRecPoint(const AliITSRecPoint &r)
426 // Add a reconstructed space point to the list
428 TClonesArray &lrecp = *fRecPoints;
429 new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
433 //____________________________________________
434 void AliITS::ResetDigits()
437 // Reset number of digits and the digits array for thE ITS detector
443 for(i=0;i<fNDetTypes;i++ ) {
444 if ((*fDtype)[i]) ((TClonesArray*)(*fDtype)[i])->Clear();
445 if (fNdtype) fNdtype[i]=0;
449 //____________________________________________
450 void AliITS::ResetDigits(Int_t i)
453 // Reset number of digits and the digits array for this branch
455 if ((*fDtype)[i]) ((TClonesArray*)(*fDtype)[i])->Clear();
456 if (fNdtype) fNdtype[i]=0;
459 //____________________________________________
460 void AliITS::ResetClusters()
463 // Reset number of clusters and the clusters array for ITS
467 for(i=0;i<fNDetTypes;i++ ) {
468 if ((*fCtype)[i]) ((TClonesArray*)(*fCtype)[i])->Clear();
469 if (fNctype) fNctype[i]=0;
473 //____________________________________________
474 void AliITS::ResetClusters(Int_t i)
477 // Reset number of clusters and the clusters array for this branch
479 if ((*fCtype)[i]) ((TClonesArray*)(*fCtype)[i])->Clear();
480 if (fNctype) fNctype[i]=0;
485 //____________________________________________
486 void AliITS::ResetRecPoints()
489 // Reset number of rec points and the rec points array
492 if (fRecPoints) fRecPoints->Clear();
496 //_____________________________________________________________________________
497 Int_t AliITS::DistancetoPrimitive(Int_t , Int_t ){
499 // Distance from mouse to ITS on the screen. Dummy routine
500 // A dummy routine used by the ROOT macro display.C to allow for the
501 // use of the mouse (pointing device) in the macro. In general this should
502 // never be called. If it is it returns the number 9999 for any value of
508 //_____________________________________________________________________________
511 // Initialise ITS after it has been built
512 // This routine initializes the AliITS class. It is intended to be called
513 // from the Init function in AliITSv?. Besides displaying a banner
514 // indicating that it has been called it initializes the array fIdSens
515 // and sets the default segmentation, response, digit and raw cluster classes
516 // Therefore it should be called after a call to CreateGeometry.
525 for(i=0;i<35;i++) printf("*");
526 printf(" ITS_INIT ");
527 for(i=0;i<35;i++) printf("*");
531 for(i=0;i<fIdN;i++) fIdSens[i] = gMC->VolId(fIdName[i]);
533 for(i=0;i<80;i++) printf("*");
537 //_____________________________________________________________________________
538 void AliITS::SetDefaults()
540 // sets the default segmentation, response, digit and raw cluster classes
543 AliITSDetType *iDetType;
547 AliITSsegmentationSPD *seg0=new AliITSsegmentationSPD(fITSgeom);
548 AliITSresponseSPD *resp0=new AliITSresponseSPD();
550 if (!iDetType->GetSegmentationModel()) SetSegmentationModel(0,seg0);
551 if (!iDetType->GetResponseModel()) SetResponseModel(0,resp0);
552 // set digit and raw cluster classes to be used
553 const char *kData=resp0->DataType();
554 if (strstr(kData,"real")) {
555 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSPD");
556 } else iDetType->ClassNames("AliITSdigitSPD","AliITSRawClusterSPD");
559 AliITSresponseSDD *resp1=new AliITSresponseSDD();
560 AliITSsegmentationSDD *seg1=new AliITSsegmentationSDD(fITSgeom,resp1);
562 if (!iDetType->GetSegmentationModel()) SetSegmentationModel(1,seg1);
563 if (!iDetType->GetResponseModel()) SetResponseModel(1,resp1);
564 kData=resp1->DataType();
565 Option_t *opt=resp1->ZeroSuppOption();
566 if ((!strstr(opt,"2D")) && (!strstr(opt,"1D")) || strstr(kData,"real") ) {
567 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSDD");
568 } else iDetType->ClassNames("AliITSdigitSDD","AliITSRawClusterSDD");
571 AliITSsegmentationSSD *seg2=new AliITSsegmentationSSD(fITSgeom);
572 AliITSresponseSSD *resp2=new AliITSresponseSSD();
574 if (!iDetType->GetSegmentationModel()) SetSegmentationModel(2,seg2);
575 if (!iDetType->GetResponseModel()) SetResponseModel(2,resp2);
576 kData=resp2->DataType();
577 if (strstr(kData,"real")) {
578 iDetType->ClassNames("AliITSdigit","AliITSRawClusterSSD");
579 } else iDetType->ClassNames("AliITSdigitSSD","AliITSRawClusterSSD");
582 Warning("SetDefaults","Only the three basic detector types are initialised!");
588 //_____________________________________________________________________________
590 void AliITS::MakeTreeC(Option_t *option)
592 // create a separate tree to store the clusters
594 char *optC = strstr(option,"C");
595 if (optC && !fTreeC) fTreeC = new TTree("TC","Clusters in ITS");
597 Int_t buffersize = 4000;
600 char *det[3] = {"SPD","SDD","SSD"};
602 // one branch for Clusters per type of detector
604 for(i=0; i<fNDetTypes ;i++) {
605 if (fNDetTypes==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
606 else sprintf(branchname,"%sClusters%d",GetName(),i+1);
607 if (fCtype && fTreeC) {
608 TreeC()->Branch(branchname,&((*fCtype)[i]), buffersize);
609 printf("Making Branch %s for Clusters of detector type %d\n",branchname,i+1);
615 //_____________________________________________________________________________
616 void AliITS::GetTreeC(Int_t event)
619 // get the clusters tree for this event and set the branch address
623 char *det[3] = {"SPD","SDD","SSD"};
630 sprintf(treeName,"TreeC%d",event);
631 fTreeC = (TTree*)gDirectory->Get(treeName);
637 for(i=0; i<fNDetTypes; i++) {
638 if (fNDetTypes==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
639 else sprintf(branchname,"%sClusters%d",GetName(),i+1);
641 branch = fTreeC->GetBranch(branchname);
642 if (branch) branch->SetAddress(&((*fCtype)[i]));
646 printf("ERROR: cannot find Clusters Tree for event:%d\n",event);
650 //_____________________________________________________________________________
651 void AliITS::MakeBranch(Option_t* option){
653 // Creates Tree branches for the ITS.
657 Int_t buffersize = 4000;
659 sprintf(branchname,"%s",GetName());
661 AliDetector::MakeBranch(option);
664 // one branch for digits per type of detector
666 char *det[3] = {"SPD","SDD","SSD"};
668 fNdtype = new Int_t[fNDetTypes];
669 fDtype = new TObjArray(fNDetTypes);
671 fNctype = new Int_t[fNDetTypes];
672 fCtype = new TObjArray(fNDetTypes);
674 const char *kDigclass;
675 const char *kClclass;
678 for(i=0; i<fNDetTypes ;i++) {
679 AliITSDetType *iDetType=DetType(i);
680 iDetType->GetClassNames(kDigclass,kClclass);
681 //printf("i, digclass, recclass %d %s %s\n",i,kDigclass,kClclass);
683 (*fDtype)[i] = new TClonesArray(kDigclass,100);
686 (*fCtype)[i] = new TClonesArray(kClclass,100);
691 for(i=0; i<fNDetTypes ;i++) {
692 if (fNDetTypes==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
693 else sprintf(branchname,"%sDigits%d",GetName(),i+1);
695 if (fDtype && gAlice->TreeD()) {
696 gAlice->TreeD()->Branch(branchname,&((*fDtype)[i]), buffersize);
697 printf("Making Branch %s for digits of type %d\n",branchname,i+1);
701 // only one branch for rec points for all detector types
702 sprintf(branchname,"%sRecPoints",GetName());
704 fRecPoints=new TClonesArray("AliITSRecPoint",1000);
706 if (fRecPoints && gAlice->TreeR()) {
707 gAlice->TreeR()->Branch(branchname,&fRecPoints, buffersize);
708 printf("Making Branch %s for reconstructed space points\n",branchname);
714 //___________________________________________
715 void AliITS::SetTreeAddress()
718 // Set branch address for the Trees.
721 AliDetector::SetTreeAddress();
723 char *det[3] = {"SPD","SDD","SSD"};
726 TTree *treeD = gAlice->TreeD();
727 TTree *treeR = gAlice->TreeR();
731 for(i=0; i<fNDetTypes; i++) {
732 if (fNDetTypes==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
733 else sprintf(branchname,"%sDigits%d",GetName(),i+1);
735 branch = treeD->GetBranch(branchname);
736 if (branch) branch->SetAddress(&((*fDtype)[i]));
743 sprintf(branchname,"%sRecPoints",GetName());
745 branch = treeR->GetBranch(branchname);
746 if (branch) branch->SetAddress(&fRecPoints);
753 //____________________________________________________________________________
754 void AliITS::InitModules(Int_t size,Int_t &nmodules){
756 //initialize the modules array
758 Int_t nl,indexMAX,index;
761 if(size<=0){ // default to using data stored in AliITSgeom
763 printf("Error in AliITS::InitModule fITSgeom not defined\n");
765 } // end if fITSgeom==0
766 nl = fITSgeom->GetNlayers();
767 indexMAX = fITSgeom->GetModuleIndex(nl,fITSgeom->GetNladders(nl),
768 fITSgeom->GetNdetectors(nl))+1;
770 fITSmodules = new TObjArray(indexMAX);
771 indSPD = fITSgeom->GetModuleIndex(2,fITSgeom->GetNladders(2),
772 fITSgeom->GetNdetectors(2));
773 indSDD = fITSgeom->GetModuleIndex(4,fITSgeom->GetNladders(4),
774 fITSgeom->GetNdetectors(4));
775 for(index=0;index<indexMAX;index++){
777 fITSmodules->AddAt( new AliITSmodule(index),index);
778 else if(index<=indSDD)
779 fITSmodules->AddAt( new AliITSmodule(index),index);
781 fITSmodules->AddAt( new AliITSmodule(index),index);
784 fITSmodules = new TObjArray(size);
789 //____________________________________________________________________________
790 void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t lastev,Int_t nmodules,Option_t *option,Text_t *filename){
792 // fill the modules with the sorted by module hits; add hits from background
795 static TTree *trH1; //Tree with background hits
796 static TClonesArray *fHits2; //List of hits for one track only
798 static Bool_t first=kTRUE;
800 char *add = strstr(option,"Add");
805 cout<<"filename "<<filename<<endl;
806 file=new TFile(filename);
807 cout<<"I have opened "<<filename<<" file "<<endl;
808 fHits2 = new TClonesArray("AliITShit",1000 );
813 // Get Hits Tree header from file
814 if(fHits2) fHits2->Clear();
815 if(trH1) delete trH1;
819 sprintf(treeName,"TreeH%d",bgrev);
820 trH1 = (TTree*)gDirectory->Get(treeName);
821 //printf("TrH1 %p of treename %s for event %d \n",trH1,treeName,bgrev);
824 printf("ERROR: cannot find Hits Tree for event:%d\n",bgrev);
826 // Set branch addresses
829 sprintf(branchname,"%s",GetName());
830 if (trH1 && fHits2) {
831 branch = trH1->GetBranch(branchname);
832 if (branch) branch->SetAddress(&fHits2);
836 //Int_t ntracks1 =(Int_t)TrH1->GetEntries();
837 //printf("background - ntracks1 - %d\n",ntracks1);
840 // store temporarily coordinates of signal particles - see where delete
841 static Int_t *signal;
843 if(!signal) signal=new Int_t[nmodules];
844 memset(signal,0,sizeof(int)*nmodules);
845 // Float_t xhit[nmodules][4];
846 Float_t **xhit = new Float_t*[nmodules];
847 for(i=0;i<nmodules;i++) xhit[i] = new Float_t[4];
848 // Float_t yhit[nmodules][4];
849 Float_t **yhit = new Float_t*[nmodules];
850 for(i=0;i<nmodules;i++) yhit[i] = new Float_t[4];
852 Int_t npart = gAlice->GetEvent(evnt);
854 TClonesArray *itsHits = this->Hits();
855 Int_t lay,lad,det,index;
859 TTree *iTH = gAlice->TreeH();
860 Int_t ntracks =(Int_t) iTH->GetEntries();
863 for(t=0; t<ntracks; t++){
866 Int_t nhits = itsHits->GetEntriesFast();
867 if (!nhits) continue;
868 // cout << nhits << " hits in track " << t << endl;
869 for(h=0; h<nhits; h++){
870 itsHit = (AliITShit *)itsHits->UncheckedAt(h);
871 itsHit->GetDetectorID(lay,lad,det);
872 index = fITSgeom->GetModuleIndex(lay,lad,det);
873 mod = this->GetModule(index);
875 xhit[index][signal[index]]=itsHit->fX;
876 yhit[index][signal[index]]=itsHit->fY;
878 if (signal[index] >4)
879 printf("index,nsignal %d %d\n",index,signal[index]);
881 if(lay == 1 || lay == 2)
882 mod->AddHit((AliITShit *) itsHit,t,h);
883 else if(lay == 3 || lay == 4)
884 mod->AddHit((AliITShit *) itsHit,t,h);
885 else if(lay == 5 || lay ==6)
886 mod->AddHit((AliITShit *)itsHit,t,h);
888 mod->AddHit(itsHit,t,h);
890 } // end loop over hits
891 } // end loop over tracks
893 // open the file with background
897 ntracks =(Int_t)trH1->GetEntries();
898 //printf("background - ntracks1 %d\n",ntracks);
899 //printf("background - Start loop over tracks \n");
902 for(track=0; track<ntracks; track++) {
904 if (fHits2) fHits2->Clear();
905 trH1->GetEvent(track);
907 for(i=0;i<fHits2->GetEntriesFast();++i) {
909 itsHit=(AliITShit*) (*fHits2)[i];
910 itsHit->GetDetectorID(lay,lad,det);
911 index = fITSgeom->GetModuleIndex(lay,lad,det);
912 mod = this->GetModule(index);
914 Float_t xbgr=itsHit->fX;
915 Float_t ybgr=itsHit->fY;
916 Float_t ebgr=itsHit->GetIonization();
919 for(isig =0; isig < signal[index]; isig++) {
921 (xbgr-xhit[index][isig])*(xbgr-xhit[index][isig])
922 +(ybgr-yhit[index][isig])*(ybgr-yhit[index][isig]);
923 if (dist<0.2&& ebgr!=0) cond=kTRUE; // check this number for ITS!
926 if(lay == 1 || lay == 2)
927 mod->AddHit((AliITShit *) itsHit,track,i);
928 else if(lay == 3 || lay == 4)
929 mod->AddHit((AliITShit *) itsHit,track,i);
930 else if(lay == 5 || lay ==6)
931 mod->AddHit((AliITShit *)itsHit,track,i);
933 mod->AddHit(itsHit,track,i);
936 } // end loop over hits
937 } // end loop over tracks
939 TTree *fAli=gAlice->TreeK();
942 if (fAli) fileAli =fAli->GetCurrentFile();
943 //printf("fAli, file %p %p\n",fAli,file);
949 for(i=0;i<nmodules;i++) delete[] xhit[i];
951 for(i=0;i<nmodules;i++) delete[] yhit[i];
953 //if (evnt==lastev) {delete [] signal; delete signal;}
955 //gObjectTable->Print();
960 //____________________________________________________________________________
961 void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t lastev,Int_t size,
962 Option_t *option,Option_t *opt,Text_t *filename)
964 // keep galice.root for signal and name differently the file for
965 // background when add! otherwise the track info for signal will be lost !
967 char *all = strstr(opt,"All");
968 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
969 //printf("Det 1 2 3 %s %s %s \n",det[0],det[1],det[2]);
972 InitModules(size,nmodules);
973 FillModules(evNumber,bgrev,lastev,nmodules,option,filename);
974 //printf("nmodules %d\n",nmodules);
977 AliITSsimulation* sim;
979 TObjArray *branches=gAlice->TreeD()->GetListOfBranches();
980 AliITSgeom *geom = GetITSgeom();
983 for(id=0;id<3;id++) {
984 //printf("id %d All %s det[id] %s \n",id,all,det[id]);
985 if (!all && !det[id]) continue;
986 branch = (TBranch*)branches->UncheckedAt(id);
987 AliITSDetType *iDetType=DetType(id);
988 sim = (AliITSsimulation*)iDetType->GetSimulationModel();
990 Error("HitsToDigits","The simulation class was not instantiated!");
992 // or SetDefaultSimulation(id,iDetType*);
994 Int_t first = geom->GetStartDet(id);
995 Int_t last = geom->GetLastDet(id);
996 //printf("det type %d first, last %d %d \n",id,first,last);
997 for(module=first;module<=last;module++) {
998 AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
999 sim->DigitiseModule(mod,module,evNumber);
1000 // fills all branches - wasted disk space
1001 gAlice->TreeD()->Fill();
1003 // try and fill only the branch
1006 } // loop over modules
1007 } // loop over detector types
1012 //Int_t nentries=(Int_t)gAlice->TreeD()->GetEntries();
1015 sprintf(hname,"TreeD%d",evNumber);
1016 gAlice->TreeD()->Write(hname);
1018 gAlice->TreeD()->Reset();
1023 //____________________________________________________________________________
1024 void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
1027 char *all = strstr(opt,"All");
1028 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
1030 static Bool_t first=kTRUE;
1038 AliITSClusterFinder* rec;
1040 TObjArray *branches=gAlice->TreeR()->GetListOfBranches();
1041 AliITSgeom *geom = GetITSgeom();
1044 for(id=0;id<fNDetTypes;id++) {
1045 if (!all && !det[id]) continue;
1046 branch = (TBranch*)branches->UncheckedAt(id);
1047 AliITSDetType *iDetType=DetType(id);
1048 rec = (AliITSClusterFinder*)iDetType->GetReconstructionModel();
1050 Error("DigitsToRecPoints","The cluster finder class was not instantiated!");
1052 // or SetDefaultClusterFinder(id,iDetType*);
1054 TClonesArray *itsDigits = this->DigitsAddress(id);
1056 Int_t first = geom->GetStartDet(id);
1057 Int_t last = geom->GetLastDet(id);
1058 //printf("det type %d first, last %d %d \n",id,first,last);
1059 for(module=first;module<=last;module++) {
1060 //printf("AliITS: module=%d\n",module);
1061 this->ResetDigits();
1062 if (all) gAlice->TreeD()->GetEvent(lastentry+module);
1063 else gAlice->TreeD()->GetEvent(lastentry+(module-first));
1064 Int_t ndigits = itsDigits->GetEntriesFast();
1065 if (ndigits) rec->FindRawClusters();
1066 gAlice->TreeR()->Fill();
1070 // try and fill only the branch
1072 //ResetRecPoints(id);
1073 } // loop over modules
1074 } // loop over detector types
1077 //Int_t nentries=(Int_t)gAlice->TreeR()->GetEntries();
1079 //Int_t ncentries=(Int_t)TC->GetEntries();
1082 sprintf(hname,"TreeR%d",evNumber);
1083 gAlice->TreeR()->Write(hname);
1085 gAlice->TreeR()->Reset();
1087 sprintf(hname,"TreeC%d",evNumber);
1093 //____________________________________________________________________________
1094 void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t lastev,Int_t size,
1095 Option_t *option,Option_t *opt,Text_t *filename)
1097 // keep galice.root for signal and name differently the file for
1098 // background when add! otherwise the track info for signal will be lost !
1100 char *all = strstr(opt,"All");
1101 char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
1104 InitModules(size,nmodules);
1105 FillModules(evNumber,bgrev,lastev,nmodules,option,filename);
1107 static AliITSsimulationFastPoints* sim=0;
1108 if (!sim) sim = new AliITSsimulationFastPoints();
1111 AliITSgeom *geom = GetITSgeom();
1114 for(id=0;id<3;id++) {
1115 if (!all && !det[id]) continue;
1116 Int_t first = geom->GetStartDet(id);
1117 Int_t last = geom->GetLastDet(id);
1118 for(module=first;module<=last;module++) {
1119 AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
1120 sim->CreateFastRecPoints(mod);
1121 gAlice->TreeR()->Fill();
1123 } // loop over modules
1124 } // loop over detector types
1129 //Int_t nentries=(Int_t)gAlice->TreeR()->GetEntries();
1132 sprintf(hname,"TreeR%d",evNumber);
1133 gAlice->TreeR()->Write(hname);
1135 gAlice->TreeR()->Reset();
1139 //____________________________________________________________________________
1140 void AliITS::Streamer(TBuffer &R__b){
1141 // Stream an object of class AliITS.
1145 AliITSresponse *response;
1146 AliITSsegmentation *seg;
1147 TClonesArray *digitsaddress;
1148 TClonesArray *clustaddr;
1151 if (R__b.IsReading()) {
1152 Version_t R__v = R__b.ReadVersion();
1154 AliDetector::Streamer(R__b);
1158 if(fIdSens!=0) delete[] fIdSens;
1159 if(fIdName!=0) delete[] fIdName;
1160 fIdSens = new Int_t[fIdN];
1161 fIdName = new char*[fIdN];
1162 for(i=0;i<fIdN;i++) R__b >> fIdSens[i];
1163 for(i=0;i<fIdN;i++){
1165 fIdName[i] = new char[l+1]; // add room for null character.
1166 for(j=0;j<l;j++) R__b >> fIdName[i][j];
1167 fIdName[i][l] = '\0'; // Null terminate this string.
1169 R__b >> fMajorVersion;
1170 R__b >> fMinorVersion;
1172 R__b.ReadArray(fNdtype);
1174 R__b.ReadArray(fNctype);
1178 R__b >> fNRecPoints;
1179 // stream out only response and segmentation
1180 for(i=0;i<fNDetTypes;i++) {
1181 det=(AliITSDetType*)(*fDetTypes)[i];
1182 det->Streamer(R__b);
1183 response=det->GetResponseModel();
1184 if (response) response->Streamer(R__b);
1185 seg=det->GetSegmentationModel();
1186 if (seg) seg->Streamer(R__b);
1187 digitsaddress=(TClonesArray*) (*fDtype)[i];
1188 digitsaddress->Streamer(R__b);
1189 clustaddr=(TClonesArray*) (*fCtype)[i];
1190 clustaddr->Streamer(R__b);
1195 R__b.WriteVersion(AliITS::IsA());
1196 AliDetector::Streamer(R__b);
1200 for(i=0;i<fIdN;i++) R__b <<fIdSens[i];
1201 for(i=0;i<fIdN;i++){
1202 l = strlen(fIdName[i]);
1204 for(j=0;j<l;j++) R__b << fIdName[i][j];
1206 R__b << fMajorVersion;
1207 R__b << fMinorVersion;
1209 R__b.WriteArray(fNdtype, fNDetTypes);
1211 R__b.WriteArray(fNctype, fNDetTypes);
1215 R__b << fNRecPoints;
1216 for(i=0;i<fNDetTypes;i++) {
1217 det=(AliITSDetType*)(*fDetTypes)[i];
1218 det->Streamer(R__b);
1219 response=det->GetResponseModel();
1220 if(response) response->Streamer(R__b);
1221 seg=det->GetSegmentationModel();
1222 if (seg) seg->Streamer(R__b);
1223 digitsaddress=(TClonesArray*) (*fDtype)[i];
1224 digitsaddress->Streamer(R__b);
1225 clustaddr=(TClonesArray*) (*fCtype)[i];
1226 clustaddr->Streamer(R__b);
1233 ClassImp(AliITSRecPoint)
1234 ClassImp(AliITSsegmentation)
1235 ClassImp(AliITSresponse)
1236 //ClassImp(AliITStrack)