]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ITS/AliITS.cxx
TPC and CONTAINERS included in the search path
[u/mrichter/AliRoot.git] / ITS / AliITS.cxx
index 9169e5b128f3a46258c11ca543b95d073fc65d41..978ef192a404bcda009e66ac249b8f1f295d188a 100644 (file)
 
 /*
 $Log$
+Revision 1.52  2001/05/30 06:04:58  hristov
+Changes made to be consitant with changes in TPC tracking classes (B.Nilsen)
+
+Revision 1.51  2001/05/16 14:57:15  alibrary
+New files for folders and Stack
+
+Revision 1.50  2001/05/11 09:15:21  barbera
+Corrected to make fast point creation working with PPR geometry
+
+Revision 1.49  2001/05/11 07:37:49  hristov
+Legacy lines commented
+
+Revision 1.48  2001/05/10 18:14:25  barbera
+A typo corrected
+
+Revision 1.47  2001/05/10 17:55:59  barbera
+Modified to create rec points also for PPR geometries
+
+Revision 1.46  2001/05/10 00:05:28  nilsen
+Allowed for HitsToDigits function to work with versions 5, 7, 8, and 9. This
+should probably be cleaned up to only check to make sure that fITSgeom has
+been properly defined.
+
+Revision 1.45  2001/05/01 22:35:48  nilsen
+Remove/commented a number of cout<< statements. and made change needed by
+SSD code.
+
+Revision 1.44  2001/04/26 22:44:01  nilsen
+Removed dependence on layer 5/6 in AliITS::HitsToDigits. This will be
+done properly in AliITSv???.cxx via SetDefaults.
+
+Revision 1.43  2001/04/26 13:22:52  barbera
+TMatrix and TVector elimininated to speed up the code
+
+Revision 1.42  2001/04/25 21:55:12  barbera
+Updated version to be compatible with actual verion of STEER and TPC
+
+Revision 1.41  2001/04/21 15:16:51  barbera
+Updated with the new SSD reconstruction code
+
+Revision 1.40  2001/03/17 15:07:06  mariana
+Update SDD response parameters
+
+Revision 1.39  2001/03/12 17:45:32  hristov
+Changes needed on Sun with CC 5.0
+
+Revision 1.38  2001/03/07 14:04:51  barbera
+Some vector dimensions increased to cope with full events
+
+Revision 1.37  2001/03/07 12:36:35  barbera
+A change added in the tracking part to manage delta rays
+
+Revision 1.36  2001/03/02 19:44:11  barbera
+ modified to taking into account new version tracking v1
+
+Revision 1.35  2001/02/28 18:16:46  mariana
+Make the code compatible with the new AliRun
+
+Revision 1.34  2001/02/11 15:51:39  mariana
+Set protection in MakeBranch
+
+Revision 1.33  2001/02/10 22:26:39  mariana
+Move the initialization of the containers for raw clusters in MakeTreeC()
+
+Revision 1.32  2001/02/08 23:55:31  nilsen
+Removed fMajor/MinorVersion variables in favor of variables in derived classes.
+Set arrays char *det[3] = {"SPD","SDD","SSD"} as const.
+
+Revision 1.31  2001/02/02 23:57:28  nilsen
+Added include file that are no londer included in AliITSgeom.h
+
+Revision 1.30  2001/01/30 09:23:13  hristov
+Streamers removed (R.Brun)
+
+Revision 1.29  2001/01/26 20:01:09  hristov
+Major upgrade of AliRoot code
+
+Revision 1.28  2000/12/18 14:02:00  barbera
+new version of the ITS tracking to take into account the new TPC track parametrization
+
+Revision 1.27  2000/12/08 13:49:27  barbera
+Hidden declaration in a for loop removed to be compliant with HP-UX compiler
+
+Revision 1.26  2000/11/27 13:12:13  barbera
+New version containing the files for tracking
+
+Revision 1.25  2000/11/12 22:38:05  barbera
+Added header file for the SPD Bari model
+
+Revision 1.24  2000/10/09 22:18:12  barbera
+Bug fixes from MAriana to le AliITStest.C run correctly
+
+Revision 1.23  2000/10/05 20:47:42  nilsen
+fixed dependencies of include files. Tryed but failed to get a root automaticly
+generates streamer function to work. Modified SetDefaults.
+
+Revision 1.9.2.15  2000/10/04 16:56:40  nilsen
+Needed to include stdlib.h
+
+=======
+Revision 1.22  2000/10/04 19:45:52  barbera
+Corrected by F. Carminati for v3.04
+
+Revision 1.21  2000/10/02 21:28:08  fca
+Removal of useless dependecies via forward declarations
+
+Revision 1.20  2000/10/02 16:31:39  barbera
+General code clean-up
+
 Revision 1.9.2.14  2000/10/02 15:43:51  barbera
 General code clean-up (e.g., printf -> cout)
 
@@ -93,18 +202,21 @@ the AliITS class.
 // futher information.
 //
 ///////////////////////////////////////////////////////////////////////////////
+#include <iostream.h>
+#include <iomanip.h>
+#include <fstream.h>
+#include <stdlib.h>
 #include <TMath.h>
 #include <TRandom.h>
 #include <TBranch.h>
 #include <TVector.h>
-#include <TObjArray.h>
+#include <TClonesArray.h>
 #include <TROOT.h>
 #include <TObjectTable.h>
 #include <TFile.h>
 #include <TTree.h>
 #include <TString.h>
-
+#include <TParticle.h>
 
 
 #include "AliRun.h"
@@ -112,9 +224,14 @@ the AliITS class.
 #include "AliITSMap.h"
 #include "AliITSDetType.h"
 #include "AliITSClusterFinder.h"
-#include "AliITSsimulation.h"
+//#include "AliITSsimulation.h"
+#include "AliITSsimulationSPD.h"
+#include "AliITSsimulationSDD.h"
+#include "AliITSsimulationSSD.h"
+#include "AliITSresponse.h"
 #include "AliITSsegmentationSPD.h"
 #include "AliITSresponseSPD.h"
+#include "AliITSresponseSPDbari.h"
 #include "AliITSsegmentationSDD.h"
 #include "AliITSresponseSDD.h"
 #include "AliITSsegmentationSSD.h"
@@ -125,8 +242,18 @@ the AliITS class.
 #include "AliITSmodule.h"
 #include "AliITSRecPoint.h"
 #include "AliITSRawCluster.h"
-
-const Int_t AliITS::fgkNTYPES=3;
+#include "AliMC.h"
+#include "stdlib.h"
+#include "AliKalmanTrack.h" 
+#include "AliMagF.h"
+
+#include "AliITStrack.h"
+#include "AliITSiotrack.h"
+#include "AliITStracking.h"
+#include "AliITSRad.h"   
+#include "AliTPC.h"
+#include "AliTPCParam.h"
+#include "AliTPCtracker.h"
 
 ClassImp(AliITS)
  
@@ -145,7 +272,7 @@ AliITS::AliITS() : AliDetector() {
   fIshunt     = 0;
   fEuclidOut  = 0;
 
-  //fNDetTypes = fgkNTYPES;
+  fNDetTypes = kNTYPES;
   fIdN        = 0;
   fIdName     = 0;
   fIdSens     = 0;
@@ -183,13 +310,13 @@ AliITS::AliITS(const char *name, const char *title):AliDetector(name,title){
   fHits       = new TClonesArray("AliITShit", 1560);
   gAlice->AddHitList(fHits);
 
-  //fNDetTypes = fgkNTYPES;
+  fNDetTypes = kNTYPES;
 
-  fNdtype = new Int_t[fgkNTYPES];
-  fDtype = new TObjArray(fgkNTYPES);
+  fNdtype = new Int_t[kNTYPES];
+  fDtype = new TObjArray(kNTYPES);
 
-  fNctype = new Int_t[fgkNTYPES];
-  fCtype = new TObjArray(fgkNTYPES);
+  fNctype = new Int_t[kNTYPES];
+  fCtype = new TObjArray(kNTYPES);
 
 
   fRecPoints = 0;
@@ -205,11 +332,11 @@ AliITS::AliITS(const char *name, const char *title):AliDetector(name,title){
   fIdName     = 0;
   fIdSens     = 0;
  
-  fDetTypes = new TObjArray(fgkNTYPES);  
+  fDetTypes = new TObjArray(kNTYPES);  
 
   Int_t i;
-  for(i=0;i<fgkNTYPES;i++) {
-    (*fDetTypes)[i]=new AliITSDetType(); 
+  for(i=0;i<kNTYPES;i++) {
+    fDetTypes->AddAt(new AliITSDetType(),i); 
     fNdtype[i]=0;
     fNctype[i]=0;
    }
@@ -256,7 +383,8 @@ AliITS::~AliITS(){
   delete fHits;
   delete fDigits;
   delete fRecPoints;
-  if(fIdName!=0) delete[] fIdName;
+//  delete fIdName;        // TObjArray of TObjStrings
+  if(fIdName!=0) delete[] fIdName;  // Array of TStrings
   if(fIdSens!=0) delete[] fIdSens;
   if(fITSmodules!=0) {
       this->ClearModules();
@@ -447,7 +575,7 @@ void AliITS::ResetDigits()
     if (!fDtype) return;
 
     Int_t i;
-    for (i=0;i<fgkNTYPES;i++ ) {
+    for (i=0;i<kNTYPES;i++ ) {
        if ((*fDtype)[i])    ((TClonesArray*)(*fDtype)[i])->Clear();
        if (fNdtype)  fNdtype[i]=0;
     }
@@ -472,7 +600,7 @@ void AliITS::ResetClusters()
     //
 
     Int_t i;
-    for (i=0;i<fgkNTYPES;i++ ) {
+    for (i=0;i<kNTYPES;i++ ) {
        if ((*fCtype)[i])    ((TClonesArray*)(*fCtype)[i])->Clear();
        if (fNctype)  fNctype[i]=0;
     }
@@ -524,18 +652,13 @@ void AliITS::Init(){
   // and sets the default segmentation, response, digit and raw cluster classes
   // Therefore it should be called after a call to CreateGeometry.
   //
+  Int_t i;
 
-
+//
   SetDefaults();
-
-  Int_t i;
-  cout << endl;
-  for(i=0;i<30;i++) cout << "*";cout << " ITS_INIT ";
-  for(i=0;i<30;i++) cout << "*";cout << endl;
+// Array of TStrings
   for(i=0;i<fIdN;i++) fIdSens[i] = gMC->VolId(fIdName[i]);
-  //
-  for(i=0;i<70;i++) cout << "*";
-  cout << endl;
+//
 }
 
 //_____________________________________________________________________________
@@ -543,50 +666,59 @@ void AliITS::SetDefaults()
 {
   // sets the default segmentation, response, digit and raw cluster classes
 
-  cout << "AliITS::SetDefaults" << endl;
+  if(fDebug) printf("%s: SetDefaults\n",ClassName());
 
   AliITSDetType *iDetType;
 
+
   //SPD 
 
-  AliITSsegmentationSPD *seg0=new AliITSsegmentationSPD(fITSgeom);
-  AliITSresponseSPD *resp0=new AliITSresponseSPD();
   iDetType=DetType(0); 
-  if (!iDetType->GetSegmentationModel()) SetSegmentationModel(0,seg0); 
-  if (!iDetType->GetResponseModel()) SetResponseModel(0,resp0); 
+  if (!iDetType->GetSegmentationModel()) {
+    AliITSsegmentationSPD *seg0=new AliITSsegmentationSPD(fITSgeom);
+    SetSegmentationModel(0,seg0); 
+  }
+  if (!iDetType->GetResponseModel()) {
+     SetResponseModel(0,new AliITSresponseSPD()); 
+  }
   // set digit and raw cluster classes to be used
-  const char *kData0=resp0->DataType();
+  
+  const char *kData0=(iDetType->GetResponseModel())->DataType();
   if (strstr(kData0,"real")) {
       iDetType->ClassNames("AliITSdigit","AliITSRawClusterSPD");
   } else iDetType->ClassNames("AliITSdigitSPD","AliITSRawClusterSPD");
 
   // SDD                                         //
-  AliITSresponseSDD *resp1=new AliITSresponseSDD();
-  AliITSsegmentationSDD *seg1=new AliITSsegmentationSDD(fITSgeom,resp1);
   iDetType=DetType(1); 
-  //printf("SetDefaults: iDetType %p\n",iDetType);
-  if (!iDetType->GetSegmentationModel()) SetSegmentationModel(1,seg1); 
-  //printf("SetDefaults: segm %p\n",iDetType->GetSegmentationModel());
-  if (!iDetType->GetResponseModel()) SetResponseModel(1,resp1); 
-  //printf("SetDefaults: resp %p\n",iDetType->GetResponseModel());
-  const char *kData1=resp1->DataType();
-  const char *kopt=resp1->ZeroSuppOption();
+  if (!iDetType->GetResponseModel()) {
+    SetResponseModel(1,new AliITSresponseSDD()); 
+  }
+  AliITSresponse *resp1=iDetType->GetResponseModel();
+  if (!iDetType->GetSegmentationModel()) {
+    AliITSsegmentationSDD *seg1=new AliITSsegmentationSDD(fITSgeom,resp1);
+    SetSegmentationModel(1,seg1); 
+  }
+  const char *kData1=(iDetType->GetResponseModel())->DataType();
+  const char *kopt=iDetType->GetResponseModel()->ZeroSuppOption();
   if ((!strstr(kopt,"2D")) && (!strstr(kopt,"1D")) || strstr(kData1,"real") ) {
       iDetType->ClassNames("AliITSdigit","AliITSRawClusterSDD");
   } else iDetType->ClassNames("AliITSdigitSDD","AliITSRawClusterSDD");
 
   // SSD
-  AliITSsegmentationSSD *seg2=new AliITSsegmentationSSD(fITSgeom);
-  AliITSresponseSSD *resp2=new AliITSresponseSSD();
   iDetType=DetType(2); 
-  if (!iDetType->GetSegmentationModel()) SetSegmentationModel(2,seg2); 
-  if (!iDetType->GetResponseModel()) SetResponseModel(2,resp2); 
-  const char *kData2=resp2->DataType();
+  if (!iDetType->GetSegmentationModel()) {
+    AliITSsegmentationSSD *seg2=new AliITSsegmentationSSD(fITSgeom);
+    SetSegmentationModel(2,seg2); 
+  }
+  if (!iDetType->GetResponseModel()) {
+    SetResponseModel(2,new AliITSresponseSSD()); 
+  }
+  const char *kData2=(iDetType->GetResponseModel())->DataType();
   if (strstr(kData2,"real")) {
       iDetType->ClassNames("AliITSdigit","AliITSRawClusterSSD");
   } else iDetType->ClassNames("AliITSdigitSSD","AliITSRawClusterSSD");
 
-  if (fgkNTYPES>3) {
+  if (kNTYPES>3) {
     Warning("SetDefaults","Only the three basic detector types are initialised!");
   } 
 
@@ -609,26 +741,33 @@ void AliITS::MakeTreeC(Option_t *option)
 {
   // create a separate tree to store the clusters
 
-  cout << "AliITS::MakeTreeC" << endl;
+//  cout << "AliITS::MakeTreeC" << endl;
 
-     char *optC = strstr(option,"C");
+     const char *optC = strstr(option,"C");
      if (optC && !fTreeC) fTreeC = new TTree("TC","Clusters in ITS");
      else return;
 
      Int_t buffersize = 4000;
      char branchname[30];
 
-     char *det[3] = {"SPD","SDD","SSD"};
+     const char *det[3] = {"SPD","SDD","SSD"};
+
+     char digclass[40];
+     char clclass[40];
 
      // one branch for Clusters per type of detector
      Int_t i;
-     for (i=0; i<fgkNTYPES ;i++) {
-        if (fgkNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
+     for (i=0; i<kNTYPES ;i++) {
+       AliITSDetType *iDetType=DetType(i); 
+       iDetType->GetClassNames(digclass,clclass);
+       // clusters
+       (*fCtype)[i] = new TClonesArray(clclass,10000); 
+        if (kNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
        else  sprintf(branchname,"%sClusters%d",GetName(),i+1);
        if (fCtype   && fTreeC) {
           TreeC()->Branch(branchname,&((*fCtype)[i]), buffersize);
-          cout << "Making Branch " << branchname;
-          cout << " for Clusters of detector type " << i+1 << endl;
+//        cout << "Making Branch " << branchname;
+//        cout << " for Clusters of detector type " << i+1 << endl;
        }       
      }
 
@@ -638,13 +777,13 @@ void AliITS::MakeTreeC(Option_t *option)
 void AliITS::GetTreeC(Int_t event)
 {
 
-  cout << "AliITS::GetTreeC" << endl;
+//  cout << "AliITS::GetTreeC" << endl;
 
   // get the clusters tree for this event and set the branch address
     char treeName[20];
     char branchname[30];
 
-    char *det[3] = {"SPD","SDD","SSD"};
+    const char *det[3] = {"SPD","SDD","SSD"};
 
     ResetClusters();
     if (fTreeC) {
@@ -654,12 +793,11 @@ void AliITS::GetTreeC(Int_t event)
     sprintf(treeName,"TreeC%d",event);
     fTreeC = (TTree*)gDirectory->Get(treeName);
 
-
     TBranch *branch;
     if (fTreeC) {
         Int_t i;
-       for (i=0; i<fgkNTYPES; i++) {
-          if (fgkNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
+       for (i=0; i<kNTYPES; i++) {
+          if (kNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
           else  sprintf(branchname,"%sClusters%d",GetName(),i+1);
           if (fCtype) {
                branch = fTreeC->GetBranch(branchname);
@@ -673,61 +811,66 @@ void AliITS::GetTreeC(Int_t event)
 
 }
 //_____________________________________________________________________________
-void AliITS::MakeBranch(Option_t* option){
+void AliITS::MakeBranch(Option_t* option, const char *file)
+{
   //
   // Creates Tree branches for the ITS.
   //
-
-
+  //
   Int_t buffersize = 4000;
   char branchname[30];
   sprintf(branchname,"%s",GetName());
 
-  AliDetector::MakeBranch(option);
+  AliDetector::MakeBranch(option,file);
 
+  const char *cD = strstr(option,"D");
+  const char *cR = strstr(option,"R");
 
-// one branch for digits per type of detector
-  
-   char *det[3] = {"SPD","SDD","SSD"};
+  if (cD) {
+  //
+  // one branch for digits per type of detector
+  //
+   const char *det[3] = {"SPD","SDD","SSD"};
 
    char digclass[40];
    char clclass[40];
 
    Int_t i;
-   for (i=0; i<fgkNTYPES ;i++) {
+   for (i=0; i<kNTYPES ;i++) {
        AliITSDetType *iDetType=DetType(i); 
        iDetType->GetClassNames(digclass,clclass);
-       //printf("i, digclass, recclass %d %s %s\n",i,digclass,clclass); 
        // digits
-       (*fDtype)[i] = new TClonesArray(digclass,10000); 
-       // clusters
-       (*fCtype)[i] = new TClonesArray(clclass,10000); 
+       if(!((*fDtype)[i])) (*fDtype)[i] = new TClonesArray(digclass,10000);
+       else ResetDigits(i);
    }
 
-
-  for (i=0; i<fgkNTYPES ;i++) {
-      if (fgkNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
-      else  sprintf(branchname,"%sDigits%d",GetName(),i+1);
-      
-      if (fDtype   && gAlice->TreeD()) {
-         gAlice->TreeD()->Branch(branchname,&((*fDtype)[i]), buffersize);
-         cout << "Making Branch " << branchname;
-         cout << " for digits of type "<< i+1 << endl;
+   for (i=0; i<kNTYPES ;i++) {
+      if (kNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
+      else  sprintf(branchname,"%sDigits%d",GetName(),i+1);      
+      if (fDtype && gAlice->TreeD()) {
+        MakeBranchInTree(gAlice->TreeD(), 
+                         branchname, &((*fDtype)[i]), buffersize, file);
+//         cout << "Making Branch " << branchname;
+//         cout << " for digits of type "<< i+1 << endl;
       }        
+    }
   }
 
+  if (cR) {
+  //
   // only one branch for rec points for all detector types
-  sprintf(branchname,"%sRecPoints",GetName());
-
-  fRecPoints=new TClonesArray("AliITSRecPoint",10000);
+  //
+    sprintf(branchname,"%sRecPoints",GetName());
+    if(!fRecPoints) fRecPoints=new TClonesArray("AliITSRecPoint",10000);
 
-  if (fRecPoints && gAlice->TreeR()) {
-    gAlice->TreeR()->Branch(branchname,&fRecPoints, buffersize);
-    cout << "Making Branch " << branchname;
-    cout << " for reconstructed space points" << endl;
+    if (fRecPoints && gAlice->TreeR()) {
+      MakeBranchInTree(gAlice->TreeR(), 
+                               branchname, &fRecPoints, buffersize, file) ;
+//      cout << "Making Branch " << branchname;
+//      cout << " for reconstructed space points" << endl;
+    }
   }    
-
-
 }
 
 //___________________________________________
@@ -739,7 +882,7 @@ void AliITS::SetTreeAddress()
   char branchname[30];
   AliDetector::SetTreeAddress();
 
-  char *det[3] = {"SPD","SDD","SSD"};
+  const char *det[3] = {"SPD","SDD","SSD"};
 
   TBranch *branch;
   TTree *treeD = gAlice->TreeD();
@@ -747,8 +890,8 @@ void AliITS::SetTreeAddress()
 
   Int_t i;
   if (treeD) {
-      for (i=0; i<fgkNTYPES; i++) {
-         if (fgkNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
+      for (i=0; i<kNTYPES; i++) {
+         if (kNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
          else  sprintf(branchname,"%sDigits%d",GetName(),i+1);
          if (fDtype) {
              branch = treeD->GetBranch(branchname);
@@ -760,10 +903,8 @@ void AliITS::SetTreeAddress()
  
   if (treeR) {
     sprintf(branchname,"%sRecPoints",GetName());
-    if (fRecPoints) {
       branch = treeR->GetBranch(branchname);
       if (branch) branch->SetAddress(&fRecPoints);
-    }
   }
   
 
@@ -817,14 +958,14 @@ void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t nmodules,Option_t *option,
 
     static Bool_t first=kTRUE;
     static TFile *file;
-    char *addBgr = strstr(option,"Add");
+    const char *addBgr = strstr(option,"Add");
 
 
     if (addBgr ) {
        if(first) {
-           cout<<"filename "<<filename<<endl;
+//         cout<<"filename "<<filename<<endl;
            file=new TFile(filename);
-           cout<<"I have opened "<<filename<<" file "<<endl;
+//         cout<<"I have opened "<<filename<<" file "<<endl;
            fHits2     = new TClonesArray("AliITShit",1000  );
        }           
        first=kFALSE;
@@ -858,8 +999,8 @@ void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t nmodules,Option_t *option,
        //printf("background - ntracks1 - %d\n",ntracks1);
    }
 
-    Int_t npart = gAlice->GetEvent(evnt);
-    if(npart<=0) return;
+    //Int_t npart = gAlice->GetEvent(evnt);
+    //if(npart<=0) return;
     TClonesArray *itsHits = this->Hits();
     Int_t lay,lad,det,index;
     AliITShit *itsHit=0;
@@ -926,6 +1067,70 @@ void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t nmodules,Option_t *option,
 
 }
 
+//____________________________________________________________________________
+
+void AliITS::SDigits2Digits()
+{
+  
+  AliITSgeom *geom = GetITSgeom();
+  
+  // SPD
+  AliITSDetType *iDetType;
+  iDetType=DetType(0);
+  AliITSsegmentationSPD *seg0=(AliITSsegmentationSPD*)iDetType->GetSegmentationModel();
+  AliITSresponseSPD *res0 = (AliITSresponseSPD*)iDetType->GetResponseModel();
+  AliITSsimulationSPD *sim0=new AliITSsimulationSPD(seg0,res0);
+  SetSimulationModel(0,sim0);
+  // test
+  // printf("SPD dimensions %f %f \n",seg0->Dx(),seg0->Dz());
+  // printf("SPD npixels %d %d \n",seg0->Npz(),seg0->Npx());
+  // printf("SPD pitches %d %d \n",seg0->Dpz(0),seg0->Dpx(0));
+  // end test
+  // 
+  // SDD
+  //Set response functions
+  // SDD compression param: 2 fDecrease, 2fTmin, 2fTmax or disable, 2 fTolerance
+
+  iDetType=DetType(1);
+  AliITSresponseSDD *res1 = (AliITSresponseSDD*)iDetType->GetResponseModel();
+  if (!res1) {
+    res1=new AliITSresponseSDD();
+    SetResponseModel(1,res1);
+  }
+  Float_t noise, baseline;
+  res1->GetNoiseParam(noise,baseline);
+  Float_t noise_after_el =  res1->GetNoiseAfterElectronics();
+  Float_t fCutAmp = baseline + 2.*noise_after_el;
+  Int_t cp[8]={0,0,(int)fCutAmp,(int)fCutAmp,0,0,0,0}; //1D
+  res1->SetCompressParam(cp);
+  AliITSsegmentationSDD *seg1=(AliITSsegmentationSDD*)iDetType->GetSegmentationModel();
+  if (!seg1) {
+    seg1 = new AliITSsegmentationSDD(geom,res1);
+    SetSegmentationModel(1,seg1);
+  }
+  AliITSsimulationSDD *sim1=new AliITSsimulationSDD(seg1,res1);
+  SetSimulationModel(1,sim1);
+
+  // SSD
+  iDetType=DetType(2);
+  AliITSsegmentationSSD *seg2=(AliITSsegmentationSSD*)iDetType->GetSegmentationModel();
+  AliITSresponseSSD *res2 = (AliITSresponseSSD*)iDetType->GetResponseModel();
+  res2->SetSigmaSpread(3.,2.);
+  AliITSsimulationSSD *sim2=new AliITSsimulationSSD(seg2,res2);
+  SetSimulationModel(2,sim2);
+  cerr<<"Digitizing ITS...\n";
+  TStopwatch timer;
+  timer.Start();
+  HitsToDigits(gAlice->GetEvNumber(),0,-1," ","All"," ");
+  timer.Stop(); timer.Print();
+
+  delete sim0;
+  delete sim1;
+  delete sim2;
+}
+
 
 //____________________________________________________________________________
 void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *option, Option_t *opt,Text_t *filename)
@@ -937,14 +1142,16 @@ void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *optio
    // initialised for all versions - for the moment it is only for v5 !
    // 7 is the SDD beam test version  
    Int_t ver = this->IsVersion(); 
-   if(ver!=5 && ver!=7) return; 
-
-   char *all = strstr(opt,"All");
-   char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
+   if(ver!=5 && ver!=7 && ver!=8 && ver!=9) return; 
 
+   const char *all = strstr(opt,"All");
+   const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
+//   cout<<" 1 AliITS "<<endl;
    Int_t nmodules;
    InitModules(size,nmodules); 
+//   cout<<" 2 AliITS "<<endl;
    FillModules(evNumber,bgrev,nmodules,option,filename);
+//   cout<<" 3 AliITS "<<endl;
 
    //TBranch *branch;
    AliITSsimulation* sim;
@@ -952,24 +1159,34 @@ void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *optio
    AliITSgeom *geom = GetITSgeom();
 
    Int_t id,module;
+//   Int_t lay, lad, detect;
    Int_t first,last;
-   for (id=0;id<fgkNTYPES;id++) {
+   for (id=0;id<kNTYPES;id++) {
         if (!all && !det[id]) continue;
        //branch = (TBranch*)branches->UncheckedAt(id);
        AliITSDetType *iDetType=DetType(id); 
        sim = (AliITSsimulation*)iDetType->GetSimulationModel();
+
        if (!sim) {
            Error("HitsToDigits","The simulation class was not instantiated!");
            exit(1);
           // or SetDefaultSimulation();
        }
+
        if(geom) {
          first = geom->GetStartDet(id);
          last = geom->GetLastDet(id);
        } else first=last=0;
-       cout << "det type " << id << " first, last "<< first << last << endl;
+//     cout << "det type " << id << " first, last "<< first << last << endl;
        for(module=first;module<=last;module++) {
            AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
+/*
+           geom->GetModuleId(module,lay, lad, detect);
+            if ( lay == 6 )
+             ((AliITSsegmentationSSD*)(((AliITSsimulationSSD*)sim)->GetSegmentation()))->SetLayer(6);
+           if ( lay == 5 )
+             ((AliITSsegmentationSSD*)(((AliITSsimulationSSD*)sim)->GetSegmentation()))->SetLayer(5);
+*/
            sim->DigitiseModule(mod,module,evNumber);
            // fills all branches - wasted disk space
            gAlice->TreeD()->Fill(); 
@@ -982,12 +1199,13 @@ void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size, Option_t *optio
 
    ClearModules();
 
-   Int_t nentries=(Int_t)gAlice->TreeD()->GetEntries();
-   cout << "nentries in TreeD" << nentries << endl;
+//   Int_t nentries=(Int_t)
+   gAlice->TreeD()->GetEntries();
+//   cout << "nentries in TreeD" << nentries << endl;
 
    char hname[30];
    sprintf(hname,"TreeD%d",evNumber);
-   gAlice->TreeD()->Write(hname);
+   gAlice->TreeD()->Write(hname,TObject::kOverwrite);
    // reset tree
    gAlice->TreeD()->Reset();
 
@@ -1003,18 +1221,19 @@ void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
    // initialised for all versions - for the moment it is only for v5 !
    // 7 is the SDD beam test version  
    Int_t ver = this->IsVersion(); 
-   if(ver!=5 && ver!=7) return; 
+   if(ver!=5 && ver!=8 && ver!=9) return;
 
-   char *all = strstr(opt,"All");
-   char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
+   const char *all = strstr(opt,"All");
+   const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
 
    static Bool_t first=kTRUE;
-   if (first) {
+   if (!TreeC() && first) {
        MakeTreeC("C");
        first=kFALSE;
    }
+
+   TTree *treeC=TreeC();
  
-   TTree *iTC=TreeC();
 
    //TBranch *branch;
    AliITSClusterFinder* rec;
@@ -1023,7 +1242,7 @@ void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
    AliITSgeom *geom = GetITSgeom();
 
    Int_t id,module;
-   for (id=0;id<fgkNTYPES;id++) {
+   for (id=0;id<kNTYPES;id++) {
         if (!all && !det[id]) continue;
        //branch = (TBranch*)branches->UncheckedAt(id);
        AliITSDetType *iDetType=DetType(id); 
@@ -1040,16 +1259,16 @@ void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
          first = geom->GetStartDet(id);
          last = geom->GetLastDet(id);
        } else first=last=0;
-       //printf("first last %d %d\n",first,last);
+       printf("first module - last module %d %d\n",first,last);
        for(module=first;module<=last;module++) {
               this->ResetDigits();
               if (all) gAlice->TreeD()->GetEvent(lastentry+module);
              else gAlice->TreeD()->GetEvent(lastentry+(module-first));
              Int_t ndigits = itsDigits->GetEntriesFast();
-             if (ndigits) rec->FindRawClusters();
+             if (ndigits) rec->FindRawClusters(module);
              gAlice->TreeR()->Fill(); 
              ResetRecPoints();
-             iTC->Fill();
+             treeC->Fill();
               ResetClusters();
              // try and fill only the branch 
              //branch->Fill();
@@ -1058,22 +1277,22 @@ void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt)
    } // loop over detector types
 
 
-   Int_t nentries=(Int_t)gAlice->TreeR()->GetEntries();
-   Int_t ncentries=(Int_t)iTC->GetEntries();
-   cout << " nentries ncentries " << nentries << ncentries <<  endl;
+//   Int_t nentries=(Int_t)
+   gAlice->TreeR()->GetEntries();
+//   Int_t ncentries=(Int_t)
+   treeC->GetEntries();
+//   cout << " nentries ncentries " << nentries << ncentries <<  endl;
 
    char hname[30];
    sprintf(hname,"TreeR%d",evNumber);
-   gAlice->TreeR()->Write(hname);
+   gAlice->TreeR()->Write(hname,TObject::kOverwrite);
    // reset tree
    gAlice->TreeR()->Reset();
 
    sprintf(hname,"TreeC%d",evNumber);
-   iTC->Write(hname);
-   iTC->Reset();
+   treeC->Write(hname,TObject::kOverwrite);
+   treeC->Reset();
 }
-
-
 //____________________________________________________________________________
 void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t size,
 Option_t *option,Option_t *opt,Text_t *filename)
@@ -1085,10 +1304,11 @@ Option_t *option,Option_t *opt,Text_t *filename)
    // the condition below will disappear when the geom class will be
    // initialised for all versions - for the moment it is only for v5 !  
    Int_t ver = this->IsVersion(); 
-   if(ver!=5) return; 
+   if(ver!=5 && ver!=8 && ver!=9) return;
+   //if(ver!=5) return; 
 
-   char *all = strstr(opt,"All");
-   char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
+   const char *all = strstr(opt,"All");
+   const char *det[3] ={strstr(opt,"SPD"),strstr(opt,"SDD"),strstr(opt,"SSD")};
 
    Int_t nmodules;
    InitModules(size,nmodules);
@@ -1111,17 +1331,23 @@ Option_t *option,Option_t *opt,Text_t *filename)
 
 
    Int_t id,module;
-   for (id=0;id<fgkNTYPES;id++) {
+   for (id=0;id<kNTYPES;id++) {
         if (!all && !det[id]) continue;
        AliITSDetType *iDetType=DetType(id); 
        sim = (AliITSsimulation*)iDetType->GetSimulationModel();
        if (!sim) {
-           Error("HitsToFastPoints","The simulation class was not instantiated!");
+           Error("HitsToFastPoints",
+                "The simulation class was not instantiated!");
            exit(1);
           // or SetDefaultSimulation();
        }
-       Int_t first = geom->GetStartDet(id);
-       Int_t last = geom->GetLastDet(id);
+
+        Int_t first,last;
+        if(geom) {
+         first = geom->GetStartDet(id);
+         last = geom->GetLastDet(id);
+        } else first=last=0;
+        printf("first module - last module %d %d\n",first,last);
        for(module=first;module<=last;module++) {
            AliITSmodule *mod = (AliITSmodule *)fITSmodules->At(module);
            sim->CreateFastRecPoints(mod,module,random);
@@ -1137,67 +1363,478 @@ Option_t *option,Option_t *opt,Text_t *filename)
 
    char hname[30];
    sprintf(hname,"TreeR%d",evNumber);
-   gAlice->TreeR()->Write(hname);
+   gAlice->TreeR()->Write(hname,TObject::kOverwrite);
    // reset tree
    gAlice->TreeR()->Reset();
 
    delete [] random;
 
 }
+//________________________________________________________________
+AliITStrack  AliITS::Tracking(AliITStrack &track, AliITStrack *reference,
+                             TObjArray *fastpoints, Int_t **vettid,
+                             Bool_t flagvert,  AliITSRad *rl  ) {
+// Origin  A. Badala' and G.S. Pappalardo:  e-mail Angela.Badala@ct.infn.it,
+// Giuseppe.S.Pappalardo@ct.infn.it
 
+  TList *list= new TList();   
 
-//____________________________________________________________________________
-void AliITS::Streamer(TBuffer &R__b){
-   // Stream an object of class AliITS.
-
-   Int_t i;
+  AliITStrack tr(track);
+  
+  list->AddLast(&tr);
+  
+  Double_t Pt=(tr).GetPt();
+//  cout << "\n Pt = " << Pt <<"\n";  //stampa
+
+  AliITStracking obj(list, reference, this, fastpoints,
+                    TMath::Abs(Pt),vettid, flagvert, rl);
+  list->Delete();
+  delete list;
+
+  Int_t itot=-1;
+  TVector VecTotLabref(18);
+  Int_t lay, k;
+  for(lay=5; lay>=0; lay--) {
+    TVector VecLabref(3); 
+    VecLabref=(*reference).GetLabTrack(lay);
+    Float_t ClustZ=(*reference).GetZclusterTrack( lay); //aggiunta il 5-3-2001
+    for(k=0; k<3; k++){
+       //{itot++; VecTotLabref(itot)=VecLabref(k);} //cambiata il 5-3-2001 
+       Int_t lpp=(Int_t)VecLabref(k);
+       if(lpp>=0) {
+           TParticle *p=(TParticle*) gAlice->Particle(lpp);
+           Int_t pcode=p->GetPdgCode();
+           if(pcode==11) VecLabref(k)=p->GetFirstMother();
+       } // end if lpp>=0    
+       itot++; VecTotLabref(itot)=VecLabref(k);
+       if(VecLabref(k)==0. && ClustZ == 0.) VecTotLabref(itot) =-3.; }  
+  } // end for lay
+  Long_t labref;
+  Int_t freq;  
+  (*reference).Search(VecTotLabref, labref, freq);
+    
+  //if(freq < 4) labref=-labref;
+  //if(freq < 6) labref=-labref;        // cinque - sei
+  if(freq < 5) labref=-labref;        // cinque - sei  
+  (*reference).SetLabel(labref);
 
-   if (R__b.IsReading()) {
-      Version_t R__v = R__b.ReadVersion();
-      if (R__v) {
-         AliDetector::Streamer(R__b);
-         R__b >> fIdN;
-         R__b.ReadArray(fIdSens); 
-         for(i=0;i<fIdN;i++) fIdName[i].Streamer(R__b);
-         R__b >> fITSgeom;
-         R__b >> fITSmodules;
-         R__b >> fEuclidOut;
-         R__b >> fMajorVersion;
-         R__b >> fMinorVersion;
-         R__b >> fDetTypes;
-         R__b >> fDtype;
-         delete []fNdtype; 
-         fNdtype = new Int_t[fgkNTYPES];
-         R__b.ReadFastArray(fNdtype,fgkNTYPES);
-         R__b >> fCtype;
-         delete []fNctype; 
-         fNctype = new Int_t[fgkNTYPES];
-         R__b.ReadFastArray(fNctype,fgkNTYPES);
-         R__b >> fRecPoints;
-         R__b >> fNRecPoints;
-         R__b >> fTreeC;
-      } // end if R__v
-   } else { // writing
-      R__b.WriteVersion(AliITS::IsA());
-      AliDetector::Streamer(R__b);
-      R__b << fIdN;
-      R__b.WriteArray(fIdSens,fIdN); 
-      for(i=0;i<fIdN;i++) fIdName[i].Streamer(R__b);
-      R__b << fITSgeom;
-      R__b << fITSmodules;
-      R__b << fEuclidOut;
-      R__b << fMajorVersion;
-      R__b << fMinorVersion;
-      R__b << fDetTypes;
-      R__b << fDtype;
-      R__b.WriteFastArray(fNdtype,fgkNTYPES);
-      R__b << fCtype;
-      R__b.WriteFastArray(fNctype,fgkNTYPES);
-      R__b << fRecPoints;
-      R__b << fNRecPoints;
-      R__b << fTreeC;
-   } // end if
+  return *reference; 
 
 }
+//________________________________________________________________------
+void AliITS::DoTracking(Int_t evNumber, Int_t min_t, Int_t max_t,
+                       TFile *file, Bool_t flagvert) {
+    //   ex macro for tracking ITS
+
+    printf("begin DoTracking - file %p\n",file);
+
+    //const char *pname="75x40_100x60";
+  
+    Int_t imax=200,jmax=450;
+    AliITSRad *rl = new AliITSRad(imax,jmax);
+    //cout<<" dopo costruttore AliITSRad\n"; getchar();
+
+    struct GoodTrack {
+       Int_t lab,code;
+       Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
+       Bool_t flag;
+    };
+
+    gAlice->GetEvent(0);
+
+    AliKalmanTrack *kkprov;
+    kkprov->SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor());  
+
+    /*  //modificato il 26-4-2001
+       AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC");
+       AliTPCParam *digp = (AliTPCParam*)file->Get(pname);
+       if (digp!=0) TPC->SetParam(digp);
+    */
+    TFile *cf=TFile::Open("AliTPCclusters.root");  
+    AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60");
+    if (!digp) { cerr<<"TPC parameters have not been found !\n"; getchar();}
+
+    AliTPCtracker *tracker = new AliTPCtracker(digp); //aggiunto il 23-5 
+    //aggiunto il 23-5 
+    // Load clusters
+    tracker->LoadInnerSectors();
+    tracker->LoadOuterSectors();
+
+    GoodTrack gt[15000];
+    Int_t ngood=0;
+    ifstream in("itsgood_tracks");
+
+    cerr<<"Reading itsgood tracks...\n";
+    while (in>>gt[ngood].lab>>gt[ngood].code
+          >>gt[ngood].px >>gt[ngood].py>>gt[ngood].pz
+          >>gt[ngood].x  >>gt[ngood].y >>gt[ngood].z
+          >>gt[ngood].pxg  >>gt[ngood].pyg >>gt[ngood].pzg
+          >>gt[ngood].ptg >>gt[ngood].flag) {
+       ngood++;
+       cerr<<ngood<<'\r';
+       if (ngood==15000) {
+           cerr<<"Too many good tracks !\n";
+           break;
+       } // end if ngood==1500
+    } // end while
+    if (!in.eof()) cerr<<"Read error (itsgood_tracks) !\n";
+
+    // Load tracks
+    // TFile *tf=TFile::Open("tpctracks.root"); //commentato 26-4-2001 
+    TFile *tf=TFile::Open("AliTPCtracks.root"); 
+    if (!tf->IsOpen()) {cerr<<"Can't open AliTPCtracks.root !\n"; return ;}
+    TObjArray tracks(200000);
+    //TTree *tracktree=(TTree*)tf->Get("TreeT");
+    TTree *tracktree=(TTree*)tf->Get("TPCf"); //aggiunto il 23-5
+    if (!tracktree) {cerr<<"Can't get a tree with TPC tracks !\n";} 
+    TBranch *tbranch=tracktree->GetBranch("tracks");
+    Int_t nentr=(Int_t)tracktree->GetEntries();
+    Int_t kk;
+    /*  commentato il 26-4-2001
+       for (kk=0; kk<nentr; kk++) {
+       AliTPCtrack *iotrack=new AliTPCtrack;
+       tbranch->SetAddress(&iotrack);
+       tracktree->GetEvent(kk);
+       tracks.AddLast(iotrack);
+       } 
+    */
+    AliTPCtrack *iotracktpc=0;    
+    for (kk=0; kk<nentr; kk++) {
+       iotracktpc=new AliTPCtrack; 
+       tbranch->SetAddress(&iotracktpc);
+       tracktree->GetEvent(kk);
+       tracker->CookLabel(iotracktpc,0.1); //aggiunto 23-5
+       tracks.AddLast(iotracktpc);    
+    } // end for kk
+    delete tracker;
+    tf->Close();
+
+    Int_t nt = tracks.GetEntriesFast();
+    cerr<<"Number of found tracks "<<nt<<endl;
+
+    TVector DataOut(9);
+    Int_t kkk=0;
+
+    Double_t ptg=0.,pxg=0.,pyg=0.,pzg=0.;
+
+    //////////////  good tracks definition in TPC  ///////////////
+
+    ofstream out1 ("AliITSTrag.out");
+    Int_t i;
+    for (i=0; i<ngood; i++) out1 << gt[i].ptg << "\n";
+    out1.close();
+
+
+    TVector vec(5);
+    TTree *TR=gAlice->TreeR();
+    Int_t nent=(Int_t)TR->GetEntries();
+    TClonesArray  *recPoints = RecPoints();
+    Int_t numbpoints;
+    Int_t totalpoints=0;
+    Int_t *np = new Int_t[nent];
+    Int_t **vettid = new Int_t* [nent];
+    Int_t mod;
+    for (mod=0; mod<nent; mod++) {
+       vettid[mod]=0;
+       this->ResetRecPoints();
+       //gAlice->TreeR()->GetEvent(mod+1); //first entry in TreeR is empty
+       gAlice->TreeR()->GetEvent(mod); //first entry in TreeR is empty
+       numbpoints = recPoints->GetEntries();
+       totalpoints+=numbpoints;
+       np[mod] = numbpoints;
+       //cout<<" mod = "<<mod<<"   numbpoints = "<<numbpoints<<"\n";getchar();
+       vettid[mod] = new Int_t[numbpoints];
+       Int_t ii;
+       for (ii=0;ii<numbpoints; ii++) *(vettid[mod]+ii)=0;
+    } // end for mod
+
+    AliTPCtrack *track=0;
+    if(min_t < 0) {min_t = 0; max_t = nt-1;}   
+
+/*
+  //////////////////// Definition of vertex end its error ////////////////////
+  ////////////// In the future it will be given by a method //////////////////
+  Double_t Vx=0.;
+  Double_t Vy=0.;
+  Double_t Vz=0.;
+
+  Float_t sigmavx=0.0050;      // 50  microns
+  Float_t sigmavy=0.0050;      // 50  microns
+  Float_t sigmavz=0.010;       // 100 microns
+
+  //Vx+=gRandom->Gaus(0,sigmavx);  Vy+=gRandom->Gaus(0,sigmavy);  Vz+=gRandom->Gaus(0,sigmavz);
+  TVector vertex(3), ervertex(3)
+  vertex(0)=Vx; vertex(1)=Vy; vertex(2)=Vz;
+  ervertex(0)=sigmavx;  ervertex(1)=sigmavy;  ervertex(2)=sigmavz;
+  //////////////////////////////////////////////////////////////////////////
+*/
+
+    TTree tracktree1("TreeT","Tree with ITS tracks");
+    AliITSiotrack *iotrack=0;
+    tracktree1.Branch("ITStracks","AliITSiotrack",&iotrack,32000,0);
+
+    ofstream out ("AliITSTra.out");
+    //ofstream outprova ("AliITSprova.out"); //commentato il 26-4-2001
+
+    Int_t j;       
+    for (j=min_t; j<=max_t; j++) {     
+       track=(AliTPCtrack*)tracks.UncheckedAt(j);
+       Int_t flaglab=0;
+       if (!track) continue;
+       ////// elimination of not good tracks ////////////       
+       Int_t ilab=TMath::Abs(track->GetLabel());
+       Int_t iii;
+       for (iii=0;iii<ngood;iii++) {
+           //cout<<" ilab, gt[iii].lab = "<<ilab<<" "<<gt[iii].lab<<"\n";getchar();
+           if (ilab==gt[iii].lab) { 
+               flaglab=1;
+               ptg=gt[iii].ptg; 
+               pxg=gt[iii].pxg;
+               pyg=gt[iii].pyg;
+               pzg=gt[iii].pzg;        
+               break;
+           } // end if ilab==
+       } // end for iii
+       //cout<<" j flaglab =  " <<j<<" "<<flaglab<<"\n";  getchar();
+       if (!flaglab) continue;  
+       //cout<<" j =  " <<j<<"\n";  getchar();
+       /*
+         ////// old propagation to the end of TPC //////////////       
+         Double_t xk=76.;
+         track->PropagateTo(xk);
+         xk-=0.11;
+         track->PropagateTo(xk,42.7,2.27); //C
+         xk-=2.6;
+         track->PropagateTo(xk,36.2,1.98e-3); //C02
+         xk-=0.051;
+         track->PropagateTo(xk,42.7,2.27); //C 
+         /////////////////////////////////////////////////// 
+       */ 
+       ////// new propagation to the end of TPC //////////////
+       Double_t xk=77.415;
+       track->PropagateTo(xk, 28.94, 1.204e-3);         //Ne
+       xk -=0.01;
+       track->PropagateTo(xk, 44.77, 1.71);     //Tedlar
+       xk -=0.04;
+       track->PropagateTo(xk, 44.86, 1.45);     //Kevlar
+       xk -=2.0;
+       track->PropagateTo(xk, 41.28, 0.029);    //Nomex         
+       xk-=16;
+       track->PropagateTo(xk,36.2,1.98e-3); //C02
+       xk -=0.01;
+       track->PropagateTo(xk, 24.01, 2.7);      //Al    
+       xk -=0.01;
+       track->PropagateTo(xk, 44.77, 1.71);     //Tedlar
+       xk -=0.04;
+       track->PropagateTo(xk, 44.86, 1.45);     //Kevlar
+       xk -=0.5;
+       track->PropagateTo(xk, 41.28, 0.029);    //Nomex
+       ///////////////////////////////////////////////////////////////
+
+       ///////////////////////////////////////////////////////////////
+       AliITStrack trackITS(*track);
+       AliITStrack result(*track);
+       AliITStrack primarytrack(*track); 
+    
+       /////////////////////////////////////////////////////////////////
+       TVector Vgeant(3);
+       Vgeant=result.GetVertex();
+
+       // Definition of Dv and Zv for vertex constraint        
+       Double_t sigmaDv=0.0050;  Double_t sigmaZv=0.010;       
+       //Double_t sigmaDv=0.0015;  Double_t sigmaZv=0.0015;
+       Double_t uniform= gRandom->Uniform();
+       Double_t signdv;
+       if(uniform<=0.5) signdv=-1.;
+       else
+           signdv=1.;
+
+       Double_t Vr=TMath::Sqrt(Vgeant(0)*Vgeant(0)+ Vgeant(1)*Vgeant(1));
+       Double_t Dv=gRandom->Gaus(signdv*Vr,(Float_t)sigmaDv); 
+       Double_t Zv=gRandom->Gaus(Vgeant(2),(Float_t)sigmaZv);
+               
+       //cout<<" Dv e Zv = "<<Dv<<" "<<Zv<<"\n";
+       trackITS.SetDv(Dv);  trackITS.SetZv(Zv);
+       trackITS.SetsigmaDv(sigmaDv); trackITS.SetsigmaZv(sigmaZv); 
+       result.SetDv(Dv);  result.SetZv(Zv);
+       result.SetsigmaDv(sigmaDv); result.SetsigmaZv(sigmaZv);
+       primarytrack.SetDv(Dv);  primarytrack.SetZv(Zv);
+       primarytrack.SetsigmaDv(sigmaDv); primarytrack.SetsigmaZv(sigmaZv);
+
+       ////////////////////////////////////////////////////////////////
+
+       primarytrack.PrimaryTrack(rl);
+       TVector  d2=primarytrack.Getd2();
+       TVector  tgl2=primarytrack.Gettgl2();
+       TVector  dtgl=primarytrack.Getdtgl();
+       trackITS.Setd2(d2); trackITS.Settgl2(tgl2);  trackITS.Setdtgl(dtgl); 
+       result.Setd2(d2); result.Settgl2(tgl2);  result.Setdtgl(dtgl);     
+       /*                       
+          trackITS.SetVertex(vertex); trackITS.SetErrorVertex(ervertex);
+          result.SetVertex(vertex);   result.SetErrorVertex(ervertex);   
+       */
+
+
+       Tracking(trackITS,&result,recPoints,vettid, flagvert,rl);
+       // cout<<" progressive track number = "<<j<<"\r";
+       // cout<<j<<"\r";
+       Int_t NumofCluster=result.GetNumClust();  
+       // cout<<" progressive track number = "<<j<<"\n";    // stampa
+       Long_t labITS=result.GetLabel();
+       //  cout << " ITS track label = " << labITS << "\n";    // stampa
+       int lab=track->GetLabel();                  
+       // cout << " TPC track label = " << lab <<"\n";      // stampa
+
+       //propagation to vertex
+
+       Double_t rbeam=3.;
+
+       result.Propagation(rbeam);
+
+       Double_t C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44;
+       result.GetCElements(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44);
+       Double_t pt=TMath::Abs(result.GetPt());
+       Double_t Dr=result.GetD();
+       Double_t Z=result.GetZ();
+       Double_t tgl=result.GetTgl();
+       Double_t C=result.GetC();
+       Double_t Cy=C/2.;
+       Double_t Dz=Z-(tgl/Cy)*TMath::ASin(result.arga(rbeam));
+       Dz-=Vgeant(2);
+
+       // cout<<" Dr e dz alla fine = "<<Dr<<" "<<Dz<<"\n"; getchar();
+       Double_t phi=result.Getphi();
+       Double_t phivertex = phi - TMath::ASin(result.argA(rbeam));
+       Double_t duepi=2.*TMath::Pi();   
+       if(phivertex>duepi) phivertex-=duepi;
+       if(phivertex<0.) phivertex+=duepi;
+       Double_t Dtot=TMath::Sqrt(Dr*Dr+Dz*Dz);
+
+       /////////////////////////////////////////////////////////////////////
+
+       Int_t idmodule,idpoint;
+       if(NumofCluster >=5)  {            // cinque - sei
+           //if(NumofCluster ==6)  {            // cinque - sei
+
+           AliITSiotrack outtrack;
+
+           iotrack=&outtrack;
+
+           iotrack->SetStatePhi(phi);
+           iotrack->SetStateZ(Z);
+           iotrack->SetStateD(Dr);
+           iotrack->SetStateTgl(tgl);
+           iotrack->SetStateC(C);
+           Double_t radius=result.Getrtrack();
+           iotrack->SetRadius(radius);
+           Int_t charge;
+           if(C>0.) charge=-1;  else charge=1;
+           iotrack->SetCharge(charge);
+
+
+           iotrack->SetCovMatrix(C00,C10,C11,C20,C21,C22,C30,C31,C32,C33,C40,C41,C42,C43,C44);  
+
+           Double_t px=pt*TMath::Cos(phivertex);
+           Double_t py=pt*TMath::Sin(phivertex);
+           Double_t pz=pt*tgl;
+
+           Double_t xtrack=Dr*TMath::Sin(phivertex);
+           Double_t ytrack=Dr*TMath::Cos(phivertex);
+           Double_t ztrack=Dz+Vgeant(2);
+
+
+           iotrack->SetPx(px);
+           iotrack->SetPy(py);
+           iotrack->SetPz(pz);
+           iotrack->SetX(xtrack);
+           iotrack->SetY(ytrack);
+           iotrack->SetZ(ztrack);
+           iotrack->SetLabel(labITS);
+
+           Int_t il;           
+           for(il=0;il<6; il++){
+               iotrack->SetIdPoint(il,result.GetIdPoint(il));
+               iotrack->SetIdModule(il,result.GetIdModule(il));
+           } // end for il
+           tracktree1.Fill();
+
+           //cout<<" labITS = "<<labITS<<"\n";
+           //cout<<" phi z Dr tgl C = "<<phi<<" "<<Z<<" "<<Dr<<" "<<tgl<<" "<<C<<"\n";  getchar();        
+
+           DataOut(kkk) = ptg; kkk++; DataOut(kkk)=labITS; kkk++; DataOut(kkk)=lab; kkk++;             
+
+           for (il=0;il<6;il++) {
+               idpoint=result.GetIdPoint(il);
+               idmodule=result.GetIdModule(il);
+               *(vettid[idmodule]+idpoint)=1; 
+               iotrack->SetIdPoint(il,idpoint);
+               iotrack->SetIdModule(il,idmodule);
+           } // end for il
+
+           //  cout<<"  +++++++++++++  pt e ptg = "<<pt<<" "<<ptg<<"  ++++++++++\n";
+           /* ///  provvisorio il 23-5-2001
+              Double_t pg=TMath::Sqrt(pxg*pxg+pyg*pyg+pzg*pzg);
+              Double_t cosl=TMath::Sqrt(1./(1.+tgl*tgl));
+              Double_t ptot=pt/cosl;
+              cout<<"ptot e pg = "<<ptot<<" "<<pg<<"\n";
+              Double_t difpt= (ptot-pg)/pg*100.;
+           */
+           ///////////////////////////////
+           Double_t difpt= (pt-ptg)/ptg*100.;  //coomentato prov il 23-5-2001
+           DataOut(kkk)=difpt; kkk++;
+           Double_t lambdag=TMath::ATan(pzg/ptg);
+           Double_t   lam=TMath::ATan(tgl);      
+           Double_t diflam = (lam - lambdag)*1000.;
+           DataOut(kkk) = diflam; kkk++;
+           Double_t phig=TMath::ATan2(pyg,pxg);
+           if(phig<0) phig=2.*TMath::Pi()+phig;       
+           Double_t phi=phivertex;
+
+           Double_t difphi = (phi - phig)*1000.;
+           DataOut(kkk)=difphi; kkk++;
+           DataOut(kkk)=Dtot*1.e4; kkk++;
+           DataOut(kkk)=Dr*1.e4; kkk++;
+           DataOut(kkk)=Dz*1.e4; kkk++; 
+           Int_t r;
+           for (r=0; r<9; r++) { out<<DataOut(r)<<" ";}
+           out<<"\n";
+           kkk=0;
+
+       } // end if on NumofCluster
+       //gObjectTable->Print();    // stampa memoria     
+    }  //  end for (int j=min_t; j<=max_t; j++)
+
+    out.close();  
+    //outprova.close();   
+    static Bool_t first=kTRUE;
+    static TFile *tfile;
+
+    if(first) {
+       tfile=new TFile("itstracks.root","RECREATE");
+       //cout<<"I have opened itstracks.root file "<<endl;
+    } // end if first      
+    first=kFALSE;
+    tfile->cd();
+    tfile->ls();
+
+    char hname[30];
+    sprintf(hname,"TreeT%d",evNumber);
 
+    tracktree1.Write(hname);
 
+
+  
+    TTree *fAli=gAlice->TreeK();
+    TFile *fileAli=0;
+           
+    if (fAli) fileAli =fAli->GetCurrentFile();
+    fileAli->cd();
+     
+    ///////////////////////////////////////////////////////////////////
+
+    printf("delete vectors\n");
+    if(np) delete [] np;
+    if(vettid) delete [] vettid;
+}