PWG/muon/AliAnalysisTaskMuonAODfromGeneral.cxx

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