Skip to content

read magnetic field from CCDB #245

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
jgrosseo merged 1 commit into from
Nov 3, 2021
Merged

read magnetic field from CCDB #245

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 10 additions & 10 deletions PWGCF/Core/PairCuts.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -61,7 +61,7 @@ class PairCuts
bool conversionCuts(T const& track1, T const& track2);

template <typename T>
bool twoTrackCut(T const& track1, T const& track2, int bSign);
bool twoTrackCut(T const& track1, T const& track2, float magField);

protected:
float mCuts[ParticlesLastEntry] = {-1};
Expand All @@ -80,7 +80,7 @@ class PairCuts
double getInvMassSquaredFast(T const& track1, double m0_1, T const& track2, double m0_2);

template <typename T>
float getDPhiStar(T const& track1, T const& track2, float radius, float bSign);
float getDPhiStar(T const& track1, T const& track2, float radius, float magField);
};

template <typename T>
Expand Down Expand Up @@ -109,28 +109,28 @@ bool PairCuts::conversionCuts(T const& track1, T const& track2)
}

template <typename T>
bool PairCuts::twoTrackCut(T const& track1, T const& track2, int bSign)
bool PairCuts::twoTrackCut(T const& track1, T const& track2, float magField)
{
// the variables & cut have been developed in Run 1 by the CF - HBT group
//
// Parameters:
// bSign: sign of B field
// magField: B field in kG

auto deta = track1.eta() - track2.eta();

// optimization
if (std::fabs(deta) < mTwoTrackDistance * 2.5 * 3) {
// check first boundaries to see if is worth to loop and find the minimum
float dphistar1 = getDPhiStar(track1, track2, mTwoTrackRadius, bSign);
float dphistar2 = getDPhiStar(track1, track2, 2.5, bSign);
float dphistar1 = getDPhiStar(track1, track2, mTwoTrackRadius, magField);
float dphistar2 = getDPhiStar(track1, track2, 2.5, magField);

const float kLimit = mTwoTrackDistance * 3;

if (std::fabs(dphistar1) < kLimit || std::fabs(dphistar2) < kLimit || dphistar1 * dphistar2 < 0) {
float dphistarminabs = 1e5;
float dphistarmin = 1e5;
for (Double_t rad = mTwoTrackRadius; rad < 2.51; rad += 0.01) {
float dphistar = getDPhiStar(track1, track2, rad, bSign);
float dphistar = getDPhiStar(track1, track2, rad, magField);

float dphistarabs = std::fabs(dphistar);

Expand All @@ -145,7 +145,7 @@ bool PairCuts::twoTrackCut(T const& track1, T const& track2, int bSign)
}

if (dphistarminabs < mTwoTrackDistance && std::fabs(deta) < mTwoTrackDistance) {
//LOGF(debug, "Removed track pair %ld %ld with %f %f %f %f %d %f %f %d %d", track1.index(), track2.index(), deta, dphistarminabs, track1.phi2(), track1.pt(), track1.sign(), track2.phi2(), track2.pt(), track2.sign(), bSign);
//LOGF(debug, "Removed track pair %ld %ld with %f %f %f %f %d %f %f %d %f", track1.index(), track2.index(), deta, dphistarminabs, track1.phi2(), track1.pt(), track1.sign(), track2.phi2(), track2.pt(), track2.sign(), magField);
return true;
}

Expand Down Expand Up @@ -308,7 +308,7 @@ double PairCuts::getInvMassSquaredFast(T const& track1, double m0_1, T const& tr
}

template <typename T>
float PairCuts::getDPhiStar(T const& track1, T const& track2, float radius, float bSign)
float PairCuts::getDPhiStar(T const& track1, T const& track2, float radius, float magField)
{
//
// calculates dphistar
Expand All @@ -322,7 +322,7 @@ float PairCuts::getDPhiStar(T const& track1, T const& track2, float radius, floa
auto pt2 = track2.pt();
auto charge2 = track2.sign();

float dphistar = phi1 - phi2 - charge1 * bSign * std::asin(0.075 * radius / pt1) + charge2 * bSign * std::asin(0.075 * radius / pt2);
float dphistar = phi1 - phi2 - charge1 * std::asin(0.015 * magField * radius / pt1) + charge2 * std::asin(0.015 * magField * radius / pt2);

if (dphistar > M_PI) {
dphistar = M_PI * 2 - dphistar;
Expand Down
50 changes: 32 additions & 18 deletions PWGCF/Tasks/correlations.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
#include "PWGCF/DataModel/CorrelationsDerived.h"
#include "PWGCF/Core/CorrelationContainer.h"
#include "PWGCF/Core/PairCuts.h"
#include "DataFormatsParameters/GRPObject.h"

#include <TH1F.h>
#include <cmath>
Expand Down Expand Up @@ -148,14 +149,30 @@ struct CorrelationTask {

// o2-ccdb-upload -p Users/jgrosseo/correlations/LHC15o -f /tmp/correction_2011_global.root -k correction

ccdb->setURL("http://ccdb-test.cern.ch:8080");
ccdb->setURL("http://alice-ccdb.cern.ch");
ccdb->setCaching(true);
ccdb->setLocalObjectValidityChecking();

long now = std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count();
ccdb->setCreatedNotAfter(now); // TODO must become global parameter from the train creation time
}

float getMagneticField(uint64_t timestamp)
{
// TODO done only once (and not per run). Will be replaced by CCDBConfigurable
static o2::parameters::GRPObject* grpo = nullptr;
if (grpo == nullptr) {
grpo = ccdb->getForTimeStamp<o2::parameters::GRPObject>("GLO/GRP/GRP", timestamp);
if (grpo == nullptr) {
LOGF(fatal, "GRP object not found for timestamp %llu", timestamp);
return 0;
}
LOGF(info, "Retrieved GRP for timestamp %llu with L3 current of %f", timestamp, grpo->getL3Current());
}
float l3current = grpo->getL3Current();
return l3current / 30000.0f * 5.0f;
}

template <typename TCollision, typename TTracks>
void fillQA(TCollision collision, float centrality, TTracks tracks)
{
Expand All @@ -180,7 +197,7 @@ struct CorrelationTask {
}

template <typename TTarget, typename TTracks>
void fillCorrelations(TTarget target, TTracks tracks1, TTracks tracks2, float centrality, float posZ, int bSign)
void fillCorrelations(TTarget target, TTracks tracks1, TTracks tracks2, float centrality, float posZ, float magField)
{
// Cache efficiency for particles (too many FindBin lookups)
float* efficiencyAssociated = nullptr;
Expand Down Expand Up @@ -228,7 +245,7 @@ struct CorrelationTask {
continue;
}

if (cfgTwoTrackCut > 0 && mPairCuts.twoTrackCut(track1, track2, bSign)) {
if (cfgTwoTrackCut > 0 && mPairCuts.twoTrackCut(track1, track2, magField)) {
continue;
}

Expand Down Expand Up @@ -277,38 +294,34 @@ struct CorrelationTask {

LOGF(info, "processSameAOD: Tracks for collision: %d | Vertex: %.1f | INT7: %d | V0M: %.1f", tracks.size(), collision.posZ(), collision.sel7(), collision.centV0M());

int bSign = 1; // TODO magnetic field from CCDB
const auto centrality = collision.centV0M();

if (fillCollisionAOD(same, collision, centrality) == false) {
return;
}
registry.fill(HIST("eventcount"), -2);
fillQA(collision, centrality, tracks);
fillCorrelations(same, tracks, tracks, centrality, collision.posZ(), bSign);
fillCorrelations(same, tracks, tracks, centrality, collision.posZ(), getMagneticField(bc.timestamp()));
}
PROCESS_SWITCH(CorrelationTask, processSameAOD, "Process same event on AOD", true);

void processSameDerived(soa::Filtered<aod::CFCollisions>::iterator const& collision, soa::Filtered<aod::CFTracks> const& tracks)
{
LOGF(info, "processSameDerived: Tracks for collision: %d | Vertex: %.1f | V0M: %.1f", tracks.size(), collision.posZ(), collision.centV0M());

int bSign = 1; // TODO magnetic field from CCDB
const auto centrality = collision.centV0M();

same->fillEvent(centrality, CorrelationContainer::kCFStepReconstructed);
registry.fill(HIST("eventcount"), -2);
fillQA(collision, centrality, tracks);
fillCorrelations(same, tracks, tracks, centrality, collision.posZ(), bSign);
fillCorrelations(same, tracks, tracks, centrality, collision.posZ(), getMagneticField(collision.timestamp()));
}
PROCESS_SWITCH(CorrelationTask, processSameDerived, "Process same event on derived data", false);

void processMixedAOD(soa::Filtered<soa::Join<aod::Collisions, aod::Hashes, aod::EvSels, aod::CentV0Ms>>& collisions, aodTracks const& tracks)
void processMixedAOD(soa::Filtered<soa::Join<aod::Collisions, aod::Hashes, aod::EvSels, aod::CentV0Ms>>& collisions, aodTracks const& tracks, aod::BCsWithTimestamps const&)
{
// TODO loading of efficiency histogram missing here, because it will happen somehow in the CCDBConfigurable

int bSign = 1; // TODO magnetic field from CCDB

collisions.bindExternalIndices(&tracks);
auto tracksTuple = std::make_tuple(tracks);
AnalysisDataProcessorBuilder::GroupSlicer slicer(collisions, tracksTuple);
Expand Down Expand Up @@ -345,10 +358,12 @@ struct CorrelationTask {
auto tracks2 = std::get<aodTracks>(it2.associatedTables());
tracks2.bindExternalIndices(&collisions);

auto bc = collision1.bc_as<aod::BCsWithTimestamps>();

// LOGF(info, "Tracks: %d and %d entries", tracks1.size(), tracks2.size());

// TODO mixed event weight missing
fillCorrelations(mixed, tracks1, tracks2, collision1.centV0M(), collision1.posZ(), bSign);
fillCorrelations(mixed, tracks1, tracks2, collision1.centV0M(), collision1.posZ(), getMagneticField(bc.timestamp()));
}
}
PROCESS_SWITCH(CorrelationTask, processMixedAOD, "Process mixed events on AOD", true);
Expand All @@ -357,8 +372,6 @@ struct CorrelationTask {
{
// TODO loading of efficiency histogram missing here, because it will happen somehow in the CCDBConfigurable

int bSign = 1; // TODO magnetic field from CCDB

collisions.bindExternalIndices(&tracks);
auto tracksTuple = std::make_tuple(tracks);
AnalysisDataProcessorBuilder::GroupSlicer slicer(collisions, tracksTuple);
Expand Down Expand Up @@ -393,19 +406,20 @@ struct CorrelationTask {

// LOGF(info, "Tracks: %d and %d entries", tracks1.size(), tracks2.size());

fillCorrelations(mixed, tracks1, tracks2, collision1.centV0M(), collision1.posZ(), bSign);
fillCorrelations(mixed, tracks1, tracks2, collision1.centV0M(), collision1.posZ(), getMagneticField(collision1.timestamp()));
}
}
PROCESS_SWITCH(CorrelationTask, processMixedDerived, "Process mixed events on derived data", false);

// Version with combinations
void processWithCombinations(soa::Join<aod::Collisions, aod::CentV0Ms>::iterator const& collision, soa::Filtered<aod::Tracks> const& tracks)
void processWithCombinations(soa::Join<aod::Collisions, aod::CentV0Ms>::iterator const& collision, aod::BCsWithTimestamps const&, soa::Filtered<aod::Tracks> const& tracks)
{
LOGF(info, "Tracks for collision (Combination run): %d", tracks.size());

const auto centrality = collision.centV0M();
// TODO will go to CCDBConfigurable
auto bc = collision.bc_as<aod::BCsWithTimestamps>();

int bSign = 1; // TODO magnetic field from CCDB
const auto centrality = collision.centV0M();

for (auto track1 = tracks.begin(); track1 != tracks.end(); ++track1) {

Expand Down Expand Up @@ -436,7 +450,7 @@ struct CorrelationTask {
continue;
}

if (cfgTwoTrackCut > 0 && mPairCuts.twoTrackCut(track1, track2, bSign)) {
if (cfgTwoTrackCut > 0 && mPairCuts.twoTrackCut(track1, track2, getMagneticField(bc.timestamp()))) {
continue;
}

Expand Down