ClassImp(AliJetEventParticles)
-/**************************************************************************/
-
AliJetEventParticles::AliJetEventParticles(Int_t size) :
fNParticles(0),
fParticles(new TClonesArray("AliJetParticle",size)),
fVertexX(0.),
fVertexY(0.),
- fVertexZ(0.)
+ fVertexZ(0.),
+ fTrials(0),
+ fNJets(0),
+ fNUQJets(0),
+ fXJet(-1),
+ fYJet(-1)
{
- //default constructor
+ // Default Constructor
+ for (Int_t i = 0; i < 4; i++) fZquench[i] = 0.;
+ for (Int_t i = 0; i < 10; i++)
+ for (Int_t j = 0; j < 4; j++) {
+ fJets[j][i]=0; // Trigger jets
+ fUQJets[j][i]=0; // Unquenched trigger jets
+ }
}
-/**************************************************************************/
-
AliJetEventParticles::AliJetEventParticles(const AliJetEventParticles& source) :
TObject(source),
fNParticles(source.fNParticles),
- fParticles(new TClonesArray("AliJetEventParticles",source.fNParticles)),
- fVertexX(0.),
- fVertexY(0.),
- fVertexZ(0.)
+ fParticles(new TClonesArray("AliJetParticle",source.fNParticles)),
+ fVertexX(source.GetVertexX()),
+ fVertexY(source.GetVertexY()),
+ fVertexZ(source.GetVertexZ()),
+ fTrials(source.Trials()),
+ fNJets(source.NTriggerJets()),
+ fNUQJets(source.NUQTriggerJets()),
+ fXJet(source.GetXJet()),
+ fYJet(source.GetXJet())
{
//copy constructor
for(Int_t i =0; i<fNParticles; i++)
{
const AliJetParticle *kjp=(const AliJetParticle *)source.fParticles->At(i);
- fParticles->AddAt(new AliJetParticle(*(kjp)),i);
+ new((*fParticles)[i]) AliJetParticle(*(kjp));
+ }
+ for (Int_t i = 0; i < 4; i++) fZquench[i] = 0.;
+ for (Int_t i = 0; i < 10; i++)
+ for (Int_t j = 0; j < 4; j++) {
+ fJets[j][i]=0; // Trigger jets
+ fUQJets[j][i]=0; // Unquenched trigger jets
}
+ source.GetZQuench(fZquench);
+ for (Int_t i = 0; i < NTriggerJets(); i++){
+ source.TriggerJet(i,fJets[0][i],fJets[1][i],fJets[2][i],fJets[3][i]);
+ }
+ for (Int_t i = 0; i < NUQTriggerJets(); i++){
+ source.UQJet(i,fUQJets[0][i],fUQJets[1][i],fUQJets[2][i],fUQJets[3][i]);
+ }
}
-/**************************************************************************/
-
AliJetEventParticles::~AliJetEventParticles()
{
//destructor
delete fParticles;
}
-
-/**************************************************************************/
-
void AliJetEventParticles::Reset(Int_t size)
{
//deletes all particles from the event
-
if(fParticles) fParticles->Clear();
if(size>=0) fParticles->Expand(size);
fNParticles = 0;
-}
-/**************************************************************************/
+ fVertexX=0.;
+ fVertexY=0.;
+ fVertexZ=0.;
+ fTrials=0;
+ fNJets=0;
+ fNUQJets=0;
+ fXJet=-1;
+ fYJet=-1;
+ for (Int_t i = 0; i < 4; i++) fZquench[i] = 0.;
+ for (Int_t i = 0; i < 10; i++)
+ for (Int_t j = 0; j < 4; j++) {
+ fJets[j][i]=0; // Trigger jets
+ fUQJets[j][i]=0; // Unquenched trigger jets
+ }
+}
void AliJetEventParticles::AddParticle(AliJetParticle* part)
{
fParticles->AddAt(part,fNParticles++);
}
-/**************************************************************************/
-
void AliJetEventParticles::AddParticle(const AliJetParticle* part)
{
//Adds new particle to the event
new((*fParticles)[fNParticles++]) AliJetParticle(*part);
}
-/**************************************************************************/
-
void AliJetEventParticles::AddParticle(const TParticle* part,Int_t idx, Int_t l)
{
//Adds new particle to the event
new((*fParticles)[fNParticles++]) AliJetParticle(part,idx,l);
}
-/**************************************************************************/
-
void AliJetEventParticles::AddParticle(Float_t px, Float_t py, Float_t pz,
Float_t etot, Int_t idx, Int_t l)
{
new((*fParticles)[fNParticles++]) AliJetParticle(px,py,pz,etot,idx,l);
}
-/**************************************************************************/
-
void AliJetEventParticles::AddParticle(Float_t px, Float_t py, Float_t pz, Float_t etot, Int_t idx, Int_t l,
Float_t pt, Float_t phi, Float_t eta)
{
new((*fParticles)[fNParticles++]) AliJetParticle(px,py,pz,etot,idx,l,pt,phi,eta);
}
-
-/**************************************************************************/
-
const AliJetParticle* AliJetEventParticles::GetParticleSafely(Int_t n)
{
//returns nth particle with range check
return (const AliJetParticle*)fParticles->At(n);
}
-/**************************************************************************/
+void AliJetEventParticles::AddJet(Float_t px, Float_t py, Float_t pz, Float_t e)
+{
+ //
+ // Add a jet
+ //
+ if (fNJets < 10) {
+ fJets[0][fNJets] = px;
+ fJets[1][fNJets] = py;
+ fJets[2][fNJets] = pz;
+ fJets[3][fNJets] = e;
+ fNJets++;
+ } else {
+ printf("\nWarning: More than 10 jets triggered !!\n");
+ }
+}
+
+void AliJetEventParticles::AddJet(Float_t p[4])
+{
+ //
+ // Add a jet
+ //
+ if (fNJets < 10) {
+ fJets[0][fNJets] = p[0];
+ fJets[1][fNJets] = p[1];
+ fJets[2][fNJets] = p[2];
+ fJets[3][fNJets] = p[3];
+ fNJets++;
+ } else {
+ printf("\nWarning: More than 10 jets triggered !!\n");
+ }
+}
+
+void AliJetEventParticles::AddUQJet(Float_t px, Float_t py, Float_t pz, Float_t e)
+{
+ //
+ // Add a jet
+ //
+ if (fNUQJets < 10) {
+ fUQJets[0][fNUQJets] = px;
+ fUQJets[1][fNUQJets] = py;
+ fUQJets[2][fNUQJets] = pz;
+ fUQJets[3][fNUQJets] = e;
+ fNUQJets++;
+ } else {
+ printf("\nWarning: More than 10 jets triggered !!\n");
+ }
+}
+
+void AliJetEventParticles::AddUQJet(Float_t p[4])
+{
+ //
+ // Add a jet
+ //
+ if (fNUQJets < 10) {
+ fUQJets[0][fNUQJets] = p[0];
+ fUQJets[1][fNUQJets] = p[1];
+ fUQJets[2][fNUQJets] = p[2];
+ fUQJets[3][fNUQJets] = p[3];
+ fNUQJets++;
+ } else {
+ printf("\nWarning: More than 10 jets triggered !!\n");
+ }
+}
+
+void AliJetEventParticles::SetZQuench(Double_t z[4])
+{
+ //
+ // Set quenching fraction
+ //
+ for (Int_t i = 0; i < 4; i++) fZquench[i] = z[i];
+}
+
+void AliJetEventParticles::GetZQuench(Double_t z[4]) const
+{
+ //
+ // Get quenching fraction
+ //
+ for (Int_t i = 0; i < 4; i++) z[i] = fZquench[i];
+}
+
+void AliJetEventParticles::TriggerJet(Int_t i, Float_t p[4]) const
+{
+ //
+ // Give back jet #i
+ //
+ if (i >= fNJets) {
+ printf("\nWarning: TriggerJet, index out of Range!!\n");
+ } else {
+ p[0] = fJets[0][i];
+ p[1] = fJets[1][i];
+ p[2] = fJets[2][i];
+ p[3] = fJets[3][i];
+ }
+}
+
+void AliJetEventParticles::TriggerJet(Int_t i, Float_t &p1, Float_t &p2, Float_t &p3, Float_t &E) const
+{
+ //
+ // Give back jet #i
+ //
+ if (i >= fNJets) {
+ printf("\nWarning: TriggerJet, index out of Range!!\n");
+ } else {
+ p1 = fJets[0][i];
+ p2 = fJets[1][i];
+ p3 = fJets[2][i];
+ E = fJets[3][i];
+ }
+}
+
+void AliJetEventParticles::UQJet(Int_t i, Float_t p[4]) const
+{
+ //
+ // Give back jet #i
+ //
+ if (i >= fNUQJets) {
+ printf("\nWarning: Unquenched Jets, index out of Range!!\n");
+ } else {
+ p[0] = fUQJets[0][i];
+ p[1] = fUQJets[1][i];
+ p[2] = fUQJets[2][i];
+ p[3] = fUQJets[3][i];
+ }
+}
+
+void AliJetEventParticles::UQJet(Int_t i, Float_t &p1, Float_t &p2, Float_t &p3, Float_t &E) const
+{
+ //
+ // Give back jet #i
+ //
+ if (i >= fNUQJets) {
+ printf("\nWarning: Unquenched Jets, index out of Range!!\n");
+ } else {
+ p1 = fUQJets[0][i];
+ p2 = fUQJets[1][i];
+ p3 = fUQJets[2][i];
+ E = fUQJets[3][i];
+ }
+}
+
+void AliJetEventParticles::SetXYJet(Double_t x, Double_t y)
+{
+ //
+ // Add jet production point
+ //
+ fXJet = x;
+ fYJet = y;
+}
void AliJetEventParticles::Print(Option_t* /*t*/) const
{
cout << "--- AliJetEventParticles ---" << endl;
if(fHeader.Length()) cout << fHeader.Data() << endl;
- cout << "no of particles: " << fNParticles << endl;
+ cout << "Particles in Event: " << fNParticles << endl;
+ if(fNUQJets){
+ cout << "Unquenched Jets: " << fNUQJets << endl;
+ for(Int_t i = 0; i<fNUQJets; i++){
+ Float_t x=fUQJets[0][i];
+ Float_t y=fUQJets[1][i];
+ Float_t z=fUQJets[2][i];
+ Float_t e=fUQJets[3][i];
+ Float_t ptj=TMath::Sqrt(x*x+y*y);
+ Float_t thj=TMath::ATan2(ptj,z);
+ Float_t etaj=-TMath::Log(TMath::Tan(thj/2));
+ Float_t phj=TMath::Pi()+TMath::ATan2(-y,-x);
+ Float_t et=e*TMath::Sin(thj);
+ cout << i << " " << et << " " << etaj << " " << phj << endl;
+ }
+ }
+ if(fNJets){
+ cout << "Triggered Jets: " << fNJets << endl;
+ for(Int_t i = 0; i<fNJets; i++){
+ Float_t x=fJets[0][i];
+ Float_t y=fJets[1][i];
+ Float_t z=fJets[2][i];
+ Float_t e=fJets[3][i];
+ Float_t ptj=TMath::Sqrt(x*x+y*y);
+ Float_t thj=TMath::ATan2(ptj,z);
+ Float_t etaj=-TMath::Log(TMath::Tan(thj/2));
+ Float_t phj=TMath::Pi()+TMath::ATan2(-y,-x);
+ Float_t et=e*TMath::Sin(thj);
+ cout << i << " " << et << " " << etaj << " " << phj << endl;
+ }
+ }
+
}
+
Float_t GetVertexY() const {return fVertexY;}
Float_t GetVertexZ() const {return fVertexZ;}
+ Int_t Trials() const {return fTrials;}
+ Int_t NTriggerJets() const {return fNJets;}
+ Int_t NUQTriggerJets() const {return fNUQJets;}
+ void TriggerJet(Int_t i, Float_t p[4]) const;
+ void UQJet(Int_t i, Float_t p[4]) const;
+ void TriggerJet(Int_t i, Float_t &p1, Float_t &p2, Float_t &p3, Float_t &E) const;
+ void UQJet(Int_t i, Float_t &p1, Float_t &p2, Float_t &p3, Float_t &E) const;
+ Double_t GetXJet() const {return fXJet;}
+ Double_t GetYJet() const {return fYJet;}
+ void GetZQuench(Double_t z[4]) const;
+ TString getHeader() const {return fHeader;}
+
+ void SetXYJet(Double_t x, Double_t y);
+ void SetZQuench(Double_t z[4]);
+ void SetTrials(Int_t trials) {fTrials = trials;}
+ void AddJet(Float_t px, Float_t py, Float_t pz, Float_t e);
+ void AddUQJet(Float_t px, Float_t py, Float_t pz, Float_t e);
+ void AddJet(Float_t p[4]);
+ void AddUQJet(Float_t p[4]);
+
void Print(Option_t *t="") const;
protected:
Int_t fNParticles; // number of particles read
TClonesArray *fParticles; //-> particles in event
- Float_t fVertexX; //vertex x
- Float_t fVertexY; //vertex y
- Float_t fVertexZ; //vertex z
+ Float_t fVertexX; // vertex x
+ Float_t fVertexY; // vertex y
+ Float_t fVertexZ; // vertex z
+
+ Int_t fTrials; // Number of trials to fulfill trigger condition
+ Int_t fNJets; // Number of triggered jets
+ Int_t fNUQJets; // Number of unquenched
+ Double_t fXJet; // Jet production point (x)
+ Double_t fYJet; // Jet production point (y)
+ Float_t fJets[4][10]; // Trigger jets
+ Float_t fUQJets[4][10]; // Unquenched trigger jets
+ Double_t fZquench[4]; // Quenching fraction
- ClassDef(AliJetEventParticles,1) //class AliJetEventParticles
+ ClassDef(AliJetEventParticles,2) //class AliJetEventParticles
};
#endif
{
retval = new TString(".");
return *retval;
- }
+ }
if ((entry>fDirs->GetEntries()) || (entry<0))
- //if out of bounds return empty string
- //note that entry==0 is accepted even if array is empty (size=0)
+ //if out of bounds return empty string
+ //note that entry==0 is accepted even if array is empty (size=0)
{
Error("GetDirName","Entry out of bounds");
retval = new TString();
virtual Int_t Next(); //call this if you want the next event
virtual void Rewind() = 0;
- void SetPtCut(Float_t ptmin=0, Float_t ptmax=1000)
+ void SetPtCut(Float_t ptmin=0, Float_t ptmax=100)
{fPtMin=ptmin;fPtMax=ptmax;}
void SetPhiCut(Float_t phi=2*TMath::Pi()){SetPhiCut(0,phi);}
void SetPhiCut(Float_t phimin, Float_t phimax)
#include <TH1.h>
#include <AliESD.h>
#include <AliESDtrack.h>
+#include <AliKalmanTrack.h>
#include <AliJetEventParticles.h>
#include "AliJetParticlesReaderESD.h"
AliJetParticlesReaderESD::~AliJetParticlesReaderESD()
{
//desctructor
- if(fFile) delete fFile;
if(fESD) delete fESD;
if(fTree) delete fTree;
if(fKeyIterator) delete fKeyIterator;
+ if(fFile) delete fFile;
}
fFile = OpenFile(fCurrentDir);
if (fFile == 0)
{
- Error("ReadNext","Cannot get fFile for dir no. %d",fCurrentDir);
+ Error("ReadNext","Can't get fFile for dir no. %d",fCurrentDir);
fCurrentDir++;
continue;
}
if(fTree)
{
+ if(AliKalmanTrack::GetConvConst()<=0.)
+ AliKalmanTrack::SetMagneticField(0.5);
if(fCurrentEvent>=fTree->GetEntries())
{
fCurrentDir++;
fESD = dynamic_cast<AliESD*>(fFile->Get(esdname));
if (fESD == 0)
{
- Info("ReadNext","Can not find AliESD object named %s",esdname.Data());
+ Info("ReadNext","Can't find AliESD object named %s",esdname.Data());
fCurrentDir++;
delete fKeyIterator;
fKeyIterator = 0;
// return kFALSE;
// }
- //Float_t mf = esd->GetMagneticField();
- //if ( (mf == 0.0) && (fNTrackPoints > 0) )
- //{
- // Error("ReadESD","Magnetic Field is 0 and Track Points demanded. Skipping to next event.");
- // return kFALSE;
- //}
+ Float_t mf = esd->GetMagneticField();
+ if (mf <= 0.0)
+ {
+ Error("ReadESD","Magnetic Field is 0. Skipping to next event.");
+ return kFALSE;
+ }
+ AliKalmanTrack::SetMagneticField(mf/10.);
Info("ReadESD","Reading Event %d",fCurrentEvent);
if((!fOwner) || (fEventParticles==0))
fEventParticles = new AliJetEventParticles();
+ const Int_t kntr = esd->GetNumberOfTracks();
+ Info("ReadESD","Found %d tracks.",kntr);
+ fEventParticles->Reset(kntr);
+
+ TString headdesc="";
+ headdesc+="Run ";
+ headdesc+=esd->GetRunNumber();
+ headdesc+=": Ev ";
+ headdesc+=esd->GetEventNumber();
+ fEventParticles->SetHeader(headdesc);
+
Double_t vertexpos[3];//vertex position
const AliESDVertex* kvertex = esd->GetVertex();
if (kvertex == 0)
}
fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
- const Int_t kntr = esd->GetNumberOfTracks();
- Info("ReadESD","Found %d tracks.",kntr);
+ //loop over tracks
for (Int_t i = 0;i<kntr; i++)
{
const AliESDtrack *kesdtrack = esd->GetTrack(i);
if (kesdtrack == 0)
{
- Error("ReadESD","Can not get track %d", i);
+ Error("ReadESD","Can't get track %d", i);
continue;
}
const TString& dirname = GetDirName(n);
if (dirname == "")
{
- Error("OpenFiles","Can not get directory name");
+ Error("OpenFiles","Can't get directory name");
return 0;
}
TString filename = dirname +"/"+ fESDFileName;
#include <AliESD.h>
#include <AliESDtrack.h>
#include <AliESDHLTtrack.h>
+#include <AliKalmanTrack.h>
#include <AliJetEventParticles.h>
#include "AliJetParticlesReaderHLT.h"
AliJetParticlesReaderHLT::AliJetParticlesReaderHLT(Bool_t bMapper, const Char_t* esdfilename) :
AliJetParticlesReaderESD(esdfilename),
- fTrackerType(bMapper)
+ fTrackerType(bMapper),
+ fMinHits(0),
+ fMinWeight(0)
{
//constructor
}
TObjArray* dirs,
const Char_t* esdfilename) :
AliJetParticlesReaderESD(dirs,esdfilename),
- fTrackerType(bMapper)
+ fTrackerType(bMapper),
+ fMinHits(0),
+ fMinWeight(0)
+
{
//constructor
}
return kFALSE;
}
+ Float_t mf = esd->GetMagneticField();
+ if (mf <= 0.0)
+ {
+ Error("ReadESD","Magnetic Field is 0. Skipping to next event.");
+ return kFALSE;
+ }
+ AliKalmanTrack::SetMagneticField(mf/10.);
+
Info("ReadESD","Reading Event %d",fCurrentEvent);
if((!fOwner) || (fEventParticles==0))
fEventParticles = new AliJetEventParticles();
+ Int_t ntr=0;
+ if(fTrackerType){
+ ntr =esd->GetNumberOfHLTHoughTracks();
+ Info("ReadESD","Found %d conformal tracks.",ntr);
+ } else {
+ ntr=esd->GetNumberOfHLTHoughTracks();
+ Info("ReadESD","Found %d hough tracks.",ntr);
+ }
+ fEventParticles->Reset(ntr);
+
+ TString headdesc="";
+ headdesc+="Run ";
+ headdesc+=esd->GetRunNumber();
+ headdesc+=": Ev ";
+ headdesc+=esd->GetEventNumber();
+ fEventParticles->SetHeader(headdesc);
+
Double_t vertexpos[3];//vertex position
const AliESDVertex* kvertex = esd->GetVertex();
if (kvertex == 0)
fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
- if(fTrackerType){
- const Int_t kntr =esd->GetNumberOfHLTHoughTracks();
- Info("ReadESD","Found %d conformal tracks.",kntr);
- for (Int_t i = 0;i<kntr; i++) {
- const AliESDHLTtrack *kesdtrack=esd->GetHLTConfMapTrack(i);
+ for (Int_t i = 0;i<ntr; i++) {
+ AliESDHLTtrack *kesdtrack;
+ if(fTrackerType){
+ kesdtrack=esd->GetHLTConfMapTrack(i);
+ } else {
+ kesdtrack=esd->GetHLTHoughTrack(i);
+ }
- if (kesdtrack == 0)
+ if (kesdtrack == 0)
{
Error("ReadESD","Can not get track %d", i);
continue;
}
- //const Float_t kpid=kesdtrack->GetPID();
- //const Int_t knhits=kesdtrack->GetNHits();
- const Float_t kpx=kesdtrack->GetPx();
- const Float_t kpy=kesdtrack->GetPy();
- const Float_t kpz=kesdtrack->GetPz();
- const Float_t kpt=kesdtrack->GetPt();
- const Float_t kp=TMath::Sqrt(kpz*kpz+kpt*kpt);
- const Float_t keta=0.5*TMath::Log((kp+kpz+1e-30)/(kp-kpz+1e-30));
- const Float_t kphi=TMath::Pi()+TMath::ATan2(-kpy,-kpx);
+ //const Float_t kpid=kesdtrack->GetPID();
+ const Int_t knhits=kesdtrack->GetNHits();
+ const Int_t kweight=kesdtrack->GetWeight();
+ if((fMinHits>0) && (knhits<fMinHits)) continue;
+ if((fMinWeight>0) && (kweight<fMinWeight)) continue;
- if(IsAcceptedParticle(kpt,kphi,keta))
- fEventParticles->AddParticle(kpx,kpy,kpz,kp,i,kesdtrack->GetMCid(),kpt,kphi,keta);
+ const Float_t kpx=kesdtrack->GetPx();
+ const Float_t kpy=kesdtrack->GetPy();
+ const Float_t kpz=kesdtrack->GetPz();
+ const Float_t kpt=kesdtrack->GetPt();
+ const Float_t kp=TMath::Sqrt(kpz*kpz+kpt*kpt);
+ const Float_t keta=0.5*TMath::Log((kp+kpz+1e-30)/(kp-kpz+1e-30));
+ const Float_t kphi=TMath::Pi()+TMath::ATan2(-kpy,-kpx);
- } //loop over conf tracks
- } else {
- const Int_t kntr=esd->GetNumberOfHLTHoughTracks();
- Info("ReadESD","Found %d hough tracks.",kntr);
-
- for (Int_t i = 0;i<kntr; i++) {
- const AliESDHLTtrack *kesdtrack=esd->GetHLTHoughTrack(i);
-
- if (kesdtrack == 0)
- {
- Error("ReadESD","Can not get track %d", i);
- continue;
- }
-
- //const Float_t kweight=kesdtrack->GetWeight();
- const Float_t kpx=kesdtrack->GetPx();
- const Float_t kpy=kesdtrack->GetPy();
- const Float_t kpz=kesdtrack->GetPz();
+ if(IsAcceptedParticle(kpt,kphi,keta))
+ fEventParticles->AddParticle(kpx,kpy,kpz,kp,i,kesdtrack->GetMCid(),kpt,kphi,keta);
+ } //loop over tracks
- const Float_t kpt=kesdtrack->GetPt();
- const Float_t kp=TMath::Sqrt(kpz*kpz+kpt*kpt);
- const Float_t keta=0.5*TMath::Log((kp+kpz+1e-30)/(kp-kpz+1e-30));
- const Float_t kphi=TMath::Pi()+TMath::ATan2(-kpy,-kpx);
-
- if(IsAcceptedParticle(kpt,kphi,keta))
- fEventParticles->AddParticle(kpx,kpy,kpz,kp,i,kesdtrack->GetMCid(),kpt,kphi,keta);
- } //loop over hough tracks
- }
return kTRUE;
}
virtual ~AliJetParticlesReaderHLT();
+ void SetMinHits(Int_t i){fMinHits=i;}
+ void SetMinWeight(Int_t i){fMinWeight=i;}
+
protected:
Int_t ReadESD(AliESD* esd); //read esd file/objects
Bool_t fTrackerType; //track type: 0==Conformal Tracker, 1==Hough Tracker
+ Int_t fMinHits; //minimum hits required
+ Int_t fMinWeight; //minimum weight required
ClassDef(AliJetParticlesReaderHLT,1) //
};
#include "AliJetParticle.h"
#include "AliJetEventParticles.h"
#include "AliJetParticlesReaderKine.h"
+#include "AliHeader.h"
+#include "AliGenPythiaEventHeader.h"
ClassImp(AliJetParticlesReaderKine)
-
-/**********************************************************/
-
AliJetParticlesReaderKine::AliJetParticlesReaderKine() :
AliJetParticlesReader(),
fFileName("galice.root"),
//constructor
}
-/**********************************************************/
-
-
AliJetParticlesReaderKine::AliJetParticlesReaderKine(TString& fname) :
AliJetParticlesReader(),
fFileName(fname),
//constructor
}
-/**********************************************************/
-
AliJetParticlesReaderKine::AliJetParticlesReaderKine(TObjArray* dirs, const Char_t *filename):
AliJetParticlesReader(dirs),
fFileName(filename),
//constructor
}
-/**********************************************************/
-
AliJetParticlesReaderKine::~AliJetParticlesReaderKine()
{
//destructor
if(fRunLoader) delete fRunLoader;
}
-/**********************************************************/
-
void AliJetParticlesReaderKine::Rewind()
{
//Rewinds to the beginning
if(fRunLoader) delete fRunLoader;
- fRunLoader = 0;
- fCurrentDir = 0;
- fNEventsRead= 0;
+ fRunLoader = 0;
+ fCurrentDir = 0;
+ fNEventsRead = 0;
}
-/**********************************************************/
-
Int_t AliJetParticlesReaderKine::ReadNext()
{
//Reads Kinematics Tree
while(fCurrentDir < GetNumberOfDirs())
{
- if (!OpenFile(fCurrentDir))
+
+ if (!OpenFile(fCurrentDir))
{
+ delete fRunLoader; //close current session
+ fRunLoader = 0; //assure pointer is null
fCurrentDir++;
continue;
}
continue;
}
- //Info("ReadNext","Reading Event %d",fCurrentEvent);
-
+ Info("ReadNext","Reading Event %d",fCurrentEvent);
fRunLoader->GetEvent(fCurrentEvent);
AliStack* stack = fRunLoader->Stack();
if (!stack)
{
- Error("ReadNext","Can not get stack for event %d",fCurrentEvent);
+ Error("ReadNext","Can't get stack for event %d",fCurrentEvent);
continue;
}
+ //clear old event
+ Int_t nprim = stack->GetNprimary();
+ Int_t npart = nprim;
+ if(fUseTracks)
+ npart = stack->GetNtrack();
+ fEventParticles->Reset(npart);
+
+ TString headdesc="";
+ AliHeader *header=fRunLoader->GetHeader();
+ if(!header) {
+ Warning("ReadNext","Header not found in event %d",fCurrentEvent);
+ } else {
+ headdesc+="Run ";
+ headdesc+=header->GetRun();
+ headdesc+=": Ev ";
+ headdesc+=header->GetEventNrInRun();
+ AliGenPythiaEventHeader *pheader=(AliGenPythiaEventHeader*)header->GenEventHeader();
+ if(!pheader) {
+ Warning("ReadNext","Pythia-Header not found in event %d",fCurrentEvent);
+ } else {
+ Int_t ntruq=pheader->NUQTriggerJets();
+ if(ntruq){
+ Double_t x0=pheader->GetXJet();
+ Double_t y0=pheader->GetYJet();
+ Double_t zquench[4];
+ pheader->GetZQuench(zquench);
+ fEventParticles->SetXYJet(x0,y0);
+ fEventParticles->SetZQuench(zquench);
+ for(Int_t j=0;j<ntruq;j++){
+ Float_t pjet[4];
+ pheader->UQJet(j,pjet);
+ fEventParticles->AddUQJet(pjet);
+ }
+ }
+ //Int_t ptyp=pheader->ProcessType();
+ Int_t ntrials=pheader->Trials();
+ headdesc+=": Tr ";
+ headdesc+=ntrials;
+ Int_t ntr=pheader->NTriggerJets();
+ if(ntr){
+ for(Int_t j=0;j<ntr;j++){
+ Float_t pjet[4];
+ pheader->TriggerJet(j,pjet);
+ fEventParticles->AddJet(pjet);
+ if(!ntruq) fEventParticles->AddUQJet(pjet);
+ }
+ }
+ }
+ }
+ fEventParticles->SetHeader(headdesc);
+
//get vertex
const TParticle *kv = stack->Particle(0);
if(kv) {
fEventParticles->SetVertex(kv->Vx(),kv->Vy(),kv->Vz());
}
- Int_t nprim = stack->GetNprimary();
- Int_t npart = nprim;
- if(fUseTracks)
- npart = stack->GetNtrack();
-
- fEventParticles->Reset(npart);
+ //loop over particles
for (Int_t i = 0;i<npart; i++)
{
TParticle *p = stack->Particle(i);
Int_t child1 = p->GetFirstDaughter();
//Int_t child2 = p->GetLastDaughter();
//Int_t mother = p->GetFirstMother();
- if((child1>=0) && (child1<nprim)) continue;
//cout << child1 << " " << child2 << " " << mother << endl;
-
+ if((child1>=0) && (child1<nprim)) continue;
+ if(p->GetStatusCode()!=1){
+ //p->Print();
+ continue;
+ }
if(IsAcceptedParticle(p)) //put particle in event
fEventParticles->AddParticle(p,i);
}
fCurrentEvent++;
- fCurrentDir++;
fNEventsRead++;
return kTRUE;
}
- //end of loop over directories specified in fDirs Obj Array
+
+ //end of loop over directories specified in fDirs ObjArray
return kFALSE;
}
-/**********************************************************/
-
Int_t AliJetParticlesReaderKine::OpenFile(Int_t n)
{
//opens file with kine tree
+ if(fRunLoader){
+ if(fCurrentEvent < fRunLoader->GetNumberOfEvents()) return kTRUE;
+ else return kFALSE;
+ }
+
const TString& dirname = GetDirName(n);
if (dirname == "")
{
- Error("OpenNextFile","Can not get directory name with index %d",n);
+ Error("OpenNextFile","Can't get directory name with index %d",n);
return kFALSE;
}
TString filename = dirname +"/"+ fFileName;
- if(fRunLoader) delete fRunLoader;
fRunLoader = AliRunLoader::Open(filename.Data());
if (fRunLoader == 0)
fRunLoader = 0;
return kFALSE;
}
-
+
if (fRunLoader->LoadKinematics())
{
Error("OpenNextFile","Error occured while loading kinematics.");
return kFALSE;
}
-
+
fCurrentEvent = 0;
return kTRUE;
}
inline Bool_t AliJetParticlesReaderKine::IsAcceptedParticle(TParticle *p) const
{
- //p->Print();
-
Int_t pcode=p->GetPdgCode();
if ((pcode==11)||(pcode==-11)||(pcode==22)) {
if(!fEM) return kFALSE;
DHDR= JetAnalysisLinkDef.h
-EINCLUDE:= STEER ANALYSIS
+EINCLUDE:= STEER ANALYSIS PYTHIA6
git://git.uio.no / u/mrichter/AliRoot.git / commitdiff