New version containing the files for tracking

diff --git a/ITS/AliITS.cxx b/ITS/AliITS.cxx
index 2b1c5d98a90e567ee7dacc31d18f3cb45eaa7b34..41b0221473b94f45c02c9946a763c8249e697211 100644 (file)
--- a/ITS/AliITS.cxx
+++ b/ITS/AliITS.cxx
@@ -15,6 +15,9 @@
 
 /*
 $Log$
+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
 
@@ -150,6 +153,12 @@ the AliITS class.
 #include "AliMC.h"
 #include "stdlib.h"
 
+#include "AliITStrack.h"
+#include "AliITSiotrack.h"
+#include "AliITStracking.h" 
+#include "../TPC/AliTPC.h"
+#include "../TPC/AliTPCParam.h"
+
 const Int_t AliITS::fgkNTYPES=3;
 
 ClassImp(AliITS)
@@ -1240,3 +1249,419 @@ void AliITS::Streamer(TBuffer &R__b){
    } // end if
 
 }
+
+//________________________________________________________________
+
+AliITStrack  AliITS::Tracking(AliITStrack &track, AliITStrack *reference,TObjArray *fastpoints, Int_t
+**vettid, Bool_t flagvert ) { 
+                                                                               
+//Origin  A. Badala' and G.S. Pappalardo:  e-mail Angela.Badala@ct.infn.it, Giuseppe.S.Pappalardo@ct.infn.it 
+
+                                                                                   
+  TList *list= new TList();   
+
+  AliITStrack tr(track);
+
+  list->AddLast(&tr);
+  
+  Double_t Pt=(tr).GetPt();
+  //cout << "\n Pt = " << Pt <<"\n";
+
+
+  AliITStracking obj(list, reference, this, fastpoints,TMath::Abs(Pt),vettid, flagvert);
+  list->Delete();
+  delete list;
+
+  Int_t itot=-1;
+  TVector VecTotLabref(18);
+  for(Int_t lay=5; lay>=0; lay--) {
+    TVector VecLabref(3); 
+    VecLabref=(*reference).GetLabTrack(lay);
+    for(Int_t k=0; k<3; k++) {itot++; VecTotLabref(itot)=VecLabref(k);}    
+  }
+  Long_t labref;
+  Int_t freq;  
+  (*reference).Search(VecTotLabref, labref, freq);
+    
+  if(freq < 5) labref=-labref;         
+  (*reference).SetLabel(labref);
+
+  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";
+  
+  struct GoodTrack {
+    Int_t lab,code;
+    Float_t px,py,pz,x,y,z,pxg,pyg,pzg,ptg;
+    Bool_t flag;
+  };
+  
+
+  gAlice->GetEvent(0);
+
+  AliTPC *TPC=(AliTPC*)gAlice->GetDetector("TPC");
+  AliTPCParam *digp = (AliTPCParam*)file->Get(pname);
+  if (digp!=0) TPC->SetParam(digp);
+  
+  GoodTrack gt[7000];
+  Int_t ngood=0;
+  ifstream in("good_tracks");
+
+  cerr<<"Reading good 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==7000) {
+      cerr<<"Too many good tracks !\n";
+      break;
+    }
+  }
+  if (!in.eof()) cerr<<"Read error (good_tracks) !\n";
+  
+  
+// Load tracks
+  TFile *tf=TFile::Open("tpctracks.root");
+  if (!tf->IsOpen()) {cerr<<"Can't open tpctracks.root !\n"; return ;}
+  TObjArray tracks(200000);
+  TTree *tracktree=(TTree*)tf->Get("TreeT");
+  TBranch *tbranch=tracktree->GetBranch("tracks");
+  Int_t nentr=(Int_t)tracktree->GetEntries();
+  for (Int_t k=0; k<nentr; k++) {
+    AliTPCtrack *iotrack=new AliTPCtrack;
+    tbranch->SetAddress(&iotrack);
+    tracktree->GetEvent(k);
+    tracks.AddLast(iotrack);
+  }   
+  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");
+  for (Int_t 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];
+  for (Int_t mod=0; mod<nent; mod++) {
+    vettid[mod]=0;
+    this->ResetRecPoints();
+    gAlice->TreeR()->GetEvent(mod+1); //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];
+    for (Int_t i=0;i<numbpoints; i++) *(vettid[mod]+i)=0;
+  }
+
+  AliTPCtrack *track;
+
+     
+  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;
+  /////////////////////////////////////////////////////////////////////////////////////////////////
+*/      
+
+  //TDirectory *savedir=gDirectory; 
+
+  TTree tracktree1("TreeT","Tree with ITS tracks");
+  AliITSiotrack *iotrack=0;
+  tracktree1.Branch("ITStracks","AliITSiotrack",&iotrack,32000,0);
+
+  ofstream out ("AliITSTra.out");
+         
+  for (int 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());
+    
+    for (Int_t i=0;i<ngood;i++) {
+        //cout<<" ilab, gt[i].lab = "<<ilab<<" "<<gt[i].lab<<"\n"; getchar();
+      if (ilab==gt[i].lab) { 
+       flaglab=1;
+       ptg=gt[i].ptg; 
+       pxg=gt[i].pxg;
+       pyg=gt[i].pyg;
+       pzg=gt[i].pzg;  
+       break;
+      }
+    }
+        //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();
+    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);  
+            
+    cout<<" progressive track number = "<<j<<"\r";
+       // cout<<" progressive track number = "<<j<<"\n";
+    Long_t labITS=result.GetLabel();
+   // cout << " ITS track label = " << labITS << "\n";                     
+    int lab=track->GetLabel();             
+    //cout << " TPC track label = " << lab <<"\n";
+    //result.GetClusters(); getchar();  //to print the cluster matrix
+        
+//propagation to vertex
+       
+    Double_t rbeam=3.;
+    result.Propagation(rbeam);   
+    TMatrix *cov;
+        cov=&result.GetCMatrix();       
+    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 NumofCluster, idmodule,idpoint;
+    NumofCluster=result.GetNumClust();
+    if(NumofCluster >=5)  {     
+
+
+      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(cov);         
+
+      Double_t px=pt*TMath::Cos(phi);
+      Double_t py=pt*TMath::Sin(phi);
+      Double_t pz=pt*tgl;
+               
+      Double_t xtrack=Dr*TMath::Sin(phi);
+      Double_t ytrack=Dr*TMath::Cos(phi);
+      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);
+               
+               for( Int_t il=0;il<6; il++){
+                 iotrack->SetIdPoint(il,result.GetIdPoint(il));
+                 iotrack->SetIdModule(il,result.GetIdModule(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 (Int_t 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);
+      }
+
+      Double_t difpt= (pt-ptg)/ptg*100.;                                       
+      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::ATan(pyg/pxg);
+      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++;
+      for (int 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();  
+  
+  
+  static Bool_t first=kTRUE;
+  static TFile *tfile;
+
+       if(first) {
+           tfile=new TFile("itstracks.root","RECREATE");
+           //cout<<"I have opened itstracks.root file "<<endl;
+       }           
+       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;
+  
+}

diff --git a/ITS/AliITS.h b/ITS/AliITS.h
index 2cc07f124b9b449553708b9543158fad213140f3..bae6046c7a1a63e66adff14dd95165c8ff00d3d8 100644 (file)
--- a/ITS/AliITS.h
+++ b/ITS/AliITS.h
@@ -15,6 +15,7 @@
 
 class TString;
 class TTree;
+class TFile;
 
 class AliITSDetType;
 class AliITSsimulation;
@@ -27,7 +28,7 @@ class AliITSdigit;
 class AliITSRecPoint;
 class AliITSRawCluster;
 class AliITSmodule;
-
+class AliITStrack;
 
 
 
@@ -133,6 +134,13 @@ class AliITS : public AliDetector {
     TTree          *TreeC() {return fTreeC;}
 
 
+    // tracking
+
+    AliITStrack Tracking(AliITStrack &track, AliITStrack *reference, TObjArray *fpoints, Int_t **vettid,
+        Bool_t flagvert );  
+
+    void DoTracking(Int_t evNumber, Int_t min_t, Int_t max_t, TFile *file, Bool_t flagvert);
+
  protected:
 
     static const Int_t fgkNTYPES;  // Number of detector types

diff --git a/ITS/ITSLinkDef.h b/ITS/ITSLinkDef.h
index f5d11d9f018481a0b012ba9237635237ca293087..cb15dba83a74d54fba1caaae7528f78e26dd0f20 100644 (file)
--- a/ITS/ITSLinkDef.h
+++ b/ITS/ITSLinkDef.h
@@ -83,6 +83,10 @@
 #pragma link C++ class  AliITSdcsSSD+;
 #pragma link C++ class  AliITSclusterSSD+;
 #pragma link C++ class  AliITSpackageSSD+;
+// Classes used for Tracking
+#pragma link C++ class  AliITStrack+;
+#pragma link C++ class  AliITStracking+;
+#pragma link C++ class  AliITSiotrack+;
 // New used for Alignment studdies
 //#pragma link C++ class  AliITSAlignmentTrack-;
 //#pragma link C++ class  AliITSAlignmentModule-;
@@ -102,8 +106,6 @@
 //#pragma link C function GetStringDialog;
 //#pragma link C function GetIntegerDialog;
 //#pragma link C function GetFloatDialog;
-// New classes used for Tracking
-//#pragma link C++ class  AliITStrack+;
 // This class will always be for ITS only
 //#pragma link C++ class  AliITSvtest-;
 

diff --git a/ITS/Makefile b/ITS/Makefile
index 7d38c0eaea0ebe2a70a3268539b387d4c8bfb603..f4a5f9309ebc61a9eb17d2e3c736ddcb5a033565 100644 (file)
--- a/ITS/Makefile
+++ b/ITS/Makefile
@@ -34,7 +34,8 @@ SRCS          = AliITS.cxx AliITSv1.cxx AliITSv3.cxx AliITSv5.cxx \
                AliITSHuffman.cxx AliITSClusterFinderSSD.cxx \
                AliITSclusterSSD.cxx AliITSpackageSSD.cxx \
                AliITSdictSSD.cxx AliITSgeomSPD300.cxx AliITSgeomSPD425.cxx \
-               AliITSstatistics.cxx AliITSstatistics2.cxx 
+               AliITSstatistics.cxx AliITSstatistics2.cxx \
+                AliITStrack.cxx AliITStracking.cxx AliITSiotrack.cxx           
 #              AliITSAlignmentTrack.cxx AliITSAlignmentModule.cxx \
 #              AliITSvtest.cxx AliITStrack.cxx