From: kleinb Date: Mon, 3 Jan 2011 16:14:39 +0000 (+0000) Subject: Added recalcualtion of rho to background subtraction (Leticia), removed static pointe... X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=commitdiff_plain;h=a18eedbb6423cda1dadbe99d36001629e9ded524 Added recalcualtion of rho to background subtraction (Leticia), removed static pointers to avoid the same pointer in different instances --- diff --git a/JETAN/AliAnalysisTaskJetBackgroundSubtract.cxx b/JETAN/AliAnalysisTaskJetBackgroundSubtract.cxx index 1b136a273c6..63c53bb3160 100644 --- a/JETAN/AliAnalysisTaskJetBackgroundSubtract.cxx +++ b/JETAN/AliAnalysisTaskJetBackgroundSubtract.cxx @@ -59,7 +59,7 @@ AliAnalysisTaskJetBackgroundSubtract::AliAnalysisTaskJetBackgroundSubtract(): fNonStdFile(""), fReplaceString1("B0"), fReplaceString2("B%d"), - fSubtraction(kArea), + fSubtraction(kRhoRecalc), fInJetArrayList(0x0), fOutJetArrayList(0x0), fHistList(0x0) @@ -78,7 +78,7 @@ AliAnalysisTaskJetBackgroundSubtract::AliAnalysisTaskJetBackgroundSubtract(const fNonStdFile(""), fReplaceString1("B0"), fReplaceString2("B%d"), - fSubtraction(kArea), + fSubtraction(kRhoRecalc), fInJetArrayList(0x0), fOutJetArrayList(0x0), fHistList(0x0) @@ -102,7 +102,8 @@ Bool_t AliAnalysisTaskJetBackgroundSubtract::Notify() for(int iJB = 0;iJBGetEntries();iJB++){ TObjString *ostr = (TObjString*)fJBArray->At(iJB); - + + TClonesArray* jarray = 0; if(!jarray&&fAODOut){ jarray = (TClonesArray*)(fAODOut->FindListObject(ostr->GetString().Data())); @@ -156,7 +157,7 @@ void AliAnalysisTaskJetBackgroundSubtract::UserCreateOutputObjects() if (fDebug > 1) printf("AnalysisTaskJetBackgroundSubtract::UserCreateOutputObjects() \n"); if(fNonStdFile.Length()!=0){ - // + // case that we have an AOD extension we need to fetch the jets from the extended output // we identifay the extension aod event by looking for the branchname AliAODHandler *aodH = dynamic_cast(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()); @@ -261,27 +262,39 @@ void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/) - static AliAODJetEventBackground* evBkg = 0; - static TClonesArray* bkgClusters = 0; - static TString bkgClusterName(fBackgroundBranch.Data()); + AliAODJetEventBackground* evBkg = 0; + TClonesArray* bkgClusters = 0; + TClonesArray* bkgClustersRC = 0; + TString bkgClusterName(fBackgroundBranch.Data()); bkgClusterName.ReplaceAll(Form("%s_",AliAODJetEventBackground::StdBranchName()),""); - if(!evBkg&&!bkgClusters&&fAODOut){ + TString bkgClusterRCName(Form("%s%s",bkgClusterName.Data(),"RandomCone")); + + if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODOut){ evBkg = (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data())); bkgClusters = (TClonesArray*)(fAODOut->FindListObject(bkgClusterName.Data())); + bkgClustersRC = (TClonesArray*)(fAODOut->FindListObject(bkgClusterRCName.Data())); + if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data()); + if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data()); if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data()); } - if(!evBkg&&!bkgClusters&&fAODExtension){ + if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODExtension){ evBkg = (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data())); bkgClusters = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(bkgClusterName.Data())); + bkgClustersRC = (TClonesArray*)(fAODExtension->GetAOD()->FindListObject(bkgClusterRCName.Data())); if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data()); + if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data()); + if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data()); } - if(!evBkg&&!bkgClusters&&fAODIn){ + if(!evBkg&&!bkgClusters&&!bkgClustersRC&&fAODIn){ evBkg = (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data())); bkgClusters = (TClonesArray*)(fAODIn->FindListObject(bkgClusterName.Data())); + bkgClustersRC = (TClonesArray*)(fAODIn->FindListObject(bkgClusterRCName.Data())); + if(fDebug&&bkgClusters)Printf("%s:%d Background cluster branch %s found",(char*)__FILE__,__LINE__,bkgClusterName.Data()); + if(fDebug&&bkgClustersRC)Printf("%s:%d Background cluster RC branch %s found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data()); if(fDebug&&evBkg)Printf("%s:%d Backgroundbranch %s found",(char*)__FILE__,__LINE__,fBackgroundBranch.Data()); } @@ -303,28 +316,44 @@ void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/) return; } - + if(!bkgClustersRC){ + if(fDebug){ + Printf("%s:%d Background cluster RC branch %s not found",(char*)__FILE__,__LINE__,bkgClusterRCName.Data()); + PrintAODContents(); + } + PostData(1,fHistList); + return; + } // LOOP over all jet branches and subtract the background - Float_t rho = 0; - if(fSubtraction==kArea)rho = evBkg->GetBackground(2); - - + Float_t rho = 0; + Double_t meanarea = 0; + if(fSubtraction==kArea)rho = evBkg->GetBackground(2); + if(fSubtraction==kRhoRecalc){ + meanarea=evBkg->GetMeanarea(2); + rho =RecalcRho(bkgClusters,meanarea); + } + if(fSubtraction==kRhoRC)rho =RhoRC(bkgClustersRC); + for(int iJB = 0;iJBGetEntries();iJB++){ TClonesArray* jarray = (TClonesArray*)fInJetArrayList->At(iJB); TClonesArray* jarrayOut = (TClonesArray*)fOutJetArrayList->At(iJB); if(!jarray||!jarrayOut){ - Printf("%s:%d Array not found %d: %p %p",(char*)__FILE__,__LINE__,iJB,jarray,jarrayOut); - continue; - } + Printf("%s:%d Array not found %d: %p %p",(char*)__FILE__,__LINE__,iJB,jarray->GetName(),jarrayOut->GetName()); + continue; + } TH2F* h2PtInOut = (TH2F*)fHistList->FindObject(Form("h2PtInPtOut_%d",iJB)); // loop over all jets Int_t nOut = 0; + + for(int i = 0;i < jarray->GetEntriesFast();i++){ AliAODJet *jet = (AliAODJet*)jarray->At(i); AliAODJet tmpNewJet(*jet); Bool_t bAdd = false; + + if(fSubtraction==kArea){ Double_t background = rho * jet->EffectiveAreaCharged(); Float_t ptSub = jet->Pt() - background; @@ -344,12 +373,43 @@ void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/) // allows to recover old p_T and rho... tmpNewJet.SetBgEnergy(background,0); }// kAREA - else if(fSubtraction==kRhoRecalc1){ - // Put a function call to calculate rho here - // * exclude edges - // * exclude clusters with small area - // * exclude areas around leading jets - } + else if(fSubtraction==kRhoRecalc){ + Double_t background = rho * jet->EffectiveAreaCharged(); + Float_t 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 = RescaleJetMomentum(&tmpNewJet,0.1); + if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),0.1); + } + else{ + bAdd = RescaleJetMomentum(&tmpNewJet,ptSub); + if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),ptSub); + } + // add background estimates to the new jet object + // allows to recover old p_T and rho... + tmpNewJet.SetBgEnergy(background,0); + + }//kRhoRecalc + else if(fSubtraction==kRhoRC){ + Double_t background = rho * jet->EffectiveAreaCharged(); + Float_t 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 = RescaleJetMomentum(&tmpNewJet,0.1); + if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),0.1); + } + else{ + bAdd = RescaleJetMomentum(&tmpNewJet,ptSub); + if(h2PtInOut)h2PtInOut->Fill(jet->Pt(),ptSub); + } + // add background estimates to the new jet object + // allows to recover old p_T and rho... + tmpNewJet.SetBgEnergy(background,0); + + }//kRhoRC if(bAdd){ AliAODJet *newJet = new ((*jarrayOut)[nOut++]) AliAODJet(tmpNewJet); @@ -357,7 +417,12 @@ void AliAnalysisTaskJetBackgroundSubtract::UserExec(Option_t */*option*/) newJet->GetRefTracks()->Clear(); } - } + + + + + + } @@ -390,6 +455,75 @@ Bool_t AliAnalysisTaskJetBackgroundSubtract::RescaleJetMomentum(AliAODJet *jet,F return kTRUE; } +Double_t AliAnalysisTaskJetBackgroundSubtract::RecalcRho(TClonesArray* bkgClusters,Double_t meanarea){ + + Double_t ptarea=0.; + Int_t count=0; + Double_t rho=0.; + const Double_t Rlimit2=0.8*0.8; //2*jet radius. + TClonesArray* jarray=0; + + for(int iJB = 0;iJBGetEntries();iJB++){ + TObjString *ostr = (TObjString*)fInJetArrayList->At(iJB); + TString jetref=ostr->GetString().Data(); + if(jetref.Contains("ANTIKT04")){ + jarray = (TClonesArray*)fInJetArrayList->At(iJB);}} + if(jarray->GetEntries()>=2){ + AliAODJet *first = (AliAODJet*)(jarray->At(0)); + AliAODJet *second= (AliAODJet*)(jarray->At(1)); + for(Int_t k=0;kGetEntriesFast();k++){ + AliAODJet *clus = (AliAODJet*)(bkgClusters->At(k)); + if(TMath::Abs(clus->Eta())>0.5) continue; + if((clus->EffectiveAreaCharged())<0.1*meanarea) continue; + Double_t distance1=(first->Eta()-clus->Eta())*(first->Eta()-clus->Eta())+ + (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((distance1Pt()/clus->EffectiveAreaCharged(); + count=count+1;} + if(count!=0) rho=ptarea/count; + } + return rho; +} + + Double_t AliAnalysisTaskJetBackgroundSubtract::RhoRC(TClonesArray* bkgClustersRC){ + + Double_t ptarea=0.; + Int_t count=0; + Double_t rho=0.; + const Double_t Rlimit2=0.8*0.8; //2*jet radius. + TClonesArray* jarray=0; + for(int iJB = 0;iJBGetEntries();iJB++){ + TObjString *ostr = (TObjString*)fInJetArrayList->At(iJB); + TString jetref=ostr->GetString().Data(); + if(jetref.Contains("ANTIKT04")){ + jarray = (TClonesArray*)fInJetArrayList->At(iJB);}} + if(jarray->GetEntries()>=2){ + AliAODJet *first = (AliAODJet*)(jarray->At(0)); + AliAODJet *second=(AliAODJet*)(jarray->At(1)); + for(Int_t k=0;kGetEntriesFast();k++){ + AliAODJet *clus = (AliAODJet*)(bkgClustersRC->At(k)); + if(TMath::Abs(clus->Eta())>0.5) continue; + Double_t distance1=(first->Eta()-clus->Eta())*(first->Eta()-clus->Eta())+ + (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((distance1Pt()/clus->EffectiveAreaCharged(); + count=count+1;} + if(count!=0) rho=ptarea/count; } + return rho; +} + + + + + + + + + void AliAnalysisTaskJetBackgroundSubtract::ResetOutJets(){ if(!fOutJetArrayList)return; for(int iJB = 0;iJBGetEntries();iJB++){ diff --git a/JETAN/AliAnalysisTaskJetBackgroundSubtract.h b/JETAN/AliAnalysisTaskJetBackgroundSubtract.h index b2bea0d80c9..212cd3c9c6a 100644 --- a/JETAN/AliAnalysisTaskJetBackgroundSubtract.h +++ b/JETAN/AliAnalysisTaskJetBackgroundSubtract.h @@ -55,9 +55,9 @@ class AliAnalysisTaskJetBackgroundSubtract : public AliAnalysisTaskSE virtual void SetToReplace(char* c){fReplaceString1 = c;} const char* GetToReplace(){return fReplaceString1.Data();} virtual void SetReplacementMask(char* c){fReplaceString2 = c;} - const char* GetReplacementMask(){fReplaceString2.Data();} + const char* GetReplacementMask(){return fReplaceString2.Data();} - enum {kNoSubtract = 0,kArea,k4Area,kRhoRecalc1}; + enum {kNoSubtract = 0,kArea,k4Area,kRhoRecalc,kRhoRC}; private: @@ -66,7 +66,8 @@ class AliAnalysisTaskJetBackgroundSubtract : public AliAnalysisTaskSE AliAnalysisTaskJetBackgroundSubtract(const AliAnalysisTaskJetBackgroundSubtract&); AliAnalysisTaskJetBackgroundSubtract& operator=(const AliAnalysisTaskJetBackgroundSubtract&); Bool_t RescaleJetMomentum(AliAODJet *jet,Float_t pT); - + Double_t RecalcRho(TClonesArray* fbkgclusters,Double_t meanarea); + Double_t RhoRC(TClonesArray* fbkgclustersRC); void ResetOutJets(); void PrintAODContents(); @@ -84,7 +85,7 @@ class AliAnalysisTaskJetBackgroundSubtract : public AliAnalysisTaskSE TList *fOutJetArrayList; //! transient list to make ease the reset of output jets TList *fHistList; //! the histograms output list - ClassDef(AliAnalysisTaskJetBackgroundSubtract, 1) + ClassDef(AliAnalysisTaskJetBackgroundSubtract, 3) }; #endif diff --git a/JETAN/AliAnalysisTaskJetCluster.cxx b/JETAN/AliAnalysisTaskJetCluster.cxx index 3bbeb03ea42..5e67a67fa0c 100644 --- a/JETAN/AliAnalysisTaskJetCluster.cxx +++ b/JETAN/AliAnalysisTaskJetCluster.cxx @@ -665,7 +665,7 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/) } } - static AliAODJetEventBackground* externalBackground = 0; + AliAODJetEventBackground* externalBackground = 0; if(!externalBackground&&fBackgroundBranch.Length()){ externalBackground = (AliAODJetEventBackground*)(AODEvent()->FindListObject(fBackgroundBranch.Data())); Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;