#pragma link C++ class AliDielectronSignalExt+;
#pragma link C++ class AliDielectronSpectrum+;
#pragma link C++ class AliDielectronDebugTree+;
+#pragma link C++ class AliDielectronTrackRotator+;
#pragma link C++ class AliDielectronPID+;
+#pragma link C++ class AliDielectronCutGroup+;
+#pragma link C++ class AliDielectronEventCuts+;
#endif
#include <AliVEvent.h>
#include <AliInputEventHandler.h>
#include <AliESDInputHandler.h>
+#include <AliAODInputHandler.h>
#include "AliDielectron.h"
#include "AliDielectronMC.h"
//_________________________________________________________________________________
AliAnalysisTaskDielectronFilter::AliAnalysisTaskDielectronFilter() :
- AliAnalysisTaskSE(),
- fDielectron(0),
- fSelectPhysics(kTRUE),
- fTriggerMask(AliVEvent::kMB),
- fEventStat(0x0)
+AliAnalysisTaskSE(),
+fDielectron(0),
+fSelectPhysics(kTRUE),
+fTriggerMask(AliVEvent::kMB),
+fEventStat(0x0)
{
//
// Constructor
//_________________________________________________________________________________
AliAnalysisTaskDielectronFilter::AliAnalysisTaskDielectronFilter(const char *name) :
- AliAnalysisTaskSE(name),
- fDielectron(0),
- fSelectPhysics(kTRUE),
- fTriggerMask(AliVEvent::kMB),
- fEventStat(0x0)
+AliAnalysisTaskSE(name),
+fDielectron(0),
+fSelectPhysics(kTRUE),
+fTriggerMask(AliVEvent::kMB),
+fEventStat(0x0)
{
//
// Constructor
fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel.");
fEventStat->GetXaxis()->SetBinLabel(3,"After Cand. Sel.");
}
-
+
PostData(2,fEventStat);
}
//load esd pid bethe bloch parameters depending on the existance of the MC handler
// yes: MC parameters
// no: data parameters
- if (!AliDielectronVarManager::GetESDpid()){
- if (AliDielectronMC::Instance()->HasMC()) {
- AliDielectronVarManager::InitESDpid();
- } else {
- AliDielectronVarManager::InitESDpid(1);
+
+ //ESD case
+ if (man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()){
+ if (!AliDielectronVarManager::GetESDpid()){
+
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitESDpid();
+ } else {
+ AliDielectronVarManager::InitESDpid(1);
+ }
+ }
+ }
+ //AOD case
+ if (man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()){
+ if (!AliDielectronVarManager::GetAODpidUtil()){
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitAODpidUtil();
+ } else {
+ AliDielectronVarManager::InitAODpidUtil(1);
+ }
}
}
}
+
// Was event selected ?
AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
UInt_t isSelected = AliVEvent::kAny;
isSelected = inputHandler->IsEventSelected();
isSelected&=fTriggerMask;
}
-
+
//Before physics selection
fEventStat->Fill(0.);
if (isSelected==0) {
fDielectron->Process(InputEvent());
if(fDielectron->HasCandidates()){
- //If input event is an AliESDevent
- // replace the references of the legs with the AOD references
- if(man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()){
- AliAODEvent *aod = ((AliAODHandler*)((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler()))->GetAOD();
- TObjArray *obj = 0x0;
- AliAODTrack *leg1 = 0x0;
- AliAODTrack *leg2 = 0x0;
- for(Int_t i=0; i < 10; i++ ){
- obj = (TObjArray*)(*(fDielectron->GetPairArraysPointer()))->UncheckedAt(i);
- if(!obj) continue;
- for(int j=0;j<obj->GetEntriesFast();j++)
- {
- AliDielectronPair *pairObj = (AliDielectronPair*)obj->UncheckedAt(j);
- Int_t id1 = ((AliESDtrack*)pairObj->GetFirstDaughter())->GetID();
- Int_t id2 = ((AliESDtrack*)pairObj->GetSecondDaughter())->GetID();
- for(Int_t it=0;it<aod->GetNumberOfTracks();it++){
- if(aod->GetTrack(it)->GetID() == id1) leg1 = aod->GetTrack(it);
- if(aod->GetTrack(it)->GetID() == id2) leg2 = aod->GetTrack(it);
- }
- if(!leg1 || !leg2) continue;
- pairObj->SetRefFirstDaughter(leg1);
- pairObj->SetRefSecondDaughter(leg2);
+ AliAODEvent *aod = ((AliAODHandler*)((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler()))->GetAOD();
+
+ //replace the references of the legs with the AOD references
+ TObjArray *obj = 0x0;
+ for(Int_t i=0; i < 10; i++ ){
+ obj = (TObjArray*)((*(fDielectron->GetPairArraysPointer()))->UncheckedAt(i));
+ if(!obj) continue;
+ for(int j=0;j<obj->GetEntriesFast();j++){
+ AliAODTrack *leg1 = 0x0;
+ AliAODTrack *leg2 = 0x0;
+ AliDielectronPair *pairObj = (AliDielectronPair*)obj->UncheckedAt(j);
+ Int_t id1 = ((AliVTrack*)pairObj->GetFirstDaughter())->GetID();
+ Int_t id2 = ((AliVTrack*)pairObj->GetSecondDaughter())->GetID();
+
+ for(Int_t it=0;it<aod->GetNumberOfTracks();it++){
+ if(aod->GetTrack(it)->GetID() == id1) leg1 = aod->GetTrack(it);
+ if(aod->GetTrack(it)->GetID() == id2) leg2 = aod->GetTrack(it);
+ }
+ if(!leg1 || !leg2) continue;
+
+ if(man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class()){
+ leg1->ResetBit(kIsReferenced);
+ leg1->SetUniqueID(0);
+ leg2->ResetBit(kIsReferenced);
+ leg2->SetUniqueID(0);
}
+ pairObj->SetRefFirstDaughter(leg1);
+ pairObj->SetRefSecondDaughter(leg2);
}
}
-
- AliAODExtension *extDielectron = dynamic_cast<AliAODHandler*>
- ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler())->GetFilteredAOD("AliAOD.Dielectron.root");
+
+ AliAODExtension *extDielectron = dynamic_cast<AliAODHandler*> ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler())->GetFilteredAOD("AliAOD.Dielectron.root");
extDielectron->SelectEvent();
//after candidate selection
fEventStat->Fill(2.);
//see if dielectron candidate branch exists, if not create is
TTree *t=extDielectron->GetTree();
+
+ if(!t->GetListOfBranches()->GetEntries() && man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class())
+ t->Branch(aod->GetList());
+
if (!t->GetBranch("dielectrons")){
t->Bronch("dielectrons","TObjArray",fDielectron->GetPairArraysPointer());
}
+
+ if(man->GetInputEventHandler()->IsA()==AliAODInputHandler::Class()) t->Fill();
}
PostData(1, const_cast<THashList*>(fDielectron->GetHistogramList()));
//load esd pid bethe bloch parameters depending on the existance of the MC handler
// yes: MC parameters
// no: data parameters
- if (!AliDielectronVarManager::GetESDpid()){
- if (AliDielectronMC::Instance()->HasMC()) {
- AliDielectronVarManager::InitESDpid();
- } else {
- AliDielectronVarManager::InitESDpid(1);
+
+ //ESD case
+ if (man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()){
+ if (!AliDielectronVarManager::GetESDpid()){
+
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitESDpid();
+ } else {
+ AliDielectronVarManager::InitESDpid(1);
+ }
+ }
+ }
+ //AOD case
+ if (man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()){
+ if (!AliDielectronVarManager::GetAODpidUtil()){
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitAODpidUtil();
+ } else {
+ AliDielectronVarManager::InitAODpidUtil(1);
+ }
}
}
}
#include <AliCFContainer.h>
#include <AliInputEventHandler.h>
#include <AliESDInputHandler.h>
+#include <AliAODHandler.h>
#include <AliAnalysisManager.h>
#include <AliVEvent.h>
+#include <AliTriggerAnalysis.h>
#include "AliDielectron.h"
#include "AliDielectronHistos.h"
fListCF(),
fSelectPhysics(kFALSE),
fTriggerMask(AliVEvent::kMB),
+ fTriggerOnV0AND(kFALSE),
+ fRejectPileup(kFALSE),
+ fTriggerAnalysis(0x0),
+ fEventFilter(0x0),
fEventStat(0x0)
{
//
fListCF(),
fSelectPhysics(kFALSE),
fTriggerMask(AliVEvent::kMB),
+ fTriggerOnV0AND(kFALSE),
+ fRejectPileup(kFALSE),
+ fTriggerAnalysis(0x0),
+ fEventFilter(0x0),
fEventStat(0x0)
{
//
DefineOutput(3, TH1D::Class());
fListHistos.SetName("Dielectron_Histos_Multi");
fListCF.SetName("Dielectron_CF_Multi");
+ fListDielectron.SetOwner();
+ fListHistos.SetOwner();
+ fListCF.SetOwner();
}
if (!fListHistos.IsEmpty()||!fListCF.IsEmpty()) return; //already initialised
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ Bool_t isESD=man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class();
+// Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODHandler::Class();
+
TIter nextDie(&fListDielectron);
AliDielectron *die=0;
while ( (die=static_cast<AliDielectron*>(nextDie())) ){
}
Int_t cuts=fListDielectron.GetEntries();
- Int_t nbins=2+2*cuts;
+ Int_t nbins=kNbinsEvent+2*cuts;
if (!fEventStat){
fEventStat=new TH1D("hEventStat","Event statistics",nbins,0,nbins);
fEventStat->GetXaxis()->SetBinLabel(1,"Before Phys. Sel.");
fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel.");
+ if (fTriggerOnV0AND&&isESD) fEventStat->GetXaxis()->SetBinLabel(3,"V0and triggers");
+ if (fEventFilter) fEventStat->GetXaxis()->SetBinLabel(4,"After Event Filter");
+ if (fRejectPileup) fEventStat->GetXaxis()->SetBinLabel(5,"After Pileup rejection");
+
for (Int_t i=0; i<cuts; ++i){
- fEventStat->GetXaxis()->SetBinLabel(3+2*i,Form("#splitline{1 candidate}{%s}",fListDielectron.At(i)->GetName()));
- fEventStat->GetXaxis()->SetBinLabel(4+2*i,Form("#splitline{With >1 candidate}{%s}",fListDielectron.At(i)->GetName()));
+ fEventStat->GetXaxis()->SetBinLabel((kNbinsEvent+1)+2*i,Form("#splitline{1 candidate}{%s}",fListDielectron.At(i)->GetName()));
+ fEventStat->GetXaxis()->SetBinLabel((kNbinsEvent+2)+2*i,Form("#splitline{With >1 candidate}{%s}",fListDielectron.At(i)->GetName()));
}
}
+
+ if (!fTriggerAnalysis) fTriggerAnalysis=new AliTriggerAnalysis;
+ fTriggerAnalysis->EnableHistograms();
+ fTriggerAnalysis->SetAnalyzeMC(AliDielectronMC::Instance()->HasMC());
PostData(1, &fListHistos);
PostData(2, &fListCF);
AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
AliESDInputHandler *esdHandler=0x0;
+ Bool_t isESD=man->GetInputEventHandler()->IsA()==AliESDInputHandler::Class();
+ Bool_t isAOD=man->GetInputEventHandler()->IsA()==AliAODHandler::Class();
if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
AliDielectronVarManager::SetESDpid(esdHandler->GetESDpid());
} else {
//load esd pid bethe bloch parameters depending on the existance of the MC handler
// yes: MC parameters
// no: data parameters
- if (!AliDielectronVarManager::GetESDpid()){
- if (AliDielectronMC::Instance()->HasMC()) {
- AliDielectronVarManager::InitESDpid();
- } else {
- AliDielectronVarManager::InitESDpid(1);
+
+ //ESD case
+ if (isESD){
+ if (!AliDielectronVarManager::GetESDpid()){
+
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitESDpid();
+ } else {
+ AliDielectronVarManager::InitESDpid(1);
+ }
+ }
+ }
+ //AOD case
+ if (isAOD){
+ if (!AliDielectronVarManager::GetAODpidUtil()){
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitAODpidUtil();
+ } else {
+ AliDielectronVarManager::InitAODpidUtil(1);
+ }
}
}
}
}
//Before physics selection
- fEventStat->Fill(0.);
+ fEventStat->Fill(kAllEvents);
if (isSelected==0) {
PostData(3,fEventStat);
return;
}
//after physics selection
- fEventStat->Fill(1.);
+ fEventStat->Fill(kSelectedEvents);
+
+ //V0and
+ if (fTriggerOnV0AND&&isESD){
+ if (!fTriggerAnalysis->IsOfflineTriggerFired(static_cast<AliESDEvent*>(InputEvent()), AliTriggerAnalysis::kV0AND)) return;
+ }
+ fEventStat->Fill(kV0andEvents);
+
+ //event filter
+ if (fEventFilter) {
+ if (!fEventFilter->IsSelected(InputEvent())) return;
+ }
+ fEventStat->Fill(kFilteredEvents);
+
+ //pileup
+ if (fRejectPileup){
+ if (InputEvent()->IsPileupFromSPD(3,0.8,3.,2.,5.)) return;
+ }
+ fEventStat->Fill(kPileupEvents);
//bz for AliKF
Double_t bz = InputEvent()->GetMagneticField();
AliKFParticle::SetField( bz );
+ AliDielectronPID::SetCorrVal((Double_t)InputEvent()->GetRunNumber());
+
//Process event in all AliDielectron instances
TIter nextDie(&fListDielectron);
AliDielectron *die=0;
die->Process(InputEvent());
if (die->HasCandidates()){
Int_t ncandidates=die->GetPairArray(1)->GetEntriesFast();
- if (ncandidates==1) fEventStat->Fill(3+2*idie);
- else if (ncandidates>1) fEventStat->Fill(4+2*idie);
+ if (ncandidates==1) fEventStat->Fill((kNbinsEvent+1)+2*idie);
+ else if (ncandidates>1) fEventStat->Fill((kNbinsEvent+2)+2*idie);
}
++idie;
}
class AliDielectron;
class TH1D;
+class AliAnalysisCuts;
+class AliTriggerAnalysis;
class AliAnalysisTaskMultiDielectron : public AliAnalysisTaskSE {
virtual void UserCreateOutputObjects();
virtual void FinishTaskOutput();
//temporary
- virtual void NotifyRun(){AliDielectronPID::SetCorrVal((Double_t)fCurrentRunNumber);}
+// virtual void NotifyRun(){AliDielectronPID::SetCorrVal((Double_t)fCurrentRunNumber);}
void UsePhysicsSelection(Bool_t phy=kTRUE) {fSelectPhysics=phy;}
void SetTriggerMask(UInt_t mask) {fTriggerMask=mask;}
UInt_t GetTriggerMask() const { return fTriggerMask; }
-
+
+ void SetEventFilter(AliAnalysisCuts * const filter) {fEventFilter=filter;}
+ void SetTriggerOnV0AND(Bool_t v0and=kTRUE) { fTriggerOnV0AND=v0and; }
+ void SetRejectPileup(Bool_t pileup=kTRUE) { fRejectPileup=pileup; }
void AddDielectron(AliDielectron * const die) { fListDielectron.Add(die); }
private:
-
+ enum {kAllEvents=0, kSelectedEvents, kV0andEvents, kFilteredEvents, kPileupEvents, kNbinsEvent};
TList fListDielectron; // List of dielectron framework instances
TList fListHistos; //! List of histogram manager lists in the framework classes
TList fListCF; //! List with CF Managers
Bool_t fSelectPhysics; // Whether to use physics selection
UInt_t fTriggerMask; // Event trigger mask
+ Bool_t fTriggerOnV0AND; // if to trigger on V0and
+ Bool_t fRejectPileup; // pileup rejection wanted
+
+ AliTriggerAnalysis *fTriggerAnalysis; //! trigger analysis class
+ AliAnalysisCuts *fEventFilter; // event filter
+
TH1D *fEventStat; //! Histogram with event statistics
AliAnalysisTaskMultiDielectron(const AliAnalysisTaskMultiDielectron &c);
8: ev2+ ev2- (same event unlike sign)
9: ev2- ev2- (same event like sign -)
-
+10: ev1+ ev1- (same event track rotation)
*/
// //
#include "AliDielectronCF.h"
#include "AliDielectronMC.h"
#include "AliDielectronVarManager.h"
+#include "AliDielectronTrackRotator.h"
#include "AliDielectronDebugTree.h"
#include "AliDielectron.h"
"ev2-"
};
-const char* AliDielectron::fgkPairClassNames[10] = {
+const char* AliDielectron::fgkPairClassNames[11] = {
"ev1+_ev1+",
"ev1+_ev1-",
"ev1-_ev1-",
"ev1+_ev2-",
"ev1-_ev2-",
"ev2+_ev2-",
- "ev2-_ev2-"
+ "ev2-_ev2-",
+ "ev1+_ev1-_TR"
};
//________________________________________________________________
fHistos(0x0),
fPairCandidates(new TObjArray(10)),
fCfManagerPair(0x0),
+ fTrackRotator(0x0),
fDebugTree(0x0)
{
//
fHistos(0x0),
fPairCandidates(new TObjArray(10)),
fCfManagerPair(0x0),
+ fTrackRotator(0x0),
fDebugTree(0x0)
{
//
// Initialise objects
//
if (fCfManagerPair) fCfManagerPair->InitialiseContainer(fPairFilter);
+ if (fTrackRotator) fTrackRotator->SetTrackArrays(&fTracks[0],&fTracks[1]);
if (fDebugTree) fDebugTree->SetDielectron(this);
}
}
}
+ //track rotation
+ if (fTrackRotator) FillPairArrayTR();
+
//in case there is a histogram manager, fill the QA histograms
if (fHistos) FillHistograms(ev1);
}
}
+//________________________________________________________________
+void AliDielectron::FillHistogramsPair(AliDielectronPair *pair)
+{
+ //
+ // Fill Histogram information for pairs and the track in the pair
+ // NOTE: in this funtion the leg information may be filled multiple
+ // times. This funtion is used in the track rotation pairing
+ // and those legs are not saved!
+ //
+ TString className,className2;
+ Double_t values[AliDielectronVarManager::kNMaxValues];
+
+ //Fill Pair information, separately for all pair candidate arrays and the legs
+ TObjArray arrLegs(100);
+ const Int_t type=pair->GetType();
+ className.Form("Pair_%s",fgkPairClassNames[type]);
+ className2.Form("Track_Legs_%s",fgkPairClassNames[type]);
+
+ Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
+ Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
+
+ //fill pair information
+ if (pairClass){
+ AliDielectronVarManager::Fill(pair, values);
+ fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
+ }
+
+ if (legClass){
+ AliVParticle *d1=pair->GetFirstDaughter();
+ AliDielectronVarManager::Fill(d1, values);
+ fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
+
+ AliVParticle *d2=pair->GetSecondDaughter();
+ AliDielectronVarManager::Fill(d2, values);
+ fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
+ }
+}
//________________________________________________________________
void AliDielectron::FillTrackArrays(AliVEvent * const ev, Int_t eventNr)
}
//________________________________________________________________
-void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2) {
+void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2)
+{
//
// Prefilter tracks from pairs
// Needed for datlitz rejections
}
//________________________________________________________________
-void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2) {
+void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2)
+{
//
// select pairs and fill pair candidate arrays
//
delete candidate;
}
+//________________________________________________________________
+void AliDielectron::FillPairArrayTR()
+{
+ //
+ // select pairs and fill pair candidate arrays
+ //
+ UInt_t selectedMask=(1<<fPairFilter.GetCuts()->GetEntries())-1;
+
+ while ( fTrackRotator->NextCombination() ){
+ AliDielectronPair candidate;
+ candidate.SetTracks(fTrackRotator->GetTrackP(), fPdgLeg1, fTrackRotator->GetTrackN(), fPdgLeg2);
+ candidate.SetType(10);
+
+ //pair cuts
+ UInt_t cutMask=fPairFilter.IsSelected(&candidate);
+
+ //CF manager for the pair
+ if (fCfManagerPair) fCfManagerPair->Fill(cutMask,&candidate);
+
+ //apply cut
+ if (cutMask==selectedMask&&fHistos) FillHistogramsPair(&candidate);
+ }
+}
+
//________________________________________________________________
void AliDielectron::FillDebugTree()
{
class THashList;
class AliDielectronCF;
class AliDielectronDebugTree;
+class AliDielectronTrackRotator;
//________________________________________________________________
class AliDielectron : public TNamed {
public:
- enum ParticleValues { kPx=0, kPy, kPz, kPt, kP, kXv, kYv, kZv, kOneOverPt,
- kPhi, kTheta, kE, kM, kEta, kY, kCharge, kNParticleValues };
- enum PairValues { kChi2NDF=kNParticleValues, kDecayLength, kR, kOpeningAngle, kMerr, kNPairValues };
+ enum EPairType { kEv1PP=0, kEv1PM, kEv1MM,
+ kEv1PEv2P, kEv1MEv2P, kEv2PP,
+ kEv1PEv2M, kEv1MEv2M, kEv2PM,
+ kEv2MM, kEv1PMRot };
+ enum ELegType { kEv1P, kEv1M, kEv2P, kEv2M };
AliDielectron();
AliDielectron(const char* name, const char* title);
void SetCFManagerPair(AliDielectronCF * const cf) { fCfManagerPair=cf; }
AliDielectronCF* GetCFManagerPair() const { return fCfManagerPair; }
+ void SetTrackRotator(AliDielectronTrackRotator * const rot) { fTrackRotator=rot; }
+ AliDielectronTrackRotator* GetTrackRotator() const { return fTrackRotator; }
+
void SetDebugTree(AliDielectronDebugTree * const tree) { fDebugTree=tree; }
static const char* TrackClassName(Int_t i) { return (i>=0&&i<4)?fgkTrackClassNames[i]:""; }
- static const char* PairClassName(Int_t i) { return (i>=0&&i<10)?fgkPairClassNames[i]:""; }
+ static const char* PairClassName(Int_t i) { return (i>=0&&i<11)?fgkPairClassNames[i]:""; }
void SaveDebugTree();
//TODO: better way to store it? TClonesArray?
AliDielectronCF *fCfManagerPair;//Correction Framework Manager for the Pair
-
+ AliDielectronTrackRotator *fTrackRotator; //Track rotator
AliDielectronDebugTree *fDebugTree; // Debug tree output
void FillTrackArrays(AliVEvent * const ev, Int_t eventNr=0);
void PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2);
void FillPairArrays(Int_t arr1, Int_t arr2);
+ void FillPairArrayTR();
Int_t GetPairIndex(Int_t arr1, Int_t arr2) const {return arr1>=arr2?arr1*(arr1+1)/2+arr2:arr2*(arr2+1)/2+arr1;}
TObjArray* PairArray(Int_t i);
static const char* fgkTrackClassNames[4]; //Names for track arrays
- static const char* fgkPairClassNames[10]; //Names for pair arrays
+ static const char* fgkPairClassNames[11]; //Names for pair arrays
void ProcessMC();
void FillHistograms(const AliVEvent *ev);
+ void FillHistogramsPair(AliDielectronPair *pair);
+
void FillDebugTree();
AliDielectron(const AliDielectron &c);
TH1* AliDielectronCFdraw::Project(Int_t ndim, Int_t *vars, Int_t slice)
{
//
- // Do an nim projection
+ // Do an ndim projection
//
switch (ndim){
case 1:
return 0x0;
}
+//________________________________________________________________
+TH1* AliDielectronCFdraw::Project(const Option_t* var, Int_t slice)
+{
+ //
+ // translate variable names and do projection
+ //
+ TObjArray *arrVars=TString(var).Tokenize(":");
+ Int_t entries=arrVars->GetEntriesFast();
+ if (entries<1||entries>3){
+ AliError("Wrong number of variables, supported are 1 - 3 dimensions");
+ delete arrVars;
+ return 0x0;
+ }
+
+ TIter next(arrVars);
+ TObjString *ostr=0x0;
+ Int_t ivar[3]={-1,-1,-1};
+ for (Int_t i=0; i<entries; ++i){
+ ostr=static_cast<TObjString*>(next());
+ if (ostr->GetString().IsDigit()){
+ ivar[i]=ostr->GetString().Atoi();
+ } else {
+ ivar[i]=fCfContainer->GetVar(ostr->GetName());
+ }
+ }
+ delete arrVars;
+ return Project(entries,ivar,slice);
+}
+
//________________________________________________________________
void AliDielectronCFdraw::DrawEfficiency(const char* varnames, const char* nominators, Int_t denominator, const char* opt)
{
TObjArray* CollectHistosProj(Int_t dim, Int_t *vars, const char* slices);
TH1* Project(Int_t ndim, Int_t *vars, Int_t slice);
-
+ TH1* Project(const Option_t* var, Int_t slice);
+
//Draw efficiencies
void DrawEfficiency(const char* varnames, const char* nominators, Int_t denominator=0, const char* opt="sameleg2");
void DrawEfficiency(Int_t var, const char* nominators, Int_t denominator=0, const char* opt="sameleg", Int_t type=0);
--- /dev/null
+/*************************************************************************
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+//////////////////////////////////////////////////////////////////////////
+// CutGroup //
+// //
+// //
+// Allow to define groups of cut conditions which are tested with //
+// an OR condition between groups and an AND within groups //
+// //
+//////////////////////////////////////////////////////////////////////////
+
+#include "AliDielectronCutGroup.h"
+
+ClassImp(AliDielectronCutGroup)
+
+AliDielectronCutGroup::AliDielectronCutGroup(Bool_t compOperator /*=kCompOR*/) :
+ AliAnalysisCuts(),
+ fCutGroupList(0x0),
+ fCompOperator(compOperator)
+{
+ //
+ // Default constructor
+ //
+}
+
+//_____________________________________________________________________
+AliDielectronCutGroup::AliDielectronCutGroup(const char* name, const char* title, Bool_t compOperator /*=kCompOR*/) :
+ AliAnalysisCuts(name, title),
+ fCutGroupList(0x0),
+ fCompOperator(compOperator)
+{
+ //
+ // Named Constructor
+ //
+}
+
+//_____________________________________________________________________
+AliDielectronCutGroup::~AliDielectronCutGroup()
+{
+ //
+ //Default Destructor
+ //
+}
+
+//_____________________________________________________________________
+Bool_t AliDielectronCutGroup::IsSelected(TObject* track)
+{
+ //
+ // Selection-finder handling different comparison operations
+ //
+
+
+ //Different init for and/or makes code shorter
+ Bool_t selectionResult=fCompOperator;
+
+ TIter listIterator(&fCutGroupList);
+ while (AliAnalysisCuts *thisCut = (AliAnalysisCuts*) listIterator()) {
+ if (fCompOperator == kCompOR) {
+ selectionResult = (selectionResult || thisCut->IsSelected(track));
+ }
+ else { //kCompAND
+ selectionResult = (selectionResult && thisCut->IsSelected(track));
+ if (selectionResult==kFALSE) break; //Save loops vs. additional check?
+ }
+
+ }
+ return selectionResult;
+}
+
+//_____________________________________________________________________
+
+void AliDielectronCutGroup::AddCut(AliAnalysisCuts* fCut)
+{
+ //
+ // Add a defined cut to the list
+ //
+
+ fCutGroupList.Add(fCut);
+}
+
+//_____________________________________________________________________
+void AliDielectronCutGroup::SetCompOperator(Bool_t compOperator)
+{
+ //
+ // Switch between AND/OR
+ //
+
+ fCompOperator = compOperator;
+}
--- /dev/null
+#ifndef ALIDIELECTRONCUTGROUP_H
+#define ALIDIELECTRONCUTGROUP_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#################################################################
+//# #
+//# Class AliDielectronCutGroup #
+//# Dielectron Group of cuts #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Christoph Baumann uni Ffm / cbaumann@ikf.uni-frankfurt.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#################################################################
+
+#include <AliAnalysisCuts.h>
+#include <TList.h>
+
+class TCollection;
+
+class AliDielectronCutGroup : public AliAnalysisCuts {
+
+public:
+ enum TruthValues {
+ kCompAND = kTRUE,
+ kCompOR = kFALSE
+ };
+
+ AliDielectronCutGroup(Bool_t compOperator=kCompOR);
+ AliDielectronCutGroup(const char*name, const char* title, Bool_t compOperator=kCompOR);
+
+ virtual ~AliDielectronCutGroup();
+
+ //Analysis cuts interface
+ //
+ virtual Bool_t IsSelected(TObject* track);
+ virtual Bool_t IsSelected(TList* /* list */ ) {return kFALSE;}
+
+ void AddCut(AliAnalysisCuts* fCut);
+ void SetCompOperator(Bool_t compOperator);
+
+private:
+ TList fCutGroupList; //for saving the different cuts
+ Bool_t fCompOperator; //determines whether the cuts are AND/OR compared
+
+ ClassDef(AliDielectronCutGroup,1) //Group of cuts
+};
+
+#endif
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Dielectron EventCuts //
+// //
+// //
+/*
+Detailed description
+
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <AliTriggerAnalysis.h>
+#include <AliESDVertex.h>
+#include <AliESDEvent.h>
+
+#include "AliDielectronEventCuts.h"
+
+ClassImp(AliDielectronEventCuts)
+
+AliDielectronEventCuts::AliDielectronEventCuts() :
+ AliAnalysisCuts(),
+ fVtxZmin(0.),
+ fVtxZmax(0.),
+ fRequireVtx(kFALSE),
+ fMinVtxContributors(0),
+ fVtxType(kVtxTracks),
+ fRequireV0and(0),
+ fTriggerAnalysis(0x0),
+ fkVertex(0x0)
+{
+ //
+ // Default Constructor
+ //
+
+}
+
+//______________________________________________
+AliDielectronEventCuts::AliDielectronEventCuts(const char* name, const char* title) :
+ AliAnalysisCuts(name, title),
+ fVtxZmin(0.),
+ fVtxZmax(0.),
+ fRequireVtx(kFALSE),
+ fMinVtxContributors(0),
+ fVtxType(kVtxTracks),
+ fRequireV0and(0),
+ fTriggerAnalysis(0x0),
+ fkVertex(0x0)
+{
+ //
+ // Named Constructor
+ //
+}
+
+//______________________________________________
+AliDielectronEventCuts::~AliDielectronEventCuts()
+{
+ //
+ // Default Destructor
+ //
+ if (fTriggerAnalysis) delete fTriggerAnalysis;
+}
+
+//______________________________________________
+Bool_t AliDielectronEventCuts::IsSelected(TObject* event)
+{
+ //
+ // check the cuts
+ //
+
+ AliESDEvent *ev=dynamic_cast<AliESDEvent*>(event);
+ if (!ev) return kFALSE;
+
+ fkVertex=0x0;
+ switch(fVtxType){
+ case kVtxTracks: fkVertex=ev->GetPrimaryVertexTracks(); break;
+ case kVtxSPD: fkVertex=ev->GetPrimaryVertexSPD(); break;
+ case kVtxTPC: fkVertex=ev->GetPrimaryVertexTPC(); break;
+ case kVtxAny: fkVertex=ev->GetPrimaryVertex(); break;
+ }
+
+ if ((fRequireVtx||fVtxZmin<fVtxZmax||fMinVtxContributors>0)&&!fkVertex) return kFALSE;
+
+ if (fVtxZmin<fVtxZmax){
+ Double_t zvtx=fkVertex->GetZv();
+ if (zvtx<fVtxZmin||zvtx>fVtxZmax) return kFALSE;
+ }
+
+ if (fMinVtxContributors>0){
+ Int_t nCtrb = fkVertex->GetNContributors();
+ if (nCtrb<fMinVtxContributors) return kFALSE;
+ }
+
+ if (fRequireV0and){
+ if (!fTriggerAnalysis) fTriggerAnalysis=new AliTriggerAnalysis;
+ Bool_t v0AND = kFALSE;
+ if (fRequireV0and==1){
+ Bool_t v0A = fTriggerAnalysis->IsOfflineTriggerFired(ev, AliTriggerAnalysis::kV0A);
+ Bool_t v0C = fTriggerAnalysis->IsOfflineTriggerFired(ev, AliTriggerAnalysis::kV0C);
+ v0AND = v0A && v0C;
+ }
+
+ if (fRequireV0and==2){
+ Bool_t v0AHW = (fTriggerAnalysis->V0Trigger(ev, AliTriggerAnalysis::kASide, kTRUE) == AliTriggerAnalysis::kV0BB);
+ Bool_t v0CHW = (fTriggerAnalysis->V0Trigger(ev, AliTriggerAnalysis::kCSide, kTRUE) == AliTriggerAnalysis::kV0BB);
+ v0AND = v0AHW && v0CHW;
+ }
+
+ if (!v0AND) return kFALSE;
+ }
+
+ return kTRUE;
+}
+
--- /dev/null
+#ifndef ALIDIELECTRONEVENTCUTS_H
+#define ALIDIELECTRONEVENTCUTS_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronEventCuts #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <AliAnalysisCuts.h>
+
+class AliTriggerAnalysis;
+class AliESDVertex;
+
+class AliDielectronEventCuts : public AliAnalysisCuts {
+public:
+ enum EVtxType { kVtxTracks=0, kVtxSPD, kVtxTPC, kVtxAny };
+
+ AliDielectronEventCuts();
+ AliDielectronEventCuts(const char*name, const char* title);
+
+ virtual ~AliDielectronEventCuts();
+
+
+ void SetVertexType(EVtxType type) { fVtxType=type; }
+ void SetVertexZ(Double_t zmin, Double_t zmax) { fVtxZmin=zmin; fVtxZmax=zmax; }
+ void SetRequireVertex(Bool_t req=kTRUE) { fRequireVtx=req; }
+ void SetRequireV0and(UChar_t type=1) { fRequireV0and=type; }
+ void SetMinVtxContributors(Int_t min=1) { fMinVtxContributors=min; }
+
+ //
+ //Analysis cuts interface
+ //
+ virtual Bool_t IsSelected(TObject* event);
+ virtual Bool_t IsSelected(TList* /* list */ ) {return kFALSE;}
+
+
+private:
+
+ Double_t fVtxZmin; // minimum z vertex position
+ Double_t fVtxZmax; // maximum z vertex position
+ Bool_t fRequireVtx; // require a vertex
+ Int_t fMinVtxContributors; // min number of vertex contributors
+ EVtxType fVtxType; // vertex type
+
+ UChar_t fRequireV0and; // use V0and triggered events only
+
+ AliTriggerAnalysis *fTriggerAnalysis; //! trigger analysis class
+ const AliESDVertex *fkVertex; //! current vertex
+
+ AliDielectronEventCuts(const AliDielectronEventCuts &c);
+ AliDielectronEventCuts &operator=(const AliDielectronEventCuts &c);
+
+
+ ClassDef(AliDielectronEventCuts,1) // Dielectron EventCuts
+};
+
+
+
+#endif
return fStack->GetNtrack();
}
+//____________________________________________________________
+Int_t AliDielectronMC::GetNPrimary()
+{
+ //
+ // return the number of primary track from MC event
+ //
+ if (!fMCEvent){ AliError("No fMCEvent"); return 0; }
+ return fMCEvent->GetNumberOfPrimaries();
+}
+
+//____________________________________________________________
+Int_t AliDielectronMC::GetNPrimaryFromStack()
+{
+ //
+ // return the number of primary track from stack
+ //
+ if (!fStack){ AliError("No fStack"); return -999; }
+ return fStack->GetNprimary();
+}
+
//____________________________________________________________
AliVParticle* AliDielectronMC::GetMCTrackFromMCEvent(Int_t itrk)
{
void Initialize(); // initialization
Int_t GetNMCTracks(); // return number of generated tracks
Int_t GetNMCTracksFromStack(); // return number of generated tracks from stack
+ Int_t GetNPrimary(); // return number of primary tracks
+ Int_t GetNPrimaryFromStack(); // return number of primary tracks from stack
Int_t GetMCPID(const AliESDtrack* _track); // return MC PID
Int_t GetMCPIDFromStack(const AliESDtrack* _track); // return MC PID
Int_t GetMotherPDG(const AliESDtrack* _track); // return mother PID from the MC stack
#include <AliVTrack.h>
#include <AliLog.h>
#include <AliESDtrack.h>
+#include <AliESDpid.h>
+#include <AliAODpidUtil.h>
#include "AliDielectronVarManager.h"
AliDielectronPID::AliDielectronPID() :
AliAnalysisCuts(),
fNcuts(0),
- fESDpid(0x0)
+ fESDpid(0x0),
+ fAODpidUtil(0x0)
{
//
// Default Constructor
AliDielectronPID::AliDielectronPID(const char* name, const char* title) :
AliAnalysisCuts(name, title),
fNcuts(0),
- fESDpid(0x0)
+ fESDpid(0x0),
+ fAODpidUtil(0x0)
{
//
// Named Constructor
Bool_t selected=kFALSE;
fESDpid=AliDielectronVarManager::GetESDpid();
+ fAODpidUtil=AliDielectronVarManager::GetAODpidUtil();
for (UChar_t icut=0; icut<fNcuts; ++icut){
Double_t pMin=fPmin[icut];
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
numberOfSigmas=fESDpid->NumberOfSigmasITS(track, fPartType[icut]);
- }else{
+ }else if(part->IsA()==AliAODTrack::Class()){
// AOD case
// FIXME: Is there a place to check whether the PID is was set in ESD???
AliAODTrack *track=static_cast<AliAODTrack*>(part);
- numberOfSigmas=NumberOfSigmasITS(track, fPartType[icut]);
+ numberOfSigmas=fAODpidUtil->NumberOfSigmasITS(track, fPartType[icut]);
}
Bool_t selected=((numberOfSigmas>=fNsigmaLow[icut])&&(numberOfSigmas<=fNsigmaUp[icut]))^fExclude[icut];
return selected;
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
numberOfSigmas=fESDpid->NumberOfSigmasTPC(track, fPartType[icut]);
- }else{
+ }else if(part->IsA()==AliAODTrack::Class()){
// AOD case
// FIXME: Is there a place to check whether the PID is was set in ESD???
AliAODTrack *track=static_cast<AliAODTrack*>(part);
- numberOfSigmas=NumberOfSigmasTPC(track, fPartType[icut]);
+ numberOfSigmas=fAODpidUtil->NumberOfSigmasTPC(track, fPartType[icut]);
+ }
+ if (fPartType[icut]==AliPID::kElectron){
+ numberOfSigmas-=fgCorr;
}
-// if (fPartType[icut]==AliPID::kElectron){
-// numberOfSigmas-=fgCorr;
-// }
Bool_t selected=((numberOfSigmas>=fNsigmaLow[icut])&&(numberOfSigmas<=fNsigmaUp[icut]))^fExclude[icut];
return selected;
}
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
numberOfSigmas=fESDpid->NumberOfSigmasTOF(track, fPartType[icut], fESDpid->GetTOFResponse().GetTimeZero());
- }else{
+ }else if(part->IsA()==AliAODTrack::Class()){
// AOD case
// FIXME: Is there a place to check whether the PID is was set in ESD???
AliAODTrack *track=static_cast<AliAODTrack*>(part);
- numberOfSigmas=NumberOfSigmasTOF(track, fPartType[icut]);
+ numberOfSigmas=fAODpidUtil->NumberOfSigmasTOF(track, fPartType[icut]);
}
Bool_t selected=((numberOfSigmas>=fNsigmaLow[icut])&&(numberOfSigmas<=fNsigmaUp[icut]))^fExclude[icut];
//#############################################################
#include <AliPID.h>
-#include <AliESDpid.h>
#include <AliAODTrack.h>
#include <AliAODPid.h>
class TList;
class AliVTrack;
class TGraph;
+class AliESDpid;
+class AliAODpidUtil;
class AliDielectronPID : public AliAnalysisCuts {
public:
UChar_t fRequirePIDbit[kNmaxPID]; //How to make use of the pid bit (see)
AliESDpid *fESDpid; //! esd pid object
-
-
+ AliAODpidUtil *fAODpidUtil; //! AOD pid object
+
static TGraph *fgFitCorr; //spline fit object to correct the nsigma deviation in the TPC electron band
static Double_t fgCorr; //!correction value for current run. Set if fgFitCorr is set and SetCorrVal(run)
// was called
Bool_t IsSelectedTRD(AliVTrack * const part, Int_t icut) const;
Bool_t IsSelectedTOF(AliVTrack * const part, Int_t icut) const;
- Float_t NumberOfSigmasITS(const AliAODTrack *track, AliPID::EParticleType type) const;
- Float_t NumberOfSigmasTPC(const AliAODTrack *track, AliPID::EParticleType type) const;
- Float_t NumberOfSigmasTOF(const AliAODTrack *track, AliPID::EParticleType type) const;
-
AliDielectronPID(const AliDielectronPID &c);
AliDielectronPID &operator=(const AliDielectronPID &c);
ClassDef(AliDielectronPID,3) // Dielectron PID
};
-
-//
-// Inline functions for AOD as long as ther is no AOD pid object we have to fake it
-//
-
-inline Float_t AliDielectronPID::NumberOfSigmasITS(const AliAODTrack *track, AliPID::EParticleType type) const {
- AliAODPid *pid=track->GetDetPid();
- if (!pid) return -1000.;
-
- return fESDpid->GetITSResponse().GetNumberOfSigmas(track->P(),pid->GetITSsignal(),type);
-}
-
-inline Float_t AliDielectronPID::NumberOfSigmasTPC(const AliAODTrack *track, AliPID::EParticleType type) const {
- AliAODPid *pid=track->GetDetPid();
- if (!pid) return -1000.;
-
- Double_t mom = pid->GetTPCmomentum();
- if (mom<0) mom=track->P();
-
- //FIXME: rough estimate of the number of clusters used for PID. Needs to be fixed!!!
- Int_t ncl=(Int_t)track->GetTPCClusterMap().CountBits();
- return fESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),ncl,type);
-}
-
-inline Float_t AliDielectronPID::NumberOfSigmasTOF(const AliAODTrack *track, AliPID::EParticleType type) const {
- AliAODPid *pid=track->GetDetPid();
- if (!pid) return -1000.;
-
- Double_t times[AliPID::kSPECIES];
- pid->GetIntegratedTimes(times);
- Double_t tofRes = fESDpid->GetTOFResponse().GetExpectedSigma(track->P(),times[type],AliPID::ParticleMass(type));
- return (pid->GetTOFsignal() - times[type])/ tofRes;
-}
-
-
-
#endif
printf("Mass res.: %.5g #pm %.5g\n", fValues(5), fErrors(5));
}
}
+
+//______________________________________________
+void AliDielectronSignalBase::ScaleHistograms(TH1* histRaw, TH1* histBackground, Double_t intMin, Double_t intMax)
+{
+ //
+ // scale histBackground to match the integral of histRaw in the interval intMin, intMax
+ //
+
+ Double_t intRaw = histRaw->Integral(histRaw->FindBin(intMin),histRaw->FindBin(intMax));
+ Double_t intBack = histBackground->Integral(histBackground->FindBin(intMin),histBackground->FindBin(intMax));
+ if (intBack>0){
+ histBackground->Sumw2();
+ histBackground->Scale(1./intBack);
+ histBackground->Scale(intRaw);
+ }
+}
enum EBackgroundMethod {
kFitted = 0,
kLikeSign,
- kEventMixing
+ kEventMixing,
+ kRotation
};
AliDielectronSignalBase();
Double_t GetMassWidth() const { return fValues(5);}
Double_t GetMassWidthError() const { return fErrors(5);}
- TH1F* GetSignalHistogram() const {return fHistSignal;}
- TH1F* GetBackgroundHistogram() const {return fHistBackground;}
-
+ TH1* GetSignalHistogram() const {return fHistSignal;}
+ TH1* GetBackgroundHistogram() const {return fHistBackground;}
+ TH1* GetUnlikeSignHistogram() const {return fHistDataPM;}
+
+ static void ScaleHistograms(TH1* histRaw, TH1* histBackground, Double_t intMin, Double_t intMax);
+
virtual void Print(Option_t *option="") const;
/**
protected:
- TH1F *fHistSignal; // histogram of pure signal
- TH1F *fHistBackground; // histogram of background (fitted=0, like-sign=1, event mixing=2)
- TH1F *fHistDataPM; // histogram of selected +- pair candidates
- TH1F *fHistDataPP; // histogram of selected ++ pair candidates
- TH1F *fHistDataMM; // histogram of selected -- pair candidates
+ TH1 *fHistSignal; // histogram of pure signal
+ TH1 *fHistBackground; // histogram of background (fitted=0, like-sign=1, event mixing=2)
+ TH1 *fHistDataPM; // histogram of selected +- pair candidates
+ TH1 *fHistDataPP; // histogram of selected ++ pair candidates
+ TH1 *fHistDataMM; // histogram of selected -- pair candidates
TVectorD fValues; // values
TVectorD fErrors; // value errors
ProcessEM(arrhist); // process event mixing method
break;
+ case kRotation:
+ ProcessRotation(arrhist);
+ break;
+
default :
AliWarning("Subtraction method not supported. Please check SetMethod() function.");
}
//
// signal subtraction
//
- fHistDataPP = (TH1F*)(arrhist->At(0))->Clone("histPP"); // ++ SE
- fHistDataPM = (TH1F*)(arrhist->At(1))->Clone("histPM"); // +- SE
- fHistDataMM = (TH1F*)(arrhist->At(2))->Clone("histMM"); // -- SE
+ fHistDataPP = (TH1*)(arrhist->At(0))->Clone("histPP"); // ++ SE
+ fHistDataPM = (TH1*)(arrhist->At(1))->Clone("histPM"); // +- SE
+ fHistDataMM = (TH1*)(arrhist->At(2))->Clone("histMM"); // -- SE
fHistDataPP->Sumw2();
fHistDataPM->Sumw2();
fHistDataMM->Sumw2();
fHistDataMM->Rebin(fRebin);
}
- fHistSignal = new TH1F("HistSignal", "Like-Sign substracted signal",
+ fHistSignal = new TH1D("HistSignal", "Like-Sign substracted signal",
fHistDataPM->GetXaxis()->GetNbins(),
fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
- fHistBackground = new TH1F("HistBackground", "Like-sign contribution",
+ fHistBackground = new TH1D("HistBackground", "Like-sign contribution",
fHistDataPM->GetXaxis()->GetNbins(),
fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
ProcessLS(arrhist);
}
+//______________________________________________
+void AliDielectronSignalExt::ProcessRotation(TObjArray* const arrhist)
+{
+ //
+ // signal subtraction
+ //
+ fHistDataPM = (TH1*)(arrhist->At(1))->Clone("histPM"); // +- SE
+ if (!fHistDataPM){
+ AliError("Unlike sign histogram not available. Cannot extract the signal.");
+ return;
+ }
+ fHistDataPM->Sumw2();
+
+ fHistBackground=(TH1*)arrhist->At(10)->Clone("histRotation");
+ if (!fHistBackground){
+ AliError("Histgram from rotation not available. Cannot extract the signal.");
+ delete fHistDataPM;
+ fHistDataPM=0x0;
+ return;
+ }
+
+ //scale histograms to match integral between 3.2 and 4. GeV
+ ScaleHistograms(fHistDataPM,fHistBackground,3.2,4.0);
+ fHistSignal=(TH1*)fHistDataPM->Clone("histSignal");
+ fHistSignal->Add(fHistBackground,-1.);
+
+ // signal
+ fValues(0) = fHistSignal->IntegralAndError(fHistSignal->FindBin(fIntMin),
+ fHistSignal->FindBin(fIntMax), fErrors(0));
+ // background
+ fValues(1) = fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
+ fHistBackground->FindBin(fIntMax),
+ fErrors(1));
+ // S/B and significance
+ SetSignificanceAndSOB();
+
+ fProcessed = kTRUE;
+
+}
+
//______________________________________________
void AliDielectronSignalExt::Draw(const Option_t* option)
{
virtual void Process(TObjArray* const arrhist);
void ProcessLS(TObjArray* const arrhist); // like-sign method
void ProcessEM(TObjArray* const arrhist); // event mixing method
-
+ void ProcessRotation(TObjArray* const arrhist); // event mixing method
+
virtual void Draw(const Option_t* option = "");
private:
#include <TPaveText.h>
#include <TList.h>
#include <TFitResult.h>
-//#include <../hist/hist/src/TF1Helper.h> //not supposed to be used!
+#include <../hist/hist/src/TF1Helper.h>
#include <AliLog.h>
// by the user in its macro
fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
TFitResultPtr pmFitPtr = fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
- //TFitResult *pmFitResult = pmFitPtr.Get(); //not used when TF1Helper out
+ TFitResult *pmFitResult = pmFitPtr.Get();
fFuncSignal->SetParameters(fFuncSigBack->GetParameters());
fFuncBackground->SetParameters(fFuncSigBack->GetParameters()+fFuncSignal->GetNpar());
Double_t epm = fHistDataPM->GetBinError(iBin);
Double_t bknd = fFuncBackground->Eval(m);
Double_t ebknd = 0;
- /* to be revised ... TF1Helper
for(Int_t iPar=fFuncSignal->GetNpar(); iPar<fFuncSigBack->GetNpar(); iPar++) {
for(Int_t jPar=iPar; jPar<fFuncSigBack->GetNpar(); jPar++) {
TF1 gradientIpar("gradientIpar",
(iPar==jPar ? 1.0 : 2.0);
}
}
- */ // TF1Helper
Double_t signal = pm-bknd;
Double_t error = TMath::Sqrt(epm*epm+ebknd);
fHistSignal->SetBinContent(iBin, signal);
// signal
fValues(0) = fFuncSignal->Integral(fIntMin, fIntMax)/fHistDataPM->GetBinWidth(1);
fErrors(0) = 0;
- /* to be revised ... TF1Helper
for(Int_t iPar=0; iPar<fFuncSignal->GetNpar(); iPar++) {
for(Int_t jPar=iPar; jPar<fFuncSignal->GetNpar(); jPar++) {
TF1 gradientIpar("gradientIpar",
(iPar==jPar ? 1.0 : 2.0);
}
}
- */ //TF1Helper
// background
fValues(1) = fFuncBackground->Integral(fIntMin, fIntMax)/fHistDataPM->GetBinWidth(1);
fErrors(1) = 0;
- /* to be revised... TF1Helper
for(Int_t iPar=fFuncSignal->GetNpar(); iPar<fFuncSigBack->GetNpar(); iPar++) {
for(Int_t jPar=iPar; jPar<fFuncSigBack->GetNpar(); jPar++) {
TF1 gradientIpar("gradientIpar",
(iPar==jPar ? 1.0 : 2.0);
}
}
- */ // TF1Helper
}
else {
// signal
TFitResultPtr mmFitPtr = fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
mmFitResult = mmFitPtr.Get();
- /* to be revised ... TF1Helper
for(Int_t iBin=1; iBin<=fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
Double_t m = fHistDataPM->GetBinCenter(iBin);
Double_t pm = fHistDataPM->GetBinContent(iBin);
(iPar==jPar ? 1.0 : 2.0);
}
}
-
Double_t emm = 0;
for(Int_t iPar=0; iPar<funcCloneMM->GetNpar(); iPar++) {
for(Int_t jPar=iPar; jPar<funcCloneMM->GetNpar(); jPar++) {
fHistBackground->SetBinContent(iBin, background);
fHistBackground->SetBinError(iBin, ebackground);
}
- */ // TF1Helper
-
+
// signal
fValues(0) = fHistSignal->IntegralAndError(fHistSignal->FindBin(fIntMin),
fHistSignal->FindBin(fIntMax), fErrors(0));
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Dielectron TrackRotator //
+// //
+// //
+/*
+Detailed description
+
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <TMath.h>
+#include <TObjArray.h>
+#include <AliAODTrack.h>
+#include <AliESDtrack.h>
+#include <TRandom3.h>
+
+#include "AliDielectronTrackRotator.h"
+
+ClassImp(AliDielectronTrackRotator)
+
+AliDielectronTrackRotator::AliDielectronTrackRotator() :
+ TNamed(),
+ fIterations(1),
+ fRotationType(kRotateBothRandom),
+ fStartAnglePhi(TMath::Pi()),
+ fConeAnglePhi(TMath::Pi()/6.),
+ fkArrTracksP(0x0),
+ fkArrTracksN(0x0),
+ fCurrentIteration(0),
+ fCurrentTackP(0),
+ fCurrentTackN(0),
+ fTrackP(0x0),
+ fTrackN(0x0)
+{
+ //
+ // Default Constructor
+ //
+ gRandom->SetSeed();
+}
+
+//______________________________________________
+AliDielectronTrackRotator::AliDielectronTrackRotator(const char* name, const char* title) :
+ TNamed(name, title),
+ fIterations(1),
+ fRotationType(kRotateBothRandom),
+ fStartAnglePhi(TMath::Pi()),
+ fConeAnglePhi(TMath::Pi()/6.),
+ fkArrTracksP(0x0),
+ fkArrTracksN(0x0),
+ fCurrentIteration(0),
+ fCurrentTackP(0),
+ fCurrentTackN(0),
+ fTrackP(0x0),
+ fTrackN(0x0)
+{
+ //
+ // Named Constructor
+ //
+ gRandom->SetSeed();
+}
+
+//______________________________________________
+AliDielectronTrackRotator::~AliDielectronTrackRotator()
+{
+ //
+ // Default Destructor
+ //
+
+}
+
+//______________________________________________
+void AliDielectronTrackRotator::Reset()
+{
+ //
+ // Reset the current iterators
+ //
+ fCurrentIteration=0;
+ fCurrentTackP=0;
+ fCurrentTackN=0;
+}
+
+//______________________________________________
+Bool_t AliDielectronTrackRotator::NextCombination()
+{
+ //
+ // Perform track rotation of the tracks in the track arrays as long as there are possible combinations
+ //
+ if (!fkArrTracksP || !fkArrTracksP) {
+ Reset();
+ return kFALSE;
+ }
+
+ Int_t nP=fkArrTracksP->GetEntriesFast();
+ Int_t nN=fkArrTracksN->GetEntriesFast();
+ if (nP==0||nN==0){
+ Reset();
+ return kFALSE;
+ }
+
+ if (fCurrentIteration==fIterations){
+ fCurrentIteration=0;
+ ++fCurrentTackP;
+ }
+
+ if (fCurrentTackP==nP){
+ ++fCurrentTackN;
+ fCurrentTackP=0;
+ }
+
+ if (fCurrentTackN==nN){
+ Reset();
+ return kFALSE;
+ }
+
+ if (!RotateTracks()){
+ Reset();
+ return kFALSE;
+ }
+
+ ++fCurrentIteration;
+ return kTRUE;
+}
+
+//______________________________________________
+Bool_t AliDielectronTrackRotator::RotateTracks()
+{
+ //
+ // Actual track rotation
+ // Find out particle type and perform the rotation
+ //
+
+ const AliVTrack *trackP=dynamic_cast<AliVTrack*>(fkArrTracksP->UncheckedAt(fCurrentTackP));
+ const AliVTrack *trackN=dynamic_cast<AliVTrack*>(fkArrTracksP->UncheckedAt(fCurrentTackN));
+ if (!trackP||!trackN) return kFALSE;
+
+
+ Double_t angle = fStartAnglePhi+(2*gRandom->Rndm()-1)*fConeAnglePhi;
+ Int_t charge = TMath::Nint(gRandom->Rndm());
+
+ if (trackP->IsA()==AliESDtrack::Class()) {
+
+ if (!fTrackP) {
+ fTrackP=new AliESDtrack;
+ fTrackN=new AliESDtrack;
+ }
+
+ trackP->Copy(*fTrackP);
+ trackN->Copy(*fTrackN);
+
+ if (fRotationType==kRotatePositive||(fRotationType==kRotateBothRandom&&charge==0)){
+ ((AliESDtrack*)fTrackP)->Rotate(angle);
+ }
+
+ if (fRotationType==kRotateNegative||(fRotationType==kRotateBothRandom&&charge==1)){
+ ((AliESDtrack*)fTrackN)->Rotate(angle);
+ }
+
+ } else if (trackP->IsA()==AliAODTrack::Class()) {
+
+ if (!fTrackP) {
+ fTrackP=new AliAODTrack;
+ fTrackN=new AliAODTrack;
+ }
+
+ (*(AliAODTrack*)fTrackP)=(*(AliAODTrack*)trackP);
+ (*(AliAODTrack*)fTrackN)=(*(AliAODTrack*)trackN);
+
+ if (fRotationType==kRotatePositive||(fRotationType==kRotateBothRandom&&charge==0)){
+ Double_t phi=fTrackP->Phi()+angle;
+ if (phi>2*TMath::Pi()) phi-=2*TMath::Pi();
+ ((AliAODTrack*)fTrackP)->SetPhi(phi);
+ }
+
+ if (fRotationType==kRotateNegative||(fRotationType==kRotateBothRandom&&charge==1)){
+ Double_t phi=fTrackN->Phi()+angle;
+ if (phi>2*TMath::Pi()) phi-=2*TMath::Pi();
+ ((AliAODTrack*)fTrackN)->SetPhi(phi);
+ }
+
+ }
+
+ return kTRUE;
+}
--- /dev/null
+#ifndef ALIDIELECTRONTRACKROTATOR_H
+#define ALIDIELECTRONTRACKROTATOR_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronTrackRotator #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TNamed.h>
+
+class TObjArray;
+class AliVTrack;
+
+class AliDielectronTrackRotator : public TNamed {
+public:
+ enum ERotationType {kRotatePositive, kRotateNegative, kRotateBothRandom};
+
+ AliDielectronTrackRotator();
+ AliDielectronTrackRotator(const char*name, const char* title);
+
+ virtual ~AliDielectronTrackRotator();
+
+ void SetTrackArrays(const TObjArray * const arrP, const TObjArray * const arrN) {fkArrTracksP=arrP;fkArrTracksN=arrN;}
+ void Reset();
+ Bool_t NextCombination();
+
+ //Setters
+ void SetIterations(UInt_t niter) { fIterations=niter; }
+ void SetRotationType(ERotationType type) { fRotationType=type; }
+ void SetStartAnglePhi(Double_t phi) { fStartAnglePhi=phi; }
+ void SetConeAnglePhi(Double_t phi) { fConeAnglePhi=phi; }
+
+ //Getters
+ Int_t GetIterations() const { return fIterations; }
+ ERotationType GetRotationType() const { return fRotationType; }
+ Double_t GetStartAnglePhi() const { return fStartAnglePhi; }
+ Double_t GetConeAnglePhi() const { return fConeAnglePhi; }
+
+
+ AliVTrack* GetTrackP() const {return fTrackP;}
+ AliVTrack* GetTrackN() const {return fTrackN;}
+
+private:
+ UInt_t fIterations; // number of iterations
+
+ ERotationType fRotationType; // which track to rotate
+
+ Double_t fStartAnglePhi; // starting angle for rotation
+ Double_t fConeAnglePhi; // opening angle in phi for multiple rotation
+
+ const TObjArray *fkArrTracksP; //! array of positive tracks
+ const TObjArray *fkArrTracksN; //! array of negative tracks
+
+ UInt_t fCurrentIteration; //! current iteration step
+ Int_t fCurrentTackP; //! current positive track in array
+ Int_t fCurrentTackN; //! current negative track in array
+
+ AliVTrack *fTrackP; //! Positive track
+ AliVTrack *fTrackN; //! Negative track
+
+ Bool_t RotateTracks();
+
+ AliDielectronTrackRotator(const AliDielectronTrackRotator &c);
+ AliDielectronTrackRotator &operator=(const AliDielectronTrackRotator &c);
+
+
+ ClassDef(AliDielectronTrackRotator,1) // Dielectron TrackRotator
+};
+
+
+
+#endif
};
AliESDpid* AliDielectronVarManager::fgESDpid = 0x0;
+AliAODpidUtil* AliDielectronVarManager::fgAODpidUtil = 0x0;
AliVEvent* AliDielectronVarManager::fgEvent = 0x0;
AliKFVertex* AliDielectronVarManager::fgKFVertex = 0x0;
//________________________________________________________________
#include <AliExternalTrackParam.h>
#include <AliESDpid.h>
+#include <AliAODpidUtil.h>
#include <AliPID.h>
#include "AliDielectronPair.h"
static void Fill(const TObject* particle, Double_t * const values);
static void InitESDpid(Int_t type=0);
+ static void InitAODpidUtil(Int_t type=0);
static void SetESDpid(AliESDpid * const pid) {fgESDpid=pid;}
static AliESDpid* GetESDpid() {return fgESDpid;}
+ static AliAODpidUtil* GetAODpidUtil() {return fgAODpidUtil;}
static void SetEvent(AliVEvent * const ev);
-
+ static Bool_t GetDCA(const AliAODTrack *track, Double_t d0z0[2]);
+
static const char* GetValueName(Int_t i) { return (i>=0&&i<kNMaxValues)?fgkParticleNames[i]:""; }
private:
static void FillVarAODEvent(const AliAODEvent *event, Double_t * const values);
static void FillVarMCEvent(const AliMCEvent *event, Double_t * const values);
- static AliESDpid* fgESDpid; // ESD pid object
- static AliVEvent* fgEvent; // current event pointer
- static AliKFVertex *fgKFVertex; // kf vertex
+ static AliESDpid *fgESDpid; // ESD pid object
+ static AliAODpidUtil *fgAODpidUtil; // AOD pid object
+ static AliVEvent *fgEvent; // current event pointer
+ static AliKFVertex *fgKFVertex; // kf vertex
AliDielectronVarManager(const AliDielectronVarManager &c);
AliDielectronVarManager &operator=(const AliDielectronVarManager &c);
values[AliDielectronVarManager::kE] = particle->E();
values[AliDielectronVarManager::kM] = particle->M();
values[AliDielectronVarManager::kCharge] = particle->Charge();
-
+
+ values[AliDielectronVarManager::kPdgCode] = particle->PdgCode();
+
if ( fgEvent ) AliDielectronVarManager::Fill(fgEvent, values);
}
// TODO: for the moment we set the bethe bloch parameters manually
// this should be changed in future!
- values[AliDielectronVarManager::kTPCnSigmaEle]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kElectron);
- //-AliDielectronPID::GetCorrVal();
+ values[AliDielectronVarManager::kTPCnSigmaEle]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kElectron)-AliDielectronPID::GetCorrVal();
values[AliDielectronVarManager::kTPCnSigmaPio]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kPion);
values[AliDielectronVarManager::kTPCnSigmaMuo]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kMuon);
values[AliDielectronVarManager::kTPCnSigmaKao]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kKaon);
values[AliDielectronVarManager::kTRDprobPio] = 0;
//TODO: This is only an approximation!!!
- values[AliDielectronVarManager::kTPCsignalN] = values[AliDielectronVarManager::kNclsTPC];
+ values[AliDielectronVarManager::kTPCsignalN] = 0;
-// Fill AliAODTrack interface information
- // ...
- values[AliDielectronVarManager::kImpactParXY] = particle->DCA();
- values[AliDielectronVarManager::kImpactParZ] = particle->ZAtDCA();
+ // Fill AliAODTrack interface information
+ //
+ Double_t d0z0[2];
+ GetDCA(particle, d0z0);
+ values[AliDielectronVarManager::kImpactParXY] = d0z0[0];
+ values[AliDielectronVarManager::kImpactParZ] = d0z0[1];
values[AliDielectronVarManager::kPIn]=0;
values[AliDielectronVarManager::kTPCsignal]=0;
AliAODPid *pid=particle->GetDetPid();
if (pid){
Double_t mom =pid->GetTPCmomentum();
- //TODO: kTPCsignalN is only an approximation (see above)!!
-
- Double_t tpcNsigmaEle=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
- TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]) ,AliPID::kElectron);
- Double_t tpcNsigmaPio=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
- TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kPion);
- Double_t tpcNsigmaMuo=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
- TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kMuon);
- Double_t tpcNsigmaKao=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
- TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kKaon);
- Double_t tpcNsigmaPro=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
- TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kProton);
+ values[AliDielectronVarManager::kTPCsignalN] = pid->GetTPCsignalN();
+ Double_t tpcNsigmaEle=fgAODpidUtil->NumberOfSigmasTPC(particle,AliPID::kElectron);
+ Double_t tpcNsigmaPio=fgAODpidUtil->NumberOfSigmasTPC(particle,AliPID::kPion);
+ Double_t tpcNsigmaMuo=fgAODpidUtil->NumberOfSigmasTPC(particle,AliPID::kMuon);
+ Double_t tpcNsigmaKao=fgAODpidUtil->NumberOfSigmasTPC(particle,AliPID::kKaon);
+ Double_t tpcNsigmaPro=fgAODpidUtil->NumberOfSigmasTPC(particle,AliPID::kProton);
values[AliDielectronVarManager::kPIn]=mom;
values[AliDielectronVarManager::kTPCsignal]=pid->GetTPCsignal();
//
// Fill pair information available for histogramming into an array
//
-
+
+ values[AliDielectronVarManager::kPdgCode]=0;
+ values[AliDielectronVarManager::kPdgCodeMother]=0;
+
// Fill common AliVParticle interface information
FillVarVParticle(pair, values);
fgESDpid->GetTPCResponse().SetSigma(3.79301e-03, 2.21280e+04);
}
+inline void AliDielectronVarManager::InitAODpidUtil(Int_t type)
+{
+ if (!fgAODpidUtil) fgAODpidUtil=new AliAODpidUtil;
+ Double_t alephParameters[5];
+ // simulation
+ alephParameters[0] = 2.15898e+00/50.;
+ alephParameters[1] = 1.75295e+01;
+ alephParameters[2] = 3.40030e-09;
+ alephParameters[3] = 1.96178e+00;
+ alephParameters[4] = 3.91720e+00;
+ fgAODpidUtil->GetTOFResponse().SetTimeResolution(80.);
+
+ // data
+ if (type==1){
+ alephParameters[0] = 0.0283086/0.97;
+ alephParameters[1] = 2.63394e+01;
+ alephParameters[2] = 5.04114e-11;
+ alephParameters[3] = 2.12543e+00;
+ alephParameters[4] = 4.88663e+00;
+ fgAODpidUtil->GetTOFResponse().SetTimeResolution(130.);
+ fgAODpidUtil->GetTPCResponse().SetMip(50.);
+ }
+
+ fgAODpidUtil->GetTPCResponse().SetBetheBlochParameters(
+ alephParameters[0],alephParameters[1],alephParameters[2],
+ alephParameters[3],alephParameters[4]);
+
+ fgAODpidUtil->GetTPCResponse().SetSigma(3.79301e-03, 2.21280e+04);
+}
+
inline void AliDielectronVarManager::SetEvent(AliVEvent * const ev)
{
fgKFVertex=0x0;
if (ev && ev->GetPrimaryVertex()) fgKFVertex=new AliKFVertex(*ev->GetPrimaryVertex());
}
+
+
+inline Bool_t AliDielectronVarManager::GetDCA(const AliAODTrack *track, Double_t d0z0[2])
+{
+ if(track->TestBit(AliAODTrack::kIsDCA)){
+ d0z0[0]=track->DCA();
+ d0z0[1]=track->ZAtDCA();
+ return kTRUE;
+ }
+
+ Double_t covd0z0[3];
+ AliAODTrack copy(*track);
+ AliAODVertex *vtx =(AliAODVertex*)(fgEvent->GetPrimaryVertex());
+ Double_t fBzkG = fgEvent->GetMagneticField(); // z componenent of field in kG
+ Bool_t ok = copy.PropagateToDCA(vtx,fBzkG,kVeryBig,d0z0,covd0z0);
+ if(!ok){
+ d0z0[0]=-999.;
+ d0z0[1]=-999.;
+ }
+ return ok;
+}
+
/*
inline void AliDielectronVarManager::FillValues(const TParticle *particle, Double_t *values)
{
//check for output aod handler
if (!mgr->GetOutputEventHandler()||mgr->GetOutputEventHandler()->IsA()!=AliAODHandler::Class()) {
- Error("AddTaskJPSIFilter","No AOD output handler available. Cannot Proceed!");
+ Warning("AddTaskJPSIFilter","No AOD output handler available. Not adding the task!");
return 0;
}
//Do we run on AOD?
Bool_t isAOD=mgr->GetInputEventHandler()->IsA()==AliAODInputHandler::Class();
+
+ if(isAOD) {
+ //add options to AliAODHandler to duplicate input event
+ AliAODHandler *aodHandler = (AliAODHandler*)mgr->GetOutputEventHandler();
+ aodHandler->SetCreateNonStandardAOD();
+ aodHandler->SetNeedsHeaderReplication();
+ aodHandler->SetNeedsTracksBranchReplication();
+ aodHandler->SetNeedsVerticesBranchReplication();
+ aodHandler->SetNeedsV0sBranchReplication();
+ aodHandler->SetNeedsCascadesBranchReplication();
+ aodHandler->SetNeedsTrackletsBranchReplication();
+ aodHandler->SetNeedsPMDClustersBranchReplication();
+ aodHandler->SetNeedsJetsBranchReplication();
+ aodHandler->SetNeedsFMDClustersBranchReplication();
+ aodHandler->SetNeedsCaloClustersBranchReplication();
+ //aodHandler->SetNeedsMCParticlesBranchReplication();
+ aodHandler->SetNeedsDimuonsBranchReplication();
+ }
//Create task and add it to the analysis manager
AliAnalysisTaskDielectronFilter *task=new AliAnalysisTaskDielectronFilter("jpsi_DielectronFilter");
// if (cutDefinition<nDie-1)
InitCFDieleData(diele, cutDefinition, isAOD);
+ AliDielectronTrackRotator *rot=new AliDielectronTrackRotator;
+ rot->SetIterations(10);
+ diele->SetTrackRotator(rot);
return diele;
}
for (Int_t i=0; i<3; ++i){
histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(i)));
}
+ //track rotation
+ histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(10)));
+ histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(10)));
+
//add histograms to event class
if (cutDefinition==0){
cf->AddVariable(AliDielectronVarManager::kY,40,-2,2);
cf->AddVariable(AliDielectronVarManager::kM,50,1.98,1.98+50*.04); //40Mev Steps
- cf->AddVariable(AliDielectronVarManager::kPairType,10,0,10);
+ cf->AddVariable(AliDielectronVarManager::kPairType,11,0,11);
//leg variables
cf->AddVariable(AliDielectronVarManager::kPt,"0.0, 0.8, 1.2, 100.0",kTRUE);
cf->AddVariable(AliDielectronVarManager::kNclsTPC,"0, 100, 120, 160",kTRUE);
Int_t baseColors[5]={kRed, kGreen+1, kAzure-4, kMagenta, kCyan+1};
Int_t sigmaColorOffset=0;
+Int_t baseColors[5]={kRed, kRed, kRed, kRed, kRed};
Double_t sigmas[5]={3,3,3,3,3};
Double_t masses[5];
alephParameters[2] = 5.04114e-11;
alephParameters[3] = 2.12543e+00;
alephParameters[4] = 4.88663e+00;
- Double_t mip=49.2;
+ Double_t mip=50;
Color_t color=kRed;
Int_t lineWidth=2;
c->SetAlias("TOFcut","abs(Leg1_TOF_nSigma_Electrons)<3&&abs(Leg2_TOF_nSigma_Electrons)<3");
c->SetAlias("TOFcut2","(Leg1_P_InnerParam<1.3&&abs(Leg1_TOF_nSigma_Electrons)<3||Leg1_P_InnerParam>=1.3)&&(Leg2_P_InnerParam<1.3&&abs(Leg2_TOF_nSigma_Electrons)<3||Leg2_P_InnerParam>=1.3)");
c->SetAlias("TPCcut","(Leg1_TPC_signal>70*(1-exp(-1*(Leg1_P_InnerParam+2))))&&(Leg2_TPC_signal>70*(1-exp(-1*(Leg2_P_InnerParam+2))))")
-c->SetAlias("NClcut","Leg1_NclsTPC>120&&Leg2_NclsTPC>120");
+c->SetAlias("NClcut","Leg1_NclsTPC>90&&Leg2_NclsTPC>90");
c->SetAlias("eleParam","Leg1_TPC_nSigma_Electrons<5&&Leg2_TPC_nSigma_Electrons<5&&Leg1_TPC_nSigma_Electrons>-2.65*exp(-0.6757*Leg1_P_InnerParam)&&Leg2_TPC_nSigma_Electrons>-2.65*exp(-0.6757*Leg2_P_InnerParam)")
c->SetAlias("cut","PairType==1&&eleParam&&Run<127719")
+c->SetAlias("NClcut","Leg1_NclsTPC>90&&Leg2_NclsTPC>90");
+c->SetAlias("Ptcut","Leg1_Pt>1&&Leg2_Pt>1")
+c->SetAlias("PairT","PairType==1");
+c->SetAlias("cut","NClcut&&Ptcut&&PairT")
+
+
+c->SetAlias("cutE","abs(Leg1_TPC_nSigma_Electrons)<3&&abs(Leg2_TPC_nSigma_Electrons)<3");
+c->SetAlias("cutPi","abs(Leg1_TPC_nSigma_Pions)>3&&abs(Leg2_TPC_nSigma_Pions)>3");
+c->SetAlias("cutP","(Leg1_TPC_nSigma_Protons)>3&&(Leg2_TPC_nSigma_Protons)>3");
+c->SetAlias("cut","NClcut&&Ptcut&&PairT&&cutE&&cutPi&&cutP")
+
+c->SetAlias("eta","abs(Leg1_Eta)<0.88&&abs(Leg1_Eta)<0.88");
+c->SetAlias("rap","abs(Y)<0.88");
+c->SetAlias("spdFirst","(Leg1_ITS_clusterMap&1)==1&&(Leg2_ITS_clusterMap&1)==1")
+c->SetAlias("cut","NClcut&&Ptcut&&PairT&&cutE&&cutPi&&cutP&&eta&&rap&&spdFirst")
+ c->Draw("M>>hM(50,1.98,1.98+50*.04)","cut","e")
+
+
+c->SetAlias("cutProL1",Form("Leg1_TPC_nSigma_Electrons>(%f*%f+(AliExternalTrackParam::BetheBlochAleph(Leg1_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)-AliExternalTrackParam::BetheBlochAleph(Leg1_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)))/%f",nSigma,resolution, AliPID::ParticleMass(AliPID::kProton), AliPID::ParticleMass(AliPID::kElectron), resolution))
+
+c->SetAlias("cutProL2",Form("Leg2_TPC_nSigma_Electrons>(%f*%f+(AliExternalTrackParam::BetheBlochAleph(Leg2_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)-AliExternalTrackParam::BetheBlochAleph(Leg2_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)))/%f",nSigma,resolution, AliPID::ParticleMass(AliPID::kProton), AliPID::ParticleMass(AliPID::kElectron), resolution))
+
+c->SetAlias("cutPioL1",Form("Leg1_TPC_nSigma_Electrons>(%f*%f+(AliExternalTrackParam::BetheBlochAleph(Leg1_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)-AliExternalTrackParam::BetheBlochAleph(Leg1_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)))/%f",nSigma,resolution, AliPID::ParticleMass(AliPID::kPion), AliPID::ParticleMass(AliPID::kElectron), resolution))
+
+c->SetAlias("cutPioL2",Form("Leg2_TPC_nSigma_Electrons>(%f*%f+(AliExternalTrackParam::BetheBlochAleph(Leg2_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)-AliExternalTrackParam::BetheBlochAleph(Leg2_P_InnerParam/%f,0.0283086/0.97,2.63394e+01,5.04114e-11,2.12543e+00,4.88663e+00)))/%f",nSigma,resolution, AliPID::ParticleMass(AliPID::kPion), AliPID::ParticleMass(AliPID::kElectron), resolution))
+
+c->SetAlias("cutPro","cutProL1&&cutProL2");
+c->SetAlias("cutPio","cutPioL1&&cutPioL2");
*/
--- /dev/null
+#include <TROOT.h>
+#include <TStyle.h>
+#include <TH1D.h>
+#include <TObjArray.h>
+#include <TCanvas.h>
+#include <TMath.h>
+#include <TAxis.h>
+
+#include "AliDielectronSignalExt.h"
+#include "AliDielectronCFdraw.h"
+#include "AliDielectron.h"
+
+/*
+gSystem->AddIncludePath("-I/data/Work/software/svngsi/dielectron/trunk/dielectron");
+.L PlotDataResults.C+g
+PlotDataResuts("/data/Work/train/V005.data/2010-10-29_2144.3552/mergedPeriods/data/7TeV/LHC10c.pass2/jpsi.root");
+PlotDataResuts("/data/Work/train/V005.data/2010-10-21_2342.3445/mergedPeriods/data/7TeV/LHC10pass2/jpsi.root")
+PlotDataResuts("/data/Work/train/V005.data/2010-10-21_2342.3445/mergedPeriods/data/7TeV/LHC10d.pass1/wiechula_jpsi.root");
+*/
+
+AliDielectronSignalBase* GetSignalLS(AliDielectronCFdraw &d, Int_t step, const char* nameAdd);
+AliDielectronSignalBase* GetSignalRot(AliDielectronCFdraw &d, Int_t step, const char* nameAdd);
+void DrawSpectra(AliDielectronSignalBase *sig, const char* cname);
+
+//_______________________________________
+void PlotDataResuts(const char* filename)
+{
+ AliDielectronCFdraw d(filename);
+
+ Int_t stepFirst=0, stepAny=1;
+
+ gStyle->SetOptStat(0);
+ //Set common Ranges
+// d.SetRangeUser("Leg1_Pt",1.,1000.);
+// d.SetRangeUser("Leg2_Pt",1.,1000.);
+
+ //============================
+ //SPD first
+ //
+
+ //--- Like sign subtraction
+ AliDielectronSignalBase *sigFirst=GetSignalLS(d,stepFirst,"First");
+ DrawSpectra(sigFirst,"cFirst");
+ //--- Rotation subtraction
+ AliDielectronSignalBase *sigFirstRot=GetSignalRot(d,stepFirst,"FirstRot");
+ DrawSpectra(sigFirstRot,"cFirstRot");
+
+ //============================
+ //SPD any
+ //
+ AliDielectronSignalBase *sigAny=GetSignalLS(d,stepAny,"Any");
+ DrawSpectra(sigAny,"cAny");
+ //--- Rotation subtraction
+ AliDielectronSignalBase *sigAnyRot=GetSignalRot(d,stepAny,"AnyRot");
+ DrawSpectra(sigAnyRot,"cAnyRot");
+
+}
+
+
+//_______________________________________
+AliDielectronSignalBase *GetSignalLS(AliDielectronCFdraw &d, Int_t step, const char* nameAdd)
+{
+ //
+ // Get Extracted signal from likesign method
+ //
+
+ TObjArray *arr=new TObjArray;
+ arr->SetOwner();
+
+ for (Int_t iType=0;iType<3;++iType){
+ d.SetRangeUser("PairType",iType,iType);
+ arr->AddAt(d.Project("M",step),iType);
+ }
+
+ AliDielectronSignalExt *sig=new AliDielectronSignalExt;
+ sig->SetIntegralRange(2.9,3.15);
+ sig->SetMethod(AliDielectronSignalBase::kLikeSign);
+ sig->Process(arr);
+
+ delete arr;
+ return sig;
+}
+
+//_______________________________________
+AliDielectronSignalBase *GetSignalRot(AliDielectronCFdraw &d, Int_t step, const char* nameAdd)
+{
+ //
+ // Get Extracted signal from likesign method
+ //
+
+ TObjArray *arr=new TObjArray;
+ arr->SetOwner();
+
+ Int_t iType=AliDielectron::kEv1PM;
+ d.SetRangeUser("PairType",iType,iType);
+ arr->AddAt(d.Project("M",step),iType);
+
+ iType=AliDielectron::kEv1PMRot;
+ d.SetRangeUser("PairType",iType,iType);
+ arr->AddAt(d.Project("M",step),iType);
+
+ AliDielectronSignalExt *sig=new AliDielectronSignalExt;
+ sig->SetIntegralRange(2.9,3.15);
+ sig->SetMethod(AliDielectronSignalBase::kRotation);
+ sig->Process(arr);
+
+ delete arr;
+ return sig;
+}
+
+//_______________________________________
+void DrawSpectra(AliDielectronSignalBase *sig, const char* cname)
+{
+ //
+ //
+ //
+ gStyle->SetOptTitle(0);
+ TCanvas *c=(TCanvas*)gROOT->FindObject(cname);
+ if (!c) c=new TCanvas(cname,cname,400,600);
+ c->Clear();
+ c->Divide(1,2,0,0);
+ c->cd(1);
+ gPad->SetTopMargin(0.01);
+ gPad->SetRightMargin(0.01);
+
+ gPad->SetBottomMargin(0);
+ gPad->SetGridx();
+ gPad->SetGridy();
+ gPad->SetTickx();
+ gPad->SetTicky();
+
+ TH1 *hUS=sig->GetUnlikeSignHistogram();
+ hUS->SetMarkerStyle(20);
+ hUS->SetMarkerSize(0.7);
+ hUS->SetMarkerColor(kRed);
+ hUS->SetLineColor(kRed);
+ hUS->SetStats(0);
+
+ TH1* hBackground=sig->GetBackgroundHistogram();
+ hBackground->SetMarkerStyle(24);
+ hBackground->SetMarkerSize(0.7);
+ hBackground->SetStats(0);
+ hBackground->SetMarkerColor(kBlue);
+ hBackground->SetLineColor(kBlue);
+
+ hUS->Draw();
+ hBackground->Draw("same");
+
+ c->cd(2);
+ gPad->SetRightMargin(0.01);
+ gPad->SetGridx();
+ gPad->SetGridy();
+ gPad->SetTickx();
+ gPad->SetTicky();
+ gPad->SetTopMargin(0);
+
+ TH1* hSignal=sig->GetSignalHistogram();
+ hSignal->SetMarkerStyle(20);
+ hSignal->SetMarkerSize(0.7);
+ hSignal->SetMarkerColor(kRed);
+ hSignal->SetLineColor(kRed);
+ hSignal->Draw();
+}
+
dielectron/AliDielectronSignalExt.cxx \
dielectron/AliDielectronSpectrum.cxx \
dielectron/AliDielectronDebugTree.cxx \
- dielectron/AliDielectronPID.cxx
+ dielectron/AliDielectronTrackRotator.cxx \
+ dielectron/AliDielectronPID.cxx \
+ dielectron/AliDielectronCutGroup.cxx \
+ dielectron/AliDielectronEventCuts.cxx
HDRS= $(SRCS:.cxx=.h)
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff