macro to be updated for newIO
authormasera <masera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 1 Aug 2003 15:43:20 +0000 (15:43 +0000)
committermasera <masera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 1 Aug 2003 15:43:20 +0000 (15:43 +0000)
38 files changed:
ITS/oldmacros/AliITSDeleteOldFiles.sh [new file with mode: 0755]
ITS/oldmacros/AliITSHits2SDigitsDubna.C [new file with mode: 0644]
ITS/oldmacros/AliITSImpactParametersPP.C [new file with mode: 0644]
ITS/oldmacros/AliITSMerge.C [new file with mode: 0644]
ITS/oldmacros/AliITSNeuralCompleteEval.C [new file with mode: 0644]
ITS/oldmacros/AliITSNeuralTracking.C [new file with mode: 0644]
ITS/oldmacros/AliITSOccupancy.C [new file with mode: 0644]
ITS/oldmacros/AliITSPlotGeom.C [new file with mode: 0644]
ITS/oldmacros/AliITSPlotTracksV1.C [new file with mode: 0644]
ITS/oldmacros/AliITSPrintDigits.C [new file with mode: 0644]
ITS/oldmacros/AliITSPrintSDigits.C [new file with mode: 0644]
ITS/oldmacros/AliITSPropagateBackV2.C [new file with mode: 0644]
ITS/oldmacros/AliITSReadPlotData.C [new file with mode: 0644]
ITS/oldmacros/AliITSSDigits2DigitsDubna.C [new file with mode: 0644]
ITS/oldmacros/AliITSSavePIDV2.C [new file with mode: 0644]
ITS/oldmacros/AliITSScanPIDV1.C [new file with mode: 0644]
ITS/oldmacros/AliITSScanPIDV2.C [new file with mode: 0644]
ITS/oldmacros/AliITSStoreFindableTracks.C [new file with mode: 0644]
ITS/oldmacros/AliITSStoreFindableTracksCompiled.C [new file with mode: 0644]
ITS/oldmacros/AliITSTrackingV1.C [new file with mode: 0644]
ITS/oldmacros/AliITSTrackingV2.C [new file with mode: 0644]
ITS/oldmacros/AliITSVertexerTracksTest.C [new file with mode: 0644]
ITS/oldmacros/AliITSVertexerTracksTest2.C [new file with mode: 0644]
ITS/oldmacros/AliITSplotSPDClusters.C [new file with mode: 0644]
ITS/oldmacros/AliITSrefitV2.C [new file with mode: 0644]
ITS/oldmacros/AliITSsd2d.C [new file with mode: 0644]
ITS/oldmacros/AliITStestAll0.sh [new file with mode: 0755]
ITS/oldmacros/AliITStestPIDV2.sh [new file with mode: 0755]
ITS/oldmacros/AliITStestV1.sh [new file with mode: 0755]
ITS/oldmacros/ConfigSPD02.C [new file with mode: 0644]
ITS/oldmacros/Config_bg.C [new file with mode: 0644]
ITS/oldmacros/Config_muon.C [new file with mode: 0644]
ITS/oldmacros/DrawS01.C [new file with mode: 0644]
ITS/oldmacros/DrawSPD.C [new file with mode: 0644]
ITS/oldmacros/MergeV1.C [new file with mode: 0644]
ITS/oldmacros/ViewITSSPD.C [new file with mode: 0644]
ITS/oldmacros/ViewS01.C [new file with mode: 0644]
ITS/oldmacros/ViewSPD.C [new file with mode: 0644]

diff --git a/ITS/oldmacros/AliITSDeleteOldFiles.sh b/ITS/oldmacros/AliITSDeleteOldFiles.sh
new file mode 100755 (executable)
index 0000000..ee37a57
--- /dev/null
@@ -0,0 +1,10 @@
+rm galice.root
+rm AliTPCclusters.root
+rm AliTPCtracks.root
+rm AliTPCtracksSorted.root
+rm AliITSclustersV2.root
+rm AliITStracksV2.root
+rm AliITStracksV2Pid.root
+rm itstracks.root
+rm good_tracks_tpc
+rm itsgood_tracks
diff --git a/ITS/oldmacros/AliITSHits2SDigitsDubna.C b/ITS/oldmacros/AliITSHits2SDigitsDubna.C
new file mode 100644 (file)
index 0000000..27c8935
--- /dev/null
@@ -0,0 +1,135 @@
+TFile* AccessFile(TString inFile="galice.root", TString acctype="R");
+void writeAR(TFile * fin, TFile *fou);
+
+Int_t AliITSHits2SDigitsDubna(Int_t evNumber1=0,Int_t evNumber2=0,
+                             TString inFile ="galice.root", 
+                             TString outFile="galiceS.root"){
+
+  // Dynamically link some shared libs
+  if (gClassTable->GetID("AliRun") < 0) {
+       gROOT->LoadMacro("loadlibs.C");
+       loadlibs();
+  } // end if
+
+  // Connect the Root Galice file containing Geometry, Kine and Hits
+
+  TFile *file;
+  if(outFile.Data() == inFile.Data()){
+    file = AccessFile(inFile,"U");
+  }
+  else {
+    file = AccessFile(inFile);
+  }
+  
+  TFile *file2 = 0;  // possible output file for TreeS
+
+  if(!(outFile.Data() == inFile.Data())){
+    // open output file and create TreeS on it
+    file2 = gAlice->InitTreeFile("S",outFile);
+  }
+
+  AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");      
+  if (!ITS) {
+       cerr<<"AliITSHits2DigitsDefault.C : AliITS object not found on file"
+           << endl;
+       return 3;
+  }  // end if !ITS
+  if(!(ITS->GetITSgeom())){
+       cerr << " AliITSgeom not found. Can't digitize with out it." << endl;
+       return 4;
+  } // end if
+
+  // For old files, must change SPD noise.
+  AliITSresponseSPDdubna *resp0 = new AliITSresponseSPDdubna();
+  if(ITS->DetType(0)->GetResponseModel() !=0){
+      delete ((AliITSresponse*)ITS->DetType(0)->GetResponseModel());
+      ITS->DetType(0)->ResponseModel(0);
+  } // end if
+  ITS->DetType(0)->ResponseModel(resp0);
+  AliITSsegmentationSPD *seg0 = (AliITSsegmentationSPD*)ITS->DetType(0)->GetSegmentationModel();
+  AliITSsimulationSPDdubna *sim0 = new AliITSsimulationSPDdubna(seg0,resp0);
+  if(ITS->DetType(0)->GetSimulationModel() !=0){
+      delete ((AliITSsimulation*)ITS->DetType(0)->GetSimulationModel());
+      ITS->DetType(0)->SimulationModel(0);
+  } // end if
+  ITS->DetType(0)->SimulationModel(sim0);
+
+  TStopwatch timer;
+  timer.Start();
+  for(Int_t nevent = evNumber1; nevent <= evNumber2; nevent++){
+    gAlice->GetEvent(nevent);
+    if(!gAlice->TreeS() && file2 == 0){ 
+      cout << "Having to create the SDigits Tree." << endl;
+      gAlice->MakeTree("S");
+    } // end if !gAlice->TreeS()
+    if(file2)gAlice->MakeTree("S",file2);
+    //    make branch
+    ITS->MakeBranch("S");
+    ITS->SetTreeAddress();
+    cout<<"Making ITS SDigits for event "<<nevent<<endl;
+    TStopwatch timer;
+    Long_t size0 = file->GetSize();
+    ITS->Hits2SDigits();
+  }
+  timer.Stop();
+  timer.Print();
+
+  // write the AliRun object to the output file
+  if(file2)writeAR(file,file2);
+
+  delete gAlice;   gAlice=0;
+  file->Close();
+}
+
+//-------------------------------------------------------------------
+TFile * AccessFile(TString FileName, TString acctype){
+
+  // Function used to open the input file and fetch the AliRun object
+
+  TFile *retfil = 0;
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(FileName);
+  if (file) {file->Close(); delete file; file = 0;}
+  if(acctype.Contains("U")){
+    file = new TFile(FileName,"update");
+  }
+  if(acctype.Contains("N") && !file){
+    file = new TFile(FileName,"recreate");
+  }
+  if(!file) file = new TFile(FileName);   // default readonly
+  if (!file->IsOpen()) {
+       cerr<<"Can't open "<<FileName<<" !" << endl;
+       return retfil;
+  } 
+
+  // Get AliRun object from file or return if not on file
+  if (gAlice) {delete gAlice; gAlice = 0;}
+  gAlice = (AliRun*)file->Get("gAlice");
+  if (!gAlice) {
+       cerr << "AliRun object not found on file"<< endl;
+       return retfil;
+  } 
+  return file;
+}
+
+//-------------------------------------------------------------------
+void writeAR(TFile * fin, TFile *fou) {
+  TDirectory *current = gDirectory;
+  TTree *Te;
+  TTree *TeNew;
+  AliHeader *alhe = new AliHeader();
+  Te = (TTree*)fin->Get("TE");
+  Te->SetBranchAddress("Header",&alhe);
+  Te->SetBranchStatus("*",1);
+  fou->cd();
+  TeNew = Te->CloneTree();
+  TeNew->Write(0,TObject::kOverwrite);
+  gAlice->Write(0,TObject::kOverwrite);
+  current->cd();
+  delete alhe;
+  cout<<"AliRun object written to file\n";
+}
+
+
+
+
+
diff --git a/ITS/oldmacros/AliITSImpactParametersPP.C b/ITS/oldmacros/AliITSImpactParametersPP.C
new file mode 100644 (file)
index 0000000..f05c54b
--- /dev/null
@@ -0,0 +1,191 @@
+/****************************************************************************
+ *                                                                          *
+ * This macro determines the impact parameter of each track in pp event     *
+ * with respect to the primaty vertex position reconstructed using all      *
+ * the other tracks in the event. This procedure is necessary to avoid      *
+ * biases in the impact parameter estimate.                                 *
+ *                                                                          *
+ * Output:  TNtuple (ntd0) on file ImpactParameters.root                    *
+ *          There is one entry per track. The elements of the ntuple are:   *
+ *          1) Event: Event number                                          *
+ *          2) nTraks: Number of tracks for this event                      *
+ *          3) pt: Transverse momentum for the current track                *
+ *          4) d0rphi: transverse impact parameter of the current track     *
+ *          5) vx:     x coordinate of the vertex                           *
+ *          6) vy:     y coordinate of the vertex                           *
+ * origin: A.Dainese    andrea.dainese@pd.infn.it                           * 
+ ****************************************************************************/ 
+#if !defined(__CINT__) || defined(__MAKECINT__)
+//-- --- standard headers------------- 
+#include <Riostream.h>
+//--------Root headers ---------------
+#include <TSystem.h>
+#include <TFile.h>
+#include <TNtuple.h>
+#include <TString.h>
+#include <TStopwatch.h>
+#include <TObject.h>
+#include <TVector3.h>
+#include <TTree.h>
+#include <TParticle.h>
+#include <TArray.h>
+//----- AliRoot headers ---------------
+#include "alles.h"
+#include "AliRun.h"
+#include "AliKalmanTrack.h"
+#include "AliITStrackV2.h"
+#include "AliHeader.h"
+#include "AliGenEventHeader.h"
+#include "AliV0vertex.h"
+#include "AliV0vertexer.h"
+#include "AliITSVertex.h"
+#include "AliITSVertexer.h"
+#include "AliITSVertexerTracks.h"
+#endif
+//-------------------------------------
+
+// field (T)
+const Double_t kBz = 0.4;
+
+// this function ckecks file existence
+Bool_t GetInputFile();
+void   MoveTracks(TTree& itsTree,TObjArray& array);
+
+
+void   AliITSImpactParametersPP(Int_t evFirst=0,Int_t evLast=0) {
+
+  const Char_t *name="AliITSImpactParametersPP";
+  cerr<<'\n'<<name<<"...\n";
+  gBenchmark->Start(name);  
+
+  AliKalmanTrack::SetConvConst(100/0.299792458/kBz);
+
+
+  // check existence of input file
+  if(!GetInputFile()) { cerr<<"No tracks file found"<<endl; return; }
+
+  TString outName("ImpactParameters.root");
+
+  // primary vertex
+  Double_t v1[3] = {0.,0.,0,};
+  
+  Int_t    ev,itsEntries,i;
+  Int_t    nTotEv=0;
+  Double_t pt,d0rphi;
+  Char_t   trksName[100];
+  AliITStrackV2 *itstrack = 0;
+
+  // output ntuple
+  TNtuple *ntd0 = new TNtuple("ntd0","ntd0","Event:nTrks:pt:d0rphi:vx:vy");
+
+  // create the AliITSVertexerTracks object
+  AliITSVertexerTracks *vertexer1 = new AliITSVertexerTracks;
+  vertexer1->SetMinTracks(3);
+  vertexer1->SetDebug(0);
+  Int_t skipped[2];
+
+  // Open file with ITS tracks
+  TFile* itstrks = TFile::Open("AliITStracksV2.root");
+
+
+  // loop on events in file
+  for(ev=evFirst; ev<=evLast; ev++) {
+    printf(" --- Processing event  %d ---\n",ev);
+    sprintf(trksName,"TreeT_ITS_%d",ev);
+
+    // tracks from ITS
+    TTree *itsTree  =(TTree*)itstrks->Get(trksName);
+    if(!itsTree) continue;
+    itsEntries = (Int_t)itsTree->GetEntries();
+    printf("+++\n+++ Number of tracks in ITS: %d\n+++\n\n",itsEntries);
+
+    TObjArray tArray(itsEntries);
+    MoveTracks(*itsTree,tArray);
+    
+    // count the total number of events
+    nTotEv++;
+
+    // loop on tracks
+    for(i=0; i<itsEntries; i++) {
+      itstrack = (AliITStrackV2*)tArray.At(i);
+
+      // pt
+      pt = 1./TMath::Abs(itstrack->Get1Pt());
+
+      // primary vertex from other tracks in event
+      Bool_t goodVtx = kFALSE;
+      vertexer1->SetVtxStart(0.,0.);
+      skipped[0] = i;
+      skipped[1] = -1;
+      vertexer1->SetSkipTracks(1,skipped);
+      AliITSVertex *vertex1 = 
+        (AliITSVertex*)vertexer1->VertexOnTheFly(*itsTree); 
+      if(vertex1->GetNContributors()>0) goodVtx = kTRUE; 
+      vertex1->GetXYZ(v1);
+      //vertex1->PrintStatus(); 
+      // impact parameter w.r.t. v1      (in microns)
+      d0rphi =  10000.*itstrack->GetD(v1[0],v1[1]);
+      delete vertex1;
+
+      // fill ntuple
+      if (goodVtx) ntd0->Fill(ev,itsEntries,pt,d0rphi,v1[0],v1[1]);
+
+      itstrack = 0;
+    }   // loop on tracks
+
+    delete itsTree;
+
+  }    // loop on events in file
+
+
+  printf("\n+++\n+++ Total number of events: %d\n+++\n",nTotEv);
+
+  delete vertexer1;
+
+  itstrks->Close();
+
+  // store ntuple in a file
+  TFile* outroot = new TFile(outName.Data(),"recreate");
+  ntd0->Write();
+  outroot->Close();
+  delete outroot;
+
+
+  gBenchmark->Stop(name);
+  gBenchmark->Show(name);
+
+  return;
+}
+//___________________________________________________________________________
+Bool_t   GetInputFile() {
+  TString itsName("AliITStracksV2.root");
+
+  if(gSystem->AccessPathName(itsName.Data(),kFileExists)) return kFALSE;
+
+  return kTRUE;
+}
+//___________________________________________________________________________
+void MoveTracks(TTree& itsTree,TObjArray& array) {
+//
+// this function creates two TObjArrays with tracks
+//
+  Int_t entr = (Int_t)itsTree.GetEntries();
+
+  // trasfer tracks from tree to array
+  for(Int_t i=0; i<entr; i++) {
+
+    AliITStrackV2 *t = new AliITStrackV2; 
+    itsTree.SetBranchAddress("tracks",&t);
+
+    itsTree.GetEvent(i);
+
+    array.AddLast(t);
+
+  }
+
+  return;
+}
+
+
+
diff --git a/ITS/oldmacros/AliITSMerge.C b/ITS/oldmacros/AliITSMerge.C
new file mode 100644 (file)
index 0000000..f2cb899
--- /dev/null
@@ -0,0 +1,192 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+
+#include "iostream.h"
+#include "TDatetime.h"
+#include "STEER/AliRun.h"
+#include "STEER/AliRunDigitizer.h"
+#include "ITS/AliITSDigitizer.h"
+#include "ITS/AliITS.h"
+#include "ITS/AliITSDetType.h"
+#include "ITS/AliITSresponseSDD.h"
+#include "TStopwatch.h"
+
+#endif
+
+TFile* AccessFile(TString inFile="galice.root", TString acctype="R");
+void writeAR(TFile * fin, TFile *fou);
+Int_t ChangeITSDefaults(TFile *hitfile,AliITS *ITS,TString opt="");
+//#define DEBUG
+Int_t AliITSMerge(TString digFile="galiceMD.root", 
+           TString sdigFileSig="galiceS_sig.root", 
+           TString sdigFileBg="",TString opt=""){
+    // Standeard ITS SDigits to Digits, with posible merging.
+
+    // Dynamically link some shared libs
+    if (gClassTable->GetID("AliRun") < 0) {
+       gROOT->LoadMacro("loadlibs.C");
+       loadlibs();
+    } // end if
+    // Connect the Root Galice file containing Geometry, Kine and Hits
+    if (gAlice) {delete gAlice; gAlice = 0;}
+
+    TFile *sdigfilesig = 0;   // pointer to signal input file.
+    TFile *sdigfilebg  = 0;   // pointer to background input file.
+    TFile *digfile = 0;  // possible output file for TreeD
+
+    // Setup to copy gAlice and the event tree to the output file.
+
+    if(digFile.CompareTo(sdigFileSig)==0){//write output to same file as input.
+       sdigfilesig = AccessFile(sdigFileSig,"U");//input file open for update.
+    }else{ // different output file then input file.
+       sdigfilesig = AccessFile(sdigFileSig,"R");//input file open read only
+       digfile = new TFile(digFile,"NEW");
+    } // end if digFile == hitFile.
+
+    AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");      
+    if (!ITS) {
+       cerr<<"AliITSHits2DigitsDefault.C : AliITS object not found on file"
+           << endl;
+       return 3;
+    }  // end if !ITS
+
+    ChangeITSDefaults(sdigfilesig,ITS,opt);
+    // write the AliRun object to the output file if different from input file.
+    if(digfile){
+       writeAR(sdigfilesig,digfile);
+       digfile->Close(); // Manager will open in update mode.
+       digfile = 0;
+    } // end if digfile
+    sdigfilesig->Close();
+    sdigfilesig = 0;
+//    delete gAlice; gAlice=0;   // there is a problem with deleting gAlice????
+
+    AliRunDigitizer *manager;
+    if(sdigFileBg.CompareTo("")==0) { // do not merge.
+       manager = new AliRunDigitizer(1,1);
+    }else{
+       manager = new AliRunDigitizer(2,1);
+       manager->SetInputStream(0,sdigFileSig.Data());
+       manager->SetInputStream(1,sdigFileBg.Data());
+    } // end if
+    if (digFile.CompareTo(sdigFileSig) !=0) {
+       manager->SetOutputFile(digFile);
+    } // end if
+//    manager->SetCopyTreesFromInput(0);
+    AliITSDigitizer *dITS = new AliITSDigitizer(manager);
+    if(opt.Contains("ROI")){
+       cout << "Region of Interest selected" << endl;
+       dITS->SetByRegionOfInterestFlag(1);
+    }else{
+       cout << "Digizing everthing" << endl;
+       dITS->SetByRegionOfInterestFlag(0);
+    } // end if
+
+    TStopwatch timer;
+    timer.Start();
+    manager->Exec("all");
+    timer.Stop(); 
+    timer.Print();
+
+    if(digfile!=0){
+       cout << digFile << " size =" << digfile->GetSize() << endl;
+    }else if(sdigfilesig!=0){
+       cout << sdigFileSig << " size =" << sdigfilesig->GetSize() << endl;
+    } // end if sdigfile!=0
+
+
+    if(digfile) digfile->Close();
+    if(sdigfilesig) sdigfilesig->Close();
+    if(sdigfilebg)  sdigfilebg->Close();
+
+    // There is a problem deleting gAlice. It is related to the destructor
+    // and the fTreeE in AliRun. So for now it is not deleted.
+//    delete gAlice; // digfile is closed by deleting gAlice if != hitfile.
+//    gAlice = 0;
+    delete manager;
+}
+//______________________________________________________________________
+TFile * AccessFile(TString FileName, TString acctype){
+    // Function used to open the input file and fetch the AliRun object
+
+    TFile *retfil = 0;
+    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(FileName);
+    if(file) {
+       file->Close();
+       delete file;
+       file = 0;
+    } // end if file
+    if(acctype.Contains("U")){
+       file = new TFile(FileName,"UPDATE");
+    } // end if open for update
+    if(acctype.Contains("N") && !file){
+       file = new TFile(FileName,"RECREATE");
+    } // end if open a new file
+    if(!file) file = new TFile(FileName,"READ");   // default readonly
+    if (!file->IsOpen()) {
+       cerr << "Can't open " << FileName << " !" << endl;
+       return retfil;
+    } // end if error opeing file
+
+    // Get AliRun object from file or return if not on file
+    if (gAlice) {delete gAlice; gAlice = 0;}
+    gAlice = (AliRun*)file->Get("gAlice");
+    if (!gAlice) {
+       cerr << "AliRun object not found on file "<< FileName << "!" << endl;
+       file->Close();  // close file and return error.
+       return retfil;
+    } // end if !gAlice
+    return file;
+}
+//______________________________________________________________________
+void writeAR(TFile * fin, TFile *fou) {
+    TDirectory *current = gDirectory;
+    TTree *TeOld;
+    TTree *TeNew;
+    AliHeader *alhe = new AliHeader();
+    TeOld = (TTree*)fin->Get("TE");
+    TeOld->SetBranchAddress("Header",&alhe);
+    TeOld->SetBranchStatus("*",1);
+    fou->cd();
+    TeNew = TeOld->CloneTree();
+    TeNew->Write(0,TObject::kOverwrite);
+    gAlice->Write(0,TObject::kOverwrite);
+    current->cd();
+    delete alhe;
+#ifdef DEBUG
+    cout << "AliRun object written to file" << endl;
+#endif
+}
+//______________________________________________________________________
+Int_t ChangeITSDefaults(TFile *hitfile,AliITS *ITS,TString opt){
+
+    TDatime *ct0 = new TDatime(2002,04,26,00,00,00);
+    TDatime ct = hitfile->GetCreationDate();
+
+    if(ct0->GetDate()>ct.GetDate()){
+       // For old files, must change SDD noise.
+       AliITSresponseSDD *resp1 = (AliITSresponseSDD*)ITS->DetType(1)->
+           GetResponseModel();
+       resp1 = new AliITSresponseSDD();
+       ITS->SetResponseModel(1,resp1);
+       cout << "Changed response class for SDD:" << endl;
+       resp1->Print();
+    } // end if
+
+    if(opt.Contains("Dubna")){
+       AliITSresponseSPDdubna *resp0 = new AliITSresponseSPDdubna();
+       if(ITS->DetType(0)->GetResponseModel() !=0){
+           delete ((AliITSresponse*)ITS->DetType(0)->GetResponseModel());
+           ITS->DetType(0)->ResponseModel(0);
+       } // end if
+       ITS->DetType(0)->ResponseModel(resp0);
+       AliITSsegmentationSPD *seg0 = (AliITSsegmentationSPD*)ITS->
+           DetType(0)->GetSegmentationModel();
+       AliITSsimulationSPDdubna *sim0 = new AliITSsimulationSPDdubna(seg0,
+                                                                     resp0);
+       if(ITS->DetType(0)->GetSimulationModel() !=0){
+           delete ((AliITSsimulation*)ITS->DetType(0)->GetSimulationModel());
+           ITS->DetType(0)->SimulationModel(0);
+       } // end if
+       ITS->DetType(0)->SimulationModel(sim0);
+    } // end if Dubna
+}
diff --git a/ITS/oldmacros/AliITSNeuralCompleteEval.C b/ITS/oldmacros/AliITSNeuralCompleteEval.C
new file mode 100644 (file)
index 0000000..f9e4640
--- /dev/null
@@ -0,0 +1,233 @@
+// Neural performance evaluation
+//
+// before using this macro, you should have make run the macro
+// 'ITSstoreFindableTracks.C' in order to have a file named following
+// the rule to append '_fnd' to the name 
+// (e.g.   galice.root --> galice_fnd.root)
+// and have in it a tree with some tracks statistics useful for evaluation
+// of performance with the same event without loosing too much time...
+
+void AliITSNeuralCompleteEval
+(Int_t nsecs, const char *filename, Bool_t low = kFALSE, 
+ Bool_t draw = kTRUE, const char *save = 0)
+{
+       Int_t N, M;
+       Float_t *edge = 0;
+       if (!low) {
+               N = 7;
+               edge = new Float_t[8];
+               edge[0] = 0.0;
+               edge[1] = 0.5;
+               edge[2] = 1.0;
+               edge[3] = 1.5;
+               edge[4] = 2.0;
+               edge[5] = 3.0;
+               edge[6] = 4.0;
+               edge[7] = 5.0;
+       }
+       else {
+               N = 5;
+               edge = new Float_t[6];
+               edge[0] = 0.0;
+               edge[1] = 0.2;
+               edge[2] = 0.4;
+               edge[3] = 0.6;
+               edge[4] = 0.8;
+               edge[5] = 1.0;
+       }
+       
+       Float_t find[2] = {0.0, 0.0}, find1[2] = {0.0, 0.0};
+       Float_t good[2] = {0.0, 0.0}, good1[2] = {0.0, 0.0};
+       Float_t fake[2] = {0.0, 0.0}, fake1[2] = {0.0, 0.0};
+       
+       // histos filled with neural results
+       TH1D *hgood[2], *hfake[2], *hfind[2];
+       TH1D *hgood1[2], *hfake1[2], *hfind1[2];
+       TH1D *hg[2], *hf[2];
+       if (draw) {
+               hgood[0] = new TH1D("hgood0", "Good found tracks", N, edge);
+               hfake[0] = new TH1D("hfake0", "Fake found tracks", N, edge);
+               hfind[0] = new TH1D("hfound0", "Findable tracks", N, edge);
+               hgood[1] = new TH1D("hgood1", "Good found tracks", N, edge);
+               hfake[1] = new TH1D("hfake1", "Fake found tracks", N, edge);
+               hfind[1] = new TH1D("hfound1", "Findable tracks", N, edge);
+       
+               hgood[0]->Sumw2(); 
+               hfake[0]->Sumw2(); 
+               hfind[0]->Sumw2(); 
+               hgood[1]->Sumw2(); 
+               hfake[1]->Sumw2(); 
+               hfind[1]->Sumw2(); 
+               
+               // histos for evaluating percentual efficiency
+               hg[0] = new TH1D("hg0", "Efficiency (%) for #geq 5 right pts.", N, edge);
+               hf[0] = new TH1D("hf0", "Fake probability (%) for #geq 5 right pts.", N, edge);
+               hg[1] = new TH1D("hg1", "Efficiency (%) for 6 right pts.", N, edge);
+               hf[1] = new TH1D("hf1", "Fake probability (%) for 6 right pts.", N, edge);
+       }
+       
+       TFile *ffind;
+       TTree *tree;
+       
+       if (low) {
+               ffind = new TFile(Form("%s.root", filename), "READ");
+               tree  = (TTree*)ffind->Get("TreeF");
+       }
+       else {
+               ffind = new TFile(Form("%s_fnd.root", filename), "READ");
+               tree  = (TTree*)ffind->Get("tree");
+       }
+       
+       TFile *ftracks = new TFile(Form("%s_%d.root", filename, nsecs), "READ");
+       
+       
+       Double_t pt;
+       Int_t i, j, count, prim, *none = 0, div;
+       Int_t entries = tree->GetEntries(), label, min[] = {5,6};
+       tree->SetBranchAddress("pt", &pt);
+       tree->SetBranchAddress("count", &count);
+       tree->SetBranchAddress("prim", &prim);
+       
+       AliITSneuralTrack *trk = 0;
+       div = low ? 50 : 500;
+       for(i = 1;;i++) {
+               trk = (AliITSneuralTrack*)ftracks->Get(Form("AliITSneuralTrack;%d", i));
+               if (i%div == 0) cout << "\rEvaluating found track " << i << flush;
+               if (!trk) break;
+               for (j = 0; j < 2; j++) {
+                       label = trk->EvaluateTrack(0, min[j], none);
+                       tree->GetEntry(abs(label));
+                       if (count >= min[j]) {
+                               if (label > 0) {
+                                       if (draw) hgood[j]->Fill(pt);
+                                       good[j]++;
+                                       if (pt >= 1.0) good1[j]++;
+                               }
+                               else {
+                                       if (draw) hfake[j]->Fill(pt);
+                                       fake[j]++;
+                                       if (pt >= 1.0) fake1[j]++;
+                               }
+                       }
+               }
+       }
+       cout << endl;
+               
+       div = low ? 200 : 20000;
+       
+       for (i = 0; i < entries; i++) {
+               if (i%div == 0) cout << "\rEvaluating findable track no. " << i << flush;
+               tree->GetEntry(i);
+               for (j = 0; j < 2; j++) {
+                       if (count >= min[j]) {
+                               find[j]++;
+                               if (draw) hfind[j]->Fill(pt);
+                               if (pt >= 1.0) find1[j]++;
+                       }
+               }
+       }
+       cout << endl;
+       cout << hgood[0]->GetEntries() << " " << hgood[1]->GetEntries() << endl;
+       cout << hfake[0]->GetEntries() << " " << hfake[1]->GetEntries() << endl;
+       cout << hfind[0]->GetEntries() << " " << hfind[1]->GetEntries() << endl << endl;
+       
+       if (draw) {
+               TCanvas *canv[2];
+               canv[0] = new TCanvas("c_0", "Tracking efficiency (soft)", 0, 0, 600, 500);
+               canv[1] = new TCanvas("c_1", "Tracking efficiency (hard)", 630, 0, 600, 500);
+               
+               TLine *line1 = new TLine(1,100.0,edge[N],100.0); line1->SetLineStyle(4);
+               TLine *line2 = new TLine(1,90,edge[N],90); line2->SetLineStyle(4);
+                       
+               Bool_t good_drawn;
+               for (i = 0; i < 2; i++) {
+                       canv[i]->cd();
+                       good_drawn = kFALSE;
+                       if (hgood[i]->GetEntries() > 0.0) {
+                               good_drawn = kTRUE;
+                               hg[i]->Divide(hgood[i], hfind[i], 100.0, 1.0);
+                               hg[i]->SetMaximum(120);
+                               hg[i]->SetMinimum(0);
+                               hg[i]->SetMarkerStyle(21);
+                               hg[i]->SetMarkerSize(1);
+                               hg[i]->SetStats(kFALSE);
+                               hg[i]->GetXaxis()->SetTitle("pt (GeV/c)");
+                               hg[i]->Draw("PE1");
+                       }
+                       if (hfake[i]->GetEntries() > 0.0) {
+                               hf[i]->Divide(hfake[i], hfind[i], 100.0, 1.0);
+                               hf[i]->SetMaximum(120);
+                               hf[i]->SetMinimum(0);
+                               hf[i]->SetMarkerStyle(25);
+                               hf[i]->SetMarkerSize(1);
+                               hf[i]->SetStats(kFALSE);
+                               if (good_drawn)
+                                       hf[i]->Draw("PE1SAME");
+                               else
+                                       hf[i]->Draw("PE1");
+                       }
+                       line1->Draw("histosame");
+                       line2->Draw("histosame");
+                       canv[i]->Update();
+               }
+               canv[0]->SaveAs(Form("%s_soft.eps", filename));
+               canv[1]->SaveAs(Form("%s_hard.eps", filename));
+               cout << endl;
+       }
+       
+       Float_t sgood[2] = {0.0, 0.0}, sgood1[2] = {0.0, 0.0};
+       Float_t sfake[2] = {0.0, 0.0}, sfake1[2] = {0.0, 0.0};
+       for (i = 0; i < 2; i++) {
+               sgood[i] = error(good[i], find[i]);
+               sgood1[i] = error(good1[i], find1[i]);
+               sfake[i] = error(fake[i], find[i]);
+               sfake1[i] = error(fake1[i], find1[i]);
+               
+               good[i] = good[i] * 100.0 / find[i];
+               fake[i] = fake[i] * 100.0 / find[i];
+               good1[i] = good1[i] * 100.0 / find1[i];
+               fake1[i] = fake1[i] * 100.0 / find1[i];
+       }
+       
+       if (save) {
+               fstream data(save, ios::app);
+               data.setf(ios::fixed);
+               data.precision(1);
+               data << good1[0] << " " << fake1[0] << " ";
+               data << good1[1] << " " << fake1[1] << endl;
+               data.close();
+       }
+       else {
+               cout.setf(ios::fixed);
+               cout.precision(1);
+               cout << "*****************************************" << endl;
+               cout << "* Tracks with at least 5 correct points *" << endl;
+               cout << "*****************************************" << endl;
+               cout << "(all particles)" << endl;
+               cout << "Efficiency: " << good[0] << " +/- " << sgood[0] << "%" << endl;
+               cout << "Fake prob.: " << fake[0] << " +/- " << sfake[0] << "%" << endl;
+               if (!low) {
+                       cout << "(pt >= 1 GeV/c)" << endl;
+                       cout << "Efficiency: " << good1[0] << " +/- " << sgood1[0] << "%" << endl;
+                       cout << "Fake prob.: " << fake1[0] << " +/- " << sfake1[0] << "%" << endl;
+               }
+               cout << endl;
+               cout << "************************************" << endl;
+               cout << "* Tracks with all 6 correct points *" << endl;
+               cout << "************************************" << endl;
+               cout << "(all particles)" << endl;
+               cout << "Efficiency: " << good[1] << " +/- " << sgood[1] << "%" << endl;
+               cout << "Fake prob.: " << fake[1] << " +/- " << sfake[1] << "%" << endl;
+               if (!low) {
+                       cout << "(pt >= 1 GeV/c)" << endl;
+                       cout << "Efficiency: " << good1[1] << " +/- " << sgood1[1] << "%" << endl;
+                       cout << "Fake prob.: " << fake1[1] << " +/- " << sfake1[1] << "%" << endl;
+               }
+               cout << endl;
+       }
+}
+
+Double_t error(Double_t x, Double_t y) {
+       Double_t radq = x + x * x / y;
+       return 100.0 * sqrt(radq) / y;
+}
diff --git a/ITS/oldmacros/AliITSNeuralTracking.C b/ITS/oldmacros/AliITSNeuralTracking.C
new file mode 100644 (file)
index 0000000..b286461
--- /dev/null
@@ -0,0 +1,111 @@
+// ARGUMENTS:
+// 1. number of azymuthal sectors (it's better not to go under 8 or over 40)
+// 2. the ROOT file to read (WITHOUT exstension)
+// 3. event number
+// 4. if specified a string, a fstream named like the argument is opened and
+//    the elapsed CPU time is stored (not useful)
+
+// the macro will save a file named, for example "galice_<nsecs>.root"
+// containing may AliITSneuralTrack objects
+
+void AliITSNeuralTracking
+(Int_t nsecs = 12, const char* rfile = "galice", Int_t event = 0, const char* save = 0)
+{
+       TStopwatch timer;
+       Double_t CONVERT = TMath::Pi() / 180.0;
+       const char* wfile = Form("%s_%d.root", rfile, nsecs);
+       cout << "Reading file " << rfile << ".root and saving in " << wfile << endl;
+       
+// ==================================
+// ==== CURVATURE CUT DEFINITION ====
+// ==================================
+
+       // These values define the curvature cuts for all steps 
+       // within a sector.
+       // For a greater clarity, the cuts are indicated in units
+       // of transverse momentum (GeV/c) but these value have no
+       // exact physical meaning, but are useful to understand
+       // well what means a choice in the value of a certain 
+       // curvature constraint
+       // NOTE: becareful to make sure that the 'ncuts' variable
+       //       have the same value of the dimension of the allocated arrays
+       
+       Int_t ncuts;
+       Double_t *p, *cut;
+       
+       ncuts = 5;
+       p = new Double_t[5];
+       cut = new Double_t[5];  
+       p[0] = 2.0;
+       p[1] = 1.0;
+       p[2] = 0.7;
+       p[3] = 0.5;
+       p[4] = 0.3;
+       
+       for (Int_t i = 0; i < ncuts; i++) cut[i] = 0.003 * 0.2 / p[i];
+       
+
+// ==========================
+// ==== OTHER PARAMETERS ====
+// ==========================
+           
+       Bool_t   flag   = kFALSE; // for now, don't change this line, please...
+       
+       Double_t diff   = 0.02;   // helicoidal cut 
+       Double_t dtheta = 1.0;    // delta-theta cut
+       Double_t temp   = 1.0;    // temperature parameter
+       Double_t var    = 0.0001; // stabilization threshold
+       
+       Double_t exp    = 7.0;    // straight-line excitator
+       Double_t gtoc   = 3.0;    // gain/cost contribution ratio
+       
+       Double_t min    = 0.4;    // minimum in random activation initialization
+       Double_t max    = 0.6;    // maximum in random activation initialization
+       
+       
+// =========================
+// ==== NEURAL TRACKING ====
+// =========================
+       
+       AliITSneuralTracker *ANN = new AliITSneuralTracker(nsecs, ncuts, cut, CONVERT*dtheta);
+               
+       TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(Form("%s.root", rfile));
+       if (!file) file = new TFile(Form("%s.root", rfile),"UPDATE");
+
+       //Double_t Xv = -0.001097;
+       //Double_t Yv = -0.00347647;
+       //Double_t Zv =  0.000631345;
+       //ANN->SetVertex(Xv, Yv, Zv);
+       // You should find the vertex with VertexMacro.C
+       // and then put by hand the found values with
+       // the above method.
+       
+       Int_t points = ANN->ReadFile(Form("%s.root", rfile), event);
+               
+       ANN->SetTemperature(temp);
+       ANN->SetVariationLimit(var);
+       ANN->SetGainToCostRatio(gtoc);
+       ANN->SetExponent(exp);
+       ANN->SetInitLimits(min, max);
+       ANN->SetDiff(diff);
+       
+       cout << points << " points found " << endl << endl;
+       
+       TStopwatch timer;
+       timer.Start();
+       
+       ANN->Go(wfile, flag);
+       
+       timer.Stop();
+       cout << endl;
+       timer.Print();
+       
+       if (save) {
+               fstream ksave(save, ios::app);
+               ksave << nsecs << " " << timer->CpuTime() << endl;
+               ksave.close();
+       }
+       
+//     delete gAlice;
+//     gAlice = 0;
+}
diff --git a/ITS/oldmacros/AliITSOccupancy.C b/ITS/oldmacros/AliITSOccupancy.C
new file mode 100644 (file)
index 0000000..a1f3db6
--- /dev/null
@@ -0,0 +1,223 @@
+/******************************************************************************
+
+  "AliITSOccupancy.C"
+  
+  this macro calculates the mean occupancy of each ITS layer, counting the
+  "fired" digits of each module, and making two overall calculations:
+     1) the calculus of the mean overall layer occupancy (as the ratio
+            between the total number of active digits and the total number of 
+                 channels in the module;
+         2) a distribution of the occupancies as a funcion of z, to obtain some
+                 kind of most significand data for this value, along the z axis
+       
+  The macro creates also the GIF and the EPS of the TCanvas where the 
+  histograms are plotted. Muliple input files are supported (see the 
+  input argument list below).
+                         
+                                                                        
+  It is also possible (and recommended) to compile this macro, 
+  to execute it faster. To do this, it is necessary to:
+  1) type, at the root prompt, the instruction 
+     gSystem->SetIncludePath("-I- -I$ALICE_ROOT/include -I$ALICE_ROOT/ITS -g")
+         to adjoin the ALICE include libraries to the main include directory
+  2) load the function via the ".L AliITSOccupancy.C++" statement 
+     (only the first time, because the first time the macro must be compiled). 
+      From the second time you use the macro, you must only load it via the 
+         ".L AliITSOccupancy.C+" instruction (note the number of '+' signs in 
+      each case
+  3) call the function with only its name, e.g.: AliITSOccupancy()
+
+ Author: Alberto Pulvirenti                                             
+******************************************************************************/
+
+
+#if !defined(__CINT__) || defined(__MAKECINT__)
+       #include <iostream.h>
+       #include <TROOT.h>
+       #include <TArrayI.h>
+       #include <TCanvas.h>
+       #include <TClassTable.h>
+       #include <TClonesArray.h>
+       #include <TFile.h>
+       #include <TObject.h>
+       #include <TObjArray.h>
+       #include <TTree.h>
+       #include <TMath.h>
+       #include <TString.h>
+       #include <TH1.h>
+       #include <AliRun.h>
+       #include <AliITS.h>
+       #include <AliITSgeom.h>
+       #include <AliITSDetType.h>
+       #include <AliITSRecPoint.h>
+       #include <AliITSdigit.h>
+       #include <AliITShit.h>
+       #include <AliITSmodule.h> 
+       #include <AliITSsegmentation.h>
+       #include <AliITSsegmentationSPD.h> 
+       #include <AliITSsegmentationSDD.h>
+       #include <AliITSsegmentationSSD.h>
+#endif
+
+TFile* AccessFile(TString inFile="galice.root", TString acctype="R");
+
+Int_t AliITSOccupancy(TString FileHits="galice.root", TString FileDigits="galice.root", Int_t evNum = 0) {       
+                       
+  // Open the Root input file containing the AliRun data and 
+  // the Root output data file
+  TFile *froot = AccessFile(FileHits);
+  froot->ls();
+
+  if(!(FileDigits.Data() == FileHits.Data())){
+    gAlice->SetTreeDFileName(FileDigits);
+  }
+  gAlice->GetEvent(evNum);
+  // Initialize the ITS and its geometry
+  AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
+  AliITSgeom *gm = ITS->GetITSgeom();
+
+  // Variables and objects definition and setting
+  Float_t zmax[] = {20.0,20.0,25.0,35.0,49.5,54.0}; // z edges for TH1F (cm)
+  Int_t nbins[]  = {40,40,50,70,99,108};             // bins number for TH1Fs
+       
+  // Histos for plotting occupancy distributions
+  TH1F *above[6], *below[6];
+       
+  for (Int_t lay = 0; lay < 6; lay++) {
+    Int_t nlads = gm->GetNladders(lay+1);
+    Int_t ndets = gm->GetNdetectors(lay+1);
+    Int_t dtype = lay / 2;
+    Int_t minID = gm->GetModuleIndex(lay+1, 1, 1);
+    Int_t maxID = gm->GetModuleIndex(lay+1, nlads, ndets);
+    Text_t ffname[20];
+    sprintf(ffname, "h_%d", lay+1);
+    below[lay] = new TH1F(ffname, "Total z distribution of digits", nbins[lay], -zmax[lay], zmax[lay]);
+    cout << "Evaluating total channels number of layer " << lay+1 << endl;
+    for (Int_t I = minID; I <= maxID; I++) {           
+      AliITSDetType *det = ITS->DetType(dtype);
+      AliITSsegmentation *seg = det->GetSegmentationModel();
+      Int_t NX = seg->Npx();
+      if(lay == 2 || lay == 3) NX = 192;
+      Int_t NZ = seg->Npz();
+      //                       cout << "ID: " << I << ", Layer: " << lay+1 << ", NX = " << NX << ", NZ = " << NZ << endl;
+      for (Int_t ix = 0; ix <= NX; ix++) {
+        for (Int_t iz = 0; iz <= NZ; iz++) {
+          Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
+          seg->DetToLocal(ix, iz, lx[0], lx[2]);
+          gm->LtoG(I, lx, gx);
+          below[lay]->Fill(gx[2]);
+        }
+      }
+    }
+    sprintf(ffname, "H_%d", lay+1);
+    above[lay] = new TH1F(ffname, "histo", nbins[lay], -zmax[lay], zmax[lay]);
+  }
+               
+  // Counting the hits, digits and recpoints contained in the ITS
+  TTree *TD = gAlice->TreeD();
+  ITS->ResetDigits();
+  Float_t mean[6];
+  Float_t average[6];
+
+  for (Int_t L = 0; L < 6; L++) {
+               
+    cout << "Layer " << L + 1 << ": " << flush;                                
+
+    // To avoid two nested loops, are calculated 
+    // the ID of the first and last module of the L
+    // (remember the L goest from 0 to 5 (not from 1 to 6)
+    Int_t first, last;
+    first = gm->GetModuleIndex(L + 1, 1, 1);
+    last = gm->GetModuleIndex(L + 1, gm->GetNladders(L + 1), gm->GetNdetectors(L + 1));
+                               
+    // Start loop on modules
+    for (Int_t ID = first; ID <= last; ID++) {
+      Int_t dtype = L / 2;
+      AliITSDetType *det = ITS->DetType(dtype);
+      AliITSsegmentation *seg = det->GetSegmentationModel();
+      if (dtype == 2) seg->SetLayer(L+1);
+                       
+      TD->GetEvent(ID);
+      TClonesArray *digits_array = ITS->DigitsAddress(dtype);
+      Int_t digits_num = digits_array->GetEntries();
+      // Get the coordinates of the module
+      for (Int_t j = 0; j < digits_num; j++) {
+        Float_t lx[] = {0.0,0.0,0.0}, gx[] = {0.0,0.0,0.0};
+        AliITSdigit *digit = (AliITSdigit*)digits_array->UncheckedAt(j);
+        Int_t iz=digit->fCoord1;  // cell number z
+        Int_t ix=digit->fCoord2;  // cell number x
+        // Get local coordinates of the element (microns)
+        seg->DetToLocal(ix, iz, lx[0], lx[2]);
+        gm->LtoG(ID, lx, gx);
+        above[L]->Fill(gx[2]);
+      }
+    }
+               
+    Float_t occupied = above[L]->GetEntries();
+    Float_t total = below[L]->GetEntries();
+    cout << "Entries: " << occupied << ", " << total << endl;
+    average[L] = 100.*occupied/total;
+    above[L]->Divide(above[L], below[L], 100.0, 1.0);
+    mean[L] = above[L]->GetSumOfWeights() / above[L]->GetNbinsX(); 
+    cout.setf(ios::fixed);
+    cout.precision(2);
+    cout << " digits occupancy = " << mean[L] << "%" << endl;
+    cout << " average digits occupancy = " << average[L] << "%" << endl;
+  }
+       
+  TCanvas *view = new TCanvas("view", "Digits occupancy distributions", 600, 900);
+  view->Divide(2, 3);
+       
+  for (Int_t I = 1; I < 7; I++) {
+    view->cd(I);
+    Text_t title[50];
+    sprintf(title, "Layer %d: %4.2f%c", I, mean[I-1], '%');
+    title[6] = (Text_t)I + '0';
+    above[I-1]->SetTitle(title);
+    above[I-1]->SetStats(kFALSE);
+    above[I-1]->SetXTitle("z (cm)");
+    above[I-1]->SetYTitle("%");
+    above[I-1]->Draw();
+    view->Update();
+  }
+
+  view->SaveAs("AliITSOccupancy_digit.gif");
+  view->SaveAs("AliITSOccupancy_digit.eps");
+  TFile *fOc = new TFile("AliITSOccupancy.root","recreate");
+  fOc->cd();
+  for (Int_t I = 0; I < 6; I++) above[I]->Write();
+  fOc->Close();
+
+  return 1;
+}
+
+//-------------------------------------------------------------------
+TFile * AccessFile(TString FileName, TString acctype){
+
+  // Function used to open the input file and fetch the AliRun object
+
+  TFile *retfil = 0;
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(FileName);
+  if (file) {file->Close(); delete file; file = 0;}
+  if(acctype.Contains("U")){
+    file = new TFile(FileName,"update");
+  }
+  if(acctype.Contains("N") && !file){
+    file = new TFile(FileName,"recreate");
+  }
+  if(!file) file = new TFile(FileName);   // default readonly
+  if (!file->IsOpen()) {
+       cerr<<"Can't open "<<FileName<<" !" << endl;
+       return retfil;
+  } 
+
+  // Get AliRun object from file or return if not on file
+  if (gAlice) {delete gAlice; gAlice = 0;}
+  gAlice = (AliRun*)file->Get("gAlice");
+  if (!gAlice) {
+       cerr << "AliRun object not found on file"<< endl;
+       return retfil;
+  } 
+  return file;
+}
diff --git a/ITS/oldmacros/AliITSPlotGeom.C b/ITS/oldmacros/AliITSPlotGeom.C
new file mode 100644 (file)
index 0000000..abda472
--- /dev/null
@@ -0,0 +1,183 @@
+void AliITSPlotGeom(const char *filename="galice.root"){
+/////////////////////////////////////////////////////////////////////////
+//   This macro displays part of the ITS sensitive volume as defined
+// in the root file.
+//   
+//     Root > .L AliITSPlotGeom.C         //this loads the macro in memory
+//     Root > AliITSPlotGeom();           //by default process first event   
+//     Root > AliITSPlotGeom("galice2.root"); //process third event from 
+//                                              galice2.root file.
+//Begin_Html
+/*
+<img src="picts/AliITSplotgeom.gif">
+*/
+//End_Html
+/////////////////////////////////////////////////////////////////////////
+    if(gAlice){
+       delete gAlice;
+       gAlice=0;
+    }else{
+       // Dynamically link some shared libs
+       if(gClassTable->GetID("AliRun") < 0) {
+           gROOT->LoadMacro("loadlibs.C");
+           loadlibs();
+       } // end if
+    } // end if gAlice
+    // Connect the Root Galice file containing Geometry, Kine and Hits
+    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+    if(!file) file = new TFile(filename);
+
+    // Get AliRun object from file or create it if not on file
+    if(!gAlice) {
+       gAlice = (AliRun*)file->Get("gAlice");
+       if(gAlice) printf("AliRun object found on file\n");
+       if(!gAlice) gAlice = new AliRun("gAlice","Alice test program");
+    } // end if !gAlice
+      
+
+// Get pointers to the ITS Alice detectors and Hits containers
+    AliITS    *ITS    = (AliITS*)    gAlice->GetDetector("ITS");
+    if(ITS==0){
+       cout << "ITS not found. Exiting." << endl;
+       return;
+    } // end if ITS==0
+    AliITSgeom *itsgeom = ITS->GetITSgeom();
+    if(itsgeom==0){
+       cout << "ITSgeom not defined. Exiting." << endl;
+       return;
+    } // end if itsgeom==0
+    Int_t is,ie,in;
+    is = itsgeom->GetStartSPD();
+    ie = itsgeom->GetLastSPD()+1;
+    in = ie-is;
+
+    // Transformations.
+//    Float_t *x = new Float_t[in];
+//    Float_t *y = new Float_t[in];
+//    Float_t *z = new Float_t[in];
+    Float_t x,y,z;
+    TRotMatrix **r = new TRotMatrix*[in];
+    Int_t i,j;
+    Double_t m[9];
+    char name[10],title[20],name2[10],title2[20];
+//    cout << in << endl;
+    for(i=0;i<in;i++){
+       itsgeom->GetRotMatrix(i+is,m);
+       sprintf(name,"ROT%d",i+is);
+       sprintf(title,"ROT%d",i+is);
+       r[i] = new TRotMatrix(name,title,m);
+//     cout << name << title << endl;
+//     itsgeom->GetTrans(i+is,x[i],y[i],z[i]);
+//     cout << i << " " << x[i] << " " << y[i] << " " << z[i] <<endl;
+    } // end for i
+
+    // Sensitive volumes.
+    if(itsgeom->IsShapeDefined(0)){
+       TBRIK *spds = (TShape*) (itsgeom->GetShape(0));
+    }else{
+       TBRIK *spds = new TBRIK("SPD","SPD","void",0.64,0.015,3.48);
+    } // end if
+//    cout << spds << endl;
+    if(itsgeom->IsShapeDefined(1)){
+       TBRIK *sdds = (TShape*) (itsgeom->GetShape(240));
+    }else{
+       TBRIK *sdds = new TBRIK("SDD","SDD","void",3.50850,0.01499,3.76320);
+    } // end if
+//    cout << sdds << endl;
+    if(itsgeom->IsShapeDefined(2)){
+       TBRIK *ssds = (TShape*) (itsgeom->GetShape(1000));
+    }else{
+       TBRIK *ssds = new TBRIK("SSD","SSD","void",3.65,0.015,2.00);
+    } // end if
+//    cout << ssds << endl;
+
+    // Set up display.
+
+    TCanvas *c1 = new TCanvas("c1","ITS geometry",10,10,700,700);
+//  TPad *p1 = new TPad("p1","p1",0.01,0.01,0.99,0.99,46,3,1);
+//  p1->Draw();
+//  p1->cd();
+//  TView *view = new TView(1);
+//  view->SetRange(-5,-5,-5,5,5,5);
+      TShape  *mother = new TBRIK("Mother","Mother Volume","void",10,10,10);
+//      TShape  *mother = new TBRIK("Mother","Mother Volume","void",30,30,30);
+//      TShape  *mother = new TBRIK("Mother","Mother Volume","void",50,50,50);
+      TNode *node = new TNode("node","Mother Node",mother);
+      node->cd();
+
+    // Set up nodes
+    TNode **itsn = new TNode*[in];
+    TPolyLine3D **pl = new TPolyLine3D*[in];
+/*    Double_t p[5*3],axis[5*3]={1.,0.,0.,
+                               0.,0.,0.,
+                              0.,1.,0.,
+                              0.,0.,0.,
+                              0.,0.,1.}; */
+    Double_t p[19*3],axis[19*3]={-0.25,0.,1.25,  // 1
+                                +0.00,0.,1.25,  // 2
+                                -0.25,0.,1.00,  // 3
+                                +0.00,0.,1.00,  // 4
+                                +0.00,0.,0.00,  // 5
+                                +1.00,0.,0.00,  // 6
+                                +1.25,0.,0.25,  // 7
+                                +1.125,0.,0.125,  // 8
+                                +1.00,0.,0.25,  // 9
+                                +1.125,0.,0.125,  // 10
+                                +1.25,0.,0.00,  // 11
+                                +1.125,0.,0.125,  // 12
+                                +1.00,0.,0.00,  // 13
+                                +0.00,0.,0.00,  // 14
+                                +0.00,1.,0.00,  // 15
+                                +0.00,1.,0.125,  // 16
+                                +0.25,1.,0.25,  // 17
+                                +0.00,1.,0.125,  // 18
+                                -0.25,1.,0.25,  // 19
+    };
+    for(i=0;i<in;i++){
+       itsgeom->GetTrans(i+is,x,y,z);
+/*     switch (itsgeom->GetGeomMatrix(i+is)->GetDetectorIndex())
+        case 0:
+           sprintf(name,"SPD%d",i+is);
+           sprintf(title,"SPD%d",i+is);
+           sprintf(name2,"BSPD%d",i+is);
+           sprintf(title2,"BSPD%d",i+is);
+           itsn[i] = new TNode(name,title,new TBRIK(name2,title2,"void",
+                      spds->GetDx(),spds->GetDy(),spds->GetDz()),x,y,z,r[i]);
+           break;
+        case 1:
+           sprintf(name,"SDD%d",i+is);
+           sprintf(title,"SDD%d",i+is);
+           sprintf(name2,"BSDD%d",i+is);
+           sprintf(title2,"BSDD%d",i+is);
+           itsn[i] = new TNode(name,title,new TBRIK(name2,title2,"void",
+                      sdds->GetDx(),sdds->GetDy(),sdds->GetDz()),x,y,z,r[i]);
+           break;
+        case 2:
+           sprintf(name,"SSD%d",i+is);
+           sprintf(title,"SSD%d",i+is);
+           sprintf(name2,"BSSD%d",i+is);
+           sprintf(title2,"BSSD%d",i+is);
+           itsn[i] = new TNode(name,title,new TBRIK(name2,title2,"void",
+                      ssds->GetDx(),ssds->GetDy(),ssds->GetDz()),x,y,z,r[i]);
+           break; */
+       for(j=0;j<19;j++) itsgeom->LtoG(i+is,(Double_t*)&axis[3*j],
+                                            (Double_t*)&p[3*j]);
+       pl[i] = new TPolyLine3D(19,p);
+    } // end for i
+
+    // display it
+    node->cd();
+    node->Draw();
+//    for(i=0;i<in;i++) itsn[i]->Draw();
+//    node->Draw();
+    for(i=0;i<in;i++) pl[i]->Draw();
+    c1->Update();
+
+    // clean up
+//    delete[] x;
+//    delete[] y;
+//    delete[] z;
+//    for(i=0;i<itsgeom->GetIndexMax();i++) delete[] r[i];
+//    delete[] r;
+}
+
diff --git a/ITS/oldmacros/AliITSPlotTracksV1.C b/ITS/oldmacros/AliITSPlotTracksV1.C
new file mode 100644 (file)
index 0000000..7c7742b
--- /dev/null
@@ -0,0 +1,140 @@
+#include <iostream.h>
+#include <fstream.h>
+
+void AliITSPlotTracksV1()
+{
+       TFile *file = new TFile("AliITSComparisonV1.root");
+       TTree *tree = (TTree*)file->Get("Eval"); 
+       
+       if (!tree) {
+               cerr << "Evaluation tree not found!" << endl;
+               return;
+       }
+       
+       // histogram definition - I (efficiency)
+       TH1D *hFindables = new TH1D("hFindables", "Findable tracks", 10, 0, 2);
+       TH1D *hGood = new TH1D("hGood", "Good found tracks", 10, 0, 2);
+       TH1D *hFake = new TH1D("hFake", "Fake found tracks", 10, 0, 2);
+       TH1D *hRatioG = new TH1D("hRatioG", "", 10, 0, 2);
+       TH1D *hRatioF = new TH1D("hRatioF", "", 10, 0, 2);
+       hGood->Sumw2(); 
+       hFake->Sumw2();  
+       hFindables->Sumw2();
+       hRatioG->SetLineColor(kBlue); 
+       hRatioG->SetLineWidth(2);
+       hRatioF->SetLineColor(kRed); 
+       hRatioF->SetLineWidth(2);
+       
+       // histograms definition - II (resolution)
+       TH1D *hPhi    = new TH1D("hPhi", "#Phi resolution", 50, -15., 15.); 
+       TH1D *hLambda = new TH1D("hLambda", "#lambda resolution", 50, -15., 15.); 
+       TH1D *hPt     = new TH1D("hPt", "Relative Pt resolution", 40, -10., 10.);
+       TH1D *hDtot   = new TH1D("hDtot", "Total impact parameter distribution", 100, 0., 2000.); 
+       TH1D *hDr     = new TH1D("hDr", "Transv. impact parameter distribution", 50, -1000., 1000.);  
+       TH1D *hDz     = new TH1D("hDz", "Long. impact parameter distribution", 50, -1000., 1000.);  
+       hPhi->SetFillColor(4);
+       hLambda->SetFillColor(4);
+       hPt->SetFillColor(2); 
+       hDtot->SetFillColor(6);
+       hDr->SetFillColor(kGreen);
+       hDz->SetFillColor(kBlue);
+       
+       // Evaluation tree settings
+       Int_t labITS, labTPC, signC;
+       Int_t i, j, tot = (Int_t)tree->GetEntries();
+       Double_t difpt, diflambda, difphi, Dz, Dr, Dtot, ptg;
+       tree->SetBranchAddress("labITS"   , &labITS   );
+       tree->SetBranchAddress("labTPC"   , &labTPC   );
+       tree->SetBranchAddress("difpt"    , &difpt    ); 
+       tree->SetBranchAddress("diflambda", &diflambda);
+       tree->SetBranchAddress("difphi"   , &difphi   );
+       tree->SetBranchAddress("Dz"       , &Dz       );
+       tree->SetBranchAddress("Dr"       , &Dr       );
+       tree->SetBranchAddress("Dtot"     , &Dtot     );
+       tree->SetBranchAddress("ptg"      , &ptg      );
+       tree->SetBranchAddress("signC"    , &signC    );
+       
+       // Filling the histogram of findable tracks (w.r. to momentum)
+       for(i = 0; i < tot; i++) {
+               tree->GetEntry(i);
+               hFindables->Fill(ptg);
+       }
+       
+       // Filling the evaluation histograms
+       Int_t neglabs = 0;
+       for(i = 0; i < tot; i++) {
+               tree->GetEntry(i);
+//             if(signC<0) continue;
+               if (labITS < 0) neglabs++;
+               if (labITS >= 0) {
+                       hGood->Fill(ptg); 
+                       hPt->Fill(difpt);
+                       hLambda->Fill(diflambda);
+                       hPhi->Fill(difphi);
+                       hDtot->Fill(Dtot);
+                       hDr->Fill(Dr);
+                       hDz->Fill(Dz);
+               }
+               else {
+                       hFake->Fill(ptg);
+                       neglabs++;
+               }
+       }
+  
+       // Drawing
+       cerr << "Findable tracks   : " << hFindables->GetEntries() << endl;
+       cerr << "Good found tracks : " << hGood->GetEntries() << endl;
+       cerr << "Fake found tracks : " << hFake->GetEntries() << endl;
+       
+       gStyle->SetOptStat(111110);
+       gStyle->SetOptFit(1);
+       
+       TCanvas *c1 = new TCanvas("c1","Parameter resolutions",0,0,700,700);
+       c1->Divide(2, 2, 0.001, 0.001);
+       c1->cd(1); hPhi->SetXTitle("(mrad)"); hPhi->Draw(); hPhi->Fit("gaus", "Q"); 
+       c1->cd(2); hLambda->SetXTitle("(mrad)"); hLambda->Draw(); hLambda->Fit("gaus", "Q");
+       c1->cd(3); hPt->SetXTitle("(%)"); hPt->Draw(); hPt->Fit("gaus", "Q");
+       c1->cd(4); hDtot->SetXTitle("(micron)"); hDtot->Draw();
+       c1->Update();
+       
+       TCanvas *c2 = new TCanvas("c2","Impact parameters resolution",100,100,700,400);
+       c2->Divide(2, 1, 0.001, 0.001);
+       c2->cd(1); hDr->SetXTitle("(micron)"); hDr->Draw(); hDr->Fit("gaus", "Q"); 
+       c2->cd(2); hDz->SetXTitle("(micron)"); hDz->Draw(); hDz->Fit("gaus", "Q"); 
+       c2->Update();
+       
+       TCanvas *c3 = new TCanvas("c3","Momentum distributions",200,200,800,500);
+       c3->Divide(3, 1, 0.001, 0.001);
+       c3->cd(1); hGood->Draw(); 
+       c3->cd(2); hFake->Draw(); 
+       c3->cd(3); hFindables->Draw();
+       c3->Update();
+       
+       TCanvas *c4 = new TCanvas("c4","Tracking efficiency",300,300,800,500);
+       hRatioG->Divide(hGood, hFindables, 1., 1.);
+       hRatioF->Divide(hFake, hFindables, 1., 1.);
+       hRatioG->SetMaximum(1.4);
+       hRatioG->SetYTitle("Tracking efficiency");
+       hRatioG->SetXTitle("Pt (GeV/c)");
+       hRatioG->Draw();  // to not plot the erro bars    hg->Draw("histo");
+       hRatioF->SetFillColor(1);
+       hRatioF->SetFillStyle(3013);
+       hRatioF->SetLineColor(2);
+       hRatioF->SetLineWidth(2);
+       hRatioF->Draw("same");  // to not plot the error bars  hRatioF->Draw("histosame");
+       // a line to mark the best efficiency
+       TLine *line1 = new TLine(0,1.0,2,1.0); line1->SetLineStyle(4);
+       line1->Draw("same");
+       TLine *line2 = new TLine(0,0.9,2,0.9); line2->SetLineStyle(4);
+       line2->Draw("histosame");
+       // a text explanation
+       TText *text = new TText(0.461176,0.248448,"Fake tracks");
+       text->SetTextSize(0.05);
+       text->Draw();
+       text = new TText(0.453919,1.11408,"Good tracks");
+       text->SetTextSize(0.05);
+       text->Draw();
+       c4->Update();
+       
+       cout << "neglabs = " << neglabs << endl;  // provvisoria
+}
diff --git a/ITS/oldmacros/AliITSPrintDigits.C b/ITS/oldmacros/AliITSPrintDigits.C
new file mode 100644 (file)
index 0000000..512d305
--- /dev/null
@@ -0,0 +1,92 @@
+#include "iostream.h"
+#include "TFile.h"
+#include "TString.h"
+#include "TClonesArray.h"
+/*
+#include "$(ALICE_ROOT)/STEER/AliRun.h"
+#include "$(ALICE_ROOT)/ITS/AliITS.h"
+#include "$(ALICE_ROOT)/ITS/AliITSgeom.h"
+#include "$(ALICE_ROOT)/ITS/AliITSdigit.h"
+*/
+void AliITSPrintDigits(TString rfn="galice.root",Int_t mod=-1,
+                         Int_t evnt=-1){
+    // Macro to print out the recpoints for all or a specific module
+
+    // Dynamically link some shared libs
+    if (gClassTable->GetID("AliRun") < 0) {
+        gROOT->LoadMacro("loadlibs.C");
+        loadlibs();
+    } // end if
+    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSstandard.C");
+
+    TFile *rf=0;
+    rf = AccessFile(rfn,"R"); // Set up to read in Data
+    AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");
+    if(!ITS){
+       cout << "Error: no ITS found. Aborting"<<endl;
+       return;
+    } // end if !ITS
+
+    Int_t evNumber1 = 0;
+    Int_t evNumber2 = gAlice->GetEventsPerRun();
+    if(evnt>=0){
+       evNumber1 = evnt;
+       evNumber2 = evnt+1;
+    } // end if evnt>=0
+    Int_t mod1 = 0;
+    Int_t mod2 = ITS->GetITSgeom()->GetIndexMax();
+    if(mod>=0){
+       mod1 = mod;
+       mod2 = mode+1;
+    } // end if mod>=0
+    TClonesArray   *pda = 0;
+    AliITSdigitSPD *dp0 = 0;
+    AliITSdigitSDD *dp1 = 0;
+    AliITSdigitSSD *dp2 = 0;
+
+    Int_t event,m,i,i2,id;
+    for(event = evNumber1; event < evNumber2; event++){
+        gAlice->GetEvent(event);
+       for(m=mod1;m<mod2;m++){
+           id = ITS->GetITSgeom()->GetModuleType(m);
+           dpa = ITS->DigitsAddress(id);
+           ITS->ResetDigits();
+           gAlice->TreeD()->GetEvent(m);
+           i2 = dpa->GetEntriesFast();
+           switch (id) {
+           case 0:
+               cout <<  "Event=" << event << " module=" << m <<
+                   " Number of SPD Digits=" << i2 <<endl;
+               for(i=0;i<i2;i++){
+                   dp0 = (AliITSdigitSPD*)(dpa->At(i));
+                   cout << i << " ";
+                   dp0->Print((ostream*)cout);
+                   cout << endl;
+               } // end for i
+               break;
+           case 1:
+               cout <<  "Event=" << event << " module=" << m <<
+                   " Number of SDD Digits=" << i2 <<endl;
+               for(i=0;i<i2;i++){
+                   dp1 = (AliITSdigitSDD*)(dpa->At(i));
+                   cout << i << " ";
+                   dp1->Print((ostream*)cout);
+                   cout << endl;
+               } // end for i
+               break;
+           case 2:
+               cout <<  "Event=" << event << " module=" << m <<
+                   " Number of SSD Digits=" << i2 <<endl;
+               for(i=0;i<i2;i++){
+                   dp2 = (AliITSdigitSSD*)(dpa->At(i));
+                   cout << i << " ";
+                   dp2->Print((ostream*)cout);
+                   cout << endl;
+               } // end for i
+               break;
+           default:
+               break;
+           } // end switch
+       } // end for m
+    } // end for event
+}
diff --git a/ITS/oldmacros/AliITSPrintSDigits.C b/ITS/oldmacros/AliITSPrintSDigits.C
new file mode 100644 (file)
index 0000000..f488d57
--- /dev/null
@@ -0,0 +1,63 @@
+#include "iostream.h"
+#include "TFile.h"
+#include "TString.h"
+#include "TClonesArray.h"
+/*
+#include "$(ALICE_ROOT)/STEER/AliRun.h"
+#include "$(ALICE_ROOT)/ITS/AliITS.h"
+#include "$(ALICE_ROOT)/ITS/AliITSgeom.h"
+#include "$(ALICE_ROOT)/ITS/AliITSpList.h"
+*/
+void AliITSPrintSDigits(TString rfn="galice.root",Int_t mod=-1,
+                         Int_t evnt=-1){
+    // Macro to print out the recpoints for all or a specific module
+
+    // Dynamically link some shared libs
+    if (gClassTable->GetID("AliRun") < 0) {
+        gROOT->LoadMacro("loadlibs.C");
+        loadlibs();
+    } // end if
+    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSstandard.C");
+
+    TFile *rf=0;
+    rf = AccessFile(rfn,"R"); // Set up to read in Data
+    AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");
+    if(!ITS){
+       cout << "Error: no ITS found. Aborting"<<endl;
+       return;
+    } // end if !ITS
+
+    Int_t evNumber1 = 0;
+    Int_t evNumber2 = gAlice->GetEventsPerRun();
+    if(evnt>=0){
+       evNumber1 = evnt;
+       evNumber2 = evnt+1;
+    } // end if evnt>=0
+    Int_t mod1 = 0;
+    Int_t mod2 = ITS->GetITSgeom()->GetIndexMax();
+    if(mod>=0){
+       mod1 = mod;
+       mod2 = mode+1;
+    } // end if mod>=0
+    TClonesArray *sdpa = ITS->GetSDigits();
+    AliITSpListItem *sdp = 0;
+
+    Int_t event,m,i,i2;
+    for(event = evNumber1; event < evNumber2; event++){
+        gAlice->GetEvent(event);
+       for(m=mod1;m<mod2;m++){
+           ITS->ResetSDigits();
+           gAlice->TreeS()->GetEvent(m);
+           i2 = sdpa->GetEntriesFast();
+           cout <<  "Event=" << event << " module=" << m <<
+               " Number of SDigits=" << i2 <<endl;
+           for(i=0;i<i2;i++){
+               sdp = (AliITSpListItem*)(sdpa->At(i));
+               cout << i << " ";
+               sdp->Print((ostream*)cout);
+               cout << endl;
+           } // end for i
+       } // end for m
+    } // end for event
+
+}
diff --git a/ITS/oldmacros/AliITSPropagateBackV2.C b/ITS/oldmacros/AliITSPropagateBackV2.C
new file mode 100644 (file)
index 0000000..4ae6e77
--- /dev/null
@@ -0,0 +1,41 @@
+#ifndef __CINT__
+  #include <Riostream.h>
+  #include "AliITSgeom.h"
+  #include "AliITStrackerV2.h"
+
+  #include "TFile.h"
+  #include "TStopwatch.h"
+#endif
+
+Int_t AliITSPropagateBackV2(Int_t nev=1) {
+   cerr<<"Propagating tracks back through the ITS...\n";
+
+   TFile *in=TFile::Open("AliITStracksV2.root");
+   if (!in->IsOpen()) {cerr<<"Can't open AliITStracksV2.root !\n"; return 1;}
+
+   TFile *out=TFile::Open("AliTPCtracks.root","update");
+   if (!out->IsOpen()) {
+      cerr<<"Can't open AliTPCtracks.root !\n"; return 2;
+   }
+   TFile *file=TFile::Open("AliITSclustersV2.root");
+   if (!file->IsOpen()) {
+      cerr<<"Can't open AliITSclustersV2.root !\n";return 3;
+   }
+   AliITSgeom *geom=(AliITSgeom*)file->Get("AliITSgeom");
+
+   Int_t rc=0;
+   TStopwatch timer;
+   AliITStrackerV2 tracker(geom);
+   for (Int_t i=0; i<nev; i++) {
+     cerr<<"Processing event number : "<<i<<endl;
+     tracker.SetEventNumber(i);
+     rc=tracker.PropagateBack(in,out);
+   }
+   timer.Stop(); timer.Print();
+
+   file->Close();
+   in->Close();
+   out->Close();
+
+   return rc;
+}
diff --git a/ITS/oldmacros/AliITSReadPlotData.C b/ITS/oldmacros/AliITSReadPlotData.C
new file mode 100644 (file)
index 0000000..9c9e195
--- /dev/null
@@ -0,0 +1,478 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include<iostream.h>
+#include<TROOT.h>
+#include<TArrayI.h>
+#include<TBranch.h>
+#include<TCanvas.h>
+#include<TClassTable.h>
+#include<TClonesArray.h>
+#include<TFile.h>
+#include<TH1.h>
+#include<TH2.h>
+#include<TLegend.h>
+#include<TObject.h>
+#include<TObjArray.h>
+#include<TTree.h>
+#include <AliRun.h>
+#include <AliITS.h>
+#include <AliITSgeom.h>
+#include <AliITSDetType.h>
+#include <AliITSRecPoint.h>
+#include <AliITSdigit.h>
+#include <AliITShit.h>
+#include <AliITSmodule.h> 
+#include <AliITSsegmentation.h>
+#include <AliITSsegmentationSPD.h> 
+#include <AliITSsegmentationSDD.h>
+#include <AliITSsegmentationSSD.h>
+#endif
+
+Int_t GetModuleHits (TObject* its, Int_t ID, Float_t*& X, Float_t*& Y, Float_t*& Z, Bool_t*& St);
+Int_t GetModuleRecPoints (TObject *its, Int_t ID, Float_t*& X, Float_t*& Z, Int_t isfastpoints);
+Int_t GetModuleDigits(TObject *its, Int_t ID, Int_t dtype, Float_t*& X, Float_t*& Z);
+
+Int_t AliITSReadPlotData( Int_t evNum = 0, char *filename = "galice.root", TString FileDigits="galice.root", TString FileRec="galice.root", Int_t isfastpoints = 0) {
+
+  /*********************************************************************
+   *                                                                   *
+   *  Macro used to read hits, digits and recpoints for a module       *
+   *  It draws 3 TH2Fs where stores the 2-D coordinates of these       *
+   *  data for a specified module (for which the macro asks when       *
+   *  it starts, with some checks to avoid wrong detector coords.      *
+   *                                                                   *
+   *  Only a little 'experimental' warning:                            *
+   *  with all the tests I have made, I wasn't able to plot the        *
+   *  digits fo the 5th layer...                                       *
+   *  I skipped this problem with an alert to beware th user, while    *
+   *  in this situation the macro plots only recpoints and hits        *
+   *                                                                   *
+   *  Author: Alberto Pulvirenti                                       *
+   *                                                                   *
+   *********************************************************************/
+
+       
+  // First of all, here are put some variable declarations
+  // that are useful in the following part:
+  Int_t nparticles; // number of particles
+  // ITS module coordinates [layer = 1, ladder = 2, det = 3] and absolute ID[0] of module [0]
+  TArrayI ID(4);
+  Int_t nmodules, dtype; // Total number of modules and module type (SSD, SPD, SDD)
+  Float_t *x = 0, *y = 0, *z = 0; // Arrays where to store read coords
+  Bool_t *St = 0; // Status of the track (hits only)
+  /*
+       // It's necessary to load the gAlice shared libs
+       // if they aren't already stored in memory...
+       if (gClassTable->GetID("AliRun") < 0) {
+       gROOT->LoadMacro("loadlibs.C");
+    loadlibs();
+       }
+    // Anyway, this macro needs to read a gAlice file, so it
+    // clears the gAlice object if there is already one in memory...
+    else {
+  */
+  if(gAlice){
+    delete gAlice;
+    gAlice = 0;
+  }
+  //   }
+
+  // Now is opened the Root input file containing Geometry, Kine and Hits
+  // by default its name must be "galice.root".
+  // When the file is opened, its contens are shown.
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+  if (!file) file = new TFile(filename);
+  file->ls();
+       
+  // Then, the macro gets the AliRun object from file.
+  // If this object is not present, an error occurs
+  // and the execution is stopped.
+  // Anyway, this operation needs some time,
+  // don't worry about an apparent absence of output and actions...
+  cout << "\nSearching in '" << filename << "' for an AliRun object ... " << flush;
+  gAlice = (AliRun*)file->Get("gAlice");
+  if (gAlice)
+    cout << "FOUND!" << endl;
+  else {
+    cout<<"NOT FOUND! The Macro can't continue!" << endl;
+    return 0;
+  }
+  if(!(FileDigits.Data() == filename)){
+    gAlice->SetTreeDFileName(FileDigits);
+  }
+  if(!(FileRec.Data() == filename)){
+    gAlice->SetTreeRFileName(FileRec);
+  }
+       
+  // Then, the macro selects the event number specified. Default is 0.
+  nparticles = gAlice->GetEvent(evNum);
+  cout << "\nNumber of particles   = " << nparticles << endl;
+  if (!nparticles) {
+    cout << "With no particles I can't do much... Goodbye!" << endl;
+    return 0;
+  }
+       
+  // The next step is filling the ITS from the AliRun object.
+  AliITS *ITS = (AliITS*)gAlice->GetModule("ITS");
+  ITS->InitModules(-1, nmodules);
+  cout << "Number of ITS modules = " << nmodules << endl;
+  cout << "\nFilling modules (it takes a while, now)..." << flush;
+  ITS->FillModules(0, 0, nmodules, " ", " ");
+  cout << "DONE!" << endl;
+  AliITSgeom *gm = ITS->GetITSgeom();
+  //   AliITSDetType *det = ITS->DetType(dtype);       
+  //   AliITSsegmentation *seg = det->GetSegmentationModel();  
+       
+  for(;;) {
+
+    // Input phase.
+    // The macro asks if the user wants to put a detector ID[0]
+    // or prefers to input layer, ladder and detector.
+    for (Int_t i = 0; i < 4; i++) ID[i] = -1;
+    Int_t answ;
+    do {
+      cout << "\nSelection modes:" << endl;
+      cout << "1) by layer - ladder - detector numbers" << endl;
+      cout << "2) by unique detector ID" << endl;
+      cout << "0) exit macro" << endl;
+      cout << "\nEnter your choice: ";
+      cin >> answ;
+    } while (answ < 0 || answ > 2);
+    switch (answ) {
+    case 0:
+      // input = 0 ---> EXIT
+      return 0;
+      break;
+    case 1:
+      // input = 1 ---> SELECTION BY COORD
+      do {
+        cout << "\nLayer number [1-6, 0 to exit]: ";
+        cin >> ID[1];
+        if (!ID[1]) return 0;
+      } while (ID[1] < 0 || ID[1] > 6);
+                               
+      // Detector type: 0 --> SPD, 1 --> SDD, 2 --> SSD.
+      // Layer 1,2 --> 0 / Layer 3,4 --> 1 / Layer 5,6 --> 2
+      dtype = (ID[1] - 1) / 2;
+                               
+      // Once fixed the layer number, the macro calculates the max number
+      // for ladder and detector from geometry, and accepts only suitable values.
+      do {
+        ID[2] = gm->GetNladders(ID[1]);
+        cout << "Ladder number [1-" << ID[2] << ", 0 to exit]: ";
+        cin >> ID[2];
+        if (!ID[2]) return 0;
+      } while (ID[2] < 0 || ID[2] > gm->GetNladders(ID[1]));
+      do {
+        ID[3] = gm->GetNdetectors(ID[1]);
+        cout << "Detector number [1-" << ID[3] << ", 0 to exit]: ";
+        cin >> ID[3];
+        if (!ID[3]) return 0;
+      } while (ID[3] < 0 || ID[3] > gm->GetNdetectors(ID[1]));
+      break;
+    case 2:
+      // input = 2 ---> SELECTION BY ID[0]
+      do {
+        ID[0] = gm->GetIndexMax();
+        cout << "\n Detector ID number [0-" << ID[0]-1 << ", -1 to exit]: ";
+        cin >> ID[0];
+        ID[0]++;
+        if (ID[0] == -1) return 0;
+      } while (ID[0] < 0 || ID[0] > gm->GetIndexMax());
+      break;
+    };
+
+    if (ID[0] == -1)
+      // If ID[0] = -1 the chioce was by coords, so it's necessary to assign the ID:
+      ID[0] = gm->GetModuleIndex(ID[1], ID[2], ID[3]);
+    else {     
+      // Else we must get the layer, ladder and detector by the ID;
+      // ...but first we must remember that the ID goest from 0 to NModules - 1!
+      ID[0]--;
+      ID[1] = ITS->GetModule(ID[0])->GetLayer();
+      ID[2] = ITS->GetModule(ID[0])->GetLadder();
+      ID[3] = ITS->GetModule(ID[0])->GetDet();
+    }
+                               
+    // Defines the histograms inside the `for' cycle, so they are destroyed at the end
+    // of every read sequqnce, in order to mek another withour segmentation faults
+    Text_t msg[250];
+    Float_t xm = 0.0, zm = 0.0;
+    dtype = (ID[1] - 1) / 2;
+    switch (dtype) {
+    case 0: xm = 1.5; zm = 7.0; break;
+    case 1: xm = 7.5; zm = 8.0; break;
+    case 2: xm = 7.5; zm = 4.5; break;
+    }
+    sprintf(msg, "Module index=%d lay=%d lad=%d det=%d", ID[0], ID[1], ID[2], ID[3]);
+    TH2F *hhits = new TH2F("hhits", msg, 500, -xm, xm, 500, -zm, zm);     // Histogram of hits
+    TH2F *hrecs = new TH2F("hrecs", msg, 500, -xm, xm, 500, -zm, zm);     // Histogram of recpoints
+    TH2F *hdigits = new TH2F("hdigits", msg, 500, -xm, xm, 500, -zm, zm); // Histogram of digits
+               
+    cout << endl;
+               
+    // Reads hits...
+    Int_t tmpint = ID[0];
+    Int_t hits = GetModuleHits(ITS, tmpint, x, y, z, St);
+    if (!hits) {
+      cout << "No hits in module!" << endl;
+      continue;
+    }
+    else
+      cout << "Hits scanned..." << endl;
+    for (Int_t i = 0; i < hits; i++) if (!St[i]) hhits->Fill(x[i], z[i]);
+               
+    // Reads recpoints...
+    Int_t recs = GetModuleRecPoints(ITS, ID[0], x, z, isfastpoints);
+    if (!recs) {
+      cout << "No recpoints in module!" << endl;
+      continue;
+    }
+    else
+      cout << "Recpoints scanned..." << endl;
+    for (Int_t i = 0; i < recs; i++) hrecs->Fill(x[i], z[i]);
+               
+    // Reads digits...
+    Int_t digits = GetModuleDigits(ITS, ID[0], dtype, x, z);
+    if (!digits) {
+      cout << "No digits in module!" << endl;
+      //continue;
+    }
+    else
+      cout << "Digits scanned..." << endl;
+    for (Int_t i = 0; i < digits; i++) hdigits->Fill(x[i], z[i]);
+
+    // Draws read data...
+    // HITS -------> red (2) crosses.
+    // DIGITS -----> green (8) boxes.
+    // REC-POINTS -> blue (4) St. Andrew's crosses.
+
+    TCanvas *viewer = new TCanvas("viewer", "Module viewer canvas", 0, 0, 800, 800);
+    viewer->cd();
+               
+    hdigits->SetMarkerStyle(25);
+    hdigits->SetMarkerColor(8);
+    hdigits->SetMarkerSize(2);
+    hdigits->SetStats(kFALSE);
+    hdigits->SetXTitle("Local X (cm)");
+    hdigits->SetYTitle("Local Z (cm)");
+    hdigits->Draw();
+               
+    hhits->SetMarkerStyle(5);
+    hhits->SetMarkerColor(2);
+    hhits->SetMarkerSize(3);
+    hhits->SetStats(kFALSE);
+    hhits->Draw("same");
+       
+    hrecs->SetMarkerStyle(2);
+    hrecs->SetMarkerColor(4);
+    hrecs->SetMarkerSize(3);
+    hrecs->SetStats(kFALSE);
+    hrecs->Draw("same");
+               
+    TLegend *legend = new TLegend(0.7, 0.8, 0.9, 0.9);
+    legend->SetMargin(0.2);
+    legend->AddEntry(hhits, "hits","P");
+    legend->AddEntry(hrecs, "recpoints","P");
+    legend->AddEntry(hdigits, "digits","P");
+    legend->SetTextSizePixels(14);
+    legend->Draw();
+               
+    viewer->Update();
+               
+               
+    Text_t fname[250],ans;
+    cout << "Do you want to save the current canvas on a file (y/n) ? ";
+    cin >> ans;
+    if(ans == 'y' || ans == 'Y') {
+      cout << "Enter filename: ";
+      cin >> fname;
+      TString *control = new TString(fname);
+      Bool_t ok=control->Contains(".C") || control->Contains(".root") || control->Contains(".ps") || control->Contains(".eps") || control->Contains(".gif");
+      if(!ok){ 
+        cout << "File extension is not recognized. The canvas will be saved as Postscript file";
+        strcat(fname,".ps");
+      }  
+      viewer->SaveAs(fname);
+    }
+  }
+       
+  cout << "Done. Goodbye" << endl;
+  return 0;
+}
+
+Int_t GetModuleHits(TObject* its, Int_t ID, Float_t*& X, Float_t*& Y, Float_t*& Z, Bool_t*& St) {      
+  // First of all, the macro selects the specified module,
+  // then gets the array of hits in it and their number.
+  AliITS *ITS = (AliITS*) its;
+  AliITSmodule *module = ITS->GetModule(ID);
+  TObjArray *hits_array = module->GetHits();
+  Int_t hits_num = hits_array->GetEntriesFast();
+       
+       // Now, if this count returns 0, there's nothing to do,
+       // while, if it doesn't, the first thing to do is dimensioning
+       // the coordinate arrays, and then the loop can start.
+  if (!hits_num)
+    return 0;
+  else {
+    if (X) delete [] X;        
+    if (Y) delete [] Y;
+    if (Z) delete [] Z;
+    if (St) delete [] St;
+    X = new Float_t[hits_num];
+    Y = new Float_t[hits_num];
+    Z = new Float_t[hits_num];
+    St = new Bool_t[hits_num];
+  }
+
+  for (Int_t j = 0; j < hits_num; j++) {
+    AliITShit *hit = (AliITShit*) hits_array->At(j);
+    X[j]  = hit->GetXL();
+    Y[j]  = hit->GetYL();
+    Z[j]  = hit->GetZL();
+    St[j] = hit->StatusEntering();
+  }
+  return hits_num;
+}
+
+Int_t GetModuleRecPoints (TObject *its, Int_t ID, Float_t*& X, Float_t*& Z, Int_t isfastpoints) {
+       
+  // First of all, the macro selects the specified module,
+  // then gets the array of recpoints in it and their number.
+  AliITS *ITS = (AliITS*) its;
+  TTree *TR = gAlice->TreeR();
+  ITS->ResetRecPoints();
+       
+  TClonesArray *recs_array = ITS->RecPoints();
+
+  TBranch *branch = 0;
+  if(TR && recs_array){
+    if(isfastpoints==1){
+      branch = gAlice->TreeR()->GetBranch("ITSRecPointsF");
+      cout<<"using fast points\n";
+    }
+    else {
+      branch = gAlice->TreeR()->GetBranch("ITSRecPoints");
+    }
+    if(branch)branch->SetAddress(&recs_array);
+  }
+        
+  branch->GetEvent(ID);
+
+  Int_t recs_num = recs_array->GetEntries();
+
+
+  // Now, if this count returns 0, there's nothing to do,
+  // while, if it doesn't, the first thing to do is dimensioning
+  // the coordinate and energy loss arrays, and then the loop can start.
+  if (!recs_num)
+    return 0;
+  else {
+    if (X) delete [] X;        
+    if (Z) delete [] Z;
+    X = new Float_t[recs_num];
+    Z = new Float_t[recs_num];
+  }
+  for(Int_t j = 0; j < recs_num; j++) {
+    AliITSRecPoint *recp = (AliITSRecPoint*)recs_array->At(j);
+    X[j] = recp->GetX();
+    Z[j] = recp->GetZ();
+  }
+  return recs_num;     
+}
+
+Int_t GetModuleDigits(TObject *its, Int_t ID, Int_t dtype, Float_t*& X, Float_t*& Z) {
+
+  // First of all, the macro selects the specified module,
+  // then gets the array of recpoints in it and their number.
+  AliITSdigit *digit;
+  AliITS *ITS = (AliITS*)its;
+  TTree *TD = gAlice->TreeD();
+  ITS->ResetDigits();
+  TD->GetEvent(ID);
+  TClonesArray *digits_array = ITS->DigitsAddress(dtype);
+  AliITSgeom *gm = ITS->GetITSgeom();  
+  AliITSDetType *det = ITS->DetType(dtype);    
+  AliITSsegmentation *seg = det->GetSegmentationModel();       
+  TArrayI ssdone(5000);  // used to match p and n side digits of strips
+  TArrayI pair(5000);    // as above   
+  Int_t digits_num = digits_array->GetEntries();
+  // Now, if this count returns 0, there's nothing to do,
+  // while, if it doesn't, the first thing to do is dimensioning
+  // the coordinate and energy loss arrays, and then the loop can start.
+
+  if(dtype==2){
+    Int_t layer, ladder, detec;
+    gm->GetModuleId(ID,layer,ladder,detec);
+    seg->SetLayer(layer);
+  }
+
+  if (!digits_num)
+    return 0;
+  else {
+    cout << "Digits to scan: " << digits_num << endl;
+    if (X) delete [] X;                        
+    if (Z) delete [] Z;
+    X = new Float_t[digits_num];               
+    Z = new Float_t[digits_num];
+  }
+       
+  // Get the coordinates of the module
+  if (dtype == 2) {
+    for (Int_t j = 0; j < digits_num; j++) {
+      ssdone[j] = 0;                   
+      pair[j] = 0;
+    }
+  }
+  //    AliITSdigit *digit;
+  for (Int_t j = 0; j < digits_num; j++) {
+    digit = (AliITSdigit*)digits_array->UncheckedAt(j);
+    Int_t iz=digit->fCoord1;  // cell number z
+    Int_t ix=digit->fCoord2;  // cell number x
+    // Get local coordinates of the element (microns)
+    // ******************************* PARTE CORRETTA ***************************************
+    if(dtype < 2) {
+      Float_t xx, zz; // aggiunta
+      seg->DetToLocal(ix, iz, xx, zz);
+      X[j] = xx; // aggiunta
+      Z[j] = zz; // aggiunta
+    }
+    // ******************************* FINE PARTE CORRETTA ***************************************
+    else {
+      // SSD: if iz==0 ---> N side; if iz==1 P side
+      if (ssdone[j] == 0) {
+        ssdone[j]=1;
+        pair[j]=-1;
+        Bool_t pside = (iz == 1);
+        Bool_t impaired = kTRUE;
+        Int_t pstrip = 0;
+        Int_t nstrip = 0;
+        if (pside) pstrip = ix; else nstrip = ix;
+        for (Int_t k = 0; k < digits_num; k++) {
+          if (ssdone[k] == 0 && impaired) {
+            AliITSdigitSSD *sdigit=(AliITSdigitSSD*)digits_array->UncheckedAt(k);
+            if (sdigit->fCoord1 != iz && sdigit->GetTracks()[0] == digit->GetTracks()[0]) {
+              ssdone[k]=2;
+              pair[j]=k;
+              if (pside) nstrip = sdigit->fCoord2; else pstrip = sdigit->fCoord2;
+              impaired=kFALSE;
+            }
+          }
+        }
+        if (!impaired) seg->GetPadCxz(pstrip, nstrip, X[j], Z[j]);
+        X[j] /= 10000.0;  // convert microns to cm
+        Z[j] /= 10000.0;  // convert microns to cm
+      }
+    }
+  }
+  return digits_num;
+} 
+
+
+
+
+
+
+
+
+
+
diff --git a/ITS/oldmacros/AliITSSDigits2DigitsDubna.C b/ITS/oldmacros/AliITSSDigits2DigitsDubna.C
new file mode 100644 (file)
index 0000000..9961c66
--- /dev/null
@@ -0,0 +1,117 @@
+TFile* AccessFile(TString inFile="galice.root", TString acctype="R");
+void writeAR(TFile * fin, TFile *fou);
+void AliITSSD2D(TString inFile, TString outFile);
+
+void AliITSSDigits2DigitsDubna(TString inFile= "galiceS.root",
+                         TString outFile = "galiceD.root"){
+    // This macro takes SDigits and produces Digits. No merging is done
+    // and only one galice.root file is used. 
+    // Dynamically link some shared libs 
+    TStopwatch timer;
+
+    if(gAlice){
+       delete gAlice;
+       gAlice = 0;
+    } // end if gAlice
+    cout << "Creating digits from summable digits for the ITS..." << endl;
+    AliITSSD2D(inFile,outFile);
+    timer.Stop(); 
+    timer.Print();
+}
+//______________________________________________________________________
+void AliITSSD2D(TString inFile, TString outFile){
+    AliRunDigitizer * manager = new AliRunDigitizer(1,1);
+    char ftmp[50];
+    sprintf(ftmp,"%s",inFile.Data());
+    TFile *file0 = AccessFile(ftmp);
+    AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");      
+    if (!ITS) {
+        cerr<<"AliITSHits2DigitsDefault.C : AliITS object not found on file"
+            << endl;
+        return 3;
+    }  // end if !ITS
+    if(!(ITS->GetITSgeom())){
+        cerr << " AliITSgeom not found. Can't digitize with out it." << endl;
+        return 4;
+    } // end if
+
+    // For old files, must change SPD noise.
+    AliITSresponseSPDdubna *resp0 = new AliITSresponseSPDdubna();
+    if(ITS->DetType(0)->GetResponseModel() !=0){
+       delete ((AliITSresponse*)ITS->DetType(0)->GetResponseModel());
+       ITS->DetType(0)->ResponseModel(0);
+    } // end if
+    ITS->DetType(0)->ResponseModel(resp0);
+    AliITSsegmentationSPD *seg0 = (AliITSsegmentationSPD*)ITS->DetType(0)->
+       GetSegmentationModel();
+    AliITSsimulationSPDdubna *sim0 = new AliITSsimulationSPDdubna(seg0,resp0);
+    if(ITS->DetType(0)->GetSimulationModel() !=0){
+       delete ((AliITSsimulation*)ITS->DetType(0)->GetSimulationModel());
+       ITS->DetType(0)->SimulationModel(0);
+    } // end if
+    ITS->DetType(0)->SimulationModel(sim0);
+    manager->SetInputStream(0,ftmp);
+    if(outFile != "")manager->SetOutputFile(outFile);
+    AliITSDigitizer *dITS  = new AliITSDigitizer(manager);
+    manager->Exec("");
+    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(inFile);
+    TFile *file2 = 0;
+    if(outFile != ""){ 
+       file2 = new TFile(outFile,"UPDATE");
+       writeAR(file,file2);
+    } // end if outFile!=""
+    delete manager;
+    if(file){
+       file->Write();
+    } // end if file
+    if(file2){
+       file2->Close();
+       delete file2;
+    } // end if file2
+}
+//______________________________________________________________________
+TFile * AccessFile(TString FileName, TString acctype){
+
+  // Function used to open the input file and fetch the AliRun object
+
+  TFile *retfil = 0;
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(FileName);
+  if (file) {file->Close(); delete file; file = 0;} 
+  if(acctype.Contains("U")){
+    file = new TFile(FileName,"update");
+  }
+  if(acctype.Contains("N") && !file){
+    file = new TFile(FileName,"recreate");
+  }
+  if(!file) file = new TFile(FileName);   // default readonly
+  if (!file->IsOpen()) {
+        cerr<<"Can't open "<<FileName<<" !" << endl;
+        return retfil;
+  } 
+
+  // Get AliRun object from file or return if not on file
+  if (gAlice) {delete gAlice; gAlice = 0;}
+  gAlice = (AliRun*)file->Get("gAlice");
+  if (!gAlice) {
+        cerr << "AliRun object not found on file"<< endl;
+        return retfil;
+  } 
+  return file;
+}
+//______________________________________________________________________
+void writeAR(TFile * fin, TFile *fou) {
+    TDirectory *current = gDirectory;
+    TTree *Te;
+    TTree *TeNew;
+    AliHeader *alhe = new AliHeader();
+    Te = (TTree*)fin->Get("TE");
+    Te->SetBranchAddress("Header",&alhe);
+    Te->SetBranchStatus("*",1);
+    fou->cd();
+    TeNew = Te->CloneTree();
+    TeNew->Write(0,TObject::kOverwrite);
+    gAlice->Write(0,TObject::kOverwrite);
+    current->cd();
+    delete alhe;
+    cout<<"AliRun object written to file\n";
+}
diff --git a/ITS/oldmacros/AliITSSavePIDV2.C b/ITS/oldmacros/AliITSSavePIDV2.C
new file mode 100644 (file)
index 0000000..50aba07
--- /dev/null
@@ -0,0 +1,99 @@
+////////////////////////////////////////////////////////////
+// Test macro for AliITSPid class and tracking version V2 // 
+// Rev. 25 July 2002 v3-08-03 + Root 3.02/07 + RedHat 6.2 //
+////////////////////////////////////////////////////////////
+void
+AliITSSavePIDV2(Int_t evNumber1=0,Int_t evNumber2=0) {
+  const char *filename="AliITStracksV2.root";
+
+   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+   if (!file) file = new TFile(filename);
+
+TFile *fpid = new TFile("AliITStracksV2Pid.root","recreate");
+
+
+ Float_t factor=0;
+ if(gAlice!=0)factor=gAlice->Field()->SolenoidField();
+ if(factor==0.){
+   cout<<" ================  WARNING ====================================\n";
+   cout<<"  The default magnetic field value of 0.4 T will be used\n";
+   cout<<" ==============================================================\n";
+   factor=4.;  // Default  mag. field = 0.4T
+ }
+ else {
+   cout<<"AliITSSavePIDV2.C:  Magnetic field is "<<factor/10.<<" T\n";
+ }
+ AliKalmanTrack::SetConvConst(1000/0.299792458/factor);
+
+
+
+ AliITSPid *pid=new AliITSPid(100);
+//
+//   Loop over events 
+//
+ for (int nev=0; nev<= evNumber2; nev++) {
+   char tname[30];
+   sprintf(tname,"TreeT_ITS_%d;1",nev);
+   TTree *tracktree=(TTree*)file->Get(tname);
+   TBranch *tbranch=tracktree->GetBranch("tracks");
+       
+   Int_t nentr=tracktree->GetEntries();
+   cout<<"Found "<<nentr<<" ITS tracks in event No "<<nev<<endl;
+
+   char tpidname[30];
+   sprintf(tpidname,"TreeT%d",nev);
+   AliITStrackV2Pid pidtmp;
+   TTree itspidTree(tpidname,"Tree with PID");
+   AliITStrackV2Pid *outpid=&pidtmp;
+   itspidTree.Branch("pids","AliITStrackV2Pid",&outpid,32000,1);
+   AliITStrackV2 *iotrack=0;
+   for (Int_t i=0; i<nentr; i++) {
+      AliITStrackV2 *iotrack=new AliITStrackV2;
+       tbranch->SetAddress(&iotrack);
+       tracktree->GetEvent(i);
+              iotrack->CookdEdx();
+              Float_t  signal=iotrack->GetdEdx();
+             iotrack->PropagateTo(3.,0.0028,65.19);
+             iotrack->PropagateToVertex();
+             Double_t xk,par[5]; iotrack->GetExternalParameters(xk,par);
+             Float_t lam=TMath::ATan(par[3]);
+             Float_t pt_1=TMath::Abs(par[4]);
+             Float_t mom=0.;
+             if( (pt_1*TMath::Cos(lam))!=0. ){ mom=1./(pt_1*TMath::Cos(lam)); }else{mom=0.;};
+             Float_t phi=TMath::ASin(par[2]) + iotrack->GetAlpha();
+             if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
+             if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
+             signal=signal/35.;
+       pidtmp.fSignal=signal;
+       pidtmp.fMom=mom;
+       pidtmp.fPhi=phi;
+       pidtmp.fLam=lam;
+       pidtmp.fGlab=TMath::Abs(iotrack->GetLabel());
+       Int_t pcode=pid->GetPcode(signal,mom);
+       pidtmp.fPcode=pcode;
+       pidtmp.fWpi=pid->GetWpi();
+       pidtmp.fWk=pid->GetWk();
+       pidtmp.fWp=pid->GetWp();
+       //cout<<" pcode,sig,mom="<<pcode<<" "<<signal<<" "<<mom<<endl;
+       //cout<<" wpi,wka,wp="<<pidtmp.fWpi<<" "<<pidtmp.fWk<<" "<<pidtmp.fWp<<endl;
+       itspidTree.Fill();
+       delete iotrack;
+   }// End for i (tracks)
+   cout<<"n ev="<<nev<<endl;
+ fpid->Write();
+ }// End for nev (events)
+ file->Close();                 
+ cout<<"File AliITStracksV2Pid.root written"<<endl;
+ delete file;
+ cout<<"end of AliITSSavePIDV2.C "<<endl;
+ //return;
+///////////////////////////////////////////////////////
+}
+
+
+
+
+
+
+
+
diff --git a/ITS/oldmacros/AliITSScanPIDV1.C b/ITS/oldmacros/AliITSScanPIDV1.C
new file mode 100644 (file)
index 0000000..4810e6a
--- /dev/null
@@ -0,0 +1,138 @@
+///////////////////////////////////////////////////////////
+// Test macro for AliITStracksV1Pid.root file            //
+// JINR Dubna Jan 2002                                   //
+///////////////////////////////////////////////////////////
+void
+AliITSScanPIDV1(Int_t evNumber1=0,Int_t evNumber2=0) {
+  //................. Prepare histogramms ................
+     TH2F *qplot =  new TH2F("Qtrm","Qtrm vs Pmom",100,0,1300,100,0,13);
+     TH2F *qplotP=  new TH2F("Qtrm","Qtrm vs Pmom",100,0,1300,100,0,13); 
+     TH2F *qplotKa= new TH2F("Qtrm","Qtrm vs Pmom",100,0,1300,100,0,13);
+     TH2F *qplotPi= new TH2F("Qtrm","Qtrm vs Pmom",100,0,1300,100,0,13);
+     TH2F *qplotE=  new TH2F("Qtrm","Qtrm vs Pmom",100,0,1300,100,0,13);
+     qplotP.SetMarkerStyle(8); qplotP.SetMarkerColor(kBlack); qplotP.SetMarkerSize(.3);
+     qplotKa.SetMarkerStyle(8); qplotKa.SetMarkerColor(kRed); qplotKa.SetMarkerSize(.3);
+     qplotPi.SetMarkerStyle(8); qplotPi.SetMarkerColor(kBlue); qplotPi.SetMarkerSize(.3);
+     qplotE.SetMarkerStyle(8); qplotE.SetMarkerColor(kGreen); qplotE.SetMarkerSize(.3);
+  //......................................................
+  TH1F *signal_mip = new TH1F("signal_mip","Signal (mips) for track",100,0.,15.);
+
+//*****************************************************************************************************************************************
+
+  const char *filename="itstracks.root";
+  
+  ///////////////// Dynamically link some shared libs ////////////////////////////////
+  
+  if (gClassTable->GetID("AliRun") < 0) {
+    gROOT->LoadMacro("loadlibs.C");
+    loadlibs();
+  } else {
+    delete gAlice;
+    gAlice=0;
+  }
+
+// Connect the Root Galice file containing Geometry, Kine and Hits
+   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+   if (!file) file = new TFile(filename);
+
+//
+//   Loop over events 
+//
+   char tname[30];
+   for (int nev=evNumber1; nev<= evNumber2; nev++) {
+
+  // for (int nev=0; nev<= evNumber2; nev++) {
+
+   sprintf(tname,"TreeT%d",nev);
+   TTree *tracktree=(TTree*)file->Get(tname);
+   TBranch *tbranch=tracktree->GetBranch("ITStracks");
+   cout<<" nev = "<<nev<<"\n";
+       //cout<<" open the file \n"; 
+       
+   Int_t nentr=tracktree->GetEntries();
+
+   TObjArray tarray(nentr);
+  // AliITSIOTrack *iotrack=0;
+   printf("nentr %d\n",nentr);
+       
+   for (Int_t i=0; i<nentr; i++) {
+      AliITSIOTrack *iotrack=new AliITSIOTrack;
+      // tarray.AddAt(new AliITSIOTrack,i);
+      // iotrack=(AliITSiotrack*)tarray.UncheckedAt(i);
+       tbranch->SetAddress(&iotrack);
+       tracktree->GetEvent(i);
+       tarray.AddLast(iotrack);
+   }
+   //file->Close();             
+       
+         AliITSIOTrack *iotrack;
+   for (Int_t i=0; i<nentr; i++) {
+         AliITSIOTrack *iotrack=new AliITSIOTrack;     
+        iotrack=(AliITSIOTrack*)tarray.UncheckedAt(i);
+        if(!iotrack) continue;
+         
+        Float_t Px=iotrack->GetPx();
+        Float_t Py=iotrack->GetPy();
+        Float_t Pz=iotrack->GetPz();
+         
+        Float_t momentum = TMath::Sqrt(Px*Px+Py*Py+Pz*Pz);
+        Float_t dEdx = iotrack->GetdEdx();
+        Int_t pcode = 211;//iotrack->GetPDG();
+
+        signal_mip->Fill(dEdx);
+        
+         if(pcode == 2212) qplotP.Fill(1000*momentum,dEdx);
+        if(pcode ==  321) qplotKa.Fill(1000*momentum,dEdx);
+        if(pcode ==  211) qplotPi.Fill(1000*momentum,dEdx);
+        if(pcode ==   11) qplotE.Fill(1000*momentum,dEdx);
+        
+    delete iotrack;             
+   }  
+
+   }   // event loop 
+   file->Close();   
+
+//*****************************************************************************************************************************************
+  
+  //...................... Draw histogramms .................
+   TCanvas *c1 = new TCanvas("PID_test","Scan PID ",200,10,900,700);
+   c1->Divide(2,1);
+  //.........................................................
+   c1->cd(1); gPad->SetFillColor(33);
+   signal_mip->Draw();
+
+   c1->cd(2); //gPad->SetFillColor(33);
+   qplot->Draw();
+   qplotP.Draw("same"); qplotKa.Draw("same"); qplotPi.Draw("same"); qplotE.Draw("same");
+   AliITSPid *pid =new AliITSPid(1000);
+    c1->Range(0,0,1300,10);
+    gStyle->SetLineColor(kRed);
+    gStyle->SetLineWidth(2);
+       TLine *lj[3],*lk[3]; 
+       for(Int_t j=0;j<3;j++){
+               Float_t x1,x2,y1,y2,xx1,xx2,yy1,yy2;
+               x1=pid->cut[j+1][0]; x2=pid->cut[j+2][0];
+               y1=y2=pid->cut[j+2][2];
+           lj[j]=new TLine(1000*x1,y1,1000*x2,y2); lj[j]->Draw();
+           if(j==0){yy1=10.;}else{yy1=lj[j-1]->GetY1();}
+           yy2=lj[j]->GetY1();
+           xx1=xx2=x1;
+           lk[j]=new TLine(1000*xx1,yy1,1000*xx2,yy2); lk[j]->Draw();
+       }
+       //Draw pions-kaons cuts.
+       TLine *mj[7],*mk[7]; 
+       for(Int_t j=0;j<7;j++){
+               Float_t x1,x2,y1,y2,xx1,xx2,yy1,yy2;
+               x1=pid->cut[j+2][0]; x2=pid->cut[j+3][0];
+               y1=y2=pid->cut[j+3][5];
+           mj[j]=new TLine(1000*x1,y1,1000*x2,y2); mj[j]->Draw();
+           if(j==0){yy1=10.;}else{yy1=mj[j-1]->GetY1();}
+           yy2=mj[j]->GetY1();
+           xx1=xx2=x1;
+           mk[j]=new TLine(1000*xx1,yy1,1000*xx2,yy2); mk[j]->Draw();
+       }
+  cout<<"End of file AliITStracksV1Pid.root "<<endl; 
+  return;
+}
+
diff --git a/ITS/oldmacros/AliITSScanPIDV2.C b/ITS/oldmacros/AliITSScanPIDV2.C
new file mode 100644 (file)
index 0000000..9e2492b
--- /dev/null
@@ -0,0 +1,76 @@
+///////////////////////////////////////////////////////////
+// Test macro for AliITStracksV2Pid.root file            //
+// JINR Dubna Jan 2002                                   //
+///////////////////////////////////////////////////////////
+void
+AliITSScanPIDV2(Int_t evNumber1=0,Int_t evNumber2=0) {
+  //................. Prepare histogramms ................
+     TH2F *qplot =  new TH2F("Qtrm","Qtrm vs Pmom",100,0,1.300,100,0,13);
+     TH2F *qplotP=  new TH2F("QtrmP","Qtrm vs Pmom",100,0,1.300,100,0,13); 
+     TH2F *qplotKa= new TH2F("QtrmKa","Qtrm vs Pmom",100,0,1.300,100,0,13);
+     TH2F *qplotPi= new TH2F("QtrmPi","Qtrm vs Pmom",100,0,1.300,100,0,13);
+     TH2F *qplotE=  new TH2F("QtrmE","Qtrm vs Pmom",100,0,1.300,100,0,13);
+     qplotP.SetMarkerStyle(8); qplotP.SetMarkerColor(kBlue); qplotP.SetMarkerSize(.3);
+     qplotKa.SetMarkerStyle(8); qplotKa.SetMarkerColor(kRed); qplotKa.SetMarkerSize(.3);
+     qplotPi.SetMarkerStyle(8); qplotPi.SetMarkerColor(kBlack); qplotPi.SetMarkerSize(.3);
+     qplotE.SetMarkerStyle(8); qplotE.SetMarkerColor(kGreen); qplotE.SetMarkerSize(.3);
+  //......................................................
+  TH1F *signal_mip = new TH1F("signal_mip","Signal (mips) for track",100,0.,15.);
+
+TFile *fpid = new TFile("AliITStracksV2Pid.root","read");
+fpid->ls();
+//
+//   Loop over events 
+//
+for (int nev=0; nev<= evNumber2; nev++) {
+  char tpidname[30];
+  sprintf(tpidname,"TreeT%d",nev);
+  TTree *tracktree=(TTree*)fpid->Get(tpidname);
+  TBranch *tbranch=tracktree->GetBranch("pids");
+       
+   Int_t nentr=tracktree->GetEntries();
+   cout<<"Found PID for "<<nentr<<" ITS V2 tracks on "<<tpidname<<endl;
+
+   AliITStrackV2Pid *iopid=0;
+for(Int_t ii=0;ii<nentr;ii++)
+  {
+      AliITStrackV2Pid *iopid=new AliITStrackV2Pid;
+      tbranch->SetAddress(&iopid);
+      tracktree->GetEvent(ii);
+
+      signal_mip->Fill(iopid->fSignal);
+
+        if(iopid->fPcode ==2212)qplotP.Fill(iopid->fMom,iopid->fSignal);
+       if(iopid->fPcode == 321)qplotKa.Fill(iopid->fMom,iopid->fSignal  );
+       if(iopid->fPcode == 211)qplotPi.Fill(iopid->fMom,iopid->fSignal  );
+       if(iopid->fPcode ==  11)qplotE.Fill(iopid->fMom,iopid->fSignal   );
+       /*
+       
+if(  (iopid->fWp<0.10)||(iopid->fWk<0.0)||(iopid->fWpi<0.0) ){
+          cout<<"PID pcode,fsignal,fmom= "<<iopid->fPcode<<","<<iopid->fSignal<<","<<iopid->fMom<<endl;
+         cout<<"wpi,wka,wp="<<iopid->fWpi<<" "<<iopid->fWk<<" "<<iopid->fWp<<endl;
+      }
+       */
+      delete iopid;
+  }// Enf for ii (tracks)
+ }// End for nev (events)
+ fpid->Close();
+  //...................... Draw histogramms .................
+   TCanvas *c1 = new TCanvas("PID_test","Scan PID ",200,10,900,700);
+   c1->Divide(2,1);
+  //.........................................................
+   c1->cd(1); gPad->SetFillColor(33);
+   signal_mip->Draw();
+
+   c1->cd(2); //gPad->SetFillColor(33);
+   qplot->Draw();
+   qplotP.Draw("same"); qplotKa.Draw("same"); qplotPi.Draw("same"); qplotE.Draw("same");
+
+   AliITSPid *pid =new AliITSPid(100);
+   fcutka.Draw("same"); fcutpr.Draw("same");
+   c1->Print("ITSPIDplot.ps");
+
+  cout<<"End of file AliITStracksV2Pid.root "<<endl; 
+  return;
+}
+
diff --git a/ITS/oldmacros/AliITSStoreFindableTracks.C b/ITS/oldmacros/AliITSStoreFindableTracks.C
new file mode 100644 (file)
index 0000000..eb3bd42
--- /dev/null
@@ -0,0 +1,8 @@
+Int_t AliITSStoreFindableTracks
+(Int_t nMinClusters = 5, const Text_t *evname = "galice", Int_t evnum = 0)
+{
+       gSystem->SetIncludePath("-I- -I$ALICE_ROOT/ITS -I$ALICE_ROOT/STEER");
+
+       gROOT->LoadMacro("$ALICE_ROOT/ITS/AliITSStoreFindableTracksCompiled.C+");
+       AliITSStoreFindableTracksCompiled(nMinClusters, evname, evnum);
+}
diff --git a/ITS/oldmacros/AliITSStoreFindableTracksCompiled.C b/ITS/oldmacros/AliITSStoreFindableTracksCompiled.C
new file mode 100644 (file)
index 0000000..74590c1
--- /dev/null
@@ -0,0 +1,185 @@
+#include <fstream.h>
+
+#include <TClassTable.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <TString.h>
+#include <TClonesArray.h>
+#include <TParticle.h>
+
+#include "AliRun.h"
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSRecPoint.h"
+
+Int_t AliITSStoreFindableTracksCompiled
+(Int_t nMinClusters = 5, const Text_t *evname = "galice", Int_t evnum = 0)
+{
+       // Make sure that ALICE objects are loaded
+       if (gAlice) {
+               delete gAlice;
+               gAlice = 0;
+       }
+       
+       // Define the names of all involved files
+       TString strEventFile(evname);
+       TString strOutputFile(evname);
+       strEventFile.Append(".root");
+       strOutputFile.Append("_tracks_");
+       strOutputFile += evnum;
+       strOutputFile.Append(".root");
+               
+       // Connect the Root Galice file containing Geometry, Kine and Hits
+       TFile *fileEvent = (TFile*)gROOT->GetListOfFiles()->FindObject(strEventFile);
+       if (!fileEvent) fileEvent = new TFile(strEventFile,"UPDATE");
+       
+       // Get AliRun object from file
+       gAlice = (AliRun*)fileEvent->Get("gAlice");
+       if (gAlice) cout << "OK, found an AliRun object in file" << endl;
+       
+       // Get ITS related objects and data
+       AliITS* ITS =(AliITS *)gAlice->GetDetector("ITS");
+       if (!ITS) {
+               cerr << "ITS object not found!" << endl;
+               return 1;
+       }
+       AliITSgeom *geometry = ITS->GetITSgeom(); 
+       if (!geometry) {
+               cerr << "ITS geometry object not found!" << endl;
+               return 2;
+       }
+       
+       // Count the number of modules per layer
+       Int_t nLadders, nDetectors, mod_min[6], mod_max[6];
+       for(Int_t i = 0; i < 6; i++) {
+               nLadders = geometry->GetNladders(i + 1);
+               nDetectors = geometry->GetNdetectors(i + 1);
+               mod_min[i] = geometry->GetModuleIndex(i + 1, 1, 1);
+               mod_max[i] = geometry->GetModuleIndex(i + 1, nLadders, nDetectors);
+       }
+       
+       // Load event and ITS recpoints
+       Int_t nParticles = gAlice->GetEvent(evnum);
+       cout << "Event number: " << evnum << endl;
+       cout << "# particles : " << nParticles <<endl;
+       if (nParticles <= 0) {
+               cerr << "Can't have <= 0 particles!" << endl;
+               return 3; 
+       }
+       AliITSRecPoint *recp = 0;
+       TClonesArray  *recPoints = ITS->RecPoints();
+       TObjArray     *particles = gAlice->Particles();
+       Int_t nTracks = gAlice->GetNtrack(); //FCA correction
+       Bool_t *hitITSLayer[6];
+       for (Int_t i = 0; i < 6; i++) {
+               hitITSLayer[i] = new Bool_t[nTracks];
+               for (Int_t j = 0; j < nTracks; j++) hitITSLayer[i][j] = kFALSE;
+       }
+       
+       // Load recpoints in event
+       TTree *TR = gAlice->TreeR();
+       if (!TR) {
+               cerr << "TreeR object not found!" << endl;
+               return 4;
+       }
+       
+       // Scan recpoints and define findable tracks
+       Int_t nModules = (Int_t)TR->GetEntries(), nPoints = 0, nEmpty = 0;
+       cout << "Found " << nModules;
+       cout << " entries in the TreeR (must be one per module!)" << endl;
+       for (Int_t layer = 1; layer <= 6; layer++) {
+               for (Int_t mod = mod_min[layer - 1]; mod <= mod_max[layer - 1]; mod++) {
+                       ITS->ResetRecPoints();
+                       TR->GetEntry(mod);
+                       nPoints = recPoints->GetEntries();
+                       if(!nPoints) {
+                               nEmpty++;
+                               continue;
+                       }
+                       for (Int_t point = 0; point < nPoints; point++) {
+                               recp = (AliITSRecPoint*)recPoints->UncheckedAt(point);
+                               for (Int_t it = 0; it < 3; it++) {
+                                       Int_t track = recp->GetLabel(it);
+                                       if(track < 0) continue;
+                                       if(track > nTracks) {
+                                               cout << "Found track index " << track;
+                                               cout << " whilw gAlice->GetNtrack() = " << nTracks << endl;
+                                               continue;
+                                       }
+                                       hitITSLayer[layer - 1][track] = kTRUE;
+                               } // loop over recpoint labels
+                       } //loop over points
+               } //loop over modules
+       } //loop over layers
+       cout << "Found " << nEmpty << " empty modules" << endl;
+       
+       // Scan the file of tracks in TPC to retrieve the findable TPC tracks
+       TString strLabelsTPC;
+       Int_t label, pdg_code, nFindablesTPC = 0;
+       Double_t dummy;
+       ifstream tpc("good_tracks_tpc");
+       while (tpc >> label >> pdg_code) {
+               for (Int_t i = 0; i < 6; i++) tpc >> dummy;
+               nFindablesTPC++;
+               strLabelsTPC.Append(Form("[%d]", label));               
+       }
+                       
+       // Define the TTree with tracks data by means of a set of variables
+       Int_t nFindablesITS = 0, nFindablesITSTPC = 0;
+       Int_t nhits, tpc_ok, mother, entry = 0;
+       Double_t vx, vy, vz;
+       Double_t px, py, pz, pt;
+
+       TTree *tree = new TTree("Tracks", "Findable tracks in ITS");
+       
+       tree->Branch("vx", &vx, "vx/D");
+       tree->Branch("vy", &vy, "vy/D");
+       tree->Branch("vz", &vz, "vz/D");
+       tree->Branch("px", &px, "px/D");
+       tree->Branch("py", &py, "py/D");
+       tree->Branch("pz", &pz, "pz/D");
+       tree->Branch("pt", &pt, "pt/D");
+       tree->Branch("label", &label, "label/I");
+       tree->Branch("entry", &entry, "entry/I");
+       tree->Branch("mother", &mother, "mother/I");
+       tree->Branch("pdg_code", &pdg_code, "pdg_code/I");
+       tree->Branch("nhits", &nhits, "nhits/I");
+       tree->Branch("tpc_ok", &tpc_ok, "tpc_ok/I");
+       
+       // Fill the tree
+       cout << endl;
+       TParticle *p = 0;
+       for (Int_t i = 0; i < nTracks; i++) {
+               nhits = 0;
+               for (Int_t j = 0; j < 6; j++) if (hitITSLayer[j][i]) nhits++;
+               if (nhits < nMinClusters) continue;
+               p = gAlice->Particle(i);
+               px = p->Px();
+               py = p->Py();
+               pz = p->Pz();
+               pt = p->Pt();
+               vx = p->Vx();
+               vy = p->Vy();
+               vz = p->Vz();
+               mother = p->GetFirstMother();
+               cout << "Track " << i << " stored\r" << flush;
+               tpc_ok = (strLabelsTPC.Contains(Form("[%d]", i)));
+               pdg_code = p->GetPdgCode();
+               label = i;
+               nFindablesITS++;
+               if (tpc_ok) nFindablesITSTPC++;
+               tree->Fill();
+               entry++;
+       }
+       
+       // Save into a file
+       TFile *fileOutput = new TFile(strOutputFile, "recreate");
+       tree->Write();
+       fileOutput->Close();
+       
+       cout << "# findable tracks in TPC     : " << nFindablesTPC << endl;
+       cout << "# findable tracks in ITS     : " << nFindablesITS << endl;
+       cout << "# findable tracks in ITS+TPC : " << nFindablesITSTPC << endl;
+       
+       return 0;
+}
diff --git a/ITS/oldmacros/AliITSTrackingV1.C b/ITS/oldmacros/AliITSTrackingV1.C
new file mode 100644 (file)
index 0000000..891153e
--- /dev/null
@@ -0,0 +1,143 @@
+#ifndef __CINT__
+#include "iostream.h"
+#endif
+
+Bool_t TPCSortTracks(Int_t event = 0)
+{
+   TFile *fileTracks   = TFile::Open("AliTPCtracks.root");
+   TFile *fileClusters = TFile::Open("AliTPCclusters.root");
+   TFile *fileEvent    = TFile::Open("galice.root");
+
+   // get TPC parameterization
+   AliTPCParam *param=(AliTPCParam *)fileEvent->Get("75x40_100x60_150x60");
+   if (!param) {
+               cerr << "(TPCSortTracks) ERROR: TPC parameters have not been found!" << endl;
+               return kFALSE;
+       }
+
+   // read and sort tracks
+   Int_t i;
+       TSortedList tracks_list;
+   AliTPCtrack *iotrack = 0;
+   TTree *tracktree = (TTree*)fileTracks->Get(Form("TreeT_TPC_%d", event));
+   Int_t nentr = (Int_t)tracktree->GetEntries();
+   for (i = 0; i < nentr; i++) {
+     iotrack = new AliTPCtrack;
+     tracktree->SetBranchAddress("tracks", &iotrack);
+     tracktree->GetEvent(i);
+     tracks_list.Add(iotrack);
+   }   
+       delete tracktree;
+       
+   // assign to each track its GEANT label
+   fileClusters->cd();
+       AliTPCtracker *tracker = new AliTPCtracker(param, event);
+   tracker->LoadInnerSectors();
+   tracker->LoadOuterSectors();
+       TListIter iter(&tracks_list);
+   for (i = 0; i < nentr; i++) {
+     iotrack = (AliTPCtrack*)iter.Next();
+         if (!iotrack) {
+            cerr << "(TPCSortTracks) WARNING: Track no. " << i << " is NULL!!!" << endl;
+                 continue;  
+         }
+     tracker->CookLabel(iotrack, 0.1);
+   }   
+   delete tracker;
+       
+   // create the new TTree of TPC tracks sorted w.r. to Pt
+       tracktree = new TTree(Form("TreeT_TPC_%d", event),"Tree with TPC tracks sorted w.r to pt");
+   tracktree->Branch("tracks", "AliTPCtrack", &iotrack, 32000, 0);
+       iter.Reset();
+   for (i = 0; i < nentr; i++) {
+     iotrack = (AliTPCtrack*)iter.Next();
+         if (!iotrack) {
+            cerr << "(TPCSortTracks) WARNING: Track no. " << i << " is NULL!!!" << endl;
+                 continue;  
+         }
+     tracktree->Fill();
+   }
+       
+       // save the new tree into new file
+       TFile *fileOutput = TFile::Open("AliTPCtracksSorted.root","recreate");
+   tracktree->Write();
+   fileOutput->Close();
+   fileEvent->Close();
+   fileClusters->Close();
+       fileTracks->Close();
+       
+       return kTRUE;
+}
+
+
+void AliITSTrackingV1(Int_t evNumber1=0,Int_t evNumber2=0, Int_t min_t=-1, Int_t max_t=0,Bool_t flagvert=1, Bool_t realmass=0, const char *filename="galice.root") {
+
+  ///////////////// Dynamically link some shared libs ////////////////////////////////
+  
+  if (gClassTable->GetID("AliRun") < 0) {
+    gROOT->LoadMacro("loadlibs.C");
+    loadlibs();
+  } else {
+    delete gAlice;
+    gAlice=0;
+  }
+  
+  // Connect the Root Galice file containing Geometry, Kine and Hits
+   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+   if (!file) file = new TFile("galice.root","UPDATE");
+   //if (!file) file = new TFile(filename);
+
+// Get AliRun object from file or create it if not on file
+   if (!gAlice) {
+      gAlice = (AliRun*)file->Get("gAlice");
+      if (gAlice) printf("AliRun object found on file\n");
+      if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
+   }
+   
+   AliKalmanTrack::SetMagneticField(gAlice->Field()->SolenoidField() / 10.0);
+
+  
+  cout << "Sorting TPC tracks w.r. to transverse momentum...";
+  Bool_t success_sorting = TPCSortTracks();
+  if (success_sorting) {
+     cout << "DONE!" << endl;
+  }
+  else {
+       cout << "Some error occurred..." << endl;
+       return 1;
+  }
+
+  AliITS* IITTSS =(AliITS *)gAlice->GetDetector("ITS");        
+  if (!IITTSS) return;                                           
+
+//
+//   Loop over events 
+//
+   Int_t Nh=0;
+   Int_t Nh1=0;
+   for (Int_t nev=0; nev<= evNumber2; nev++) {
+       AliITSTrackerV1 *ITStracker = new AliITSTrackerV1(IITTSS,nev,flagvert);
+     Int_t nparticles = gAlice->GetEvent(nev);
+     cout << "nev         " << nev <<endl;
+     cout << "nparticles  " << nparticles <<endl;
+     if (nev < evNumber1) continue;
+     if (nparticles <= 0) return;
+
+     TTree *TR=gAlice->TreeR();
+     Int_t nent=TR->GetEntries();
+     //printf("Found %d entries in the TreeR (must be one per module per event!)\n",nent);
+
+
+     TStopwatch timer;
+         
+         timer.Start();
+     ITStracker->DoTracking(nev,min_t,max_t,file,realmass);    
+     timer.Stop(); timer.Print();
+        AliITSgeom *g1 = IITTSS->GetITSgeom();
+    Int_t NumOfModules = g1->GetIndexMax();       
+         ITStracker->DelMatrix(NumOfModules);
+         delete ITStracker;  
+   }   // event loop  
+   file->Close();   
+}
+
diff --git a/ITS/oldmacros/AliITSTrackingV2.C b/ITS/oldmacros/AliITSTrackingV2.C
new file mode 100644 (file)
index 0000000..aed8dfe
--- /dev/null
@@ -0,0 +1,171 @@
+////////////////////////////////////////////////////////////////////////
+//
+// AliITSTrackingV2.C 
+//
+// date: 18.03.2003
+// author: Thomas Kuhr based on AliTPCTracking.C, AliITSFindClustersV2.C
+//         and AliITSFindTracksV2.C
+// version: 1.0
+// description: 
+//      reconstructs of tracks in ITS in the following steps:
+//         ITS cluster finding
+//         ITS track finding
+// input parameters: 
+//        Int_t nEvents      ... nr of events to process
+//        Int_t firstEventNr ... first event number (starts from 0)
+//        char* fileNameHits ... name of file with hits
+//        char* fileNameITSDigits .. name of file with ITS digits
+//        char* fileNameITSTracks .. name of file with TPC tracks
+//        char* fileNameITSClusters .. name of file with ITS clusters (output)
+//        char* fileNameITSTracks .. name of file with ITS tracks (output)
+//
+////////////////////////////////////////////////////////////////////////
+
+
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include "Riostream.h"
+#include "TFile.h"
+#include "TBenchmark.h"
+#include "AliITS.h"
+#include "AliITSgeom.h"
+#include "AliITSclustererV2.h"
+#include "AliITStrackerV2.h"
+#include "AliRun.h"
+#endif
+
+Int_t gDEBUG = 2;
+
+Int_t AliITSTrackingV2(Int_t nEvents=1, Int_t firstEvent=0,
+                      const char* fileNameHits="galice.root",
+                      const char* fileNameITSDigits="its.digits.root",
+                      const char* fileNameTPCTracks="tpc.tracks.root",
+                      const char* fileNameITSClusters="its.clusters.root",
+                      const char* fileNameITSTracks="its.tracks.root");
+            
+Int_t ITSFindClustersV2(const char* fileNameITSDigits, 
+                       const char* fileNameITSClusters, 
+                       Int_t n, Int_t first);
+Int_t ITSFindTracksV2(const char* fileNameTPCTracks, 
+                     const char* fileNameITSClusters, 
+                     const char* fileNameITSTracks, 
+                     Int_t n, Int_t first);
+
+////////////////////////////////////////////////////////////////////////
+Int_t AliITSTrackingV2(Int_t nEvents, Int_t firstEvent,
+                      const char* fileNameHits,
+                      const char* fileNameITSDigits,
+                      const char* fileNameTPCTracks,
+                      const char* fileNameITSClusters,
+                      const char* fileNameITSTracks) {
+
+  AliTracker::SetFieldFactor(fileNameHits,kFALSE);
+
+// find clusters
+
+  if (fileNameITSDigits && fileNameITSClusters) {
+    if(ITSFindClustersV2(fileNameITSDigits,fileNameITSClusters,nEvents,firstEvent)) {
+      cerr << "ITS clustering failed \n";
+      return 1;
+    }
+  }
+
+// find tracks
+
+  if (fileNameTPCTracks && fileNameITSClusters && fileNameITSTracks) {
+    if(ITSFindTracksV2(fileNameTPCTracks,fileNameITSClusters,fileNameITSTracks,nEvents,firstEvent)) {
+      cerr << "ITS tracking failed \n";
+      return 2;
+    }
+  }
+
+  return 0;
+}
+////////////////////////////////////////////////////////////////////////
+Int_t ITSFindClustersV2(const char* fileNameITSDigits, const char* fileNameITSClusters, Int_t nEvents, Int_t firstEvent) {
+//
+// create ITS clusters, store them in the file fileNameITSClusters
+// gAlice object must be in memory
+
+  const char *name="ITSFindClustersV2";
+  if (gDEBUG>1) cout<<name<<" starts...\n";
+  if (gDEBUG>1) gBenchmark->Start(name);
+
+  TFile *in =TFile::Open(fileNameITSDigits);
+  if (!in->IsOpen()) {
+    cerr<<"Can't open file "<<fileNameITSDigits<<endl; 
+    return 1;
+  }
+  gAlice->SetTreeDFileName(fileNameITSDigits);
+  TFile *out=TFile::Open(fileNameITSClusters,"recreate");
+  if (!out->IsOpen()) {
+    cerr<<"Can't open file "<<fileNameITSClusters<<endl; 
+    return 1;
+  }
+
+  AliITS *ITS  = (AliITS*)gAlice->GetModule("ITS");
+  if (!ITS) { cerr<<"Can't find the ITS !\n"; return 3; }
+
+  AliITSgeom *geom=ITS->GetITSgeom();
+  out->cd();
+  geom->Write();
+  gROOT->cd();
+
+  AliITSclustererV2 clusterer(geom);
+  for (Int_t iEvent = firstEvent; iEvent<firstEvent+nEvents; iEvent++){
+    cout<<"ITSFindClusters: processing event "<<iEvent<<endl;
+    gAlice->GetEvent(iEvent);
+    in->cd();
+    clusterer.SetEvent(iEvent);
+    clusterer.Digits2Clusters(in,out); 
+  }
+
+  out->Close();
+  in->Close();
+
+  if (gDEBUG>1) gBenchmark->Show(name);
+  return 0;
+}
+////////////////////////////////////////////////////////////////////////
+Int_t ITSFindTracksV2(const char* fileNameTPCTracks, 
+                     const char* fileNameITSClusters, 
+                     const char* fileNameITSTracks, 
+                     Int_t nEvents, Int_t first) {
+  Int_t rc=0;
+  const char *name="ITSFindTracksV2";
+  if (gDEBUG>1) cout<<name<<" starts...\n";
+  if (gDEBUG>1) gBenchmark->Start(name);
+
+  TFile *out=TFile::Open(fileNameITSTracks,"recreate");
+  if (!out->IsOpen()) {
+    cerr<<"Can't open file "<<fileNameITSTracks<<endl; 
+    return 1;
+  }
+  TFile *in =TFile::Open(fileNameTPCTracks);
+  if (!in->IsOpen()) {
+    cerr<<"Can't open file "<<fileNameTPCTracks<<endl; 
+    return 1;
+  }
+  TFile *in2 =TFile::Open(fileNameITSClusters);
+  if (!in2->IsOpen()) {
+    cerr<<"Can't open file "<<fileNameITSClusters<<endl; 
+    return 1;
+  }
+
+  AliITSgeom *geom=(AliITSgeom*)in2->Get("AliITSgeom");
+  if (!geom) { cerr<<"can't get ITS geometry !\n"; return 1;}
+
+  AliITStrackerV2 tracker(geom);
+  for (Int_t iEvent=first;iEvent<first+nEvents;iEvent++){
+    cout<<"ITSFindTracks -- event "<<iEvent<<endl;
+    tracker.SetEventNumber(iEvent);
+    rc=tracker.Clusters2Tracks(in,out);
+  }
+
+  in->Close();
+  in2->Close();
+  out->Close();
+
+  if (gDEBUG>1) gBenchmark->Show(name);
+  return rc;
+}
+
diff --git a/ITS/oldmacros/AliITSVertexerTracksTest.C b/ITS/oldmacros/AliITSVertexerTracksTest.C
new file mode 100644 (file)
index 0000000..5ca857e
--- /dev/null
@@ -0,0 +1,79 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+//-- --- standard headers------------- 
+#include <Riostream.h>
+//--------Root headers ---------------
+#include <TSystem.h>
+#include <TFile.h>
+#include <TStopwatch.h>
+#include <TObject.h>
+#include <TTree.h>
+//----- AliRoot headers ---------------
+#include "alles.h"
+#include "AliRun.h"
+#include "AliMagF.h"
+#include "AliKalmanTrack.h"
+#include "AliITSVertex.h"
+#include "AliITSVertexer.h"
+#include "AliITSVertexerTracks.h"
+//-------------------------------------
+#endif
+void AliITSVertexerTracksTest(Int_t evFirst=0,Int_t evLast=0,
+                             const Char_t *galiceName="galice.root",
+                             const Char_t *trksName="AliITStracksV2.root",
+                             const Char_t *vtxName="AliITSVertices.root") {
+  /*******************************************************************
+   *                                                                 *
+   * Test macro for vertexing in pp using tracks.                    *
+   * Input file must contain trees with AliITStrackV2 objects.       *
+   * Output file can be the same file with the tracks                *
+   * or another file.                                                *
+   * If the file galice.root is available, B is taken from there,    *
+   * otherwise is can be set here "by hand".                         *
+   *                                                                 *
+   * Origin: A.Dainese, Padova  andrea.dainese@pd.infn.it            *
+   *******************************************************************/     
+
+  // Look for field value in galice.root
+  Double_t field = 0.4;
+  if(!gSystem->AccessPathName(galiceName,kFileExists)) {
+    TFile *galice = new TFile(galiceName);
+    gAlice = (AliRun*)galice->Get("gAlice");
+    AliMagF *fiel = (AliMagF*)gAlice->Field();
+    field=(Double_t)fiel->SolenoidField()/10.;
+    AliKalmanTrack::SetConvConst(100/0.299792458/field);
+    printf(" B = %3.1f read from gAlice and set\n",field);
+    delete gAlice;
+    gAlice = 0;
+    galice->Close();
+    delete galice;
+  } else {
+    printf(" File galice.root not found: default 0.4 T being used!\n");
+  }
+
+  // Open input and output files
+  TFile *inFile = TFile::Open(trksName);
+  TFile *outFile = TFile::Open(vtxName,"recreate");
+
+  // Create vertexer
+  AliITSVertexerTracks *vertexer = 
+    new AliITSVertexerTracks(inFile,outFile,field);
+  vertexer->SetFirstEvent(evFirst);
+  vertexer->SetLastEvent(evLast);
+  vertexer->SetDebug(0);
+  vertexer->SetUseThrustFrame(0);
+  vertexer->PrintStatus();
+  // Find vertices
+  vertexer->FindVertices();
+
+
+  delete vertexer;
+
+
+
+  inFile->Close();
+  outFile->Close();
+  delete inFile;
+  delete outFile;
+
+  return;
+}
diff --git a/ITS/oldmacros/AliITSVertexerTracksTest2.C b/ITS/oldmacros/AliITSVertexerTracksTest2.C
new file mode 100644 (file)
index 0000000..470909e
--- /dev/null
@@ -0,0 +1,108 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+//-- --- standard headers------------- 
+#include <Riostream.h>
+//--------Root headers ---------------
+#include <TSystem.h>
+#include <TFile.h>
+#include <TStopwatch.h>
+#include <TObject.h>
+#include <TTree.h>
+//----- AliRoot headers ---------------
+#include "alles.h"
+#include "AliRun.h"
+#include "AliMagF.h"
+#include "AliKalmanTrack.h"
+#include "AliITSVertex.h"
+#include "AliITSVertexer.h"
+#include "AliITSVertexerTracks.h"
+#include <AliHeader.h>
+#include <AliGenEventHeader.h>
+//-------------------------------------
+#endif
+void AliITSVertexerTracksTest2(Int_t evFirst=0,Int_t evLast=0,
+                             const Char_t *galiceName="galice.root",
+                             const Char_t *trksName="AliITStracksV2.root",
+                             const Char_t *vtxName="AliITSVertices.root") {
+  /*******************************************************************
+   *                                                                 *
+   * Test macro for vertexing in pp using tracks.                    *
+   * Input file must contain trees with AliITStrackV2 objects.       *
+   * Output file can be the same file with the tracks                *
+   * or another file.                                                *
+   * If the file galice.root is available, B is taken from there,    *
+   * otherwise is can be set here "by hand".                         *
+   *                                                                 *
+   * Origin: A.Dainese, Padova  andrea.dainese@pd.infn.it            *
+   *******************************************************************/     
+
+  // Look for field value in galice.root
+  Double_t field = 0.4;
+  Int_t kDebug = 0;
+  TFile *galice = 0;
+  if(!gSystem->AccessPathName(galiceName,kFileExists)) {
+    galice = new TFile(galiceName);
+    gAlice = (AliRun*)galice->Get("gAlice");
+    AliMagF *fiel = (AliMagF*)gAlice->Field();
+    field=(Double_t)fiel->SolenoidField()/10.;
+    AliKalmanTrack::SetConvConst(100/0.299792458/field);
+    printf(" B = %3.1f read from gAlice and set\n",field);
+  } else {
+    printf(" File galice.root not found: default 0.4 T being used!\n");
+  }
+
+  // Open input and output files
+  TFile *inFile = TFile::Open(trksName);
+  TFile *outFile = TFile::Open(vtxName,"recreate");
+
+  // Create vertexer
+  AliITSVertexerTracks *vertexer = 
+    new AliITSVertexerTracks(inFile,outFile,field,0.,0.);
+  vertexer->SetDebug(0);
+  vertexer->SetUseThrustFrame(0);
+  vertexer->PrintStatus();
+
+  AliITSVertex *vert = 0;
+  // Find vertices
+
+  for(Int_t i=evFirst; i<=evLast; i++){
+    if(i%100==0)cout<<"processing event "<<i<<endl;
+    gAlice->GetEvent(i);
+    // The true Z coord. is fetched for comparison
+    AliHeader *header = gAlice->GetHeader();
+    AliGenEventHeader* genEventHeader = header->GenEventHeader();
+    TArrayF primaryVertex(3);
+    genEventHeader->PrimaryVertex(primaryVertex);
+    vert = vertexer->FindVertexForCurrentEvent(i);
+    if(kDebug>0){
+      // Prints the results
+      cout <<"========================================================\n";
+      cout << "Event number: "<<i<<")  Z Vertex:"<<endl;
+      if(vert){
+       cout<<"FOUND: "<<vert->GetZv()<<"; "<<vert->GetZRes()<<endl;
+       cout <<" True Z position "<<primaryVertex[2]<<endl;
+       cout<<", diff= "<<(primaryVertex[2]-vert->GetZv())*10000.<<endl;
+      }
+      else {
+       cout<<"NOT FOUND "<<endl;
+      }
+    }
+    if(vert){
+      Double_t pos[3];
+      for(Int_t kk=0;kk<3;kk++)pos[kk]=(Double_t)primaryVertex[kk];
+      vert->SetTruePos(pos);
+      vertexer->WriteCurrentVertex();
+    }
+  }
+
+  delete vertexer;
+
+
+
+  inFile->Close();
+  outFile->Close();
+  delete inFile;
+  delete outFile;
+  galice->Close();
+  delete galice;
+  return;
+}
diff --git a/ITS/oldmacros/AliITSplotSPDClusters.C b/ITS/oldmacros/AliITSplotSPDClusters.C
new file mode 100644 (file)
index 0000000..1b95d72
--- /dev/null
@@ -0,0 +1,302 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+
+// Standard includes
+#include <iostream.h>
+// Root includes
+#include <TFile.h>
+#include <TObjArray.h>
+#include <TClonesArray.h>
+#include <TTree.h>
+#include <TBranch.h>
+#include <TH1D.h>
+// AliRoot includes
+#include "STEER/AliRun.h"
+#include "ITS/AliITS.h"
+#include "ITS/AliITSRawCluster.h"
+
+#endif
+
+void AliITSplotSPDClusters(const char* filename="galice_80.root"){
+    //------------------------------------------------------------------
+    // This macro will read the TreeC produced during reconstruction and
+    // plot the number and type of clusters found in the SPD.
+    // root [0] .L AliITSplotSPDClusters.C  //this loads the macro in memory
+    // root [1] AliITSplotSPDClusters();    //by default process first event
+    // or
+    // root [1] AliITSplotSPDClusters("galice2.root"); // process all events
+    //                                                    from the file
+    //                                                    galice2.root.
+    // or
+    // root [0] .x AliITSplotSPDClusters.C("galice2.root") // will all of the
+    //                                                        events from the
+    //                                                        galice2.root
+    //------------------------------------------------------------------
+
+    // delete the current gAlice object, the one from input file will be used
+    if(gAlice){
+       delete gAlice;
+       gAlice = 0;
+    }else{ // Dynamically link some shared libs
+        if(gClassTable->GetID("AliRun") < 0) {
+            gROOT->LoadMacro("loadlibs.C");
+            loadlibs();
+        } // end if
+    } // end if gAlice
+
+    // Connect the Root Galice file containing Geometry, Kine and Hits
+    TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
+    if(!file) file = new TFile(filename);
+    // Get AliRun object from file or create it if not on file
+    if(!gAlice) {
+        gAlice = (AliRun*)file->Get("gAlice");
+        if(gAlice) cout << "AliRun object found on file" << endl;
+        if(!gAlice) gAlice = new AliRun("gAlice","SPD Clulster analysis");
+    } // end if !gAlice
+    Int_t nevents = gAlice->GetEventsPerRun();
+
+    // Get pointers to the ITS Alice detectors.
+    AliITS *ITS = (AliITS*) gAlice->GetDetector("ITS");
+    if(ITS==0){
+        cout << "ITS not found. Exiting." << endl;
+        return;
+    } // end if ITS==0
+    AliITSresponseSPDdubna *resp0 = ITS->DetType(0)->GetResponseModel();
+    Float_t diffCoeff= resp0->DistanceOverVoltage();//Get Value of Diffusion Coefficient parameter d/v.
+    diffCoeff *= resp0->Temperature();
+    diffCoeff = TMath::Sqrt(diffCoeff*8.6173376e-5);// K_Boltzman/e_Coulumb
+    char    aDiffCoeff[100];            //Character array for sprintf
+
+    sprintf(aDiffCoeff,"Number of SPD Clusters in x, layer 1, DiffCoeff=%f #sqrt{cm}",
+           diffCoeff);
+    TH1D *hclx1 = new TH1D("hclx1",aDiffCoeff,15,0.5,15.5);
+    sprintf(aDiffCoeff,"Number of SPD Clusters in z, layer 1, DiffCoeff=%f #sqrt{cm}",
+           diffCoeff);
+    TH1D *hclz1 = new TH1D("hclz1",aDiffCoeff,5,0.5,5.5);
+    sprintf(aDiffCoeff,"Number of SPD Clusters in x, layer 2, DiffCoeff=%f #sqrt{cm}",
+           diffCoeff);
+    TH1D *hclx2 = new TH1D("hclx2",aDiffCoeff,15,0.5,15.5);
+    sprintf(aDiffCoeff,"Number of SPD Clusters in z, layer 2, DiffCoeff=%f #sqrt{cm}",
+           diffCoeff);
+    TH1D *hclz2 = new TH1D("hclz2",aDiffCoeff,5,0.5,5.5);
+    // Create Arrays with clusters from:  Data, Ba/Sa Model, old version of
+    // Dubna
+    Float_t dataX1[9] = {0.493, 0.493, 0.0273, 0.00617, 0.00112, 0.000257,
+                        0.000257,  0.000257, 0.000257};
+                        //clusters from the data in the x-dimension with
+                        //standard paramaters
+    Float_t baSaX1[11] = {0.942, 0.0180, 0.0112, 0.00389, 0.00203, 0.000257,
+                         0.001, 0.000257, 0.000257, 0.000257, 0.001};
+                        //same for Ba/Sa model
+    Float_t dubnaX1[7] =  {0.889, 0.0789, 0.0151, 0.00389, 0.001, 0.001,
+                          0.001};
+                        //same for old version of Dubna model
+    Float_t dataX2[9] = {0.482, 0.482, 0.0325, 0.00791, 0.00157, 0.000275,
+                         0.000275, 0.000275, 0.000275};
+                        //clusters from the data in the x-dimension with
+                        //optimized paramaters
+    Float_t baSaX2[11] = {0.946, 0.0175, 0.0245, 0.00482, 0.001, 0.001,
+                         0.000275,  0.001, 0.000275, 0.000275, 0.001};
+                        //same for Ba/Sa model
+    Float_t dubnaX2[8] = {0.638, 0.325, 0.0275, 0.01, 0.00431, 0.000275,
+                         0.001, 0.001};
+                        //same for old version of Dubna model
+    Float_t dataZ1[4] = {0.889, 0.0624, 0.000257, 0.000257};
+                        //clusters from the data in the z-dimension with
+                        //standard paramaters
+    Float_t baSaZ1[3] = {1., 0.0160, 0.000257}; //same for Ba/Sa model
+    Float_t dubnaZ1[3] = {0.889, 0.0180, 0.000257};
+                        //same for old version of Dubna model
+    Float_t dataZ2[4] = {0.889, 0.0624, 0.000257, 0.000257};
+                        //clusters from the data in the z-dimension with
+                        //optimized paramaters
+    Float_t baSaZ2[4] = {1., 0.0160, 0.00215, 0.000257}; //same for Ba/Sa model
+    Float_t dubnaZ2[3] = {0.889, 0.0412, 0.000257};
+                        //same for old version of Dubna model
+    
+    TH1F *hdataX1  = new TH1F("hdataX1","Number of SPD Clusters in x, layer 1",
+                             9, 0.5, 9.5);
+    hdataX1->SetMarkerColor(2);
+    hdataX1->SetMarkerStyle(20);
+    hdataX1->SetMarkerSize(0.7);
+    TH1F *hbaSaX1  = new TH1F("hbaSaX1", "Ba/Sa", 11, 0.5, 11.5);
+    TH1F *hdubnaX1 = new TH1F("hdubnaX1", "Old Dubna", 7, 0.5, 7.5);
+    TH1F *hdataX2  = new TH1F("hdataX2","Number of SPD Clusters in x, layer 2",
+                             9, 0.5, 9.5);
+    hdataX2->SetMarkerColor(2);
+    hdataX2->SetMarkerStyle(20);
+    hdataX2->SetMarkerSize(0.7);
+    TH1F *hbaSaX2  = new TH1F("hbaSaX2", "Ba/Sa", 11, 0.5, 11.5);
+    TH1F *hdubnaX2 = new TH1F("hdubnaX2", "Old Dubna", 8, 0.5, 8.5);
+    TH1F *hdataZ1  = new TH1F("hdataZ1","Number of SPD Clusters in z, layer 1",
+                             4, 0.5, 4.5);
+    hdataZ1->SetMarkerColor(2);
+    hdataZ1->SetMarkerStyle(20);
+    hdataZ1->SetMarkerSize(0.7);
+    TH1F *hbaSaZ1  = new TH1F("hbaSaZ1", "Ba/Sa", 3, 0.5, 3.5);
+    TH1F *hdubnaZ1 = new TH1F("hdubnaZ1", "Old Dubna", 3, 0.5, 3.5);
+    TH1F *hdataZ2  = new TH1F("hdataZ2","Number of SPD Clusters in z, layer 2",
+                             4, 0.5, 4.5);
+    hdataZ2->SetMarkerColor(2);
+    hdataZ2->SetMarkerStyle(20);
+    hdataZ2->SetMarkerSize(0.7);
+    TH1F *hbaSaZ2  = new TH1F("hbaSaZ2", "Ba/Sa", 4, 0.5, 4.5);
+    TH1F *hdubnaZ2 = new TH1F("hdubnaZ2", "Old Dubna", 3, 0.5, 3.5);
+
+    Int_t j = 0; //j will be a counter for the loops to fill the histograms
+                 //with the values for clusters for the Data, the Ba/Sa model
+                 //and the old Dubna model
+    for(j=0; j<9; j++){
+       hdataX1->SetBinContent((j+1), (Double_t)dataX1[j]);
+       hdataX1->SetBinError((j+1), 0.00001);
+    }
+    for(j=0; j<11; j++) hbaSaX1->Fill((Float_t)j +0.5, baSaX1[j]);
+    for(j=0; j<7; j++)  hdubnaX1->Fill((Float_t)j +0.5, dubnaX1[j]);
+    for(j=0; j<9; j++){
+       hdataX2->SetBinContent((j+1), (Double_t)dataX2[j]);
+       hdataX2->SetBinError((j+1), 0.00001);
+    }
+    for(j=0; j<11; j++) hbaSaX2->Fill((Float_t)j +0.5, baSaX2[j]);
+    for(j=0; j<8; j++)  hdubnaX2->Fill((Float_t)j +0.5, dubnaX2[j]);
+    for(j=0; j<4; j++){
+       hdataZ1->SetBinContent((j+1), (Double_t)dataZ1[j]);
+       hdataZ1->SetBinError((j+1), 0.00001);
+    }
+    for(j=0; j<3; j++)  hbaSaZ1->Fill((Float_t)j +0.5, baSaZ1[j]);
+    for(j=0; j<3; j++)  hdubnaZ1->Fill((Float_t)j +0.5, dubnaZ1[j]);
+    for(j=0; j<4; j++){
+       hdataZ2->SetBinContent((j+1), (Double_t)dataZ2[j]);
+       hdataZ2->SetBinError((j+1), 0.00001);
+    }
+    for(j=0; j<4; j++)  hbaSaZ2->Fill((Float_t)j +0.5, baSaZ2[j]);
+    for(j=0; j<3; j++)  hdubnaZ2->Fill((Float_t)j +0.5, dubnaZ2[j]);
+
+    TTree *tc = 0;
+    TBranch *spdBranch = 0;
+    TClonesArray *spdClustcA = new TClonesArray("AliITSRawClusterSPD",1000);
+    AliITSRawClusterSPD *spdClust = 0;
+    char tn[20];
+    Int_t evnt,i,n,nc;
+    Float_t nclx = 0,nclz = 0;
+    for(evnt=0;evnt<nevents;evnt++){
+       sprintf(tn,"TreeC%d",evnt);
+       tc = (TTree*) file->Get(tn);
+       spdBranch = tc->GetBranch("ITSClustersSPD");
+       spdBranch->SetAddress(&spdClustcA);
+       n = (Int_t) tc->GetEntries();
+       for(i=0;i<n;i++){
+           tc->GetEvent(i);
+           nc = spdClustcA->GetEntries();
+           if(nc>240) nc = 240;
+           for(j=0;j<nc;j++){
+               spdClust = (AliITSRawClusterSPD*) spdClustcA->At(j);
+               nclx = spdClust->NclX();
+               nclz = spdClust->NclZ();
+               if(spdClust->Module()<80){
+                   hclx1->Fill(nclx-0.5);
+                   hclz1->Fill(nclz-0.5);
+               }else if(spdClust->Module()<240){
+                   hclx2->Fill(nclx-0.5);
+                   hclz2->Fill(nclz-0.5);
+               } //endif
+           } // end for j
+       } // end for i
+       spdClustcA->Clear();
+       delete spdBranch; spdBranch = 0;
+       delete tc; tc = 0;
+    } // end for evnt
+
+    //Begin Process of normalizing histograms
+    Double_t integral = (Double_t) hclx1->Integral();
+    if(integral>0.0) hclx1->Scale(1./integral);
+    integral = hclz1->Integral();
+    if(integral>0.0) hclz1->Scale(1./integral);
+    integral = hclx2->Integral();
+    if(integral>0.0) hclx2->Scale(1./integral);
+    integral = hclz2->Integral();
+    if(integral>0.0) hclz2->Scale(1./integral);
+
+    hclx1->SetMinimum(0.000257);
+    hclx1->SetMaximum(1.01);
+    hclz1->SetMinimum(0.000274);
+    hclz1->SetMaximum(1.01);
+    hclx2->SetMinimum(0.000257);
+    hclx2->SetMaximum(1.01);
+    hclz2->SetMinimum(0.000257);
+    hclz2->SetMaximum(1.01);
+
+    sprintf(aDiffCoeff,"SPD Clusters with Diffusion Coefficent=%f",diffCoeff);
+    TCanvas *cSPDclusters = new TCanvas("cSPDclusters",aDiffCoeff,
+                                      400,10,600,776);
+    cSPDclusters->Divide(2, 2);
+
+    cSPDclusters->cd(1);
+    cSPDclusters->SetLogy(1);
+    hclx1->Draw("hist");
+    hdataX1->Draw("same e1p");
+    hbaSaX1->SetLineColor(4);
+    hbaSaX1->SetLineStyle(2);    
+    hbaSaX1->Draw("same");
+    hdubnaX1->SetLineColor(3);
+    hdubnaX1->SetLineStyle(4);
+    hdubnaX1->Draw("same");
+    TLegend *l1 = new TLegend(0.55,0.65,0.76,0.82);
+    l1->AddEntry(hclx1,"New simulation","l");
+    l1->AddEntry(hdataX1,"Test Beam Results","p");
+    l1->AddEntry(hbaSaX1,"Bari/Selero simulation","l");
+    l1->AddEntry(hdubnaX1,"Dubna simulation","l");
+    l1->Draw();
+
+    cSPDclusters->cd(2);
+    cSPDclusters->SetLogy(1);
+    hclz1->Draw("hist");
+    hdataZ1->Draw("same e1p");
+    hbaSaZ2->SetLineColor(4);
+    hbaSaZ1->SetLineStyle(2);
+    hbaSaZ1->Draw("same");
+    hdubnaZ1->SetLineColor(3);
+    hdubnaZ1->SetLineStyle(4);
+    hdubnaZ1->Draw("same");
+    TLegend *l2 = new TLegend(0.55,0.65,0.76,0.82);
+    l2->AddEntry(hclz1,"New simulation","l");
+    l2->AddEntry(hdataZ1,"Test Beam Results","p");
+    l2->AddEntry(hbaSaZ1,"Bari/Selero simulation","l");
+    l2->AddEntry(hdubnaZ1,"Dubna simulation","l");
+    l2->Draw();
+
+    cSPDclusters->cd(3);
+    cSPDclusters->SetLogy(1);
+    hclx2->Draw("hist");
+    TLegend *l3 = new TLegend(0.55,0.65,0.76,0.82);
+    l3->AddEntry(hclx2,"New simulation","l");
+    l3->Draw();
+/*
+    hdataX2->Draw("e1p");
+    hbaSaX2->SetLineColor(4);
+    hbaSaX2->SetLineStyle(2);
+    hbaSaX2->Draw("same");
+    hdubnaX2->SetLineColor(3);
+    hdubnaX2->SetLineStyle(4);
+    hdubnaX2->Draw("same");
+*/
+    cSPDclusters->cd(4);
+    cSPDclusters->SetLogy(1);
+    hclz2->Draw("hist");
+    TLegend *l4 = new TLegend(0.55,0.65,0.76,0.82);
+    l4->AddEntry(hclz2,"New simulation","l");
+    l4->Draw();
+/*
+    hdataZ2->Draw("e1p");
+    hbaSaZ2->SetLineColor(4);
+    hbaSaZ2->SetLineStyle(2);
+    hbaSaZ2->Draw("same");
+    hdubnaZ2->SetLineColor(3);
+    hdubnaZ2->SetLineStyle(4);
+    hdubnaZ2->Draw("same");
+*/
+    cSPDclusters->Update();
+
+    if(gROOT->IsBatch()){
+       cSPDclusters->Print("SPDClusters.eps");
+    } // end if gROOT->IsBatch()
+
+}
diff --git a/ITS/oldmacros/AliITSrefitV2.C b/ITS/oldmacros/AliITSrefitV2.C
new file mode 100644 (file)
index 0000000..6b752ef
--- /dev/null
@@ -0,0 +1,43 @@
+#ifndef __CINT__
+  #include <Riostream.h>
+  #include "AliITSgeom.h"
+  #include "AliITStrackerV2.h"
+
+  #include "TFile.h"
+  #include "TStopwatch.h"
+#endif
+
+Int_t AliITSrefitV2(Int_t nev=1) {
+   cerr<<"Propagating tracks inward through the ITS...\n";
+
+   TFile *in=TFile::Open("AliTPCrefited.root");
+   if (!in->IsOpen()) {cerr<<"Can't open AliTPCrefited.root !\n"; return 1;}
+
+   TFile *out=TFile::Open("AliITStracksV2.root","update");
+   if (!out->IsOpen()) {
+      cerr<<"Can't open AliITStracksV2.root !\n"; return 2;
+   }
+   TFile *file=TFile::Open("AliITSclustersV2.root");
+   if (!file->IsOpen()) {
+      cerr<<"Can't open AliITSclustersV2.root !\n";return 3;
+   }
+   AliITSgeom *geom=(AliITSgeom*)file->Get("AliITSgeom");
+
+   Int_t rc=0;
+   TStopwatch timer;
+   AliITStrackerV2 tracker(geom);
+   for (Int_t i=0; i<nev; i++) {
+     cerr<<"Processing event number : "<<i<<endl;
+     tracker.SetEventNumber(i);
+     rc=tracker.RefitInward(in,out);
+   }
+   timer.Stop(); timer.Print();
+
+   delete geom;
+
+   file->Close();
+   in->Close();
+   out->Close();
+
+   return rc;
+}
diff --git a/ITS/oldmacros/AliITSsd2d.C b/ITS/oldmacros/AliITSsd2d.C
new file mode 100644 (file)
index 0000000..ffba047
--- /dev/null
@@ -0,0 +1,84 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+
+#include "iostream.h"
+#include "TDatime.h"
+#include "TFile.h"
+#include "TString.h"
+#include "../STEER/AliRun.h"
+#include "../STEER/AliRunDigitizer.h"
+#include "ITS/AliITSDigitizer.h"
+#include "ITS/AliITS.h"
+#include "ITS/AliITSDetType.h"
+#include "ITS/AliITSresponseSDD.h"
+#include "TStopwatch.h"
+
+Bool_t GaliceITSok();
+TFile* AccessFile(TString inFile="galice.root", TString acctype="R");
+void writeAR(TFile * fin, TFile *fou);
+Int_t ChangeITSDefaults(TFile *hitfile,AliITS *ITS,TString opt="");
+
+#endif
+
+//#define DEBUG
+
+Int_t AliITSsd2d(TString df="galice.root",TString sf1="galice.root",
+                TString sf2="",TString opt=""){
+    // Produce ITS Digits from SDigits, with psible merging.
+
+    // Dynamically link some shared libs
+    if (gClassTable->GetID("AliRun") < 0) {
+       gROOT->LoadMacro("loadlibs.C");
+       loadlibs();
+    } // end if
+    gROOT->LoadMacro("$(ALICE_ROOT)/ITS/AliITSstandard.C");
+
+    TFile *sfp1 = 0, *sfp2 = 0, *dfp = 0;
+    if(!GaliceITSok()){
+       // gAlice not define. Must open a file and read it in.
+       if(df.CompareTo(sf1) == 0 || df.CompareTo("") == 0) {
+           // Input file = output file
+           sfp1 = AccessFile(sf1,"U");  // input file open for update.
+       }else{ // Input file different from output file.
+           sfp1 = AccessFile(sf1,"R"); // input file open as read only
+           // open output file and create TreeR on it
+           dfp = gAlice->InitTreeFile("S",df);
+       } // end if
+    } // end if !GALICEITSOK()
+    AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS"); 
+
+    ChangeITSDefaults(sfp1,ITS,opt);
+    // write the AliRun object to the output file if different from input file.
+    if(dfp) {
+       writeAR(sfp1,dfp);
+       dfp->Close(); // Manager will open in update mode.
+    } // end if
+
+    AliRunDigitizer *manager;
+    if(sf2.CompareTo("")==0) { // do not merge.
+       manager = new AliRunDigitizer(1,1);
+       manager->SetInputStream(0,sf1);
+    }else{
+       manager = new AliRunDigitizer(2,1);
+       manager->SetInputStream(0,sf1.Data());
+       manager->SetInputStream(1,sf2.Data());
+    } // end if
+    if (df.CompareTo(sf1) !=0) {
+       manager->SetOutputFile(df);
+    } // end if
+    AliITSDigitizer *dITS = new AliITSDigitizer(manager);
+    if(opt.Contains("ROI")==0) dITS->SetByRegionOfInterestFlag(1);
+
+    TStopwatch timer;
+    timer.Start();
+    manager->Exec("all");
+    timer.Stop(); 
+    timer.Print();
+    delete manager;
+
+    if(dfp!=0){
+       cout << df << " size =" << df->GetSize() << endl;
+    }else if(sfp1!=0){
+       cout << sf1 << " size =" << sf1->GetSize() << endl;
+    } // end if sf1!=0
+    return 0;
+}
diff --git a/ITS/oldmacros/AliITStestAll0.sh b/ITS/oldmacros/AliITStestAll0.sh
new file mode 100755 (executable)
index 0000000..41baf4b
--- /dev/null
@@ -0,0 +1,29 @@
+#
+set echo on
+#
+aliroot -b -q grun.C\(2,\"$ALICE_ROOT/ITS/Config_muon.C\"\) >& aliroot_muon.log 
+#
+mv galice.root galice_muon.root
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSHits2SDigits.C\(\"galice_muon.root\"\,\"galice_muonS.root\"\) >& aliroot_muonS.log
+#
+aliroot -b -q grun.C\(2,\"$ALICE_ROOT/ITS/Config_bg.C\"\) >& aliroot_bg.log
+#
+cp galice.root galice_all.root
+mv galice.root galice_bg.root
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSHits2SDigits.C\(\"galice_bg.root\"\,\"galice_bgS.root\"\) >& aliroot_bgS.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSHits2Digits.C\(\"galice_bg.root\"\,\"galice_D.root\"\) >& aliroot_D.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSHits2Digits.C\(\"galice_all.root\"\,\"galice_all.root\"\) >& aliroot_allD.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSMerge.C\(\"galice_mD.root\",\"galice_muonS.root\",\"galice_bgS.root\"\) >& aliroot_m.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSDigits2RecPoints.C\(\"galice_mD.root\",\"galice_mR.root\"\) >& aliroot_mRec.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSDigits2RecPoints.C\(\"galice_D.root\",\"galice_R.root\"\) >& aliroot_Rec.log
+#
+aliroot -q -b $ALICE_ROOT/ITS/AliITSDigits2RecPoints.C\(\"galice_all.root\",\"galice_all.root\"\) >& aliroot_allRec.log
+#
+
diff --git a/ITS/oldmacros/AliITStestPIDV2.sh b/ITS/oldmacros/AliITStestPIDV2.sh
new file mode 100755 (executable)
index 0000000..6f3a944
--- /dev/null
@@ -0,0 +1,69 @@
+#!/bin/sh
+# 12 May 2003,Dubna
+# This script processes the following steps for n events (AliRoot v3-09-08,
+# root v3.05/04)
+# (use the parameter n):
+# - the TPC and ITS slow simulation (hits+digits+clusters),
+# - the TPC+ITS tracking,
+# - the TPC and ITS PID
+# - the AliITStracksPidV2.root is created with the following information:
+#    fGlab - track number
+#    fPcode - particle code after the TPC PID
+#    fMom   - particle momentum      (from the AliTPCtracks.root)
+#    fLam   - particle lambder angle (from the AliTPCtracks.root)
+#    fPhi   - particle phi angle     (from the AliTPCtracks.root)
+#    fSignal- TPC signal (mips)
+#    fWpi   - the PID weight for the identified pion
+#    fWk    - the PID weight for the identified kaon
+#    fWp    - the PID weight for the identified proton
+#       (the weghts are tuned for the HIJING events)
+# and the AliITSScanPIDV2.C is the example of macros to read the PID
+# information.
+# See ITSPIDplot.ps after run of this script.
+
+if [ $# = 0   ]; then nev=1; else nev=$1; fi
+if [ $nev = 0 ]; then nev=1; fi
+#
+# delete eventual old files from the last run
+echo "Start simulation for " $nev " event(s)"
+$ALICE_ROOT/ITS/AliITSDeleteOldFiles.sh
+#
+# run the hit generation
+aliroot -q -b "$ALICE_ROOT/macros/grun.C($nev)"  
+# digitize TPC
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCHits2Digits.C($nev)" 
+# TPC tracking
+ln -s galice.root digits.root
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCFindClustersMI.C($nev)" 
+aliroot -b <<EOI
+.L $ALICE_ROOT/TPC/AliTPCFindTracksMI.C
+AliTPCFindTracks($nev);
+.q
+EOI
+#
+# digitize ITS
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSHits2SDigits.C"
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSSDigits2Digits.C"
+# ITS-TPC tracking
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSFindClustersV2.C('s',$nev)"
+aliroot -b <<EOI
+TFile *inkin = TFile::Open("galice.root");
+if (!inkin->IsOpen()) cerr<<"Can't open galice.root !\n";               
+if (gAlice) {delete gAlice; gAlice=0;}
+                                  
+gAlice = (AliRun*)inkin->Get("gAlice");
+cout<<"AliRun object found on file "<<gAlice<<endl;
+cout<<"!!!! field ="<<gAlice->Field()->SolenoidField()<<endl;
+AliKalmanTrack::SetConvConst(1000/0.299792458/gAlice->Field()->SolenoidField());
+inkin->Close();
+.x $ALICE_ROOT/ITS/AliITSFindTracksV2.C($nev);
+.q
+EOI
+#
+# Do the PID procedure for the ITS. The output file AliITStrackV2Pid.root is
+# created with the information for each track:
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSSavePIDV2.C(0,$[$nev-1])"
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSScanPIDV2.C(0,$[$nev-1])"
+
+
+
diff --git a/ITS/oldmacros/AliITStestV1.sh b/ITS/oldmacros/AliITStestV1.sh
new file mode 100755 (executable)
index 0000000..e50179b
--- /dev/null
@@ -0,0 +1,31 @@
+#!/bin/sh
+
+# delete eventual old files from the last run
+./AliITSDeleteOldFiles.sh
+
+# do everything for the TPC
+aliroot -q -b "$ALICE_ROOT/TPC/AliTPCtest.C"
+
+# create summable digits for the ITS
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSHits2SDigits.C"
+
+# go from summable digits to digits for the ITS
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSSDigits2Digits.C"
+
+# create reconstructed points for the ITS
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSDigits2RecPoints.C"
+
+# do the tracking V1
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSTrackingV1.C"
+
+# prepare results of tracking V1 for comparison
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSStoreFindableTracks.C"
+
+# do ITS tracking V1 comparison
+aliroot -q -b "$ALICE_ROOT/ITS/AliITSComparisonV1.C"
+
+#
+# after all of the above you can run AliITSPlotTracksV1.C macro under 
+# aliroot to see plots of the efficiency and track parameter resolution
+# or AliITSDisplayTracksV1.C to display found tracks on top of the
+# ITS detailed geometry
diff --git a/ITS/oldmacros/ConfigSPD02.C b/ITS/oldmacros/ConfigSPD02.C
new file mode 100644 (file)
index 0000000..46a8584
--- /dev/null
@@ -0,0 +1,102 @@
+void Config(){
+    // Set Random Number seed
+    // gRandom->SetSeed(12345);
+    // libraries required by geant321
+    gSystem->Load("libgeant321");
+    new TGeant3("C++ Interface to Geant3");
+    if (!gSystem->Getenv("CONFIG_FILE")){
+        cout<<"Config.C: Creating Run Loader ..."<<endl;
+        AliRunLoader *rl = AliRunLoader::Open("galice.root",
+                      AliConfig::fgkDefaultEventFolderName,"recreate");
+        if (rl == 0x0){
+           gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
+           return;
+       } // end if rl==0x0
+        rl->SetCompressionLevel(2);
+        rl->SetNumberOfEventsPerFile(1000);
+        gAlice->SetRunLoader(rl);
+    } // end if !gSystem
+    TGeant3 *geant3 = (TGeant3 *) gMC;
+    // Set External decayer
+    AliDecayer *decayer = new AliDecayerPythia();
+    decayer->SetForceDecay(kAll);
+    decayer->Init();
+    gMC->SetExternalDecayer(decayer);
+    //=======================================================================
+    // ******* GEANT STEERING parameters FOR ALICE SIMULATION *******
+    geant3->SetTRIG(1);         //Number of events to be processed 
+    geant3->SetSWIT(4, 10);
+    geant3->SetDEBU(0, 0, 1);
+    //geant3->SetSWIT(2,2);
+    geant3->SetDCAY(1);
+    geant3->SetPAIR(1);
+    geant3->SetCOMP(1);
+    geant3->SetPHOT(1);
+    geant3->SetPFIS(0);
+    geant3->SetDRAY(0);
+    geant3->SetANNI(1);
+    geant3->SetBREM(1);
+    geant3->SetMUNU(1);
+    geant3->SetCKOV(1);
+    geant3->SetHADR(1);//Select pure GEANH (HADR 1) or GEANH/NUCRIN (HADR 3)
+    geant3->SetLOSS(2);
+    geant3->SetMULS(1);
+    geant3->SetRAYL(1);
+    geant3->SetAUTO(1);//Select automatic STMIN etc... calc. (AUTO 1) or manual (AUTO 0)
+    geant3->SetABAN(0);//Restore 3.16 behaviour for abandoned tracks
+    geant3->SetOPTI(2);//Select optimisation level for GEANT geometry searches (0,1,2)
+    geant3->SetERAN(5.e-7);
+    Float_t cut = 1.e-3;        // 1MeV cut by default
+    Float_t tofmax = 1.e10;
+    //             GAM ELEC NHAD CHAD MUON EBREM MUHAB EDEL MUDEL MUPA TOFMAX
+    geant3->SetCUTS(cut, cut, cut, cut, cut, cut, cut, cut, cut, cut,
+                    tofmax);
+    //=======================================================================
+    // ************* STEERING parameters FOR ALICE SIMULATION **************
+    // --- Specify event type to be tracked through the ALICE setup
+    // --- All positions are in cm, angles in degrees, and P and E in GeV
+    if (gSystem->Getenv("CONFIG_NPARTICLES")){
+        int     nParticles = atoi(gSystem->Getenv("CONFIG_NPARTICLES"));
+    }else{
+        int     nParticles = 1;
+    } // end if
+    //*********************************************
+    // Example for Moving Particle Gun            *
+    //*********************************************
+    AliGenBox *gener = new AliGenBox(nParticles);
+    gener->SetMomentumRange(100.,300.);
+    gener->SetPhiRange(0,0.1);
+    gener->SetThetaRange(0.0, .1);
+    gener->SetOrigin(0.,0.,-50.);
+    //vertex position
+    gener->SetSigma(0.1,0.1,0.0); //Sigma in (X,Y,Z) (cm) on IP position
+    gener->SetPart(211);                //GEANT particle type
+    gener->Init();
+    // Activate this line if you want the vertex smearing to happen
+    // track by track
+    //
+    //gener->SetVertexSmear(perTrack); 
+    // Field (L3 0.4 T)
+    rootfile->cd();
+    //gAlice->SetField(field);
+
+    Int_t   iHALL  =  0;
+    Int_t   iITS   =  1;
+    rl->CdGAFile();
+    //=================== Alice BODY parameters =============================
+    AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
+
+    if (iHALL){
+        //=================== HALL parameters ============================
+        AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
+    } // end if
+    if(iITS) {
+       //=================== ITS parameters ============================
+       AliITSvSPD02 *ITS  = new AliITSvSPD02("SPD test beam 2002");
+    }
+    return;
+}
+
+Float_t EtaToTheta(Float_t arg){
+  return (180./TMath::Pi())*2.*atan(exp(-arg));
+}
diff --git a/ITS/oldmacros/Config_bg.C b/ITS/oldmacros/Config_bg.C
new file mode 100644 (file)
index 0000000..34fb9cd
--- /dev/null
@@ -0,0 +1,260 @@
+void Config(){
+    // 7-DEC-2000 09:00
+    // Switch on Transition Radiation simulation. 6/12/00 18:00
+    // iZDC=1  7/12/00 09:00
+    // ThetaRange is (0., 180.). It was (0.28,179.72) 7/12/00 09:00
+    // Theta range given through pseudorapidity limits 22/6/2001
+
+    // Set Random Number seed
+    // gRandom->SetSeed(12345);
+
+    // libraries required by geant321
+    gSystem->Load("libgeant321");
+
+    new     TGeant3("C++ Interface to Geant3");
+
+    if (!gSystem->Getenv("CONFIG_FILE"))
+    {
+        TFile  *rootfile = new TFile("galice.root", "recreate");
+
+        rootfile->SetCompressionLevel(2);
+    }
+
+    TGeant3 *geant3 = (TGeant3 *) gMC;
+
+    //
+    // Set External decayer
+    AliDecayer *decayer = new AliDecayerPythia();
+
+    decayer->SetForceDecay(kAll);
+    decayer->Init();
+    gMC->SetExternalDecayer(decayer);
+    //
+    //
+    //=======================================================================
+    // ******* GEANT STEERING parameters FOR ALICE SIMULATION *******
+    geant3->SetTRIG(1);         //Number of events to be processed 
+    geant3->SetSWIT(4, 10);
+    geant3->SetDEBU(0, 0, 1);
+    //geant3->SetSWIT(2,2);
+    geant3->SetDCAY(1);
+    geant3->SetPAIR(1);
+    geant3->SetCOMP(1);
+    geant3->SetPHOT(1);
+    geant3->SetPFIS(0);
+    geant3->SetDRAY(0);
+    geant3->SetANNI(1);
+    geant3->SetBREM(1);
+    geant3->SetMUNU(1);
+    geant3->SetCKOV(1);
+    geant3->SetHADR(1); //Select pure GEANH (HADR 1) or GEANH/NUCRIN (HADR 3)
+    geant3->SetLOSS(2);
+    geant3->SetMULS(1);
+    geant3->SetRAYL(1);
+    geant3->SetAUTO(1); //Select automatic STMIN etc... calc. (AUTO 1) or manual (AUTO 0)
+    geant3->SetABAN(0); //Restore 3.16 behaviour for abandoned tracks
+    geant3->SetOPTI(2); //Select optimisation level for GEANT geometry searches (0,1,2)
+    geant3->SetERAN(5.e-7);
+
+    Float_t cut = 1.e-3;        // 1MeV cut by default
+    Float_t tofmax = 1.e10;
+
+    //             GAM ELEC NHAD CHAD MUON EBREM MUHAB EDEL MUDEL MUPA TOFMAX
+    geant3->SetCUTS(cut, cut, cut, cut, cut, cut, cut, cut, cut, cut,
+                    tofmax);
+    //
+    //=======================================================================
+    // ************* STEERING parameters FOR ALICE SIMULATION **************
+    // --- Specify event type to be tracked through the ALICE setup
+    // --- All positions are in cm, angles in degrees, and P and E in GeV
+    if (gSystem->Getenv("CONFIG_NPARTICLES")){
+        int     nParticles = atoi(gSystem->Getenv("CONFIG_NPARTICLES"));
+    } else{
+        int     nParticles = 500;
+    } // end if
+    AliGenHIJINGpara *gener = new AliGenHIJINGpara(nParticles);
+
+    gener->SetMomentumRange(0, 999);
+    gener->SetPhiRange(0, 360);
+    // Set pseudorapidity range from -8 to 8.
+    Float_t thmin = EtaToTheta(8);   // theta min. <---> eta max
+    Float_t thmax = EtaToTheta(-8);  // theta max. <---> eta min 
+    gener->SetThetaRange(thmin,thmax);
+    gener->SetOrigin(0, 0, 0);  //vertex position
+    gener->SetSigma(0, 0, 0);   //Sigma in (X,Y,Z) (cm) on IP position
+    gener->Init();
+    // 
+    // Activate this line if you want the vertex smearing to happen
+    // track by track
+    //
+    //gener->SetVertexSmear(perTrack); 
+
+    gAlice->SetField(-999, 2);  //Specify maximum magnetic field in Tesla (neg. ==> default field)
+
+    Int_t   iABSO  = 1;
+    Int_t   iDIPO  = 1;
+    Int_t   iFMD   = 1;
+    Int_t   iFRAME = 1;
+    Int_t   iHALL  = 1;
+    Int_t   iITS   = 1;
+    Int_t   iMAG   = 1;
+    Int_t   iMUON  = 1;
+    Int_t   iPHOS  = 1;
+    Int_t   iPIPE  = 1;
+    Int_t   iPMD   = 1;
+    Int_t   iRICH  = 1;
+    Int_t   iSHIL  = 1;
+    Int_t   iSTART = 1;
+    Int_t   iTOF   = 1;
+    Int_t   iTPC   = 1;
+    Int_t   iTRD   = 1;
+    Int_t   iZDC   = 1;
+    Int_t   iEMCAL = 1;
+
+    //=================== Alice BODY parameters =============================
+    AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
+
+    if (iMAG){
+        //=================== MAG parameters ============================
+        // --- Start with Magnet since detector layouts may be depending ---
+        // --- on the selected Magnet dimensions ---
+        AliMAG *MAG = new AliMAG("MAG", "Magnet");
+    }
+    if (iABSO){
+        //=================== ABSO parameters ============================
+        AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
+    }
+    if (iDIPO){
+        //=================== DIPO parameters ============================
+        AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
+    }
+    if (iHALL){
+        //=================== HALL parameters ============================
+        AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
+    }
+    if (iFRAME){
+        //=================== FRAME parameters ============================
+        AliFRAME *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
+    }
+    if (iSHIL){
+        //=================== SHIL parameters ============================
+        AliSHIL *SHIL = new AliSHILvF("SHIL", "Shielding");
+    }
+    if (iPIPE){
+        //=================== PIPE parameters ============================
+        AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
+    }
+    if(iITS) {
+       //=================== ITS parameters ============================
+       //
+       // As the innermost detector in ALICE, the Inner Tracking System 
+       // "impacts" on almost all other detectors. This involves the fact
+       // that the ITS geometry still has several options to be followed
+       //in parallel in order to determine the best set-up which minimizes
+       // the induced background. All the geometries available to date are
+       // described in the following. Read carefully the comments and use
+       // the default version (the only one uncommented) unless you are making
+       // comparisons and you know what you are doing. In this case just
+       // uncomment the ITS geometry you want to use and run Aliroot.
+       // Detailed geometries:
+       AliITSvPPRasymm *ITS  = new AliITSvPPRasymm("ITS","New ITS PPR detailed version with asymmetric services");
+       ITS->SetMinorVersion(2); // don't touch this parameter if you're not an ITS developer
+       ITS->SetReadDet(kFALSE); // don't touch this parameter if you're not an ITS developer
+       ITS->SetThicknessDet1(200.); // detector thickness on layer 1 must be in the range [100,300]
+       ITS->SetThicknessDet2(200.); // detector thickness on layer 2 must be in the range [100,300]
+       ITS->SetThicknessChip1(200.); // chip thickness on layer 1 must be in the range [150,300]
+       ITS->SetThicknessChip2(200.);// chip thickness on layer 2 must be in the range [150,300]
+       ITS->SetRails(1);            // 1 --> rails in ; 0 --> rails out
+       ITS->SetCoolingFluid(1);   // 1 --> water ; 0 --> freon
+       // Geant3 <-> EUCLID conversion
+       // ============================
+       // SetEUCLID is a flag to output (=1) or not to output (=0) both
+       // geometry and media to two ASCII files (called by default
+       // ITSgeometry.euc and ITSgeometry.tme) in a format understandable
+       // to the CAD system EUCLID.  The default (=0) means that you dont
+       // want to use this facility.
+       ITS->SetEUCLID(0);  
+    }
+    if (iTPC){
+        //============================ TPC parameters =========================
+        // This allows the user to specify sectors for the SLOW (TPC geometry
+       // 2) Simulator. SecAL (SecAU) <0 means that ALL lower (upper)
+        // sectors are specified, any value other than that requires at least
+       // one sector (lower or upper)to be specified!
+        // Reminder: sectors  1-24 are lower sectors  ( 1-12->z>0,13-24->z<0)
+        //           sectors 25-72 are the upper ones (25-48->z>0,49-72->z<0)
+        // SecLows - number of lower sectors specified (up to 6)
+        // SecUps - number of upper sectors specified (up to 12)
+        // Sens - sensitive strips for the Slow Simulator !!!
+        // This does NOT work if all S or L-sectors are specified, i.e.
+        // if SecAL or SecAU < 0
+        //---------------------------------------------------------------------
+        //  gROOT->LoadMacro("SetTPCParam.C");
+        //  AliTPCParam *param = SetTPCParam();
+        AliTPC *TPC = new AliTPCv2("TPC", "Default");
+        // All sectors included 
+        TPC->SetSecAL(-1);
+        TPC->SetSecAU(-1);
+    }
+    if (iTOF){
+        //=================== TOF parameters ============================
+        AliTOF *TOF = new AliTOFv2("TOF", "normal TOF");
+    }
+    if (iRICH){
+        //=================== RICH parameters ===========================
+        AliRICH *RICH = new AliRICHv3("RICH", "normal RICH");
+    }
+    if (iZDC){
+        //=================== ZDC parameters ============================
+        AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
+    }
+    if (iTRD){
+        //=================== TRD parameters ============================
+        AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
+        // Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 90% Xe +
+       // 10% CO2)
+        TRD->SetGasMix(1);
+        // With hole in front of PHOS
+        TRD->SetPHOShole();
+        // With hole in front of RICH
+        TRD->SetRICHhole();
+        // Switch on TR
+        AliTRDsim *TRDsim = TRD->CreateTR();
+    }
+    if (iFMD){
+        //=================== FMD parameters ============================
+        AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
+        FMD->SetRingsSi1(256);
+        FMD->SetRingsSi2(64);
+        FMD->SetSectorsSi1(20);
+        FMD->SetSectorsSi2(24);
+    }
+    if (iMUON){
+        //=================== MUON parameters ===========================
+        AliMUON *MUON = new AliMUONv1("MUON", "default");
+    }
+    if (iPHOS){
+       //=================== PHOS parameters ===========================
+        AliPHOS *PHOS = new AliPHOSv1("PHOS", "GPS2");
+    }
+    if (iPMD){
+        //=================== PMD parameters ============================
+        AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
+        PMD->SetPAR(1., 1., 0.8, 0.02);
+        PMD->SetIN(6., 18., -580., 27., 27.);
+        PMD->SetGEO(0.0, 0.2, 4.);
+        PMD->SetPadSize(0.8, 1.0, 1.0, 1.5);
+    }
+    if (iEMCAL){
+        //=================== EMCAL parameters ============================
+        AliEMCAL *EMCAL = new AliEMCALv1("EMCAL", "G56_2_55_19_104_14");
+    }
+    if (iSTART){
+        //=================== START parameters ============================
+        AliSTART *START = new AliSTARTv1("START", "START Detector");
+    }
+}
+//----------------------------------------------------------------------
+Float_t EtaToTheta(Float_t arg){
+  return (180./TMath::Pi())*2.*atan(exp(-arg));
+}
diff --git a/ITS/oldmacros/Config_muon.C b/ITS/oldmacros/Config_muon.C
new file mode 100644 (file)
index 0000000..126ac32
--- /dev/null
@@ -0,0 +1,259 @@
+void Config(){
+    // 7-DEC-2000 09:00
+    // Switch on Transition Radiation simulation. 6/12/00 18:00  
+    // iZDC=1  7/12/00 09:00
+    // ThetaRange is (0., 180.). It was (0.28,179.72) 7/12/00 09:00
+    // Theta range given through pseudorapidity limits 22/6/2001
+
+    // Set Random Number seed
+    // gRandom->SetSeed(12345);
+
+    // libraries required by geant321
+    gSystem->Load("libgeant321");
+
+    new     TGeant3("C++ Interface to Geant3");
+
+    if (!gSystem->Getenv("CONFIG_FILE"))
+    {
+        TFile  *rootfile = new TFile("galice.root", "recreate");
+
+        rootfile->SetCompressionLevel(2);
+    }
+
+    TGeant3 *geant3 = (TGeant3 *) gMC;
+
+    //
+    // Set External decayer
+    AliDecayer *decayer = new AliDecayerPythia();
+
+    decayer->SetForceDecay(kAll);
+    decayer->Init();
+    gMC->SetExternalDecayer(decayer);
+    //
+    //
+    //=======================================================================
+    // ******* GEANT STEERING parameters FOR ALICE SIMULATION *******
+    geant3->SetTRIG(1);         //Number of events to be processed 
+    geant3->SetSWIT(4, 10);
+    geant3->SetDEBU(0, 0, 1);
+    //geant3->SetSWIT(2,2);
+    geant3->SetDCAY(1);
+    geant3->SetPAIR(1);
+    geant3->SetCOMP(1);
+    geant3->SetPHOT(1);
+    geant3->SetPFIS(0);
+    geant3->SetDRAY(0);
+    geant3->SetANNI(1);
+    geant3->SetBREM(1);
+    geant3->SetMUNU(1);
+    geant3->SetCKOV(1);
+    geant3->SetHADR(1);         //Select pure GEANH (HADR 1) or GEANH/NUCRIN (HADR 3)
+    geant3->SetLOSS(2);
+    geant3->SetMULS(1);
+    geant3->SetRAYL(1);
+    geant3->SetAUTO(1);         //Select automatic STMIN etc... calc. (AUTO 1) or manual (AUTO 0)
+    geant3->SetABAN(0);         //Restore 3.16 behaviour for abandoned tracks
+    geant3->SetOPTI(2);         //Select optimisation level for GEANT geometry searches (0,1,2)
+    geant3->SetERAN(5.e-7);
+
+    Float_t cut = 1.e-3;        // 1MeV cut by default
+    Float_t tofmax = 1.e10;
+
+    //             GAM ELEC NHAD CHAD MUON EBREM MUHAB EDEL MUDEL MUPA TOFMAX
+    geant3->SetCUTS(cut, cut, cut, cut, cut, cut, cut, cut, cut, cut,
+                    tofmax);
+    //
+    //=======================================================================
+    // ************* STEERING parameters FOR ALICE SIMULATION **************
+    // --- Specify event type to be tracked through the ALICE setup
+    // --- All positions are in cm, angles in degrees, and P and E in GeV
+    if (gSystem->Getenv("CONFIG_NPARTICLES")){
+        int     nParticles = atoi(gSystem->Getenv("CONFIG_NPARTICLES"));
+    }else{
+        int     nParticles = 2;
+    } // end if
+//*********************************************
+// Example for Fixed Particle Gun
+//*********************************************
+     AliGenFixed *gener = new AliGenFixed(nParticles);
+     gener->SetMomentum(50);
+     gener->SetPhi(180.);
+     gener->SetTheta(95.);
+     gener->SetOrigin(0,0,0);          //vertex position
+     gener->SetPart(13);                //GEANT particle type
+    gener->Init();
+    // 
+    // Activate this line if you want the vertex smearing to happen
+    // track by track
+    //
+    //gener->SetVertexSmear(perTrack); 
+
+    gAlice->SetField(-999, 2);  //Specify maximum magnetic field in Tesla (neg. ==> default field)
+
+    Int_t   iABSO  = 1;
+    Int_t   iDIPO  = 1;
+    Int_t   iFMD   = 1;
+    Int_t   iFRAME = 1;
+    Int_t   iHALL  = 1;
+    Int_t   iITS   = 1;
+    Int_t   iMAG   = 1;
+    Int_t   iMUON  = 1;
+    Int_t   iPHOS  = 1;
+    Int_t   iPIPE  = 1;
+    Int_t   iPMD   = 1;
+    Int_t   iRICH  = 1;
+    Int_t   iSHIL  = 1;
+    Int_t   iSTART = 1;
+    Int_t   iTOF   = 1;
+    Int_t   iTPC   = 1;
+    Int_t   iTRD   = 1;
+    Int_t   iZDC   = 1;
+    Int_t   iEMCAL = 1;
+
+    //=================== Alice BODY parameters =============================
+    AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
+
+    if (iMAG){
+        //=================== MAG parameters ============================
+        // --- Start with Magnet since detector layouts may be depending ---
+        // --- on the selected Magnet dimensions ---
+        AliMAG *MAG = new AliMAG("MAG", "Magnet");
+    }
+    if (iABSO){
+        //=================== ABSO parameters ============================
+        AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
+    }
+    if (iDIPO){
+        //=================== DIPO parameters ============================
+        AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
+    }
+    if (iHALL){
+        //=================== HALL parameters ============================
+        AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
+    }
+    if (iFRAME){
+        //=================== FRAME parameters ============================
+        AliFRAME *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
+    }
+    if (iSHIL){
+        //=================== SHIL parameters ============================
+        AliSHIL *SHIL = new AliSHILvF("SHIL", "Shielding");
+    }
+    if (iPIPE){
+        //=================== PIPE parameters ============================
+        AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
+    }
+    if(iITS) {
+       //=================== ITS parameters ============================
+       //
+       // As the innermost detector in ALICE, the Inner Tracking System 
+       // "impacts" on almost all other detectors. This involves the fact
+       // that the ITS geometry still has several options to be followed
+       //in parallel in order to determine the best set-up which minimizes
+       // the induced background. All the geometries available to date are
+       // described in the following. Read carefully the comments and use
+       // the default version (the only one uncommented) unless you are making
+       // comparisons and you know what you are doing. In this case just
+       // uncomment the ITS geometry you want to use and run Aliroot.
+       // Detailed geometries:
+       AliITSvPPRasymm *ITS  = new AliITSvPPRasymm("ITS","New ITS PPR detailed version with asymmetric services");
+       ITS->SetMinorVersion(2); // don't touch this parameter if you're not an ITS developer
+       ITS->SetReadDet(kFALSE); // don't touch this parameter if you're not an ITS developer
+       ITS->SetThicknessDet1(200.); // detector thickness on layer 1 must be in the range [100,300]
+       ITS->SetThicknessDet2(200.); // detector thickness on layer 2 must be in the range [100,300]
+       ITS->SetThicknessChip1(200.); // chip thickness on layer 1 must be in the range [150,300]
+       ITS->SetThicknessChip2(200.);// chip thickness on layer 2 must be in the range [150,300]
+       ITS->SetRails(1);            // 1 --> rails in ; 0 --> rails out
+       ITS->SetCoolingFluid(1);   // 1 --> water ; 0 --> freon
+       // Geant3 <-> EUCLID conversion
+       // ============================
+       // SetEUCLID is a flag to output (=1) or not to output (=0) both
+       // geometry and media to two ASCII files (called by default
+       // ITSgeometry.euc and ITSgeometry.tme) in a format understandable
+       // to the CAD system EUCLID.  The default (=0) means that you dont
+       // want to use this facility.
+       ITS->SetEUCLID(0);  
+    }
+    if (iTPC){
+        //============================ TPC parameters =========================
+        // This allows the user to specify sectors for the SLOW (TPC geometry
+       // 2) Simulator. SecAL (SecAU) <0 means that ALL lower (upper)
+        // sectors are specified, any value other than that requires at least
+       // one sector (lower or upper)to be specified!
+        // Reminder: sectors  1-24 are lower sectors  ( 1-12->z>0,13-24->z<0)
+        //           sectors 25-72 are the upper ones (25-48->z>0,49-72->z<0)
+        // SecLows - number of lower sectors specified (up to 6)
+        // SecUps - number of upper sectors specified (up to 12)
+        // Sens - sensitive strips for the Slow Simulator !!!
+        // This does NOT work if all S or L-sectors are specified, i.e.
+        // if SecAL or SecAU < 0
+        //---------------------------------------------------------------------
+        //  gROOT->LoadMacro("SetTPCParam.C");
+        //  AliTPCParam *param = SetTPCParam();
+        AliTPC *TPC = new AliTPCv2("TPC", "Default");
+        // All sectors included 
+        TPC->SetSecAL(-1);
+        TPC->SetSecAU(-1);
+    }
+    if (iTOF){
+        //=================== TOF parameters ============================
+        AliTOF *TOF = new AliTOFv2("TOF", "normal TOF");
+    }
+    if (iRICH){
+        //=================== RICH parameters ===========================
+        AliRICH *RICH = new AliRICHv3("RICH", "normal RICH");
+    }
+    if (iZDC){
+        //=================== ZDC parameters ============================
+        AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
+    }
+    if (iTRD){
+        //=================== TRD parameters ============================
+        AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
+        // Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 90% Xe +
+       // 10% CO2)
+        TRD->SetGasMix(1);
+        // With hole in front of PHOS
+        TRD->SetPHOShole();
+        // With hole in front of RICH
+        TRD->SetRICHhole();
+        // Switch on TR
+        AliTRDsim *TRDsim = TRD->CreateTR();
+    }
+    if (iFMD){
+        //=================== FMD parameters ============================
+        AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
+        FMD->SetRingsSi1(256);
+        FMD->SetRingsSi2(64);
+        FMD->SetSectorsSi1(20);
+        FMD->SetSectorsSi2(24);
+    }
+    if (iMUON){
+        //=================== MUON parameters ===========================
+        AliMUON *MUON = new AliMUONv1("MUON", "default");
+    }
+    if (iPHOS){
+       //=================== PHOS parameters ===========================
+        AliPHOS *PHOS = new AliPHOSv1("PHOS", "GPS2");
+    }
+    if (iPMD){
+        //=================== PMD parameters ============================
+        AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
+        PMD->SetPAR(1., 1., 0.8, 0.02);
+        PMD->SetIN(6., 18., -580., 27., 27.);
+        PMD->SetGEO(0.0, 0.2, 4.);
+        PMD->SetPadSize(0.8, 1.0, 1.0, 1.5);
+    }
+    if (iEMCAL){
+        //=================== EMCAL parameters ============================
+        AliEMCAL *EMCAL = new AliEMCALv1("EMCAL", "G56_2_55_19_104_14");
+    }
+    if (iSTART){
+        //=================== START parameters ============================
+        AliSTART *START = new AliSTARTv1("START", "START Detector");
+    }
+}
+//----------------------------------------------------------------------
+Float_t EtaToTheta(Float_t arg){
+  return (180./TMath::Pi())*2.*atan(exp(-arg));
+}
diff --git a/ITS/oldmacros/DrawS01.C b/ITS/oldmacros/DrawS01.C
new file mode 100644 (file)
index 0000000..be45e53
--- /dev/null
@@ -0,0 +1,16 @@
+void DrawS01()
+{
+   gMC->Gsatt("*", "seen", -1);
+   gMC->Gsatt("alic", "seen", 0);
+   gROOT->LoadMacro("ViewS01.C");
+   gInterpreter->ProcessLine("ViewS01()");
+   gMC->Gdopt("hide", "on");
+   gMC->Gdopt("shad", "on");
+   gMC->Gsatt("*", "fill", 7);
+   gMC->SetClipBox(".");
+   gMC->SetClipBox("*", 0, 100, -1000, 1000, -1000, 1000);
+   gMC->DefaultRange();
+   gMC->Gdraw("alic", 40, 30, 0, 10, 10, .1, .1);
+   gMC->Gdhead(1111, "Inner Tracking System");
+   gMC->Gdman(16, 6, "MAN");
+}
diff --git a/ITS/oldmacros/DrawSPD.C b/ITS/oldmacros/DrawSPD.C
new file mode 100644 (file)
index 0000000..cbe4d49
--- /dev/null
@@ -0,0 +1,45 @@
+void DrawSPD(TString view="Prospective" ){
+    AliITS *ITS = (AliITS*) gAlice->GetModule("ITS");
+    if(!ITS) {
+       cout << "You must initilize AliRoot with an ITS before plotting it"<<endl;
+       return;
+    } // end if
+    Int_t version1 = ITS->GetMajorVersion();
+    Int_t version2 = ITS->GetMinorVersion();
+    gMC->Gsatt("*", "seen", -1);
+    gMC->Gsatt("*", "fill", 7);
+    gMC->Gsatt("alic", "seen", 0);
+    gROOT->LoadMacro("ViewITSSPD.C");
+    ViewITSSPD(version1,version2);
+    //gInterpreter->ProcessLine("ViewSPD()");
+    gMC->Gdopt("hide", "on");
+    gMC->Gdopt("shad", "on");
+    gMC->SetClipBox(".");
+    if(view.Contains("Prospective")){
+       //SetClipBox(volume name,xmin,xmax,ymin,ymax,zmin,zmax)
+       gMC->SetClipBox("*",0.0,100.0,-1000.0,1000.0,0.0,1000.0);
+       gMC->DefaultRange();
+       //Gdraw(volume name,theta,phi,psi,axis originX,axis originY,scaleX,scaleY) angles in degrees
+       gMC->Gdraw("alic", 40.0, 30.0, 0.0, 11, 10, .4, .4);
+       gMC->Gdhead(1111, "Inner Tracking System");
+       //gMC->Gdman(16,6,"MAN");
+    } // end if
+    if(view.Contains("EndView")){
+       //SetClipBox(name,xmin,xmax,ymin,ymax,zmin,zmax)
+       gMC->SetClipBox("*",-100.0,100.0,-1000.0,1000.0,1.0,1000.0);
+       gMC->DefaultRange();
+       //Gdraw(volume name,theta,phi,psi,axis originX,axis originY,scaleX,scaleY) angles in degrees
+       gMC->Gdraw("alic",0.0,0.0,0.0,10.0,10.0,1.0,1.0);
+       gMC->Gdhead(1111, "Inner Tracking System");
+       //gMC->Gdman(16, 6, "MAN");
+    } // end if
+    if(view.Contains("SideView")){
+       //SetClipBox(name,xmin,xmax,ymin,ymax,zmin,zmax)
+       gMC->SetClipBox("*",0.0,100.0,-1000.0,1000.0,0.0,1000.0);
+       gMC->DefaultRange();
+       //Gdraw(volume name,theta,phi,psi,axis originX,axis originY,scaleX,scaleY) angles in degrees
+       gMC->Gdraw("alic",90.0,0.0,0.0,10.0,10.0,0.35,0.35);
+       gMC->Gdhead(1111, "Inner Tracking System");
+       //gMC->Gdman(16, 6, "MAN");
+    } // end if
+}
diff --git a/ITS/oldmacros/MergeV1.C b/ITS/oldmacros/MergeV1.C
new file mode 100644 (file)
index 0000000..6c1f80e
--- /dev/null
@@ -0,0 +1,152 @@
+////////////////////////////////////////////////////////////////////////
+//
+// name: MergeV1
+// date: 11.4.2002
+// last update: 11.4.2002
+// author: Jiri Chudoba
+// version: 1.0
+//
+// description: 
+//       creates digits from sdigits for several detectors
+//       stores sdigits in separate file (or in the source file
+//       with sdigits). Stores gAlice object and copies TE to the
+//       file with digits
+//       ITS region of Interest is set
+//       test
+//
+// input:
+//       TString fileNameSDigits ... input file with sdigits
+//       TString fileNameDigits ... output file with digits
+//       Int_t nEvents  ... how many events to process
+//       Int_t ITS, TPC, ...   many flags for diff. detectors
+//
+// History:
+//
+// 04.04.02 - first version
+// 
+////////////////////////////////////////////////////////////////////////
+
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include "iostream.h"
+#include "TDatetime.h"
+#include "STEER/AliRun.h"
+#include "STEER/AliRunDigitizer.h"
+#include "ITS/AliITSDigitizer.h"
+#include "ITS/AliITS.h"
+#include "ITS/AliITSDetType.h"
+#include "ITS/AliITSresponseSDD.h"
+#include "TPC/AliTPCDigitizer.h"
+#include "TRD/AliTRDdigitizer.h"
+#include "PHOS/AliPHOSDigitizer.h"
+#include "MUON/AliMUONDigitizer.h"
+#include "RICH/AliRICHDigitizer.h"
+#include "TStopwatch.h"
+#endif
+
+// #include "AliHits2SDigits.C"
+
+// void AliCopyN(TString inputFile, TString outputFile);
+
+Int_t MergeV1(TString fileNameDigits="digits.root", 
+           TString fileNameSDigitsSig="sig.sdigits.root", 
+           TString fileNameSDigitsBgr="bgr.sdigits.root", 
+           Int_t nEvents = 1, Int_t iITS = 2, Int_t iTPC = 0,
+           Int_t iTRD = 0,  Int_t iPHOS = 0, Int_t iMUON = 0,
+           Int_t iRICH = 0, Int_t iCopy = 1)
+{
+// delete the current gAlice object, the one from input file
+//  will be used
+
+  if(gAlice){
+    delete gAlice;
+    gAlice = 0;
+  } // end if gAlice
+
+  // Connect the Root Galice file containing Geometry, Kine and Hits
+  TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(fileNameSDigitsSig.Data());
+  if(!file) file = new TFile(fileNameSDigitsSig.Data());
+  TDatime *ct0 = new TDatime(2002,04,26,00,00,00), ct = file->GetCreationDate();
+  
+  // Get AliRun object from file or create it if not on file
+  if(!gAlice) {
+      gAlice = (AliRun*)file->Get("gAlice");
+      if(gAlice) printf("AliRun object found on file\n");
+      if(!gAlice) gAlice = new AliRun("gAlice","Alice test program");
+  } // end if !gAlice
+
+  AliRunDigitizer * manager = new AliRunDigitizer(2,1);
+  manager->SetInputStream(0,fileNameSDigitsSig.Data());
+  manager->SetInputStream(1,fileNameSDigitsBgr.Data());
+  if (fileNameDigits != "") {
+//    if (iCopy) {
+//      AliCopyN(fileNameSDigitsSig,fileNameDigits);
+//    }
+    manager->SetOutputFile(fileNameDigits);
+  }
+  manager->SetNrOfEventsToWrite(nEvents);
+  
+  if (iITS) {
+    AliITSDigitizer *dITS  = new AliITSDigitizer(manager);
+    if (iITS == 2) dITS->SetByRegionOfInterestFlag(1);
+    if(ct0->GetDate()>ct.GetDate()){
+       // For old files, must change SDD noise.
+       AliITS *ITS = (AliITS*) gAlice->GetDetector("ITS");
+       AliITSresponseSDD *resp1 = ITS->DetType(1)->GetResponseModel();
+       resp1->SetNoiseParam();
+       resp1->SetNoiseAfterElectronics();
+       Float_t n,b;
+       Int_t cPar[8];
+       resp1->GetNoiseParam(n,b);
+       n = resp1->GetNoiseAfterElectronics();
+       cPar[0]=0;
+       cPar[1]=0;
+       cPar[2]=(Int_t)(b + 2.*n + 0.5);
+       cPar[3]=(Int_t)(b + 2.*n + 0.5);
+       cPar[4]=0;
+       cPar[5]=0;
+       cPar[6]=0;
+       cPar[7]=0;
+       resp1->SetCompressParam(cPar);
+    } // end if
+  }
+  if (iTPC) AliTPCDigitizer *dTPC  = new AliTPCDigitizer(manager);
+  if (iTRD) AliTRDdigitizer *dTRD  = new AliTRDdigitizer(manager);
+  if (iPHOS) AliPHOSDigitizer *dPHOS  = new AliPHOSDigitizer(manager);
+  if (iMUON) AliMUONDigitizer *dMUON  = new AliMUONDigitizer(manager);
+  if (iRICH) AliRICHDigitizer *dRICH  = new AliRICHDigitizer(manager);
+  TStopwatch timer;
+  timer.Start();
+  manager->Exec("deb all");
+  timer.Stop(); 
+  timer.Print();
+//  delete gAlice;  
+//  gAlice = 0;
+  delete manager;
+}
+
+
+/*
+////////////////////////////////////////////////////////////////////////
+void AliCopyN(TString inputFileName, TString outputFileName) {
+// copy some objects
+
+  TFile *inputFile = OpenFile(inputFileName);
+  if (!inputFile) return;
+
+  TFile *outputFile = TFile::Open(outputFileName.Data(),"update");
+  if (!outputFile->IsOpen()) {
+    cerr<<"Can't open "<<outputFileName.Data()<<" !\n";
+    return;
+  }
+  if (!ImportgAlice(inputFile)) return;
+  AliCopy(inputFile, outputFile);
+  delete gAlice;
+  gAlice=0;
+  inputFile->Close();
+  delete inputFile;
+  outputFile->Close();
+  delete outputFile;
+}
+////////////////////////////////////////////////////////////////////////
+*/
diff --git a/ITS/oldmacros/ViewITSSPD.C b/ITS/oldmacros/ViewITSSPD.C
new file mode 100644 (file)
index 0000000..3800c2c
--- /dev/null
@@ -0,0 +1,70 @@
+void ViewITSSPD(Int_t version1,Int_t version2){
+    gMC->Gsatt("ITSV","seen",0);// Air
+      gMC->Gsatt("ITSD","seen",0);// Air
+        gMC->Gsatt("IT12","seen",0);// Air
+       if(version1==10){gMC->Gsatt("I651","seen",1);gMC->Gsatt("I651","colo",7);} // Services
+          gMC->Gsatt("I12B","seen",0);// Air
+            gMC->Gsatt("I10B","seen",0);// Air
+              gMC->Gsatt("I105","seen",1);gMC->Gsatt("I105","colo",7);// SPD End ladder
+              gMC->Gsatt("I108","seen",1);gMC->Gsatt("I108","colo",7);// SPD bus
+              gMC->Gsatt("I109","seen",1);gMC->Gsatt("I109","colo",7);// SPD bus
+              gMC->Gsatt("I107","seen",0);//Air
+                gMC->Gsatt("I106","seen",1);gMC->Gsatt("I106","colo",6);// Silicon
+                gMC->Gsatt("I101","seen",1);gMC->Gsatt("I101","colo",6);// Silicon
+                  gMC->Gsatt("ITS1","seen",1);gMC->Gsatt("ITS1","colo",6);// Silicon
+            gMC->Gsatt("I20B","seen",0);// Air
+              gMC->Gsatt("I105","seen",1);gMC->Gsatt("I105","colo",7);// SPD End ladder
+              gMC->Gsatt("I109","seen",1);gMC->Gsatt("I109","colo",7);// SPD bus
+              gMC->Gsatt("I108","seen",1);gMC->Gsatt("I108","colo",7);// SPD bus
+              gMC->Gsatt("I107","seen",0);// Air
+                gMC->Gsatt("I1D6","seen",1);gMC->Gsatt("I1D6","colo",6);// Silicon
+                gMC->Gsatt("I1D1","seen",1);gMC->Gsatt("I1D1","colo",6);// Silicon
+                  gMC->Gsatt("ITS2","seen",1);gMC->Gsatt("ITS2","colo",6);// Silicon
+            gMC->Gsatt("I123","seen",1);gMC->Gsatt("I123","colo",1);// carbon
+            gMC->Gsatt("I121","seen",1);gMC->Gsatt("I121","colo",1);// carbon
+            gMC->Gsatt("I122","seen",1);gMC->Gsatt("I122","colo",1);// carbon
+            gMC->Gsatt("I120","seen",1);gMC->Gsatt("I120","colo",1);// carbon
+            gMC->Gsatt("I144","seen",1);gMC->Gsatt("I144","colo",1);// carbon
+            gMC->Gsatt("I143","seen",1);gMC->Gsatt("I143","colo",1);// carbon
+            gMC->Gsatt("I142","seen",1);gMC->Gsatt("I142","colo",1);// carbon
+            gMC->Gsatt("I141","seen",1);gMC->Gsatt("I141","colo",1);// carbon
+            gMC->Gsatt("I140","seen",1);gMC->Gsatt("I140","colo",1);// carbon
+            gMC->Gsatt("I139","seen",1);gMC->Gsatt("I139","colo",1);// carbon
+            gMC->Gsatt("I138","seen",1);gMC->Gsatt("I138","colo",1);// carbon
+            gMC->Gsatt("I137","seen",1);gMC->Gsatt("I137","colo",1);// carbon
+            gMC->Gsatt("I136","seen",1);gMC->Gsatt("I136","colo",1);// carbon
+            gMC->Gsatt("I135","seen",1);gMC->Gsatt("I135","colo",1);// carbon
+            gMC->Gsatt("I134","seen",1);gMC->Gsatt("I134","colo",1);// carbon
+            gMC->Gsatt("I133","seen",1);gMC->Gsatt("I133","colo",1);// carbon
+            gMC->Gsatt("I132","seen",1);gMC->Gsatt("I132","colo",1);// carbon
+            gMC->Gsatt("I131","seen",1);gMC->Gsatt("I131","colo",1);// carbon
+            gMC->Gsatt("I130","seen",1);gMC->Gsatt("I130","colo",1);// carbon
+            gMC->Gsatt("I129","seen",1);gMC->Gsatt("I129","colo",1);// carbon
+            gMC->Gsatt("I128","seen",1);gMC->Gsatt("I128","colo",1);// carbon
+            gMC->Gsatt("I126","seen",1);gMC->Gsatt("I126","colo",1);// carbon
+            gMC->Gsatt("I125","seen",1);gMC->Gsatt("I125","colo",1);// carbon
+            gMC->Gsatt("I124","seen",1);gMC->Gsatt("I124","colo",1);// carbon
+            gMC->Gsatt("I113","seen",0);//Air
+              gMC->Gsatt("I112","seen",1);gMC->Gsatt("I112","colo",1);// carbon
+              gMC->Gsatt("I111","seen",1);gMC->Gsatt("I111","colo",1);// carbon
+              gMC->Gsatt("I118","seen",1);gMC->Gsatt("I118","colo",2);// glue
+              gMC->Gsatt("I110","seen",1);gMC->Gsatt("I110","colo",1);// carbon
+              gMC->Gsatt("I114","seen",1);gMC->Gsatt("I114","colo",3);// INOX
+              gMC->Gsatt("I115","seen",1);gMC->Gsatt("I115","colo",4);// water
+              gMC->Gsatt("I116","seen",1);gMC->Gsatt("I116","colo",3);// INOX
+                gMC->Gsatt("I117","seen",1);gMC->Gsatt("I117","colo",4);// Water
+          gMC->Gsatt("I650","seen",0); //air
+           gMC->Gsatt("I666","seen",1);gMC->Gsatt("I666","colo",5);//Aluminum
+            gMC->Gsatt("I667","seen",1);gMC->Gsatt("I667","colo",5);//Aluminum
+              gMC->Gsatt("I668","seen",1);gMC->Gsatt("I668","colo",4);// water
+            gMC->Gsatt("I669","seen",1);gMC->Gsatt("I669","colo",3);// INOX
+              gMC->Gsatt("I670","seen",1);gMC->Gsatt("I670","colo",4);// water
+            gMC->Gsatt("I671","seen",1);gMC->Gsatt("I671","colo",6);// Silicon
+              gMC->Gsatt("I672","seen",1);gMC->Gsatt("I672","colo",4);// water
+            gMC->Gsatt("I673","seen",1);gMC->Gsatt("I673","colo",6);// Silicon
+              gMC->Gsatt("I674","seen",1);gMC->Gsatt("I674","colo",4);// water
+              gMC->Gsatt("I675","seen",1);gMC->Gsatt("I675","colo",5);//Aluminum
+            gMC->Gsatt("I676","seen",1);gMC->Gsatt("I676","colo",6);// Silicon
+              gMC->Gsatt("I677","seen",1);gMC->Gsatt("I677","colo",4);// water
+              gMC->Gsatt("I678","seen",1);gMC->Gsatt("I678","colo",5);//Aluminum
+}
diff --git a/ITS/oldmacros/ViewS01.C b/ITS/oldmacros/ViewS01.C
new file mode 100644 (file)
index 0000000..40d5e81
--- /dev/null
@@ -0,0 +1,14 @@
+void ViewS01(){
+    gMC->Gsatt("ITSV","seen",0); // Air "Gen Air"
+      gMC->Gsatt("ITSD","seen",0); // Air "Gen Air"
+        gMC->Gsatt("IS01","seen",0); //  Air "Air"
+         gMC->Gsatt("ICY2","seen",1); // Carbon Fiber "SDD C (M55J)
+         gMC->Gsatt("I212","seen",1); // Carbon "ITS SandW C"
+         gMC->Gsatt("I211","seen",1); // PC board "G10FR4"
+         gMC->Gsatt("I217","seen",1); // Carbon Fiber "SDD/SSD rings"
+         gMC->Gsatt("I219","seen",1); // Carbon Fiber "SDD/SSD rings"
+         gMC->Gsatt("I214","seen",1); // PC board "G10FR4"
+         gMC->Gsatt("I213","seen",1); // PC board "G10FR4"
+         gMC->Gsatt("I215","seen",2); // Air "Air"
+         gMC->Gsatt("I216","seen",1); // Carbon "ITS Sandw C"
+}
diff --git a/ITS/oldmacros/ViewSPD.C b/ITS/oldmacros/ViewSPD.C
new file mode 100644 (file)
index 0000000..4abf2c9
--- /dev/null
@@ -0,0 +1,55 @@
+void ViewSPD(){
+    gMC->Gsatt("ITSV","seen",0);
+      gMC->Gsatt("ITSD","seen",0);
+        gMC->Gsatt("IT12","seen",0);
+          gMC->Gsatt("I12B","seen",0);
+            gMC->Gsatt("I10B","seen",0);
+              gMC->Gsatt("I105","seen",1);
+              gMC->Gsatt("I108","seen",1);
+              gMC->Gsatt("I109","seen",1);
+              gMC->Gsatt("I107","seen",0);
+                gMC->Gsatt("I106","seen",1);
+                gMC->Gsatt("I101","seen",0);
+                  gMC->Gsatt("ITS1","seen",1);
+            gMC->Gsatt("I20B","seen",0);
+              gMC->Gsatt("I105","seen",1);
+              gMC->Gsatt("I109","seen",1);
+              gMC->Gsatt("I108","seen",1);
+              gMC->Gsatt("I107","seen",0);
+                gMC->Gsatt("I1D6","seen",1);
+                gMC->Gsatt("I1D1","seen",0);
+                  gMC->Gsatt("ITS2","seen",1);
+            gMC->Gsatt("I123","seen",1);
+            gMC->Gsatt("I121","seen",1);
+            gMC->Gsatt("I122","seen",1);
+            gMC->Gsatt("I120","seen",1);
+            gMC->Gsatt("I144","seen",1);
+            gMC->Gsatt("I143","seen",1);
+            gMC->Gsatt("I142","seen",1);
+            gMC->Gsatt("I141","seen",1);
+            gMC->Gsatt("I140","seen",1);
+            gMC->Gsatt("I139","seen",1);
+            gMC->Gsatt("I138","seen",1);
+            gMC->Gsatt("I137","seen",1);
+            gMC->Gsatt("I136","seen",1);
+            gMC->Gsatt("I135","seen",1);
+            gMC->Gsatt("I134","seen",1);
+            gMC->Gsatt("I133","seen",1);
+            gMC->Gsatt("I132","seen",1);
+            gMC->Gsatt("I131","seen",1);
+            gMC->Gsatt("I130","seen",1);
+            gMC->Gsatt("I129","seen",1);
+            gMC->Gsatt("I128","seen",1);
+            gMC->Gsatt("I126","seen",1);
+            gMC->Gsatt("I125","seen",1);
+            gMC->Gsatt("I124","seen",1);
+            gMC->Gsatt("I113","seen",0);
+              gMC->Gsatt("I112","seen",1);
+              gMC->Gsatt("I111","seen",1);
+              gMC->Gsatt("I118","seen",1);
+              gMC->Gsatt("I110","seen",1);
+              gMC->Gsatt("I114","seen",1);
+              gMC->Gsatt("I115","seen",1);
+              gMC->Gsatt("I116","seen",0);
+                gMC->Gsatt("I117","seen",1);
+}