#include <TH2F.h>
#include <THnSparse.h>
#include <TSystem.h>
+#include <TObjArray.h>
#include <TInterpreter.h>
#include <TList.h>
#include <TLorentzVector.h>
ClassImp(AliAnalysisTaskJetBackgroundSubtract)
AliAnalysisTaskJetBackgroundSubtract::~AliAnalysisTaskJetBackgroundSubtract(){
+ //
+ // destructor
+ //
delete fJBArray;
delete fOutJetArrayList;
delete fInJetArrayList;
fAODOut(0x0),
fAODIn(0x0),
fAODExtension(0x0),
- fJBArray(new TObjArray()),
+ fJBArray(0),
fBackgroundBranch(""),
fNonStdFile(""),
fReplaceString1("B0"),
fReplaceString2("B%d"),
fSubtraction(k4Area),
+ fKeepJets(kFALSE),
+ fExclude2Leading(kTRUE),
fInJetArrayList(0x0),
fOutJetArrayList(0x0),
fh2CentvsRho(0x0),
fAODOut(0x0),
fAODIn(0x0),
fAODExtension(0x0),
- fJBArray(new TObjArray()),
+ fJBArray(0),
fBackgroundBranch(""),
fNonStdFile(""),
fReplaceString1("B0"),
fReplaceString2("B%d"),
fSubtraction(k4Area),
+ fKeepJets(kFALSE),
+ fExclude2Leading(kTRUE),
fInJetArrayList(0x0),
fOutJetArrayList(0x0),
fh2CentvsRho(0x0),
Bool_t AliAnalysisTaskJetBackgroundSubtract::Notify()
{
+ //
+ // exec with every new file
//
fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
fInJetArrayList->Clear();
fOutJetArrayList->Clear();
- for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
+ for(int iJB = 0;iJB<(fJBArray?fJBArray->GetEntries():0);iJB++){
TObjString *ostr = (TObjString*)fJBArray->At(iJB);
if(!fInJetArrayList)fInJetArrayList =new TList();
if(!fOutJetArrayList)fOutJetArrayList =new TList();
- for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
+ for(int iJB = 0;iJB<(fJBArray?fJBArray->GetEntries():0);iJB++){
TObjString *ostr = (TObjString*)fJBArray->At(iJB);
TString newName(ostr->GetString().Data());
if(!newName.Contains(fReplaceString1.Data())){
const Double_t xminDelta[nSparseBinsDelta] = {-120.5, 0, 0, 0};
const Double_t xmaxDelta[nSparseBinsDelta] = { 120.5, 1.0, 100,5000};
- for(int iJB = 0;iJB<fJBArray->GetEntries();iJB++){
+ for(int iJB = 0;iJB<(fJBArray?fJBArray->GetEntries():0);iJB++){
TObjString *ostr = (TObjString*)fJBArray->At(iJB);
TString oldName(ostr->GetString().Data());
TString newName(ostr->GetString().Data());
if(fBackgroundBranch.Length()==0)
AliError(Form("%s:%d No BackgroundBranch defined",(char*)__FILE__,__LINE__));
- if(fJBArray->GetEntries()==0)
+ if((fJBArray?fJBArray->GetEntries():0)==0)
AliError(Form("%s:%d No Jet Branches defined defined",(char*)__FILE__,__LINE__));
}
void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/)
{
+ //
+ // Execute for every selected event
+ //
if (fDebug > 1) printf("AnalysisTaskJetBackgroundSubtract::UserExec() \n");
ResetOutJets();
- if(fBackgroundBranch.Length()==0||fJBArray->GetEntries()==0){
+ if(fBackgroundBranch.Length()==0||(fJBArray?fJBArray->GetEntries():0)==0){
if(fDebug)Printf("%s:%d No background subtraction done",(char*)__FILE__,__LINE__);
PostData(1,fHistList);
}
TLorentzVector backgroundv;
Float_t cent=0.;
- if(fAODOut)cent = fAODOut->GetHeader()->GetCentrality();
- if(fAODIn) cent = fAODIn->GetHeader()->GetCentrality();
+ if(fAODOut)cent = ((AliVAODHeader*)fAODOut->GetHeader())->GetCentrality();
+ if(fAODIn) cent = ((AliVAODHeader*)fAODIn->GetHeader() )->GetCentrality();
if(evBkg)sigma=evBkg->GetSigma(1);
+
+ if(fSubtraction==kArea) {rho = evBkg->GetBackground(1);
+ if(fExclude2Leading==kFALSE){
+ rho=evBkg->GetBackground(3);
+ sigma=evBkg->GetSigma(3);}}
- if(fSubtraction==kArea) rho = evBkg->GetBackground(1);
if(fSubtraction==k4Area){
+
rho = evBkg->GetBackground(0);
sigma=evBkg->GetSigma(0);
}
+
+
if(fSubtraction==kRhoRecalc){
meanarea=evBkg->GetMeanarea(1);
if(fDebug>2){
Printf("%s:%d Jet %d %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub);
}
- if(ptSub<0){
+ if(ptSub<=0){
// optionally rescale it and keep??
- bAdd = false; // RescaleJetMomentum(&tmpNewJet,0.1);
+ if(fKeepJets){
+ bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
+ }
+ else{
+ bAdd = false;
+ }
}
else{
bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
// add background estimates to the new jet object
// allows to recover old p_T and rho...
tmpNewJet.SetBgEnergy(background,0);
+ tmpNewJet.SetPtSubtracted(ptSub,0);
}// kAREA
else if(fSubtraction==kRhoRecalc){
Double_t background = rho * jet->EffectiveAreaCharged();
ptSub = jet->Pt() - background;
if(fDebug>2){
Printf("%s:%d Jet %d %3.3f %3.3f %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub,background,rho);}
- if(ptSub<0){
- // optionally rescale it and keep??
- bAdd = false;// RescaleJetMomentum(&tmpNewJet,0.1);
+ if(ptSub<=0){
+ // optionally rescale it and keep
+ if(fKeepJets){
+ bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
+ }
+ else{
+ bAdd = false;
+ }
}
else{
bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
// add background estimates to the new jet object
// allows to recover old p_T and rho...
tmpNewJet.SetBgEnergy(background,0);
+ tmpNewJet.SetPtSubtracted(ptSub,0);
}//kRhoRecalc
else if(fSubtraction==kRhoRC){
Double_t background = rho * jet->EffectiveAreaCharged();
ptSub = jet->Pt() - background;
if(fDebug>2){ Printf("%s:%d Jet %d %3.3f %3.3f %3.3f %3.3f",(char*)__FILE__,__LINE__,i,jet->Pt(),ptSub,background,rho);}
- if(ptSub<0){
- // optionally rescale it and keep??
- bAdd = false; // RescaleJetMomentum(&tmpNewJet,0.1);
+ if(ptSub<=0){
+ if(fKeepJets){
+ bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
+ }
+ else{
+ bAdd = false;
+ }
}
else{
bAdd = RescaleJetMomentum(&tmpNewJet,ptSub);
// add background estimates to the new jet object
// allows to recover old p_T and rho...
tmpNewJet.SetBgEnergy(background,0);
-
+ tmpNewJet.SetPtSubtracted(ptSub,0);
}//kRhoRC
else if(fSubtraction==k4Area&&jet->VectorAreaCharged()){
backgroundv.SetPxPyPzE(rho*(jet->VectorAreaCharged())->Px(),rho*(jet->VectorAreaCharged())->Py(),rho*(jet->VectorAreaCharged())->Pz(),rho*(jet->VectorAreaCharged())->E());
ptSub = jet->Pt()-backgroundv.Pt();
- if((backgroundv.E()>jet->E())&&(backgroundv.Pt()>jet->Pt())){
- // optionally rescale it and keep??
- bAdd = false; // RescaleJetMomentum(&tmpNewJet,0.1);
+ if((backgroundv.E()>=jet->E())||(backgroundv.Pt()>=jet->Pt())){
+ if(fKeepJets){
+ bAdd = RescaleJetMomentum(&tmpNewJet,0.1);
+ }
+ else{
+ bAdd = false;
+ }
}
else{
bAdd = RescaleJet4vector(&tmpNewJet,backgroundv);
}
// add background estimates to the new jet object
// allows to recover old p_T and rho...
- tmpNewJet.SetBgEnergy(backgroundv.P(),0);
+ tmpNewJet.SetBgEnergy(backgroundv.Pt(),0);
+ tmpNewJet.SetPtSubtracted(ptSub,0);
}//kArea4vector
-
-
-
-
-
if(bAdd){
AliAODJet *newJet = new ((*jarrayOut)[nOut++]) AliAODJet(tmpNewJet);
// what about track references, clear for now...
Double_t AliAnalysisTaskJetBackgroundSubtract::RecalcRho(TClonesArray* bkgClusters,Double_t meanarea){
+ //
+ // recalc rhoo
+ //
+
Double_t ptarea=0.;
Int_t count=0;
Double_t rho=0.;
- const Double_t Rlimit2=0.8*0.8; //2*jet radius.
+ const Double_t rLimit2=0.8*0.8; //2*jet radius.
TClonesArray* jarray=0;
for(int iJB = 0;iJB<fInJetArrayList->GetEntries();iJB++){
(first->Phi()-clus->Phi())*(first->Phi()-clus->Phi());
Double_t distance2= (second->Eta()-clus->Eta())*(second->Eta()-clus->Eta())+
(second->Phi()-clus->Phi())*(second->Phi()-clus->Phi());
- if((distance1<Rlimit2)||(distance2<Rlimit2)) continue;
+ if((distance1<rLimit2)||(distance2<rLimit2)) continue;
ptarea=ptarea+clus->Pt()/clus->EffectiveAreaCharged();
count=count+1;}
if(count!=0) rho=ptarea/count;
return rho;
}
- Double_t AliAnalysisTaskJetBackgroundSubtract::RhoRC(TClonesArray* bkgClustersRC){
+Double_t AliAnalysisTaskJetBackgroundSubtract::RhoRC(TClonesArray* bkgClustersRC){
+
+ //
+ // calc rho from random cones
+ //
Double_t ptarea=0.;
Int_t count=0;
Double_t rho=0.;
- const Double_t Rlimit2=0.8*0.8; //2*jet radius.
+ const Double_t rLimit2=0.8*0.8; //2*jet radius.
TClonesArray* jarray=0;
for(int iJB = 0;iJB<fInJetArrayList->GetEntries();iJB++){
TObjString *ostr = (TObjString*)fInJetArrayList->At(iJB);
(first->Phi()-clus->Phi())*(first->Phi()-clus->Phi());
Double_t distance2= (second->Eta()-clus->Eta())*(second->Eta()-clus->Eta())+
(second->Phi()-clus->Phi())*(second->Phi()-clus->Phi());
- if((distance1<Rlimit2)||(distance2<Rlimit2)) continue;
+ if((distance1<rLimit2)||(distance2<rLimit2)) continue;
ptarea=ptarea+clus->Pt()/clus->EffectiveAreaCharged();
count=count+1;}
if(count!=0) rho=ptarea/count; }
void AliAnalysisTaskJetBackgroundSubtract::ResetOutJets(){
+ //
+ // Reset the output jets
+ //
+
if(!fOutJetArrayList)return;
for(int iJB = 0;iJB<fOutJetArrayList->GetEntries();iJB++){
TClonesArray* jarray = (TClonesArray*)fOutJetArrayList->At(iJB);
void AliAnalysisTaskJetBackgroundSubtract::PrintAODContents(){
+
+ //
+ // guess from the name what this function does
+ //
+
if(fAODIn){
Printf("%s:%d >>>>>> Input",(char*)__FILE__,__LINE__);
fAODIn->Print();
}
Int_t AliAnalysisTaskJetBackgroundSubtract::MultFromJetRefs(TClonesArray* jets){
+ //
+ // calculate multiplicty based on jet references
+ //
+
if(!jets)return 0;
Int_t refMult = 0;
return refMult;
}
+
+void AliAnalysisTaskJetBackgroundSubtract::AddJetBranch(const char* c){
+ if(!fJBArray)fJBArray = new TObjArray();
+ fJBArray->Add(new TObjString(c));
+}