]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Adding TPRofiles for mean and sigma(pT), clean up of EP from tracks, added EP from...
authorkleinb <kleinb@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 17 Aug 2011 14:32:40 +0000 (14:32 +0000)
committerkleinb <kleinb@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 17 Aug 2011 14:32:40 +0000 (14:32 +0000)
PWG4/JetTasks/AliAnalysisTaskJetServices.cxx
PWG4/JetTasks/AliAnalysisTaskJetServices.h
PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.cxx
PWG4/JetTasks/AliAnalysisTaskJetSpectrum2.h
PWG4/macros/AddTaskJetSpectrum2.C
PWG4/macros/AnalysisTrainPWG4Jets.C
PWG4/macros/ConfigTrainGrid.C

index a5f2bb39974e6c2642a6ea866753ca99de59ba39..feb5f47bcbe0b9046a59f07200f5c2ee5850b157 100644 (file)
@@ -53,6 +53,7 @@
 #include "AliAODEvent.h"
 #include "AliAODHandler.h"
 #include "AliAODTrack.h"
+#include "AliAODVZERO.h"
 #include "AliAODJet.h"
 #include "AliAODMCParticle.h"
 #include "AliMCEventHandler.h"
@@ -69,6 +70,7 @@
 
 ClassImp(AliAnalysisTaskJetServices)
 
+AliAODVZERO*  AliAnalysisTaskJetServices::fgAODVZERO = NULL;
 AliAODHeader*  AliAnalysisTaskJetServices::fgAODHeader = NULL;
 TClonesArray*   AliAnalysisTaskJetServices::fgAODVertices = NULL;
 
@@ -82,7 +84,7 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices():
   fSelectionInfoESD(0),
   fEventCutInfoESD(0),
   fFilterMask(0),
-  fRPSubeventMethod(0),
+  fRPMethod(0),
   fCollisionType(kPbPb),
   fAvgTrials(1),
   fVtxXMean(0),
@@ -97,6 +99,8 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices():
   fTrackRecEtaWindow(0.9),
   fMinTrackPt(0.15),
   fRPAngle(0),
+  fPsiVZEROA(0),
+  fPsiVZEROC(0),
   fRandomizer(0),
   fNonStdFile(""),
   fh1Xsec(0x0),
@@ -115,21 +119,22 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices():
   fh2ESDTriggerRun(0x0),
   fh2VtxXY(0x0),
   fh1NCosmicsPerEvent(0x0),
-  fh2RPSubevents(0x0),
+  fp1RPXA(0x0),
+  fp1RPYA(0x0),
+  fp1RPXC(0x0),
+  fp1RPYC(0x0),
+  fh2RPAC(0x0),
+  fh2RPAT(0x0),
+  fh2RPCT(0x0),
+  fh2XYA(0x0),
+  fh2XYC(0x0),
   fh2RPCentrality(0x0),
-  fh2RPDeltaRP(0x0),
-  fh2RPQxQy(0x0),
-  fh2RPCosDeltaRP(0x0),
-  fh3PhiWeights(0x0),
-  fh3RPPhiTracks(0x0),
+  fh2RPACentrality(0x0),
+  fh2RPCCentrality(0x0),
   fTriggerAnalysis(0x0),
   fHistList(0x0)  
 {
   fRunRange[0] = fRunRange[1] = 0; 
-  fFlatA[0] =   fFlatA[1] = 0;
-  fFlatB[0] =   fFlatB[1] = 0;
-  fDeltaQxy[0] =   fDeltaQxy[1] = 0; 
-
 }
 
 AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
@@ -142,7 +147,7 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
   fSelectionInfoESD(0),
   fEventCutInfoESD(0),
   fFilterMask(0),
-  fRPSubeventMethod(0),
+  fRPMethod(0),
   fCollisionType(kPbPb),
   fAvgTrials(1),
   fVtxXMean(0),
@@ -157,6 +162,8 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
   fTrackRecEtaWindow(0.9),
   fMinTrackPt(0.15),
   fRPAngle(0),
+  fPsiVZEROA(0),
+  fPsiVZEROC(0),
   fRandomizer(0),
   fNonStdFile(""),
   fh1Xsec(0x0),
@@ -175,21 +182,22 @@ AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
   fh2ESDTriggerRun(0x0),
   fh2VtxXY(0x0),
   fh1NCosmicsPerEvent(0x0),
-  fh2RPSubevents(0x0),
+  fp1RPXA(0x0),
+  fp1RPYA(0x0),
+  fp1RPXC(0x0),
+  fp1RPYC(0x0),
+  fh2RPAC(0x0),
+  fh2RPAT(0x0),
+  fh2RPCT(0x0),
+  fh2XYA(0x0),
+  fh2XYC(0x0),
   fh2RPCentrality(0x0),
-  fh2RPDeltaRP(0x0),
-  fh2RPQxQy(0x0),
-  fh2RPCosDeltaRP(0x0),
-  fh3PhiWeights(0x0),
-  fh3RPPhiTracks(0x0),
-
+  fh2RPACentrality(0x0),
+  fh2RPCCentrality(0x0),
   fTriggerAnalysis(0x0),
   fHistList(0x0)  
 {
   fRunRange[0] = fRunRange[1] = 0; 
-  fFlatA[0] =   fFlatA[1] = 0;
-  fFlatB[0] =   fFlatB[1] = 0;
-  fDeltaQxy[0] =   fDeltaQxy[1] = 0; 
   DefineOutput(1,TList::Class());
 }
 
@@ -317,28 +325,46 @@ void AliAnalysisTaskJetServices::UserCreateOutputObjects()
   }
 
   fh1NCosmicsPerEvent = new TH1F("fh1NCosmicsPerEvent","Number of cosmic candidates per event",10,0.,10.);
+  fHistList->Add(fh1NCosmicsPerEvent),
 
-  
-  fh2RPSubevents = new TH2F("fh2RPSubevents" ,"Reaction Plane Angle" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
-  fHistList->Add( fh2RPSubevents);
 
-  fh2RPCentrality = new TH2F("fh2RPCentrality" ,"Reaction Plane Angle" , 20, 0.,100., 180, 0, TMath::Pi());
+  fh2RPCentrality = new TH2F("fh2RPCentrality" ,"Reaction Plane Angle from tracks" , 20, 0.,100., 180, 0, TMath::Pi());
   fHistList->Add(fh2RPCentrality);
 
-  fh2RPDeltaRP   = new TH2F("fh2DeltaRP" ,"Delta Reaction Plane Angle" , 100, -TMath::Pi()/2, TMath::Pi()/2,20,0.,100.0);
-  fHistList->Add(fh2RPDeltaRP);
 
-  fh2RPQxQy      = new TH2F("fh2RPQxQy" ,"" , 100, -100,100,100,-100,100);
-  fHistList->Add(fh2RPQxQy);
+  fh2RPACentrality = new TH2F("fh2RPACentrality" ,"Reaction Plane Angle from vzero A" , 20, 0.,100., 180, 0, TMath::Pi());
+  fHistList->Add(fh2RPACentrality);
 
-  fh2RPCosDeltaRP = new TH2F("fh2RPCosDeltaRP" ,"" , 20, 0.001,100.001,100,-1,1);
-  fHistList->Add(fh2RPCosDeltaRP);
+  fh2RPCCentrality = new TH2F("fh2RPCCentrality" ,"Reaction Plane Angle from vzero C" , 20, 0.,100., 180, 0, TMath::Pi());
+  fHistList->Add(fh2RPCCentrality);
 
-  fh3RPPhiTracks = new TH3F("fh3RPPhiTracks","Phi Tracks Pt Centrality", 10, 0.,100.,20,-5,5,180, 0, 2*TMath::Pi());
-  fHistList->Add(fh3RPPhiTracks);
-  
+  fh2RPAC = new TH2F("fh2RPAC" ,"Reaction Plane Angle vzero a vs c" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
+  fHistList->Add( fh2RPAC);
 
-  fHistList->Add(fh1NCosmicsPerEvent),
+  fh2RPAT = new TH2F("fh2RPAT" ,"Reaction Plane Angle vzero a vs tracks" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
+  fHistList->Add( fh2RPAT);
+
+  fh2RPCT = new TH2F("fh2RPCT" ,"Reaction Plane Angle vzero c vs tracks" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
+  fHistList->Add( fh2RPCT);
+
+  fh2XYA = new TH2F("fh2XYA" ,"XY vzeroa ;X;Y;" ,100,-0.3,0.3,100,-0.3,0.3);
+  fHistList->Add(fh2XYA);
+  fh2XYC = new TH2F("fh2XYC" ,"XY vzeroc ;X;Y;" ,100,-0.3,0.3,100,-0.3,0.3);
+  fHistList->Add(fh2XYC);
+
+  // profiles for mean 
+  fp1RPXA = new TProfile("fp1RPXA","mean vzeroa x vs run number;run;x",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
+  fHistList->Add(fp1RPXA);
+
+  fp1RPYA = new TProfile("fp1RPYA","mean vzeroa y vs run number;run;y",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
+  fHistList->Add(fp1RPYA);
+
+
+  fp1RPXC = new TProfile("fp1RPXC","mean vzeroc x vs run number;run;x",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
+  fHistList->Add(fp1RPXC);
+
+  fp1RPYC = new TProfile("fp1RPYC","mean vzeroa y vs run number;run;y",(Int_t)(1+fRunRange[1]-fRunRange[0]),fRunRange[0]-0.5,fRunRange[1]+0.5);
+  fHistList->Add(fp1RPYC);
 
 
   TH1::AddDirectory(oldStatus);
@@ -348,6 +374,9 @@ void AliAnalysisTaskJetServices::UserCreateOutputObjects()
      if (fDebug > 1) AliInfo("Replicating header");
      fgAODHeader = new AliAODHeader;
      AddAODBranch("AliAODHeader",&fgAODHeader,fNonStdFile.Data());
+     if (fDebug > 1) AliInfo("Replicating vzeros");
+     fgAODVZERO = new AliAODVZERO;
+     AddAODBranch("AliAODVZERO",&fgAODVZERO,fNonStdFile.Data());
      if (fDebug > 1) AliInfo("Replicating primary vertices");
      fgAODVertices = new TClonesArray("AliAODVertex",3);
      fgAODVertices->SetName("vertices");
@@ -485,6 +514,9 @@ void AliAnalysisTaskJetServices::UserExec(Option_t */*option*/)
   }
   fCentrality = cent;
   fRPAngle = 0;
+  fPsiVZEROA = 0;
+  fPsiVZEROC = 0;
+
 
   if(esd){
     const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
@@ -598,12 +630,19 @@ void AliAnalysisTaskJetServices::UserExec(Option_t */*option*/)
 
       TList recTracks;
       GetListOfTracks(&recTracks);
-      //      CalculateReactionPlaneAngle(&recTracks);
+      CalculateReactionPlaneAngleVZERO(aod);
       fRPAngle = aod->GetHeader()->GetEventplane();
       fh1RP->Fill(fRPAngle);
       fh2RPCentrality->Fill(fCentrality,fRPAngle);
+      fh2RPACentrality->Fill(fCentrality,fPsiVZEROA);
+      fh2RPCCentrality->Fill(fCentrality,fPsiVZEROC);
+      fh2RPAC->Fill(fPsiVZEROA,fPsiVZEROC);
+      fh2RPAT->Fill(fPsiVZEROA,fRPAngle);
+      fh2RPCT->Fill(fPsiVZEROC,fRPAngle);
+      if(fRPMethod==kRPTracks)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fRPAngle); // set slection to false
+      else if(fRPMethod==kRPVZEROA)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fPsiVZEROA); // set slection to false
+      else if(fRPMethod==kRPVZEROC)AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fPsiVZEROA); // set slection to false
 
-      AliAnalysisHelperJetTasks::ReactionPlane(kTRUE,fRPAngle); // set slection to false
       if(fUseAODInput&&fCentrality<=80){
        if(fFilterAODCollisions&&aod){
          aodH->SetFillAOD(kTRUE);
@@ -647,8 +686,11 @@ void AliAnalysisTaskJetServices::UserExec(Option_t */*option*/)
   if(fNonStdFile.Length()&&aod){
     if (fgAODHeader){
       *fgAODHeader =  *(dynamic_cast<AliAODHeader*>(aod->GetHeader()));
-      Double_t q[2] = {fRPAngle,fRPAngle};
-      fgAODHeader->SetQTheta(q,2);
+      Double_t q[kRPMethods] = {fRPAngle,fPsiVZEROA,fPsiVZEROC};
+      fgAODHeader->SetQTheta(q,kRPMethods);
+    }
+    if (fgAODVZERO){
+      *fgAODVZERO =  *(dynamic_cast<AliAODVZERO*>(aod->GetVZEROData()));
     }
     if(fgAODVertices){
       fgAODVertices->Delete();
@@ -945,134 +987,80 @@ void AliAnalysisTaskJetServices::Terminate(Option_t */*option*/)
   //
 }
 
-Bool_t AliAnalysisTaskJetServices::CalculateReactionPlaneAngle(const TList *trackList)
-{
-
-  if(!trackList)return kFALSE;
-  fRPAngle=0;
-
-  // need to get this info from elsewhere??
-
-  Double_t fPsiRP =0,fDeltaPsiRP = 0;
-   
-   
-    
-  TVector2 mQ,mQ1,mQ2;
-  Float_t mQx= fDeltaQxy[0], mQy=fDeltaQxy[1];
-  
-  Float_t mQx1=fDeltaQxy[0], mQy1=fDeltaQxy[1];
-  Float_t mQx2=fDeltaQxy[0], mQy2=fDeltaQxy[1];
-  
-  AliVParticle *track=0x0;
-  Int_t count[3]={0,0,0};
-  
+Bool_t AliAnalysisTaskJetServices::CalculateReactionPlaneAngleVZERO(AliAODEvent *aod){
 
-  for (Int_t iter=0;iter<trackList->GetEntries();iter++){
-
-    track=(AliVParticle*)trackList->At(iter);
-    
-    //cuts already applied before
-    // Comment DCA not correctly implemented yet for AOD tracks
-    
-    Double_t momentum;
-    if(track->Charge()>0){momentum=track->Pt();}
-    else{momentum=-track->Pt();}
-
-       
-
-    // For Weighting
-    fh3RPPhiTracks->Fill(fCentrality,momentum,track->Phi());
-    count[0]++;
-
-    Double_t phiweight=GetPhiWeight(track->Phi(),momentum);
-    //    Double_t phiweight=1; 
-    Double_t weight=2;
-    if(track->Pt()<2){weight=track->Pt();}
-    
-
-    mQx += (cos(2*track->Phi()))*weight*phiweight;
-    mQy += (sin(2*track->Phi()))*weight*phiweight;
-
-    // Make random Subevents
-
-    if(fRPSubeventMethod==0){
-      if(fRandomizer->Binomial(1,0.5)){
-       mQx1 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy1 += (sin(2*track->Phi()))*weight*phiweight;
-       count[1]++;}
-      else{
-       mQx2 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy2 += (sin(2*track->Phi()))*weight*phiweight;
-       count[2]++;}
-    }
-    else if(fRPSubeventMethod==1){
-      // Make eta dependent subevents
-      if(track->Eta()>0){
-       mQx1 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy1 += (sin(2*track->Phi()))*weight*phiweight;
-       count[1]++;}
-      else{
-       mQx2 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy2 += (sin(2*track->Phi()))*weight*phiweight;
-       count[2]++;}
+  //  const Double_t arr_eta[11]={-3.7, -3.2, -2.7, -2.2, -1.7,0, 2.8, 3.4, 3.9, 4.5,5.1};
+  if(!aod)return kFALSE;
+  static bool bFirst = true;
+  static Double_t v0phi[64] = {0,};
+
+  if(bFirst){
+    int is=0;
+    for(int iArm = 0; iArm<2; iArm++){
+      for(int iRing = 0; iRing<4 ; iRing++){
+       for(int iSec = 0; iSec<8 ; iSec++){
+         v0phi[is] = 22.5 + 45. * iSec;
+         v0phi[is] *= TMath::Pi()/180; 
+         // cout<< is <<" "<< v0phi[is]<<endl;
+         is++;
+       }
+      }
     }
-
+    bFirst = false;
   }
 
-
-
-  //If no track passes the cuts, the ,Q.Phi() will return Pi and a peak at Pi/2 in the RP Angular Distribution will appear
-  if(count[0]==0||count[1]==0||count[2]==0){
-    return kFALSE;
-  }
-
-  mQ.Set(mQx,mQy);
-  mQ1.Set(mQx1,mQy1);
-  mQ2.Set(mQx2,mQy2);
-
-  // cout<<"MQ"<<mQx<<" " <<mQy<<" psi"<<endl;
-
-  fPsiRP=mQ.Phi()/2;
-    
-  //Correction
-  fPsiRP+=fFlatA[0]*TMath::Cos(2*fPsiRP)+fFlatB[0]*TMath::Sin(2*fPsiRP)+fFlatA[1]*TMath::Cos(4*fPsiRP)+fFlatB[1]*TMath::Sin(4*fPsiRP);
-
-  Double_t fPsiRP1=mQ1.Phi()/2;
-  fPsiRP1+=fFlatA[0]*TMath::Cos(2*fPsiRP1)+fFlatB[0]*TMath::Sin(2*fPsiRP1)+fFlatA[1]*TMath::Cos(4*fPsiRP1)+fFlatB[1]*TMath::Sin(4*fPsiRP1);
-  Double_t fPsiRP2=mQ2.Phi()/2;
-  fPsiRP2+=fFlatA[0]*TMath::Cos(2*fPsiRP2)+fFlatB[0]*TMath::Sin(2*fPsiRP2)+fFlatA[1]*TMath::Cos(4*fPsiRP2)+fFlatB[1]*TMath::Sin(4*fPsiRP2);
-  fDeltaPsiRP=fPsiRP1-fPsiRP2;
-  
-  if(fPsiRP>TMath::Pi()){fPsiRP-=TMath::Pi();}
-  if(fPsiRP<0){fPsiRP+=TMath::Pi();}
-  
-  // reactionplaneangle + Pi() is the same angle
-  if(TMath::Abs(fDeltaPsiRP)>TMath::Pi()/2){
-    if(fDeltaPsiRP>0)fDeltaPsiRP-=TMath::Pi();
-    else fDeltaPsiRP+=TMath::Pi();
-  }
+  // 
+  const AliAODVZERO *aodVZERO = aod->GetVZEROData();
+  Double_t numYZNA = 0,numXZNA = 0,sumZNA = 0;
+  Double_t meanXA = 0,meanYA = 0;
+
+  Double_t numYZNC = 0,numXZNC = 0,sumZNC = 0;
+  Double_t meanXC = 0,meanYC = 0;
+
+  for (int i=0; i<64; i++) {  
+    Double_t mult = aodVZERO->GetMultiplicity(i);
+    Double_t phi = v0phi[i];
+    if (mult>0) {
+      if (i<32) { //C-side
+        Double_t wZNC= mult;
+       numYZNC += sin(2.*phi)*wZNC; 
+        numXZNC += cos(2.*phi)*wZNC;
+       sumZNC+=wZNC;
+      }
+      else if(i>31){ //A-side
+       Double_t wZNA=mult;
+       numYZNA += sin(2.*phi)*wZNA; 
+        numXZNA += cos(2.*phi)*wZNA;
+       sumZNA+=wZNA; 
+      } 
+    }// mult>0
+  }// 64 sectors
+
+  Double_t   XC = numXZNC/sumZNC; 
+  Double_t   YC = numYZNC/sumZNC; 
   
-  Double_t cos2deltaRP=TMath::Cos(2*fDeltaPsiRP);
+  Double_t   XA = numXZNA/sumZNA;
+  Double_t   YA = numYZNA/sumZNA;
   
-  // FillHistograms
-  fh2RPSubevents->Fill(fPsiRP1,fPsiRP2);
-  fh1RP->Fill(fPsiRP);
-  fh2RPCentrality->Fill(fCentrality,fPsiRP);
-  fh2RPDeltaRP->Fill(fDeltaPsiRP,fCentrality);
-  fh2RPQxQy->Fill(mQx,mQy);
-  fh2RPCosDeltaRP->Fill(fCentrality,cos2deltaRP);
+
+  fPsiVZEROA = 0.5*TMath::ATan2(YA-meanYA, XA-meanXA);
+  if(fPsiVZEROA>TMath::Pi()){fPsiVZEROA-=TMath::Pi();}
+  if(fPsiVZEROA<0){fPsiVZEROA+=TMath::Pi();}
+
+  fPsiVZEROC = 0.5*TMath::ATan2(YC-meanYC, XA-meanXC);
+  if(fPsiVZEROC>TMath::Pi()){fPsiVZEROC-=TMath::Pi();}
+  if(fPsiVZEROC<0){fPsiVZEROC+=TMath::Pi();}
   
-  fRPAngle=fPsiRP;  
+  fh2XYA->Fill(XA,YA);
+  fp1RPXA->Fill(aod->GetRunNumber(),XA);
+  fp1RPYA->Fill(aod->GetRunNumber(),YA);
+  fh2XYC->Fill(XC,YC);
+  fp1RPXC->Fill(aod->GetRunNumber(),XC);
+  fp1RPYC->Fill(aod->GetRunNumber(),YC);
   return kTRUE;
-}
 
-Double_t AliAnalysisTaskJetServices::GetPhiWeight(Double_t phi,Double_t signedpt){
-  if(!fh3PhiWeights)return 1;
-  else return fh3PhiWeights->GetBinContent(fh3PhiWeights->GetXaxis()->FindBin(fCentrality),fh3PhiWeights->GetYaxis()->FindBin(signedpt),fh3PhiWeights->GetZaxis()->FindBin(phi));
 }
 
- //________________________________________________________________________
-
 Int_t  AliAnalysisTaskJetServices::GetListOfTracks(TList *list){
   Int_t iCount = 0;
   AliAODEvent *aod = 0;
index 958d1c79ab1b71564d9a069b634ec2d08b88de79..6ac81eef7a2fd643903eabeabd53a1699d1a8748 100644 (file)
@@ -18,6 +18,7 @@ class AliESDEvent;
 class AliESDVertex;
 class AliAODEvent;
 class AliAODVertex;
+class AliAODVZERO;
 class AliAODJet;
 class AliGenPythiaEventHeader;
 class AliCFManager;
@@ -80,22 +81,14 @@ class AliAnalysisTaskJetServices : public AliAnalysisTaskSE
     virtual void SetFilterMask(UInt_t i){fFilterMask = i;}
     virtual void SetMinTrackPt(Float_t f){fMinTrackPt = f;}
     virtual void SetTrackEtaWindow(Float_t f){fTrackRecEtaWindow = f;}
+    virtual void SetRPMethod(Int_t i){fRPMethod = i;}
 
-    virtual void SetPhiWeights(TH3F *phiw){fh3PhiWeights = phiw;}
-    virtual void SetFlatteningCoeff(Float_t *fA,Float_t *fB){
-      fFlatA[0] = fA[0];fFlatA[1] = fA[1];
-      fFlatA[0] = fB[0];fFlatB[1] = fB[1];
-    }
-    virtual void SetDeltaQxy(Float_t *fD){
-      fDeltaQxy[0] = fD[0];
-      fDeltaQxy[1] = fD[1];
-    }
-
-    Bool_t   CalculateReactionPlaneAngle(const TList *trackList);
-    Double_t GetPhiWeight(Double_t phi,Double_t signedpt);
+    Bool_t   CalculateReactionPlaneAngleVZERO(AliAODEvent *aod);
     Int_t   GetListOfTracks(TList *list);
 
+
     enum { kAllTriggered = 0,kTriggeredVertex,kTriggeredVertexIn,kSelectedALICE,kSelectedALICEVertexValid,kSelectedALICEVertexIn,kSelected,kConstraints};
+    enum { kRPTracks = 0, kRPVZEROA,kRPVZEROC,kRPMethods};
 
     enum { kNoEventCut=1<<0,
           kPhysicsSelectionCut=1<<1,
@@ -124,7 +117,7 @@ class AliAnalysisTaskJetServices : public AliAnalysisTaskSE
     UInt_t        fSelectionInfoESD;   // slection info bit mask
     UInt_t        fEventCutInfoESD;   // event selection info of what is cutted after physics selection
     UInt_t        fFilterMask;         // filter bit for slecected tracks
-    Int_t         fRPSubeventMethod;   // method for subevent calculation
+    Int_t         fRPMethod;           // method for subevent calculation
     Int_t         fCollisionType;           // type of collisions
     Float_t       fAvgTrials;          // Average number of trials
     Float_t       fVtxXMean;           // mean x for cuts
@@ -140,9 +133,8 @@ class AliAnalysisTaskJetServices : public AliAnalysisTaskSE
     Float_t       fTrackRecEtaWindow;     // eta window for rec tracks
     Float_t       fMinTrackPt;            // limits the track p_T 
     Float_t       fRPAngle;               // ! RP angle of the reaction plane
-    Float_t       fFlatA[2];              // flattening for RP
-    Float_t       fFlatB[2];              // flattening for RP
-    Float_t       fDeltaQxy[2];           // centering of QX QY
+    Float_t       fPsiVZEROA;              // ! RP angle from vzeroa
+    Float_t       fPsiVZEROC;              // ! RP angle from vzeroc
 
     TRandom3      *fRandomizer;        // ! randomizer
 
@@ -162,23 +154,28 @@ class AliAnalysisTaskJetServices : public AliAnalysisTaskSE
     TH2F*         fh2ESDTriggerVtx;    //! vtx. position vs. trigger decision 
     TH2F*         fh2ESDTriggerRun;    //! fired triggers vs. run number
     TH2F*         fh2VtxXY;            //! XY position of VTX were available
-    TH1F*         fh1NCosmicsPerEvent;  //! Number of coscmic candidates found in event
-    TH2F*         fh2RPSubevents;   //! subevent RP 
-    TH2F*         fh2RPCentrality;   //! RP vs centrality
-    TH2F*         fh2RPDeltaRP;     //! centrality vs. RP dela  
-    TH2F*         fh2RPQxQy;        //! QX QY moments
-    TH2F*         fh2RPCosDeltaRP;  //! RP resolution
-
-    TH3F*         fh3PhiWeights;  // RP phi weights, need to be set externally
-    TH3F*         fh3RPPhiTracks; //! RP angle
+    TH1F*         fh1NCosmicsPerEvent; //! Number of coscmic candidates found in event
+    TProfile*     fp1RPXA;              //! mean XA vs run
+    TProfile*     fp1RPYA;              //! mean YA vs run
+    TProfile*     fp1RPXC;              //! mean XA vs run
+    TProfile*     fp1RPYC;              //! mean YA vs run
+    TH2F*         fh2RPAC;              //! RP A vs C 
+    TH2F*         fh2RPAT;              //! RP A vs tracks 
+    TH2F*         fh2RPCT;              //! RP C vs tracks 
+    TH2F*         fh2XYA;               //! XY correlations VZERO C
+    TH2F*         fh2XYC;               //! XY correlations VZERO C
+    TH2F*         fh2RPCentrality;      //! RP vs centrality
+    TH2F*         fh2RPACentrality;     //! RP vs centrality
+    TH2F*         fh2RPCCentrality;     //! RP vs centrality
 
     AliTriggerAnalysis *fTriggerAnalysis; //! Trigger Analysis to get the background rates etc.
     TList *fHistList; //! Output list
 
         // Provisions for replication
     static AliAODHeader*    fgAODHeader;        //! Header for replication
+    static AliAODVZERO*    fgAODVZERO;        //! vzero for replication
     static TClonesArray*  fgAODVertices;        //! primary vertex for replication
-    ClassDef(AliAnalysisTaskJetServices,13)
+    ClassDef(AliAnalysisTaskJetServices,14)
 };
  
 #endif
index 51357d7a54d0ae2158396229fe8a347b1cb29386..49fcf9007a4ff851891b6aa1f2583ae0a5fb845b 100644 (file)
@@ -32,6 +32,7 @@
 #include <TH2F.h>
 #include <TH3F.h>
 #include <TProfile.h>
+#include <TProfile2D.h>
 #include <TList.h>
 #include <TLorentzVector.h>
 #include <TClonesArray.h>
@@ -98,7 +99,7 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2():
   fTrackTypeRec(kTrackUndef),
   fTrackTypeGen(kTrackUndef),
   fEventClass(0),
-  fRPSubeventMethod(0),
+  fRPMethod(0),
   fAvgTrials(1),
   fExternalWeight(1),    
   fJetRecEtaWindow(0.5),
@@ -121,24 +122,12 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2():
   fh1TmpRho(0x0),
   fh2MultRec(0x0),
   fh2MultGen(0x0),
-  fh2RPSubevents(0x0),
   fh2RPCentrality(0x0),
-  fh2RPDeltaRP(0x0),
-  fh2RPQxQy(0x0),
-  fh2RPCosDeltaRP(0x0),
   fh2PtFGen(0x0),
   fh2RelPtFGen(0x0),
-  fh3PhiWeights(0x0),
-  fh3RPPhiTracks(0x0),
   fHistList(0x0)  
 {
 
-  fFlatA[0] =   fFlatA[1] = 0;
-  fFlatB[0] =   fFlatB[1] = 0;
-  fDeltaQxy[0] =   fDeltaQxy[1] = 0; 
-
-
   for(int ij = 0;ij <kJetTypes;++ij){    
     fFlagJetType[ij] = 1; // default = on
     fh1NJets[ij] = 0;
@@ -148,6 +137,8 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2():
     fh1PtTracksInLow[ij] = 0;
     fh2NJetsPt[ij]  = 0;
     fh2NTracksPt[ij]  = 0;
+    fp2MultRPPhiTrackPt[ij] = 0;
+    fp2CentRPPhiTrackPt[ij] = 0;
     fhnJetPt[ij] = 0;
     fhnJetPtQA[ij] = 0;
     fhnTrackPt[ij] = 0;
@@ -198,7 +189,7 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name):
   fTrackTypeRec(kTrackUndef),
   fTrackTypeGen(kTrackUndef),
   fEventClass(0),
-  fRPSubeventMethod(0),
+  fRPMethod(0),
   fAvgTrials(1),
   fExternalWeight(1),    
   fJetRecEtaWindow(0.5),
@@ -221,22 +212,12 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name):
   fh1TmpRho(0x0),
   fh2MultRec(0x0),
   fh2MultGen(0x0),
-  fh2RPSubevents(0x0),
   fh2RPCentrality(0x0),
-  fh2RPDeltaRP(0x0),
-  fh2RPQxQy(0x0),
-  fh2RPCosDeltaRP(0x0),
   fh2PtFGen(0x0),
   fh2RelPtFGen(0x0),
-  fh3PhiWeights(0x0),
-  fh3RPPhiTracks(0x0),
   fHistList(0x0)
 {
 
-  fFlatA[0] =   fFlatA[1] = 0;
-  fFlatB[0] =   fFlatB[1] = 0;
-  fDeltaQxy[0] =   fDeltaQxy[1] = 0; 
-
   for(int ij = 0;ij <kJetTypes;++ij){    
     fFlagJetType[ij] = 1; // default = on
     fh1NJets[ij] = 0;
@@ -246,6 +227,8 @@ AliAnalysisTaskJetSpectrum2::AliAnalysisTaskJetSpectrum2(const char* name):
     fh1PtTracksInLow[ij] = 0;
     fh2NJetsPt[ij]  = 0;
     fh2NTracksPt[ij]  = 0;
+    fp2MultRPPhiTrackPt[ij] = 0;
+    fp2CentRPPhiTrackPt[ij] = 0;
     fhnJetPt[ij] = 0;
     fhnJetPtQA[ij] = 0;
     fhnTrackPt[ij] = 0;
@@ -411,6 +394,7 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
   fh1ZVtx = new TH1F("fh1ZVtx","z vtx;z_{vtx} (cm)",400,-20,20);
   fHistList->Add(fh1ZVtx);
 
+
   fh1RP = new TH1F("fh1RP","RP;#Psi",440, -1.*TMath::Pi(), 2.*TMath::Pi());
   fHistList->Add(fh1RP);
 
@@ -422,30 +406,15 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
   fh2MultGen = new TH2F("fh2MultGen","multiplicity gen;# tracks;# jetrefs",400,0,5000,500,0.,5000);
   fHistList->Add(fh2MultGen);
 
-  fh2RPSubevents = new TH2F("fh2RPSubevents" ,"Reaction Plane Angle" , 180, 0, TMath::Pi(), 180, 0, TMath::Pi());
-  fHistList->Add( fh2RPSubevents);
 
   fh2RPCentrality = new TH2F("fh2RPCentrality" ,"Reaction Plane Angle" , 20, 0.,100., 180, 0, TMath::Pi());
   fHistList->Add(fh2RPCentrality);
 
-  fh2RPDeltaRP   = new TH2F("fh2DeltaRP" ,"Delta Reaction Plane Angle" , 100, -TMath::Pi()/2, TMath::Pi()/2,20,0.,100.0);
-  fHistList->Add(fh2RPDeltaRP);
-
-  fh2RPQxQy      = new TH2F("fh2RPQxQy" ,"" , 100, -100,100,100,-100,100);
-  fHistList->Add(fh2RPQxQy);
-
-  fh2RPCosDeltaRP = new TH2F("fh2RPCosDeltaRP" ,"" , 20, 0.001,100.001,100,-1,1);
-  fHistList->Add(fh2RPCosDeltaRP);
-
-
   fh2PtFGen = new TH2F("fh2PtFGen",Form("%s vs. %s;p_{T,gen};p_{T,rec}",fBranchRec.Data(),fBranchGen.Data()),nBinPt,binLimitsPt,nBinPt,binLimitsPt);
   fHistList->Add(fh2PtFGen);
 
   fh2RelPtFGen = new TH2F("fh2RelPtFGen",";p_{T,gen};p_{T,rec}-p_{T,gen}/p_{T,Gen}",nBinPt,binLimitsPt,241,-2.41,2.41);
   fHistList->Add(fh2RelPtFGen);
-
-  fh3RPPhiTracks = new TH3F("fh3RPPhiTracks","Phi Tracks Pt Centrality", 10, 0.,100.,20,-5,5,200, 0, 2*TMath::Pi());
-  fHistList->Add(fh3RPPhiTracks);
   
   for(int ij = 0;ij <kJetTypes;++ij){    
     TString cAdd = "";
@@ -479,7 +448,7 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
     fh1PtTracksIn[ij] = new TH1F(Form("fh1PtTracks%sIn",cAdd.Data()),Form("%s track p_T;p_{T} (GeV/c)",cAdd.Data()),nBinPt,binLimitsPt);
     fHistList->Add(fh1PtTracksIn[ij]);
     
-    fh1PtTracksInLow[ij] = new TH1F(Form("fh1PtTracks%sInLow",cAdd.Data()),Form("%s track p_T;p_{T} (GeV/c)",cAdd.Data()),50,0.,5.);
+    fh1PtTracksInLow[ij] = new TH1F(Form("fh1PtTracks%sInLow",cAdd.Data()),Form("%s track p_T;p_{T} (GeV/c)",cAdd.Data()),100,0.,5.);
     fHistList->Add(fh1PtTracksInLow[ij]);
     
     fh1SumPtTrack[ij] = new TH1F(Form("fh1SumPtTrack%s",cAdd.Data()),Form("Sum %s track p_T;p_{T} (GeV/c)",cAdd.Data()),1000,0.,3000.);
@@ -490,19 +459,24 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
     
     fh2NTracksPt[ij]  = new TH2F(Form("fh2N%sTracksPt",cAdd.Data()),Form("Number of %s tracks above threshhold;p_{T,cut} (GeV/c);N_{tracks}",cAdd.Data()),nBinPt,binLimitsPt,1000,0.,5000);
     fHistList->Add(fh2NTracksPt[ij]);
-    
-
-      // Bins:  Jet number: pTJet, cent, mult, RP, Area.   total bins = 4.5M
-      const Int_t nBinsSparse1 = 6;
-      const Int_t nBins1[nBinsSparse1] = {     kMaxJets+1,100, 10,  25,    fNRPBins, 10};
-      const Double_t xmin1[nBinsSparse1]  = {        -0.5,  0,  0,   0,        -0.5, 0.};
-      const Double_t xmax1[nBinsSparse1]  = {kMaxJets+0.5,250,100,5000,fNRPBins-0.5,1.0};
-
-      fhnJetPt[ij] = new THnSparseF(Form("fhnJetPt%s",cAdd.Data()),";jet number;p_{T,jet};cent;# tracks;RP;area",nBinsSparse1,nBins1,xmin1,xmax1);
-      fHistList->Add(fhnJetPt[ij]);
 
-      // Bins:  Jet number: pTJet, cent, eta, phi, Area.   total bins = 9.72 M
-      const Int_t nBinsSparse2 = 6;
+    
+    fp2MultRPPhiTrackPt[ij] = new TProfile2D(Form("fp2MultRPPhiTrackPt%s",cAdd.Data()),"RP phi vs Number of tracks;# tracks;#Delta#phi_{RP}; <p_{T}>",20,0,4000,181,-1./180.*TMath::Pi(),TMath::Pi(),"S");
+    fHistList->Add(fp2MultRPPhiTrackPt[ij]);
+    fp2CentRPPhiTrackPt[ij] = new TProfile2D(Form("fp2CentRPPhiTrackPt%s",cAdd.Data()),"RP phi vs cent;# cent;#Delta#phi_{RP}; <p_{T}>",10,0,100,181,-1./180.*TMath::Pi(),TMath::Pi(),"S");
+    fHistList->Add(fp2CentRPPhiTrackPt[ij]);    
+
+    // Bins:  Jet number: pTJet, cent, mult, RP, Area.   total bins = 4.5M
+    const Int_t nBinsSparse1 = 6;
+    const Int_t nBins1[nBinsSparse1] = {     kMaxJets+1,100, 10,  25,    fNRPBins, 10};
+    const Double_t xmin1[nBinsSparse1]  = {        -0.5,  0,  0,   0,        -0.5, 0.};
+    const Double_t xmax1[nBinsSparse1]  = {kMaxJets+0.5,250,100,5000,fNRPBins-0.5,1.0};
+    
+    fhnJetPt[ij] = new THnSparseF(Form("fhnJetPt%s",cAdd.Data()),";jet number;p_{T,jet};cent;# tracks;RP;area",nBinsSparse1,nBins1,xmin1,xmax1);
+    fHistList->Add(fhnJetPt[ij]);
+    
+    // Bins:  Jet number: pTJet, cent, eta, phi, Area.   total bins = 9.72 M
+    const Int_t nBinsSparse2 = 6;
       const Int_t nBins2[nBinsSparse2] = {     kMaxJets+1, 25,   5,  18,             360, 10};
       const Double_t xmin2[nBinsSparse2]  = {        -0.5,  0,   0,-0.9,              0,  0.};
       const Double_t xmax2[nBinsSparse2]  = {kMaxJets+0.5,250, 100, 0.9, 2.*TMath::Pi(),1.0};
@@ -511,19 +485,20 @@ void AliAnalysisTaskJetSpectrum2::UserCreateOutputObjects()
 
       // Bins:track number  pTtrack, cent, mult, RP.   total bins = 224 k
       const Int_t nBinsSparse3 = 5;
-      const Int_t nBins3[nBinsSparse3] = {       2,  75,     10,  25,    fNRPBins};
-      const Double_t xmin3[nBinsSparse3]  = { -0.5,   0,   0,      0,        -0.5};
-      const Double_t xmax3[nBinsSparse3]  = { 1.5,  200, 100,   5000,fNRPBins-0.5};  
+      const Int_t nBins3[nBinsSparse3] = {       2,    100,     10,  20,    fNRPBins};
+      const Double_t xmin3[nBinsSparse3]  = { -0.5,     0,   0,      0,        -0.5};
+      const Double_t xmax3[nBinsSparse3]  = { 1.5,    200, 100,   4000,fNRPBins-0.5};  
 
       // change the binning ot the pT axis:
       Double_t *xPt3 = new Double_t[nBins3[1]+1];
       xPt3[0] = 0.;
       for(int i = 1; i<=nBins3[1];i++){
-       if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.1;
-       else if(xPt3[i-1]<10)xPt3[i] = xPt3[i-1] + 0.5;
-       else if(xPt3[i-1]<20)xPt3[i] = xPt3[i-1] +  1.;
-       else if(xPt3[i-1]<30)xPt3[i] = xPt3[i-1] + 2.5;
-       else xPt3[i] = xPt3[i-1] + 5.;
+       if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.05; // 1 - 40
+       else if(xPt3[i-1]<4)xPt3[i] = xPt3[i-1] + 0.2; // 41 - 50
+       else if(xPt3[i-1]<10)xPt3[i] = xPt3[i-1] + 0.5; // 50 - 62
+       else if(xPt3[i-1]<20)xPt3[i] = xPt3[i-1] +  1.; // 62 - 72
+       else if(xPt3[i-1]<30)xPt3[i] = xPt3[i-1] + 2.5; // 74 - 78
+       else xPt3[i] = xPt3[i-1] + 5.; // 78 - 100 = 140 
       }
       
       fhnTrackPt[ij] = new THnSparseF(Form("fhnTrackPt%s",cAdd.Data()),";track number;p_{T};cent;#tracks;RP",nBinsSparse3,nBins3,xmin3,xmax3);
@@ -813,7 +788,14 @@ void AliAnalysisTaskJetSpectrum2::UserExec(Option_t */*option*/){
   if(fDebug>2)Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__,nT,genParticles.GetEntries());
 
   //  CalculateReactionPlaneAngle(&recParticles);
-  fRPAngle = aod->GetHeader()->GetEventplane();
+  fRPAngle = 0;
+  
+  if(fRPMethod==0)fRPAngle = aod->GetHeader()->GetEventplane();
+  else if(fRPMethod==1||fRPMethod==2){
+    fRPAngle = aod->GetHeader()->GetQTheta(fRPMethod);
+  }
+  fh1RP->Fill(fRPAngle);
+  fh2RPCentrality->Fill(fCentrality,fRPAngle);
   // Event control and counting ...  
   // MC
   fh1PtHard->Fill(ptHard,eventW);
@@ -1082,6 +1064,20 @@ void AliAnalysisTaskJetSpectrum2::FillTrackHistos(TList &particlesList,int iType
       sumPt += tmpPt;
       Float_t tmpPhi = tmpTrack->Phi();
       if(tmpPhi<0)tmpPhi+=TMath::Pi()*2.;    
+
+      
+      Float_t phiRP = tmpPhi-fRPAngle;
+      if(phiRP>TMath::Pi())phiRP -= TMath::Pi();
+      if(phiRP<0)phiRP += TMath::Pi();
+      if(phiRP<0)phiRP += TMath::Pi();
+      const float allPhi = -1./180.*TMath::Pi();
+
+      fp2MultRPPhiTrackPt[iType]->Fill(refMult,phiRP,tmpPt);
+      fp2MultRPPhiTrackPt[iType]->Fill(refMult,allPhi,tmpPt);
+
+      fp2CentRPPhiTrackPt[iType]->Fill(fCentrality,phiRP,tmpPt);
+      fp2CentRPPhiTrackPt[iType]->Fill(fCentrality,allPhi,tmpPt);
+
       Int_t phiBin = GetPhiBin(tmpPhi-fRPAngle);
       var3[0] = 1;
       var3[1] = tmpPt;
@@ -1485,128 +1481,6 @@ AliVParticle *AliAnalysisTaskJetSpectrum2::LeadingTrackInCone(AliAODJet* jet,TLi
   return leading;
 }
 
-Bool_t AliAnalysisTaskJetSpectrum2::CalculateReactionPlaneAngle(const TList *trackList)
-{
-
-  if(!trackList)return kFALSE;
-  fRPAngle=0;
-
-  // need to get this info from elsewhere??
-
-  Double_t fPsiRP =0,fDeltaPsiRP = 0;
-   
-   
-    
-  TVector2 mQ,mQ1,mQ2;
-  Float_t mQx= fDeltaQxy[0], mQy=fDeltaQxy[1];
-  
-  Float_t mQx1=fDeltaQxy[0], mQy1=fDeltaQxy[1];
-  Float_t mQx2=fDeltaQxy[0], mQy2=fDeltaQxy[1];
-  
-  AliVParticle *track=0x0;
-  Int_t count[3]={0,0,0};
-  
-
-  for (Int_t iter=0;iter<trackList->GetEntries();iter++){
-
-    track=(AliVParticle*)trackList->At(iter);
-    
-    //cuts already applied before
-    // Comment DCA not correctly implemented yet for AOD tracks
-    
-    Double_t momentum;
-    if(track->Charge()>0){momentum=track->Pt();}
-    else{momentum=-track->Pt();}
-
-       
-
-    // For Weighting
-    fh3RPPhiTracks->Fill(fCentrality,momentum,track->Phi());
-    count[0]++;
-
-    Double_t phiweight=GetPhiWeight(track->Phi(),momentum);
-    //    Double_t phiweight=1; 
-    Double_t weight=2;
-    if(track->Pt()<2){weight=track->Pt();}
-    
-
-    mQx += (cos(2*track->Phi()))*weight*phiweight;
-    mQy += (sin(2*track->Phi()))*weight*phiweight;
-
-    // Make random Subevents
-
-    if(fRPSubeventMethod==0){
-      if(fRandomizer->Binomial(1,0.5)){
-       mQx1 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy1 += (sin(2*track->Phi()))*weight*phiweight;
-       count[1]++;}
-      else{
-       mQx2 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy2 += (sin(2*track->Phi()))*weight*phiweight;
-       count[2]++;}
-    }
-    else if(fRPSubeventMethod==1){
-      // Make eta dependent subevents
-      if(track->Eta()>0){
-       mQx1 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy1 += (sin(2*track->Phi()))*weight*phiweight;
-       count[1]++;}
-      else{
-       mQx2 += (cos(2*track->Phi()))*weight*phiweight;
-       mQy2 += (sin(2*track->Phi()))*weight*phiweight;
-       count[2]++;}
-    }
-
-  }
-
-
-
-  //If no track passes the cuts, the ,Q.Phi() will return Pi and a peak at Pi/2 in the RP Angular Distribution will appear
-  if(count[0]==0||count[1]==0||count[2]==0){
-    return kFALSE;
-  }
-
-  mQ.Set(mQx,mQy);
-  mQ1.Set(mQx1,mQy1);
-  mQ2.Set(mQx2,mQy2);
-
-  // cout<<"MQ"<<mQx<<" " <<mQy<<" psi"<<endl;
-
-  fPsiRP=mQ.Phi()/2;
-    
-  //Correction
-  fPsiRP+=fFlatA[0]*TMath::Cos(2*fPsiRP)+fFlatB[0]*TMath::Sin(2*fPsiRP)+fFlatA[1]*TMath::Cos(4*fPsiRP)+fFlatB[1]*TMath::Sin(4*fPsiRP);
-
-  Double_t fPsiRP1=mQ1.Phi()/2;
-  fPsiRP1+=fFlatA[0]*TMath::Cos(2*fPsiRP1)+fFlatB[0]*TMath::Sin(2*fPsiRP1)+fFlatA[1]*TMath::Cos(4*fPsiRP1)+fFlatB[1]*TMath::Sin(4*fPsiRP1);
-  Double_t fPsiRP2=mQ2.Phi()/2;
-  fPsiRP2+=fFlatA[0]*TMath::Cos(2*fPsiRP2)+fFlatB[0]*TMath::Sin(2*fPsiRP2)+fFlatA[1]*TMath::Cos(4*fPsiRP2)+fFlatB[1]*TMath::Sin(4*fPsiRP2);
-  fDeltaPsiRP=fPsiRP1-fPsiRP2;
-  
-  if(fPsiRP>TMath::Pi()){fPsiRP-=TMath::Pi();}
-  if(fPsiRP<0){fPsiRP+=TMath::Pi();}
-  
-  // reactionplaneangle + Pi() is the same angle
-  if(TMath::Abs(fDeltaPsiRP)>TMath::Pi()/2){
-    if(fDeltaPsiRP>0)fDeltaPsiRP-=TMath::Pi();
-    else fDeltaPsiRP+=TMath::Pi();
-  }
-  
-  Double_t cos2deltaRP=TMath::Cos(2*fDeltaPsiRP);
-  
-  // FillHistograms
-  fh2RPSubevents->Fill(fPsiRP1,fPsiRP2);
-  fh1RP->Fill(fPsiRP);
-  fh2RPCentrality->Fill(fCentrality,fPsiRP);
-  fh2RPDeltaRP->Fill(fDeltaPsiRP,fCentrality);
-  fh2RPQxQy->Fill(mQx,mQy);
-  fh2RPCosDeltaRP->Fill(fCentrality,cos2deltaRP);
-  
-  fRPAngle=fPsiRP;
-  
-  return kTRUE;
-}
-
 Int_t AliAnalysisTaskJetSpectrum2::GetPhiBin(Double_t phi)
 {
     Int_t phibin=-1;
@@ -1617,11 +1491,3 @@ Int_t AliAnalysisTaskJetSpectrum2::GetPhiBin(Double_t phi)
     return phibin;
 }
 
-
-
-Double_t AliAnalysisTaskJetSpectrum2::GetPhiWeight(Double_t phi,Double_t signedpt){
-  if(!fh3PhiWeights)return 1;
-  else return fh3PhiWeights->GetBinContent(fh3PhiWeights->GetXaxis()->FindBin(fCentrality),fh3PhiWeights->GetYaxis()->FindBin(signedpt),fh3PhiWeights->GetZaxis()->FindBin(phi));
-}
-
- //________________________________________________________________________
index ab5dc4791136fc1df3572ace3fce6c1e2dddff54..9330fdf6cb7bedc6a1f1d693a311c57fcffa5e78 100644 (file)
@@ -32,6 +32,7 @@ class TH2F;
 class TH3F;
 class TRandom3;
 class TProfile;
+class TProfile2D;
 class TSTring;
 
 
@@ -74,16 +75,9 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     virtual void SetTrackTypeRec(Int_t i){fTrackTypeRec = i;}
     virtual void SetFilterMask(UInt_t i){fFilterMask = i;}
     virtual void SetMatching(Bool_t b = kTRUE){fDoMatching = b;}
+    virtual void SetRPMethod(Int_t i){fRPMethod = i;}
     virtual void SetEventSelectionMask(UInt_t i){fEventSelectionMask = i;}
-    virtual void SetPhiWeights(TH3F *phiw){fh3PhiWeights = phiw;}
-    virtual void SetFlatteningCoeff(Float_t *fA,Float_t *fB){
-      fFlatA[0] = fA[0];fFlatA[1] = fA[1];
-      fFlatA[0] = fB[0];fFlatB[1] = fB[1];
-    }
-    virtual void SetDeltaQxy(Float_t *fD){
-      fDeltaQxy[0] = fD[0];
-      fDeltaQxy[1] = fD[1];
-    }
+
     virtual void SetNonStdFile(char* c){fNonStdFile = c;} 
 
 
@@ -166,7 +160,7 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     Int_t         fTrackTypeGen;          // type of tracks used for FF 
     Int_t         fFlagJetType[kJetTypes]; // disable the filling and booking of certain JetType histos
     Int_t         fEventClass;            // event class to be looked at for this instance of the task
-    Int_t         fRPSubeventMethod;      // method for subevent calculation
+    Int_t         fRPMethod;              // method for subevent calculation
     Float_t       fAvgTrials;             // Average nimber of trials
     Float_t       fExternalWeight;        // external weight
     Float_t       fJetRecEtaWindow;       // eta window for rec jets
@@ -176,9 +170,6 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     Float_t       fDeltaPhiWindow;        // minium angle between dijets
     Float_t       fCentrality;            // ! centrality
     Float_t       fRPAngle;               // ! RP angle of the reaction plane
-    Float_t       fFlatA[2];              // flattening for RP
-    Float_t       fFlatB[2];              // flattening for RP
-    Float_t       fDeltaQxy[2];           // centering of QX QY
     Int_t         fMultRec;               // ! reconstructed track multiplicity
     Int_t         fMultGen;               // ! generated track multiplicity
     
@@ -194,17 +185,11 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     TH1F*         fh1TmpRho;        //! just temporary histo for calculation    
     TH2F*         fh2MultRec;       //! reconstructed track multiplicity   
     TH2F*         fh2MultGen;       //! generated track multiplicity   
-    TH2F*         fh2RPSubevents;   //! subevent RP 
     TH2F*         fh2RPCentrality;   //! RP vs centrality
-    TH2F*         fh2RPDeltaRP;     //! centrality vs. RP dela  
-    TH2F*         fh2RPQxQy;        //! QX QY moments
-    TH2F*         fh2RPCosDeltaRP;  //! RP resolution
 
     TH2F*         fh2PtFGen;                //! found vs generated 
     TH2F*         fh2RelPtFGen;             //! relative difference between generated and found 
 
-    TH3F*         fh3PhiWeights;  // RP phi weights, need to be set externally
-    TH3F*         fh3RPPhiTracks; //! RP angle
     
 
     // Jet histos second go
@@ -219,7 +204,8 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     
     TH2F*         fh2NJetsPt[kJetTypes];    //! Number of found jets above threshold
     TH2F*         fh2NTracksPt[kJetTypes];  //! Number of tracks above threshold
-
+    TProfile2D    *fp2MultRPPhiTrackPt[kJetTypes];         //! for mean pT vs RP   
+    TProfile2D    *fp2CentRPPhiTrackPt[kJetTypes];         //! for mean pT vs RP   
     THnSparseF    *fhnJetPt[kJetTypes];                  //! jet pt information for analysis
     THnSparseF    *fhnJetPtQA[kJetTypes];                //! jet pt information for QA
     THnSparseF    *fhnTrackPt[kJetTypes];                //! track pt information for analysis
@@ -237,7 +223,7 @@ class AliAnalysisTaskJetSpectrum2 : public AliAnalysisTaskSE
     TList *fHistList;                  //! Output list
    
 
-    ClassDef(AliAnalysisTaskJetSpectrum2, 16) // Analysis task for standard jet analysis
+    ClassDef(AliAnalysisTaskJetSpectrum2, 17) // Analysis task for standard jet analysis
 };
  
 #endif
index 51a17ade2931e12639dfb22d66405773fc162ed8..11ad5e88e6d14ae003c4f7e3f1393b2a501d43b0 100644 (file)
@@ -94,10 +94,6 @@ AliAnalysisTaskJetSpectrum2 *AddTaskJetSpectrum2(const char* bRec,const char* bG
    pwg4spec->SetMinJetPt(5.);\r
    pwg4spec->SetJetEtaWindow(0.4);\r
 \r
-   Float_t fDQxy[2] = {-0.398,-0.379};\r
-   pwg4spec->SetDeltaQxy(fDQxy);\r
-\r
-   \r
 \r
 \r
    if(type == "AOD"){\r
index edfd341aa04ab298445fa11f1614009787339aac..d10ee2a9dae5cd64036edaea9196fca01b9b89a9 100644 (file)
@@ -91,9 +91,11 @@ Int_t       kJetMapOffset[3] = {10000,100,1};
 TString     kDefaultJetBranch     = "";      // is currently set when filled (iJETAN or clusters) or from config macro 
 TString     kDefaultJetBackgroundBranch            = "";      // is currently set when filled (jet clsuters  
 TString     kDefaultJetBackgroundBranchCut1        = "";      // is currently set when filled (jet clsuters  
+TString     kDefaultJetBackgroundBranchCut2        = "";      // is currently set when filled (jet clsuters  
 TString     kDefaultJetBackgroundBranch_extra     = "";      // is currently set when filled (jet clsuters) 
 TString     kJetSubtractBranches     = "";      // is currently set when filled (jet clsuters  
 TString     kJetSubtractBranchesCut1     = "";      // is currently set when filled (jet clsuters  
+TString     kJetSubtractBranchesCut2     = "";      // is currently set when filled (jet clsuters  
 TString     kJetSubtractBranches_extra     = "";      // is currently set when filled (jet clsuters  
 
 TString     kDefaultJetBranchMC     = "";      // is currently set when filled (iJETAN or clusters) or from config macro 
@@ -605,7 +607,7 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
         kJetSubtractBranches += Form("%s%s ",taskCl->GetJetOutputBranch(),"RandomConeSkip00");
         kJetSubtractBranches += Form("%s%s ",taskCl->GetJetOutputBranch(),"RandomCone_random");
 
-        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"KT",0.4,0,1, kDeltaAODJetName.Data(),2.0,fTrackEtaWindow,0); // this one is for the background and random jets, random cones with no skip
+        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"KT",0.4,0,1, kDeltaAODJetName.Data(),1.0,fTrackEtaWindow,0); // this one is for the background and random jets, random cones with no skip
         taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
         taskCl->SetNRandomCones(1);
         taskCl->SetBackgroundCalc(kTRUE);
@@ -619,6 +621,21 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
         kJetListSpectrum.Add(new TObjString(Form("%sRandomConeSkip%02d",taskCl->GetJetOutputBranch(),0)));
         kJetListSpectrum.Add(new TObjString(Form("%sRandomCone_random",taskCl->GetJetOutputBranch())));
 
+        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"KT",0.4,0,1, kDeltaAODJetName.Data(),2.0,fTrackEtaWindow,0); // this one is for the background and random jets, random cones with no skip
+        taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
+        taskCl->SetNRandomCones(1);
+        taskCl->SetBackgroundCalc(kTRUE);
+        taskCl->SetCentralityCut(fCenLo,fCenUp);
+        taskCl->SetGhostEtamax(fTrackEtaWindow);
+        if(iAODanalysis==2)taskCl->SetAODTrackInput(kTRUE);
+        kDefaultJetBackgroundBranchCut2 = Form("%s_%s",AliAODJetEventBackground::StdBranchName(),taskCl->GetJetOutputBranch());
+        kJetSubtractBranchesCut2 += Form("%s%s ",taskCl->GetJetOutputBranch(),"RandomConeSkip00");
+        kJetSubtractBranchesCut2 += Form("%s%s ",taskCl->GetJetOutputBranch(),"RandomCone_random");
+        kJetListSpectrum.Add(new TObjString(taskCl->GetJetOutputBranch()));
+        kJetListSpectrum.Add(new TObjString(Form("%sRandomConeSkip%02d",taskCl->GetJetOutputBranch(),0)));
+        kJetListSpectrum.Add(new TObjString(Form("%sRandomCone_random",taskCl->GetJetOutputBranch())));
+
+
          if (iPWG4FastEmbedding) {
            AliAnalysisTaskJetCluster *taskClEmb = 0;
            taskClEmb = AddTaskJetCluster("AODextra","",kHighPtFilterMask,iPhysicsSelectionFlag,"KT",0.4,0,1, kDeltaAODJetName.Data(),0.15,fTrackEtaWindow); // this one is for the background and random jets
@@ -696,7 +713,7 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
        kJetListSpectrum.Add(new TObjString(Form("%sRandomConeSkip%02d",taskCl->GetJetOutputBranch(),2)));
        kJetListSpectrum.Add(new TObjString(Form("%sRandomCone_random",taskCl->GetJetOutputBranch())));
 
-       taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.4,2,1,kDeltaAODJetName.Data(),2.0);
+       taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.4,2,1,kDeltaAODJetName.Data(),1.0);
        taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
        taskCl->SetCentralityCut(fCenLo,fCenUp);
        if(iAODanalysis==2)taskCl->SetAODTrackInput(kTRUE);
@@ -704,12 +721,23 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
        kJetSubtractBranchesCut1 += Form("%s ",taskCl->GetJetOutputBranch());
        kJetListSpectrum.Add(new TObjString(taskCl->GetJetOutputBranch()));
 
+
+       taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.4,2,1,kDeltaAODJetName.Data(),2.0);
+       taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
+       taskCl->SetCentralityCut(fCenLo,fCenUp);
+       if(iAODanalysis==2)taskCl->SetAODTrackInput(kTRUE);
+       if(kIsPbPb)taskCl->SetBackgroundBranch(kDefaultJetBackgroundBranchCut2.Data());
+       kJetSubtractBranchesCut2 += Form("%s ",taskCl->GetJetOutputBranch());
+       kJetListSpectrum.Add(new TObjString(taskCl->GetJetOutputBranch()));
+
+
        // tmp track qa...
+       /*
        taskCl = AddTaskJetCluster("AOD","",1<<8,iPhysicsSelectionFlag,"ANTIKT",0.4,2,1,kDeltaAODJetName.Data(),2.0);
        taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
        taskCl->SetCentralityCut(fCenLo,fCenUp);
        taskCl->SetFilterMask(1<<4|1<<8,1);
-
+       */
        taskCl = AddTaskJetCluster("AOD","",kHighPtFilterMask,iPhysicsSelectionFlag,"ANTIKT",0.2,0,1,kDeltaAODJetName.Data(),0.15);
        taskCl->SetCentralityCut(fCenLo,fCenUp);
        taskCl->SetEventSelection(kTRUE); // saves some computing time, not all vertices are processed
@@ -811,7 +839,7 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
        //
        // cut1
        Int_t iB = 2;
-       taskSubtract = AddTaskJetBackgroundSubtract(kJetSubtractBranchesCut1,iB,kJetSubtractMask1.Data(),kJetSubtractMask2.Data(),"Cut2000");
+       taskSubtract = AddTaskJetBackgroundSubtract(kJetSubtractBranchesCut1,iB,kJetSubtractMask1.Data(),kJetSubtractMask2.Data(),"Cut1000");
        taskSubtract->SetBackgroundBranch(kDefaultJetBackgroundBranchCut1.Data());              
        taskSubtract->SelectCollisionCandidates(iPhysicsSelectionFlag);
        if(kDeltaAODJetName.Length()>0)taskSubtract->SetNonStdOutputFile(kDeltaAODJetName.Data());
@@ -824,6 +852,18 @@ void AnalysisTrainPWG4Jets(const char *analysis_mode="local",
         kJetListSpectrum.Add(new TObjString(cTmp.Data()));
        }
 
+       taskSubtract = AddTaskJetBackgroundSubtract(kJetSubtractBranchesCut2,iB,kJetSubtractMask1.Data(),kJetSubtractMask2.Data(),"Cut2000");
+       taskSubtract->SetBackgroundBranch(kDefaultJetBackgroundBranchCut2.Data());              
+       taskSubtract->SelectCollisionCandidates(iPhysicsSelectionFlag);
+       if(kDeltaAODJetName.Length()>0)taskSubtract->SetNonStdOutputFile(kDeltaAODJetName.Data());
+
+       objArr = kJetSubtractBranchesCut2.Tokenize(" ");
+       for(int iJB = 0;iJB<objArr->GetEntries();iJB++){
+        TObjString *ostr = (TObjString*)objArr->At(iJB);
+        cTmp = ostr->GetString().Data();
+        cTmp.ReplaceAll(kJetSubtractMask1.Data(),Form(kJetSubtractMask2.Data(),iB));
+        kJetListSpectrum.Add(new TObjString(cTmp.Data()));
+       }
 
      }
      if(kJetSubtractBranches_extra.Length()){
index 451061a0a4c5201537ada1687c3abe67a45b655f..34a6b4a46599a0d2648e788393f32e222eb18e34 100644 (file)
@@ -39,7 +39,7 @@
   // bextra == 0 4 plus
   // bextra == 1 large pass1 split..
   // bextra == 2 3 plus
-  Int_t bRun = 802; Int_t bExtra = 0;  char* cDate = "110728a";
+  Int_t bRun = 802; Int_t bExtra = 0;  char* cDate = "110816a";
   //  Int_t bRun = 8102; Int_t bExtra = 1;  char* cDate = "110725a";
   iAODanalysis = 0; 
   // 1 == Read Jets and tracks form the input AOD
       30: clustersAOD_KT04_B2_Filter00016_Cut02000_Skip00RandomCone_random
       31: clustersAOD_ANTIKT04_B2_Filter00016_Cut02000_Skip02
             */
-      // UA1
+      /*
+       ############# Possible jet branches ###################
+       1: jetsAOD_UA104_B0_Filter00272_Cut01000
+       2: jetsAOD_UA104_B0_Filter00272_Cut02000
+       3: jetsAOD_UA104_B2_Filter00272_Cut01000
+       4: jetsAOD_UA104_B2_Filter00272_Cut02000
+       5: clustersAOD_KT04_B0_Filter00272_Cut00150_Skip00
+       6: clustersAOD_KT04_B0_Filter00272_Cut00150_Skip00RandomConeSkip00
+       7: clustersAOD_KT04_B0_Filter00272_Cut00150_Skip00RandomCone_random
+       8: clustersAOD_KT04_B0_Filter00272_Cut01000_Skip00
+       9: clustersAOD_KT04_B0_Filter00272_Cut01000_Skip00RandomConeSkip00
+       10: clustersAOD_KT04_B0_Filter00272_Cut01000_Skip00RandomCone_random
+       11: clustersAOD_KT04_B0_Filter00272_Cut02000_Skip00
+       12: clustersAOD_KT04_B0_Filter00272_Cut02000_Skip00RandomConeSkip00
+       13: clustersAOD_KT04_B0_Filter00272_Cut02000_Skip00RandomCone_random
+       14: clustersAOD_KT02_B0_Filter00272_Cut00150_Skip00
+       15: clustersAOD_ANTIKT04_B0_Filter00272_Cut00150_Skip02
+       16: clustersAOD_ANTIKT04_B0_Filter00272_Cut00150_Skip02RandomConeSkip02
+       17: clustersAOD_ANTIKT04_B0_Filter00272_Cut00150_Skip02RandomCone_random
+       18: clustersAOD_ANTIKT04_B0_Filter00272_Cut01000_Skip02
+       19: clustersAOD_ANTIKT04_B0_Filter00272_Cut02000_Skip02
+       20: clustersAOD_ANTIKT02_B0_Filter00272_Cut00150_Skip00
+       21: clustersAOD_KT04_B1_Filter00272_Cut00150_Skip00RandomConeSkip00
+       22: clustersAOD_KT04_B1_Filter00272_Cut00150_Skip00RandomCone_random
+       23: clustersAOD_ANTIKT04_B1_Filter00272_Cut00150_Skip02
+       24: clustersAOD_ANTIKT04_B1_Filter00272_Cut00150_Skip02RandomConeSkip02
+       25: clustersAOD_ANTIKT04_B1_Filter00272_Cut00150_Skip02RandomCone_random
+       26: clustersAOD_ANTIKT02_B1_Filter00272_Cut00150_Skip00
+       27: clustersAOD_KT04_B2_Filter00272_Cut00150_Skip00RandomConeSkip00
+       28: clustersAOD_KT04_B2_Filter00272_Cut00150_Skip00RandomCone_random
+       29: clustersAOD_ANTIKT04_B2_Filter00272_Cut00150_Skip02
+       30: clustersAOD_ANTIKT04_B2_Filter00272_Cut00150_Skip02RandomConeSkip02
+       31: clustersAOD_ANTIKT04_B2_Filter00272_Cut00150_Skip02RandomCone_random
+       32: clustersAOD_ANTIKT02_B2_Filter00272_Cut00150_Skip00
+       33: clustersAOD_KT04_B2_Filter00272_Cut01000_Skip00RandomConeSkip00
+       34: clustersAOD_KT04_B2_Filter00272_Cut01000_Skip00RandomCone_random
+       35: clustersAOD_ANTIKT04_B2_Filter00272_Cut01000_Skip02
+       36: clustersAOD_KT04_B2_Filter00272_Cut02000_Skip00RandomConeSkip00
+       37: clustersAOD_KT04_B2_Filter00272_Cut02000_Skip00RandomCone_random
+       38: clustersAOD_ANTIKT04_B2_Filter00272_Cut02000_Skip02
+      */
+      
+      // in the first map we fill the correlations we want to plot
+      // in the jet back map we associated the branche used for background calculation
+      // to fetch the multiplicity
+
+      // UA1 
       kJetMapSpectrum.Add(4,2);
       kJetBackMapSpectrum.Add(4,8);
 
       // anti kT 150 MeV
-      kJetMapSpectrum.Add(25,12);
-      kJetBackMapSpectrum.Add(25,5);
-      kJetBackMapSpectrum.Add(13,5);
+      kJetMapSpectrum.Add(29,15);
+      kJetBackMapSpectrum.Add(29,5);
+      kJetBackMapSpectrum.Add(15,5);
+
+      // anti kT 1000 MeV
+      kJetMapSpectrum.Add(35,18);
+      kJetBackMapSpectrum.Add(35,8);      
+      kJetBackMapSpectrum.Add(35,8);
 
       // anti kT 2000 MeV
-      kJetMapSpectrum.Add(31,15);
-      kJetBackMapSpectrum.Add(31,8);      
-      kJetBackMapSpectrum.Add(15,8);
+      kJetMapSpectrum.Add(38,19);
+      kJetBackMapSpectrum.Add(35,11);      
+      kJetBackMapSpectrum.Add(19,11);
+
 
       // anti kT 0.2
-      kJetMapSpectrum.Add(28,16);
-      kJetBackMapSpectrum.Add(28,5);      
-      kJetBackMapSpectrum.Add(17,5);
+      kJetMapSpectrum.Add(32,14);
+      kJetBackMapSpectrum.Add(32,5);      
+      kJetBackMapSpectrum.Add(14,5);
 
       // random cones
       kJetMapSpectrum.Add(6,7);