Ernesto's changes for train running

diff --git a/PWG4/JetTasks/AliAnalysisTaskUE.cxx b/PWG4/JetTasks/AliAnalysisTaskUE.cxx
index 1e668feb366484d347cab10151c0c8235fbf3d4a..addbd842cdad6bc2be81be1ad97849f7b600e23d 100644 (file)
--- a/PWG4/JetTasks/AliAnalysisTaskUE.cxx
+++ b/PWG4/JetTasks/AliAnalysisTaskUE.cxx
@@ -101,6 +101,7 @@ fConeRadius(0.7),
 fConePosition(1),
 fAreaReg(1.5393), // Pi*0.7*0.7
 fUseChPartJet(kFALSE),
+fUseChargeHadrons(kFALSE),
 fUseSingleCharge(kFALSE),
 fUsePositiveCharge(kTRUE),
 fOrdering(1),
@@ -170,6 +171,7 @@ Bool_t AliAnalysisTaskUE::Notify()
     }
     AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials); 
     fh1Xsec->Fill("<#sigma>",xsection);
+
    // construct average trials
    Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
    if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
@@ -255,6 +257,7 @@ void  AliAnalysisTaskUE::Exec(Option_t */*option*/)
        if (fDebug > 1) AliInfo(" Trigger Selection: event REJECTED ... ");
        return;
   }
+  
   //Event selection (vertex) *****************************************
   AliKFVertex primVtx(*(fAOD->GetPrimaryVertex()));
   Int_t nTracksPrim=primVtx.GetNContributors();
@@ -265,7 +268,6 @@ void  AliAnalysisTaskUE::Exec(Option_t */*option*/)
        }
   if (fDebug > 1) AliInfo(" Primary-vertex Selection: event ACCEPTED ...");
  
-  
   // Execute analysis for current event
   //
   if ( fDebug > 3 ) AliInfo( " Processing event..." );
@@ -336,6 +338,7 @@ void  AliAnalysisTaskUE::AnalyseUE()
       jet = fAODjets->GetJet(i);
       }
       Double_t jetPt = jet->Pt();//*1.666; // FIXME Jet Pt Correction ?????!!!
+ 
       if( jetPt > maxPtJet1 ) {
        maxPtJet3 = maxPtJet2; index3 = index2;
        maxPtJet2 = maxPtJet1; index2 = index1;
@@ -347,6 +350,7 @@ void  AliAnalysisTaskUE::AnalyseUE()
        maxPtJet3 = jetPt; index3 = i;
       }
     }
+
     if( index1 != -1 ) {
       AliAODJet *jet = 0;
       if (fDeltaAOD) {
@@ -390,13 +394,13 @@ void  AliAnalysisTaskUE::AnalyseUE()
       if( nJets > 1 ) {
        index2 = 1;
        AliAODJet* jet = (AliAODJet*)jets->At(1);
-       maxPtJet1 = jet->Pt();
+        maxPtJet2 = jet->Pt();
        jetVect[1].SetXYZ(jet->Px(), jet->Py(), jet->Pz());
       }
       if( nJets > 2 ) {
        index3 = 2;
        AliAODJet* jet = (AliAODJet*)jets->At(2);
-       maxPtJet1 = jet->Pt();
+        maxPtJet3 = jet->Pt();
        jetVect[2].SetXYZ(jet->Px(), jet->Py(), jet->Pz());
       }
       
@@ -465,7 +469,10 @@ fhNJets->Fill(nJets);
   Double_t maxPartPtRegion  = 0.;
   Int_t    nTrackRegionPosit = 0;
   Int_t    nTrackRegionNegat = 0;
-  static const Double_t k270rad = 270.*TMath::Pi()/180.;
+  static Double_t const  kPI     = TMath::Pi();
+  static Double_t const  kTWOPI  = 2.*kPI;
+  static Double_t const  k270rad = 270.*kPI/180.;
+
   
   if (!fUseMCParticleBranch){
     fhEleadingPt->Fill( maxPtJet1 );
@@ -480,7 +487,7 @@ fhNJets->Fill(nJets);
       Bool_t isHadron = part->GetMostProbablePID()==AliAODTrack::kPion || 
                         part->GetMostProbablePID()==AliAODTrack::kKaon || 
                         part->GetMostProbablePID()==AliAODTrack::kProton;
-      if ( !isHadron ) continue;
+      if ( fUseChargeHadrons && !isHadron ) continue;
       
       if ( !part->Charge() ) continue; //Only charged
       if ( fUseSingleCharge ) { // Charge selection
@@ -494,8 +501,8 @@ fhNJets->Fill(nJets);
       TVector3 partVect(part->Px(), part->Py(), part->Pz());
       
       Double_t deltaPhi = jetVect[0].DeltaPhi(partVect)+k270rad;
-      if( deltaPhi > 2.*TMath::Pi() )  deltaPhi-= 2.*TMath::Pi();
-      fhdNdEtaPhiDist->Fill( deltaPhi );
+      if( deltaPhi > kTWOPI )  deltaPhi-= kTWOPI;
+      fhdNdEtaPhiDist->Fill( deltaPhi, maxPtJet1 );
       
       fhFullRegPartPtDistVsEt->Fill( part->Pt(), maxPtJet1 );
       
@@ -619,13 +626,13 @@ fhNJets->Fill(nJets);
                           TMath::Abs(mctrk->GetPdgCode())==2212 ||
                           TMath::Abs(mctrk->GetPdgCode())==321;
         
-        if (!isHadron) continue;
+        if ( fUseChargeHadrons && !isHadron ) continue;
         
         TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
 
         Double_t deltaPhi = jetVectnew[index].DeltaPhi(partVect)+k270rad;
         if( deltaPhi > 2.*TMath::Pi() )  deltaPhi-= 2.*TMath::Pi();
-        fhdNdEtaPhiDist->Fill( deltaPhi );
+        fhdNdEtaPhiDist->Fill( deltaPhi, maxPtJet1 );
         
         fhFullRegPartPtDistVsEt->Fill( mctrk->Pt(), maxPtJet1 );
         
@@ -671,13 +678,13 @@ fhNJets->Fill(nJets);
                           TMath::Abs(mctrk->GetPdgCode())==2212 ||
                           TMath::Abs(mctrk->GetPdgCode())==321;
         
-        if (!isHadron) continue;
+        if ( fUseChargeHadrons && !isHadron ) continue;
         
         TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
 
         Double_t deltaPhi = jetVectnew[index].DeltaPhi(partVect)+k270rad;
         if( deltaPhi > 2.*TMath::Pi() )  deltaPhi-= 2.*TMath::Pi();
-        fhdNdEtaPhiDist->Fill( deltaPhi );
+        fhdNdEtaPhiDist->Fill( deltaPhi, maxPtJet1 );
 
         fhFullRegPartPtDistVsEt->Fill( mctrk->Pt(), maxPtJet1 );
         
@@ -848,14 +855,14 @@ TObjArray*  AliAnalysisTaskUE::FindChargedParticleJets()
 {
   // Return a TObjArray of "charged particle jets"
   //
-  // Charged particle jet deï¬\81nition from reference:
+  // Charged particle jet deï¬Â\81nition from reference:
   //
   // "Charged jet evolution and the underlying event
   //  in proton-antiproton collisions at 1.8 TeV"
   //  PHYSICAL REVIEW D 65 092002, CDF Collaboration
   //
-  // We deï¬\81ne "jets" as circular regions in eta-phi space with
-  // radius deï¬\81ned by R = sqrt( (eta-eta0)^2 +(phi-phi0)^2 ).
+  // We deï¬Â\81ne "jets" as circular regions in eta-phi space with
+  // radius deï¬Â\81ned by R = sqrt( (eta-eta0)^2 +(phi-phi0)^2 ).
   // Our jet algorithm is as follows:
   //   1- Order all charged particles according to their pT .
   //   2- Start with the highest pT particle and include in the jet all
@@ -866,7 +873,7 @@ TObjArray*  AliAnalysisTaskUE::FindChargedParticleJets()
   //      a jet and add to the jet all particles not already included in
   //      a jet within R=0.7.
   //   4- Continue until all particles are in a jet.
-  // We deï¬\81ne the transverse momentum of the jet to be
+  // We deï¬Â\81ne the transverse momentum of the jet to be
   // the scalar pT sum of all the particles within the jet, where pT
   // is measured with respect to the beam axis
   
@@ -886,6 +893,7 @@ TObjArray*  AliAnalysisTaskUE::FindChargedParticleJets()
   
   nTracks = tracks.GetEntriesFast();
   if( !nTracks ) return 0;
+
   TObjArray *jets = new TObjArray(nTracks);
   TIter itrack(&tracks);
   while( nTracks ) {
@@ -900,17 +908,21 @@ TObjArray*  AliAnalysisTaskUE::FindChargedParticleJets()
     jets->AddLast( new TLorentzVector(px, py, pz, pt) );
     tracks.Remove( track );
     TLorentzVector* jet = (TLorentzVector*)jets->Last();
+    jet->SetPtEtaPhiE( 1., jet->Eta(), jet->Phi(), pt );
     // 3- Go to the next highest pT particle not already included...
     AliAODTrack* track1;
+    Double_t fPt = jet->E();
     while ( (track1  = (AliAODTrack*)(itrack.Next())) ) {
-      Double_t dphi = TVector2::Phi_mpi_pi(jet->Phi()-track1->Phi());
+      Double_t tphi = track1->Phi(); // here Phi is from 0 <-> 2Pi
+      if (tphi > TMath::Pi()) tphi -= 2. * TMath::Pi(); // convert to  -Pi <-> Pi
+      Double_t dphi = TVector2::Phi_mpi_pi(jet->Phi()-tphi);
       Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
                                dphi*dphi );
       if( r < fConeRadius ) {
-        Double_t fPt   = jet->E()+track1->Pt();  // Scalar sum of Pt
+        fPt   = jet->E()+track1->Pt();  // Scalar sum of Pt
         // recalculating the centroid
         Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
-        Double_t phi = jet->Phi()*jet->E()/fPt + track1->Phi()*track1->Pt()/fPt;
+        Double_t phi = jet->Phi()*jet->E()/fPt + tphi*track1->Pt()/fPt;
         jet->SetPtEtaPhiE( 1., eta, phi, fPt );
         tracks.Remove( track1 );
       }
@@ -934,6 +946,7 @@ TObjArray*  AliAnalysisTaskUE::FindChargedParticleJets()
     py = jet->E() * TMath::Sin(jet->Phi());  // Pt * sin(phi)
     pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
     en = TMath::Sqrt(px * px + py * py + pz * pz);
+
     aodjets->AddLast( new AliAODJet(px, py, pz, en) );
   }
   // Order jets according to their pT .
@@ -1010,7 +1023,7 @@ void  AliAnalysisTaskUE::CreateHistos()
   fListOfHistos = new TList();
   
   
-  fhNJets = new TH1F("fhNJets", "n Jet",  10, 0, 10);
+  fhNJets = new TH1F("fhNJets", "n Jet",  20, 0, 20);
   fhNJets->SetXTitle("# of jets");
   fhNJets->Sumw2();
   fListOfHistos->Add( fhNJets );                 // At(0) 
@@ -1055,23 +1068,24 @@ void  AliAnalysisTaskUE::CreateHistos()
   fhMinRegSumPtvsMult->Sumw2();
   fListOfHistos->Add( fhMinRegSumPtvsMult );     // At(7);
   
-  fhdNdEtaPhiDist  = new TH1F("hdNdEtaPhiDist",   "Charge particle density |#eta|< 1 vs #Delta#phi",  120, 0.,   2.*TMath::Pi());
+  fhdNdEtaPhiDist  = new TH2F("hdNdEtaPhiDist", Form("Charge particle density |#eta|<%3.1f vs #Delta#phi", fTrackEtaCut),  
+                               62, 0.,   2.*TMath::Pi(), fBinsPtInHist, fMinJetPtInHist, fMaxJetPtInHist);
   fhdNdEtaPhiDist->SetXTitle("#Delta#phi");
-  fhdNdEtaPhiDist->SetYTitle("dN_{ch}/d#etad#phi");
+  fhdNdEtaPhiDist->SetYTitle("Leading Jet P_{T}");
   fhdNdEtaPhiDist->Sumw2();
   fListOfHistos->Add( fhdNdEtaPhiDist );        // At(8)
   
   // Can be use to get part pt distribution for differente Jet Pt bins
-  fhFullRegPartPtDistVsEt = new TH2F("hFullRegPartPtDistVsEt", "dN/dP_{T} |#eta|<1 vs Leading Jet P_{T}",  
-                                     50,0.,50., fBinsPtInHist, fMinJetPtInHist,   fMaxJetPtInHist);
+  fhFullRegPartPtDistVsEt = new TH2F("hFullRegPartPtDistVsEt", Form( "dN/dP_{T} |#eta|<%3.1f vs Leading Jet P_{T}", fTrackEtaCut),
+                                     100,0.,50., fBinsPtInHist, fMinJetPtInHist, fMaxJetPtInHist);
   fhFullRegPartPtDistVsEt->SetYTitle("Leading Jet P_{T}");
   fhFullRegPartPtDistVsEt->SetXTitle("p_{T}");
   fhFullRegPartPtDistVsEt->Sumw2();
   fListOfHistos->Add( fhFullRegPartPtDistVsEt );  // At(9) 
   
   // Can be use to get part pt distribution for differente Jet Pt bins
-  fhTransRegPartPtDistVsEt = new TH2F("hTransRegPartPtDistVsEt", "dN/dP_{T} in tranvese regions |#eta|<1 vs Leading Jet P_{T}",  
-                                      50,0.,50., fBinsPtInHist, fMinJetPtInHist,   fMaxJetPtInHist);
+  fhTransRegPartPtDistVsEt = new TH2F("hTransRegPartPtDistVsEt", Form( "dN/dP_{T} in tranvese regions |#eta|<%3.1f vs Leading Jet P_{T}", fTrackEtaCut),  
+                                     100,0.,50., fBinsPtInHist, fMinJetPtInHist,   fMaxJetPtInHist);
   fhTransRegPartPtDistVsEt->SetYTitle("Leading Jet P_{T}");
   fhTransRegPartPtDistVsEt->SetXTitle("p_{T}");
   fhTransRegPartPtDistVsEt->Sumw2();
@@ -1152,6 +1166,7 @@ void  AliAnalysisTaskUE::CreateHistos()
   fSettingsTree->Branch("fRegionType", &fRegionType,"Reg/I");
   fSettingsTree->Branch("fOrdering", &fOrdering,"OrderMeth/I");
   fSettingsTree->Branch("fUseChPartJet", &fUseChPartJet,"UseChPart/O");
+  fSettingsTree->Branch("fUseChargeHadrons", &fUseChargeHadrons,"UseChHadrons/O");
   fSettingsTree->Branch("fUseSingleCharge", &fUseSingleCharge,"UseSingleCh/O");
   fSettingsTree->Branch("fUsePositiveCharge", &fUsePositiveCharge,"UsePositiveCh/O");
   fSettingsTree->Fill();
@@ -1191,16 +1206,15 @@ void  AliAnalysisTaskUE::Terminate(Option_t */*option*/)
       AliError("Histogram List is not available");
       return;
     }
+    fhNJets              = (TH1F*)fListOfHistos->At(0);
     fhEleadingPt         = (TH1F*)fListOfHistos->At(1);
-    fhdNdEtaPhiDist     = (TH1F*)fListOfHistos->At(8);
+    fhdNdEtaPhiDist      = (TH2F*)fListOfHistos->At(8);
     fhRegionSumPtMaxVsEt = (TH1F*)fListOfHistos->At(11);
     fhRegionSumPtMinVsEt = (TH1F*)fListOfHistos->At(12);
     fhRegionMultMaxVsEt  = (TH1F*)fListOfHistos->At(14);
     fhRegionMultMinVsEt  = (TH1F*)fListOfHistos->At(15);
     fhRegionAveSumPtVsEt = (TH1F*)fListOfHistos->At(16);
     
-    fhNJets              = (TH1F*)fListOfHistos->At(0);
-    
     //fhValidRegion  = (TH2F*)fListOfHistos->At(21);
     
     // Canvas
@@ -1269,15 +1283,18 @@ void  AliAnalysisTaskUE::Terminate(Option_t */*option*/)
     c2->cd(1);
     fhEleadingPt->SetMarkerStyle(20);
     fhEleadingPt->SetMarkerColor(2);
-    fhEleadingPt->Scale(normFactor);
+    if( normFactor > 0.) fhEleadingPt->Scale(normFactor);
     //fhEleadingPt->Draw("p");
     fhEleadingPt->DrawCopy("p");
     gPad->SetLogy();
     
     c2->cd(2);
-    fhdNdEtaPhiDist->SetMarkerStyle(20);
-    fhdNdEtaPhiDist->SetMarkerColor(2);
-    fhdNdEtaPhiDist->DrawCopy("p");
+    Int_t xbin1 = fhdNdEtaPhiDist->GetYaxis()->FindFixBin(fMinJetPtInHist);
+    Int_t xbin2 = fhdNdEtaPhiDist->GetYaxis()->FindFixBin(fMaxJetPtInHist);
+    TH1D* dNdEtaPhiDistAllJets = fhdNdEtaPhiDist->ProjectionX("dNdEtaPhiDistAllJets",xbin1,xbin2);
+    dNdEtaPhiDistAllJets->SetMarkerStyle(20);
+    dNdEtaPhiDistAllJets->SetMarkerColor(2);
+    dNdEtaPhiDistAllJets->DrawCopy("p");
     gPad->SetLogy();
     
     c2->cd(3);      

diff --git a/PWG4/JetTasks/AliAnalysisTaskUE.h b/PWG4/JetTasks/AliAnalysisTaskUE.h
index a0adc6d83898c9d4a9340ab2a23d650eea5e2a0e..bbd2b943f2b0d24a2eebccf14c76a715a82b2887 100644 (file)
--- a/PWG4/JetTasks/AliAnalysisTaskUE.h
+++ b/PWG4/JetTasks/AliAnalysisTaskUE.h
@@ -52,6 +52,7 @@ class  AliAnalysisTaskUE : public AliAnalysisTask
     void   SetAnaTopology( Int_t val )    { fAnaType = val;    }
     void   SetRegionType( Int_t val )     { fRegionType = val; }
     void   SetUseChPartJet( Int_t val )   { fUseChPartJet = val; }
+    void   SetUseChargeHadrons( Bool_t val ){ fUseChargeHadrons = val; }
     void   SetPtSumOrdering( Bool_t val ) { fOrdering = val;   }
     void   SetFilterBit( UInt_t val )     { fFilterBit = val;  }
     void   SetJetsOnFly( Bool_t val )     { fJetsOnFly = val;  }
@@ -134,13 +135,13 @@ class  AliAnalysisTaskUE : public AliAnalysisTask
     // if fRegionType = 2 not always it is included within eta range
     Bool_t   fUseChPartJet;     // Use "Charged Particle Jet" instead of jets from AOD
     // see FindChargedParticleJets()
+    Bool_t     fUseChargeHadrons;   // Only use charge hadrons
     
     // Theoreticians ask for tools charge-aware
     // especially those which are related to multiplicity and MC-tunings
     // see arXiv:hep-ph/0507008v3
     Bool_t   fUseSingleCharge;     //Make analysis for a single type of charge (=kFALSE default)
     Bool_t   fUsePositiveCharge;   //If Single type of charge used then set which one (=kTRUE default positive)
-    
     Int_t   fOrdering;         //  Pt and multiplicity summation ordering:
     //     1=CDF-like -independent sorting according quantity to be scored: Double sorting- (default)
     //       if Pt summation will be scored take Pt minimum between both zones and 
@@ -179,7 +180,7 @@ class  AliAnalysisTaskUE : public AliAnalysisTask
     TH1F*  fhMinRegMaxPtPart;        //!
     TH1F*  fhMinRegSumPtvsMult;      //!
     
-    TH1F*  fhdNdEtaPhiDist;         //!
+    TH2F*  fhdNdEtaPhiDist;          //!
     TH2F*  fhFullRegPartPtDistVsEt;  //!
     TH2F*  fhTransRegPartPtDistVsEt; //!
     
@@ -204,7 +205,7 @@ class  AliAnalysisTaskUE : public AliAnalysisTask
     
     
     
-    ClassDef( AliAnalysisTaskUE, 3); // Analysis task for Underlying Event analysis
+    ClassDef( AliAnalysisTaskUE, 4); // Analysis task for Underlying Event analysis
   };
 
 #endif