X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PWGJE%2FAliPWG4HighPtSpectra.cxx;h=5d440599280679e0c1cad3b56394a1d06a945435;hb=cba7b46e9dd5f74cc7faeae6d3d7a399c31127ba;hp=631be43cf315568651f0c495d85edcc0dc84e1e7;hpb=adf3920d00ae9a524206d0973b4e67b04e336b2d;p=u%2Fmrichter%2FAliRoot.git

diff --git a/PWGJE/AliPWG4HighPtSpectra.cxx b/PWGJE/AliPWG4HighPtSpectra.cxx
index 631be43cf31..5d440599280 100644
--- a/PWGJE/AliPWG4HighPtSpectra.cxx
+++ b/PWGJE/AliPWG4HighPtSpectra.cxx
@@ -54,6 +54,7 @@
 #include "AliMCEventHandler.h"
 #include "AliCFContainer.h"
 #include "AliGenPythiaEventHeader.h"
+#include "AliGenHijingEventHeader.h"
 #include "AliGenCocktailEventHeader.h"
 //#include "$ALICE_ROOT/PWG4/JetTasks/AliAnalysisHelperJetTasks.h"
 
@@ -71,11 +72,13 @@ AliPWG4HighPtSpectra::AliPWG4HighPtSpectra() : AliAnalysisTask("AliPWG4HighPtSpe
   fMC(0x0),
   fStack(0x0),
   fVtx(0x0),
+  fTriggerMask(AliVEvent::kMB),
   fIsPbPb(0),
   fCentClass(10),
   fTrackType(0),
   fTrackCuts(0x0),
   fTrackCutsReject(0x0),
+  fbSelectHIJING(kFALSE),
   fSigmaConstrainedMax(100.),
   fAvgTrials(1),
   fHistList(0),
@@ -104,11 +107,13 @@ AliPWG4HighPtSpectra::AliPWG4HighPtSpectra(const Char_t* name) :
   fMC(0x0),
   fStack(0x0),
   fVtx(0x0),
+  fTriggerMask(AliVEvent::kMB),
   fIsPbPb(0),
   fCentClass(10),
   fTrackType(0),
   fTrackCuts(0x0),
   fTrackCutsReject(0x0),
+  fbSelectHIJING(kFALSE),
   fSigmaConstrainedMax(100.),
   fAvgTrials(1),
   fHistList(0),
@@ -203,7 +208,7 @@ Bool_t AliPWG4HighPtSpectra::SelectEvent() {
 
   //Trigger
   UInt_t isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
-  if(!(isSelected&AliVEvent::kMB)) { //Select collison candidates
+  if(fTriggerMask != AliVEvent::kAny && !(isSelected&fTriggerMask)) { //Select collison candidates
     AliDebug(2,Form(" Trigger Selection: event REJECTED ... "));
     fNEventReject->Fill("Trigger",1);
     selectEvent = kFALSE;
@@ -307,7 +312,6 @@ void AliPWG4HighPtSpectra::Exec(Option_t *)
   fNEventAll->Fill(0.);
 
   if(!SelectEvent()) {
-    fNEventReject->Fill("NTracks<2",1);
     // Post output data
     PostData(0,fHistList);
     PostData(1,fCFManagerPos->GetParticleContainer());
@@ -347,117 +351,179 @@ void AliPWG4HighPtSpectra::Exec(Option_t *)
   //Now go to rec level
   for (Int_t iTrack = 0; iTrack<nTracks; iTrack++) 
     {   
+   
+
       //Get track for analysis
       AliESDtrack *track = 0x0;
       AliESDtrack *esdtrack = fESD->GetTrack(iTrack);
-      if(!esdtrack) continue;
+      if(!esdtrack)
+	continue;
+      
 
-      if(fTrackType==1) {
-	track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
-	if(!track) continue;
+      if(fTrackType==4) {
+	if (!(fTrackCuts->AcceptTrack(esdtrack)))
+	  continue;
       }
-      else if(fTrackType==2) {
-	track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
-	if(!track) continue;
 
+      if(fTrackType==1)
+	track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
+      else if(fTrackType==2 || fTrackType==4) {
+	track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(fESD),esdtrack->GetID());
+	if(!track)
+	  continue;
+	
 	AliExternalTrackParam exParam;
 	Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
 	if( !relate ) {
-	  delete track;
+	  if(track) delete track;
 	  continue;
 	}
 	track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
       }
+      else if(fTrackType==5 || fTrackType==6) {
+	if(fTrackCuts->AcceptTrack(esdtrack)) {
+	  continue;
+	}
+	else {
+	  if( !(fTrackCutsReject->AcceptTrack(esdtrack)) && fTrackCuts->AcceptTrack(esdtrack) ) {
+
+	    if(fTrackType==5) {
+	      //use TPConly constrained track
+	      track = AliESDtrackCuts::GetTPCOnlyTrack(fESD,esdtrack->GetID());
+	      if(!track)
+		continue;
+	      
+	      AliExternalTrackParam exParam;
+	      Bool_t relate = track->RelateToVertexTPC(fVtx,fESD->GetMagneticField(),kVeryBig,&exParam);
+	      if( !relate ) {
+		if(track) delete track;
+		continue;
+	      }
+	      track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
+	    }
+	    else if(fTrackType==6) {
+	      //use global constrained track
+	      track = new AliESDtrack(*esdtrack);
+	      track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
+
+	    }
+	  }
+	}
+      }
       else if(fTrackType==7) {
 	//use global constrained track
 	track = new AliESDtrack(*esdtrack);
-	//	track = esdtrack;
-	//	track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
-
       }
       else
 	track = esdtrack;
     
-      if(!track) continue;
- 
-      if(fTrackType==2) {
+      if(!track)
+	continue;
+      
+
+      if(fTrackType==2 || fTrackType==4 || fTrackType==5) {
 	//Cut on chi2 of constrained fit
-	if(track->GetConstrainedChi2TPC() > fSigmaConstrainedMax*fSigmaConstrainedMax) {
-	  delete track;
+	if(track->GetConstrainedChi2TPC() > fSigmaConstrainedMax*fSigmaConstrainedMax && fSigmaConstrainedMax>0.) {
+	  if(track) delete track;
 	  continue;
 	}
       }
+      if (!(fTrackCuts->AcceptTrack(track)) && fTrackType!=4 && fTrackType!=5 && fTrackType!=6) {
+	if(fTrackType==1 || fTrackType==2 || fTrackType==7) {
+	  if(track) delete track;
+	}
+	continue;
+      }
 
-      if (fTrackCuts->AcceptTrack(track)) {
-
-	if(fTrackType==7) {
-	  if(fTrackCutsReject ) {
-	    if(fTrackCutsReject->AcceptTrack(track) ) {
-	      if(track) delete track;
-	      continue;
-	    }
+      if(fTrackType==7) {
+	if(fTrackCutsReject ) {
+	  if(fTrackCutsReject->AcceptTrack(track) ) {
+	    if(track) delete track;
+	    continue;
 	  }
-	  
-	  if(esdtrack->GetConstrainedParam()) 
-	    track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
 	}
+      
+	if(esdtrack->GetConstrainedParam()) 
+	  track->Set(esdtrack->GetConstrainedParam()->GetX(),esdtrack->GetConstrainedParam()->GetAlpha(),esdtrack->GetConstrainedParam()->GetParameter(),esdtrack->GetConstrainedParam()->GetCovariance());
+      }
 
-	//fill the container
-	containerInputRec[0] = track->Pt();
-	containerInputRec[1] = track->Phi();
-	containerInputRec[2] = track->Eta();
-	containerInputRec[3] = track->GetTPCNclsIter1();
+      if(!track) {
+	if(fTrackType==1 || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
+	  if(track) delete track;
+	}
+	continue;
+      }
 
-	if(track->GetSign()>0.) fCFManagerPos->GetParticleContainer()->Fill(containerInputRec,kStepReconstructed);
-	if(track->GetSign()<0.) fCFManagerNeg->GetParticleContainer()->Fill(containerInputRec,kStepReconstructed);
+      //fill the container
+      containerInputRec[0] = track->Pt();
+      containerInputRec[1] = track->Phi();
+      containerInputRec[2] = track->Eta();
+      containerInputRec[3] = track->GetTPCNclsIter1();
+
+      if(track->GetSign()>0.) fCFManagerPos->GetParticleContainer()->Fill(containerInputRec,kStepReconstructed);
+      if(track->GetSign()<0.) fCFManagerNeg->GetParticleContainer()->Fill(containerInputRec,kStepReconstructed);
 
   	
-	//Only fill the MC containers if MC information is available
-	if(fMC) {
-	  Int_t label = TMath::Abs(track->GetLabel());
-	  TParticle *particle = fStack->Particle(label) ;
-	  if(!particle) {
+      //Only fill the MC containers if MC information is available
+      if(fMC) {
+	Int_t label = TMath::Abs(track->GetLabel());
+	if(label>fStack->GetNtrack()) {
+	  if(fTrackType==1 || fTrackType==2 || fTrackType==7)
+	    delete track;
+	  continue;
+	}
+	//Only select particles generated by HIJING if requested
+	if(fbSelectHIJING) {
+	  if(!IsHIJINGParticle(label)) {
 	    if(fTrackType==1 || fTrackType==2 || fTrackType==7)
 	      delete track;
 	    continue;
 	  }
-	  containerInputRecMC[0] = particle->Pt();      
-	  containerInputRecMC[1] = particle->Phi();      
-	  containerInputRecMC[2] = particle->Eta();  
-	  containerInputRecMC[3] = track->GetTPCNclsIter1();
-
-	  //Container with primaries
-	  if(fStack->IsPhysicalPrimary(label)) {
-	    if(particle->GetPDG()->Charge()>0.) {
-	      fCFManagerPos->GetParticleContainer()->Fill(containerInputRecMC,kStepReconstructedMC);
-	    }
-	    if(particle->GetPDG()->Charge()<0.) {
-	      fCFManagerNeg->GetParticleContainer()->Fill(containerInputRecMC,kStepReconstructedMC);
-	    }
+	}
+	TParticle *particle = fStack->Particle(label) ;
+	if(!particle) {
+	  if(fTrackType==1 || fTrackType==2 || fTrackType==7)
+	    delete track;
+	  continue;
+	}
+	containerInputRecMC[0] = particle->Pt();      
+	containerInputRecMC[1] = particle->Phi();      
+	containerInputRecMC[2] = particle->Eta();  
+	containerInputRecMC[3] = track->GetTPCNclsIter1();
+
+	//Container with primaries
+	if(fStack->IsPhysicalPrimary(label)) {
+	  if(particle->GetPDG()->Charge()>0.) {
+	    fCFManagerPos->GetParticleContainer()->Fill(containerInputRecMC,kStepReconstructedMC);
+	  }
+	  if(particle->GetPDG()->Charge()<0.) {
+	    fCFManagerNeg->GetParticleContainer()->Fill(containerInputRecMC,kStepReconstructedMC);
+	  }
 
-	    //Fill pT resolution plots for primaries
-	    fPtRelUncertainty1PtPrim->Fill(containerInputRec[0],containerInputRec[0]*TMath::Sqrt(track->GetSigma1Pt2()));
+	  //Fill pT resolution plots for primaries
+	  fPtRelUncertainty1PtPrim->Fill(containerInputRec[0],containerInputRec[0]*TMath::Sqrt(track->GetSigma1Pt2()));
 
-	  }
+	}
 
-	  //Container with secondaries
-	  if (!fStack->IsPhysicalPrimary(label) ) {
-	    if(particle->GetPDG()->Charge()>0.) {
-	      fCFManagerPos->GetParticleContainer()->Fill(containerInputMC,kStepSecondaries);
-	    }
-	    if(particle->GetPDG()->Charge()<0.) {
-	      fCFManagerNeg->GetParticleContainer()->Fill(containerInputMC,kStepSecondaries);
-	    }
-	    //Fill pT resolution plots for primaries
-	    fPtRelUncertainty1PtSec->Fill(containerInputRec[0],containerInputRec[0]*TMath::Sqrt(track->GetSigma1Pt2()));
+	//Container with secondaries
+	if (!fStack->IsPhysicalPrimary(label) ) {
+	  if(particle->GetPDG()->Charge()>0.) {
+	    fCFManagerPos->GetParticleContainer()->Fill(containerInputMC,kStepSecondaries);
+	  }
+	  if(particle->GetPDG()->Charge()<0.) {
+	    fCFManagerNeg->GetParticleContainer()->Fill(containerInputMC,kStepSecondaries);
 	  }
+	  //Fill pT resolution plots for primaries
+	  fPtRelUncertainty1PtSec->Fill(containerInputRec[0],containerInputRec[0]*TMath::Sqrt(track->GetSigma1Pt2()));
 	}
+      }
 	
-      }//trackCuts global tracks
 
-      if(fTrackType==1 || fTrackType==2 || fTrackType==7) {
+      if(fTrackType==1  || fTrackType==2 || fTrackType==4 || fTrackType==5 || fTrackType==6 || fTrackType==7) {
 	if(track) delete track;
       }
+
+
     }//track loop
   
 
@@ -467,6 +533,12 @@ void AliPWG4HighPtSpectra::Exec(Option_t *)
       AliMCParticle *mcPart  = (AliMCParticle*)fMC->GetTrack(iPart);
       if(!mcPart) continue;
 
+      //Only select particles generated by HIJING if requested
+      if(fbSelectHIJING) {
+	if(!IsHIJINGParticle(iPart))
+	  continue;
+      }
+
       Int_t pdg = mcPart->PdgCode();
       
       // select charged pions, protons, kaons , electrons, muons
@@ -521,7 +593,7 @@ Bool_t AliPWG4HighPtSpectra::PythiaInfoFromFile(const char* currFile,Float_t &fX
     // next trial fetch the histgram file
     fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
     if(!fxsec){
-	// not a severe condition but inciate that we have no information
+      // not a severe condition but indicates that we have no information
       return kFALSE;
     }
     else{
@@ -618,6 +690,36 @@ AliGenPythiaEventHeader*  AliPWG4HighPtSpectra::GetPythiaEventHeader(AliMCEvent
 
 }
 
+//________________________________________________________________________
+AliGenHijingEventHeader*  AliPWG4HighPtSpectra::GetHijingEventHeader(AliMCEvent *mcEvent){
+  
+  if(!mcEvent)return 0;
+  AliGenEventHeader* genHeader = mcEvent->GenEventHeader();
+  AliGenHijingEventHeader* hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
+  if(!hijingGenHeader){
+    // cocktail ??
+    AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
+    
+    if (!genCocktailHeader) {
+      //      AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Unknown header type (not Pythia or Cocktail)");
+      //      AliWarning(Form("%s %d: Unknown header type (not Pythia or Cocktail)",(char*)__FILE__,__LINE__));
+      return 0;
+    }
+    TList* headerList = genCocktailHeader->GetHeaders();
+    for (Int_t i=0; i<headerList->GetEntries(); i++) {
+      hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(headerList->At(i));
+      if (hijingGenHeader)
+        break;
+    }
+    if(!hijingGenHeader){
+      AliWarningGeneral(Form(" %s:%d",(char*)__FILE__,__LINE__),"Hijing event header not found");
+      return 0;
+    }
+  }
+  return hijingGenHeader;
+
+}
+
 
 //___________________________________________________________________________
 void AliPWG4HighPtSpectra::Terminate(Option_t*)
@@ -706,4 +808,31 @@ void AliPWG4HighPtSpectra::CreateOutputObjects() {
 
 }
 
+//________________________________________________________________________
+Bool_t AliPWG4HighPtSpectra::IsHIJINGParticle(Int_t label) {
+  //
+  // Return kTRUE in case particle is from HIJING event
+  //
+
+  AliGenHijingEventHeader*  hijingHeader = GetHijingEventHeader(fMC);
+
+  Int_t nproduced = hijingHeader->NProduced();
+
+  TParticle * mom = fStack->Particle(label);
+  Int_t iMom = label;
+  Int_t iParent = mom->GetFirstMother();
+  while(iParent!=-1){
+    iMom = iParent;
+    mom = fStack->Particle(iMom);
+    iParent = mom->GetFirstMother();
+  }
+
+  if(iMom<nproduced) {
+    return kTRUE;
+  } else {
+    return kFALSE;
+  }
+
+}
+
 #endif