]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ITS/AliITSTrackerV1.cxx
Coding conventions
[u/mrichter/AliRoot.git] / ITS / AliITSTrackerV1.cxx
index 310787bb8f94b738ca680f29c84c844267f2271c..30b44ac2910fc1cca3bc5672dbd5d7594193aad8 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
  
-/*
-$Log$
-Revision 1.10.2.1  2001/10/24 07:26:04  hristov
-All the changes from the head are merged with the release
-
-Revision 1.14  2001/10/24 07:19:57  hristov
-Some pointer correctly initialised in one of the constructors
-
-Revision 1.13  2001/10/21 19:17:12  hristov
-Several pointers were set to zero in the default constructors to avoid memory management problems
-
-Revision 1.12  2001/10/19 21:32:35  nilsen
-Minor changes to remove compliation warning on gcc 2.92.2 compiler, and
-cleanded up a little bit of code.
+/* $Id$ */
 
-*/
-//     The purpose of this class is to permorm the ITS tracking. The 
+// The purpose of this class is to permorm the ITS tracking. The 
 // constructor has the task to inizialize some private members. The method 
 // DoTracking is written to be called by a macro. It gets the event number,
 // the minimum and maximum order number of TPC tracks that are to be tracked
@@ -43,104 +29,127 @@ cleanded up a little bit of code.
 // imposition respectively. The authors thank Mariana Bondila to have help 
 // them to resolve some problems.  July-2000
 
-#include <fstream.h>
+#include <Riostream.h>
 #include <TMath.h>
 #include <TBranch.h>
 #include <TVector.h>
 #include <TFile.h>
+#include <TRandom.h>
 #include <TTree.h>
-#include <TStopwatch.h>
 
 #include "TParticle.h"
 #include "AliRun.h"
-#include "AliITS.h"
-#include "AliITSsegmentationSSD.h"
+#include "AliITSDetTypeRec.h"
 #include "AliITSgeomSPD.h"
 #include "AliITSgeomSDD.h"
 #include "AliITSgeomSSD.h"
 #include "AliITSgeom.h"
 #include "AliITSRecPoint.h"
-#include "stdlib.h"
 #include "AliKalmanTrack.h" 
 #include "AliMagF.h"
 #include "AliITSTrackV1.h"
 #include "AliITSIOTrack.h"
 #include "AliITSRad.h"   
-#include "../TPC/AliTPCtracker.h"
+#include "AliTPCtracker.h"
 #include "AliITSTrackerV1.h"
-#include "AliITSVertex.h"
+#include "AliITSPid.h"
+#include "AliMC.h"
 
 ClassImp(AliITSTrackerV1)
  //______________________________________________________________________
 AliITSTrackerV1::AliITSTrackerV1() {
   //Default constructor
-  fITS = 0;
+  //fITS = 0;
   fresult = 0;
+  fPtref=0.;
+  fChi2max=0.;
+  //fepsphi=0.;
+  //fepsz=0.;
   frecPoints = 0;
   fvettid = 0;
+  fflagvert=0;
   frl = 0;
+  Int_t ia;
+  for(ia=0; ia<6; ia++) {
+  fNlad[ia]=0;
+  fNdet[ia]=0;
+  fAvrad[ia]=0.;
+  fDetx[ia]=0.;
+  fDetz[ia]=0.; 
+  } // end for ia  
   fzmin = 0;
   fzmax = 0;
   fphimin = 0;
   fphimax = 0;
   fphidet = 0;
+  fNRecPoints=0;
+  fRecCylR=0;
+  fRecCylPhi=0;
+  fRecCylZ=0;
+  fFieldFactor=0;
+  fDetTypeRec = 0;
 }
 //______________________________________________________________________
-AliITSTrackerV1::AliITSTrackerV1(AliITS* IITTSS, Bool_t flag) {
+AliITSTrackerV1::AliITSTrackerV1(Int_t evnumber, Bool_t flag) {
     //Origin   A. Badala' and G.S. Pappalardo:  
     // e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
     // Class constructor. It does some initializations.
 
   //PH Initialisation taken from the default constructor
-    fITS      = IITTSS;
-    fresult = 0;
-    frecPoints = 0;
-    fvettid = 0;
-    frl = 0;
-    fzmin = 0;
-    fzmax = 0;
-    fphimin = 0;
-    fphimax = 0;
-    fphidet = 0;
-
-    fPtref    = 0.;
-    fChi2max  =0.;     
-    fflagvert =flag;
-    Int_t imax = 200,jmax = 450;
-    frl       = new AliITSRad(imax,jmax);
-
-    //////////  gets information on geometry /////////////////////////////
-    AliITSgeom *g1 = ((AliITS*)gAlice->GetDetector("ITS"))->GetITSgeom();
-    // Why not AliITS *g1 = fITS->GetITSgeom(); // ?? BSN
-    Int_t ll=1, dd=1;
-    TVector det(9);
-
-    //cout<<" nlad ed ndet \n";
-    Int_t ia;
-    for(ia=0; ia<6; ia++) {
-       fNlad[ia]=g1->GetNladders(ia+1);
-       fNdet[ia]=g1->GetNdetectors(ia+1);
-       //cout<<fNlad[i]<<" "<<fNdet[i]<<"\n"; 
-    } // end for ia
-    //getchar();
-
-    //cout<<" raggio medio = ";
-    Int_t ib;
-    for(ib=0; ib<6; ib++) {  
-       g1->GetCenterThetaPhi(ib+1,ll,dd,det);
-       Double_t r1=TMath::Sqrt(det(0)*det(0)+det(1)*det(1));
-       g1->GetCenterThetaPhi(ib+1,ll,dd+1,det);
-       Double_t r2=TMath::Sqrt(det(0)*det(0)+det(1)*det(1));
-       fAvrad[ib]=(r1+r2)/2.;
-       //cout<<fAvrad[ib]<<" ";
-    } // end for ib
+  //fITS      = IITTSS;
+  fresult = 0;
+  fPtref    = 0.;
+  fChi2max  =0.; 
+  frecPoints = 0;               
+  fvettid = 0;
+  fflagvert = flag;     
+  frl = 0;
+  fzmin = 0;
+  fzmax = 0;
+  fphimin = 0;
+  fphimax = 0;
+  fphidet = 0;
+  
+  Int_t imax = 200,jmax = 450;
+  frl       = new AliITSRad(imax,jmax);
+  
+  //////////  gets information on geometry /////////////////////////////
+  AliRunLoader* rl = AliRunLoader::Open("galice.root");
+  AliITSLoader* loader = (AliITSLoader*)runLoader->GetLoader("ITSLoader");
+  if (!loader) {
+    Error("AliITSTrackerV1", "ITS loader not found");
+    return;
+  }
+  fDetTypeRec = new AliITSDetTypeRec(loader);
+  AliITSgeom* g1 = loader->GetITSgeom();
+
+  Int_t ll=1, dd=1;
+  TVector det(9);
+  
+  Int_t ia;
+  for(ia=0; ia<6; ia++) {
+    fNlad[ia]=g1->GetNladders(ia+1);
+    fNdet[ia]=g1->GetNdetectors(ia+1);
+    //cout<<fNlad[i]<<" "<<fNdet[i]<<"\n"; 
+  } // end for ia
+  
+    //cout<<" mean radius = ";
+  Int_t ib;
+  for(ib=0; ib<6; ib++) {  
+    g1->GetCenterThetaPhi(ib+1,ll,dd,det);
+    Double_t r1=TMath::Sqrt(det(0)*det(0)+det(1)*det(1));
+    g1->GetCenterThetaPhi(ib+1,ll,dd+1,det);
+    Double_t r2=TMath::Sqrt(det(0)*det(0)+det(1)*det(1));
+    fAvrad[ib]=(r1+r2)/2.;
+    //cout<<fAvrad[ib]<<" ";
+  } // end for ib
     //cout<<"\n"; getchar();
-
-    fDetx[0] = ((AliITSgeomSPD*)(g1->GetShape(1, ll, dd)))->GetDx();
-    fDetz[0] = ((AliITSgeomSPD*)(g1->GetShape(1, ll, dd)))->GetDz();
-
-    fDetx[1] = ((AliITSgeomSPD*)(g1->GetShape(2, ll, dd)))->GetDx();
-    fDetz[1] = ((AliITSgeomSPD*)(g1->GetShape(2, ll, dd)))->GetDz();
+  
+  fDetx[0] = ((AliITSgeomSPD*)(g1->GetShape(1, ll, dd)))->GetDx();
+  fDetz[0] = ((AliITSgeomSPD*)(g1->GetShape(1, ll, dd)))->GetDz();
+  
+  fDetx[1] = ((AliITSgeomSPD*)(g1->GetShape(2, ll, dd)))->GetDx();
+  fDetz[1] = ((AliITSgeomSPD*)(g1->GetShape(2, ll, dd)))->GetDz();
 
     fDetx[2] = ((AliITSgeomSDD*)(g1->GetShape(3, ll, dd)))->GetDx();
     fDetz[2] = ((AliITSgeomSDD*)(g1->GetShape(3, ll, dd)))->GetDz();
@@ -157,6 +166,7 @@ AliITSTrackerV1::AliITSTrackerV1(AliITS* IITTSS, Bool_t flag) {
     //for(Int_t la=0; la<6; la++) cout<<"    "<<fDetx[la]<<"     "<<
     //                                 fDetz[la]<<endl;
     //getchar();
+        
     // allocate memory and define matrices fzmin, fzmax, fphimin and fphimax //
     Double_t epsz=1.2;
     Double_t epszdrift=0.05;
@@ -198,7 +208,7 @@ AliITSTrackerV1::AliITSTrackerV1(AliITS* IITTSS, Bool_t flag) {
        fphidet[im1] = new Double_t[im2max];
     } // end for im1
 
-    Float_t global[3],local[3];
+       Double_t global[3],local[3];
     Double_t pigre=TMath::Pi();
     Double_t xmin,ymin,xmax,ymax;
 
@@ -227,20 +237,64 @@ AliITSTrackerV1::AliITSTrackerV1(AliITS* IITTSS, Bool_t flag) {
            if(fphimax[im1][im2]<0.) fphimax[im1][im2]+=2.*pigre;
        } // end for im2
     } // end for im1
+//////////////////////////////////////////////////////////////////////////////////////////////////////////
+/////////////// allocate memory and define vector fNRecPoints and matrices fRecCylR, fRecCylPhi, fRecCylZ /////////////
+       gAlice->GetEvent(evnumber);
+  Int_t numOfModules = g1->GetIndexMax();
+  fRecCylR = new Double_t *[numOfModules];
+  fRecCylPhi = new Double_t *[numOfModules]; 
+  fRecCylZ = new Double_t *[numOfModules];
+  AliITSRecPoint *recp;
+  fNRecPoints = new Int_t[numOfModules];
+   
+                for(Int_t module=0; module<numOfModules; module++) {                           
+                 fDetTypeRec->ResetRecPoints();                     
+                 gAlice->TreeR()->GetEvent(module);              
+                 frecPoints=fDetTypeRec->RecPoints();
+                 Int_t nRecPoints=fNRecPoints[module]=frecPoints->GetEntries();
+                 fRecCylR[module] = new Double_t[nRecPoints];
+                 fRecCylPhi[module] = new Double_t[nRecPoints];
+                 fRecCylZ[module] = new  Double_t[nRecPoints];           
+                 Int_t ind;
+                 for(ind=0; ind<fNRecPoints[module]; ind++) {    
+                   recp=(AliITSRecPoint*)frecPoints->UncheckedAt(ind);                                          
+                       // Float_t global[3], local[3];
+                        Double_t global[3], local[3];
+              local[0]=recp->GetX();
+              local[1]=0.;
+              local[2]= recp->GetZ();                                  
+                        g1->LtoG(module,local,global);
+                                                                                        
+                        Double_t r = TMath::Sqrt(global[0]*global[0]+global[1]*global[1]);                     // r hit
+                        Double_t phi = TMath::ATan2(global[1],global[0]); if(phi<0.) phi+=2.*TMath::Pi();      // phi hit                      
+                Double_t z = global[2];                                                                // z hit
+                                                                                                                                                                                                                                    
+                        fRecCylR[module][ind]=r;
+                        fRecCylPhi[module][ind]=phi;
+                        fRecCylZ[module][ind]=z;                        
+                 }             
+               }        
+        //}  
+  //}
+////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+        
 
     ////////// gets magnetic field factor //////////////////////////////
 
     AliMagF * fieldPointer = gAlice->Field();
-    fFieldFactor = (Double_t)fieldPointer->Factor();
-    //cout<< " field factor = "<<fFieldFactor<<"\n"; getchar();
+   // fFieldFactor = (Double_t)fieldPointer->Factor();
+    fFieldFactor =(Double_t)fieldPointer-> SolenoidField()/10/.2;
+   // cout<< " field factor = "<<fFieldFactor<<"\n"; getchar();
+    delete rl;
+    
 }
 //______________________________________________________________________
-AliITSTrackerV1::AliITSTrackerV1(const AliITSTrackerV1 &cobj) {
+AliITSTrackerV1::AliITSTrackerV1(const AliITSTrackerV1 &cobj) : TObject(cobj) {
     // Origin  A. Badala' and G.S. Pappalardo:
     // e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it
     // copy constructor
-
-    *fITS = *cobj.fITS;
+        
+  //    *fITS = *cobj.fITS;
     *fresult = *cobj.fresult;
     fPtref = cobj.fPtref;
     fChi2max = cobj.fChi2max;    
@@ -291,23 +345,85 @@ AliITSTrackerV1::AliITSTrackerV1(const AliITSTrackerV1 &cobj) {
            fphidet[im1][im2]=cobj.fphidet[im1][im2];  
        } // end for im2
     } // end for im2
+
+
+    AliRunLoader* rl = AliRunLoader::Open("galice.root");
+    rl->CdGAFile();
+    AliITSgeom* g1 = (AliITSgeom*)gDirectory->Get("AliITSgeom");
+
+   Int_t numOfModules = g1->GetIndexMax();
+       /*
+  fRecCylR = new Float_t *[numOfModules];
+  fRecCylPhi = new Float_t *[numOfModules]; 
+  fRecCylZ = new Float_t *[numOfModules];
+  */
+  fRecCylR = new Double_t *[numOfModules];
+  fRecCylPhi = new Double_t *[numOfModules]; 
+  fRecCylZ = new Double_t *[numOfModules];  
+  fNRecPoints = new Int_t[numOfModules];       
+               for(Int_t module=0; module<numOfModules; module++) {            
+                 Int_t nRecPoints=fNRecPoints[module]=cobj.fNRecPoints[module];
+                 /*
+                 fRecCylR[module] = new Float_t[nRecPoints];
+                 fRecCylPhi[module] = new Float_t[nRecPoints];
+                 fRecCylZ[module] = new Float_t[nRecPoints];
+                 */
+                 fRecCylR[module] = new Double_t[nRecPoints];
+                 fRecCylPhi[module] = new Double_t[nRecPoints];
+                 fRecCylZ[module] = new Double_t[nRecPoints];            
+                 Int_t ind;    
+                 for(ind=0; ind<nRecPoints; ind++) {       
+                        fRecCylR[module][ind]=cobj.fRecCylR[module][ind];
+                        fRecCylPhi[module][ind]=cobj.fRecCylPhi[module][ind];
+                        fRecCylZ[module][ind]=cobj.fRecCylZ[module][ind];                       
+                 }             
+               }        
+               delete rl;
+}
+/*
+//______________________________________________________________________
+void AliITSTrackerV1::DelMatrix(Int_t numOfModules) {
+  // cleanup of some data members
+  for(Int_t mod=0; mod<numOfModules; mod++) {
+    delete fRecCylR[mod];
+        delete fRecCylPhi[mod];
+        delete fRecCylZ[mod];
+  }
+    delete fRecCylR;
+        delete fRecCylPhi;
+        delete fRecCylZ;
 }
+*/
 //______________________________________________________________________
 AliITSTrackerV1::~AliITSTrackerV1(){
     // Origin  A. Badala' and G.S. Pappalardo:
     // e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it  
-    // class destructor
-
-    //delete fTimerKalman;
-    //delete fTimerIntersection;
+    // class destructor         
+  delete frl;
+  delete fNRecPoints;
+  for(Int_t i=0; i<6; i++) {
+    delete fzmin[i];
+        delete fzmax[i];
+        delete fphimin[i];
+        delete fphimax[i];
+        delete fphidet[i];
+  }
+
+  delete fzmin;
+  delete fzmax;
+  delete fphimin;
+  delete fphimax;
+  delete fphidet;
+        
 }
 //______________________________________________________________________
-AliITSTrackerV1 &AliITSTrackerV1::operator=(AliITSTrackerV1 obj) {
+AliITSTrackerV1 &AliITSTrackerV1::operator=(const AliITSTrackerV1 &obj) {
     // Origin  A. Badala' and G.S. Pappalardo:
     // e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it  
     // assignement operator
 
-    *fITS = *obj.fITS;
+  //    *fITS = *obj.fITS;
     *fresult = *obj.fresult;
     fPtref = obj.fPtref;
     fChi2max = obj.fChi2max;      
@@ -360,6 +476,30 @@ AliITSTrackerV1 &AliITSTrackerV1::operator=(AliITSTrackerV1 obj) {
        } // end for im2
     } // end for im1
 
+    AliRunLoader* rl = AliRunLoader::Open("galice.root");
+    rl->CdGAFile();
+    AliITSgeom* g1 = (AliITSgeom*)gDirectory->Get("AliITSgeom");
+    // AliITSgeom *g1 = fITS->GetITSgeom();  
+   Int_t numOfModules = g1->GetIndexMax();
+  fRecCylR = new Double_t *[numOfModules];
+  fRecCylPhi = new Double_t *[numOfModules]; 
+  fRecCylZ = new Double_t *[numOfModules];  
+  fNRecPoints = new Int_t[numOfModules];  
+         for(Int_t module=0; module<numOfModules; module++) {            
+                 Int_t nRecPoints=fNRecPoints[module]=obj.fNRecPoints[module];
+                 fRecCylR[module] = new Double_t[nRecPoints];
+                 fRecCylPhi[module] = new Double_t[nRecPoints];
+                 fRecCylZ[module] = new Double_t[nRecPoints];            
+                 Int_t ind;
+                 for(ind=0; ind<nRecPoints; ind++) {     
+                        fRecCylR[module][ind]=obj.fRecCylR[module][ind];
+                        fRecCylPhi[module][ind]=obj.fRecCylPhi[module][ind];
+                        fRecCylZ[module][ind]=obj.fRecCylZ[module][ind];                        
+                 }             
+               }        
+         
+        
+         delete rl;
     return *this;
 }
 //______________________________________________________________________
@@ -374,26 +514,36 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
     // The method can be called by a macro. It preforms the tracking for
     // all good TPC tracks
 
-    printf("begin DoTracking - file %p\n",file);
+    printf("begin DoTracking - file %p\n",(void*)file);
 
     gAlice->GetEvent(evNumber);  //modificato per gestire hbt
  
-    AliKalmanTrack *kkprov;
-    kkprov->SetConvConst(100/0.299792458/0.2/fFieldFactor);
+    AliKalmanTrack::SetFieldMap(gAlice->Field());
+   // cout<<" field = "<<gAlice->Field()->SolenoidField()<<endl;
+
 
     TFile *cf=TFile::Open("AliTPCclusters.root");  
-    AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60");
+    AliTPCParam *digp= (AliTPCParam*)cf->Get("75x40_100x60_150x60");
     if (!digp) { cerr<<"TPC parameters have not been found !\n"; getchar();}
 
     cf->cd();
-    AliTPCtracker *tracker = new AliTPCtracker(digp,evNumber);  
+    AliRunLoader* rl = AliRunLoader::Open("galice.root");
+    AliITSLoader* itsl = (AliITSLoader*)rl->GetLoader("ITSLoader");
+    TString foldname(itsl->GetEventFolder()->GetName());
+    //TString foldname(fITS->GetLoader()->GetEventFolder()->GetName());
+    delete rl;   
+    printf("This method is not converted to the NewIO !\n");  //I.B.
+    return; //I.B.
+    AliTPCtracker *tracker = new AliTPCtracker(digp);            //I.B.
+    //PH    tracker->SetEventNumber(evNumber);                           //I.B.
 
     // Load clusters
-    tracker->LoadInnerSectors();
-    tracker->LoadOuterSectors();
+    printf("This method is not converted to the NewIO !\n");  //I.B.
+    return; //I.B.
+    tracker->LoadClusters(0);       //I.B.
 
     // Load tracks
-    TFile *tf=TFile::Open("AliTPCtracksSorted.root");  //modificato per hbt
+    TFile *tf=TFile::Open("AliTPCtracksSorted.root");  
     if (!tf->IsOpen()) {
        cerr<<"Can't open AliTPCtracksSorted.root !\n";
        return ;
@@ -408,7 +558,8 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
     Int_t nentr=(Int_t)tracktree->GetEntries();
     Int_t kk;
 
-    AliTPCtrack *ioTrackTPC=0;    
+    AliITSRecPoint *recp;    // oggi
+       AliTPCtrack *ioTrackTPC=0;    
     for (kk=0; kk<nentr; kk++) {
        ioTrackTPC=new AliTPCtrack; 
        tbranch->SetAddress(&ioTrackTPC);
@@ -425,7 +576,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
     TVector vec(5);
     TTree *tr=gAlice->TreeR();
     Int_t nent=(Int_t)tr->GetEntries();  
-    frecPoints = fITS->RecPoints();
+    frecPoints = fDetTypeRec->RecPoints();
 
     Int_t numbpoints;
     Int_t totalpoints=0;
@@ -435,7 +586,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
 
     for (mod=0; mod<nent; mod++) {
        fvettid[mod]=0;
-       fITS->ResetRecPoints();  
+       fDetTypeRec->ResetRecPoints();  
        gAlice->TreeR()->GetEvent(mod); 
        numbpoints = frecPoints->GetEntries();
        totalpoints+=numbpoints;
@@ -454,9 +605,11 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
 
     TTree tracktree1("TreeT","Tree with ITS tracks");
     AliITSIOTrack *ioTrack=0;
+       AliITSPid *pid=new AliITSPid(1000);  // oggi
+
     tracktree1.Branch("ITStracks","AliITSIOTrack",&ioTrack,32000,0);
   
-   TDatabasePDG * db = new TDatabasePDG;   
+    TDatabasePDG * db = new TDatabasePDG;   
   
     Int_t j;       
     for (j=minTr; j<=maxTr; j++) {     
@@ -466,41 +619,56 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
        ///   mass definition ////////////////////////
     Double_t mass=0.13956995;
         Int_t pcode=211;  // a pion by default
+        
         if(realmass) {
-    Int_t TPClabel=TMath::Abs( track->GetLabel() );   
-        TParticle *p = (TParticle*)gAlice->Particle(TPClabel);
-         pcode=p->GetPdgCode();
-        // Int_t mothercode=p->GetFirstMother();
-        //if(mothercode>0 ) numofsecondaries++; else numofprimaries++;
+       if(TMath::Abs(pcode)<20443) mass=db->GetParticle(pcode)->Mass();
         }
-        //if(!pcode) pcode=211;         
-        if(TMath::Abs(pcode)<20443) mass=db->GetParticle(pcode)->Mass();       
+        else {
+                mass = track->GetMass();
+//              cout << "Mass = " << mass << endl;
+        }
+        
        
        
-       ///////////////////////////////////////////////
-
-       //////   propagation to the end of TPC //////////////
-       Double_t xk=77.415;
-       track->PropagateTo(xk, 28.94, 1.204e-3,mass);    //Ne    
-       xk -=0.01;
-       track->PropagateTo(xk, 44.77, 1.71,mass);        //Tedlar
-       xk -=0.04;
-       track->PropagateTo(xk, 44.86, 1.45,mass);        //kevlar
-       xk -=2.0;
-       track->PropagateTo(xk, 41.28, 0.029,mass);       //Nomex         
-       xk-=16;
-       track->PropagateTo(xk,36.2,1.98e-3,mass); //C02
-       xk -=0.01;
-       track->PropagateTo(xk, 24.01, 2.7,mass);         //Al    
-       xk -=0.01;
-       track->PropagateTo(xk, 44.77, 1.71,mass);        //Tedlar
-       xk -=0.04;
-       track->PropagateTo(xk, 44.86, 1.45,mass);        //kevlar
-       xk -=0.5;
-       track->PropagateTo(xk, 41.28, 0.029,mass);       //Nomex
-    ////////////////////////////////////////////////////////////////////
-
-       AliITSTrackV1 trackITS(*track);     
+        //   new propagation to the end of TPC
+    Double_t xk=80.;
+   // track->PropagateTo(xk,0.,0.); //Ne if it's still there   //attenzione funziona solo se modifica in TPC
+       // Double_t xk=77.415;   
+        track->PropagateTo(xk, 28.94, 1.204e-3);
+    xk-=0.005;
+    track->PropagateTo(xk, 44.77,1.71); //Tedlar        
+    xk-=0.02;
+    track->PropagateTo(xk, 44.86, 1.45);   //Kevlar
+    xk-=2.0;
+    track->PropagateTo(xk, 41.28, 0.029);//Nomex
+    xk-=0.02;
+    track->PropagateTo(xk, 44.86, 1.45);   //Kevlar
+    xk-=0.005;
+    track->PropagateTo(xk, 44.77, 1.71); //Tedlar
+
+    xk=61.;
+   // track->PropagateTo(xk,0.,0.); //C02
+        track->PropagateTo(xk,36.2,1.98e-3); //C02        //attenzione funziona solo se modifica in TPC
+
+    xk -=0.005;
+    track->PropagateTo(xk, 24.01, 2.7);    //Al    
+    xk -=0.005;
+    track->PropagateTo(xk, 44.77, 1.71);  //Tedlar
+    xk -=0.02;
+    track->PropagateTo(xk, 44.86, 1.45);    //Kevlar
+    xk -=0.5;
+    track->PropagateTo(xk, 41.28, 0.029);  //Nomex    
+    xk -=0.02;
+    track->PropagateTo(xk, 44.86, 1.45);    //Kevlar
+    xk -=0.005;
+    track->PropagateTo(xk, 44.77, 1.71);  //Tedlar
+    xk -=0.005;
+    track->PropagateTo(xk, 24.01, 2.7);    //Al 
+    
+       ////////////////////////////////////////////////////////////////////////////////////////////////////////        
+       //AliITSTrackV1 trackITS(*track);
+       AliITSTrackV1 trackITS(*track, fFieldFactor);
+       //cout<<" fFieldFactor = "<<fFieldFactor<<"\n";      
        trackITS.PutMass(mass);   //new to add mass to track
        if(fresult){ delete fresult; fresult=0;}         
        fresult = new AliITSTrackV1(trackITS);   
@@ -533,7 +701,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
        primaryTrack.SetZv(zv);
        primaryTrack.SetsigmaDv(sigmaDv);
        primaryTrack.SetsigmaZv(sigmaZv);
-       primaryTrack.PrimaryTrack(frl);
+       primaryTrack.PrimaryTrack();
        TVector  d2=primaryTrack.Getd2();
        TVector  tgl2=primaryTrack.Gettgl2();
        TVector  dtgl=primaryTrack.Getdtgl();
@@ -550,6 +718,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
        list->AddLast(&trackITS);
   
        fPtref=TMath::Abs( (trackITS).GetPt() );
+       //cout<<" fPtref = " <<fPtref<<"\n";
        if(fPtref>1.0) fChi2max=40.;         
        if(fPtref<=1.0) fChi2max=20.;
        if(fPtref<0.4 ) fChi2max=100.;
@@ -573,12 +742,12 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
            for(k=0; k<3; k++){  
                Int_t lpp=(Int_t)vecLabRef(k);
                if(lpp>=0) {
-                   TParticle *p=(TParticle*) gAlice->Particle(lpp);
+                   TParticle *p=(TParticle*) gAlice->GetMCApp()->Particle(lpp);
                    Int_t pcode=p->GetPdgCode();
                    if(pcode==11) vecLabRef(k)=p->GetFirstMother();
                } // end if
                itot++; vecTotLabRef(itot)=vecLabRef(k);
-               if(vecLabRef(k)==0. && clustZ == 0.) vecTotLabRef(itot) =-3.;
+               if(vecLabRef(k)==0. && clustZ == -1.) vecTotLabRef(itot) =-3.;
            } // end for k
        } // end for lay
        Long_t labref;
@@ -600,7 +769,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
        //propagation to vertex
 
        Double_t rbeam=3.;
-    if((*fresult).DoNotCross(rbeam)) continue;  //no intersection with beampipe        
+       if((*fresult).DoNotCross(rbeam)) continue;  //no intersection with beampipe     
        (*fresult).Propagation(rbeam);
        Double_t c00,c10,c11,c20,c21,c22,c30,c31,c32,c33,c40,c41,c42,c43,c44;
        (*fresult).GetCElements(c00,
@@ -639,6 +808,8 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
            Int_t charge;
            if(c>0.) charge=-1;  else charge=1;
            ioTrack->SetCharge(charge);
+               Double_t trackmass=(*fresult).GetMass();         // oggi
+               ioTrack->SetMass(trackmass);                     // oggi
            ioTrack->SetCovMatrix(c00,
                                  c10,c11,
                                  c20,c21,c22,
@@ -660,11 +831,15 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
            ioTrack->SetTPCLabel(lab);
                 ioTrack->SetDz(dz);
            Int_t il;           
-           for(il=0;il<6; il++){
+           /*
+               for(il=0;il<6; il++){
                ioTrack->SetIdPoint(il,(*fresult).GetIdPoint(il));
                ioTrack->SetIdModule(il,(*fresult).GetIdModule(il));
            } // end for il
-           tracktree1.Fill();
+           */
+           //tracktree1.Fill();
+               Float_t q[4]={-1.,-1.,-1.,-1.};
+           Float_t  globaldedx=0.;        
            for (il=0;il<6;il++) {
                idpoint=(*fresult).GetIdPoint(il);
                idmodule=(*fresult).GetIdModule(il);
@@ -672,7 +847,53 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
                                                               
                ioTrack->SetIdPoint(il,idpoint);
                ioTrack->SetIdModule(il,idmodule);
-           } // end for il         
+               ////  for q definition
+               if(il>1){
+                 if(idmodule>0.){                      
+                   fDetTypeRec->ResetRecPoints();
+                   gAlice->TreeR()->GetEvent(idmodule);
+                   recp=(AliITSRecPoint*)frecPoints->UncheckedAt(idpoint);
+                   q[il-2]=recp->GetQ()*(*fresult).Getfcor(il-2);
+                 }
+               }                               
+           } // end for il      
+               q[0]/=280.; q[1]/=280.;
+           q[2]/=38.; q[3]/=38.;
+
+ // cout<<" q prima = "<<q[0]<<" "<<q[1]<<" "<<q[2]<<" "<<q[3]<<"\n"; getchar(); 
+     
+     Int_t swap;
+  do{
+    swap=0;   
+    for (il=0; il<3; il++) {
+      if (q[il]<=q[il+1]) continue;
+      Float_t tmp=q[il];
+      q[il]=q[il+1]; q[il+1]=tmp;
+      swap++;
+    }
+  } while(swap); 
+
+  // cout<<" q dopo = "<<q[0]<<" "<<q[1]<<" "<<q[2]<<" "<<q[3]<<"\n"; getchar();
+      
+    if(q[0]<0.) {
+      q[0]=q[1];
+      q[1]=q[2];
+      q[2]=q[3];
+      q[3]=-1.;            
+    } 
+  // cout<<" q  dopo if = "<<q[0]<<" "<<q[1]<<" "<<q[2]<<" "<<q[3]<<"\n"; getchar(); 
+     
+    globaldedx=(q[0]+q[1])/2.;
+    
+   // if(q[3]> 0.) globaldedx=(q[0]+q[1]+q[2]+q[3])/4.;
+   //      else    globaldedx=(q[0]+q[1]+q[2])/3.; 
+   
+    ioTrack->SetdEdx(globaldedx);
+    ioTrack->SetPid(pid->GetPcode(ioTrack));
+           
+           tracktree1.Fill();         
        } // end if on numOfCluster
        //gObjectTable->Print();    // stampa memoria     
     }  //  end for (int j=minTr; j<=maxTr; j++)
@@ -688,6 +909,7 @@ void AliITSTrackerV1::DoTracking(Int_t evNumber,Int_t minTr,Int_t maxTr,
     tfile->ls();
     char hname[30];
     sprintf(hname,"TreeT%d",evNumber);
+       cout << "Number of saved ITS tracks " << tracktree1.GetEntries() << endl;
     tracktree1.Write(hname);
   
     TTree *fAli=gAlice->TreeK();
@@ -729,7 +951,10 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
     sigmazl[5]=0.6889; 
     ///////////////////////////////////////////////////////////
     Int_t index; 
-    AliITSgeom *g1 = fITS->GetITSgeom();
+    AliRunLoader* rl = AliRunLoader::Open("galice.root");
+    rl->CdGAFile();
+    AliITSgeom* g1 = (AliITSgeom*)gDirectory->Get("AliITSgeom");
+
     AliITSRecPoint *recp;        
     for(index =0; index<trackITSlist->GetSize(); index++) {
        AliITSTrackV1 *trackITS = (AliITSTrackV1 *) trackITSlist->At(index);
@@ -768,7 +993,7 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
 
            if(numClustNow<numClustRef && chi2Now>fresult->GetChi2()) continue;
            //cout<<" chi2Now =  "<<chi2Now<<"\n";   
-           //  commentato il 30-7-2001   
+             
            chi2Now/=numClustNow;
            if(fPtref > 1.0 && chi2Now > 30.) continue; 
            if((fPtref >= 0.6 && fPtref<=1.0) && chi2Now > 40.) continue;
@@ -785,13 +1010,7 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
        Int_t ladp, ladm, detp,detm,ladinters,detinters;        
        Int_t layerfin=layerInit-1;
        // cout<<"Prima di intersection \n";
-       //if(!fTimerIntersection) fTimerIntersection = new TStopwatch();
-       // timer
-       //fTimerIntersection->Continue();
-       // timer 
        Int_t  outinters=Intersection(*trackITS,layerfin,ladinters,detinters);
-       //fTimerIntersection->Stop();
-       // timer
        // cout<<" outinters = "<<outinters<<"\n";
        //  cout<<" Layer ladder detector intersection ="
        //      <<layerfin<<" "<<ladinters<<" "<<detinters<<"\n";
@@ -799,6 +1018,8 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
        //       << " "<<(*trackITS)(1)<<"\n"; getchar();
        if(outinters==-1) continue;
        Int_t flaghit=0;
+       (*trackITS).SetLayer(layerfin);  // oggi
+       (*trackITS).Setfcor();           // oggi
        if(outinters==0){
            TVector toucLad(9), toucDet(9);      
            Int_t lycur=layerfin;
@@ -830,7 +1051,7 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
            */
            Float_t epsphi=5.0, epsz=5.0;                  
            if(fPtref<0.2) {epsphi=3.; epsz=3.;}     
-           // nuova definizione idetot e toucLad e toucDet to be
+           // new definition of idetot e toucLad e toucDet to be
            // transformed in a method
            // these values could be modified
            Float_t pigre=TMath::Pi();
@@ -863,11 +1084,12 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
            Float_t distz = 0.0;
            Float_t phicm, phicp, distphim, distphip;
            phicm=phinters;
-           if(phinters>fphimax[layerfin-1][ladm]) phicm=phinters-2*pigre;
-           distphim=TMath::Abs(phicm-fphimax[layerfin-1][ladm]);
+                if(phinters>fphimax[layerfin-1][ladm-1]) phicm=phinters-2*pigre;  //corretto il 20-11-2001
+                distphim=TMath::Abs(phicm-fphimax[layerfin-1][ladm-1]);  //corretto il 20-11-2001
            phicp=phinters;
-           if(phinters>fphimin[layerfin-1][ladp]) phicp=phinters-2.*pigre;
-           distphip=TMath::Abs(phicp-fphimin[layerfin-1][ladp]);
+                //cout<<" fNlad[layerfin-1] e ladp = "<<fNlad[layerfin-1]<<" "<<ladp<<endl;
+                if(phinters>fphimin[layerfin-1][ladp-1]) phicp=phinters-2.*pigre;   //corretto il 20-11-2001
+                distphip=TMath::Abs(phicp-fphimin[layerfin-1][ladp-1]);      //corretto il 20-11-2001
            Int_t flagzmin=0;
            Int_t flagzmax=0;
            idetot=1;
@@ -936,74 +1158,39 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
                ///////////////////////////////////////////////////////
                /*** Rec points sorted by module *****/
                /**************************************/
-               Int_t indexmod;       //mod ott
-               // AliITSRecPoint *recp;
+               Int_t indexmod;       
                indexmod = g1->GetModuleIndex(lycur,(Int_t)toucLad(iriv),
-                                             (Int_t)toucDet(iriv));
-               //mod ott 
-               fITS->ResetRecPoints();   
+                                             (Int_t)toucDet(iriv)); 
+               fDetTypeRec->ResetRecPoints();   
                gAlice->TreeR()->GetEvent(indexmod); 
                Int_t npoints=frecPoints->GetEntries();
-               /* mod ott
-                  Int_t *indlist=new Int_t[npoints+1];
-                  Int_t counter=0;
-                  Int_t ind;
-                  for (ind=0; ind<=npoints; ind++) {
-                       indlist[ind]=-1;
-                       if (*(fvettid[index]+ind)==0) {
-                            indlist[counter]=ind;
-                            counter++;
-                       } // end if
-                  } // end for ind
-                  ind=-1;
-                  for(;;) { 
-                        ind++;
-                        if(indlist[ind] < 0) recp=0;
-                        else 
-                             recp = (AliITSRecPoint*)frecPoints->
-                                                    UncheckedAt(indlist[ind]);
-                       if((!recp)  )  break;
-                  } // end for ;;
-               */
-               ///////////////////////// new //////////////////////////       
+       
                Int_t indnew;
                for(indnew=0; indnew<npoints; indnew++){
                    if (*(fvettid[indexmod]+indnew)==0)
                        recp =(AliITSRecPoint*)frecPoints->UncheckedAt(indnew);
                    else
                        continue;
-                   ////////////////////////////////////////////////////
                    TVector cluster(3),vecclust(9);
-                   vecclust(6)=vecclust(7)=vecclust(8)=-1.;
-                   Double_t sigma[2];               
-                   // set veclust in global
-                   Float_t global[3], local[3];
-                   local[0]=recp->GetX();
-                   local[1]=0.;
-                   local[2]= recp->GetZ();
-                   Int_t play = lycur;
-                   Int_t plad = TMath::Nint(toucLad(iriv));   
-                   Int_t pdet = TMath::Nint(toucDet(iriv));            
-                   g1->LtoG(play,plad,pdet,local,global);
-                   vecclust(0)=global[0];
-                   vecclust(1)=global[1];
-                   vecclust(2)=global[2];
-
+                   //vecclust(6)=vecclust(7)=vecclust(8)=-1.;
+                   Double_t sigma[2];
+  // now vecclust is with cylindrical cohordinates
+              vecclust(0)=(Float_t)fRecCylR[indexmod][indnew];     
+              vecclust(1)=(Float_t)fRecCylPhi[indexmod][indnew];
+              vecclust(2)=(Float_t)fRecCylZ[indexmod][indnew];                          
                    vecclust(3) = (Float_t)recp->fTracks[0]; 
-                   //vecclust(4) = (Float_t)indlist[ind];
                    vecclust(4) = (Float_t)indnew;                       
-                   vecclust(5) = (Float_t)indexmod;    //mod ott
+                   vecclust(5) = (Float_t)indexmod;    
                    vecclust(6) = (Float_t)recp->fTracks[0];
                    vecclust(7) = (Float_t)recp->fTracks[1];
                    vecclust(8) = (Float_t)recp->fTracks[2];
                    sigma[0] = (Double_t)  recp->GetSigmaX2();     
                    sigma[1] = (Double_t) recp->GetSigmaZ2();
-                   //now we are in r,phi,z in global
-                   cluster(0) = TMath::Sqrt(vecclust(0)*vecclust(0)+
-                                            vecclust(1)*vecclust(1));//r hit
-                   cluster(1) = TMath::ATan2(vecclust(1),vecclust(0));
-                   if(cluster(1)<0.) cluster(1)+=2.*TMath::Pi();
-                   cluster(2) = vecclust(2);                   // z hit
+                        
+                        cluster(0)=fRecCylR[indexmod][indnew];
+          cluster(1)=fRecCylPhi[indexmod][indnew];
+                        cluster(2)=fRecCylZ[indexmod][indnew];
+                        
                    // cout<<" layer = "<<play<<"\n";
                    // cout<<" cluster prima = "<<vecclust(0)<<" "
                    //     <<vecclust(1)<<" "
@@ -1030,7 +1217,7 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
                        continue;
                    // cout<<" supero sigmaphi \n";      
                    AliITSTrackV1 *newTrack = new AliITSTrackV1((*trackITS));
-                   (*newTrack).SetLayer((*trackITS).GetLayer()-1);
+                   //(*newTrack).SetLayer((*trackITS).GetLayer()-1);
                    if (TMath::Abs(rTrack-cluster(0))/rTrack>1e-6) 
                        (*newTrack).Correct(Double_t(cluster(0)));      
                    //cout<<" cluster(2) e(*newTrack).GetZ()="<<cluster(2)<<" "
@@ -1049,34 +1236,30 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
                    // cout<<" chi2pred = "<<chi2pred<<"\n";
                    // if(chi2pred>fChi2max) continue; //aggiunto il 30-7-2001
                    if(iriv == 0) flaghit=1;
-                   (*newTrack).AddMS(frl);  // add the multiple scattering 
+                   (*newTrack).AddMS();  // add the multiple scattering 
                                             //matrix to the covariance matrix 
-                   (*newTrack).AddEL(frl,1.,0);
+                   (*newTrack).AddEL(1.,0);
 
-                   //if(!fTimerKalman) fTimerKalman = new TStopwatch();//timer
-                   //fTimerKalman->Continue();                        // timer
                    if(fflagvert){
                        KalmanFilterVert(newTrack,cluster,sigmanew);
                        //KalmanFilterVert(newTrack,cluster,sigmanew,chi2pred);
                    }else{
                        KalmanFilter(newTrack,cluster,sigmanew);
                    } // end if
-                   //fTimerKalman->Stop();                         // timer
                    (*newTrack).PutCluster(layernew, vecclust);
                    newTrack->AddClustInTrack();
                    listoftrack.AddLast(newTrack);
                }   // end for indnew
-               // delete [] indlist;  //mod ott
            }  // end of for on detectors (iriv)
        }//end if(outinters==0)
 
        if(flaghit==0 || outinters==-2) {
            AliITSTrackV1 *newTrack = new AliITSTrackV1(*trackITS);
            (*newTrack).Setfnoclust();           
-           (*newTrack).SetLayer((*trackITS).GetLayer()-1); 
-           (*newTrack).AddMS(frl);  // add the multiple scattering matrix
+           //(*newTrack).SetLayer((*trackITS).GetLayer()-1); 
+           (*newTrack).AddMS();  // add the multiple scattering matrix
                                     // to the covariance matrix  
-           (*newTrack).AddEL(frl,1.,0);
+           (*newTrack).AddEL(1.,0);
            listoftrack.AddLast(newTrack);
        } // end if
 
@@ -1087,6 +1270,8 @@ void AliITSTrackerV1::RecursiveTracking(TList *trackITSlist) {
     } // end of for on tracks (index)
 
     //gObjectTable->Print();   // stampa memoria
+     delete rl;
+    
 }
 //______________________________________________________________________
 Int_t AliITSTrackerV1::Intersection(AliITSTrackV1 &track,Int_t layer,
@@ -1534,3 +1719,4 @@ void AliITSTrackerV1::KalmanFilterVert(AliITSTrackV1 *newTrack,
     newTrack->SetChi2(newTrack->GetChi2()+chi2);   
     //   newTrack->SetChi2(newTrack->GetChi2()+chi2pred);
 }
+