#include <TCanvas.h>
#include <TParticle.h>
#include <TDatabasePDG.h>
+#include <TProfile.h>
#include <TH1I.h>
#include <TStyle.h>
#include <TFile.h>
fCFManager(0x0),
fHistEventsProcessed(0x0),
fCorrelation(0x0),
+ fListProfiles(0),
fCountMC(0),
fCountAcc(0),
fCountVertex(0),
fUseFlatPtWeight(kFALSE),
fUseZWeight(kFALSE),
fUseNchWeight(kFALSE),
+ fUseTrackletsWeight(kFALSE),
fNvar(0),
fPartName(""),
fDauNames(""),
fHistoMCNch(0x0),
fResonantDecay(0),
fLctoV0bachelorOption(1),
- fGenLctoV0bachelorOption(0)
+ fGenLctoV0bachelorOption(0),
+ fUseSelectionBit(kTRUE),
+ fPDGcode(0),
+ fMultiplicityEstimator(kNtrk10),
+ fRefMult(9.26),
+ fZvtxCorrectedNtrkEstimator(kFALSE),
+ fIsPPData(kFALSE),
+ fIsPPbData(kFALSE)
{
//
//Default ctor
//
+ for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
}
//___________________________________________________________________________
AliCFTaskVertexingHF::AliCFTaskVertexingHF(const Char_t* name, AliRDHFCuts* cuts, TF1* func) :
fCFManager(0x0),
fHistEventsProcessed(0x0),
fCorrelation(0x0),
+ fListProfiles(0),
fCountMC(0),
fCountAcc(0),
fCountVertex(0),
fUseFlatPtWeight(kFALSE),
fUseZWeight(kFALSE),
fUseNchWeight(kFALSE),
+ fUseTrackletsWeight(kFALSE),
fNvar(0),
fPartName(""),
fDauNames(""),
fHistoMCNch(0x0),
fResonantDecay(0),
fLctoV0bachelorOption(1),
- fGenLctoV0bachelorOption(0)
+ fGenLctoV0bachelorOption(0),
+ fUseSelectionBit(kTRUE),
+ fPDGcode(0),
+ fMultiplicityEstimator(kNtrk10),
+ fRefMult(9.26),
+ fZvtxCorrectedNtrkEstimator(kFALSE),
+ fIsPPData(kFALSE),
+ fIsPPbData(kFALSE)
{
//
// Constructor. Initialization of Inputs and Outputs
DefineOutput(2,AliCFContainer::Class());
DefineOutput(3,THnSparseD::Class());
DefineOutput(4,AliRDHFCuts::Class());
+ for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
+ DefineOutput(5,TList::Class()); // slot #5 keeps the zvtx Ntrakclets correction profiles
fCuts->PrintAll();
}
fFuncWeight = c.fFuncWeight;
fHistoMeasNch = c.fHistoMeasNch;
fHistoMCNch = c.fHistoMCNch;
+ for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
}
return *this;
}
fCFManager(c.fCFManager),
fHistEventsProcessed(c.fHistEventsProcessed),
fCorrelation(c.fCorrelation),
+ fListProfiles(c.fListProfiles),
fCountMC(c.fCountMC),
fCountAcc(c.fCountAcc),
fCountVertex(c.fCountVertex),
fUseFlatPtWeight(c.fUseFlatPtWeight),
fUseZWeight(c.fUseZWeight),
fUseNchWeight(c.fUseNchWeight),
+ fUseTrackletsWeight(c.fUseTrackletsWeight),
fNvar(c.fNvar),
fPartName(c.fPartName),
fDauNames(c.fDauNames),
fHistoMCNch(c.fHistoMCNch),
fResonantDecay(c.fResonantDecay),
fLctoV0bachelorOption(c.fLctoV0bachelorOption),
- fGenLctoV0bachelorOption(c.fGenLctoV0bachelorOption)
+ fGenLctoV0bachelorOption(c.fGenLctoV0bachelorOption),
+ fUseSelectionBit(c.fUseSelectionBit),
+ fPDGcode(c.fPDGcode),
+ fMultiplicityEstimator(c.fMultiplicityEstimator),
+ fRefMult(c.fRefMult),
+ fZvtxCorrectedNtrkEstimator(c.fZvtxCorrectedNtrkEstimator),
+ fIsPPData(c.fIsPPData),
+ fIsPPbData(c.fIsPPbData)
{
//
// Copy Constructor
//
+ for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
}
//___________________________________________________________________________
if (fCFManager) delete fCFManager ;
if (fHistEventsProcessed) delete fHistEventsProcessed ;
if (fCorrelation) delete fCorrelation ;
+ if (fListProfiles) delete fListProfiles;
if (fCuts) delete fCuts;
if (fFuncWeight) delete fFuncWeight;
if (fHistoMeasNch) delete fHistoMeasNch;
if (fHistoMCNch) delete fHistoMCNch;
+ for(Int_t i=0; i<4; i++) { if(fMultEstimatorAvg[i]) delete fMultEstimatorAvg[i]; }
}
//_________________________________________________________________________-
if (fDebug>1) printf("AliCFTaskVertexingHF::Init()");
if(fUseWeight && fUseZWeight) { AliFatal("Can not use at the same time pt and z-vtx weights, please choose"); return; }
- if(fUseWeight && fUseNchWeight) { AliFatal("Can not use at the same time pt and Nch weights, please choose"); return; }
+ if(fUseWeight && fUseNchWeight) { AliInfo("Beware, using at the same time pt and Nch weights, please check"); }
if(fUseNchWeight && !fHistoMCNch) { AliFatal("Need to pass the MC Nch distribution to use Nch weights"); return; }
- if(fUseNchWeight) CreateMeasuredNchHisto();
+ if(fUseNchWeight && !fHistoMeasNch) CreateMeasuredNchHisto();
AliRDHFCuts *copyfCuts = 0x0;
if (!fCuts){
switch (fDecayChannel){
case 2:{
+ fPDGcode = 421;
copyfCuts = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 21:{
+ fPDGcode = 413;
copyfCuts = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 22:{
+ fPDGcode = 4122;
copyfCuts = new AliRDHFCutsLctoV0(*(static_cast<AliRDHFCutsLctoV0*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 31:{
+ fPDGcode = 411;
copyfCuts = new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 32:{
+ fPDGcode = 4122;
copyfCuts = new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 33:{
+ fPDGcode = 431;
copyfCuts = new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
break;
}
case 4:{
+ fPDGcode = 421;
copyfCuts = new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
switch (fConfiguration) {
case kSnail: // slow configuration: all variables in
}
else{
AliFatal("Failing initializing AliRDHFCuts object - Exiting...");
- }
+ }
+
+ fListProfiles = new TList();
+ fListProfiles->SetOwner();
+ TString period[4];
+ Int_t nProfiles=4;
+
+ if (fIsPPbData) { //if pPb, use only two estimator histos
+ period[0] = "LHC13b"; period[1] = "LHC13c";
+ nProfiles = 2;
+ } else { // else assume pp (four histos for LHC10)
+ period[0] = "LHC10b"; period[1] = "LHC10c"; period[2] = "LHC10d"; period[3] = "LHC10e";
+ nProfiles = 4;
+ }
+
+ for(Int_t i=0; i<nProfiles; i++){
+ if(fMultEstimatorAvg[i]){
+ TProfile* hprof=new TProfile(*fMultEstimatorAvg[i]);
+ hprof->SetName(Form("ProfileTrkVsZvtx%s\n",period[i].Data()));
+ fListProfiles->Add(hprof);
+ }
+ }
+ PostData(5,fListProfiles);
return;
}
return;
}
+ fHistEventsProcessed->Fill(0.5);
+
Double_t* containerInput = new Double_t[fNvar];
Double_t* containerInputMC = new Double_t[fNvar];
Double_t zPrimVertex = aodVtx ->GetZ();
Double_t zMCVertex = mcHeader->GetVtxZ();
Int_t runnumber = aodEvent->GetRunNumber();
+
+ // Multiplicity definition with tracklets
+ Double_t nTracklets = 0;
+ Int_t nTrackletsEta10 = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.,1.);
+ Int_t nTrackletsEta16 = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.6,1.6);
+ nTracklets = (Double_t)nTrackletsEta10;
+ if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16 - nTrackletsEta10); }
+
+ // Apply the Ntracklets z-vtx data driven correction
+ if(fZvtxCorrectedNtrkEstimator) {
+ TProfile* estimatorAvg = GetEstimatorHistogram(aodEvent);
+ if(estimatorAvg) {
+ Int_t nTrackletsEta10Corr = AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta10,zPrimVertex,fRefMult);
+ Int_t nTrackletsEta16Corr = AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta16,zPrimVertex,fRefMult);
+ nTracklets = (Double_t)nTrackletsEta10Corr;
+ if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16Corr - nTrackletsEta10Corr); }
+ }
+ }
+
+
fWeight=1.;
if(fUseZWeight) fWeight *= GetZWeight(zMCVertex,runnumber);
if(fUseNchWeight){
Int_t nChargedMCPhysicalPrimary=AliVertexingHFUtils::GetGeneratedPhysicalPrimariesInEtaRange(mcArray,-1.0,1.0);
- fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
- AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
+ if(!fUseTrackletsWeight) fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
+ else fWeight *= GetNchWeight(nTracklets);
+ AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
}
+ Double_t eventWeight=fWeight;
if (TMath::Abs(zMCVertex) > fCuts->GetMaxVtxZ()){
AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fCuts->GetMaxVtxZ(), fCuts->GetMaxVtxZ()));
delete[] containerInput;
delete[] containerInputMC;
+ delete cfVtxHF;
+ return;
+ }
+
+ if(aodEvent->GetTriggerMask()==0 &&
+ (runnumber>=195344 && runnumber<=195677)){
+ AliDebug(3,"Event rejected because of null trigger mask");
+ delete[] containerInput;
+ delete[] containerInputMC;
+ delete cfVtxHF;
return;
}
delete[] containerInput;
delete[] containerInputMC;
delete [] trackCuts;
+ delete cfVtxHF;
return;
}
} else { // keep all centralities
fCuts->SetMaxCentrality(100.);
}
- Float_t centValue = fCuts->GetCentrality(aodEvent);
+ Float_t centValue = 0.;
+ if(!fIsPPData) centValue = fCuts->GetCentrality(aodEvent);
cfVtxHF->SetCentralityValue(centValue);
- // number of tracklets - multiplicity estimator
- Double_t multiplicity = (Double_t)(aodEvent->GetTracklets()->GetNumberOfTracklets()); // casted to double because the CF is filled with doubles
+ // multiplicity estimator with VZERO
+ Double_t vzeroMult=0;
+ AliAODVZERO *vzeroAOD = (AliAODVZERO*)aodEvent->GetVZEROData();
+ if(vzeroAOD) vzeroMult = vzeroAOD->GetMTotV0A() + vzeroAOD->GetMTotV0C();
+
+ Double_t multiplicity = nTracklets; // set to the Ntracklet estimator
+ if(fMultiplicityEstimator==kVZERO) { multiplicity = vzeroMult; }
+
+
cfVtxHF->SetMultiplicity(multiplicity);
+
+ // printf("Multiplicity estimator %d, value %2.2f\n",fMultiplicityEstimator,multiplicity);
for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
AliError("Failed casting particle from MC array!, Skipping particle");
continue;
}
+
+ //counting c quarks
+ cquarks += cfVtxHF->MCcquarkCounting(mcPart);
+
// check the MC-level cuts, must be the desidered particle
if (!fCFManager->CheckParticleCuts(0, mcPart)) {
AliDebug(2,"Check the MC-level cuts - not desidered particle");
}
cfVtxHF->SetMCCandidateParam(iPart);
- //counting c quarks
- cquarks += cfVtxHF->MCcquarkCounting(mcPart);
if (!(cfVtxHF->SetLabelArray())){
AliDebug(2,Form("Impossible to set the label array (decaychannel = %d)",fDecayChannel));
continue;
}
else{
- AliInfo(Form("Check on the family OK!!! (decaychannel = %d)",fDecayChannel));
+ AliDebug(2,Form("Check on the family OK!!! (decaychannel = %d)",fDecayChannel));
}
//Fill the MC container
if (fUseWeight){
if (fFuncWeight){ // using user-defined function
AliDebug(2,"Using function");
- fWeight = fFuncWeight->Eval(containerInputMC[0]);
+ fWeight = eventWeight*fFuncWeight->Eval(containerInputMC[0]);
}
else{ // using FONLL
AliDebug(2,"Using FONLL");
- fWeight = GetWeight(containerInputMC[0]);
+ fWeight = eventWeight*GetWeight(containerInputMC[0]);
}
AliDebug(2,Form("pt = %f, weight = %f",containerInputMC[0], fWeight));
}
if (fUseWeight){
if (fFuncWeight){ // using user-defined function
AliDebug(2, "Using function");
- fWeight = fFuncWeight->Eval(containerInput[0]);
+ fWeight = eventWeight*fFuncWeight->Eval(containerInput[0]);
}
else{ // using FONLL
AliDebug(2, "Using FONLL");
- fWeight = GetWeight(containerInput[0]);
+ fWeight = eventWeight*GetWeight(containerInput[0]);
}
AliDebug(2, Form("pt = %f, weight = %f",containerInput[0], fWeight));
}
Bool_t vtxCheck = fCuts->IsEventSelected(aodEvent);
+ // Selection on the filtering bit
+ Bool_t isBitSelected = kTRUE;
+ if(fDecayChannel==2) {
+ if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)) isBitSelected = kFALSE;
+ }else if(fDecayChannel==31){
+ if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDplusCuts)) isBitSelected = kFALSE;
+ }else if(fDecayChannel==33){
+ if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDsCuts)) isBitSelected = kFALSE;
+ }
+ if(!isBitSelected) continue;
+
+
+
if (recoStep && recoContFilled && vtxCheck){
fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed, fWeight) ;
icountReco++;
if(fAcceptanceUnf){
Double_t fill[4]; //fill response matrix
- Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
+ Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
if (bUnfolding) fCorrelation->Fill(fill);
}
Bool_t keepDs=ProcessDs(recoPidSelection);
if(keepDs) recoPidSelection=3;
} else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
- Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoAnalysisCuts);
- if (keepLctoV0bachelor) recoAnalysisCuts=3;
+ Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoPidSelection);
+ if (keepLctoV0bachelor) recoPidSelection=3;
}
Bool_t tempPid=(recoPidSelection == 3 || recoPidSelection == isPartOrAntipart);
AliDebug(3,"Reco PID cuts passed and container filled \n");
if(!fAcceptanceUnf){
Double_t fill[4]; //fill response matrix
- Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fill);
+ Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
if (bUnfolding) fCorrelation->Fill(fill);
}
}
fCountRecoPPR+= icountRecoPPR;
fCountRecoPID+= icountRecoPID;
- fHistEventsProcessed->Fill(0);
+ fHistEventsProcessed->Fill(1.5);
delete[] containerInput;
delete[] containerInputMC;
} else if (fDecayChannel==22) {
//nvarToPlot = 16;
titles = new TString[nvarToPlot];
- titles[0]="pT_Lc (GeV/c)";
- titles[1]="rapidity";
- titles[2]="phi (rad)";
- titles[3]="cosPAV0";
- titles[4]="onTheFlyStatusV0";
+ titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
+ titles[1]="y(#Lambda_{c})";
+ titles[2]="#varphi(#Lambda_{c}) [rad]";
+ titles[3]="onTheFlyStatusV0";
+ titles[4]="z_{vtx} [cm]";
titles[5]="centrality";
titles[6]="fake";
titles[7]="multiplicity";
- titles[8]="pT_bachelor (GeV/c)";
- titles[9]="pT_V0pos (GeV/c)";
- titles[10]="pT_V0neg (GeV/c)";
- titles[11]="invMassV0 (GeV/c2)";
- titles[12]="dcaV0 (nSigma)";
- titles[13]="c#tauV0 (#mum)";
- titles[14]="c#tau (#mum)";
- titles[15]="cosPA";
+ //titles[8]="pT(bachelor) [GeV/c]";
+ titles[8]="p(bachelor) [GeV/c]";
+ titles[9]="p_{T}(V0) [GeV/c]";
+ titles[10]="y(V0)";
+ titles[11]="#varphi(V0) [rad]";
+ titles[12]="m_{inv}(#pi^{+}#pi^{+}) [GeV/c^{2}]";
+ titles[13]="dcaV0 (nSigma)";
+ titles[14]="cosine pointing angle (V0)";
+ titles[15]="cosine pointing angle (#Lambda_{c})";
+ //titles[16]="c#tauV0 (#mum)";
+ //titles[17]="c#tau (#mum)";
} else {
//nvarToPlot = 12;
titles = new TString[nvarToPlot];
//nvarToPlot = 8;
titles = new TString[nvarToPlot];
if (fDecayChannel==22) {
- titles[0]="pT_candidate (GeV/c)";
- titles[1]="rapidity";
- titles[2]="phi (rad)";
- titles[3]="cosPAV0";
- titles[4]="onTheFlyStatusV0";
+ titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
+ titles[1]="y(#Lambda_{c})";
+ titles[2]="#varphi(#Lambda_{c}) [rad]";
+ titles[3]="onTheFlyStatusV0";
+ titles[4]="z_{vtx} [cm]";
titles[5]="centrality";
titles[6]="fake";
titles[7]="multiplicity";
file_projection->Close();
for (Int_t ih = 0; ih<3; ih++) delete [] h[ih];
delete [] h;
-
+ delete [] titles;
}
//TO BE SET BEFORE THE EXECUTION OF THE TASK
//
Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
-
- AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
- AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
- AliPIDResponse *localPIDResponse = (AliPIDResponse*)inputHandler->GetPIDResponse();
-
- if (fCuts->GetIsUsePID() && fDecayChannel==22) {
-
- fCuts->GetPidHF()->SetPidResponse(localPIDResponse);
- (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0pos()->SetPidResponse(localPIDResponse);
- (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0neg()->SetPidResponse(localPIDResponse);
- fCuts->GetPidHF()->SetOldPid(kFALSE);
- (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0pos()->SetOldPid(kFALSE);
- (dynamic_cast<AliRDHFCutsLctoV0*>(fCuts))->GetPidV0neg()->SetOldPid(kFALSE);
- }
//slot #1
OpenFile(1);
const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
- fHistEventsProcessed = new TH1I(nameoutput,"",1,0,1) ;
+ fHistEventsProcessed = new TH1I(nameoutput,"",2,0,2) ;
+ fHistEventsProcessed->GetXaxis()->SetBinLabel(1,"Events processed (all)");
+ fHistEventsProcessed->GetXaxis()->SetBinLabel(2,"Events analyzed (after selection)");
PostData(1,fHistEventsProcessed) ;
PostData(2,fCFManager->GetParticleContainer()) ;
}
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL276overLHC12a17a(){
+ // ad-hoc weight function from ratio of
+ // pt spectra from FONLL 2.76 TeV and
+ // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","[0]+[1]*TMath::Exp(-[2]*x)",0.,50.);
+ fFuncWeight->SetParameter(0,4.63891e-02);
+ fFuncWeight->SetParameter(1,1.51674e+01);
+ fFuncWeight->SetParameter(2,4.09941e-01);
+ fUseWeight=kTRUE;
+}
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromDataPbPb276overLHC12a17a(){
+ // ad-hoc weight function from ratio of
+ // D0 pt spectra in PbPb 2011 0-10% centrality and
+ // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","[0]+[1]/TMath::Power(x,[2])",0.05,50.);
+ fFuncWeight->SetParameter(0,1.43116e-02);
+ fFuncWeight->SetParameter(1,4.37758e+02);
+ fFuncWeight->SetParameter(2,3.08583);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL276overLHC12a17b(){
+ // weight function from the ratio of the LHC12a17b MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,50.);
+ fFuncWeight->SetParameters(1.92381e+01, 5.05055e+00, 1.05314e+01, 2.5, 1.88214e-03, 3.44871e+00, -9.74325e-02, 1.97671e+00, -3.21278e-01);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL276andBAMPSoverLHC12a17b(){
+ // weight function from the ratio of the LHC12a17b MC
+ // and FONLL calculations for pp data
+ // corrected by the BAMPS Raa calculation for 30-50% CC
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,50.);
+ fFuncWeight->SetParameters(6.10443e+00, 1.53487e+00, 1.99474e+00, 2.5, 5.51172e-03, 5.86590e+00, -5.46963e-01, 9.41201e-02, -1.64323e-01);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL5overLHC13d3(){
+ // weight function from the ratio of the LHC13d3 MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,30.);
+ fFuncWeight->SetParameters(2.94999e+00,3.47032e+00,2.81278e+00,2.5,1.93370e-02,3.86865e+00,-1.54113e-01,8.86944e-02,2.56267e-02);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL7overLHC10f6a(){
+ // weight function from the ratio of the LHC10f6a MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,40.);
+ fFuncWeight->SetParameters(2.41522e+01,4.92146e+00,6.72495e+00,2.5,6.15361e-03,4.78995e+00,-4.29135e-01,3.99421e-01,-1.57220e-02);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL7overLHC12a12(){
+ // weight function from the ratio of the LHC12a12 MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
+ fFuncWeight->SetParameters(1.31497e+01,3.30503e+00,3.45594e+00,2.5,2.28642e-02,1.42372e+00,2.32892e-04,5.21986e-02,-2.14236e-01,3.86200e+00);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL7overLHC12a12bis(){
+ // weight function from the ratio of the LHC12a12bis MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
+ fFuncWeight->SetParameters(6.54519e+00,2.74007e+00,2.48325e+00,2.5,1.61113e-01,-5.32546e-01,-3.75916e-04,2.38189e-01,-2.17561e-01,2.35975e+00);
+ fUseWeight=kTRUE;
+}
+
+//_________________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL7overLHC13e2fix(){
+ // weight function from the ratio of the LHC13e2fix MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
+ fFuncWeight->SetParameters(1.85862e+01,2.48171e+00,3.39356e+00,2.5,1.70426e-02,2.28453e+00,-4.57749e-02,5.84585e-02,-3.19719e-01,4.16789e+00);
+ fUseWeight=kTRUE;
+}
+
+//________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL5overLHC10f6a(){
+ // weight function from the ratio of the LHC10f6a MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,40.);
+ fFuncWeight->SetParameters(2.77730e+01,4.78942e+00,7.45378e+00,2.5,9.86255e-02,2.30120e+00,-4.16435e-01,3.43770e-01,-2.29380e-02);
+ fUseWeight=kTRUE;
+}
+
+//________________________________________________________________
+void AliCFTaskVertexingHF::SetPtWeightsFromFONLL276overLHC10f6a(){
+ // weight function from the ratio of the LHC10f6a MC
+ // and FONLL calculations for pp data
+ if(fFuncWeight) delete fFuncWeight;
+ fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,40.);
+ fFuncWeight->SetParameters(1.34412e+01,3.20068e+00,5.14481e+00,2.5,7.59405e-04,7.51821e+00,-3.93811e-01,2.16849e-02,-3.37768e-02,2.40308e+00);
+ fUseWeight=kTRUE;
+}
+
//_________________________________________________________________________
Double_t AliCFTaskVertexingHF::GetWeight(Float_t pt)
{
// calculates the Nch weight using the measured and generateed Nch distributions
//
if(nch<=0) return 0.;
+ if(!fHistoMeasNch || !fHistoMCNch) { AliError("Input histos to evaluate Nch weights missing"); return 0.; }
Double_t pMeas=fHistoMeasNch->GetBinContent(fHistoMeasNch->FindBin(nch));
Double_t pMC=fHistoMCNch->GetBinContent(fHistoMCNch->FindBin(nch));
- return pMeas/pMC;
+ Double_t weight = pMC>0 ? pMeas/pMC : 0.;
+ if(fUseTrackletsWeight) weight = pMC;
+ return weight;
}
//__________________________________________________________________________________________________
void AliCFTaskVertexingHF::CreateMeasuredNchHisto(){
// creates historgam with measured multiplcity distribution in pp 7 TeV collisions (from Eur. Phys. J. C (2010) 68: 345–354)
- Double_t nchbins[66]={0.50,1.50,2.50,3.50,4.50,5.50,6.50,7.50,8.50,9.50,
+ //
+ // for Nch > 70 the points were obtained with a double NBD distribution fit
+ // TF1 *fit1 = new TF1("fit1","[0]*(TMath::Gamma(x+[1])/(TMath::Gamma(x+1)*TMath::Gamma([1])))*(TMath::Power(([2]/[1]),x))*(TMath::Power((1+([2]/[1])),-x-[1]))"); fit1->SetParameter(0,1.);// normalization constant
+ // fit1->SetParameter(1,1.63); // k parameter
+ // fit1->SetParameter(2,12.8); // mean multiplicity
+ //
+ Double_t nchbins[82]={0.50,1.50,2.50,3.50,4.50,5.50,6.50,7.50,8.50,9.50,
10.50,11.50,12.50,13.50,14.50,15.50,16.50,17.50,18.50,19.50,
20.50,21.50,22.50,23.50,24.50,25.50,26.50,27.50,28.50,29.50,
30.50,31.50,32.50,33.50,34.50,35.50,36.50,37.50,38.50,39.50,
40.50,41.50,42.50,43.50,44.50,45.50,46.50,47.50,48.50,49.50,
50.50,51.50,52.50,53.50,54.50,55.50,56.50,57.50,58.50,59.50,
- 60.50,62.50,64.50,66.50,68.50,70.50};
- Double_t pch[65]={0.062011,0.072943,0.070771,0.067245,0.062834,0.057383,0.051499,0.04591,0.041109,0.036954,
+ 60.50,62.50,64.50,66.50,68.50,70.50,72.50,74.50,76.50,78.50,
+ 80.50,82.50,84.50,86.50,88.50,90.50,92.50,94.50,96.50,98.50,
+ 100.50,102.50};
+ Double_t pch[81]={0.062011,0.072943,0.070771,0.067245,0.062834,0.057383,0.051499,0.04591,0.041109,0.036954,
0.03359,0.030729,0.028539,0.026575,0.024653,0.0229,0.021325,0.019768,0.018561,0.017187,
0.01604,0.014836,0.013726,0.012576,0.011481,0.010393,0.009502,0.008776,0.008024,0.007452,
0.006851,0.006428,0.00594,0.005515,0.005102,0.00469,0.004162,0.003811,0.003389,0.003071,
0.002708,0.002422,0.002184,0.001968,0.00186,0.00165,0.001577,0.001387,0.001254,0.001118,
0.001037,0.000942,0.000823,0.000736,0.000654,0.000579,0.000512,0.00049,0.00045,0.000355,
- 0.000296,0.000265,0.000193,0.00016,0.000126};
+ 0.000296,0.000265,0.000193,0.00016,0.000126,0.0000851, 0.0000676,0.0000537,0.0000426, 0.0000338,
+ 0.0000268,0.0000213,0.0000166,0.0000133,0.0000106,0.00000837,0.00000662, 0.00000524,0.00000414, 0.00000327,
+ 0.00000258};
if(fHistoMeasNch) delete fHistoMeasNch;
- fHistoMeasNch=new TH1F("hMeaseNch","",65,nchbins);
- for(Int_t i=0; i<65; i++){
+ fHistoMeasNch=new TH1F("hMeaseNch","",81,nchbins);
+ for(Int_t i=0; i<81; i++){
fHistoMeasNch->SetBinContent(i+1,pch[i]);
fHistoMeasNch->SetBinError(i+1,0.);
}
}
return keep;
}
+
+//____________________________________________________________________________
+TProfile* AliCFTaskVertexingHF::GetEstimatorHistogram(const AliVEvent* event){
+ // Get Estimator Histogram from period event->GetRunNumber();
+ //
+ // If you select SPD tracklets in |eta|<1 you should use type == 1
+ //
+
+ Int_t runNo = event->GetRunNumber();
+ Int_t period = -1; // pp: 0-LHC10b, 1-LHC10c, 2-LHC10d, 3-LHC10e
+ // pPb: 0-LHC13b, 1-LHC13c
+
+ if (fIsPPbData) { // setting run numbers for LHC13 if pPb
+ if (runNo>195343 && runNo<195484) period = 0;
+ if (runNo>195528 && runNo<195678) period = 1;
+ if (period<0 || period>1) return 0;
+ } else { //else assume pp 2010
+ if(runNo>114930 && runNo<117223) period = 0;
+ if(runNo>119158 && runNo<120830) period = 1;
+ if(runNo>122373 && runNo<126438) period = 2;
+ if(runNo>127711 && runNo<130841) period = 3;
+ if(period<0 || period>3) return 0;
+ }
+
+ return fMultEstimatorAvg[period];
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff