#include <AliRunLoader.h>
#include <AliStack.h>
#include <AliESDtrack.h>
+#include <AliKalmanTrack.h>
#include <AliESD.h>
#include "AliAnalysis.h"
#include "AliAODRun.h"
#include "AliAOD.h"
-#include "AliAODStdParticle.h"
+#include "AliAODParticle.h"
#include "AliAODParticleCut.h"
#include "AliTrackPoints.h"
#include "AliClusterMap.h"
fNTrackPoints(0),
fdR(0.0),
fClusterMap(kFALSE),
+ fITSTrackPoints(kFALSE),
+ fMustTPC(kFALSE),
fNTPCClustMin(0),
fNTPCClustMax(150),
fTPCChi2PerClustMin(0.0),
fNTrackPoints(0),
fdR(0.0),
fClusterMap(kFALSE),
+ fITSTrackPoints(kFALSE),
+ fMustTPC(kFALSE),
fNTPCClustMin(0),
fNTPCClustMax(150),
fTPCChi2PerClustMin(0.0),
//reads next event from fFile
//fRunLoader is for reading Kine
- if (AliAODParticle::GetDebug())
+ if (AliVAODParticle::GetDebug())
Info("ReadNext","Entered");
if (fEventSim == 0x0) fEventSim = new AliAOD();
TKey* key = (TKey*)fKeyIterator->Next();
if (key == 0x0)
{
- if (AliAODParticle::GetDebug() > 2 )
+ if (AliVAODParticle::GetDebug() > 2 )
{
Info("ReadNext","No more keys.");
}
// TObject* esdobj = key->ReadObj();
// if (esdobj == 0x0)
// {
-// if (AliAODParticle::GetDebug() > 2 )
+// if (AliVAODParticle::GetDebug() > 2 )
// {
// Info("ReadNext","Key read NULL. Key Name is %s",key->GetName());
// key->Dump();
AliESD* esd = dynamic_cast<AliESD*>(fFile->Get(esdname));
if (esd == 0x0)
{
-// if (AliAODParticle::GetDebug() > 2 )
+// if (AliVAODParticle::GetDebug() > 2 )
// {
// Info("ReadNext","This key is not an AliESD object %s",key->GetName());
// }
- if (AliAODParticle::GetDebug() > 2 )
+ if (AliVAODParticle::GetDebug() > 2 )
{
Info("ReadNext","Can not find AliESD object named %s",esdname.Data());
}
Error("ReadESD","ESD is NULL");
return 1;
}
-
+
TDatabasePDG* pdgdb = TDatabasePDG::Instance();
if (pdgdb == 0x0)
{
Float_t mf = esd->GetMagneticField();
- if ( (mf == 0.0) && (fNTrackPoints > 0) )
+ if ( (mf == 0.0) && ((fNTrackPoints > 0) || fITSTrackPoints) )
{
Error("ReadESD","Magnetic Field is 0 and Track Points Demended. Skipping to next event.");
return 1;
}
+ if (fITSTrackPoints)
+ {
+ Info("ReadESD","Magnetic Field is %f",mf/10.);
+ AliKalmanTrack::SetMagneticField(mf/10.);
+ }
+
AliStack* stack = 0x0;
if (fReadSim && fRunLoader)
{
vertex->GetXYZ(vertexpos);
}
- if (AliAODParticle::GetDebug() > 0)
+ if (AliVAODParticle::GetDebug() > 0)
{
Info("ReadESD","Primary Vertex is (%f,%f,%f)",vertexpos[0],vertexpos[1],vertexpos[2]);
}
//if (esdtrack->HasVertexParameters() == kFALSE)
if ((esdtrack->GetStatus() & AliESDtrack::kITSrefit) == kFALSE)
{
- if (AliAODParticle::GetDebug() > 2)
+ if (AliVAODParticle::GetDebug() > 2)
Info("ReadNext","Particle skipped: Data at vertex not available.");
continue;
}
+ if (fMustTPC)
+ {
+ if ((esdtrack->GetStatus() & AliESDtrack::kTPCin) == kFALSE)
+ {
+ if (AliVAODParticle::GetDebug() > 2)
+ Info("ReadNext","Particle skipped: Was not reconstructed in TPC.");
+ continue;
+ }
+ }
+
if ((esdtrack->GetStatus() & AliESDtrack::kESDpid) == kFALSE)
{
- if (AliAODParticle::GetDebug() > 2)
+ if (AliVAODParticle::GetDebug() > 2)
Info("ReadNext","Particle skipped: PID BIT is not set.");
continue;
}
esdtrack->GetConstrainedExternalParameters(extx,extp);
if (extp[4] == 0.0)
{
- if (AliAODParticle::GetDebug() > 2)
+ if (AliVAODParticle::GetDebug() > 2)
Info("ReadNext","Track has 0 contrianed curvature -> Probobly parameters never updated. Skipping.");
continue;
}
Int_t charge = (extp[4] > 0)?1:-1;//if curvature=charg/Pt is positive charge is positive
//Particle from kinematics
- AliAODStdParticle* particle = 0;
+ AliAODParticle* particle = 0;
Bool_t keeppart = kFALSE;
if ( fReadSim && stack )
{
Error("ReadNext","Can not find track with such label.");
continue;
}
- if(Pass(p->GetPdgCode()))
+
+ if (fCheckParticlePID)
{
- if ( AliAODParticle::GetDebug() > 5 )
- Info("ReadNext","Simulated Particle PID (%d) did not pass the cut.",p->GetPdgCode());
- continue; //check if we are intersted with particles of this type
- }
+ if(Rejected(p->GetPdgCode()))
+ {
+ if ( AliVAODParticle::GetDebug() > 5 )
+ Info("ReadNext","Simulated Particle PID (%d) did not pass the cut.",p->GetPdgCode());
+ continue; //check if we are intersted with particles of this type
+ }
+ }
// if(p->GetPdgCode()<0) charge = -1;
- particle = new AliAODStdParticle(*p,i);
+ particle = new AliAODParticle(*p,i);
}
for (Int_t s=0; s<AliESDtrack::kSPECIES; s++) rc+=concentr[s]*pidtable[s];
if (rc==0.0)
{
- if (AliAODParticle::GetDebug() > 2)
+ if (AliVAODParticle::GetDebug() > 2)
Info("ReadNext","Particle rejected since total bayessian PID probab. is zero.");
continue;
}
for (Int_t s=0; s<AliESDtrack::kSPECIES; s++) w[s]=concentr[s]*pidtable[s]/rc;
- if (AliAODParticle::GetDebug() > 4)
+ if (AliVAODParticle::GetDebug() > 4)
{
Info("ReadNext","###########################################################################");
Info("ReadNext","Momentum: %f %f %f",mom[0],mom[1],mom[2]);
msg+=")";
}
Info("ReadNext","%s",msg.Data());
- }//if (AliAODParticle::GetDebug()>4)
+ }//if (AliVAODParticle::GetDebug()>4)
AliTrackPoints* tpts = 0x0;
if (fNTrackPoints > 0)
{
tpts = new AliTrackPoints(fNTrackPoints,esdtrack,mf,fdR);
+// tpts->Move(-vertexpos[0],-vertexpos[1],-vertexpos[2]);
+ }
+
+ AliTrackPoints* itstpts = 0x0;
+ if (fITSTrackPoints)
+ {
+ itstpts = new AliTrackPoints(AliTrackPoints::kITS,esdtrack);
+// itstpts->Move(-vertexpos[0],-vertexpos[1],-vertexpos[2]);
}
AliClusterMap* cmap = 0x0;
Float_t pp = w[s];
if (pp == 0.0)
{
- if ( AliAODParticle::GetDebug() > 5 )
+ if ( AliVAODParticle::GetDebug() > 5 )
Info("ReadNext","Probability of being PID %d is zero. Continuing.",pdgcode);
continue;
}
- if(Pass(pdgcode))
+ if(Rejected(pdgcode))
{
- if ( AliAODParticle::GetDebug() > 5 )
+ if ( AliVAODParticle::GetDebug() > 5 )
Info("ReadNext","PID (%d) did not pass the cut.",pdgcode);
continue; //check if we are intersted with particles of this type
}
Double_t mass = pdgdb->GetParticle(pdgcode)->Mass();
Double_t tEtot = TMath::Sqrt( mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2] + mass*mass);//total energy of the track
- AliAODStdParticle* track = new AliAODStdParticle(pdgcode, w[s],i,
+ AliAODParticle* track = new AliAODParticle(pdgcode, w[s],i,
mom[0], mom[1], mom[2], tEtot,
pos[0], pos[1], pos[2], 0.);
//copy probabilitis of other species (if not zero)
track->SetPIDprobability(charge*GetSpeciesPdgCode( (ESpecies)k ),w[k]);
}
- if(Pass(track))//check if meets all criteria of any of our cuts
+ if(Rejected(track))//check if meets all criteria of any of our cuts
//if it does not delete it and take next good track
{
- if ( AliAODParticle::GetDebug() > 4 )
+ if ( AliVAODParticle::GetDebug() > 4 )
Info("ReadNext","Track did not pass the cut");
delete track;
continue;
//Single Particle cuts on cluster map and track points rather do not have sense
if (tpts)
{
- track->SetTrackPoints(tpts);
+ track->SetTPCTrackPoints(tpts);
+ }
+
+ if (itstpts)
+ {
+ track->SetITSTrackPoints(itstpts);
}
if (cmap)
if (particle) fEventSim->AddParticle(particle);
keeppart = kTRUE;
- if (AliAODParticle::GetDebug() > 4 )
+ if (AliVAODParticle::GetDebug() > 4 )
{
Info("ReadNext","\n\nAdding Particle with incarnation %d",pdgcode);
track->Print();
{
delete particle;//particle was not stored in event
delete tpts;
+ delete itstpts;
delete cmap;
}
+ else
+ {
+ if ( fReadSim && stack )
+ {
+ if (particle->P() < 0.00001)
+ {
+ Info("ReadNext","###################################");
+ Info("ReadNext","###################################");
+ Info("ReadNext","Track Label %d",esdtrack->GetLabel());
+ TParticle *p = stack->Particle(esdtrack->GetLabel());
+ Info("ReadNext","");
+ p->Print();
+ Info("ReadNext","");
+ particle->Print();
+ }
+ }
+ }
}//for (Int_t i = 0;i<ntr; i++) -- loop over tracks
fEventRec->GetNumberOfParticles(), fEventSim->GetNumberOfParticles(),
fNEventsRead,fCurrentEvent,fCurrentDir);
fTrackCounter->Fill(fEventRec->GetNumberOfParticles());
+
+ /******************************************************/
+ /****** Setting glevet properties *************/
+ /******************************************************/
+
+ if (fEventRec->GetNumberOfParticles() > 0)
+ {
+ fEventRec->SetPrimaryVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
+ }
return 0;
}