#include "TMath.h"
#include "TObjArray.h"
#include "TObjString.h"
+#include <cassert>
ClassImp(AliAnalysisMuMuBinning::Range)
ClassImp(AliAnalysisMuMuBinning)
/// add one bin
AddBin(bin.Particle().Data(),bin.Type().Data(),
bin.Xmin(),bin.Xmax(),
- bin.Ymin(),bin.Ymax());
+ bin.Ymin(),bin.Ymax(),bin.Flavour());
}
-
+
//______________________________________________________________________________
void AliAnalysisMuMuBinning::AddBin(const char* particle, const char* type,
Double_t xmin, Double_t xmax,
Double_t ymin, Double_t ymax,
- Bool_t warn)
+ const char* flavour)
{
/// Add a bin
/// Note that particle and type are not case sensitive.
fBins->Add(new TObjString(sparticle),b);
}
- Range* r = new Range(sparticle.Data(),stype.Data(),xmin,xmax,ymin,ymax);
+ Range* r = new Range(sparticle.Data(),stype.Data(),xmin,xmax,ymin,ymax,flavour);
if ( b->FindObject(r) )
{
- if (warn)
- {
- AliWarning("Trying to add an already existing bin. Not doing it.");
- }
+ AliDebug(1,"Trying to add an already existing bin. Not doing it.");
delete r;
}
else
//______________________________________________________________________________
-TObjArray* AliAnalysisMuMuBinning::CreateBinObjArray(const char* particle, const char* type) const
+TObjArray* AliAnalysisMuMuBinning::CreateBinObjArray(const char* particle, const char* type, const char* flavour) const
{
/// Get the list of bins for a given particle and given type
/// The returned array must be deleted by the user
TString stype(type);
stype.ToUpper();
+ TString sflavour(flavour);
+
TObjArray* types = stype.Tokenize(",");
TObjString* onetype;
TIter nextType(types);
nextType.Reset();
while ( ( onetype = static_cast<TObjString*>(nextType()) ) )
{
- if ( r->Type() == onetype->String() )
+ if ( r->Type() == onetype->String() &&
+ ( ( sflavour.Length() > 0 && r->Flavour() == sflavour.Data() ) || sflavour.Length()==0 ) )
{
a->Add(r->Clone());
}
//______________________________________________________________________________
void AliAnalysisMuMuBinning::CreateMesh(const char* particle,
const char* type1, const char* type2,
+ const char* flavour,
Bool_t remove12)
{
/// Create 2D bins from existing 1d ones of type1 and type2
- TObjArray* a1 = CreateBinObjArray(particle,type1);
+ TObjArray* a1 = CreateBinObjArray(particle,type1,flavour);
if (!a1)
{
AliError(Form("No bin for type %s. Done nothing.",type1));
return;
}
- TObjArray* a2 = CreateBinObjArray(particle,type2);
+ TObjArray* a2 = CreateBinObjArray(particle,type2,flavour);
if (!a2)
{
AliError(Form("No bin for type %s. Done nothing.",type2));
return;
}
- TString meshType(Form("%s VS %s",type1,type2));
+ TString meshType(Form("%s VS %s - %s",type1,type2,flavour));
for ( Int_t i1 = 0; i1 <= a1->GetLast(); ++i1 )
{
{
Range* r2 = static_cast<Range*>(a2->At(i2));
- AddBin(particle,meshType,r2->Xmin(),r2->Xmax(),r1->Xmin(),r1->Xmax(),kFALSE);
+ AddBin(particle,meshType,r2->Xmin(),r2->Xmax(),r1->Xmin(),r1->Xmax(),Form("%s VS %s",r1->Flavour().Data(),r2->Flavour().Data()));
}
}
//______________________________________________________________________________
AliAnalysisMuMuBinning*
-AliAnalysisMuMuBinning::Project(const char* particle, const char* type) const
+AliAnalysisMuMuBinning::Project(const char* particle, const char* type, const char* flavour) const
{
/// Create a sub-binning object with only the bins pertaining to (particle,type)
- TObjArray* bins = CreateBinObjArray(particle,type);
+ TObjArray* bins = CreateBinObjArray(particle,type,flavour);
if (!bins) return 0x0;
AliAnalysisMuMuBinning* p = new AliAnalysisMuMuBinning;
TIter next(bins);
while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(next())) )
{
- if (bin->Type()==stype)
+ assert (bin->Type()==stype && bin->Flavour()==flavour);
{
- p->AddBin(particle,bin->Type(),bin->Xmin(),bin->Xmax(),bin->Ymin(),bin->Ymax());
+ p->AddBin(particle,bin->Type(),bin->Xmin(),bin->Xmax(),bin->Ymin(),bin->Ymax(),bin->Flavour().Data());
}
}
//______________________________________________________________________________
AliAnalysisMuMuBinning::Range::Range(const char* particle, const char* type,
Double_t xmin, Double_t xmax,
- Double_t ymin, Double_t ymax)
-: TObject(), fParticle(particle), fType(type), fXmin(xmin), fXmax(xmax), fYmin(ymin), fYmax(ymax)
+ Double_t ymin, Double_t ymax,
+ const char* flavour)
+: TObject(), fParticle(particle), fType(type),
+ fXmin(xmin), fXmax(xmax), fYmin(ymin), fYmax(ymax),
+ fFlavour(flavour)
{
/// ctor
fParticle.ToUpper();
TString s;
- s.Form("%s_%05.2f_%05.2f",Type().Data(),Xmin(),Xmax());
+ if ( fFlavour.Length() > 0 )
+ {
+ s.Form("%s_%s_%05.2f_%05.2f",Type().Data(),Flavour().Data(),Xmin(),Xmax());
+ }
+ else
+ {
+ s.Form("%s_%05.2f_%05.2f",Type().Data(),Xmin(),Xmax());
+ }
if (Is2D())
{
if (s) return s;
+ s = strcmp(Flavour().Data(),other->Flavour().Data());
+
+ if (s) return s;
+
+ if ( IsNullObject() ) return 0;
+
if ( Xmin() < other->Xmin() )
{
return -1;
std::cout << Form(" ; %5.2f : %5.2f",Ymin(),Ymax());
}
+ if (Flavour().Length()>0)
+ {
+ std::cout << " - " << Flavour().Data();
+ }
+
std::cout << "->" << AsString().Data() << std::endl;
}
Double_t xmax=TMath::Limits<Double_t>::Max(),
Double_t ymin=TMath::Limits<Double_t>::Max(),
Double_t ymax=TMath::Limits<Double_t>::Max(),
- Bool_t warn=kTRUE);
+ const char* flavour="");
TObjArray* CreateParticleArray() const;
TObjArray* CreateBinObjArray() const;
TObjArray* CreateBinObjArray(const char* particle) const;
- TObjArray* CreateBinObjArray(const char* particle, const char* type) const;
+ TObjArray* CreateBinObjArray(const char* particle, const char* type, const char* flavour) const;
- AliAnalysisMuMuBinning* Project(const char* particle, const char* type) const;
+ AliAnalysisMuMuBinning* Project(const char* particle, const char* type, const char* flavour="") const;
virtual void Print(Option_t* opt="") const;
- void CreateMesh(const char* particle, const char* type1, const char* type2, Bool_t remove12=kFALSE);
+ void CreateMesh(const char* particle, const char* type1, const char* type2, const char* flavour="", Bool_t remove12=kFALSE);
Long64_t Merge(TCollection* list);
Double_t xmin=TMath::Limits<Double_t>::Max(),
Double_t xmax=TMath::Limits<Double_t>::Max(),
Double_t ymin=TMath::Limits<Double_t>::Max(),
- Double_t ymax=TMath::Limits<Double_t>::Max());
+ Double_t ymax=TMath::Limits<Double_t>::Max(),
+ const char* version="");
virtual Int_t Compare(const TObject* obj) const;
Bool_t IsEqual(const TObject* obj) const { return Compare(obj)==0; }
Bool_t IsInRange(Double_t x, Double_t y=TMath::Limits<Double_t>::Max()) const;
+ TString Flavour() const { return fFlavour; }
+
private:
TString fParticle; // particle
TString fType; // binning type (e.g. pt, y, phi)
Double_t fXmax; // x-max of the range
Double_t fYmin; // x-min of the range
Double_t fYmax; // x-max of the range
+ TString fFlavour; // flavour (if any) this range, e.g. coarse, fine, etc...
- ClassDef(AliAnalysisMuMuBinning::Range,1)
+ ClassDef(AliAnalysisMuMuBinning::Range,2)
};
TMap* fBins; // list of bins (particle -> list of bins) = (TObjString -> TObjArray)
- ClassDef(AliAnalysisMuMuBinning,1)
+ ClassDef(AliAnalysisMuMuBinning,1) // custom binning for MuMu analysis
};
#endif
//_____________________________________________________________________________
void AliAnalysisTaskMuMu::AddBin(const char* particle, const char* type,
Double_t xmin, Double_t xmax,
- Double_t ymin, Double_t ymax)
+ Double_t ymin, Double_t ymax,
+ const char* flavour)
{
/// Add one bin
- fBinning->AddBin(particle,type,xmin,xmax,ymin,ymax);
+ fBinning->AddBin(particle,type,xmin,xmax,ymin,ymax,flavour);
// invalidate cached bins, if any
if (fBinArray)
}
//_____________________________________________________________________________
-void AliAnalysisTaskMuMu::CreateMesh(const char* particle, const char* type1, const char* type2, Bool_t remove12)
+void AliAnalysisTaskMuMu::CreateMesh(const char* particle, const char* type1, const char* type2, const char* flavour, Bool_t remove12)
{
/// Create a 2d binning from 2 1d binning.
/// WARNING : not fully tested yet
AliError("Cannot create a mesh as I have no bin at all !");
return;
}
- fBinning->CreateMesh(particle,type1,type2,remove12);
+ fBinning->CreateMesh(particle,type1,type2,flavour,remove12);
}
//_____________________________________________________________________________
Int_t nMCMinvBins = GetNbins(minvMin,minvMax,0.1);
- TObjArray* bins = fBinning->CreateBinObjArray("psi","pt vs y,pt,y,phi");
+ TObjArray* bins = fBinning->CreateBinObjArray("psi","y vs pt,integrated,pt,y,phi","");
CreatePairHisto(physics,triggerClassName,"Pt","#mu+#mu- Pt distribution",
200,0,20);
void AliAnalysisTaskMuMu::DefineDefaultBinning()
{
fBinning = new AliAnalysisMuMuBinning("BIN");
- fBinning->AddBin("psi","pt vs y");
+ fBinning->AddBin("psi","integrated");
}
//_____________________________________________________________________________
if (!fBinArray)
{
- fBinArray = fBinning->CreateBinObjArray("psi","pt vs y,pt,y");
+ fBinArray = fBinning->CreateBinObjArray("psi","y vs pt,integrated,pt,y","");
}
TIter nextBin(fBinArray);
{
Bool_t ok(kFALSE);
- if ( r->Is2D() || r->IsNullObject() )
+ if ( r->IsNullObject() )
{
- ok = r->IsInRange(part->Y(),part->Pt());
+ ok = kTRUE;
+ }
+ else if ( r->Is2D() )
+ {
+ if ( r->AsString().BeginsWith("PTVSY") )
+ {
+ ok = r->IsInRange(part->Y(),part->Pt());
+ }
+ else if ( r->AsString().BeginsWith("YVSPT") )
+ {
+ ok = r->IsInRange(part->Pt(),part->Y());
+ }
+ else
+ {
+ AliError(Form("Don't know how to deal with 2D bin %s",r->AsString().Data()));
+ }
}
else
{
{
ok = r->IsInRange(part->Y());
}
+ else if ( r->Type() == "PHI" )
+ {
+ ok = r->IsInRange(part->Phi());
+ }
}
if ( ok )
{
// Fill event-wise histograms
+// AliCodeTimerAuto("",0);
+
if (!IsHistogramDisabled("BCX"))
{
Histo(physics,triggerClassName,centrality,"BCX")->Fill(1.0*Event()->GetBunchCrossNumber());
}
+
//_____________________________________________________________________________
void AliAnalysisTaskMuMu::FillHistogramCollection(const char* physics, const char* triggerClassName)
{
{
/// Fill histograms for one track
+// AliCodeTimerAuto("",0);
+
TLorentzVector p(track.Px(),track.Py(),track.Pz(),
TMath::Sqrt(MuonMass2()+track.P()*track.P()));
{
/// Fill histograms for /physics/triggerClassName/centrality
+// AliCodeTimerAuto("",0);
+
FillEventHistos(physics,triggerClassName,centrality);
// Track loop
if (!fBinArray)
{
- fBinArray = fBinning->CreateBinObjArray("psi","pt vs y,pt,y,phi");
+ fBinArray = fBinning->CreateBinObjArray("psi","y vs pt,integrated,pt,y,phi","");
}
Int_t nMuonTracks = AliAnalysisMuonUtility::GetNTracks(Event());
while ( ( r = static_cast<AliAnalysisMuMuBinning::Range*>(nextBin()) ) )
{
- AliDebug(2,Form("bin %s pt %e rad %e = %d",r->AsString().Data(),
- pj.Pt(),pj.Rapidity(),r->IsInRange(pj.Pt(),pj.Rapidity())));
+// AliInfo(Form("%s y %e pt %e ok1 %d ok2 %d",
+// r->AsString().Data(),
+// pj.Rapidity(),
+// pj.Pt(),
+// r->IsInRange(pj.Rapidity(),pj.Pt()),
+// r->IsInRange(pj.Pt(),pj.Rapidity())));
+
Bool_t ok(kFALSE);
- if ( r->Is2D() || r->IsNullObject() )
+ if ( r->IsNullObject() )
{
- ok = r->IsInRange(pj.Rapidity(),pj.Pt());
+ ok = kTRUE;
+ }
+ else if ( r->Is2D() )
+ {
+ if ( r->AsString().BeginsWith("PTVSY") )
+ {
+ ok = r->IsInRange(pj.Rapidity(),pj.Pt());
+ }
+ else if ( r->AsString().BeginsWith("YVSPT") )
+ {
+ ok = r->IsInRange(pj.Pt(),pj.Rapidity());
+ }
+ else
+ {
+ AliError(Form("Don't know how to deal with 2D bin %s",r->AsString().Data()));
+ }
}
else
{
{
AliError(Form("Could not get %s",hname.Data()));
}
- else
- {
- AliDebug(1,Form("filling %s",hname.Data()));
- }
h->Fill(pj.M());
}
}
*/
+// AliCodeTimerAuto("",0);
+
UInt_t eMask = EventCuts()->GetSelectionMask(fInputHandler);
UInt_t m(AliAnalysisTaskMuMu::kEventAll);
m |= AliAnalysisTaskMuMu::kEventZ7;
}
- AliVVZERO* vzero = Event()->GetVZEROData();
-
- if (vzero)
- {
- Float_t v0a = vzero->GetV0ATime();
- Float_t v0c = vzero->GetV0CTime();
-
- Float_t x = v0a-v0c;
- Float_t y = v0a+v0c;
-
- if ( ( x > 6 && x < 10 ) && y > 20 )
- {
- m |= AliAnalysisTaskMuMu::kEventV0UP;
- }
- }
+// AliVVZERO* vzero = Event()->GetVZEROData();
+//
+// if (vzero)
+// {
+// Float_t v0a = vzero->GetV0ATime();
+// Float_t v0c = vzero->GetV0CTime();
+//
+// Float_t x = v0a-v0c;
+// Float_t y = v0a+v0c;
+//
+// if ( ( x > 6 && x < 10 ) && y > 20 )
+// {
+// m |= AliAnalysisTaskMuMu::kEventV0UP;
+// }
+// }
Bool_t backgroundFlag(kFALSE);
Bool_t pileupFlag(kFALSE);
str->Print();
}
}
+
+ if ( !fEventCutNames )
+ {
+ cout << "No event cut defined yet" << endl;
+ }
+ else
+ {
+ TIter nextEventCutName(fEventCutNames);
+ TObjString* eventCutName;
+
+ while ( ( eventCutName = static_cast<TObjString*>(nextEventCutName()) ) )
+ {
+ cout << Form("EVENT CUT %s MASK %x",eventCutName->String().Data(),eventCutName->GetUniqueID()) << endl;
+ }
+ }
-// if ( !fTriggerClasses || !fTriggerClasses->First() )
-// {
-// cout << "No trigger classes defined yet" << endl;
-// }
-// else
-// {
-// cout << "Trigger classes that will be considered:" << endl;
-// TIter next(fTriggerClasses);
-// TObjString* s;
-//
-// while ( ( s = static_cast<TObjString*>(next()) ) )
-// {
-// cout << s->String().Data() << endl;
-// }
-// }
-//
if ( fBinning )
{
cout << "Binning for Minv plots" << endl;
{
/// Executed at each event
+// AliCodeTimerAuto("",0);
+
fHasMC = (MCEvent()!=0x0);
TString firedTriggerClasses(AliAnalysisMuonUtility::GetFiredTriggerClasses(Event()));
fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "HASMC", fCurrentRunNumber));
}
-// fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), eventType.Data(), fCurrentRunNumber));
-
if ( firedTriggerClasses == "" )
{
fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "EMPTY", fCurrentRunNumber));
{
fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEMBTRIGGER", fCurrentRunNumber));
}
-
- if ( AtLeastOneEmcalTrigger(firedTriggerClasses) )
- {
- fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEEMCALTRIGGER", fCurrentRunNumber));
-
- if ( AtLeastOneMBTrigger(firedTriggerClasses) )
- {
- fEventCounters->Count(Form("event:%s/trigger:%s/run:%d", et->String().Data(), "ATLEASTONEEMCALORMBTRIGGER", fCurrentRunNumber));
- }
- }
-
}
}
// second loop to count only the triggers we're interested in
-
+
TIter next(EventCuts()->GetSelectedTrigClassesInEvent(Event()));
TObjString* tname;
# include "TArrayI.h"
#endif
+#ifndef ROOT_TMath
+# include "TMath.h"
+#endif
+
class AliAnalysisMuMuBinning;
class AliCounterCollection;
class AliMergeableCollection;
AliAnalysisTaskMuMu(Bool_t fromESD, TList* triggerClassesToConsider, const char* beamYear=0x0, TArrayF* centralities=0x0);
AliAnalysisTaskMuMu(Bool_t fromESD, const char* beamYear, TArrayF* centralities=0x0);
virtual ~AliAnalysisTaskMuMu();
-
+
void AddBin(const char* particle, const char* type,
Double_t xmin, Double_t xmax,
- Double_t ymin=0.0, Double_t ymax=0.0);
+ Double_t ymin,
+ Double_t ymax,
+ const char* flavour="");
- void CreateMesh(const char* particle, const char* type1, const char* type2, Bool_t remove12=kFALSE);
+ void AddBin(const char* particle, const char* type,
+ Double_t xmin, Double_t xmax,
+ const char* flavour="") { AddBin(particle,type,xmin,xmax,TMath::Limits<Double_t>::Max(),TMath::Limits<Double_t>::Max(),flavour); }
+
+ void CreateMesh(const char* particle, const char* type1, const char* type2, const char* flavour="", Bool_t remove12=kFALSE);
virtual void AddEventCut(const char* cutName, UInt_t mask);
//--- Functions ---
class AliGenPythia;
-AliGenerator* CreateGenerator();
void Config()
{
//=========================//
// Generator Configuration //
//=========================//
- AliGenerator* gener = CreateGenerator();
+ //AliGenerator* gener = CreateGenerator();
+
+ std::cout << "VAR_GENERATOR settings " << std::endl;
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/EVGEN");
+ gROOT->LoadMacro("VAR_GENERATOR.C+");
+ AliGenerator* gener = VAR_GENERATOR();
+
gener->SetOrigin(0., 0., 0.); // Taken from OCDB
gener->SetSigma(0., 0., 0.); // Sigma in (X,Y,Z) (cm) on IP position, sigmaz taken from OCDB
// gener->SetVertexSmear(kPerEvent);
AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
}
}
-
-//
-// generators
-//
-AliGenerator* CreateGenerator()
-{
- AliGenMC* generator = 0x0;
- TString genType("VAR_GEN_TYPE");
- genType.ToUpper();
- if ( genType == "GENLIB" ) {
- generator = new AliGenParam(1, VAR_GENLIB_TYPE,VAR_GENLIB_PARNAME);
- generator->SetMomentumRange(0,999);
- generator->SetPtRange(0,50.);
- generator->SetYRange(-4.1,-2.4);
- generator->SetPhiRange(0., 360.);
- generator->SetChildThetaRange(0.,180.);
- generator->SetForceDecay(kDiMuon);
- }
- else if ( genType == "GENCORR" ) {
- generator = new AliGenCorrHF(1, VAR_GENCORR_QUARK, VAR_GENCORR_ENERGY);
- generator->SetMomentumRange(0,9999);
- generator->SetChildThetaRange(160.0,180.0);
- generator->SetForceDecay(kSemiMuonic);
- }
-
- generator->SetCutOnChild(1);
- generator->SetChildPhiRange(0.,360.);
- generator->SetTrackingFlag(1);
-
- return generator;
-}
// MUON Tracker Residual Alignment
reco.SetSpecificStorage("MUON/Align/Data","alien://folder=/alice/simulation/2008/v4-15-Release/Residual");
- reco.Run();
+ reco.Run();
}
#include "TParameter.h"
#include "TPaveText.h"
#include "TROOT.h"
-#include "TStyle.h"
+#include "TStopwatch.h"
#include "TStyle.h"
#include "TSystem.h"
#include <cassert>
TString AliAnalysisMuMu::fgOCDBPath("raw://");
-TString AliAnalysisMuMu::fgDefaultDimuonTriggers("CMUL7-B-NOPF-MUON,CMUL7-S-NOPF-ALLNOTRD,CMUL7-S-NOPF-MUON,CMUL8-S-NOPF-MUON,CMUL7-B-NOPF-ALLNOTRD,CMUU7-B-NOPF-ALLNOTRD,CMUU7-B-NOPF-MUON,CPBI1MUL-B-NOPF-MUON,CMULLO-B-NOPF-MUON");
+TString AliAnalysisMuMu::fgDefaultDimuonTriggers("CMUL7-B-NOPF-MUON");
+
+//,CMUL7-S-NOPF-ALLNOTRD,CMUL7-S-NOPF-MUON,CMUL8-S-NOPF-MUON,CMUL7-B-NOPF-ALLNOTRD,CMUU7-B-NOPF-ALLNOTRD,CMUU7-B-NOPF-MUON,CPBI1MUL-B-NOPF-MUON,CMULLO-B-NOPF-MUON");
+
+TString AliAnalysisMuMu::fgDefaultMuonTriggers("CMSL7-S-NOPF-MUON");
-TString AliAnalysisMuMu::fgDefaultMuonTriggers("CMSL7-S-NOPF-MUON,CMSL7-S-NOPF-ALLNOTRD,CMSL8-S-NOPF-MUON,CMSL8-S-NOPF-ALLNOTRD,CMSL7-B-NOPF-MUON,CMUS1-B-NOPF-MUON,CMUS7-B-NOPF-MUON,CMSNGL-B-NOPF-MUON");
+//,CMSL7-S-NOPF-ALLNOTRD,CMSL8-S-NOPF-MUON,CMSL8-S-NOPF-ALLNOTRD,CMSL7-B-NOPF-MUON,CMUS1-B-NOPF-MUON,CMUS7-B-NOPF-MUON,CMSNGL-B-NOPF-MUON");
-TString AliAnalysisMuMu::fgDefaultMinbiasTriggers("CINT7-B-NOPF-ALLNOTRD,CINT7-S-NOPF-ALLNOTRD,CINT8-B-NOPF-ALLNOTRD,CINT8-S-NOPF-ALLNOTRD,CINT1-B-NOPF-ALLNOTRD,CPBI2_B1-B-NOPF-ALLNOTRD");
+TString AliAnalysisMuMu::fgDefaultMinbiasTriggers("CINT7-B-NOPF-ALLNOTRD");
-TString AliAnalysisMuMu::fgDefaultEventSelectionList("PSALL");
+//,CINT7-S-NOPF-ALLNOTRD,CINT8-B-NOPF-ALLNOTRD,CINT8-S-NOPF-ALLNOTRD,CINT1-B-NOPF-ALLNOTRD,CPBI2_B1-B-NOPF-ALLNOTRD");
+
+TString AliAnalysisMuMu::fgDefaultEventSelectionList("PSALL"); // for real data, for simulation see AliAnalysisMuMu ctor
TString AliAnalysisMuMu::fgDefaultPairSelectionList("pMATCHLOWRABSBOTH");
TString AliAnalysisMuMu::fgDefaultCentralitySelectionList("PP");
+//TString AliAnalysisMuMu::fgDefaultFitTypeList("PSILOW:2,PSILOWalphaLow0.984nLow5.839alphaUp1.972nUp3.444:2,PSILOWMCTAILS:2");
+TString AliAnalysisMuMu::fgDefaultFitTypeList("PSILOWalphaLow0.984nLow5.839alphaUp1.972nUp3.444:2,PSILOWMCTAILS:2");
+
+TString AliAnalysisMuMu::fgDefaultEventSelectionForSimulations("ALL");
+TString AliAnalysisMuMu::fgDefaultDimuonTriggerForSimulations("CMULLO-B-NOPF-MUON");
+
Bool_t AliAnalysisMuMu::fgIsCompactGraphs(kFALSE);
+//_____________________________________________________________________________
+TString First(const TString s)
+{
+ TString rv;
+
+ TObjArray* tokens = s.Tokenize(",");
+
+ if (!tokens) return rv;
+
+ rv = static_cast<TObjString*>(tokens->First())->String();
+
+ delete tokens;
+
+ return rv;
+}
+
//_____________________________________________________________________________
TString FindTrigger(const AliMergeableCollection& mc,
const char* base,
}
//_____________________________________________________________________________
-AliAnalysisMuMu::AliAnalysisMuMu(const char* filename) : TObject(),
+AliAnalysisMuMu::AliAnalysisMuMu(const char* filename, const char* associatedSimFileName) : TObject(),
fFilename(filename),
fCounterCollection(0x0),
fDimuonTriggers(fgDefaultDimuonTriggers),
fEventSelectionList(fgDefaultEventSelectionList),
fPairSelectionList(fgDefaultPairSelectionList),
fCentralitySelectionList(fgDefaultCentralitySelectionList),
+fFitTypeList(fgDefaultFitTypeList),
fBinning(0x0),
fMergeableCollection(0x0),
fRunNumbers(),
-fCorrectionPerRun(0x0)
+fCorrectionPerRun(0x0),
+fAssociatedSimulation(0x0)
{
// ctor
GetCollections(fFilename,fMergeableCollection,fCounterCollection,fBinning,fRunNumbers);
+
+ if (IsSimulation())
+ {
+ SetEventSelectionList("ALL");
+ SetDimuonTriggerList("CMULLO-B-NOPF-MUON");
+ SetFitTypeList("PSI1:1,COUNTJPSI:1");
+ }
+
+ if ( strlen(associatedSimFileName) )
+ {
+ fAssociatedSimulation = new AliAnalysisMuMu(associatedSimFileName);
+ }
}
//_____________________________________________________________________________
delete fBinning;
delete fMergeableCollection;
delete fCorrectionPerRun;
+ delete fAssociatedSimulation;
}
//_____________________________________________________________________________
}
+//_____________________________________________________________________________
+void AliAnalysisMuMu::CleanAllSpectra()
+{
+ /// Delete all the spectra we may have
+
+ MC()->RemoveByType("AliAnalysisMuMuSpectra");
+ Update();
+}
//_____________________________________________________________________________
void AliAnalysisMuMu::Compact(TGraph& g)
const char* trigger,
const char* eventType,
const char* pairCut,
- const char* centrality) const
+ const char* centrality,
+ const char* subresultname,
+ const char* flavour) const
{
/// Draw minv spectra for binning of given type
if (!MC() || !BIN()) return;
- TObjArray* bins = BIN()->CreateBinObjArray(particle,type);
+ TObjArray* bins = BIN()->CreateBinObjArray(particle,type,flavour);
if (!bins)
{
AliError(Form("Could not get %s bins",type));
return;
}
+ Double_t xmin(-1);
+ Double_t xmax(-1);
+
+ TString sparticle(particle);
+ if ( sparticle=="PSI" )
+ {
+ xmin = 2;
+ xmax = 6;
+ }
+
+ Int_t nx(1);
+ Int_t ny(1);
+
+ Int_t n = bins->GetEntries();
+
+ if ( n == 2 )
+ {
+ nx = 2;
+ }
+ else if ( n > 2 )
+ {
+ ny = TMath::Nint(TMath::Sqrt(n));
+ nx = n/ny;
+ }
+
+ TString stype(type);
+ stype.ToUpper();
+
+ TString spectraName(Form("/%s/%s/%s/%s/%s-%s-%s",eventType,trigger,centrality,pairCut,particle,stype.Data(),flavour));
+
+ AliAnalysisMuMuSpectra* spectra = static_cast<AliAnalysisMuMuSpectra*>(MC()->GetObject(spectraName.Data()));
+
+ AliDebug(1,Form("spectraName=%s spectra=%p",spectraName.Data(),spectra));
+
+ TObjArray* spectraBins(0x0);
+ if ( spectra )
+ {
+ spectraBins = spectra->Bins();
+ }
+
+ TCanvas* c = new TCanvas;
+ c->Divide(nx,ny);
+ c->Draw();
+ gStyle->SetOptFit(1112);
+
+ c->Connect("ProcessedEvent(Int_t,Int_t,Int_t,TObject*)", "AliAnalysisMuMu",
+ (void*)this, "ExecuteCanvasEvent(Int_t,Int_t,Int_t,TObject*)");
+
+
TIter next(bins);
AliAnalysisMuMuBinning::Range* r;
+ Int_t ci(0);
while ( ( r = static_cast<AliAnalysisMuMuBinning::Range*>(next())) )
{
TString name(Form("/%s/%s/%s/%s/MinvUS%s",eventType,trigger,centrality,pairCut,r->AsString().Data()));
AliDebug(1,name.Data());
-
+
+ AliAnalysisMuMuResult* spectraBin(0x0);
+
+ if ( spectraBins )
+ {
+ AliAnalysisMuMuResult* sr = static_cast<AliAnalysisMuMuResult*>(spectraBins->At(ci));
+
+ spectraBin = sr->SubResult(subresultname);
+
+ AliDebug(1,Form("spectraBin(%s)=%p",subresultname,spectraBin));
+ }
+
TH1* h = MC()->Histo(name.Data());
+
+ if ( spectraBin )
+ {
+ h = spectraBin->Minv();
+ }
+
if (h)
{
- TString cname(name);
- cname.ReplaceAll("/","_");
- TCanvas* c = new TCanvas(cname.Data(),cname.Data());
- c->SetLogy();
- h->Draw("histe");
+ ++ci;
+ c->cd(ci);
+ gPad->SetLogy();
+ if (xmin>0)
+ {
+ h->GetXaxis()->SetRangeUser(xmin,xmax);
+ }
+ h->Draw("histes");
+
+ TObject* f1 = h->GetListOfFunctions()->FindObject("fitTotal");
+ if (f1)
+ {
+ f1->Draw("same");
+ }
+
+ gPad->Modified();
+ gPad->Update();
+
+ TObject* stats = h->FindObject("stats");
+ if (stats)
+ {
+ stats->Draw("same");
+ }
}
}
}
//_____________________________________________________________________________
-void AliAnalysisMuMu::DrawMinv(const char* type, const char* particle) const
+void AliAnalysisMuMu::DrawMinv(const char* type, const char* particle, const char* flavour, const char* subresultname) const
{
/// Draw minv spectra for binning of given type
First(DimuonTriggerList()).Data(),
First(EventSelectionList()).Data(),
First(PairSelectionList()).Data(),
- First(CentralitySelectionList()).Data());
+ First(CentralitySelectionList()).Data(),
+ subresultname,
+ flavour);
+}
+
+//___________________________________________________________________
+void AliAnalysisMuMu::ExecuteCanvasEvent(Int_t event, Int_t /*px*/, Int_t /*py*/, TObject *sel)
+{
+ // Actions in reponse to mouse button events.
+
+ TCanvas* c = static_cast<TCanvas*>(gTQSender);
+ TPad* pad = static_cast<TPad*>(c->GetSelectedPad());
+ if (!pad) return;
+
+// if ((event == kButton1Down) ||
+ if (event == kButton1Double)
+ {
+
+// Float_t x = pad->AbsPixeltoX(px);
+// Float_t y = pad->AbsPixeltoY(py);
+// x = pad->PadtoX(x);
+// y = pad->PadtoY(y);
+
+// std::cout << "event=" << event << " px=" << px << " py=" << py << " ";
+
+ if ( sel && sel->InheritsFrom("TH1") )
+ {
+ TCanvas* clocal = new TCanvas;
+ clocal->SetLogy();
+ clocal->Draw();
+ sel->Draw();
+ }
+ else
+ {
+ TList* list = pad->GetListOfPrimitives();
+ TIter next(list);
+ TObject* h;
+
+ while ( ( h = next() ) )
+ {
+ if ( h->InheritsFrom("TH1") )
+ {
+ TCanvas* clocal = new TCanvas;
+ clocal->SetLogy();
+ clocal->Draw();
+ h->Draw();
+ break;
+ }
+ }
+
+ }
+
+// std::cout << std::endl;
+
+ pad->Modified();
+ }
+
}
//_____________________________________________________________________________
return 0x0;
}
-// AliDebug(1,Form("will project fMergeableCollection=%p...",fMergeableCollection));
-//
-// AliMergeableCollection* hp = fMergeableCollection->Project(Form("/%s/%s/%s/%s/",eventType,trigger,centrality,pairCut));
-//
-// AliDebug(1,Form("projection=hp=%p",hp));
-
-// if (!hp)
-// {
-// AliError(Form("projection %s/%s/%s/%s failed",eventType,trigger,centrality,pairCut));
-// delete bins;
-// return 0x0;
-// }
-
- binning.Print();
+// binning.Print();
AliAnalysisMuMuSpectra* spectra(0x0);
AliAnalysisMuMuBinning::Range* bin;
TIter next(bins);
- Int_t nrebin=1;
+
+ TObjArray* fitTypeArray = fFitTypeList.Tokenize(",");
+ TIter nextFitType(fitTypeArray);
+ TObjString* fitType;
+ TString flavour;
while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(next())) )
{
TString hname(Form("MinvUS%s",bin->AsString().Data()));
-// TH1* hminv = hp->Histo(hname.Data());
TH1* hminv = fMergeableCollection->Histo(Form("/%s/%s/%s/%s",eventType,trigger,centrality,pairCut),hname.Data());
if (!hminv)
r->SetNofTriggers(ntrigger);
-// r->AddFit("PSILOW",nrebin);
- if (IsSimulation())
- {
- r->AddFit("PSI1",1);
- r->AddFit("COUNTPSI",1);
- }
- else
+ nextFitType.Reset();
+
+ while ( ( fitType = static_cast<TObjString*>(nextFitType())) )
{
- r->AddFit("PSILOW",nrebin+1);
+ AliDebug(1,Form("<<<<<< fitType=%s bin=%s",fitType->String().Data(),bin->Flavour().Data()));
+
+ if ( fitType->String().BeginsWith("PSILOWMCTAILS") )
+ {
+ std::vector<Double_t> par;
+ par = GetMCCB2Tails(*bin);
+ if (!par.empty())
+ {
+ r->AddFit(fitType->String().Data(),par.size(),&par[0]);
+ }
+ }
+ else
+ {
+ r->AddFit(fitType->String());
+ }
}
-// r->AddFit("PSILOW",nrebin+2);
-// r->AddFit("PSILOW",nrebin+3);
-
-// r->Print();
+ flavour = bin->Flavour();
+
if (!spectra)
{
- spectra = new AliAnalysisMuMuSpectra(binning.GetName());
+ TString spectraName(binning.GetName());
+ if ( flavour.Length() > 0 )
+ {
+ spectraName += "-";
+ spectraName += flavour;
+ }
+ spectra = new AliAnalysisMuMuSpectra(spectraName.Data());
}
spectra->AdoptResult(*bin,r);
} // loop on bins
+ delete fitTypeArray;
delete bins;
return spectra;
}
+//_____________________________________________________________________________
+std::vector<Double_t>
+AliAnalysisMuMu::GetMCCB2Tails(const AliAnalysisMuMuBinning::Range& bin) const
+{
+ /// Get the tails from the associated simulation
+
+ std::vector<Double_t> par;
+
+ if (!SIM())
+ {
+ AliError("Cannot get MC tails without an associated simulation file !");
+ return par;
+ }
+
+
+ AliAnalysisMuMuSpectra* s = static_cast<AliAnalysisMuMuSpectra*>(SIM()->GetSpectra(bin.Type().Data(),bin.Flavour().Data()));
+
+ if (!s)
+ {
+ AliError(Form("Could not find spectra %s,%s for associated simulation",bin.Type().Data(),bin.Flavour().Data()));
+ fAssociatedSimulation->MC()->Print("*:Ali*");
+ return par;
+ }
+ else
+ {
+ AliDebug(1,Form("AliAnalysisMuMuSpectra* s = reinterpret_cast<AliAnalysisMuMuSpectra*>(%p)",s));
+ }
+
+ AliAnalysisMuMuResult* r = s->GetResultForBin(bin);
+
+ if ( r )
+ {
+ AliAnalysisMuMuResult* r1 = r->SubResult("JPSI:1");
+ if (r1)
+ {
+ TF1* func = static_cast<TF1*>(r1->Minv()->GetListOfFunctions()->FindObject("fitTotal"));
+ if (func)
+ {
+ par.push_back(func->GetParameter("alphaLow"));
+ par.push_back(func->GetParameter("nLow"));
+ par.push_back(func->GetParameter("alphaUp"));
+ par.push_back(func->GetParameter("nUp"));
+ }
+ }
+ }
+
+ return par;
+}
+
//_____________________________________________________________________________
ULong64_t AliAnalysisMuMu::GetTriggerScalerCount(const char* triggerList, Int_t runNumber)
{
return n;
}
+//_____________________________________________________________________________
+AliAnalysisMuMuSpectra* AliAnalysisMuMu::GetSpectra(const char* what, const char* flavour) const
+{
+ /// get a given spectra
+
+ TString swhat(what);
+ TString sflavour(flavour);
+ swhat.ToUpper();
+ sflavour.ToUpper();
+
+ TString spectraName(Form("/PSALL/%s/PP/%s/PSI-%s",
+ First(fgDefaultDimuonTriggers).Data(),
+ First(fgDefaultPairSelectionList).Data(),
+ swhat.Data()));
+
+ if (sflavour.Length()>0)
+ {
+ spectraName += "-";
+ spectraName += sflavour.Data();
+ }
+
+ if (IsSimulation())
+ {
+ spectraName.ReplaceAll("PSALL",fgDefaultEventSelectionForSimulations.Data());
+ spectraName.ReplaceAll(First(fgDefaultDimuonTriggers).Data(),fgDefaultDimuonTriggerForSimulations.Data());
+ }
+
+ return SPECTRA(spectraName.Data());
+}
+
//_____________________________________________________________________________
UInt_t AliAnalysisMuMu::GetSum(AliCounterCollection& cc, const char* triggerList,
const char* eventSelection, Int_t runNumber)
{
// whether or not we have MC information
- return ( fMergeableCollection->Histo("/INPUT/ALL/MinvUS") != 0x0 );
+ return ( fMergeableCollection->Histo(Form("/INPUT/%s/MinvUS",fgDefaultEventSelectionForSimulations.Data())) != 0x0 );
}
//_____________________________________________________________________________
Int_t
-AliAnalysisMuMu::Jpsi(const char* what)
+AliAnalysisMuMu::Jpsi(const char* what, const char* binningFlavour)
{
// Fit the J/psi (and psiprime) peaks for the triggers in fDimuonTriggers list
- // what="" => fit only fully integrated MinvUS
+ // what="integrated" => fit only fully integrated MinvUS
// what="pt" => fit MinvUS in pt bins
// what="y" => fit MinvUS in y bins
// what="pt,y" => fit MinvUS in (pt,y) bins
+ TStopwatch timer;
+
if (!fMergeableCollection)
{
AliError("No mergeable collection. Consider Upgrade()");
TObjArray* pairCutArray = fPairSelectionList.Tokenize(",");
TObjArray* whatArray = TString(what).Tokenize(",");
- if ( whatArray->GetEntries() <= 0 )
- {
- whatArray->Add(new TObjString(""));
- }
-
TIter nextTrigger(triggerArray);
TIter nextEventType(eventTypeArray);
TIter nextPairCut(pairCutArray);
if ( fBinning && swhat->String().Length() > 0 )
{
- binning = fBinning->Project("psi",swhat->String().Data());
+ binning = fBinning->Project("psi",swhat->String().Data(),binningFlavour);
}
else
{
delete triggerArray;
delete eventTypeArray;
delete pairCutArray;
-
+
+ timer.Print();
+
if (nfits)
{
- ReOpen(fFilename,"UPDATE");
- fMergeableCollection->Write("MC",TObjArray::kOverwrite);// | TObjArray::kSingleKey);
- ReOpen(fFilename,"READ");
+ Update();
+// ReOpen(fFilename,"UPDATE");
+// fMergeableCollection->Write("MC",TObjArray::kOverwrite);// | TObjArray::kSingleKey);
+// ReOpen(fFilename,"READ");
}
+
+
return nfits;
}
{
Double_t xcorr,ycorr;
fCorrectionPerRun->GetPoint(i,xcorr,ycorr); // note that the fact that xcorr==runNumber has been checked by the SetCorrectionPerRun method
- if ( TMath::Abs(ycorr) > 1E-12 )
- {
- y /= ycorr;
- }
- else
- {
- y = 0.0;
- }
+ y *= ycorr;
+ // FIXME: should get the correction error here
}
if ( asRejection ) y = 100*(1.0 - y);
if (f)
{
- f->Close();
delete f;
}
- f = TFile::Open(filename,"update");
+ f = TFile::Open(filename,mode);
if ( !f || !f->IsOpen() )
{
}
//_____________________________________________________________________________
-Bool_t AliAnalysisMuMu::SetCorrectionPerRun(const TGraph& corr)
+Bool_t AliAnalysisMuMu::SetCorrectionPerRun(const TGraph& corr, const char* formula)
{
/// Sets the graph used to correct values per run
delete fCorrectionPerRun;
++i;
}
- fCorrectionPerRun = static_cast<TGraph*>(corr.Clone());
+ fCorrectionPerRun = new TGraphErrors(corr.GetN());
+
+ TFormula* tformula(0x0);
+ if ( strlen(formula) > 0 )
+ {
+ tformula = new TFormula("SetCorrectionPerRunFormula",formula);
+ }
+
+ i = 0;
+
+ for ( std::set<int>::const_iterator it = RunNumbers().begin(); it != RunNumbers().end(); ++it )
+ {
+ Double_t y = corr.GetY()[i];
+
+ if ( tformula )
+ {
+ y = tformula->Eval(y);
+ }
+ fCorrectionPerRun->SetPoint(i,corr.GetX()[i],y);
+ ++i;
+ }
+ delete formula;
+
return kTRUE;
}
if (!hinput)
{
AliError(Form("Got a simulation file where I did not find histogram /INPUT/INYRANGE/%s",hname.Data()));
-
+
}
else
{
}
}
+//_____________________________________________________________________________
+AliAnalysisMuMuSpectra* AliAnalysisMuMu::SPECTRA(const char* fullpath) const
+{
+ /// Shortcut method to get to a spectra
+ if (!MC()) return 0x0;
+
+ return static_cast<AliAnalysisMuMuSpectra*>(MC()->GetObject(fullpath));
+}
+
//_____________________________________________________________________________
void AliAnalysisMuMu::TriggerCountCoverage(const char* triggerList,
Bool_t compact,
total += n;
totalExpected += expected;
- msg += TString::Format("%30s %9lld expected %9lld ",strigger->String().Data(),n,expected);
+ msg += TString::Format("%30s %9lld expected %9lld [%s] ",strigger->String().Data(),n,expected,
+ (n>expected ? "!" : " "));
if ( expected > 0 ) {
msg += TString::Format("fraction %5.1f %%",n*100.0/expected);
}
+
if (!compact)
{
msg += "\n";
{
++n;
current += it->first;
- std::cout << Form("%10lld",it->first) << " " << it->second << " percentage of total = " << Form("%7.2f %% %3d",current*100.0/total,n ) << std::endl;
+ Double_t percent = ( total > 0.0 ? current*100.0/total : 0.0);
+ std::cout << Form("%10lld",it->first) << " " << it->second << " percentage of total = " << Form("%7.2f %% %3d",percent,n ) << std::endl;
}
std::cout << Form("--- TOTAL %lld expected %lld fraction %5.1f %%",
fCorrectionPerRun=0x0;
}
+//_____________________________________________________________________________
+void AliAnalysisMuMu::Update()
+{
+ /// update the current file with memory
+
+ ReOpen(fFilename,"UPDATE");
+
+ MC()->Write("MC",TObject::kSingleKey);
+
+ ReOpen(fFilename,"READ");
+}
+
//_____________________________________________________________________________
Bool_t AliAnalysisMuMu::Upgrade(const char* filename)
{
}
//_____________________________________________________________________________
-AliAnalysisMuMuSpectra* AliAnalysisMuMu::CorrectSpectra(const char* realFile, const char* simFile, const char* particle, const char* type)
+AliAnalysisMuMuSpectra* AliAnalysisMuMu::CorrectSpectra(const char* type, const char* flavour)
{
+ /// Correct one spectra
- AliAnalysisMuMu real(realFile);
- AliAnalysisMuMu sim(simFile);
+ if (!SIM())
+ {
+ AliError("Cannot compute corrected yield without associated MC file !");
+ return 0x0;
+ }
+
+ const char* accEffSubResultName="COUNTJPSI:1";
+ AliAnalysisMuMuSpectra* realSpectra = GetSpectra(type,flavour);
+ AliAnalysisMuMuSpectra* simSpectra = SIM()->GetSpectra(type,flavour);
- AliAnalysisMuMuSpectra* realpt = static_cast<AliAnalysisMuMuSpectra*>(real.MC()->GetObject(Form("/PSALL/CMUL7-B-NOPF-MUON/PP/pMATCHLOWRABSBOTH/%s",type)));
- AliAnalysisMuMuSpectra* simpt = static_cast<AliAnalysisMuMuSpectra*>(sim.MC()->GetObject(Form("/ALL/CMULLO-B-NOPF-MUON/PP/pMATCHLOWRABSBOTH/%s",type)));
+ if ( !realSpectra )
+ {
+ AliError("could not get real spectra");
+ return 0x0;
+ }
- if ( !realpt )
+ if ( !simSpectra)
{
- AliErrorClass("could not get real spectra");
+ AliError("could not get sim spectra");
return 0x0;
}
- if ( !simpt )
+ realSpectra->Correct(*simSpectra,"Jpsi",accEffSubResultName);
+
+ Update();
+
+ return realSpectra;
+}
+
+//_____________________________________________________________________________
+AliAnalysisMuMuSpectra* AliAnalysisMuMu::ComputeYield(const char* type, const char* flavour)
+{
+ if (!SIM())
+ {
+ AliError("Cannot compute corrected yield without associated MC file !");
+ return 0x0;
+ }
+
+ const char* accEffSubResultName="COUNTJPSI:1";
+
+ AliAnalysisMuMuSpectra* realSpectra = GetSpectra(type,flavour);
+
+ if ( !realSpectra )
+ {
+ AliError("could not get real spectra");
+ return 0x0;
+ }
+
+ if (!realSpectra->HasValue("CoffNofJpsi"))
+ {
+ if (!CorrectSpectra(type,flavour))
+ {
+ AliError("Could not get corrected spectra");
+ return 0x0;
+ }
+ }
+
+ AliAnalysisMuMuSpectra* simSpectra = SIM()->GetSpectra(type,flavour);
+
+ if ( !simSpectra)
{
AliErrorClass("could not get sim spectra");
return 0x0;
}
- AliInfoClass("REAL >>>");
- realpt->Print("4");
- AliInfoClass("SIM >>>");
- simpt->Print("1");
+ Double_t nofCMUL7 = CC()->GetSum(Form("trigger:CMUL7-B-NOPF-MUON/event:PSALL"));
+ Double_t nofCINT7 = CC()->GetSum(Form("trigger:CINT7-B-NOPF-ALLNOTRD/event:PSALL"));
+ Double_t nofCINT7w0MUL = CC()->GetSum(Form("trigger:CINT7-B-NOPF-ALLNOTRD&0MUL/event:PSALL"));
- AliInfoClass("CORRECTED >>>");
-
- AliAnalysisMuMuSpectra* spectra = static_cast<AliAnalysisMuMuSpectra*>(realpt->Clone());
- spectra->Correct(*simpt,particle);
+ AliAnalysisMuMuBinning* binning = realSpectra->Binning();
+ TObjArray* bins = binning->CreateBinObjArray();
+ TIter nextBin(bins);
+ AliAnalysisMuMuBinning::Range* bin;
+ Int_t i(0);
+ AliAnalysisMuMuResult* r;
+
+ while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(nextBin()) ) )
+ {
+ r = static_cast<AliAnalysisMuMuResult*>(realSpectra->Bins()->At(i));
+
+ StdoutToAliDebug(1,std::cout << "bin=";r->Print(););
+
+ AliAnalysisMuMuResult* rsim = static_cast<AliAnalysisMuMuResult*>(simSpectra->Bins()->At(i));
+
+ Double_t mbeq = nofCINT7w0MUL / ( nofCINT7 * nofCMUL7);
+ Double_t mbeqError = mbeq * AliAnalysisMuMuResult::ErrorABC( nofCINT7w0MUL, TMath::Sqrt(nofCINT7w0MUL),
+ nofCINT7,TMath::Sqrt(nofCINT7),
+ nofCMUL7,TMath::Sqrt(nofCMUL7));
+
+ r->Set("MBR",nofCINT7/nofCINT7w0MUL,(nofCINT7/nofCINT7w0MUL)*AliAnalysisMuMuResult::ErrorAB( nofCINT7w0MUL, TMath::Sqrt(nofCINT7w0MUL),
+ nofCINT7,TMath::Sqrt(nofCINT7)));
+
+ Double_t yield = r->GetValue("CorrNofJpsi") * mbeq;
+
+ Double_t yieldError = yield * AliAnalysisMuMuResult::ErrorAB( r->GetValue("CorrNofJpsi"), r->GetErrorStat("CorrNofJpsi"),
+ mbeq,mbeqError);
+
+ r->Set("YJpsi",yield,yieldError);
+
+ r->Set("NofInputJpsi",rsim->GetValue("NofInputJpsi",accEffSubResultName),rsim->GetErrorStat("NofInputJpsi",accEffSubResultName));
+ r->Set("AccEffJpsi",rsim->GetValue("AccEffJpsi",accEffSubResultName),rsim->GetErrorStat("AccEffJpsi",accEffSubResultName));
+
+ ++i;
+ }
- spectra->Print("1");
+ delete bins;
- return spectra;
+ Update();
+
+ return realSpectra;
}
+////_____________________________________________________________________________
+//AliAnalysisMuMuSpectra* AliAnalysisMuMu::ComputeYield(const char* realFile, const char* simFile,
+// const char* type)
+//{
+// const char* accEffSubResultName="COUNTJPSI-1";
+//
+// AliAnalysisMuMu real(realFile);
+// AliAnalysisMuMu sim(simFile);
+//
+//
+// AliAnalysisMuMuSpectra* realSpectra = static_cast<AliAnalysisMuMuSpectra*>(real.MC()->GetObject(Form("/PSALL/CMUL7-B-NOPF-MUON/PP/pMATCHLOWRABSBOTH/PSI-%s",type)));
+// AliAnalysisMuMuSpectra* simSpectra = static_cast<AliAnalysisMuMuSpectra*>(sim.MC()->GetObject(Form("/ALL/CMULLO-B-NOPF-MUON/PP/pMATCHLOWRABSBOTH/PSI-%s",type)));
+//
+// if ( !realSpectra )
+// {
+// AliErrorClass("could not get real spectra");
+// return 0x0;
+// }
+//
+// if ( !simSpectra)
+// {
+// AliErrorClass("could not get sim spectra");
+// return 0x0;
+// }
+//
+// AliAnalysisMuMuSpectra* corrSpectra = static_cast<AliAnalysisMuMuSpectra*>(realSpectra->Clone());
+// corrSpectra->Correct(*simSpectra,"Jpsi",accEffSubResultName);
+//
+// Double_t nofCMUL7 = real.CC()->GetSum(Form("trigger:CMUL7-B-NOPF-MUON/event:PSALL"));
+// Double_t nofCINT7 = real.CC()->GetSum(Form("trigger:CINT7-B-NOPF-ALLNOTRD/event:PSALL"));
+// Double_t nofCINT7w0MUL = real.CC()->GetSum(Form("trigger:CINT7-B-NOPF-ALLNOTRD&0MUL/event:PSALL"));
+//
+// AliAnalysisMuMuBinning* binning = corrSpectra->Binning();
+// TObjArray* bins = binning->CreateBinObjArray();
+// TIter nextBin(bins);
+// AliAnalysisMuMuBinning::Range* bin;
+// Int_t i(0);
+// AliAnalysisMuMuResult* r;
+//
+// while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(nextBin()) ) )
+// {
+// r = static_cast<AliAnalysisMuMuResult*>(corrSpectra->Bins()->At(i));
+//
+// AliAnalysisMuMuResult* rsim = static_cast<AliAnalysisMuMuResult*>(simSpectra->Bins()->At(i));
+//
+// Double_t mbeq = nofCINT7w0MUL / ( nofCINT7 * nofCMUL7);
+// Double_t mbeqError = mbeq * AliAnalysisMuMuResult::ErrorABC( nofCINT7w0MUL, TMath::Sqrt(nofCINT7w0MUL),
+// nofCINT7,TMath::Sqrt(nofCINT7),
+// nofCMUL7,TMath::Sqrt(nofCMUL7));
+//
+// r->Set("MBR",nofCINT7/nofCINT7w0MUL,(nofCINT7/nofCINT7w0MUL)*AliAnalysisMuMuResult::ErrorAB( nofCINT7w0MUL, TMath::Sqrt(nofCINT7w0MUL),
+// nofCINT7,TMath::Sqrt(nofCINT7)));
+//
+// Double_t yield = r->GetValue("CorrNofJpsi") * mbeq;
+//
+// Double_t yieldError = yield * AliAnalysisMuMuResult::ErrorAB( r->GetValue("CorrNofJpsi"), r->GetErrorStat("CorrNofJpsi"),
+// mbeq,mbeqError);
+//
+// r->Set("YJpsi",yield,yieldError);
+//
+// r->Set("NofInputJpsi",rsim->GetValue("NofInputJpsi",accEffSubResultName),rsim->GetErrorStat("NofInputJpsi",accEffSubResultName));
+// r->Set("AccEffJpsi",rsim->GetValue("AccEffJpsi",accEffSubResultName),rsim->GetErrorStat("AccEffJpsi",accEffSubResultName));
+//
+// ++i;
+// }
+//
+// delete bins;
+//
+// return corrSpectra;
+//}
+
//_____________________________________________________________________________
AliAnalysisMuMuSpectra* AliAnalysisMuMu::RABy(const char* realFile, const char* simFile, const char* type,
const char* direction)
const Double_t ymax=TMath::Log(sqrts*1000.0/3.096916);
const Double_t tab = 0.093e-6; // nb^-1
const Double_t tabError = 0.0035E-6; // nb^-1
- const char* accEffSubResultName="COUNTJPSI-1";
+ const char* accEffSubResultName="COUNTJPSI:1";
TF1 ydist("ydist","[0]*TMath::Exp(-(x*x)/(2.0*0.39*0.39))",0.,0.5);
ydist.SetParameter(0,1.);
#include <string>
#include <TString.h>
#include <vector>
+#include "RQ_OBJECT.h"
class AliAnalysisMuMuResult;
class AliAnalysisMuMuSpectra;
class AliAnalysisMuMu : public TObject
{
+
+ RQ_OBJECT("AliAnalysisMuMu")
+
public:
- enum EColor {
+ enum EColor
+ {
kBlue=1500,
kOrange=1501,
kGreen=1502
};
- AliAnalysisMuMu(const char* filename="LHC12c_muons_AOD000_000179687.saf.root");
+ AliAnalysisMuMu(const char* filename="LHC12c_muons_AOD000_000179687.saf.root",
+ const char* associatedSimFileName="");
+
virtual ~AliAnalysisMuMu();
/* Basic checks */
const char* centrality,
const AliAnalysisMuMuBinning& binning);
+ AliAnalysisMuMuSpectra* CorrectSpectra(const char* type, const char* flavour="");
+
+ AliAnalysisMuMuSpectra* ComputeYield(const char* type, const char* flavour="");
+
+ void CleanAllSpectra();
+
///------
- static AliAnalysisMuMuSpectra* CorrectSpectra(const char* realFile="ds.list.saf.root", const char* simFile="ds.sim.list.saf.root", const char* particle="Jpsi", const char* type="PSI-PT");
+// static AliAnalysisMuMuSpectra* ComputeYield(const char* realFile="ds.list.saf.root",
+// const char* simFile="ds.sim.list.saf.root",
+// const char* type="PSI-Y VS PT");
static AliAnalysisMuMuSpectra* RABy(const char* realFile="ds.list.saf.root", const char* simFile="ds.sim.list.saf.root", const char* type="", const char* direction="pPb");
AliAnalysisMuMuBinning*& bin,
std::set<int>& runnumbers);
+ AliAnalysisMuMuSpectra* GetSpectra(const char* what, const char* flavour="") const;
+
static UInt_t GetSum(AliCounterCollection& cc, const char* triggerList, const char* eventSelection, Int_t runNumber=-1);
static ULong64_t GetTriggerScalerCount(const char* triggerList, Int_t runNumber);
- Int_t Jpsi(const char* what="");
+ Int_t Jpsi(const char* what="integrated", const char* binningFlavour="");
Bool_t IsSimulation() const;
TString CentralitySelectionList() const { return fCentralitySelectionList; }
void SetCentralitySelectionList(const char* centralitySelectionList) { fCentralitySelectionList = centralitySelectionList; }
+ TString FitTypeList() const { return fFitTypeList; }
+ void SetFitTypeList(const char* fitTypelist) { fFitTypeList = fitTypelist; }
+
static void SetDefaultDimuonTriggerList(const char* dimuonTriggerList) { fgDefaultDimuonTriggers = dimuonTriggerList; }
static void SetDefaultMuonTriggerList(const char* muonTriggerList) { fgDefaultMuonTriggers = muonTriggerList; }
static void SetDefaultMinbiasTriggerList(const char* minbiasTriggerList) { fgDefaultMinbiasTriggers = minbiasTriggerList; }
static void SetDefaultEventSelectionList(const char* eventSelectionList) { fgDefaultEventSelectionList = eventSelectionList; }
static void SetDefaultPairSelectionList(const char* pairSelectionList) { fgDefaultPairSelectionList = pairSelectionList; }
static void SetDefaultCentralitySelectionList(const char* centralitySelectionList) { fgDefaultCentralitySelectionList = centralitySelectionList; }
+ static void SetDefaultFitTypes(const char* fitTypes) { fgDefaultFitTypeList = fitTypes; }
-
+ static void SetDefaultEventSelectionForSimulations(const char* value) { fgDefaultEventSelectionForSimulations = value; }
+ static void SetDefaultDimuonTriggerForSimulations(const char* value) { fgDefaultDimuonTriggerForSimulations = value; }
+
AliMergeableCollection* MC() const { return fMergeableCollection; }
AliCounterCollection* CC() const { return fCounterCollection; }
AliAnalysisMuMuBinning* BIN() const { return fBinning; }
const char* trigger,
const char* eventType,
const char* pairCut,
- const char* centrality) const;
+ const char* centrality,
+ const char* subresultname="",
+ const char* flavour="") const;
- void DrawMinv(const char* type="PT", const char* particle="PSI") const;
+ void DrawMinv(const char* type="PT", const char* particle="PSI", const char* flavour="", const char* subresultname="") const;
- Bool_t SetCorrectionPerRun(const TGraph& corr);
+ Bool_t SetCorrectionPerRun(const TGraph& corr, const char* formula="");
void UnsetCorrectionPerRun();
+ void ExecuteCanvasEvent(Int_t event, Int_t px, Int_t py, TObject *sel);
+
+ std::vector<Double_t> GetMCCB2Tails(const AliAnalysisMuMuBinning::Range& bin) const;
+
+ AliAnalysisMuMu* SIM() const { return fAssociatedSimulation; }
+
+ AliAnalysisMuMuSpectra* SPECTRA(const char* fullpath) const;
+
+ void Update();
+
private:
AliAnalysisMuMu(const AliAnalysisMuMu& rhs); // not implemented on purpose
AliAnalysisMuMu& operator=(const AliAnalysisMuMu& rhs); // not implemented on purpose
static TString ExpandPathName(const char* file);
+
TString fFilename; // file containing the result collections (of objects and counters) from AliAnalysisTaskMuMu
AliCounterCollection* fCounterCollection; // collection of counters in file
TString fDimuonTriggers; // list of dimuon triggers to consider
TString fEventSelectionList; // list of event types to consider
TString fPairSelectionList; // list of pair cuts to consider
TString fCentralitySelectionList; // list of centrality cuts to consider
-
+ TString fFitTypeList; // list of fit types to perform
+
AliAnalysisMuMuBinning* fBinning; // binning
static TString fgOCDBPath; // OCDB to be used (raw:// by default)
static TString fgDefaultEventSelectionList; // default list of event selections
static TString fgDefaultPairSelectionList; // default list of pair selections
static TString fgDefaultCentralitySelectionList; // default list of centrality selections
+ static TString fgDefaultFitTypeList; // default list of fit types
+ static TString fgDefaultEventSelectionForSimulations; // default event selection (simulations)
+ static TString fgDefaultDimuonTriggerForSimulations; // default dimuon trigger (simulations)
static Bool_t fgIsCompactGraphs; // whether the graph produced should be compact
TGraph* fCorrectionPerRun; // correction factor per run
- ClassDef(AliAnalysisMuMu,10) // class to analysis results from AliAnalysisTaskMuMuXXX tasks
+ AliAnalysisMuMu* fAssociatedSimulation; // associated simulations (if any)
+
+ ClassDef(AliAnalysisMuMu,11) // class to analysis results from AliAnalysisTaskMuMuXXX tasks
};
#endif
const std::map<std::string,Double_t>& MassMap()
{
+ /// a simple map of masses...
static std::map<std::string,Double_t> massMap;
// not super elegant, but that way we do not depend on TDatabasePDG and thus
// can decide on our particle naming
}
-Double_t ErrorAB(Double_t a, Double_t aerr, Double_t b, Double_t berr)
+Double_t funcCB(Double_t* xx, Double_t* par)
{
- Double_t e(0.0);
-
- if ( TMath::Abs(a) > 1E-12 )
- {
- e += (aerr*aerr)/(a*a);
- }
+ /// Crystal ball
- if ( TMath::Abs(b) > 1E-12 )
- {
- e += (berr*berr)/(b*b);
- }
-
- return TMath::Sqrt(e);
-}
-
-Double_t funcCB(Double_t* xx, Double_t* par)
-{
- Double_t N = par[0];
+ Double_t norm = par[0];
Double_t alpha = par[1];
Double_t n = par[2];
Double_t mean = par[3];
Double_t x = xx[0];
- Double_t A = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
- Double_t B = n/TMath::Abs(alpha) - TMath::Abs(alpha);
+ Double_t a = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
+ Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
Double_t y = ( TMath::Abs(sigma) > 1E-12 ? (x-mean)/sigma : 0 );
if ( y > alpha*-1.0 )
{
- return N*TMath::Exp(-0.5*y*y);
+ return norm*TMath::Exp(-0.5*y*y);
}
else
{
- return N*A*TMath::Power(B-y,-n);
+ return norm*a*TMath::Power(b-y,-n);
}
}
Double_t funcJpsiGCBE(Double_t* xx, Double_t* par)
{
+ /// crystal ball + expo + gaussian
Double_t x = xx[0];
Double_t g = par[0]*TMath::Gaus(x,par[1],par[2]);
}
Double_t funcJpsiPsiPrimeCustom(Double_t* xx, Double_t* par)
-{
- Double_t N = par[0];
+{
+ // custom fit for jpsi + psi prime
+
+ Double_t norm = par[0];
Double_t alpha = par[1];
Double_t n = par[2];
Double_t mean = par[3];
Double_t x = xx[0];
- Double_t A = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
- Double_t B = n/TMath::Abs(alpha) - TMath::Abs(alpha);
- Double_t C = TMath::Power(nprime/TMath::Abs(alphaprime),nprime)*TMath::Exp(-0.5*alphaprime*alphaprime);
- Double_t D = nprime/TMath::Abs(alphaprime) - TMath::Abs(alphaprime);
+ Double_t a = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
+ Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
+ Double_t c = TMath::Power(nprime/TMath::Abs(alphaprime),nprime)*TMath::Exp(-0.5*alphaprime*alphaprime);
+ Double_t d = nprime/TMath::Abs(alphaprime) - TMath::Abs(alphaprime);
Double_t y = ( TMath::Abs(sigma) > 1E-12 ? (x-mean)/sigma : 0 );
if ( y > alphaprime )
{
- cb = N*C*TMath::Power(D+y,-nprime);
+ cb = norm*c*TMath::Power(d+y,-nprime);
}
else if ( y > alpha*-1.0 )
{
- cb = N*TMath::Exp(-0.5*y*y);
+ cb = norm*TMath::Exp(-0.5*y*y);
}
else
{
- cb = N*A*TMath::Power(B-y,-n);
+ cb = norm*a*TMath::Power(b-y,-n);
}
if ( x < mean )
Double_t funcCB2(Double_t* xx, Double_t* par)
-{
- Double_t N = par[0];
+{
+ /// CB2 = extended crystal ball
+
+ Double_t norm = par[0];
Double_t alpha = par[1];
Double_t n = par[2];
Double_t mean = par[3];
Double_t x = xx[0];
- Double_t A = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
- Double_t B = n/TMath::Abs(alpha) - TMath::Abs(alpha);
- Double_t C = TMath::Power(nprime/TMath::Abs(alphaprime),nprime)*TMath::Exp(-0.5*alphaprime*alphaprime);
- Double_t D = nprime/TMath::Abs(alphaprime) - TMath::Abs(alphaprime);
+ Double_t a = TMath::Power(n/TMath::Abs(alpha),n)*TMath::Exp(-0.5*alpha*alpha);
+ Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
+ Double_t c = TMath::Power(nprime/TMath::Abs(alphaprime),nprime)*TMath::Exp(-0.5*alphaprime*alphaprime);
+ Double_t d = nprime/TMath::Abs(alphaprime) - TMath::Abs(alphaprime);
Double_t y = ( TMath::Abs(sigma) > 1E-12 ? (x-mean)/sigma : 0 );
if ( y > alphaprime )
{
- return N*C*TMath::Power(D+y,-nprime);
+ return norm*c*TMath::Power(d+y,-nprime);
}
else if ( y > alpha*-1.0 )
{
- return N*TMath::Exp(-0.5*y*y);
+ return norm*TMath::Exp(-0.5*y*y);
}
else
{
- return N*A*TMath::Power(B-y,-n);
+ return norm*a*TMath::Power(b-y,-n);
}
}
Double_t funcJpsiPsiPrime(Double_t* xx, Double_t* par)
{
+ /// CB + CB2
+
Double_t jpsi = funcCB(xx,par);
Double_t psiprime = funcCB2(xx,par+5);
Double_t funcJpsiPCBE(Double_t* xx, Double_t* par)
{
+ // CB + expo + pol2
+
Double_t x = xx[0];
Double_t pol2 = par[0] + par[1]*x + par[2]*x*x;
return e1+e2+jpsi;
}
+Double_t funcJpsiNA48(Double_t*x, Double_t* par)
+{
+ /// Fit function from e.q. 4.8 of Ruben's PhD.
+ Double_t c1 = par[0];
+ Double_t c2 = par[1];
+ Double_t d1 = par[2];
+ Double_t d2 = par[3];
+ Double_t g1 = par[4];
+ Double_t g2 = par[5];
+ Double_t m0 = par[6];
+ Double_t sigma1 = par[7];
+ Double_t sigma2 = par[8];
+ Double_t b1 = par[9];
+ Double_t b2 = par[10];
+ Double_t norm = par[11];
+
+ Double_t m = x[0];
+
+ Double_t rv(0);
+
+ if ( m <= c1*m0 )
+ {
+ Double_t e = d1-g1*TMath::Sqrt(c1*m0-m);
+ rv = TMath::Power(b1*(c1*m0-m),e);
+ rv += sigma1;
+ }
+ else if( m >= c1*m0 && m <= m0 )
+ {
+ rv = sigma1;
+ }
+ else if ( m >= m0 && m < c2*m0 )
+ {
+ rv = sigma2;
+ }
+ else if( m >= c2*m0 )
+ {
+ Double_t e = d2-g2*TMath::Sqrt(m-c2*m0);
+ rv = TMath::Power(b2*(m-c2*m0),e);
+ rv += sigma2;
+ }
+
+ return norm*TMath::Exp(-(m-m0)*(m-m0)/(2.0*rv*rv));
+}
+
const char* NormalizeName(const char* name, const char* suffix)
{
+ /// Remove - and / from the name, and adds _suffix
+
TString str(Form("%s_%s",name,suffix));
str.ReplaceAll("-","_");
//------------------------------------------------------------------------------
Double_t func2CB2VWG(Double_t *x, Double_t *par)
{
+ /// 2 extended crystal balls + variable width gaussian
+ /// width of the second CB related to the first (free) one.
+
Double_t par2[7] = {
par[11],
3.68609+(par[5]-3.096916)/3.096916*3.68609,
fMap(0x0),
fMother(0x0),
fKeys(0x0),
-fWeight(0.0)
+fWeight(0.0),
+fRebin(0),
+fTriggerClass(),
+fEventSelection(),
+fPairSelection(),
+fCentralitySelection(),
+fAlias()
{
}
fMap(0x0),
fMother(0x0),
fKeys(0x0),
-fWeight(0.0)
-
+fWeight(0.0),
+fRebin(0),
+fTriggerClass(),
+fEventSelection(),
+fPairSelection(),
+fCentralitySelection(),
+fAlias()
{
SetMinv(hminv);
}
const char* fitType,
Int_t nrebin)
:
-TNamed(Form("%s-%d",fitType,nrebin),""),
+TNamed(Form("%s:%d",fitType,nrebin),""),
fNofRuns(1),
fNofTriggers(-1),
fMinv(0x0),
fMap(0x0),
fMother(0x0),
fKeys(0x0),
-fWeight(0.0)
-
+fWeight(0.0),
+fRebin(nrebin),
+fTriggerClass(),
+fEventSelection(),
+fPairSelection(),
+fCentralitySelection(),
+fAlias()
{
SetMinv(hminv);
}
fMap(0x0),
fMother(0x0),
fKeys(0x0),
-fWeight(0.0)
-
+fWeight(0.0),
+fRebin(1),
+fTriggerClass(triggerName),
+fEventSelection(eventSelection),
+fPairSelection(pairSelection),
+fCentralitySelection(centSelection),
+fAlias()
{
SetMinv(hminv);
}
fMap(0x0),
fMother(0x0),
fKeys(0x0),
-fWeight(rhs.fWeight)
+fWeight(rhs.fWeight),
+fRebin(rhs.fRebin),
+fAlias()
{
/// copy ctor
/// Note that the mother is lost
fMap = static_cast<TMap*>(rhs.fMap->Clone());
}
+ if ( rhs.fAlias.Length() > 0 )
+ {
+ fAlias = rhs.fAlias;
+ }
}
//_____________________________________________________________________________
AliAnalysisMuMuResult& AliAnalysisMuMuResult::operator=(const AliAnalysisMuMuResult& rhs)
{
+ /// Assignment operator
+
if (this!=&rhs)
{
delete fMinv;
fNofTriggers = rhs.NofTriggers();
fBin = rhs.Bin();
fWeight = rhs.fWeight;
+ fRebin = rhs.fRebin;
+
+ fAlias="";
+
+ if ( rhs.fAlias.Length() > 0 )
+ {
+ fAlias = rhs.fAlias;
+ }
}
return *this;
//_____________________________________________________________________________
AliAnalysisMuMuResult::~AliAnalysisMuMuResult()
{
+ // dtor
delete fMap;
delete fMinv;
delete fSubResults;
//_____________________________________________________________________________
const AliAnalysisMuMuBinning::Range& AliAnalysisMuMuResult::Bin() const
{
+ /// Get the bin of this result
if ( !Mother() ) return fBin;
else return Mother()->Bin();
}
//_____________________________________________________________________________
TObject* AliAnalysisMuMuResult::Clone(const char* /*newname*/) const
{
+ /// Clone this result
return new AliAnalysisMuMuResult(*this);
}
return kTRUE;
}
-
+ else
+ {
+ AliError(Form("Result does not have Nof%s : cannot correct it !",particle));
+ }
return kFALSE;
}
//_____________________________________________________________________________
AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsiGCBE(TH1& h)
{
+ /// Fit Jpsi spectra with crystal ball + gaussian + exponential
+
std::cout << "Fit with jpsi alone (gaus + CB + expo)" << std::endl;
Int_t nrebin = fMinv->GetXaxis()->GetNbins() / h.GetXaxis()->GetNbins();
//_____________________________________________________________________________
AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsi(TH1& h)
{
+ /// Fit Jpsi spectra using extended crystall ball (CB2) with free tails
+
std::cout << "Fit with jpsi alone" << std::endl;
Int_t nrebin = fMinv->GetXaxis()->GetNbins() / h.GetXaxis()->GetNbins();
const Double_t xmax(8.0);
TF1* fitTotal = new TF1("fitTotal",funcCB2,xmin,xmax,7);
- fitTotal->SetParNames("N","alpha","n","mean","sigma","alphaprime","nprime");
- fitTotal->SetParameters(h.GetMaximum(),1,5,3.0,0.07,1.5,3);
+ fitTotal->SetParNames("N","alphaLow","nLow","mean","sigma","alphaUp","nUp");
+ fitTotal->SetParameters(h.GetMaximum(),1,5,3.1,0.07,1.5,3);
fitTotal->SetParLimits(0,0,h.GetMaximum()*2); // N
fitTotal->SetParLimits(1,0,10); // alpha
- fitTotal->SetParLimits(2,1,10); // n
- fitTotal->SetParLimits(3,1,4); // mean
+ fitTotal->SetParLimits(2,0.1,10); // n
+ fitTotal->SetParLimits(3,3,3.1); // mean
fitTotal->SetParLimits(4,0.01,1); // sigma
fitTotal->SetParLimits(5,0,10); // alpha
- fitTotal->SetParLimits(6,1,10); // n
+ fitTotal->SetParLimits(6,0.1,10); // n
hfit->Fit(fitTotal,"+","",2,5);
}
//_____________________________________________________________________________
-AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsi2CB2VWG(const TH1& h)
+AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsiCB2VWG(const TH1& h)
{
- std::cout << "Fit with jpsi + psiprime VWG" << std::endl;
+ /// Fit Jpsi spectra using extended crystal ball (CB2) + variable width gaussian (VWG)
+
+ std::cout << "Fit with jpsi VWG" << std::endl;
Int_t nrebin = fMinv->GetXaxis()->GetNbins() / h.GetXaxis()->GetNbins();
- AliAnalysisMuMuResult* r = new AliAnalysisMuMuResult(h,"JPSI2CB2VWG",nrebin);
+ AliAnalysisMuMuResult* r = new AliAnalysisMuMuResult(h,"JPSICB2VWG",nrebin);
+ TH1* hfit = r->Minv();
+
+ const Double_t xmin(2.0);
+ const Double_t xmax(5.0);
+
+// // gaussian variable width
+// Double_t sigma = par[2]+par[3]*((x[0]-par[1])/par[1]);
+// return par[0]*TMath::Exp(-(x[0]-par[1])*(x[0]-par[1])/(2.*sigma*sigma));
+// Double_t CrystalBallExtended(Double_t *x,Double_t *par)
+// //par[0] = Normalization
+// //par[1] = mean
+// //par[2] = sigma
+// //par[3] = alpha
+// //par[4] = n
+// //par[5] = alpha'
+// //par[6] = n'
+
+ TF1* fitTotal = new TF1("fitTotal",fitFunctionCB2VWG,xmin,xmax,11);
+ fitTotal->SetParNames("kVWG","mVWG","sVWG1","sVWG2","norm","mean","sigma","alpha","n","alpha'","n'");
+
+ fitTotal->SetParameter(0, 10000.); // kVWG
+ fitTotal->SetParameter(1, 1.9); // mVWG
+
+ fitTotal->SetParameter(2, 0.5); // sVWG1
+ fitTotal->SetParLimits(2, 0., 100.);
+
+ fitTotal->SetParameter(3, 0.3); // sVWG2
+ fitTotal->SetParLimits(3, 0., 100.);
+
+ fitTotal->SetParameter(4, h.GetMaximum()); // norm
+
+ fitTotal->SetParameter(5, 3.1); // mean
+ fitTotal->SetParLimits(5, 3.0, 3.2);
+
+ fitTotal->SetParameter(6, 0.08); // sigma
+ fitTotal->SetParLimits(6, 0.04, 0.20);
+
+ fitTotal->SetParameter(7,1.0); // alpha
+ fitTotal->SetParameter(8,5); // n
+ fitTotal->SetParameter(9,2.0); // alpha'
+ fitTotal->SetParameter(10,4); // n'
+
+// fitTotal->FixParameter(7, 0.93);
+// fitTotal->FixParameter(8, 5.59);
+// fitTotal->FixParameter(9, 2.32);
+// fitTotal->FixParameter(10, 3.39);
+// fitTotal->SetParameter(11, 10.);
+
+ const char* fitOption = "SIER"; //+";
+
+ TFitResultPtr fitResult = hfit->Fit(fitTotal,fitOption,"");
+
+ r->Set("MeanJpsi",fitTotal->GetParameter(5),fitTotal->GetParError(5));
+ r->Set("SigmaJpsi",fitTotal->GetParameter(6),fitTotal->GetParError(6));
+
+ double m = r->GetValue("MeanJpsi");
+ double s = r->GetValue("SigmaJpsi");
+ double n = 3.0;
+
+ TF1* fcb = new TF1("fcb",CrystalBallExtended,xmin,xmax,7);
+ fcb->SetParameters(fitTotal->GetParameter(4),
+ fitTotal->GetParameter(5),
+ fitTotal->GetParameter(6),
+ fitTotal->GetParameter(7),
+ fitTotal->GetParameter(8),
+ fitTotal->GetParameter(9),
+ fitTotal->GetParameter(10));
+
+
+ fcb->SetLineColor(1);
+ fcb->SetNpx(1000);
+ TLine* l1 = new TLine(m-n*s,0,m-n*s,fitTotal->GetParameter(4));
+ TLine* l2 = new TLine(m+n*s,0,m+n*s,fitTotal->GetParameter(4));
+ l1->SetLineColor(6);
+ l2->SetLineColor(6);
+ hfit->GetListOfFunctions()->Add(fcb);
+ hfit->GetListOfFunctions()->Add(l1);
+ hfit->GetListOfFunctions()->Add(l2);
+
+ Double_t cbParameters[7];
+ Double_t covarianceMatrix[7][7];
+
+ for ( int ix = 0; ix < 7; ++ix )
+ {
+ cbParameters[ix] = fitTotal->GetParameter(ix+4);
+ }
+
+ for ( int iy = 0; iy < 5; ++iy )
+ {
+ for ( int ix = 0; ix < 5; ++ix )
+ {
+ covarianceMatrix[ix][iy] = (fitResult->GetCovarianceMatrix())(ix+4,iy+4);
+ }
+ }
+
+ double njpsi = fcb->Integral(m-n*s,m+n*s)/h.GetBinWidth(1);
+ double nerr = fcb->IntegralError(m-n*s,m+n*s,&cbParameters[0],&covarianceMatrix[0][0])/h.GetBinWidth(1);
+
+ r->Set("NofJpsi",njpsi,nerr);
+
+ return r;
+}
+
+//_____________________________________________________________________________
+AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsi2CB2VWG(const TH1& h,
+ Double_t alphaLow,
+ Double_t nLow,
+ Double_t alphaUp,
+ Double_t nUp)
+{
+ /// Fit using extended crystal ball + variable width gaussian
+
+ std::cout << Form("Fit with jpsi + psiprime VWG alphaLow=%5.2f nLow=%5.2f alphaUp=%5.2f nUp=%5.2f",
+ alphaLow,nLow,alphaUp,nUp) << std::endl;
+
+ Int_t nrebin = fMinv->GetXaxis()->GetNbins() / h.GetXaxis()->GetNbins();
+
+ TString resultName("JPSI2CB2VWG");
+ if ( alphaLow > 0 )
+ {
+ resultName += TString::Format("alphaLow=%5.2f",alphaLow);
+ }
+ if ( nLow > 0 )
+ {
+ resultName += TString::Format("nLow=%5.2f",nLow);
+ }
+ if ( alphaUp > 0 )
+ {
+ resultName += TString::Format("alphaUp=%5.2f",alphaUp);
+ }
+ if ( nUp > 0 )
+ {
+ resultName += TString::Format("nUp=%5.2f",nUp);
+ }
+ resultName.ReplaceAll(" ","");
+
+ AliAnalysisMuMuResult* r = new AliAnalysisMuMuResult(h,resultName.Data(),nrebin);
+
TH1* hfit = r->Minv();
const Double_t xmin(2.2);
fitTotal->SetParameter(3, 0.3);
fitTotal->SetParLimits(3, 0., 100.);
fitTotal->SetParameter(4, 100.);
- fitTotal->SetParameter(5, 3.13);
+ fitTotal->SetParameter(5, 3.1);
fitTotal->SetParLimits(5, 3.08, 3.2);
fitTotal->SetParameter(6, 0.08);
fitTotal->SetParLimits(6, 0.05, 0.15);
- fitTotal->FixParameter(7, 0.93);
- fitTotal->FixParameter(8, 5.59);
- fitTotal->FixParameter(9, 2.32);
- fitTotal->FixParameter(10, 3.39);
+
+// r = FitJpsi2CB2VWG(*hminv,0.93,5.59,2.32,3.39);
+
+ if ( alphaLow > 0 )
+ {
+ fitTotal->FixParameter(7, alphaLow);
+ }
+ else
+ {
+ fitTotal->SetParameter(7,0.9);
+ fitTotal->SetParLimits(7,0.1,10.0);
+ }
+
+ if ( nLow > 0 )
+ {
+ fitTotal->FixParameter(8, nLow);
+ }
+ else
+ {
+ fitTotal->SetParameter(8,5.0);
+ fitTotal->SetParLimits(8,0.0,10.0);
+ }
+
+ if ( alphaUp > 0 )
+ {
+ fitTotal->FixParameter(9, alphaUp);
+ }
+ else
+ {
+ fitTotal->SetParameter(9, 2.0);
+ fitTotal->SetParLimits(9,0.1,10.0);
+ }
+
+ if ( nUp > 0 )
+ {
+ fitTotal->FixParameter(10, nUp);
+ }
+ else
+ {
+ fitTotal->SetParameter(10,3.0);
+ fitTotal->SetParLimits(10,0.0,10.0);
+ }
+
fitTotal->SetParameter(11, 10.);
- const char* fitOption = "SIER"; //+";
+ const char* fitOption = "SER"; //+";
TFitResultPtr fitResult = hfit->Fit(fitTotal,fitOption,"");
return r;
}
+//_____________________________________________________________________________
+AliAnalysisMuMuResult* AliAnalysisMuMuResult::FitJpsiNA48(const TH1& h)
+{
+ /// fit using functional form from Ruben Shahoyan's thesis (2001) (eq. 4.8.)
+
+ std::cout << "Fit with jpsi NA50 Ruben eq. 4.8" << std::endl;
+
+ Int_t nrebin = fMinv->GetXaxis()->GetNbins() / h.GetXaxis()->GetNbins();
+
+ AliAnalysisMuMuResult* r = new AliAnalysisMuMuResult(h,"JPSINA",nrebin);
+
+ TH1* hfit = r->Minv();
+
+ const Double_t xmin(2.0);
+ const Double_t xmax(5.0);
+
+ TF1* fitTotal = new TF1("fitTotal",funcJpsiNA48,xmin,xmax,12);
+
+ fitTotal->SetParName( 0, "c1");
+ fitTotal->FixParameter(0,0.97);
+
+ fitTotal->SetParName( 1, "c2");
+ fitTotal->FixParameter(1,1.05);
+
+ fitTotal->SetParName( 2, "d1");
+ fitTotal->SetParameter(2,0.0);
+ fitTotal->SetParLimits(2,0,1);
+
+ fitTotal->SetParName( 3, "d2");
+ fitTotal->SetParameter(3,0.0);
+ fitTotal->SetParLimits(3,0,1);
+
+ fitTotal->SetParName( 4, "g1");
+ fitTotal->SetParameter(4,0.0);
+ fitTotal->SetParLimits(4,0.0,1);
+
+ fitTotal->SetParName( 5, "g2");
+ fitTotal->SetParameter(5,0.0);
+ fitTotal->SetParLimits(5,0.0,1);
+
+ fitTotal->SetParName( 6, "m0");
+ fitTotal->SetParameter(6,3.1);
+ fitTotal->SetParLimits(6,2.8,3.4);
+
+ fitTotal->SetParName( 7, "sigma1");
+ fitTotal->SetParameter(7,0.05);
+ fitTotal->SetParLimits(7,0.01,0.2);
+
+ fitTotal->SetParName( 8, "sigma2");
+ fitTotal->SetParameter(8,0.05);
+ fitTotal->SetParLimits(8,0.01,0.2);
+
+ fitTotal->SetParName( 9, "b1");
+ fitTotal->SetParameter(9,1.0);
+ fitTotal->SetParLimits(9,0,1);
+
+ fitTotal->SetParName(10, "b2");
+ fitTotal->SetParameter(10,1.0);
+ fitTotal->SetParLimits(10,0,1);
+
+ fitTotal->SetParName(11, "norm");
+ fitTotal->SetParameter(11,h.GetMaximum());
+
+ const char* fitOption = "SIER"; //+";
+
+ TFitResultPtr fitResult = hfit->Fit(fitTotal,fitOption,"");
+
+ r->Set("MeanJpsi",fitTotal->GetParameter(6),fitTotal->GetParError(6));
+ r->Set("SigmaJpsi",
+ fitTotal->GetParameter(7)+fitTotal->GetParameter(8),
+ 0.0);
+
+ double m = r->GetValue("MeanJpsi");
+ double s = r->GetValue("SigmaJpsi");
+ double n = 3.0;
+
+ TLine* l1 = new TLine(m-n*s,0,m-n*s,fitTotal->GetParameter(11));
+ TLine* l2 = new TLine(m+n*s,0,m+n*s,fitTotal->GetParameter(11));
+ l1->SetLineColor(6);
+ l2->SetLineColor(6);
+ hfit->GetListOfFunctions()->Add(l1);
+ hfit->GetListOfFunctions()->Add(l2);
+
+ double njpsi = fitTotal->Integral(m-n*s,m+n*s)/h.GetBinWidth(1);
+ double nerr = fitTotal->IntegralError(m-n*s,m+n*s)/h.GetBinWidth(1);
+
+ r->Set("NofJpsi",njpsi,nerr);
+
+ return r;
+}
+
/*
//_____________________________________________________________________________
void AliAnalysisMuMuResult::FitUpsilon(TH1& h)
//_____________________________________________________________________________
Double_t AliAnalysisMuMuResult::ErrorAB(Double_t a, Double_t aerr, Double_t b, Double_t berr)
{
+ /// Compute the quadratic sum of 2 errors
+
Double_t e(0.0);
if ( TMath::Abs(a) > 1E-12 )
//_____________________________________________________________________________
Double_t AliAnalysisMuMuResult::ErrorABC(Double_t a, Double_t aerr, Double_t b, Double_t berr, Double_t c, Double_t cerror)
{
+ /// Compute the quadratic sum of 3 errors
+
Double_t e(0.0);
if ( TMath::Abs(a) > 1E-12 )
//_____________________________________________________________________________
-Bool_t AliAnalysisMuMuResult::AddFit(const char* fitType, Int_t nrebin)
+Bool_t AliAnalysisMuMuResult::AddFit(const char* fitType, Int_t npar, Double_t* par)
{
- // new TCanvas;
+ // Add a fit to this result
+
+ TString msg(Form("fitType=%s npar=%d par[]=",fitType,npar));
+
+ for ( Int_t i = 0; i < npar; ++i )
+ {
+ msg += TString::Format("%e,",par[i]);
+ }
+
+ msg += TString::Format(" minv=%p %d",fMinv,fMinv?TMath::Nint(fMinv->GetEntries()):0);
if ( !fMinv ) return kFALSE;
TString sFitType(fitType);
sFitType.ToUpper();
+ Int_t nrebin(1);
+
+ if (sFitType.CountChar(':'))
+ {
+ Int_t index = sFitType.Index(":");
+ nrebin = TString(sFitType(index+1,sFitType.Length()-index-1)).Atoi();
+ sFitType = sFitType(0,index);
+ }
+
+ msg += TString::Format(" nrebin=%d",nrebin);
+ AliDebug(1,msg.Data());
+
+
if ( fMinv->GetEntries()<100 && !sFitType.Contains("COUNT")) return kFALSE;
TH1* hminv = static_cast<TH1*>(fMinv->Clone());
AliAnalysisMuMuResult* r(0x0);
-// if ( sFitType == "PSI12" )
-// {
-// r = FitJpsiPsiPrimeCustom(*hminv);
-// }
-//
if ( sFitType=="PSI1")
- {//
+ {
r = FitJpsi(*hminv);
-// r = FitJpsiGCBE(*hminv);
}
-
- if ( sFitType == "PSILOW")
+ else if ( sFitType == "PSILOW")
{
- r = FitJpsi2CB2VWG(*hminv);
+ r = FitJpsi2CB2VWG(*hminv,-1,-1,-1,-1); // free tails
}
-
- if ( sFitType == "COUNTPSI" )
+ else if ( sFitType == "PSILOWMCTAILS" )
+ {
+ if ( npar!= 4 )
+ {
+ AliError("Cannot use PSILOWMCTAILS without being given the MC tails !");
+ delete hminv;
+ return kFALSE;
+ }
+ r = FitJpsi2CB2VWG(*hminv,par[0],par[1],par[2],par[3]);
+ if (r)
+ {
+ r->SetAlias("MCTAILS");
+ }
+ }
+ else if ( sFitType.BeginsWith("PSILOWALPHA") )
+ {
+ Float_t lpar[] = { 0.0, 0.0, 0.0, 0.0 };
+
+ AliDebug(1,Form("sFitType=%s",sFitType.Data()));
+
+ sscanf(sFitType.Data(),"PSILOWALPHALOW%fNLOW%fALPHAUP%fNUP%f",
+ &lpar[0],&lpar[1],&lpar[2],&lpar[3]);
+
+ AliDebug(1,Form("PSILOW ALPHALOW=%f NLOW=%f ALPHAUP=%f NUP=%f",lpar[0],lpar[1],lpar[2],lpar[3]));
+
+ if ( lpar[0] == 0.0 && lpar[1] == 0.0 && lpar[0] == 0.0 && lpar[1] == 0.0 )
+ {
+ AliError("Cannot work with zero tails !");
+ }
+ else
+ {
+ r = FitJpsi2CB2VWG(*hminv,lpar[0],lpar[1],lpar[2],lpar[3]);
+ }
+ }
+ else if ( sFitType == "COUNTJPSI" )
{
r = new AliAnalysisMuMuResult(*hminv,"COUNTJPSI",1);
Double_t n = CountParticle(*hminv,"Jpsi");
r->Set("NofJpsi",n,TMath::Sqrt(n));
}
-// if (r)
-// {
-// r->SetNofInputParticles("Jpsi",GetValue("NofInputJpsi"));
-// }
-
- if (!fSubResults)
- {
- fSubResults = new TObjArray;
- fSubResults->SetOwner(kTRUE);
- }
if ( r )
{
r->SetBin(Bin());
r->SetNofTriggers(NofTriggers());
r->SetNofRuns(NofRuns());
+
+ if (!fSubResults)
+ {
+ fSubResults = new TObjArray;
+ fSubResults->SetOwner(kTRUE);
+ }
+
fSubResults->Add(r);
}
//_____________________________________________________________________________
Bool_t AliAnalysisMuMuResult::HasValue(const char* name, const char* subResultName) const
{
+ /// Whether this result (or subresult if subResultName is provided) has a property
+ /// named "name"
+
if ( strlen(subResultName) > 0 )
{
if ( !fSubResults)
{
AliAnalysisMuMuResult* result = dynamic_cast<AliAnalysisMuMuResult*>(currObj);
+ if (!result)
+ {
+ continue;
+ }
+
if (!result->HasValue(key->String()))
{
AliError(Form("Did not find key %s in of the result to merge",key->String().Data()));
//_____________________________________________________________________________
Int_t AliAnalysisMuMuResult::NofRuns() const
{
+ /// Get the number of runs
if ( !Mother() ) return fNofRuns;
else return Mother()->NofRuns();
}
//_____________________________________________________________________________
Int_t AliAnalysisMuMuResult::NofTriggers() const
{
+ /// Get the number of triggers
+
if ( !Mother() ) return fNofTriggers;
else return Mother()->NofTriggers();
}
//_____________________________________________________________________________
void AliAnalysisMuMuResult::Print(Option_t* opt) const
{
+ /// printout
+
TString sopt(opt);
sopt.ToUpper();
pot.ReplaceAll("ALL","");
pot.ReplaceAll("FULL","");
- std::cout << pot.Data() << Form("%50s- NRUNS %d - NTRIGGER %10d - %s%s",
+ std::cout << pot.Data();
+
+ if ( fAlias.Length() > 0 )
+ {
+ std::cout << Form("%s - ",fAlias.Data());
+ }
+
+ std::cout << Form("%50s - NRUNS %d - NTRIGGER %10d - %s%s",
GetName(),
NofRuns(),
NofTriggers(),
std::cout << " (" << fSubResults->GetEntries()-1 << " subresults)";
}
+ if (!fSubResults )
+ {
+ std::cout << Form(" - REBIN %d",fRebin) << std::endl;
+ }
+
std::cout << std::endl;
if ( sopt.Contains("DUMP") )
std::cout << Form("\t\tMBR %e +- %e",GetValue("MBR"),GetErrorStat("MBR")) << std::endl;
}
- if ( fSubResults && fSubResults->GetEntries() > 1 && ( sopt.Contains("ALL") || sopt.Contains("FULL") ) )
+ if ( fSubResults /* && fSubResults->GetEntries() > 1 */ && ( sopt.Contains("ALL") || sopt.Contains("FULL") ) )
{
std::cout << pot.Data() << "\t===== sub results ===== " << std::endl;
//_____________________________________________________________________________
void AliAnalysisMuMuResult::PrintParticle(const char* particle, const char* opt) const
{
+ /// Print all information about one particule type
+
Double_t npart = GetValue(Form("Nof%s",particle));
if (npart<=0) return;
}
//_____________________________________________________________________________
-void AliAnalysisMuMuResult::PrintValue(const char* key, const char* opt, Double_t value, Double_t errorStat) const
+void AliAnalysisMuMuResult::PrintValue(const char* key, const char* opt, Double_t value, Double_t errorStat)
{
// print one value and its associated error
{
std::cout << opt << Form("\t\t%20s %9.2f +- %5.2f MeV/c^2",key,value*1E3,1E3*errorStat) << std::endl;
}
- else if ( TString(key).BeginsWith("Nof") )
+ else if ( TString(key).Contains("Nof") )
{
std::cout << opt << Form("\t\t%20s %9.2f +- %5.2f",key,value,errorStat) << std::endl;
}
//_____________________________________________________________________________
void AliAnalysisMuMuResult::Set(const char* name, Double_t value, Double_t errorStat)
{
+ /// Set a (value,error) pair with a given name
+
if ( !fMap )
{
fMap = new TMap;
//_____________________________________________________________________________
void AliAnalysisMuMuResult::SetBin(const AliAnalysisMuMuBinning::Range& bin)
{
+ /// Set the bin
+
if (!Mother()) fBin = bin;
else Mother()->SetBin(bin);
}
//_____________________________________________________________________________
void AliAnalysisMuMuResult::SetNofInputParticles(const TH1& hminv)
{
+ /// Set the number of input particle from the invariant mass spectra
+
static const char* particleNames[] = { "Jpsi" , "PsiPrime", "Upsilon","UpsilonPrime" };
for ( Int_t i = 0; i < 4; ++i )
fMinv->SetDirectory(0);
}
+//_____________________________________________________________________________
+AliAnalysisMuMuResult*
+AliAnalysisMuMuResult::SubResult(const char* subResultName) const
+{
+ /// get a given subresult
+ if (!fSubResults)
+ {
+ return 0x0;
+ }
+ TIter next(fSubResults);
+ AliAnalysisMuMuResult* r;
+ while ( ( r = static_cast<AliAnalysisMuMuResult*>(next())) )
+ {
+ if ( r->Alias() == subResultName )
+ {
+ return r;
+ }
+ }
+ return 0x0;
+}
+
const AliAnalysisMuMuBinning::Range& bin);
AliAnalysisMuMuResult(const AliAnalysisMuMuResult& rhs);
- AliAnalysisMuMuResult& operator=(const AliAnalysisMuMuResult&);
+ AliAnalysisMuMuResult& operator=(const AliAnalysisMuMuResult& rhs);
virtual ~AliAnalysisMuMuResult();
void Print(Option_t* opt="") const;
- Bool_t AddFit(const char* fitType, Int_t rebin=1);
+ Bool_t AddFit(const char* fitType, Int_t npar=0, Double_t* par=0x0);
AliAnalysisMuMuResult* CountJpsi(TH1& h);
-// void FitJpsiPsiPrimeCB(TH1& h);
AliAnalysisMuMuResult* FitJpsi(TH1& h);
-// void FitJpsiCBE(TH1& h);
-// void FitJpsiECBE(TH1& h);
-// void FitJpsiPCBE(TH1& h);
-// void FitUpsilon(TH1& h);
- AliAnalysisMuMuResult* FitJpsi2CB2VWG(const TH1& h);
- AliAnalysisMuMuResult* FitJpsiGCBE(TH1& h);
+ AliAnalysisMuMuResult* FitJpsiNA48(const TH1& h);
+ AliAnalysisMuMuResult* FitJpsiCB2VWG(const TH1& h);
+ AliAnalysisMuMuResult* FitJpsi2CB2VWG(const TH1& h, Double_t alphaLow=-1.0, Double_t nLow=-1.0, Double_t alphaUp=-1.0, Double_t nUp=-1.0);
-// SubResult* FitJpsiPsiPrimeCustom(TH1& h);
+ AliAnalysisMuMuResult* FitJpsiGCBE(TH1& h);
Int_t NofRuns() const;
void SetMinv(const TH1& hminv);
+ AliAnalysisMuMuResult* SubResult(const char* subResultName) const;
+
TObjArray* SubResults() const { return fSubResults; }
- static Double_t CountParticle(const TH1& hminv, const char* particle, Double_t sigma=-1);
-
- static Double_t ErrorAB(Double_t a, Double_t aerr, Double_t b, Double_t berr);
-
- static Double_t ErrorABC(Double_t a, Double_t aerr, Double_t b, Double_t berr, Double_t c, Double_t cerror);
-
Long64_t Merge(TCollection* list);
AliAnalysisMuMuResult* Mother() const { return fMother; }
Double_t Weight() const { return fWeight > 0 ? fWeight : fNofTriggers; }
void SetWeight(Double_t w) { fWeight=w; }
+
+ static Double_t CountParticle(const TH1& hminv, const char* particle, Double_t sigma=-1.0);
+
+ static Double_t ErrorAB(Double_t a, Double_t aerr, Double_t b, Double_t berr);
+
+ static Double_t ErrorABC(Double_t a, Double_t aerr, Double_t b, Double_t berr, Double_t c, Double_t cerror);
+
+ static void PrintValue(const char* key, const char* opt, Double_t value, Double_t errorStat);
+
+ void SetAlias(const char* alias) { fAlias = alias; }
+
+ TString Alias() const { if ( fAlias.Length()>0) return fAlias; else return GetName(); }
private:
kErrorStat=1
};
-
void PrintParticle(const char* particle, const char* opt) const;
- void PrintValue(const char* key, const char* opt, Double_t value, Double_t errorStat) const;
private:
Int_t fNofRuns; // number of runs used to get this result
AliAnalysisMuMuResult* fMother; // mother result
mutable THashList* fKeys; //! keys we have in our internal map (or the one of our subresults)
Double_t fWeight; // weight of this result (default 1.0)
+ Int_t fRebin; // rebin level of minv spectra
+
+ TString fTriggerClass; // trigger class for this result
+ TString fEventSelection; // event selection for this result
+ TString fPairSelection; // pair selection for this result
+ TString fCentralitySelection; // centrality selection for this result
+
+ TString fAlias; // alias name
- ClassDef(AliAnalysisMuMuResult,5) // a class to hold invariant mass analysis results (counts, yields, AccxEff, R_AB, etc...)
+ ClassDef(AliAnalysisMuMuResult,8) // a class to hold invariant mass analysis results (counts, yields, AccxEff, R_AB, etc...)
};
#endif
+/**************************************************************************
+ * Copyright(c) 1998-1999, 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. *
+ **************************************************************************/
+
+
+// AliAnalysisMuMuSpectra : a class to encapsulate results from MuMu analysis
+//
+// Spectra can be merged and converted into histograms
+//
+// author: L. Aphecetche (Subatech)
+//
+
#include "AliAnalysisMuMuSpectra.h"
#include "AliLog.h"
fBinning(0x0),
fBins(0x0)
{
+ // default ctor
+}
+
+//______________________________________________________________________________
+AliAnalysisMuMuSpectra::AliAnalysisMuMuSpectra(const AliAnalysisMuMuSpectra& rhs)
+: TNamed(rhs.GetName(),rhs.GetTitle()),
+fBinning(0x0),
+fBins(0x0)
+{
+ // copy ctor
+
+ if ( rhs.fBinning )
+ {
+ fBinning = new AliAnalysisMuMuBinning(*rhs.fBinning);
+ }
+
+ TIter next(rhs.fBins);
+ AliAnalysisMuMuBinning::Range* bin;
+
+ while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(next()) ) )
+ {
+ if (!fBins)
+ {
+ fBins = new TObjArray;
+ fBins->SetOwner(kTRUE);
+ }
+ fBins->Add(bin);
+ }
+}
+
+//______________________________________________________________________________
+AliAnalysisMuMuSpectra&
+AliAnalysisMuMuSpectra::operator=(const AliAnalysisMuMuSpectra& rhs)
+{
+ // assignment operator
+
+ if (this==&rhs) return *this;
+
+ delete fBinning;
+ fBinning = 0x0;
+ delete fBins;
+ fBins = 0x0;
+
+ if ( rhs.fBinning )
+ {
+ fBinning = new AliAnalysisMuMuBinning(*rhs.fBinning);
+ }
+
+ TIter next(rhs.fBins);
+ AliAnalysisMuMuBinning::Range* bin;
+ while ( ( bin = static_cast<AliAnalysisMuMuBinning::Range*>(next()) ) )
+ {
+ if (!fBins)
+ {
+ fBins = new TObjArray;
+ fBins->SetOwner(kTRUE);
+ }
+ fBins->Add(bin);
+ }
+
+ return *this;
}
//______________________________________________________________________________
AliAnalysisMuMuSpectra::~AliAnalysisMuMuSpectra()
{
+ // dtor
delete fBinning;
delete fBins;
}
void AliAnalysisMuMuSpectra::AdoptResult(const AliAnalysisMuMuBinning::Range& bin,
AliAnalysisMuMuResult* result)
{
+ // adopt (i.e. we are becoming the owner) a result for a given bin
if (!fBinning)
{
fBinning = new AliAnalysisMuMuBinning;
{
AliAnalysisMuMuResult* thisResult = static_cast<AliAnalysisMuMuResult*>(fBins->At(i));
AliAnalysisMuMuResult* accResult = static_cast<AliAnalysisMuMuResult*>(bins->At(i));
- AliInfoClass(Form("i=%d",i ));
- StdoutToAliInfoClass(thisResult->Print("full");
- std::cout << "----" << std::endl;
- accResult->Print("full"));
+// AliInfoClass(Form("i=%d",i ));
+// StdoutToAliInfoClass(thisResult->Print("full");
+// std::cout << "----" << std::endl;
+// accResult->Print("full"));
thisResult->Correct(*accResult,particle,subResultName);
}
return kTRUE;
}
+//______________________________________________________________________________
+AliAnalysisMuMuResult*
+AliAnalysisMuMuSpectra::GetResultForBin(const AliAnalysisMuMuBinning::Range& bin) const
+{
+ /// Get result for a given bin
+ /// Warning: this requires a loop on bins
+
+ if ( IsEmpty() ) return 0x0;
+
+ TObjArray* bins = fBinning->CreateBinObjArray();
+
+ Int_t j(-1);
+
+ StdoutToAliDebug(1,std::cout << "searching for "; bin.Print());
+
+ for ( Int_t i = 0; i <= bins->GetLast() && j < 0 ; ++i )
+ {
+ AliAnalysisMuMuBinning::Range* b = static_cast<AliAnalysisMuMuBinning::Range*>(bins->At(i));
+
+ StdoutToAliDebug(1,b->Print(););
+
+ if ( bin == *b )
+ {
+ j = i;
+ }
+ }
+
+ delete bins;
+
+ if (j>=0)
+ {
+ return static_cast<AliAnalysisMuMuResult*>(fBins->At(j));
+ }
+ else
+ {
+ StdoutToAliDebug(1,std::cout << "Could not find result for bin:" << std::endl; bin.Print(););
+ }
+ return 0x0;
+}
+
+//______________________________________________________________________________
+Bool_t AliAnalysisMuMuSpectra::HasValue(const char* what) const
+{
+ // whether or not our result(s) has a given property
+ if ( IsEmpty() ) return kFALSE;
+
+ AliAnalysisMuMuResult* r = static_cast<AliAnalysisMuMuResult*>(fBins->First());
+
+ return r->HasValue(what);
+}
+
//______________________________________________________________________________
Bool_t AliAnalysisMuMuSpectra::IsEmpty() const
{
+ // whether this spectra is empty or not
return ( fBins==0x0 || fBins->GetEntries()<=0 );
}
}
//_____________________________________________________________________________
-TH1* AliAnalysisMuMuSpectra::Plot(const char* what, const char* subresult) const
+TH1* AliAnalysisMuMuSpectra::Plot(const char* what, const char* subresult, Bool_t divideByBinWidth) const
{
+ // Convert the spectra into an histogram
+
TString swhat(what);
swhat.ToUpper();
TH1* h(0x0);
- AliDebugClass(1,Form("nbins=%d nresults=%d",binArray->GetEntries(),fBins->GetEntries()));
+ AliDebug(1,Form("nbins=%d nresults=%d",binArray->GetEntries(),fBins->GetEntries()));
for ( Int_t j = 0; j < TMath::Min(binArray->GetEntries(),fBins->GetEntries()); ++j )
{
AliAnalysisMuMuResult* r = static_cast<AliAnalysisMuMuResult*>(fBins->At(j));
- if ( strlen(subresult) > 0 )
+ if ( strlen(subresult) > 0 && r->SubResults() )
{
- TObjArray* sr = r->SubResults();
- if (!sr) continue;
TString sub(subresult);
sub.ToUpper();
- r = static_cast<AliAnalysisMuMuResult*>(sr->FindObject(sub.Data()));
+ r = r->SubResult(sub.Data());
if (!r) continue;
}
const AliAnalysisMuMuBinning::Range& b = r->Bin();
- r->Print();
- b.Print();
-
if (!h)
{
h = new TH1F(r->GetName(),r->GetName(),binArray->GetEntries(),bins);
Double_t y = r->GetValue(what);
Double_t yerr = r->GetErrorStat(what);
-
-// if ( !swhat.BeginsWith("ACC") )
- if ( swhat.Contains("NOF") )
+
+ if ( divideByBinWidth && b.WidthX()>0 )
{
y /= (b.WidthX());
yerr /= (b.WidthX());
}
+ std::cout << b.AsString();
+ AliAnalysisMuMuResult::PrintValue(swhat.Data(),"",y,yerr);
+
h->SetBinContent(j+1,y);
h->SetBinError(j+1,yerr);
- AliInfoClass(Form("%e +- %e",y,yerr));
}
delete binArray;
//______________________________________________________________________________
void AliAnalysisMuMuSpectra::Print(Option_t* opt) const
{
+ // printout
+
if (!IsEmpty())
{
TString sopt(opt);
- sopt.ToUpper();
- if ( sopt.Contains("BINNING") )
- {
- fBinning->Print(opt);
- }
Int_t nmax = sopt.Atoi();
if ( nmax <= 0 ) nmax = fBins->GetEntries();
for ( Int_t i = 0; i < nmax; ++i )
{
- AliAnalysisMuMuBinning::Range* r = static_cast<AliAnalysisMuMuBinning::Range*>(fBins->At(i));
+ AliAnalysisMuMuResult* r = static_cast<AliAnalysisMuMuResult*>(fBins->At(i));
if (r) r->Print(opt);
}
-
}
}
{
public:
AliAnalysisMuMuSpectra(const char* name="", const char* title="");
+ AliAnalysisMuMuSpectra(const AliAnalysisMuMuSpectra& rhs);
+ AliAnalysisMuMuSpectra& operator=(const AliAnalysisMuMuSpectra& rhs);
+
virtual ~AliAnalysisMuMuSpectra();
void AdoptResult(const AliAnalysisMuMuBinning::Range& bin, AliAnalysisMuMuResult* result);
Long64_t Merge(TCollection* list);
- TH1* Plot(const char* what="NofJpsi", const char* subresult="") const;
+ TH1* Plot(const char* what="NofJpsi", const char* subresult="", Bool_t divideByBinWidth=kTRUE) const;
void Print(Option_t* opt="") const;
Bool_t Correct(const AliAnalysisMuMuSpectra& accEff, const char* particle, const char* subResultName="");
+ AliAnalysisMuMuResult* GetResultForBin(const AliAnalysisMuMuBinning::Range& bin) const;
+
+ Bool_t HasValue(const char* what="NofJpsi") const;
+
private:
AliAnalysisMuMuBinning* fBinning; // internal binning
TObjArray* fBins; // the results (bin by bin)
AliTriggerConfiguration* tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
AliTriggerRunScalers* trs = static_cast<AliTriggerRunScalers*>(GetOCDBObject("GRP/CTP/Scalers",runNumber));
AliGRPObject* grp = static_cast<AliGRPObject*>(GetOCDBObject("GRP/GRP/Data",runNumber));
+
+ if (!tc || !trs || !grp) return 0x0;
TString diCurrent(Form("L3:%5.0f;DIP:%5.0f [L3:%d;DIP:%d]",
grp->GetL3Current((AliGRPObject::Stats)0),
grp->GetL3Polarity(),
grp->GetDipolePolarity()));
- if (!tc || !trs || !grp) return 0x0;
-
const TObjArray& trClasses = tc->GetClasses();
const TObjArray* scalers = trs->GetScalersRecords();
if (!lumiB)
{
AliError(Form("Did not find lumiTrigger %s for run %09d",lumiTriggerClassName.Data(),runNumber));
+ continue;
}
Float_t pacCorrection(1.0);
/// - mu ( = -TMath::Log( 1 - P(0) ) where P(0) is the proba to have zero collisions in a bunch crossing)
/// - pileupfactor = mu/(1-exp(-mu)) : the factor to apply to correct the cint1b count rate
/// - vsnb = L0B/(NumberOfInteractingBunches*11245)
+ /// - NBCX = NumberOfInteractingBunches
TString swhat(what);
swhat.ToUpper();
UInt_t l2a = scaler->GetL2CA() - refl2a;
UInt_t timelapse = seconds - reft;
- if ( l0b <= 2 || ( l0a <= 2 && l0a != 0 ) || timelapse <= 9 ) continue;
+ if ( removeZeros && ( l0b <= 2 || ( l0a <= 2 && l0a != 0 ) || timelapse <= 9 ) )
+ {
+ AliInfo(Form("Skipping point for %s l0b %d l0a %d timelapse %d ts %s",
+ triggerClassName,l0b,l0a,timelapse,ts.AsString()));
+ continue;
+ }
reft = seconds;
refl0b = scaler->GetLOCB();
value = l0b/(11245.0*numberOfInteractingBunches);
error = -1.0; // FIXME
}
+ else if ( swhat.Contains("NBCX"))
+ {
+ value = numberOfInteractingBunches;
+ error = 1.0; // FIXME
+ }
else
{
value = timelapse;
if ( ! swhat.Contains("OVER") && ! swhat.Contains("RATIO") &&
! swhat.Contains("MU") && ! swhat.Contains("PILEUPFACTOR") &&
- ! swhat.Contains("RAW") )
+ ! swhat.Contains("RAW") & ! swhat.Contains("NBCX") )
{
value /= timelapse;
}
{
/// Plots the evolution of one trigger ratio
- TGraph* g1 = PlotTriggerEvolution(triggerClassName1,what1,kFALSE);
- TGraph* g2 = PlotTriggerEvolution(triggerClassName2,what2,kFALSE);
+ Bool_t draw(kFALSE);
+ Bool_t removeZeros(kFALSE);
+
+ TGraph* g1 = PlotTriggerEvolution(triggerClassName1,what1,draw,0x0,removeZeros);
+ TGraph* g2 = PlotTriggerEvolution(triggerClassName2,what2,draw,0x0,removeZeros);
if (!g1 || !g2) return 0x0;
+/**************************************************************************
+ * Copyright(c) 1998-1999, 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. *
+ **************************************************************************/
+
+//
+// AliMuonAccEffSubmitter : a class to help submit Acc x Eff simulations
+// anchored to real runs for J/psi, upsilon, single muons, etc...
+//
+// This class is dealing with 3 different directories :
+//
+// - template directory ($ALICE_ROOT/PWG/muondep/AccEffTemplates) containing the
+// basic template files to be used for a simuation. A template can contain
+// some variables that will be replaced during during the copy from template
+// to local dir
+//
+// - local directory, where the files from the template directory, are copied
+// once the class has been configured properly (i.e. using the various Set, Use,
+// etc... methods). Some other files (e.g. JDL ones) are generated from
+// scratch and also copied into this directory.
+// At this point one could(should) check the files, as they are the ones
+// to be copied to the remote directory for the production
+//
+// - remote directory, the alien directory where the files will be copied
+// (from the local directory) before the actual submission
+//
+// ==========================================================
+//
+// Basic usage
+//
+// AliMuonAccEffSubmitter a;
+// a.UseOCDBSnapshots(kFALSE);
+// a.SetRemoteDir("/alice/cern.ch/user/l/laphecet/Analysis/LHC13d/simjpsi/pp503z0");
+// a.ShouldOverwriteFiles(true);
+// a.MakeNofEventsPropToTriggerCount("CMUL7-B-NOPF-MUON");
+// a.SetVar("VAR_GENLIB_PARNAME","\"pp 5.03\"");
+// a.SetRunList(195682);
+// a.Print();
+// a.Run("test"); // will do everything but the submit
+// a.Submit(false); // actual submission
+//
+//
+// author: Laurent Aphecetche (Subatech
+//
+
#include "AliMuonAccEffSubmitter.h"
#include "AliAnalysisTriggerScalers.h"
#include "TMap.h"
#include "TMath.h"
#include "TObjString.h"
+#include "TROOT.h"
#include "TString.h"
#include "TSystem.h"
#include <vector>
}
//______________________________________________________________________________
-AliMuonAccEffSubmitter::AliMuonAccEffSubmitter()
+AliMuonAccEffSubmitter::AliMuonAccEffSubmitter(const char* generator)
: TObject(),
fScalers(0x0),
fRemoteDir(""),
fReferenceTrigger(""),
fRatio(1.0),
+fFixedNofEvents(10000),
fMaxEventsPerChunk(5000),
fLocalDir(gSystem->pwd()),
fOCDBPath("raw://"),
fIsValid = kFALSE;
}
- SetPackages("VO_ALICE@AliRoot::v5-03-Rev-09","VO_ALICE@GEANT3::v1-14-6","VO_ALICE@ROOT::v5-34-02-1");
+ SetPackages("VO_ALICE@AliRoot::v5-03-Rev-18","VO_ALICE@GEANT3::v1-14-8","VO_ALICE@ROOT::v5-34-05-1");
-// SetVar("VAR_ENERGY","5.03");
- SetVar("VAR_GENLIB_TYPE","AliGenMUONlib::kJpsi");
- SetVar("VAR_GENLIB_PARNAME","\"pPb 5.03\"");
SetVar("VAR_OCDB_PATH","\"raw://\"");
+
+ SetVar("VAR_GENPARAM_GENLIB_TYPE","AliGenMUONlib::kJpsi");
+ SetVar("VAR_GENPARAM_GENLIB_PARNAME","\"pPb 5.03\"");
+
+ SetVar("VAR_GENCORRHF_QUARK","5");
+ SetVar("VAR_GENCORRHF_ENERGY","5");
+
+ SetVar("VAR_GENPARAMCUSTOM_PDGPARTICLECODE","443");
+
+ // default values below are from J/psi p+Pb (from muon_calo pass)
+ SetVar("VAR_GENPARAMCUSTOM_Y_P0","4.08E5");
+ SetVar("VAR_GENPARAMCUSTOM_Y_P1","7.1E4");
+
+ SetVar("VAR_GENPARAMCUSTOM_PT_P0","1.13E9");
+ SetVar("VAR_GENPARAMCUSTOM_PT_P1","18.05");
+ SetVar("VAR_GENPARAMCUSTOM_PT_P2","2.05");
+ SetVar("VAR_GENPARAMCUSTOM_PT_P3","3.34");
+
UseOCDBSnapshots(kTRUE);
+
+ SetGenerator(generator);
}
//______________________________________________________________________________
delete fVars;
}
+//______________________________________________________________________________
+Bool_t AliMuonAccEffSubmitter::CheckCompilation(const char* file) const
+{
+ /// Check whether file can be compiled or not
+ /// FIXME: use gSystem->TempFileName for tmpfile !
+
+ Bool_t rv(kTRUE);
+
+ TString sfile(gSystem->BaseName(file));
+ TString tmpfile(Form("tmpfile_%s",sfile.Data()));
+
+ gSystem->Exec(Form("cp %s %s",file,tmpfile.Data()));
+
+ ReplaceVars(tmpfile.Data());
+
+ gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/EVGEN");
+
+ if (gROOT->LoadMacro(Form("%s++",tmpfile.Data())))
+ {
+ AliError(Form("macro %s can not be compiled. Please check.",file));
+ rv = kFALSE;
+ }
+
+ gSystem->Exec(Form("rm %s",tmpfile.Data()));
+
+ return rv;
+}
+
+
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::CheckLocal() const
{
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::CopyFile(const char* localFile)
{
+ /// copy a local file to remote destination
TString local;
if ( gSystem->IsAbsoluteFileName(localFile) )
if ( ok )
{
+ AliDebug(1,Form("cp %s alien://%s",local.Data(),remote.Data()));
return TFile::Cp(local.Data(),Form("alien://%s",remote.Data()));
}
else
if (!IsValid()) return kFALSE;
+ AliDebug(1,"");
+
if ( CheckRemoteDir() )
{
TString sdir(gSystem->ExpandPathName(LocalDir()));
if (!IsValid()) return kFALSE;
+ AliDebug(1,"");
+
TIter next(TemplateFileList());
TObjString* file;
//______________________________________________________________________________
std::ostream* AliMuonAccEffSubmitter::CreateJDLFile(const char* name) const
{
+ /// Create a JDL file
+ AliDebug(1,"");
+
TString jdl(Form("%s/%s",fLocalDir.Data(),name));
if ( !ShouldOverwriteFiles() && !gSystem->AccessPathName(jdl.Data()) )
///______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::GenerateMergeJDL(const char* name)
{
+ /// Create the JDL for merging jobs
+ /// FIXME: not checked !
+
+ AliDebug(1,"");
+
std::ostream* os = CreateJDLFile(name);
if (!os)
/// Generate (locally) the JDL to perform the simulation+reco+aod filtering
/// (to be then copied to the grid and finally submitted)
+ AliDebug(1,"");
+
std::ostream* os = CreateJDLFile(name);
if (!os)
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::GetLastStage(const char* remoteDir) const
{
+ /// Get the last staging phase already performed
+ /// FIXME : not checked !
+
Int_t n = 0, lastStage = 0;
gSystem->Exec(Form("alien_ls -F %s | grep Stage_.*/ > __stage__", remoteDir));
ifstream f("__stage__");
return lastStage;
}
+//______________________________________________________________________________
+TObjArray* AliMuonAccEffSubmitter::GetVariables(const char* file) const
+{
+ /// Find the variables in the file
+
+ std::ifstream in(file);
+ char line[1024];
+ TObjArray* variables(0x0);
+
+ while ( in.getline(line,1023,'\n') )
+ {
+ TString sline(line);
+ while (sline.Contains("VAR_") && !sline.BeginsWith("//") )
+ {
+ Int_t i1 = sline.Index("VAR_");
+ Int_t i2(i1);
+
+ while ( ( i2 < sline.Length() ) && ( isalnum(sline[i2]) || sline[i2]=='_' ) ) ++i2;
+
+ if (!variables)
+ {
+ variables = new TObjArray;
+ variables->SetOwner(kTRUE);
+ }
+
+ TString var = sline(i1,i2-i1);
+ if ( !variables->FindObject(var) )
+ {
+ variables->Add(new TObjString(var));
+ }
+ sline.ReplaceAll(var,"");
+ }
+ }
+
+ in.close();
+
+ return variables;
+}
+
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::HasVars(const char* file) const
{
if (!fScalers) return kFALSE;
+ AliDebug(1,"");
+
const std::vector<int>& runs = fScalers->GetRunList();
Bool_t ok(kTRUE);
TString ocdbSim(Form("%s/OCDB/%d/OCDB_sim.root",SnapshotDir().Data(),runNumber));
TString ocdbRec(Form("%s/OCDB/%d/OCDB_rec.root",SnapshotDir().Data(),runNumber));
- if ( !gSystem->AccessPathName(ocdbSim.Data()) &&
+ if ( !gSystem->AccessPathName(ocdbSim.Data()) &&
!gSystem->AccessPathName(ocdbRec.Data()) )
{
AliWarning(Form("Local OCDB snapshots already there for run %d. Will not redo them. If you want to force them, delete them by hand !",runNumber));
+ continue;
}
else
{
void AliMuonAccEffSubmitter::Output(std::ostream& out, const char* key,
const TObjArray& values) const
{
+ /// output to ostream of key,{values} pair
+
out << key << " = ";
Int_t n = values.GetEntries();
const char* v5, const char* v6, const char* v7,
const char* v8, const char* v9) const
{
+ /// output to ostream
+
TObjArray values;
values.SetOwner(kTRUE);
//______________________________________________________________________________
void AliMuonAccEffSubmitter::Print(Option_t* /*opt*/) const
{
+ /// Printout
+
if (!IsValid())
{
std::cout << std::string(80,'*') << std::endl;
std::cout << std::string(80,'*') << std::endl;
}
- std::cout << "Template directory = " << fTemplateDir.Data() << std::endl;
- std::cout << "Local directory = " << fLocalDir.Data() << std::endl;
- std::cout << "Remote directory = " << fRemoteDir.Data() << std::endl;
+ std::cout << "Template directory = " << fTemplateDir.Data() << std::endl;
+ std::cout << "Local directory = " << fLocalDir.Data() << std::endl;
+ std::cout << "Remote directory = " << fRemoteDir.Data() << std::endl;
+
+ if ( fSnapshotDir != fLocalDir )
+ {
+ std::cout << "Snapshots directory = " << fSnapshotDir.Data() << std::endl;
+ }
std::cout << "OCDB path = " << fOCDBPath.Data() << std::endl;
}
//______________________________________________________________________________
-Bool_t AliMuonAccEffSubmitter::ReplaceVars(const char* file)
+Bool_t AliMuonAccEffSubmitter::ReplaceVars(const char* file) const
{
+ /// Replace the variables (i.e. things starting by VAR_) found in file
+
std::ifstream in(file);
char line[1024];
TObjArray lines;
const char* geant3,
const char* api)
{
+ /// Set the packages to be used by the jobs
+ /// Must be a valid combination, see http://alimonitor.cern.ch/packages/
+ ///
fPackageAliroot = aliroot;
fPackageRoot = root;
fPackageGeant3 = geant3;
return dir;
}
+//______________________________________________________________________________
+Bool_t AliMuonAccEffSubmitter::SetGenerator(const char* generator)
+{
+ // set the variable to select the generator macro in Config.C
+
+ gSystem->Load("libEVGEN");
+
+ fIsValid = kFALSE;
+
+ TString generatorFile(Form("%s/%s.C",fTemplateDir.Data(),generator));
+
+ Int_t nofMissingVariables(0);
+
+ // first check we indeed have such a macro
+ if (!gSystem->AccessPathName(generatorFile.Data()))
+ {
+ TObjArray* variables = GetVariables(generatorFile.Data());
+
+ TIter next(variables);
+ TObjString* var;
+
+ while ( ( var = static_cast<TObjString*>(next())) )
+ {
+ if ( !fVars->GetValue(var->String()) )
+ {
+ ++nofMissingVariables;
+ AliError(Form("file %s expect the variable %s to be defined, but we've not defined it !",generatorFile.Data(),var->String().Data()));
+ }
+ }
+
+ delete variables;
+
+ if ( !nofMissingVariables )
+ {
+ if (CheckCompilation(generatorFile.Data()))
+ {
+ fIsValid = kTRUE;
+ SetVar("VAR_GENERATOR",Form("%s",generator));
+ TemplateFileList()->Add(new TObjString(Form("%s.C",generator)));
+ return kTRUE;
+ }
+ }
+ else
+ {
+ return kFALSE;
+ }
+ }
+ else
+ {
+ AliError(Form("Can not work with the macro %s",generatorFile.Data()));
+ }
+ return kFALSE;
+}
+
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::SetMergedDir(const char* dir, Bool_t create)
{
+ // Set the merged directory to be used
fMergedDir = GetRemoteDir(dir,create);
return (fMergedDir.Length()>0);
}
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::SetRemoteDir(const char* dir, Bool_t create)
{
+ // Set the remote directory to be used
fRemoteDir = GetRemoteDir(dir,create);
return (fIsValid = (fRemoteDir.Length()>0));
}
UpdateLocalFileList(kTRUE);
}
+//______________________________________________________________________________
+void AliMuonAccEffSubmitter::SetOCDBPath(const char* ocdbPath)
+{
+ /// Sets the OCDB path to be used
+
+ fOCDBPath = ocdbPath;
+
+ if (fScalers)
+ {
+ // redefine trigger scalers to use the new ocdb path
+ AliAnalysisTriggerScalers* ts = new AliAnalysisTriggerScalers(fScalers->GetRunList(),
+ fOCDBPath.Data());
+
+ delete fScalers;
+ fScalers = ts;
+ }
+}
+
+
+//______________________________________________________________________________
+void AliMuonAccEffSubmitter::SetOCDBSnapshotDir(const char* dir)
+{
+ // change the directory used for snapshot
+
+ if (gSystem->AccessPathName(Form("%s/OCDB",dir)))
+ {
+ AliError(Form("Snapshot top directory (%s) should contain an OCDB subdir with runnumbers in there",dir));
+ }
+ else
+ {
+ fSnapshotDir = dir;
+ }
+}
+
//______________________________________________________________________________
Bool_t AliMuonAccEffSubmitter::SetVar(const char* varname, const char* value)
{
+ /// Set a variable
+
TString s(varname);
s.ToUpper();
if (!s.BeginsWith("VAR_"))
if (!IsValid()) return 0;
+ AliDebug(1,"");
+
gGrid->Cd(RemoteDir());
if (!RemoteFileExists(RunJDLName()))
//______________________________________________________________________________
void AliMuonAccEffSubmitter::UpdateLocalFileList(Bool_t clearSnapshots)
{
+ /// Update the list of local files
+
if (!fScalers) return;
if ( clearSnapshots )
//______________________________________________________________________________
void AliMuonAccEffSubmitter::UseOCDBSnapshots(Bool_t flag)
{
+ /// Whether or not to use OCDB snapshots
+ /// Using OCDB snapshots will speed-up both the sim and reco initialization
+ /// phases on each worker node, but takes time to produce...
+ /// So using them is not always a win-win...
+
fUseOCDBSnapshots = flag;
if ( flag )
{
#ifndef ALIMUONACCEFFSUBMITTER_H
#define ALIMUONACCEFFSUBMITTER_H
+//
+// AliMuonAccEffSubmitter : a class to help submit Acc x Eff simulations
+// anchored to real runs for J/psi, upsilon, single muons, etc...
+//
+// author: Laurent Aphecetche (Subatech)
+//
+
#include "TObject.h"
#include "TString.h"
#include "Riostream.h"
class AliMuonAccEffSubmitter : public TObject
{
public:
- AliMuonAccEffSubmitter();
+ AliMuonAccEffSubmitter(const char* generator="GenParamCustom");
virtual ~AliMuonAccEffSubmitter();
Bool_t SetRemoteDir(const char* dir, Bool_t create = kTRUE);
- void SetLocalDir(const char* LocalDir) { fLocalDir = LocalDir; }
+ void SetLocalDir(const char* localdir) { fLocalDir = localdir; }
Bool_t SetMergedDir(const char* dir, Bool_t create = kTRUE);
void UseOCDBSnapshots(Bool_t flag);
UInt_t NofRuns() const;
void SetRatio(Float_t ratio) { fRatio = ratio; }
- void SetOCDBPath(const char* ocdbPath) { fOCDBPath = ocdbPath; }
-
+ void SetOCDBPath(const char* ocdbPath);
+
void MakeNofEventsPropToTriggerCount(const char* trigger="CMUL8-S-NOPF-MUON", Float_t ratio=1.0) { fReferenceTrigger = trigger; fRatio = ratio; }
void MakeNofEventsFixed(Int_t nevents) { fFixedNofEvents = nevents; fReferenceTrigger=""; fRatio=0.0; }
void SetOCDBSnapshotDir(const char* dir);
+ Bool_t SetGenerator(const char* generator);
+
+ TObjArray* GetVariables(const char* file) const;
+
+ Bool_t IsValid() const { return fIsValid; }
+
private:
+ TString SnapshotDir() const { return fSnapshotDir; }
+
TString GetRemoteDir(const char* dir, Bool_t create);
std::ostream* CreateJDLFile(const char* name) const;
void Output(std::ostream& out, const char* key, const TObjArray& values) const;
- Bool_t ReplaceVars(const char* file);
+ Bool_t ReplaceVars(const char* file) const;
TObjArray* TemplateFileList() const;
TObjArray* LocalFileList() const;
- Bool_t IsValid() const { return fIsValid; }
-
Bool_t HasVars(const char* localFile) const;
void UpdateLocalFileList(Bool_t clearSnapshot=kFALSE);
- TString SnapshotDir() const { return fSnapshotDir; }
+ Bool_t CheckCompilation(const char* file) const;
+
+private:
+ AliMuonAccEffSubmitter(const AliMuonAccEffSubmitter& rhs);
+ AliMuonAccEffSubmitter& operator=(const AliMuonAccEffSubmitter& rhs);
private:
- AliAnalysisTriggerScalers* fScalers;
- TString fRemoteDir;
- TString fReferenceTrigger;
- Float_t fRatio;
- Int_t fFixedNofEvents;
- Int_t fMaxEventsPerChunk;
- TString fLocalDir;
- TString fOCDBPath;
- TString fTemplateDir;
- TString fPackageAliroot;
- TString fPackageGeant3;
- TString fPackageRoot;
- TString fPackageApi;
- TString fMergedDir;
- Int_t fSplitMaxInputFileNumber;
- Int_t fCompactMode;
- Bool_t fShouldOverwriteFiles;
- TMap* fVars;
- TString fExternalConfig;
- Bool_t fUseOCDBSnapshots;
- Bool_t fIsValid;
- mutable TObjArray* fTemplateFileList;
- mutable TObjArray* fLocalFileList;
- TString fSnapshotDir;
+ AliAnalysisTriggerScalers* fScalers; // helper class used to handle the runlist and the scalers
+ TString fRemoteDir; // remote directory to used in alien
+ TString fReferenceTrigger; // reference trigger (if any) to be used to get the number of events to be used per run
+ Float_t fRatio; // ratio simulated events vs real events
+ Int_t fFixedNofEvents; // fixed number of events to be used per run
+ Int_t fMaxEventsPerChunk; // max events to generate per subjob
+ TString fLocalDir; // local directory
+ TString fOCDBPath; // OCDB path
+ TString fTemplateDir; // template directory
+ TString fPackageAliroot; // which aliroot package to use
+ TString fPackageGeant3; // which geant3 package to use
+ TString fPackageRoot; // which root package to use (for valid root,geant3,aliroot combinations see http://alimonitor.cern.ch/packages/)
+ TString fPackageApi; // which API package to use
+ TString fMergedDir; // merge directory
+ Int_t fSplitMaxInputFileNumber; // used for merging jdl
+ Int_t fCompactMode; // controls which outputs are kept (0=everything, 1=only aods)
+ Bool_t fShouldOverwriteFiles; // whether any copy (of template to local) is allowed to overwrite existing files
+ TMap* fVars; // map of the variables we can replace in template files
+ TString fExternalConfig; // path to an (optional) external config file
+ Bool_t fUseOCDBSnapshots; // whether to use OCDB snapshots or not
+ Bool_t fIsValid; // whether this object is valid (i.e. properly configured)
+ mutable TObjArray* fTemplateFileList; // list of template files
+ mutable TObjArray* fLocalFileList; // list of local files
+ TString fSnapshotDir; // directory for OCDB snapshots
ClassDef(AliMuonAccEffSubmitter,1) // Helper class to submit AccxEff single particle simulations
};
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff