AliAnalysisTaskMyTask.cxx 9.83 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* AliAnaysisTaskMyTask
 *
 * empty task which can serve as a starting point for building an analysis
 * as an example, one histogram is filled
 */

#include "TChain.h"
#include "TH1F.h"
#include "TList.h"
#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"
#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
#include "AliAnalysisTaskMyTask.h"
#include "AliESDtrack.h"
Rafael Pezzi's avatar
Rafael Pezzi committed
31
#include "AliESDVertex.h"
32
33
#include "AliVertex.h"
#include "Riostream.h"
34
35


36
Int_t esd_event_id = 0; // global variable to store unique event id
37
38
39
40
41
42
43
44

class AliAnalysisTaskMyTask;    // your analysis class

using namespace std;            // std namespace: so you can do things like 'cout'

ClassImp(AliAnalysisTaskMyTask) // classimp: necessary for root

AliAnalysisTaskMyTask::AliAnalysisTaskMyTask() : AliAnalysisTaskSE(),
45
    fESD(0), fOutputList(0), fHistPt(0),    fHistEvents(0), fHistMass(0)
46
47
48
49
50
51
{
    // default constructor, don't allocate memory here!
    // this is used by root for IO purposes, it needs to remain empty
}
//_____________________________________________________________________________
AliAnalysisTaskMyTask::AliAnalysisTaskMyTask(const char* name) : AliAnalysisTaskSE(name),
52
    fESD(0), fOutputList(0), fHistPt(0),  fHistEvents(0), fHistMass(0)
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
{
    // constructor
    DefineInput(0, TChain::Class());    // define the input of the analysis: in this case we take a 'chain' of events
                                        // this chain is created by the analysis manager, so no need to worry about it,
                                        // it does its work automatically
    DefineOutput(1, TList::Class());    // define the ouptut of the analysis: in this case it's a list of histograms
                                        // you can add more output objects by calling DefineOutput(2, classname::Class())
                                        // if you add more output objects, make sure to call PostData for all of them, and to
                                        // make changes to your AddTask macro!
}
//_____________________________________________________________________________
AliAnalysisTaskMyTask::~AliAnalysisTaskMyTask()
{
    // destructor
    if(fOutputList) {
        delete fOutputList;     // at the end of your task, it is deleted from memory by calling this function
    }
}
//_____________________________________________________________________________
void AliAnalysisTaskMyTask::UserCreateOutputObjects()
{
    // create output objects
    //
    // this function is called ONCE at the start of your analysis (RUNTIME)
    // here you ceate the histograms that you want to use
    //
    // the histograms are in this case added to a tlist, this list is in the end saved
    // to an output file
    //
    fOutputList = new TList();          // this is a list which will contain all of your histograms
                                        // at the end of the analysis, the contents of this list are written
                                        // to the output file
    fOutputList->SetOwner(kTRUE);       // memory stuff: the list is owner of all objects it contains and will delete them
                                        // if requested (dont worry about this now)

    // example of a histogram
89
    fHistPt = new TH1F("fHistPt", "fHistPt", 100, 0, 10);       // create your histogram
90
91
92
    fOutputList->Add(fHistPt);          // don't forget to add it to the list! the list will be written to file, so if you want
                                        // your histogram in the output file, add it to the list!

93
94
95
96
97
98
    // my particle mass histogram
    Double_t fHistMassEdges[12] = {0.0,0.0005,0.0405,0.08,0.12,0.13,0.17,0.48,0.52,0.92,0.96,1.0}; // 11 bins =>> has 11+1 = 12 edges

    fHistMass = new TH1F("fHistMass","Particle Histogram;M_{particle}", 11, fHistMassEdges);
    fOutputList->Add(fHistMass);

99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119

    // Histograms for dimuons
    fHistEvents = new TH1F("fHistEvents","fHistEvents;N_{events}",100,0.,10000.);

    fOutputList->Add(fHistEvents);

    PostData(1, fOutputList);           // postdata will notify the analysis manager of changes / updates to the
                                        // fOutputList object. the manager will in the end take care of writing your output to file
                                        // so it needs to know what's in the output
}
//_____________________________________________________________________________
void AliAnalysisTaskMyTask::UserExec(Option_t *)
{
    // user exec
    // this function is called once for each event
    // the manager will take care of reading the events from file, and with the static function InputEvent() you
    // have access to the current event.
    // once you return from the UserExec function, the manager will retrieve the next event from the chain

    Int_t Event=0;
    Int_t TrigEvent=0;
120
    ofstream summaryfile, eventsdetail;
121
122
123
    summaryfile.open ("events_summary.txt",std::ofstream::app);
    eventsdetail.open ("events_detail.txt",std::ofstream::app);

124
125
126
127


    fESD = dynamic_cast<AliESDEvent*>(InputEvent());    // get an event (called fESD) from the input file

128
                                                        // there's another event format (ESD) which works in a similar way
129

130
                                                  // but is more cpu/memory unfriendly. for now, we'll stick with aod's
131
132
133
    if(!fESD) return;                                   // if the pointer to the event is empty (getting it failed) skip this event
        // example part: i'll show how to loop over the tracks in an event
        // and extract some information from them which we'll store in a histogram
134

135
    Int_t iTracks(fESD->GetNumberOfTracks());           // see how many tracks there are in the event
136

137
138
139
140
141
    Double_t Vx = fESD->GetPrimaryVertex()->GetX();	// gets vertexes from individual events
    Double_t Vy = fESD->GetPrimaryVertex()->GetY();
    Double_t Vz = fESD->GetPrimaryVertex()->GetZ();
    Double_t MagneticField = fESD->GetMagneticField();

142
143
144
145
146
147
148
    // Add Event details to summaryfile and eventsdetail
    summaryfile << " Event " << esd_event_id << ": iTracks = "  << iTracks << endl;
    summaryfile << " Magnetic Field: " << MagneticField << endl;
    summaryfile << " Vertex: " << endl << " X = " << Vx << endl << " Y = " << Vy << endl << " Z = " << Vz << endl << endl;

    eventsdetail << "\n" << "********** Event " << esd_event_id << ": iTracks = "  << iTracks << " **********" << endl << endl;
    eventsdetail << "Event\tTrack\tMass\t\tEnergy\t\tPx\t\tPy\t\tPt(file)\tPt(calculado)\tPz\t\tCharge" << endl << endl;
149
150

/*
Rafael Pezzi's avatar
Rafael Pezzi committed
151

152
153
154
155
Assumed Units: Mass (GeV/c^2)[CONFIRMED] || Energy (GeV) || Momentum (GeV/c) || Charge (* 1.6*10^-19 C)

*/

156

157
    for(Int_t i(0); i < iTracks; i++) {                 // loop over all these tracks
158

159
        AliESDtrack* track = static_cast<AliESDtrack*>(fESD->GetTrack(i));         // get a track (type AliESDtrack) from the event
160

161
        if(!track) continue;                            // if we failed, skip this track
162
163
164
165
166
167
168
169
170
171
	
	Double_t Mass = track->M(); // returns the pion mass, if the particle can't be identified properly
	Double_t Energy = track->E(); // Returns the energy of the particle given its assumed mass, but assumes the pion mass if the particle can't be identified properly.
	Double_t Px = track->Px();
	Double_t Py = track->Py();
	Double_t Pt = track->Pt();
	Double_t PtCal = sqrt(Px*Px+Py*Py); // here, we calculate transversal momentum to make sure of its meaning (then compare it to Pt)
	Double_t Pz = track->Pz();
	Double_t Charge = track->Charge();

172
173
174
175
        // Add EVENT, TRACK NUMBER, MASS, ENERGY, corresponding MOMENTUM (x, y, transversal, z) and CHARGE to eventsdetail txt file
	eventsdetail << esd_event_id << "\t" << i << "\t";
	eventsdetail << fixed << Mass << "\t" << Energy << "\t";
	eventsdetail << fixed << Px << "\t" << Py << " \t " << Pt << " \t " << PtCal << " \t " << Pz << "\t" << Charge << endl;
176
177
178

	//'fixed' fixes the number of decimal places so numbers are vertically aligned

179

180
        fHistPt->Fill(track->Pt());                     // plot the pt value of the track in a histogram
181
182
	fHistMass->Fill(Mass);

183
184
185
    }

    Event++;
186
    esd_event_id++; // Increment global esd_event_id
187
    fHistEvents->Fill(Event);
188
189
    summaryfile.close();
    eventsdetail.close();
190
191
192
193
194
195
196
197
198
199
200
201
202

                                                       // continue until all the tracks are processed
    PostData(1, fOutputList);                           // stream the results the analysis of this event to
                                                        // the output manager which will take care of writing
                                                        // it to a file
}
//_____________________________________________________________________________
void AliAnalysisTaskMyTask::Terminate(Option_t *)
{
    // terminate
    // called at the END of the analysis (when all events are processed)
}
//_____________________________________________________________________________