//==============================================================
#include "TH2.h"
#include "TH1.h"
+#include "TClonesArray.h"
#include "TFile.h"
+#include "TTree.h"
#include "TRandom.h"
#include "TObjectTable.h"
-#include "AliPHOSIndexToObject.h"
+#include "AliRun.h"
+#include "AliPHOSGetter.h"
#include "AliPHOSRecParticle.h"
#include "TLorentzVector.h"
#include "TGraphErrors.h"
#include "TF1.h"
TObjectTable * gObjectTable;
-TRandom * gRandom;
+TRandom * gRandom;
+AliRun * gAlice;
-void AnaMinv(char * filename)
+Bool_t PhotonId(AliPHOSRecParticle * RecParticle, char *opt)
+{
+ Bool_t PhotonId = kTRUE ;
+ if (strstr(opt,"narrow"))
+ {
+ if ( RecParticle->GetType() & 1) PhotonId= kFALSE ;
+ }
+ if (strstr(opt,"cpv"))
+ {
+ if ( (RecParticle->GetType() & 1<<2) >>2 ) PhotonId= kFALSE ;
+ }
+ if (strstr(opt,"time"))
+ {
+
+ }
+
+ return PhotonId;
+}
+
+void AnaMinv(char * filename,char * opt)
{
TH2F * h_Minv_lowpT = new TH2F("h_Minv_lowpT","Minv vs pT low",500,0.0,1.0,40,0.,10.);
TH2F * h_Minv_highpT = new TH2F("h_Minv_highpT","Minv vs pT high",500,0.0,1.0,50,0.,100.);
TH2F * h_Peta_Pt = new TH2F("h_Peta_Pt","Pseudo vs pT",40,0.,10.,50,-1.0,1.0);
TH1F * h_Phi = new TH1F("h_Phi","Phi photons",400,-4.,4.);
TH2F * h_Peta_Phi = new TH2F("h_Peta_Phi","Pseudo vs Phi",200,-4,4,200,-1.0,1.0);
+ TH1F * h_Dispersion= new TH1F("h_Dispersion","Dispersion",50,0.,10.);
+ TH1F * h_Type = new TH1F("h_Type","Particle Type",9,-1.5,7.5);
- TH1F * h_DeltaR = new TH1F("h_DeltaR","Delta R",400,0,4);
+ TH1F * h_DeltaR = new TH1F("h_DeltaR","Delta R",400,0.,2.);
TH1F * h_Asymmetry= new TH1F("h_Asymmetry","Asymmetry",400, -2., 2.);
- AliPHOSIndexToObject * RecData = AliPHOSIndexToObject::GetInstance(filename) ;
+ AliPHOSGetter * RecData = AliPHOSGetter::GetInstance(filename,"Gines") ;
AliPHOSRecParticle * RecParticle1;
AliPHOSRecParticle * RecParticle2;
-
- Float_t RelativeRCut = 0.0001 ;
+ Float_t RelativeRCut = 0.00001 ;
Float_t AsymmetryCut = 0.7 ;
Float_t Asymmetry;
+ Float_t Type;
Int_t iEvent, iRecParticle1, iRecParticle2;
Int_t nRecParticle;
Float_t invariant_mass, invariant_mass_mixed;
- TLorentzVector P_photon1, P_photon2, P_photonMixed1, P_photonMixed2 ;
+ TClonesArray * RecParticles;
+
Float_t average_multiplicity = 0.;
- for(iEvent=0; iEvent<RecData->GetMaxEvent(); iEvent++)
- // for(iEvent=0; iEvent<1000; iEvent++)
+ for(iEvent=0; iEvent<gAlice->TreeE()->GetEntries(); iEvent++)
{
- // if (iEvent==2) gObjectTable->Print();
- //if (iEvent==15) gObjectTable->Print();
- RecData->GetEvent(iEvent);
+ TLorentzVector P_photon1, P_photon2, P_photonMixed1, P_photonMixed2 ;
+
+ RecData->Event(iEvent,"R");
printf(">>> Event %d \n",iEvent);
- nRecParticle=RecData->GimeNRecParticles();
- average_multiplicity += ((Float_t) (nRecParticle) ) / ( (Float_t)RecData->GetMaxEvent() ) ;
+ RecParticles = RecData->RecParticles();
+ nRecParticle = RecParticles->GetEntries();
+ printf(">>> >> Recparticle multilicity is %d \n",nRecParticle);
+ average_multiplicity += ((Float_t) (nRecParticle) ) / ( (Float_t)gAlice->TreeE()->GetEntries() ) ;
// Construction de la masse invariante des pairs
if (nRecParticle > 1)
{
for(iRecParticle1=0; iRecParticle1<nRecParticle; iRecParticle1++)
{
- RecParticle1 = (AliPHOSRecParticle *) RecData->GimeRecParticle(iRecParticle1);
+ RecParticle1 = (AliPHOSRecParticle *) RecParticles->At(iRecParticle1);
RecParticle1->Momentum(P_photon1);
-
- h_Pseudoeta->Fill(P_photon1.PseudoRapidity());
- h_Pt->Fill(P_photon1.Pt());
- h_Phi->Fill(P_photon1.Phi());
- h_Peta_Pt->Fill(P_photon1.Pt(), P_photon1.PseudoRapidity());
- h_Peta_Phi->Fill(P_photon1.Phi(), P_photon1.PseudoRapidity() );
-
- for(iRecParticle2=iRecParticle1+1; iRecParticle2<nRecParticle; iRecParticle2++)
- {
- RecParticle2 = (AliPHOSRecParticle *) RecData->GimeRecParticle(iRecParticle2);
- RecParticle2->Momentum(P_photon2);
- Asymmetry = TMath::Abs((P_photon1.E()-P_photon2.E())/(P_photon1.E()+P_photon2.E()));
- if ( (P_photon1 != P_photon2) &&
- (P_photon1.DeltaR(P_photon2) > RelativeRCut) &&
- (Asymmetry < AsymmetryCut) )
- {
- h_DeltaR->Fill(P_photon1.DeltaR(P_photon2));
- h_Asymmetry->Fill( Asymmetry );
-
- // printf("A. p1 es %f \n",P_photon1->E());
- invariant_mass = (P_photon1 + P_photon2).M();
- // printf("B. p1 es %f \n",P_photon1->E());
- h_Minv_lowpT->Fill(invariant_mass, (P_photon1 + P_photon2).Pt() );
- h_Minv_highpT->Fill(invariant_mass,(P_photon1 + P_photon2).Pt() );
- }
+ printf(">>> >> Photon momentum is %f \n", P_photon1.Pt() );
+ Type = RecParticle1->GetType();
+ h_Type->Fill(Type);
+ // Photon Id check
+ if ( PhotonId(RecParticle1, opt) )
+ {
+
+ h_Pseudoeta->Fill(P_photon1.PseudoRapidity());
+ h_Pt->Fill(P_photon1.Pt());
+ h_Phi->Fill(P_photon1.Phi());
+ h_Peta_Pt->Fill(P_photon1.Pt(), P_photon1.PseudoRapidity());
+ h_Peta_Phi->Fill(P_photon1.Phi(), P_photon1.PseudoRapidity() );
+
+ for(iRecParticle2=iRecParticle1+1; iRecParticle2<nRecParticle; iRecParticle2++)
+ {
+ RecParticle2 = (AliPHOSRecParticle *) RecParticles->At(iRecParticle2);
+ RecParticle2->Momentum(P_photon2);
+ Asymmetry = TMath::Abs((P_photon1.E()-P_photon2.E())/(P_photon1.E()+P_photon2.E()));
+ //Photon Id Check
+ if ( PhotonId(RecParticle2, opt) )
+ {
+ if ( (P_photon1 != P_photon2) &&
+ (P_photon1.DeltaR(P_photon2) > RelativeRCut) &&
+ (Asymmetry < AsymmetryCut) )
+ {
+ h_DeltaR->Fill(P_photon1.DeltaR(P_photon2));
+ h_Asymmetry->Fill( Asymmetry );
+
+ // printf("A. p1 es %f \n",P_photon1->E());
+ invariant_mass = (P_photon1 + P_photon2).M();
+ // printf("B. p1 es %f \n",P_photon1->E());
+ h_Minv_lowpT->Fill(invariant_mass, (P_photon1 + P_photon2).Pt() );
+ h_Minv_highpT->Fill(invariant_mass,(P_photon1 + P_photon2).Pt() );
+ }
+ }
+ }
}
}
}
}
printf(">>> Average Multiplicity is %f \n",average_multiplicity);
- Int_t Background = (Int_t) (RecData->GetMaxEvent() * average_multiplicity * (average_multiplicity-1.)/2.) ;
+ Int_t Background = (Int_t) (gAlice->TreeE()->GetEntries() * average_multiplicity * (average_multiplicity-1.)/2.) ;
printf(">>> Background is %d \n",Background);
Double_t Pt_Mixed1, Pt_Mixed2;
for(iEvent=0; iEvent<Background; iEvent++)
{
- // printf(">>> Background Event %d \n",iEvent);
+ TLorentzVector P_photon1, P_photon2, P_photonMixed1, P_photonMixed2 ;
+ //printf(">>> Background Event %d \n",iEvent);
Pt_Mixed1 = h_Pt->GetRandom();
Pt_Mixed2 = h_Pt->GetRandom();
h_Peta_Phi->GetRandom2(Phi_Mixed1, Y_Mixed1);
h_Minv_highpT_back->Fill(invariant_mass_mixed,(P_photonMixed1 + P_photonMixed2).Pt() );
}
}
-
+ printf (">>> Background calculation finished ! \n");
char outputname[80];
- sprintf(outputname,"%s.Minv",filename);
+ sprintf(outputname,"%s.Minv%s",filename,opt);
TFile output(outputname,"recreate");
h_Minv_lowpT->Write();
h_Minv_highpT->Write();
h_Peta_Pt->Write();
h_Phi->Write();
h_Peta_Phi->Write();
+ h_Dispersion->Write();
+ h_Type->Write();
h_Asymmetry->Write();
h_DeltaR->Write();
Int_t Norma_1 = 100; Float_t Norma_minv_1 = 0.2;
Int_t Norma_2 = 200; Float_t Norma_minv_2 = 0.4;
- Int_t Minv_1 = 56;
+ Int_t Minv_1 = 46;
Int_t Minv_2 = 76;
if (strstr(particle,"eta"))
{
- Minv_1 = 234;
+ Minv_1 = 214;
Minv_2 = 314;
}
{
signal = h_Minv_pT->ProjectionX("signal", NumberOfPtBins*iHistos+1,NumberOfPtBins*(iHistos+1));
background = h_Minv_pT_back->ProjectionX("background",NumberOfPtBins*iHistos+1,NumberOfPtBins*(iHistos+1));
- //signal->Rebin();
- //background->Rebin();
+
ratio = new TH1D(*signal);
ratio->Sumw2();
ratio->Add(background,-1.0);
else
Renorma = signal->Integral(Norma_1,Norma_2)/background->Integral(Norma_1,Norma_2);
difference->Add(background,(-1.)*Renorma);
-
- //ratio->Draw();
- // background->Draw("same");
- // difference->Draw();
+
Ntota = signal->Integral(Minv_1,Minv_2);
Nback = background->Integral(Minv_1,Minv_2);
NmesonsError[iHistos] = TMath::Sqrt( Ntota + Nback*Renorma*Renorma + Nback*Nback*NormaError*NormaError );
Pt[iHistos] = (iHistos+0.5)*NumberOfPtBins*PtCalibration;
PtError[iHistos] = NumberOfPtBins*PtCalibration/2.;
- // ratio->Delete("");
- //difference->Delete("");
+
}
- // in->Close();
-
char filenameout[80];
sprintf(filenameout,"%s.PtSpectrum_%d_%s_%s",filename, NumberPerPtBin, particle, opt);
PtSpectrum->SetMarkerStyle(27);
PtSpectrum->Draw("P");
+
}