PWG3/muon/AliAnalysisTaskMuonAODfromGeneral.cxx

   1 #define AliAnalysisTaskMuonAODfromGeneral_cxx
   2
   3 // 19 Nov 2007
   4 // Class implementation for the specific muon AOD generation
   5 // Extracts only muon tracks from a general AOD and builds dimuons
   6 // Livio Bianchi, Universita' di Torino
   7
   8
   9 #include "TTree.h"
  10 #include "TROOT.h"
  11 #include "TH1.h"
  12 #include "TCanvas.h"
  13 #include "TSystem.h"
  14 #include "TRandom.h"
  15
  16 #include "AliAODEvent.h"
  17 #include "AliAnalysisTask.h"
  18 #include "AliAnalysisTaskMuonAODfromGeneral.h"
  19 #include "AliAODHandler.h"
  20 #include "AliAnalysisManager.h"
  21 #include "AliAODEventInfo.h"
  22 #include "AliAODDimuon.h"
  23
  24 ClassImp(AliAnalysisTaskMuonAODfromGeneral)
  25
  26 //________________________________________________________________________
  27 AliAnalysisTaskMuonAODfromGeneral::AliAnalysisTaskMuonAODfromGeneral(const char *name, Double_t BeamEnergy):
  28         AliAnalysisTask(name, "AnalysisTaskMuonAODfromGeneral"),
  29         fInfos(0),
  30         fDimuons(0),
  31         fChain(0),
  32         fOrgAOD(0),
  33         fNewAOD(0),
  34         ft(0),
  35         fBeamEnergy(0)
  36 {  // Constructor.
  37   // Input slot #0 works with an Ntuple
  38   SetBeamEnergy(BeamEnergy);
  39   DefineInput(0, TChain::Class());
  40   // Output slot #0 writes into a TTree container
  41   DefineOutput(0, TTree::Class());
  42 }
  43
  44 //________________________________________________________________________
  45 void AliAnalysisTaskMuonAODfromGeneral::ConnectInputData(Option_t *) {
  46   printf("   ConnectInputData %s\n", GetName());
  47   // New connection to the AOD
  48   fChain = (TChain*)GetInputData(0);
  49   fOrgAOD = new AliAODEvent();
  50   fOrgAOD->ReadFromTree(fChain); // This checks also if the branch address is already set
  51 }
  52
  53 //________________________________________________________________________
  54 void AliAnalysisTaskMuonAODfromGeneral::CreateOutputObjects() {
  55   printf("Creating output objects\n");
  56   if(!fNewAOD)fNewAOD=new AliAODEvent();
  57   fNewAOD->CreateStdContent();
  58
  59   if (!ft) {
  60     ft = new TTree("AOD","Muon AOD tree",0);
  61     ft->BranchRef();
  62     ft->Branch(fNewAOD->GetList());
  63     //ft->Print();
  64  }
  65
  66   //Add new Objects: MuonInfos & Dimuons
  67    fInfos = new AliAODEventInfo();
  68    fNewAOD->AddObject(fInfos);
  69    const char *nameInfos = "MuonInfos";
  70    fInfos->SetName(nameInfos);
  71    ft->Branch(nameInfos, &fInfos);
  72    fDimuons = new TClonesArray("AliAODDimuon",0);
  73    fNewAOD->AddObject(fDimuons);
  74    const char *nameDimuons = "Dimuons";
  75    fDimuons->SetName(nameDimuons);
  76    ft->Branch(nameDimuons, &fDimuons);
  77
  78 }
  79
  80 //________________________________________________________________________
  81 void AliAnalysisTaskMuonAODfromGeneral::Exec(Option_t *) {
  82   static int ncall=0;
  83   // Task making a Muon AOD
  84   // Get input data
  85   TChain *chain = (TChain*)GetInputData(0);
  86   Long64_t ientry = chain->GetReadEntry();
  87   ientry=ientry;
  88
  89   if (!fOrgAOD) return;
  90
  91   Int_t nTracks=fOrgAOD->GetNumberOfTracks();
  92   Int_t nMuTracks=0;
  93   Int_t nPosTracks = 0;
  94   Int_t mutrNumb[10]; for (Int_t i=0; i<10; mutrNumb[i++]=0);
  95   for (Int_t iTrack=0; iTrack<nTracks; iTrack++){
  96        const Double_t *trackpid=fOrgAOD->GetTrack(iTrack)->PID();
  97        if (trackpid[AliAODTrack::kMuon]==1.) {
  98                 mutrNumb[nMuTracks]=iTrack;
  99                 nMuTracks++;
 100                 }
 101        if (fOrgAOD->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
 102        }
 103
 104   Bool_t ExistMuon=kFALSE;
 105   if (nMuTracks>0) ExistMuon=kTRUE;
 106
 107   //--------------------------------------------------------------------
 108   // set arrays and pointers
 109   Double_t pos[3];
 110   Double_t covVtx[6];
 111
 112     // Access to the header
 113     AliAODHeader *header = fNewAOD->GetHeader();
 114
 115   // fill the header
 116   header->SetRunNumber       (fOrgAOD->GetRunNumber()       );
 117   header->SetBunchCrossNumber(fOrgAOD->GetBunchCrossNumber());
 118   header->SetOrbitNumber     (fOrgAOD->GetOrbitNumber()     );
 119   header->SetPeriodNumber    (fOrgAOD->GetPeriodNumber()    );
 120   header->SetTriggerMask     (fOrgAOD->GetTriggerMask()     );
 121   header->SetTriggerCluster  (fOrgAOD->GetTriggerCluster()  );
 122   header->SetEventType       (fOrgAOD->GetEventType()       );
 123   header->SetMagneticField   (fOrgAOD->GetMagneticField()   );
 124   header->SetZDCN1Energy     (fOrgAOD->GetZDCN1Energy()     );
 125   header->SetZDCP1Energy     (fOrgAOD->GetZDCP1Energy()     );
 126   header->SetZDCN2Energy     (fOrgAOD->GetZDCN2Energy()     );
 127   header->SetZDCP2Energy     (fOrgAOD->GetZDCP2Energy()     );
 128   header->SetZDCEMEnergy     (0,fOrgAOD->GetZDCEMEnergy(0)  );
 129   header->SetZDCEMEnergy     (1,fOrgAOD->GetZDCEMEnergy(1)  );
 130   header->SetRefMultiplicity   (nTracks);
 131   header->SetRefMultiplicityPos(nPosTracks);
 132   header->SetRefMultiplicityNeg(nTracks - nPosTracks);
 133   header->SetMuonMagFieldScale(-999.); // FIXME
 134   header->SetCentrality(-999.);        // FIXME
 135
 136   Int_t nVertices = 0;
 137   const AliAODVertex *vtx = fOrgAOD->GetPrimaryVertex();
 138   if(vtx)nVertices = 1;
 139
 140   Int_t nDimuons=nMuTracks*(nMuTracks-1)/2;
 141
 142   // Access to the AOD container of vertices
 143   TClonesArray &vertices = *(fNewAOD->GetVertices());
 144   Int_t jVertices=0;
 145
 146   // Access to the AOD container of tracks
 147    TClonesArray &tracks = *(fNewAOD->GetTracks());
 148    tracks.Delete();
 149    Int_t jTracks=0;
 150
 151
 152   vtx->GetCovMatrix(covVtx); //covariance matrix
 153
 154   AliAODVertex * primary = new(vertices[jVertices++])
 155       AliAODVertex(pos, covVtx, vtx->GetChi2perNDF(), NULL, AliAODVertex::kPrimary);
 156
 157   static int ncal=0;
 158   static int numtracks=0;
 159   for (Int_t iTrack=0; iTrack<nMuTracks; iTrack++) {
 160         AliAODTrack *aodTrack;
 161         primary->AddDaughter( aodTrack=new(tracks[jTracks++])
 162         AliAODTrack((*(fOrgAOD->GetTrack(mutrNumb[iTrack])))));
 163         ncal++;
 164   }
 165   numtracks+=jTracks;
 166
 167   fDimuons->Delete();
 168   fDimuons->Expand(nDimuons);
 169   TClonesArray &rDimuons = *fDimuons;
 170
 171   Int_t jDimuons=0;
 172   for(Int_t i=0; i<nMuTracks; i++){
 173      for(Int_t j=i+1; j<nMuTracks; j++){
 174        new(rDimuons[jDimuons++]) AliAODDimuon(tracks[i],tracks[j],fInfos);
 175      }
 176   }
 177
 178   fInfos->SetBeamEnergy(fBeamEnergy);
 179   fInfos->SetEv(fNewAOD);
 180   fInfos->SetEi(fInfos);
 181   fInfos->SetHe(header);
 182   fInfos->SetTr(fNewAOD->GetTracks());
 183   fInfos->SetDi(fDimuons);
 184   fInfos->SelectTriggerBits(0,1,2,3,4,5);
 185   if(ExistMuon) ft->Fill();
 186   ncall++;
 187
 188   // Post final data. It will be written to a file with option "RECREATE"
 189   PostData(0, ft);
 190 }
 191
 192 //________________________________________________________________________
 193 void AliAnalysisTaskMuonAODfromGeneral::Terminate(Option_t *) {
 194   // Draw some histogram at the end.
 195   ft->GetUserInfo()->Add(fNewAOD);
 196   ft->Write("",TObject::kOverwrite);
 197   ft->ResetBranchAddresses();
 198   if (!gROOT->IsBatch()) {
 199     TCanvas *c1 = new TCanvas("c1","Pt",10,10,310,310);
 200     c1->SetFillColor(10);
 201     c1->SetHighLightColor(10);
 202
 203     c1->cd(1)->SetLeftMargin(0.15);
 204     c1->cd(1)->SetBottomMargin(0.15);
 205     c1->cd(1)->SetLogy();
 206     ft->Draw("tracks.Pt()");
 207     //ft->Draw("dimuons.M()");
 208   }
 209 }