AliESDtrack *track=static_cast<AliESDtrack*>(particle);
track->SetESDEvent(static_cast<AliESDEvent*>(ev)); //only in trunk...
}
-
//apply track cuts
if (fTrackFilter.IsSelected(particle)!=selectedMask) continue;
// remove contribution of all tracks to the Q-vector that are in invariant mass window
//
AliEventplane *evplane = const_cast<AliVEvent *>(ev)->GetEventplane();
+// AliEventplane *evplane = ev->GetEventplane();
if(!evplane) return;
// do not change these vectors qref
- TVector2 *qref = evplane->GetQVector();
+ TVector2 * const qref = evplane->GetQVector();
if(!qref) return;
// random subevents
TVector2 *qrsub1 = evplane->GetQsub1();
TVector2 *qrsub2 = evplane->GetQsub2();
// copy references
- TVector2 *qstd = dynamic_cast<TVector2 *>(qref->Clone());
- TVector2 *qssub1 = dynamic_cast<TVector2 *>(qrsub1->Clone());
- TVector2 *qssub2 = dynamic_cast<TVector2 *>(qrsub2->Clone());
+ TVector2 *qcorr = new TVector2(*qref);
+ TVector2 *qcsub1 = 0x0;
+ TVector2 *qcsub2 = 0x0;
+ // printf("qrsub1 %p %f \n",qrsub1,qrsub1->X());
+
+
+ // eta gap ?
+ Bool_t etagap = kFALSE;
+ for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
+ TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
+ if(cutName.Contains("eta") || cutName.Contains("Eta")) etagap=kTRUE;
+ }
- // subevents by charge separation sub1+ sub2-
- if(fLikeSignSubEvents) {
- qssub1 = dynamic_cast<TVector2 *>(qstd->Clone());
- qssub2 = dynamic_cast<TVector2 *>(qstd->Clone());
+ // subevent separation
+ if(fLikeSignSubEvents || etagap) {
+ qcsub1 = new TVector2(*qcorr);
+ qcsub2 = new TVector2(*qcorr);
Int_t ntracks=ev->GetNumberOfTracks();
Double_t cQX = evplane->GetQContributionX(track);
Double_t cQY = evplane->GetQContributionY(track);
- Short_t charge=track->Charge();
- if (charge<0) qssub1->Set(qssub1->X()-cQX, qssub1->Y()-cQY);
- if (charge>0) qssub2->Set(qssub2->X()-cQX, qssub2->Y()-cQY);
- }
- }
-
- // subevents eta division sub1+ sub2-
- Bool_t etagap = kFALSE;
- for (Int_t iCut=0; iCut<fEventPlanePreFilter.GetCuts()->GetEntries();++iCut) {
- TString cutName=fEventPlanePreFilter.GetCuts()->At(iCut)->GetName();
- if(cutName.Contains("eta") || cutName.Contains("Eta")) {
- etagap=kTRUE;
-
- qssub1 = dynamic_cast<TVector2 *>(qstd->Clone());
- qssub2 = dynamic_cast<TVector2 *>(qstd->Clone());
-
- Int_t ntracks=ev->GetNumberOfTracks();
-
- // track removals
- for (Int_t itrack=0; itrack<ntracks; ++itrack){
- AliVParticle *particle=ev->GetTrack(itrack);
- AliVTrack *track= static_cast<AliVTrack*>(particle);
- if (!track) continue;
-
- Double_t cQX = evplane->GetQContributionX(track);
- Double_t cQY = evplane->GetQContributionY(track);
-
+ // by charge sub1+ sub2-
+ if(fLikeSignSubEvents) {
+ Short_t charge=track->Charge();
+ if (charge<0) qcsub1->Set(qcsub1->X()-cQX, qcsub1->Y()-cQY);
+ if (charge>0) qcsub2->Set(qcsub2->X()-cQX, qcsub2->Y()-cQY);
+ }
+ // by eta sub1+ sub2-
+ if(etagap) {
Double_t eta=track->Eta();
- if (eta<0.0) qssub1->Set(qssub1->X()-cQX, qssub1->Y()-cQY);
- if (eta>0.0) qssub2->Set(qssub2->X()-cQX, qssub2->Y()-cQY);
+ if (eta<0.0) qcsub1->Set(qcsub1->X()-cQX, qcsub1->Y()-cQY);
+ if (eta>0.0) qcsub2->Set(qcsub2->X()-cQX, qcsub2->Y()-cQY);
}
}
}
+ else {
+ // by a random
+ qcsub1 = new TVector2(*qrsub1);
+ qcsub2 = new TVector2(*qrsub2);
+ }
// apply cuts, e.g. etagap
if(fEventPlanePreFilter.GetCuts()->GetEntries()) {
Double_t cQYsub2 = evplane->GetQContributionYsub2(track);
// update Q vectors
- qstd->Set(qstd->X()-cQX, qstd->Y()-cQY);
- qssub1->Set(qssub1->X()-cQXsub1, qssub1->Y()-cQYsub1);
- qssub2->Set(qssub2->X()-cQXsub2, qssub2->Y()-cQYsub2);
+ qcorr->Set(qcorr->X()-cQX, qcorr->Y()-cQY);
+ qcsub1->Set(qcsub1->X()-cQXsub1, qcsub1->Y()-cQYsub1);
+ qcsub2->Set(qcsub2->X()-cQXsub2, qcsub2->Y()-cQYsub2);
}
}
- if(!qstd || !qssub1 || !qssub2) return;
-
- TVector2 *qcorr = dynamic_cast<TVector2 *>(qstd->Clone());
- TVector2 *qcsub1 = dynamic_cast<TVector2 *>(qssub1->Clone());
- TVector2 *qcsub2 = dynamic_cast<TVector2 *>(qssub2->Clone());
-
-
// POI (particle of interest) rejection
Int_t pairIndex=GetPairIndex(arr1,arr2);
qcsub2->Set(qcsub2->X()-cQXsub2, qcsub2->Y()-cQYsub2);
}
+ // printf("qrsub1 %p %f \t qcsub1 %p %f \n",qrsub1,qrsub1->X(),qcsub1,qcsub1->X());
// set AliEventplane with corrected values
cevplane->SetQVector(qcorr);
cevplane->SetQsub(qcsub1, qcsub2);
kFitted,
kLikeSign,
kLikeSignArithm,
+ kLikeSignRcorr,
+ kLikeSignArithmRcorr,
kEventMixing,
kRotation
};
void SetFitRange(Double_t min, Double_t max) {fFitMin=min; fFitMax=max;}
void SetRebin(Int_t factor) {fRebin = factor;}
void SetMethod(EBackgroundMethod method) {fMethod = method;}
+ Int_t GetMethod() const { return (Int_t)fMethod; }
const TVectorD& GetValues() const {return fValues;}
const TVectorD& GetErrors() const {return fErrors;}
switch ( fMethod ){
case kLikeSign :
case kLikeSignArithm :
+ case kLikeSignRcorr:
+ case kLikeSignArithmRcorr:
ProcessLS(arrhist); // process like-sign subtraction method
break;
Float_t background = 2*TMath::Sqrt(pp*mm);
Float_t ebackground = TMath::Sqrt(mm+pp);
- if (fMethod==kLikeSignArithm){
+ if (fMethod==kLikeSignArithm || fMethod==kLikeSignArithmRcorr ){
//Arithmetic mean instead of geometric
background=(pp+mm);
ebackground=TMath::Sqrt(pp+mm);
}
//correct LS spectrum bin-by-bin with R factor obtained in mixed events
- if(fMixingCorr) {
+ if(fMixingCorr || fMethod==kLikeSignRcorr || fMethod==kLikeSignArithmRcorr) {
TH1* histMixPP = (TH1*)(arrhist->At(AliDielectron::kEv1PEv2P))->Clone("mixPP"); // ++ ME
TH1* histMixMM = (TH1*)(arrhist->At(AliDielectron::kEv1MEv2M))->Clone("mixMM"); // -- ME
- TH1* histMixPM = (TH1*)(arrhist->At(AliDielectron::kEv1PEv2M))->Clone("mixPM"); // +- ME
+
+ TH1* histMixPM = 0x0;
+ if (arrhist->At(AliDielectron::kEv1MEv2P)){
+ histMixPM = (TH1*)(arrhist->At(AliDielectron::kEv1MEv2P))->Clone("mixPM"); // -+ ME
+ }
+
+ if (arrhist->At(AliDielectron::kEv1PEv2M)){
+ TH1 *htmp=(TH1*)(arrhist->At(AliDielectron::kEv1PEv2M));
+ if (!histMixPM) fHistDataME = (TH1*)htmp->Clone("mixPM"); // +- ME
+ else histMixPM->Add(htmp);
+ }
+
histMixPM->Sumw2();
- histMixPM->Add((TH1*)arrhist->At(AliDielectron::kEv1MEv2P)); // -+ ME
// rebin the histograms
if (fRebin>1) {
//
// event mixing of +- and -+
//
- fHistDataPM = (TH1*)(arrhist->At(AliDielectron::kEv1PM))->Clone("histPMSE"); // +- SE
- fHistDataME = (TH1*)(arrhist->At(AliDielectron::kEv1MEv2P))->Clone("histMPME"); // -+ ME
- fHistDataPM->Sumw2();
- fHistDataME->Sumw2();
- fHistDataPM->SetDirectory(0);
- fHistDataME->SetDirectory(0);
-
- fHistDataME->Add((TH1*)arrhist->At(AliDielectron::kEv1PEv2M)); // +- ME
- // rebin the histograms
- if (fRebin>1) {
- fHistDataPM->Rebin(fRebin);
- fHistDataME->Rebin(fRebin);
+ if (!arrhist->At(AliDielectron::kEv1PM) || !(arrhist->At(AliDielectron::kEv1MEv2P) || arrhist->At(AliDielectron::kEv1PEv2M)) ){
+ AliError("Either OS or mixed histogram missing");
+ return;
}
- fHistSignal = new TH1D("HistSignal", "Mixed events background substracted signal",
- fHistDataPM->GetXaxis()->GetNbins(),
- fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
- fHistSignal->SetDirectory(0);
- fHistBackground = new TH1D("HistBackground", "background contribution from mixed events",
- fHistDataPM->GetXaxis()->GetNbins(),
- fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
- fHistBackground->SetDirectory(0);
+ delete fHistDataPM; fHistDataPM=0x0;
+ delete fHistDataME; fHistDataME=0x0;
+ delete fHistBackground; fHistBackground=0x0;
+
+ fHistDataPM = (TH1*)(arrhist->At(AliDielectron::kEv1PM))->Clone("histPM"); // +- SE
+ fHistDataPM->Sumw2();
+ fHistDataPM->SetDirectory(0x0);
- // fill out background and subtracted histogram
- for(Int_t ibin=1; ibin<=fHistDataPM->GetXaxis()->GetNbins(); ibin++) {
- Float_t pm = fHistDataPM->GetBinContent(ibin);
- Float_t epm = fHistDataPM->GetBinError(ibin);
- Float_t background = fHistDataME->GetBinContent(ibin);
- Float_t ebackground = fHistDataME->GetBinError(ibin);
+ if (arrhist->At(AliDielectron::kEv1MEv2P)){
+ fHistDataME = (TH1*)(arrhist->At(AliDielectron::kEv1MEv2P))->Clone("histMPME"); // -+ ME
+ }
+
+ if (arrhist->At(AliDielectron::kEv1PEv2M)){
+ TH1 *htmp=(TH1*)(arrhist->At(AliDielectron::kEv1PEv2M));
+ if (!fHistDataME) fHistDataME = (TH1*)htmp->Clone("histMPME"); // -+ ME
+ else fHistDataME->Add(htmp);
+ }
- fHistSignal->SetBinContent(ibin, pm);
- fHistSignal->SetBinError(ibin, epm);
- fHistBackground->SetBinContent(ibin, background);
- fHistBackground->SetBinError(ibin, ebackground);
+ fHistBackground = (TH1*)fHistDataME->Clone("ME_Background");
+ fHistBackground->SetDirectory(0x0);
+ fHistBackground->Sumw2();
+
+ // rebin the histograms
+ if (fRebin>1) {
+ fHistDataPM->Rebin(fRebin);
+ fHistDataME->Rebin(fRebin);
+ fHistBackground->Rebin(fRebin);
}
+ //scale histograms to match integral between fScaleMin and fScaleMax
+ // or if fScaleMax < fScaleMin use fScaleMin as scale factor
if (fScaleMax>fScaleMin) fScaleFactor=ScaleHistograms(fHistDataPM,fHistBackground,fScaleMin,fScaleMax);
- else if (fScaleMin>0.){
+ else if (fScaleMin>0.){
fScaleFactor=fScaleMin;
fHistBackground->Scale(fScaleFactor);
}
- //subract background
- fHistSignal->Add(fHistBackground,-1);
+ fHistSignal=(TH1*)fHistDataPM->Clone("Signal");
+ fHistSignal->Add(fHistBackground,-1.);
- // signal
+ // signal
fValues(0) = fHistSignal->IntegralAndError(fHistSignal->FindBin(fIntMin),
- fHistSignal->FindBin(fIntMax), fErrors(0));
+ fHistSignal->FindBin(fIntMax), fErrors(0));
// background
fValues(1) = fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
- fHistBackground->FindBin(fIntMax),
- fErrors(1));
+ fHistBackground->FindBin(fIntMax),
+ fErrors(1));
// S/B and significance
SetSignificanceAndSOB();
fActiveCuts[i]=0;
fCutMin[i]=0;
fCutMax[i]=0;
+ fCutExclude[i]=kFALSE;
}
}