HBTAN/hbtanalysis.C

   1 void hbtanalysis(Option_t* datatype, Option_t* processopt="TracksAndParticles",
   2                 Int_t first = -1,Int_t last = -1,
   3                 char *outfile = "hbtanalysis.root")
   4  {
   5 //HBT Anlysis Macro
   6 //Anlyzes TPC recontructed tracks and simulated particles that corresponds to them
   7
   8 //datatype defines type of data to be read
   9 //  Kine  - analyzes Kine Tree: simulated particles
  10 //  TPC   - analyzes TPC   tracking + particles corresponding to these tracks
  11 //  ITSv1 - analyzes ITSv1 ----------------------//--------------------------
  12 //  ITSv2 - analyzes ITSv2 ----------------------//--------------------------
  13
  14 //processopt defines option passed to AliHBTAnlysis::Process method
  15 // default: TracksAndParticles - process both recontructed tracks and sim. particles corresponding to them
  16 //          Tracks - process only recontructed tracks
  17 //          Particles - process only simulated particles
  18
  19 //Reads data from diroctories from first to last(including)
  20 // For examples if first=3 and last=5 it reads from
  21 //  ./3/
  22 //  ./4/
  23 //  ./5/
  24 //if first or last is negative (or both), it reads from current directory
  25 //
  26 //these names I use when analysis is done directly from CASTOR files via RFIO
  27 //  const char* basedir="rfio:/castor/cern.ch/user/s/skowron/";
  28 //  const char* serie="standard";
  29 //  const char* field = "0.4";
  30   const char* basedir=".";
  31   const char* serie="";
  32   const char* field = "";
  33
  34   // Dynamically link some shared libs
  35   gROOT->LoadMacro("loadlibs.C");
  36   loadlibs();
  37
  38   gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libHBTAN");
  39
  40
  41   /***********************************************************/
  42   //Create analysis and reader
  43   AliHBTAnalysis * analysis = new AliHBTAnalysis();
  44
  45   AliHBTReader* reader;
  46   Int_t kine = strcmp(datatype,"Kine");
  47   Int_t TPC = strcmp(datatype,"TPC");
  48   Int_t ITSv1 = strcmp(datatype,"ITSv1");
  49   Int_t ITSv2 = strcmp(datatype,"ITSv2");
  50   Int_t intern = strcmp(datatype,"Intern");
  51
  52   if(!kine)
  53    {
  54     reader = new AliHBTReaderKineTree();
  55     processopt="Particles"; //this reader by definition reads only simulated particles
  56    }
  57   else if(!TPC)
  58    {
  59     reader = new AliHBTReaderTPC();
  60    }
  61   else if(!ITSv1)
  62    {
  63     reader = new AliHBTReaderITSv1();
  64    }
  65   else if(!ITSv2)
  66    {
  67     reader = new AliHBTReaderITSv2();
  68    }
  69   else if(!intern)
  70    {
  71     reader = new AliHBTReaderInternal("Kine.sum.root");
  72    }
  73   else
  74    {
  75     cerr<<"Option "<<datatype<<"  not recognized. Exiting"<<endl;
  76     return;
  77    }
  78
  79   TObjArray* dirs=0;
  80   if ( (first >= 0) && (last>=0) && ( (last-first)>=0 ) )
  81    {//read from many dirs dirs
  82      char buff[50];
  83      dirs = new TObjArray(last-first+1);
  84      for (Int_t i = first; i<=last; i++)
  85       {
  86         sprintf(buff,"%s/%s/%s/%d",basedir,field,serie,i);
  87         TObjString *odir= new TObjString(buff);
  88         dirs->Add(odir);
  89       }
  90     }
  91    reader->SetDirs(dirs);
  92
  93   /***********************************************************/
  94
  95   //we want have only low pt pi+ so set a cut to reader
  96   AliHBTParticleCut* readerpartcut= new AliHBTParticleCut();
  97   readerpartcut->SetPtRange(0.0,1.5);
  98   readerpartcut->SetPID(kPiPlus);
  99   reader->AddParticleCut(readerpartcut);//read this particle type with this cut
 100
 101   analysis->SetReader(reader);
 102   /************************************************************/
 103
 104   AliHBTPairCut *paircut = new AliHBTPairCut();
 105   paircut->SetQInvRange(0.0,0.20);
 106   analysis->SetGlobalPairCut(paircut);
 107
 108   AliHBTQInvCorrelFctn * qinvcfT= new AliHBTQInvCorrelFctn();
 109   AliHBTQInvCorrelFctn * qinvcfP= new AliHBTQInvCorrelFctn();
 110
 111   analysis->AddTrackFunction(qinvcfT);
 112   analysis->AddParticleFunction(qinvcfP);
 113
 114
 115   qinvcfP->Rename("qinvcfP","Particle (simulated) Qinv CF \\pi^{+} \\pi^{+}");
 116   qinvcfT->Rename("qinvcfT","Track (recontructed) Qinv CF \\pi^{+} \\pi^{+}");
 117
 118   AliHBTQOutCMSLCCorrelFctn* qoutP = new AliHBTQOutCMSLCCorrelFctn();
 119   qoutP->Rename("qoutP","Particle (simulated) Q_{out} CF \\pi^{+} \\pi^{+}");
 120   AliHBTQOutCMSLCCorrelFctn* qoutT = new AliHBTQOutCMSLCCorrelFctn();
 121   qoutT->Rename("qoutT","Track (recontructed) Q_{out} CF \\pi^{+} \\pi^{+}");
 122
 123   AliHBTPairCut *outPairCut = new AliHBTPairCut();
 124   outPairCut->SetQOutCMSLRange(0.0,0.15);
 125   outPairCut->SetQSideCMSLRange(0.0,0.02);
 126   outPairCut->SetQLongCMSLRange(0.0,0.02);
 127   qoutP->SetPairCut(outPairCut);
 128   qoutT->SetPairCut(outPairCut);
 129
 130   AliHBTQSideCMSLCCorrelFctn* qsideP = new AliHBTQSideCMSLCCorrelFctn();
 131   qsideP->Rename("qsideP","Particle (simulated) Q_{side} CF \\pi^{+} \\pi^{+}");
 132   AliHBTQSideCMSLCCorrelFctn* qsideT = new AliHBTQSideCMSLCCorrelFctn();
 133   qsideT->Rename("qsideT","Track (recontructed) Q_{side} CF \\pi^{+} \\pi^{+}");
 134
 135   AliHBTPairCut *sidePairCut = new AliHBTPairCut();
 136   sidePairCut->SetQOutCMSLRange(0.0,0.02);
 137   sidePairCut->SetQSideCMSLRange(0.0,0.15);
 138   sidePairCut->SetQLongCMSLRange(0.0,0.02);
 139   qsideP->SetPairCut(sidePairCut);
 140   qsideT->SetPairCut(sidePairCut);
 141
 142
 143   AliHBTQLongCMSLCCorrelFctn* qlongP = new AliHBTQLongCMSLCCorrelFctn();
 144   qlongP->Rename("qlongP","Particle (simulated) Q_{long} CF \\pi^{+} \\pi^{+}");
 145   AliHBTQLongCMSLCCorrelFctn* qlongT = new AliHBTQLongCMSLCCorrelFctn();
 146   qlongT->Rename("qlongT","Track (recontructed) Q_{long} CF \\pi^{+} \\pi^{+}");
 147
 148   AliHBTPairCut *longPairCut = new AliHBTPairCut();
 149   longPairCut->SetQOutCMSLRange(0.0,0.02);
 150   longPairCut->SetQSideCMSLRange(0.0,0.02);
 151   longPairCut->SetQLongCMSLRange(0.0,0.15);
 152   qlongP->SetPairCut(longPairCut);
 153   qlongT->SetPairCut(longPairCut);
 154
 155   analysis->AddTrackFunction(qoutP);
 156   analysis->AddParticleFunction(qoutT);
 157
 158   analysis->AddTrackFunction(qsideP);
 159   analysis->AddParticleFunction(qsideT);
 160
 161   analysis->AddTrackFunction(qlongT);
 162   analysis->AddParticleFunction(qlongT);
 163
 164
 165   analysis->Process(processopt);
 166
 167   TFile histoOutput(outfile,"recreate");
 168   analysis->WriteFunctions();
 169   histoOutput.Close();
 170
 171   delete qinvcfP;
 172   delete qinvcfT;
 173   delete paircut;
 174   delete readerpartcut;
 175   if (dirs)
 176    {
 177     dirs->SetOwner();
 178     delete dirs;
 179    }
 180   delete reader;
 181   delete analysis;
 182
 183  }
 184