#include <TRandom.h>
// From AliRoot ...
#include "AliAODJet.h"
+#include "AliPDG.h"
#include "AliJetKineReaderHeader.h"
#include "AliJetKineReader.h"
#include "AliMCEventHandler.h"
AliJetKineReader::AliJetKineReader():
AliJetReader(),
fAliHeader(0),
- fMCEventHandler(0),
+ fMCEvent(0),
fGenJets(0)
{
// Default constructor
AliJetKineReader::~AliJetKineReader()
{
// Destructor
- delete fAliHeader;
+ if (fAliHeader) {
+ delete fAliHeader;
+ fAliHeader = 0;
+ }
}
//____________________________________________________________________________
ClearArray();
// Get the stack
- AliStack* stack = fMCEventHandler->MCEvent()->Stack();
+ AliStack* stack = fMCEvent->Stack();
// Number of primaries
Int_t nt = stack->GetNprimary();
Double_t ptMin = ((AliJetKineReaderHeader*) fReaderHeader)->GetPtCut();
Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
- fAliHeader = fMCEventHandler->MCEvent()->Header();
+ fAliHeader = fMCEvent->Header();
TLorentzVector p4;
Float_t p0 = p;
Float_t eta = part->Eta();
Float_t phi = part->Phi();
+
+
+ if (((AliJetKineReaderHeader*)fReaderHeader)->ChargedOnly()) {
+ // Charged particles only
+ Float_t charge =
+ TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
+ if (charge == 0) continue;
+ } // End charged only
+
// Fast simulation of TPC if requested
if (((AliJetKineReaderHeader*)fReaderHeader)->FastSimTPC()) {
Float_t px = r * part->Px();
Float_t py = r * part->Py();
Float_t pz = r * part->Pz();
- Float_t m = part->GetMass();
+ TParticlePDG *pPDG = part->GetPDG();
+ Float_t m = 0;
+ if(!pPDG){
+ // this is very rare...
+ // Is avoided by AliPDG::AddParticlesToPdgDataBase();
+ // but this should be called only once... (how in proof?)
+ // Calucluate mass with unsmeared momentum values
+ m = part->Energy()*part->Energy() -
+ (px * px + py * py + pz * pz)/r/r;
+ if(m>0)m = TMath::Sqrt(m);
+ else m = 0;
+ AliInfo(Form("Unknown Particle using %d calculated mass m = %3.3f",part->GetPdgCode(),m));
+
+ }
+ else{
+ m = pPDG->Mass();
+ }
Float_t e = TMath::Sqrt(px * px + py * py + pz * pz + m * m);
p4 = TLorentzVector(px, py, pz, e);
if ( (p4.Eta()>etaMax) || (p4.Eta()<etaMin)) continue;
// set the signal flags
fSignalFlag.Set(goodTrack, sflag);
fCutFlag.Set(goodTrack, cflag);
- AliInfo(Form("\n Number of good tracks %d \n", goodTrack));
+ AliInfo(Form(" Number of good tracks %d \n", goodTrack));
delete[] sflag;
delete[] cflag;
return fGenJets;
}
-void AliJetKineReader::SetInputEvent(TObject* /*esd*/, TObject* /*aod*/, TObject* mc)
+void AliJetKineReader::SetInputEvent(const TObject* /*esd*/, const TObject* /*aod*/, const TObject* mc)
{
// Connect the input event
- fMCEventHandler = (AliMCEventHandler*) mc;
+ fMCEvent = (AliMCEvent*) mc;
}