#include "AliJetAODReaderHeader.h"
#include "AliAODEvent.h"
#include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
#include "AliJetDummyGeo.h"
#include "AliJetAODFillUnitArrayTracks.h"
#include "AliJetAODFillUnitArrayEMCalDigits.h"
-#include "AliJetHadronCorrection.h"
+#include "AliJetHadronCorrectionv1.h"
#include "AliJetUnitArray.h"
ClassImp(AliJetAODReader)
fGrid2(0),
fGrid3(0),
fGrid4(0),
- fPtCut(0),
fApplyElectronCorrection(kFALSE),
fApplyMIPCorrection(kTRUE),
fApplyFractionHadronicCorrection(kFALSE),
if (a>=naod) continue;
if (strstr(name,pattern)){
- char path[256];
- sprintf(path,"%s/%s/aod.root",dirName,name);
- fChain->AddFile(path);
+ fChain->AddFile(Form("%s/%s/aod.root",dirName,name));
a++;
}
}
}
}
+
+Bool_t AliJetAODReader::AcceptAODMCParticle(AliAODMCParticle *mcP,Short_t flag){
+
+ //
+ // filter for charge and for charged and neutral, no detector
+ // response yet
+ // physical priamries already selected
+
+ Int_t pdg = TMath::Abs(mcP->GetPdgCode());
+
+ // exclude neutrinos anyway
+ if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
+
+ if(flag==AliJetAODReaderHeader::kAllMC)return kTRUE;
+ if(flag==AliJetAODReaderHeader::kChargedMC){
+ if(mcP->Charge()==0)return kFALSE;
+ return kTRUE;
+ }
+
+ return kTRUE;
+
+}
+
+
+Bool_t AliJetAODReader::FillMomentumArrayMC(){
+
+ //
+ // This routine fetches the MC particles from the AOD
+ // Depending on the flag all particles except neurinos are use
+ // or only charged particles
+ //
+
+ TClonesArray *mcArray = (TClonesArray*)fAOD->FindListObject(AliAODMCParticle::StdBranchName());
+ if(!mcArray){
+ Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
+ return kFALSE;
+ }
+
+ Int_t nMC = mcArray->GetEntriesFast();
+
+ // get number of tracks in event (for the loop)
+ if(fDebug>0)printf("AOD MC tracks: %5d \t", nMC);
+
+ // temporary storage of signal and pt cut flag
+ Int_t* sflag = new Int_t[nMC];
+ Int_t* cflag = new Int_t[nMC];
+
+ // get cuts set by user
+ Float_t ptMin = fReaderHeader->GetPtCut();
+ Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
+ Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
+ Int_t mcTrack = 0;
+ Float_t pt, eta;
+ TVector3 p3;
+
+ Short_t readerFlag = ((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC();
+
+
+ for (Int_t it = 0; it < nMC; it++) {
+ AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
+ if(!track->IsPhysicalPrimary())continue;
+
+ p3.SetXYZ(track->Px(),track->Py(),track->Pz());
+ eta = p3.Eta();
+ if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
+ if(!AcceptAODMCParticle(track,readerFlag))continue;
+ pt = p3.Pt();
+
+
+
+
+ if (mcTrack == 0){
+ fRef->Delete(); // make sure to delete before placement new...
+ new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
+ }
+ new ((*fMomentumArray)[mcTrack]) TLorentzVector(p3,p3.Mag());
+ sflag[mcTrack] = 1;
+ cflag[mcTrack] = ( pt > ptMin ) ? 1: 0;
+ fRef->Add(track);
+ mcTrack++;
+ }
+ if(fDebug>0)printf("Used MC tracks: %5d \n", mcTrack);
+ // set the signal flags
+ fSignalFlag.Set(mcTrack,sflag);
+ fCutFlag.Set(mcTrack,cflag);
+
+ delete [] sflag;
+ delete [] cflag;
+
+ return kTRUE;
+}
+
//____________________________________________________________________________
Bool_t AliJetAODReader::FillMomentumArray()
ClearArray();
fRef->Clear();
fDebug = fReaderHeader->GetDebug();
-
+
if (!fAOD) {
return kFALSE;
}
-
+
+ if(((AliJetAODReaderHeader*)fReaderHeader)->GetReadAODMC()){
+ return FillMomentumArrayMC();
+ }
+
// get number of tracks in event (for the loop)
Int_t nt = fAOD->GetNTracks();
- printf("AOD tracks: %5d \t", nt);
+ if(fDebug>0)printf("AOD tracks: %5d \t", nt);
// temporary storage of signal and pt cut flag
Int_t* sflag = new Int_t[nt];
if((filterMask>0)&&!(track->TestFilterBit(filterMask)))continue;
if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
+ if (aodTrack == 0){
+ fRef->Delete(); // make sure to delete before placement new...
+ new(fRef) TRefArray(TProcessID::GetProcessWithUID(track));
+ }
new ((*fMomentumArray)[aodTrack]) TLorentzVector(p3,p3.Mag());
sflag[aodTrack] = (TMath::Abs(track->GetLabel()) < 10000) ? 1 : 0;
cflag[aodTrack] = ( pt > ptMin ) ? 1: 0;
aodTrack++;
fRef->Add(track);
}
- printf("Used AOD tracks: %5d \n", aodTrack);
+ if(fDebug>0)printf("Used AOD tracks: %5d \n", aodTrack);
// set the signal flags
fSignalFlag.Set(aodTrack,sflag);
fCutFlag.Set(aodTrack,cflag);
//____________________________________________________________________________
void AliJetAODReader::CreateTasks(TChain* tree)
{
+ //
+ // For reader task initialization
+ //
+
fDebug = fReaderHeader->GetDebug();
fDZ = fReaderHeader->GetDZ();
fTree = tree;
fFillUAFromEMCalDigits->SetApplyFractionHadronicCorrection(fApplyFractionHadronicCorrection);
fFillUAFromEMCalDigits->SetFractionHadronicCorrection(fFractionHadronicCorrection);
fFillUAFromEMCalDigits->SetApplyElectronCorrection(fApplyElectronCorrection);
-
+ // Add the task to global task
fFillUnitArray->Add(fFillUAFromTracks);
fFillUnitArray->Add(fFillUAFromEMCalDigits);
fFillUAFromTracks->SetActive(kFALSE);
// clear momentum array
ClearArray();
+ fRef->Clear();
fDebug = fReaderHeader->GetDebug();
fOpt = fReaderHeader->GetDetector();
fFillUAFromTracks->SetActive(kTRUE);
fFillUAFromTracks->SetUnitArray(fUnitArray);
fFillUAFromTracks->SetRefArray(fRefArray);
+ fFillUAFromTracks->SetReferences(fRef);
+ fFillUAFromTracks->SetSignalFlag(fSignalFlag);
+ fFillUAFromTracks->SetCutFlag(fCutFlag);
fFillUAFromTracks->SetProcId(fProcId);
// fFillUAFromTracks->ExecuteTask("tpc"); // => Temporarily changed
fFillUAFromTracks->Exec("tpc");
fNumCandidate = fFillUAFromTracks->GetMult();
fNumCandidateCut = fFillUAFromTracks->GetMultCut();
}
+ fSignalFlag = fFillUAFromTracks->GetSignalFlag();
+ fCutFlag = fFillUAFromTracks->GetCutFlag();
}
// Digits only or Digits+TPC
{
// Initialise parameters
fOpt = fReaderHeader->GetDetector();
- // fHCorrection = 0; // For hadron correction
- fHadCorr = 0; // For hadron correction
- if(fEFlag==kFALSE){
- if(fOpt==0 || fOpt==1)
- fECorrection = 0; // For electron correction
- else fECorrection = 1; // For electron correction
- }
fNumUnits = fGeom->GetNCells(); // Number of cells in EMCAL
if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
}
} // end grid type == 1
fArrayInitialised = 1;
}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff