+
+ fPsi2sAODTracks->Clear("C");
+ if(nGoodTracks == 4){
+
+ TDatabasePDG *pdgdat = TDatabasePDG::Instance();
+ TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
+ Double_t muonMass = partMuon->Mass();
+ TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
+ Double_t electronMass = partElectron->Mass();
+ TParticlePDG *partPion = pdgdat->GetParticle( 211 );
+ Double_t pionMass = partPion->Mass();
+
+ Double_t KFcov[21];
+ Double_t KFpar[6];
+ Double_t KFmass = pionMass;
+ Double_t fRecTPCsignal;
+ AliKFParticle *KFpart[4];
+ AliKFVertex *KFvtx = new AliKFVertex();
+ KFvtx->SetField(aod->GetMagneticField());
+
+ 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((*fPsi2sAODTracks)[i]) AliAODTrack(*trk);
+ ((AliAODTrack*)((*fPsi2sAODTracks)[i]))->SetDCA(dca[0],dca[1]);//to get DCAxy trk->DCA(); to get DCAz trk->ZAtDCA();
+
+
+ fPIDTPCMuon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kMuon);
+ fPIDTPCElectron[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kElectron);
+ fPIDTPCPion[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kPion);
+ fPIDTPCKaon[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kKaon);
+ fPIDTPCProton[i] = fPIDResponse->NumberOfSigmasTPC(trk,AliPID::kProton);
+
+ fPIDTOFMuon[i] = fPIDResponse->NumberOfSigmasTOF(trk,AliPID::kMuon);
+ fPIDTOFElectron[i] = fPIDResponse->NumberOfSigmasTOF(trk,AliPID::kElectron);
+ fPIDTOFPion[i] = fPIDResponse->NumberOfSigmasTOF(trk,AliPID::kPion);
+ fPIDTOFKaon[i] = fPIDResponse->NumberOfSigmasTOF(trk,AliPID::kKaon);
+ fPIDTOFProton[i] = fPIDResponse->NumberOfSigmasTOF(trk,AliPID::kProton);
+
+ trk->GetPosition(KFpar);
+ trk->PxPyPz(KFpar+3);
+ trk->GetCovMatrix(KFcov);
+
+ if(trk->Pt() > 1){
+ fRecTPCsignal = trk->GetTPCsignal();
+ if(fRecTPCsignal > 40 && fRecTPCsignal < 70) KFmass = muonMass;
+ if(fRecTPCsignal > 70 && fRecTPCsignal < 100)KFmass = electronMass;
+ }
+ else KFmass = pionMass;
+
+ KFpart[i] = new AliKFParticle();
+ KFpart[i]->SetField(aod->GetMagneticField());
+ KFpart[i]->AliKFParticleBase::Initialize(KFpar,KFcov,(Int_t) trk->Charge(), KFmass);
+ KFvtx->AddDaughter(*KFpart[i]);
+
+ Double_t pos[3]={0,0,0};
+ AliExternalTrackParam *parTrk = new AliExternalTrackParam();
+ parTrk->CopyFromVTrack((AliVTrack*) trk);
+ if(!parTrk->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());
+ }
+ delete parTrk;
+ }
+ fKfVtxPos[0]= KFvtx->GetX();
+ fKfVtxPos[1]= KFvtx->GetY();
+ fKfVtxPos[2]= KFvtx->GetZ();
+ for(UInt_t i=0; i<4; i++)delete KFpart[i];
+ delete KFvtx;
+ if(!isMC) fPsi2sTree ->Fill();
+ }
+
+ if(isMC){
+ fJPsiTree ->Fill();
+ fPsi2sTree ->Fill();
+ }
+
+ PostData(1, fJPsiTree);
+ PostData(2, fPsi2sTree);
+
+}//RunAOD
+
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcPsi2s::RunAODMC(AliAODEvent *aod)
+{
+
+ fGenPart->Clear("C");
+
+ TClonesArray *arrayMC = (TClonesArray*) aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
+ if(!arrayMC) return;
+
+ Int_t nmc=0;
+ //loop over mc particles
+ for(Int_t imc=0; imc<arrayMC->GetEntriesFast(); imc++) {
+ AliAODMCParticle *mcPart = (AliAODMCParticle*) arrayMC->At(imc);
+ if(!mcPart) continue;
+
+ if(mcPart->GetMother() >= 0) continue;
+
+ TParticle *part = (TParticle*) fGenPart->ConstructedAt(nmc++);
+ part->SetMomentum(mcPart->Px(), mcPart->Py(), mcPart->Pz(), mcPart->E());
+ part->SetPdgCode(mcPart->GetPdgCode());
+ part->SetUniqueID(imc);
+ }//loop over mc particles
+
+}//RunAODMC
+
+
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcPsi2s::RunESDtrig()
+{
+
+ //input event
+ AliESDEvent *esd = (AliESDEvent*) InputEvent();
+ if(!esd) return;
+
+ fRunNum = esd ->GetRunNumber();
+ //Trigger
+ TString trigger = esd->GetFiredTriggerClasses();
+
+ if(trigger.Contains("CCUP4-B")) fHistCcup4TriggersPerRun->Fill(fRunNum); //CCUP4 triggers
+ if(trigger.Contains("CCUP7-B")) fHistCcup7TriggersPerRun->Fill(fRunNum); //CCUP7 triggers
+ if(trigger.Contains("CCUP2-B")) fHistCcup2TriggersPerRun->Fill(fRunNum); //CCUP2 triggers
+
+ if(trigger.Contains("CVLN_B2-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - synchronously downscaled
+ if(trigger.Contains("CVLN_R1-B")) fHistCvlnTriggersPerRun->Fill(fRunNum); //CVLN triggers - randomly downscaled
+
+ if(esd->GetHeader()->IsTriggerInputFired("1ZED")) fHistZedTriggersPerRun->Fill(fRunNum); //1ZED trigger inputs
+
+ //MB, Central and SemiCentral triggers
+ AliCentrality *centrality = esd->GetCentrality();
+ UInt_t selectionMask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
+
+ //Double_t percentile = centrality->GetCentralityPercentile("V0M");
+ Double_t percentile = centrality->GetCentralityPercentileUnchecked("V0M");
+
+ if(((selectionMask & AliVEvent::kMB) == AliVEvent::kMB) && percentile<=80 && percentile>=0) fHistMBTriggersPerRun->Fill(fRunNum);
+
+ if(((selectionMask & AliVEvent::kCentral) == AliVEvent::kCentral) && percentile<=6 && percentile>=0 && (trigger.Contains("CVHN_R2-B"))) fHistCentralTriggersPerRun->Fill(fRunNum);
+
+ if(((selectionMask & AliVEvent::kSemiCentral) == AliVEvent::kSemiCentral) && percentile<=50 && percentile>=15) fHistSemiCentralTriggersPerRun->Fill(fRunNum);
+
+
+PostData(3, fListTrig);
+
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskUpcPsi2s::RunESDhist()
+{
+
+ TDatabasePDG *pdgdat = TDatabasePDG::Instance();
+
+ TParticlePDG *partMuon = pdgdat->GetParticle( 13 );
+ Double_t muonMass = partMuon->Mass();
+
+ TParticlePDG *partElectron = pdgdat->GetParticle( 11 );
+ Double_t electronMass = partElectron->Mass();
+
+ TParticlePDG *partPion = pdgdat->GetParticle( 211 );
+ Double_t pionMass = partPion->Mass();
+
+ //input event
+ AliESDEvent *esd = (AliESDEvent*) InputEvent();
+ if(!esd) return;
+
+ fHistNeventsJPsi->Fill(1);
+ fHistNeventsPsi2s->Fill(1);
+
+ //Trigger
+ TString trigger = esd->GetFiredTriggerClasses();
+
+ if(!isMC && !trigger.Contains("CCUP") ) return;
+
+ fHistNeventsJPsi->Fill(2);
+ fHistNeventsPsi2s->Fill(2);
+
+ //primary vertex
+ AliESDVertex *fESDVertex = (AliESDVertex*) esd->GetPrimaryVertex();
+ fVtxContrib = fESDVertex->GetNContributors();
+ if(fVtxContrib < 2) return;
+
+ fHistNeventsJPsi->Fill(3);
+ fHistNeventsPsi2s->Fill(3);
+
+ //VZERO, ZDC
+ AliESDVZERO *fV0data = esd->GetVZEROData();
+ AliESDZDC *fZDCdata = esd->GetESDZDC();
+
+ fV0Adecision = fV0data->GetV0ADecision();
+ fV0Cdecision = fV0data->GetV0CDecision();
+ if(fV0Adecision != AliESDVZERO::kV0Empty || fV0Cdecision != AliESDVZERO::kV0Empty) return;
+
+ fZDCAenergy = fZDCdata->GetZN2TowerEnergy()[0];
+ fZDCCenergy = fZDCdata->GetZN1TowerEnergy()[0];
+ if( fZDCAenergy > 8200 || fZDCCenergy > 8200) return;
+
+ fHistNeventsJPsi->Fill(4);
+ fHistNeventsPsi2s->Fill(4);
+
+ Int_t nGoodTracks=0;
+ //Two tracks loop
+ Int_t TrackIndex[5] = {-1,-1,-1,-1,-1};
+
+ TLorentzVector vLepton[4], vPion[4], vCandidate, vDilepton;
+ Short_t qLepton[4], qPion[4];
+ UInt_t nLepton=0, nPion=0, nHighPt=0;
+ Double_t fRecTPCsignal[5];