New task for K0sK0s and small changes in Psi2s task
authorMichal Broz <Michal.Broz@cern.ch>
Mon, 20 Jan 2014 15:21:06 +0000 (16:21 +0100)
committerMichal Broz <Michal.Broz@cern.ch>
Mon, 20 Jan 2014 15:21:06 +0000 (16:21 +0100)
PWGUD/CMakelibPWGUDupc.pkg
PWGUD/PWGUDupcLinkDef.h
PWGUD/UPC/AddTaskUpcK0sK0s.C [new file with mode: 0644]
PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.cxx [new file with mode: 0644]
PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.h [new file with mode: 0644]
PWGUD/UPC/AliAnalysisTaskUpcPsi2s.cxx
PWGUD/UPC/AliAnalysisTaskUpcPsi2s.h

index 279f13f..ca606ab 100644 (file)
@@ -27,7 +27,8 @@
 #--------------------------------------------------------------------------------#
 
 set ( SRCS
-  UPC/AliAnalysisTaskUpcPsi2s.cxx)
+  UPC/AliAnalysisTaskUpcPsi2s.cxx
+  UPC/AliAnalysisTaskUpcK0sK0s.cxx)
 
 string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" )
 
index fc19197..f5881f5 100644 (file)
@@ -5,6 +5,7 @@
 #pragma link off all functions;
 
 #pragma link C++ class AliAnalysisTaskUpcPsi2s+;
+#pragma link C++ class AliAnalysisTaskUpcK0sK0s+;
 
 
 
diff --git a/PWGUD/UPC/AddTaskUpcK0sK0s.C b/PWGUD/UPC/AddTaskUpcK0sK0s.C
new file mode 100644 (file)
index 0000000..e81d52b
--- /dev/null
@@ -0,0 +1,40 @@
+AliAnalysisTaskUpcK0sK0s *AddTaskUpcK0sK0s(Bool_t runTree = kTRUE,Bool_t runHist = kTRUE){
+
+  
+  //--- get the current analysis manager ---//
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+      Error("AddTask_UpcPsi2s", "No analysis manager found.");
+      return 0;
+   }
+  
+  // Check the analysis type using the event handlers connected to the analysis manager.
+  //==============================================================================
+  if (!mgr->GetInputEventHandler()) {
+    Error("AddTask_UpcPsi2s", "This task requires an input event handler");
+    return 0;
+  }
+       
+  TString inputDataType = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
+  
+  // Create tasks
+  AliAnalysisTaskUpcK0sK0s *task = new AliAnalysisTaskUpcK0sK0s(inputDataType.Data());
+  task->SetRunTree(runTree);
+  task->SetRunHist(runHist);
+  mgr->AddTask(task);
+
+
+   // Create containers for input/output
+  AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+  AliAnalysisDataContainer *coutput = mgr->CreateContainer("K0sTree", TTree::Class(), AliAnalysisManager::kOutputContainer,Form("%s:K0sK0s", AliAnalysisManager::GetCommonFileName()));
+  AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("Counter", TH1I::Class(), AliAnalysisManager::kOutputContainer, Form("%s:K0sK0s", AliAnalysisManager::GetCommonFileName()));
+  AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("ListHist", TList::Class(), AliAnalysisManager::kOutputContainer, Form("%s:K0sK0s", AliAnalysisManager::GetCommonFileName()));
+
+  // Connect input/output
+  mgr->ConnectInput(task, 0, cinput);
+  mgr->ConnectOutput(task, 1, coutput);
+  mgr->ConnectOutput(task, 2, coutput2);
+  mgr->ConnectOutput(task, 3, coutput3);
+
+return task;
+}
diff --git a/PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.cxx b/PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.cxx
new file mode 100644 (file)
index 0000000..33f8586
--- /dev/null
@@ -0,0 +1,609 @@
+/*************************************************************************
+* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
+*                                                                        *
+* Author: The ALICE Off-line Project.                                    *
+* Contributors are mentioned in the code where appropriate.              *
+*                                                                        *
+* Permission to use, copy, modify and distribute this software and its   *
+* documentation strictly for non-commercial purposes is hereby granted   *
+* without fee, provided that the above copyright notice appears in all   *
+* copies and that both the copyright notice and this permission notice   *
+* appear in the supporting documentation. The authors make no claims     *
+* about the suitability of this software for any purpose. It is          *
+* provided "as is" without express or implied warranty.                  * 
+**************************************************************************/
+
+// c++ headers
+#include <iostream>
+#include <string.h>
+
+// root headers
+#include "TH1I.h"
+#include "TTree.h"
+#include "TClonesArray.h"
+#include "TParticle.h"
+#include "TObjString.h"
+#include "TFile.h"
+#include "TDatabasePDG.h"
+#include "TLorentzVector.h"
+
+// aliroot headers
+#include "AliAnalysisManager.h"
+#include "AliInputEventHandler.h"
+#include "AliESDEvent.h"
+#include "AliAODEvent.h"
+#include "AliMCEvent.h"
+#include "AliAODVZERO.h"
+#include "AliAODZDC.h"
+#include "AliESDVZERO.h"
+#include "AliESDZDC.h"
+#include "AliPIDResponse.h"
+#include "AliAODTrack.h"
+#include "AliAODPid.h"
+#include "AliAODVertex.h"
+#include "AliESDVertex.h"
+#include "AliMultiplicity.h"
+#include "AliESDtrack.h"
+#include "AliESDMuonTrack.h"
+#include "AliAODMCParticle.h"
+#include "AliMCParticle.h"
+
+// my headers
+#include "AliAnalysisTaskUpcK0sK0s.h"
+
+ClassImp(AliAnalysisTaskUpcK0sK0s);
+
+using std::cout;
+using std::endl;
+
+//trees for UPC analysis,
+// michal.broz@cern.ch
+
+//_____________________________________________________________________________
+AliAnalysisTaskUpcK0sK0s::AliAnalysisTaskUpcK0sK0s() 
+  : AliAnalysisTaskSE(),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),hCounter(0),fK0sTree(0),
+    fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),fVtxContrib(0),fBCrossNum(0),fNtracklets(0),
+    fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
+    fDataFilnam(0),fRecoPass(0),fEvtNum(0),
+    fK0sAODv0s(0),fK0sESDv0s(0),fK0sAODTracks(0),fK0sESDTracks(0),
+    fListHist(0),fHistTriggersPerRun(0),fHistK0sMass(0)
+
+{
+
+//Dummy constructor
+
+}//AliAnalysisTaskUpcK0sK0s
+
+
+//_____________________________________________________________________________
+AliAnalysisTaskUpcK0sK0s::AliAnalysisTaskUpcK0sK0s(const char *name) 
+  : AliAnalysisTaskSE(name),fType(0),fRunTree(kTRUE),fRunHist(kTRUE),hCounter(0),fK0sTree(0),
+    fRunNum(0),fPerNum(0),fOrbNum(0),fL0inputs(0),fL1inputs(0),fVtxContrib(0),fBCrossNum(0),fNtracklets(0),
+    fZDCAenergy(0),fZDCCenergy(0),fV0Adecision(0),fV0Cdecision(0),
+    fDataFilnam(0),fRecoPass(0),fEvtNum(0),
+    fK0sAODv0s(0),fK0sESDv0s(0),fK0sAODTracks(0),fK0sESDTracks(0),
+    fListHist(0),fHistTriggersPerRun(0),fHistK0sMass(0)
+
+{
+
+  // Constructor
+  if( strstr(name,"ESD") ) fType = 0;
+  if( strstr(name,"AOD") ) fType = 1;
+  
+  Init();
+
+  DefineOutput(1, TTree::Class());
+  DefineOutput(2, TH1I::Class());
+  DefineOutput(3, TList::Class());
+
+}//AliAnalysisTaskUpcK0sK0s
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::Init()
+{
+  
+  for(Int_t i=0; i<ntrg; i++) fTrigger[i] = kFALSE;
+
+}//Init
+
+//_____________________________________________________________________________
+AliAnalysisTaskUpcK0sK0s::~AliAnalysisTaskUpcK0sK0s() 
+{
+  // Destructor
+  if(fK0sTree){
+     delete fK0sTree;
+     fK0sTree = 0x0;
+  }
+  if(hCounter){
+     delete hCounter;
+     hCounter = 0x0;
+  }
+  if(fListHist){
+     delete fListHist;
+     fListHist = 0x0;
+  }
+
+}//~AliAnalysisTaskUpcK0sK0s
+
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::UserCreateOutputObjects()
+{
+   hCounter = new TH1I("hCounter", "hCounter", 34000, 1., 34001.);
+
+  //input file
+  fDataFilnam = new TObjString();
+  fDataFilnam->SetString("");
+
+    //vertices
+  fK0sAODv0s = new TClonesArray("AliAODv0", 1000);
+  fK0sESDv0s = new TClonesArray("AliESDv0", 1000);
+  
+    //tracks
+  fK0sAODTracks = new TClonesArray("AliAODTrack", 1000);
+  fK0sESDTracks = new TClonesArray("AliESDtrack", 1000);
+
+  //output tree with K0s candidate events
+  fK0sTree = new TTree("fK0sTree", "fK0sTree");
+  fK0sTree ->Branch("fRunNum", &fRunNum, "fRunNum/I");
+  fK0sTree ->Branch("fPerNum", &fPerNum, "fPerNum/i");
+  fK0sTree ->Branch("fOrbNum", &fOrbNum, "fOrbNum/i");
+  
+  fK0sTree ->Branch("fBCrossNum", &fBCrossNum, "fBCrossNum/s");
+  fK0sTree ->Branch("fTrigger", &fTrigger[0], Form("fTrigger[%i]/O", ntrg));
+  fK0sTree ->Branch("fL0inputs", &fL0inputs, "fL0inputs/i");
+  fK0sTree ->Branch("fL1inputs", &fL1inputs, "fL1inputs/i");
+  fK0sTree ->Branch("fNtracklets", &fNtracklets, "fNtracklets/s");
+  fK0sTree ->Branch("fVtxContrib", &fVtxContrib, "fVtxContrib/I");
+  fK0sTree ->Branch("fZDCAenergy", &fZDCAenergy, "fZDCAenergy/D");
+  fK0sTree ->Branch("fZDCCenergy", &fZDCCenergy, "fZDCCenergy/D");
+  fK0sTree ->Branch("fV0Adecision", &fV0Adecision, "fV0Adecision/I");
+  fK0sTree ->Branch("fV0Cdecision", &fV0Cdecision, "fV0Cdecision/I");  
+  fK0sTree ->Branch("fDataFilnam", &fDataFilnam);
+  fK0sTree ->Branch("fRecoPass", &fRecoPass, "fRecoPass/S");
+  fK0sTree ->Branch("fEvtNum", &fEvtNum, "fEvtNum/L");                        
+  if( fType == 0 ) {
+    fK0sTree ->Branch("fK0sESDv0s", &fK0sESDv0s);
+    fK0sTree ->Branch("fK0sESDTracks", &fK0sESDTracks);
+  }
+  if( fType == 1 ) {
+    fK0sTree ->Branch("fK0sAODv0s", &fK0sAODv0s);
+    fK0sTree ->Branch("fK0sAODTracks", &fK0sAODTracks);
+  }
+  
+  fListHist = new TList();
+  fListHist ->SetOwner();
+  
+  fHistTriggersPerRun = new TH1D("fHistTriggersPerRun", "fHistTriggersPerRun", 3000, 167000.5, 170000.5);
+  fListHist->Add(fHistTriggersPerRun);
+    
+  fHistK0sMass = new TH1D("fHistK0sMass","fHistK0sMass",200,0.4,0.6);
+  fListHist->Add(fHistK0sMass);
+  
+  PostData(1, fK0sTree);
+  PostData(2, hCounter);
+  PostData(3, fListHist);
+
+}//UserCreateOutputObjects
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::UserExec(Option_t *) 
+{
+
+  //cout<<"#################### Next event ##################"<<endl;
+
+  if( fType == 0 ) RunESD();
+  if( fType == 1 ){ 
+       if(fRunHist) RunAODhist();
+       if(fRunTree) RunAODtree();
+       }
+
+}//UserExec
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::RunAODhist()
+{
+
+  //input event
+  AliAODEvent *aod = (AliAODEvent*) InputEvent();
+  if(!aod) return;
+
+
+  //Trigger
+  TString trigger = aod->GetFiredTriggerClasses();
+  
+  if( !trigger.Contains("CCUP4-B") ) return;
+  
+  fRunNum = aod ->GetRunNumber();
+  fHistTriggersPerRun->Fill(fRunNum);
+
+
+  //primary vertex
+  AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
+  fVtxContrib = fAODVertex->GetNContributors();
+  if(fVtxContrib < 2) return;
+
+
+  //VZERO, ZDC
+  AliAODVZERO *fV0data = aod ->GetVZEROData();
+  
+  fV0Adecision = fV0data->GetV0ADecision();
+  fV0Cdecision = fV0data->GetV0CDecision();
+  if(fV0Adecision != AliAODVZERO::kV0Empty || fV0Cdecision != AliAODVZERO::kV0Empty) return;
+
+  AliAODZDC *fZDCdata = aod->GetZDCData();
+  
+  fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
+  fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
+  if( fZDCAenergy > 6000 || fZDCCenergy > 6000) return;
+  
+  Int_t nGoodV0s=0;
+  Int_t V0Index[3] = {-1,-1,-1};
+  Int_t V0TrackID[6] = {-1,-1,-1,-1,-1,-1};
+
+  Int_t nGoodTracks=0;
+  Int_t TrackID[5] = {-1,-1,-1,-1,-1};
+
+
+  //V0s loop
+  for(Int_t iV0=0; iV0<aod ->GetNumberOfV0s(); iV0++) {
+    AliAODv0 *v0 = aod->GetV0(iV0);
+    if( !v0 ) continue;
+    Bool_t lOnFlyStatus = v0->GetOnFlyStatus();
+    if (lOnFlyStatus) continue;
+    
+    AliAODTrack *pTrack=(AliAODTrack *)v0->GetDaughter(0); //0->Positive Daughter
+    AliAODTrack *nTrack=(AliAODTrack *)v0->GetDaughter(1); //1->Negative Daughter
+    if (!pTrack || !nTrack) continue;
+
+    if ( pTrack->Charge() == nTrack->Charge())continue;
+
+      if(!(pTrack->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(!(nTrack->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(pTrack->GetTPCNcls() < 50)continue;
+      if(nTrack->GetTPCNcls() < 50)continue;
+      if(pTrack->Chi2perNDF() > 4)continue;
+      if(nTrack->Chi2perNDF() > 4)continue;
+      
+      Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+      if(!pTrack->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+      if(!nTrack->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+      
+      V0Index[nGoodV0s] = iV0;
+      V0TrackID[2*nGoodV0s] = pTrack->GetID();
+      V0TrackID[2*nGoodV0s+1] = nTrack->GetID();
+      nGoodV0s++; 
+      if(nGoodV0s > 2) break;
+  }//V0s loop
+  
+  //Track loop
+  for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
+    AliAODTrack *trk = aod->GetTrack(itr);
+    if( !trk ) continue;
+
+      if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+      if(trk->GetTPCNcls() < 50)continue;
+      if(trk->Chi2perNDF() > 4)continue;
+      Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+      if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+     
+      TrackID[nGoodTracks] = trk->GetID();
+      nGoodTracks++;
+                                 
+      if(nGoodTracks > 4) break;  
+  }//Track loop
+
+  if(nGoodV0s == 2 && nGoodTracks == 4){
+       //SortArray(TrackID);
+       //SortArray(V0TrackID);
+       //for(Int_t i=0; i<4; i++) if (TrackID[i] != V0TrackID[i]) return;
+       for(Int_t i=0; i<2; i++){
+               AliAODv0 *v0 = aod->GetV0(V0Index[i]);                          
+               fHistK0sMass->Fill(v0->MassK0Short());
+               }
+  }
+
+  
+  PostData(3, fListHist);
+
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::RunAODtree()
+{
+  //input event
+  AliAODEvent *aod = (AliAODEvent*) InputEvent();
+  if(!aod) return;
+
+  //input data
+  const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
+  fDataFilnam->Clear();
+  fDataFilnam->SetString(filnam);
+  fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
+  fRunNum = aod ->GetRunNumber();
+
+  hCounter->Fill( 1 );
+
+  //Trigger
+  TString trigger = aod->GetFiredTriggerClasses();
+  
+  if( !trigger.Contains("CCUP4-B") ) return;
+
+  Bool_t isTRG = kFALSE;
+  for(Int_t i=1; i<ntrg; i++) {
+    if( fTrigger[i] ) {isTRG = kTRUE; hCounter->Fill( fRunNum - 167806 + 1 + i*2000 );}
+  }
+  if( !isTRG ) {PostData(2, hCounter); return;}
+
+  hCounter->Fill( 2 );
+
+  //trigger inputs
+  fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
+  fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
+
+  //Event identification
+  fPerNum = aod ->GetPeriodNumber();
+  fOrbNum = aod ->GetOrbitNumber();
+  fBCrossNum = aod ->GetBunchCrossNumber();
+
+  //primary vertex
+  AliAODVertex *fAODVertex = aod->GetPrimaryVertex();
+  fVtxContrib = fAODVertex->GetNContributors();
+
+  //Tracklets
+  fNtracklets = aod->GetTracklets()->GetNumberOfTracklets();
+
+  //VZERO, ZDC
+  AliAODVZERO *fV0data = aod ->GetVZEROData();
+  AliAODZDC *fZDCdata = aod->GetZDCData();
+  
+  fV0Adecision = fV0data->GetV0ADecision();
+  fV0Cdecision = fV0data->GetV0CDecision();
+  fZDCAenergy = fZDCdata->GetZNATowerEnergy()[0];
+  fZDCCenergy = fZDCdata->GetZNCTowerEnergy()[0];
+  
+  Int_t nGoodV0s=0;
+  Int_t V0Index[3] = {-1,-1,-1};
+  Int_t V0TrackID[6] = {-1,-1,-1,-1,-1,-1};
+
+  Int_t nGoodTracks=0;
+  Int_t TrackID[5] = {-1,-1,-1,-1,-1};
+
+  //V0s loop
+  for(Int_t iV0=0; iV0<aod ->GetNumberOfV0s(); iV0++) {
+    AliAODv0 *v0 = aod->GetV0(iV0);
+    if( !v0 ) continue;
+    Bool_t lOnFlyStatus = v0->GetOnFlyStatus();
+    if (lOnFlyStatus) continue;
+    
+    AliAODTrack *pTrack=(AliAODTrack *)v0->GetDaughter(0); //0->Positive Daughter
+    AliAODTrack *nTrack=(AliAODTrack *)v0->GetDaughter(1); //1->Negative Daughter
+    if (!pTrack || !nTrack) continue;
+
+    if ( pTrack->Charge() == nTrack->Charge())continue;
+
+      if(!(pTrack->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(!(nTrack->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(pTrack->GetTPCNcls() < 50)continue;
+      if(nTrack->GetTPCNcls() < 50)continue;
+      if(pTrack->Chi2perNDF() > 4)continue;
+      if(nTrack->Chi2perNDF() > 4)continue;
+      
+      Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+      if(!pTrack->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+      if(!nTrack->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+      
+      V0Index[nGoodV0s] = iV0;
+      V0TrackID[2*nGoodV0s] = pTrack->GetID();
+      V0TrackID[2*nGoodV0s+1] = nTrack->GetID();
+      nGoodV0s++; 
+      if(nGoodV0s > 2) break;
+  }//V0s loop
+  
+  //Track loop
+  for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
+    AliAODTrack *trk = aod->GetTrack(itr);
+    if( !trk ) continue;
+
+      if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+      if(trk->GetTPCNcls() < 50)continue;
+      if(trk->Chi2perNDF() > 4)continue;
+      Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
+      if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
+      if(TMath::Abs(dca[1]) > 2) continue;
+     
+      TrackID[nGoodTracks] = trk->GetID();
+      nGoodTracks++;
+                                 
+      if(nGoodTracks > 4) break;  
+  }//Track loop
+    
+  if(nGoodV0s == 2 && nGoodTracks == 4){
+       //SortArray(TrackID);
+       //SortArray(V0TrackID);
+       //for{Int_t i=0; i<4; i++} if (TrackID[i] != V0TrackID[i]) {PostData(2, hCounter); return;}
+       for(Int_t i=0; i<2; i++){
+               AliAODv0 *v0 = aod->GetV0(V0Index[i]);
+               AliAODTrack *pTrack=(AliAODTrack *)v0->GetDaughter(0); //0->Positive Daughter
+               AliAODTrack *nTrack=(AliAODTrack *)v0->GetDaughter(1); //1->Negative Daughter
+                               
+               new((*fK0sAODv0s)[i]) AliAODv0(*v0);
+               new((*fK0sAODTracks)[2*i]) AliAODTrack(*pTrack);
+               new((*fK0sAODTracks)[2*i+1]) AliAODTrack(*nTrack);
+               }
+  fK0sTree ->Fill();
+  
+  }
+  
+  PostData(1, fK0sTree);   
+  PostData(2, hCounter);
+
+}//RunAOD
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::SortArray(Double_t *dArray) {
+    for (Int_t i = 3; i > 0; --i) {
+        for (Int_t j = 0; j < i; ++j) {
+            if (dArray[j] > dArray[j+1]) {
+                Double_t dTemp = dArray[j];
+                dArray[j] = dArray[j+1];
+                dArray[j+1] = dTemp;
+            }
+        }
+    }
+}
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::RunESD()
+{
+
+  /*/input event
+  AliESDEvent *esd = (AliESDEvent*) InputEvent();
+  if(!esd) return;
+
+  //input data
+  const char *filnam = ((TTree*) GetInputData(0))->GetCurrentFile()->GetName();
+  fDataFilnam->Clear();
+  fDataFilnam->SetString(filnam);
+  fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
+  fRunNum = esd->GetRunNumber();
+
+  hCounter->Fill( 1 );
+
+  //Trigger
+  TString trigger = esd->GetFiredTriggerClasses();
+  
+  fTrigger[0]   = trigger.Contains("CINT7-B");
+  fTrigger[1]   = trigger.Contains("CCUP4-B"); // CE 
+  fTrigger[2]   = trigger.Contains("CCUP4-E"); // CE 
+
+  Bool_t isTRG = kFALSE;
+  for(Int_t i=1; i<ntrg; i++) {
+    if( fTrigger[i] ) {isTRG = kTRUE; hCounter->Fill( fRunNum - 167806 + 1 + i*2000 );}
+  }
+  if( !isTRG ) {PostData(3, hCounter); return;}
+
+  hCounter->Fill( 2 );
+
+  //trigger inputs
+  fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
+  fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
+
+  //Event identification
+  fPerNum = esd->GetPeriodNumber();
+  fOrbNum = esd->GetOrbitNumber();
+  fBCrossNum = esd->GetBunchCrossNumber();
+
+  //primary vertex
+  AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
+  fVtxContrib = fESDVertex->GetNContributors();
+
+  //Tracklets
+  fNtracklets = esd->GetMultiplicity()->GetNumberOfTracklets();
+
+  //VZERO, ZDC
+  AliESDVZERO *fV0data = esd->GetVZEROData();
+  AliESDZDC *fZDCdata = esd->GetESDZDC();
+  
+  fV0Adecision = fV0data->GetV0ADecision();
+  fV0Cdecision = fV0data->GetV0CDecision();
+  fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
+  fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
+  
+  Int_t nGoodTracks=0;
+  Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
+  
+  //Track loop
+  for(Int_t itr=0; itr<esd ->GetNumberOfTracks(); itr++) {
+    AliESDtrack *trk = esd->GetTrack(itr);
+    if( !trk ) continue;
+
+      if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
+      if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
+      if(trk->GetTPCNcls() < 50)continue;
+      if(trk->GetTPCchi2()/trk->GetTPCNcls() > 4)continue;
+      Float_t dca[2] = {0.0,0.0}; AliExternalTrackParam cParam;
+      if(!trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam)) continue;
+      trk->GetImpactParameters(dca[0],dca[1]);
+      if(TMath::Abs(dca[1]) > 2) continue;
+      
+      TrackIndex[nGoodTracks] = itr;
+      nGoodTracks++;
+      if(nGoodTracks > 4) break;   
+  }//Track loop
+
+  if(nGoodTracks == 2){
+         for(Int_t i=0; i<2; i++){
+               AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
+               
+               AliExternalTrackParam cParam;
+               trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
+                               
+               new((*fJPsiESDTracks)[i]) AliESDtrack(*trk); 
+               
+               }
+  fJPsiTree ->Fill();
+  PostData(1, fJPsiTree);
+  }
+  
+  if(nGoodTracks == 4){
+         for(Int_t i=0; i<4; i++){
+               AliESDtrack *trk = esd->GetTrack(TrackIndex[i]);
+               
+               AliExternalTrackParam cParam;
+               trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
+
+               new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk); 
+               
+               }
+  fPsi2sTree ->Fill();
+  PostData(2, fPsi2sTree);
+  }
+    
+  PostData(3, hCounter);
+/*/
+}//RunESD
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcK0sK0s::Terminate(Option_t *) 
+{
+
+  cout<<"Analysis complete."<<endl;
+}//Terminate
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.h b/PWGUD/UPC/AliAnalysisTaskUpcK0sK0s.h
new file mode 100644 (file)
index 0000000..55db3ad
--- /dev/null
@@ -0,0 +1,86 @@
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+* See cxx source for full Copyright notice */
+/* $Id$ */
+
+#ifndef ALIANALYSISTASKUPCK0SK0S_H
+#define ALIANALYSISTASKUPCK0SK0S_H
+
+class TClonesArray;
+class TTree;
+class TH1;
+class TH2;
+class TList;
+
+#define ntrg 17
+#include "AliAnalysisTaskSE.h"
+
+class AliAnalysisTaskUpcK0sK0s : public AliAnalysisTaskSE {
+ public:
+  AliAnalysisTaskUpcK0sK0s();
+  AliAnalysisTaskUpcK0sK0s(const char *name);
+  virtual ~AliAnalysisTaskUpcK0sK0s();
+
+  virtual void Init();
+  virtual void UserCreateOutputObjects();
+  virtual void UserExec(Option_t *option);
+  virtual void RunAODhist();
+  virtual void RunAODtree();
+  virtual void RunESD();
+  virtual void Terminate(Option_t *);
+  void SetRunTree(Bool_t runTree){fRunTree = runTree;}
+  void SetRunHist(Bool_t runHist){fRunHist = runHist;}
+  void SortArray(Double_t *dArray);
+
+ private:
+  Int_t fType; // 0 - ESD, 1 - AOD
+  Bool_t fRunTree; 
+  Bool_t fRunHist;
+  
+  //counter
+  TH1I *hCounter; //!
+
+  //event tree
+  TTree *fK0sTree;
+  //tree variables
+  Int_t fRunNum;
+  UInt_t fPerNum, fOrbNum;
+  //trigger
+  Bool_t fTrigger[ntrg];
+  UInt_t fL0inputs, fL1inputs;
+  Int_t fVtxContrib;
+  UShort_t fBCrossNum, fNtracklets;
+  //vzero, zdc
+  Double_t fZDCAenergy, fZDCCenergy;
+  Int_t fV0Adecision, fV0Cdecision;
+  //input data
+  TObjString *fDataFilnam;
+  Short_t fRecoPass;
+  Long64_t fEvtNum;
+  //vertices
+  TClonesArray *fK0sAODv0s;
+  TClonesArray *fK0sESDv0s;
+  //tracks
+  TClonesArray *fK0sAODTracks;
+  TClonesArray *fK0sESDTracks; 
+  
+  TList *fListHist;
+  TH1D *fHistTriggersPerRun;
+  TH1D *fHistK0sMass;
+  
+  AliAnalysisTaskUpcK0sK0s(const AliAnalysisTaskUpcK0sK0s&); //not implemented
+  AliAnalysisTaskUpcK0sK0s& operator =(const AliAnalysisTaskUpcK0sK0s&); //not implemented
+  
+  ClassDef(AliAnalysisTaskUpcK0sK0s, 1); 
+};
+
+#endif
+
+
+
+
+
+
+
+
+
index 998a969..3429c67 100644 (file)
@@ -117,6 +117,7 @@ void AliAnalysisTaskUpcPsi2s::Init()
 {
   
   for(Int_t i=0; i<ntrg; i++) fTrigger[i] = kFALSE;
+  for(Int_t i=0; i<4; i++) fTOFphi[i] = -666;
 
 }//Init
 
@@ -172,6 +173,7 @@ void AliAnalysisTaskUpcPsi2s::UserCreateOutputObjects()
   
   fJPsiTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
   fJPsiTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
+  fJPsiTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
   
   fJPsiTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
   fJPsiTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
@@ -211,6 +213,7 @@ void AliAnalysisTaskUpcPsi2s::UserCreateOutputObjects()
   
   fPsi2sTree ->Branch("fTOFtrig1", &fTOFtrig1, "fTOFtrig1/O");
   fPsi2sTree ->Branch("fTOFtrig2", &fTOFtrig2, "fTOFtrig2/O");
+  fPsi2sTree ->Branch("fTOFphi", &fTOFphi[0], "fTOFphi[4]/D");
   
   fPsi2sTree ->Branch("fVtxPosX", &fVtxPosX, "fVtxPosX/D");
   fPsi2sTree ->Branch("fVtxPosY", &fVtxPosY, "fVtxPosY/D");
@@ -382,7 +385,7 @@ void AliAnalysisTaskUpcPsi2s::RunAODhist()
   //Trigger
   TString trigger = aod->GetFiredTriggerClasses();
   
-  if( !trigger.Contains("CCUP4-B") ) return;
+  if( !trigger.Contains("CCUP") ) return;
   
   fHistNeventsJPsi->Fill(2);
   fHistNeventsPsi2s->Fill(2);
@@ -585,14 +588,19 @@ void AliAnalysisTaskUpcPsi2s::RunAODtree()
   //Trigger
   TString trigger = aod->GetFiredTriggerClasses();
   
-  if( !trigger.Contains("CCUP4-B") ) return;
+  fTrigger[0]  = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
+  fTrigger[1]  = trigger.Contains("CCUP2-B"); // Double gap
+  fTrigger[2]  = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
+  
+  Bool_t isTriggered = kFALSE;
+  for(Int_t i=0; i<ntrg; i++) {
+    if( fTrigger[i] ) isTriggered = kTRUE;
+  }
+  if( !isTriggered ) return;
 
   //trigger inputs
   fL0inputs = aod->GetHeader()->GetL0TriggerInputs();
-  fL1inputs = aod->GetHeader()->GetL1TriggerInputs();
-  
-  //TOF trigger info (0OMU)
-  
+  fL1inputs = aod->GetHeader()->GetL1TriggerInputs();  
 
   //Event identification
   fPerNum = aod ->GetPeriodNumber();
@@ -628,85 +636,95 @@ void AliAnalysisTaskUpcPsi2s::RunAODtree()
   Int_t nGoodTracks=0;
   Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
   
-   //Four track loop
+   //Two track loop
   for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
     AliAODTrack *trk = aod->GetTrack(itr);
     if( !trk ) continue;
 
       if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
       if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
-      if(trk->GetTPCNcls() < 50)continue;
+      if(trk->GetTPCNcls() < 70)continue;
       if(trk->Chi2perNDF() > 4)continue;
+      if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
       Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
       AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
       if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
       delete trk_clone;
-      if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
       if(TMath::Abs(dca[1]) > 2) continue;
-
+      if(TMath::Abs(dca[0]) > 0.2) continue;
+     
       TrackIndex[nGoodTracks] = itr;
       nGoodTracks++;
                                  
-      if(nGoodTracks > 4) break;  
+      if(nGoodTracks > 2) break;  
   }//Track loop
-  
-  
-  if(nGoodTracks == 4){
-         for(Int_t i=0; i<4; i++){
+
+  if(nGoodTracks == 2){
+         for(Int_t i=0; i<2; i++){
                AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
                
                Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
                AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
                if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
                delete trk_clone;
+                               
+               new((*fJPsiAODTracks)[i]) AliAODTrack(*trk); 
+               ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
                
-               new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
-               ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
-                               
+               Double_t pos[3]={0,0,0};
+               if(!trk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
+               else {
+                    fTOFphi[i] =  TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+                    if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+                    }          
                }
-  fPsi2sTree ->Fill();
+  fJPsiTree ->Fill();
   }
- //
-  nGoodTracks = 0;
-  //Two track loop
+  
+   nGoodTracks = 0;
+   //Four track loop
   for(Int_t itr=0; itr<aod ->GetNumberOfTracks(); itr++) {
     AliAODTrack *trk = aod->GetTrack(itr);
     if( !trk ) continue;
 
       if(!(trk->GetStatus() & AliESDtrack::kTPCrefit) ) continue;
       if(!(trk->GetStatus() & AliESDtrack::kITSrefit) ) continue;
-      if(trk->GetTPCNcls() < 70)continue;
+      if(trk->GetTPCNcls() < 50)continue;
       if(trk->Chi2perNDF() > 4)continue;
-      if((!trk->HasPointOnITSLayer(0))&&(!trk->HasPointOnITSLayer(1))) continue;
       Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
       AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
       if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
       delete trk_clone;
+      if(!trk->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
       if(TMath::Abs(dca[1]) > 2) continue;
-      if(TMath::Abs(dca[0]) > 0.2) continue;
-     
+
       TrackIndex[nGoodTracks] = itr;
       nGoodTracks++;
                                  
-      if(nGoodTracks > 2) break;  
+      if(nGoodTracks > 4) break;  
   }//Track loop
-
-  if(nGoodTracks == 2){
-         for(Int_t i=0; i<2; i++){
+    
+  if(nGoodTracks == 4){
+         for(Int_t i=0; i<4; i++){
                AliAODTrack *trk = aod->GetTrack(TrackIndex[i]);
                
                Double_t dca[2] = {0.0,0.0}, cov[3] = {0.0,0.0,0.0};
                AliAODTrack* trk_clone=(AliAODTrack*)trk->Clone("trk_clone");
                if(!trk_clone->PropagateToDCA(fAODVertex,aod->GetMagneticField(),300.,dca,cov)) continue;
                delete trk_clone;
-                               
-               new((*fJPsiAODTracks)[i]) AliAODTrack(*trk); 
-               ((AliAODTrack*)((*fJPsiAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
                
+               new((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
+               ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
+                               
+               Double_t pos[3]={0,0,0};
+               if(!trk->GetXYZAt(378,aod->GetMagneticField(),pos)) fTOFphi[i] = -666;
+               else {
+                    fTOFphi[i] =  TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+                    if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+                    }          
                }
-  fJPsiTree ->Fill();
+  fPsi2sTree ->Fill();
   }
-
   
   PostData(1, fJPsiTree);
   PostData(2, fPsi2sTree);
@@ -774,7 +792,7 @@ void AliAnalysisTaskUpcPsi2s::RunESDhist()
   //Trigger
   TString trigger = esd->GetFiredTriggerClasses();
   
-  if( !trigger.Contains("CCUP4-B") ) return;
+  if( !trigger.Contains("CCUP") ) return;
   
   fHistNeventsJPsi->Fill(2);
   fHistNeventsPsi2s->Fill(2);
@@ -963,11 +981,19 @@ void AliAnalysisTaskUpcPsi2s::RunESDtree()
   fEvtNum = ((TTree*) GetInputData(0))->GetTree()->GetReadEntry();
   fRunNum = esd->GetRunNumber();
 
-  //Trigger
+   //Trigger
   TString trigger = esd->GetFiredTriggerClasses();
   
-  if( !trigger.Contains("CCUP4-B") ) return;
-
+  fTrigger[0]  = trigger.Contains("CCUP4-B"); // Central UPC Pb-Pb 2011
+  fTrigger[1]  = trigger.Contains("CCUP2-B"); // Double gap
+  fTrigger[2]  = trigger.Contains("CCUP7-B"); // Central UPC p-Pb 2013
+  
+  Bool_t isTriggered = kFALSE;
+  for(Int_t i=0; i<ntrg; i++) {
+    if( fTrigger[i] ) isTriggered = kTRUE;
+  }
+  if( !isTriggered ) return;
+  
   //trigger inputs
   fL0inputs = esd->GetHeader()->GetL0TriggerInputs();
   fL1inputs = esd->GetHeader()->GetL1TriggerInputs();
@@ -1040,6 +1066,12 @@ void AliAnalysisTaskUpcPsi2s::RunESDtree()
                                
                new((*fJPsiESDTracks)[i]) AliESDtrack(*trk); 
                
+               Double_t pos[3]={0,0,0};
+               if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
+               else {
+                    fTOFphi[i] =  TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+                    if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+                    }          
                }
   fJPsiTree ->Fill();
   }
@@ -1071,8 +1103,14 @@ void AliAnalysisTaskUpcPsi2s::RunESDtree()
                AliExternalTrackParam cParam;
                trk->RelateToVertex(fESDVertex, esd->GetMagneticField(),300.,&cParam);// to get trk->GetImpactParameters(DCAxy,DCAz);
 
-               new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk); 
-               
+               new((*fPsi2sESDTracks)[i]) AliESDtrack(*trk);
+                
+               Double_t pos[3]={0,0,0};
+               if(!trk->GetXYZAt(378,esd->GetMagneticField(),pos)) fTOFphi[i] = -666;
+               else {
+                    fTOFphi[i] =  TMath::ATan2(pos[1],pos[0])*TMath::RadToDeg();
+                    if(fTOFphi[i] < 0) fTOFphi[i]+=(2*TMath::Pi()*TMath::RadToDeg());
+                    }          
                }
   fPsi2sTree ->Fill();
   }
@@ -1089,32 +1127,3 @@ void AliAnalysisTaskUpcPsi2s::Terminate(Option_t *)
   cout<<"Analysis complete."<<endl;
 }//Terminate
 
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
index ae87cf4..8c6762f 100644 (file)
@@ -48,6 +48,7 @@ class AliAnalysisTaskUpcPsi2s : public AliAnalysisTaskSE {
   Bool_t fTrigger[ntrg];
   UInt_t fL0inputs, fL1inputs;
   Bool_t fTOFtrig1, fTOFtrig2;
+  Double_t fTOFphi[4];
   Int_t fVtxContrib;
   Double_t fVtxPosX,fVtxPosY,fVtxPosZ;
   Double_t fVtxErrX,fVtxErrY,fVtxErrZ;