Skip to content

Updated ITSTrackSimulation to run with O2DPG #1115

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
iravasen merged 2 commits into from
Mar 29, 2022
Merged

Updated ITSTrackSimulation to run with O2DPG #1115

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
ad48553
Updated ITSTrackSimulation to run with O2DPG
JianLIUhep Mar 24, 2022
88cc537
returned task parameters, added AngularDistribution of track histogram
JianLIUhep Mar 29, 2022
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
10 changes: 7 additions & 3 deletions Modules/ITS/include/ITS/ITSTrackSimTask.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -105,10 +105,14 @@ class ITSTrackSimTask : public TaskInterface
TH1F* hTrackImpactTransvFake;
TH1F* hTrackImpactTransvValid;

std::string mRunNumber;
std::string mRunNumberPath;
std::string mMCKinePath;
TH1D* hPrimaryReco_pt;
TH1D* hPrimaryGen_pt;

TH2D* hAngularDistribution;

int mRunNumber = 0;
std::string mO2GrpPath;
std::string mCollisionsContextPath;

o2::its::GeometryTGeo* mGeom;

Expand Down
7 changes: 3 additions & 4 deletions Modules/ITS/itsTrackSim.json
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,10 +37,9 @@
},
"location" : "remote",
"taskParameters" : {
"runNumberPath" : "/home/its/QC/workdir/infiles/RunNumber.dat",
"MCKinePath" : "./o2sim_Kine.root",
"O2GrpPath" : "./o2sim_grp.root"
}
"O2GrpPath" : "./o2sim_grp.root",
"collisionsContextPath": "./collisioncontext.root"
}

}
}
Expand Down
152 changes: 87 additions & 65 deletions Modules/ITS/src/ITSTrackSimTask.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -87,15 +87,18 @@ ITSTrackSimTask::~ITSTrackSimTask()

delete hTrackImpactTransvFake;
delete hTrackImpactTransvValid;

delete hPrimaryReco_pt;
delete hPrimaryGen_pt;

delete hAngularDistribution;
}

void ITSTrackSimTask::initialize(o2::framework::InitContext& /*ctx*/)
{
ILOG(Info, Support) << "initialize ITSTrackSimTask" << ENDM;

mRunNumberPath = mCustomParameters["runNumberPath"];
mMCKinePath = mCustomParameters["MCKinePath"];
mO2GrpPath = mCustomParameters["o2GrpPath"];
mCollisionsContextPath = mCustomParameters["collisionsContextPath"];

createAllHistos();
publishHistos();
Expand All @@ -108,37 +111,25 @@ void ITSTrackSimTask::initialize(o2::framework::InitContext& /*ctx*/)
mGeom = o2::its::GeometryTGeo::Instance();
}

void ITSTrackSimTask::startOfActivity(Activity& /*activity*/)
void ITSTrackSimTask::startOfActivity(Activity& activity)
{
mRunNumber = activity.mId;
ILOG(Info, Support) << "startOfActivity" << ENDM;
}

void ITSTrackSimTask::startOfCycle()
{
TFile* file = new TFile("o2sim_Kine.root");
ILOG(Info, Support) << "startOfCycle" << ENDM;
}

void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)
{
ILOG(Info, Support) << "START DOING QC General" << ENDM;

TFile* file = new TFile(mMCKinePath.c_str()); // MC Kinematics files is used to get particle level information
TTree* mcTree = (TTree*)file->Get("o2sim");
mcTree->SetBranchStatus("MCTrack*", 1);
mcTree->SetBranchStatus("MCEventHeader.*", 1);

auto mcArr = new std::vector<o2::MCTrack>;
auto mcHeader = new o2::dataformats::MCEventHeader;

mcTree->SetBranchAddress("MCTrack", &mcArr);
mcTree->SetBranchAddress("MCEventHeader.", &mcHeader);

info.resize(mcTree->GetEntriesFast());
for (int i = 0; i < mcTree->GetEntriesFast(); ++i) {
if (!mcTree->GetEvent(i))
continue;
info[i].resize(mcArr->size());
o2::steer::MCKinematicsReader reader(mCollisionsContextPath.c_str());
info.resize(reader.getNEvents(0));
for (int i = 0; i < reader.getNEvents(0); ++i) {
std::vector<MCTrack> const& mcArr = reader.getTracks(i);
info[i].resize(mcArr.size());
}

auto clusArr = ctx.inputs().get<gsl::span<o2::itsmft::CompClusterExt>>("compclus"); // used to get hit information
Expand All @@ -152,7 +143,7 @@ void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)

int TrackID = lab.getTrackID();

if (TrackID < 0 || TrackID >= mcTree->GetEntriesFast()) {
if (TrackID < 0) {
continue;
}

Expand All @@ -178,13 +169,13 @@ void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)
ok |= 0b1000000;
}

for (int i = 0; i < mcTree->GetEntriesFast(); ++i) { // filling denominator
for (int i = 0; i < reader.getNEvents(0); ++i) {
std::vector<MCTrack> const& mcArr = reader.getTracks(i);
auto mcHeader = reader.getMCEventHeader(0, i); // SourceID=0 for ITS

if (!mcTree->GetEvent(i))
continue;
for (int mc = 0; mc < mcArr->size(); mc++) {
for (int mc = 0; mc < mcArr.size(); mc++) {

const auto& mcTrack = (*mcArr)[mc];
const auto& mcTrack = (mcArr)[mc];

info[i][mc].isFilled = false;
if (mcTrack.Vx() * mcTrack.Vx() + mcTrack.Vy() * mcTrack.Vy() > 1)
Expand All @@ -196,14 +187,17 @@ void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)
if (info[i][mc].clusters != 0b1111111)
continue;

Double_t distance = sqrt(pow(mcHeader->GetX() - mcTrack.Vx(), 2) + pow(mcHeader->GetY() - mcTrack.Vy(), 2) + pow(mcHeader->GetZ() - mcTrack.Vz(), 2));
Double_t distance = sqrt(pow(mcHeader.GetX() - mcTrack.Vx(), 2) + pow(mcHeader.GetY() - mcTrack.Vy(), 2) + pow(mcHeader.GetZ() - mcTrack.Vz(), 2));
info[i][mc].isFilled = true;
info[i][mc].r = distance;
info[i][mc].pt = mcTrack.GetPt();
info[i][mc].eta = mcTrack.GetEta();
info[i][mc].phi = mcTrack.GetPhi();
info[i][mc].z = mcTrack.Vz();
info[i][mc].isPrimary = mcTrack.isPrimary();
if (mcTrack.isPrimary()) {
hPrimaryGen_pt->Fill(mcTrack.GetPt());
}

hDenTrue_r->Fill(distance);
hDenTrue_pt->Fill(mcTrack.GetPt());
Expand All @@ -227,6 +221,8 @@ void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)
Float_t vx = 0., vy = 0., vz = 0.; // Assumed primary vertex at 0,0,0
track.getImpactParams(vx, vy, vz, bz, ip);

hAngularDistribution->Fill(track.getEta(), track.getPhi());

if (info[MCinfo.getEventID()][MCinfo.getTrackID()].isFilled) {
if (MCinfo.isFake()) {

Expand All @@ -235,51 +231,46 @@ void ITSTrackSimTask::monitorData(o2::framework::ProcessingContext& ctx)
hNumRecoFake_eta->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].eta);
hNumRecoFake_z->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].z);
hNumRecoFake_r->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].r);
if (info[MCinfo.getEventID()][MCinfo.getTrackID()].isPrimary)
if (info[MCinfo.getEventID()][MCinfo.getTrackID()].isPrimary) {
hTrackImpactTransvFake->Fill(ip[0]);
hPrimaryReco_pt->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].pt);
}
} else {
hNumRecoValid_pt->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].pt);
hNumRecoValid_phi->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].phi);
hNumRecoValid_eta->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].eta);
hNumRecoValid_z->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].z);
hNumRecoValid_r->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].r);
if (info[MCinfo.getEventID()][MCinfo.getTrackID()].isPrimary)
if (info[MCinfo.getEventID()][MCinfo.getTrackID()].isPrimary) {
hTrackImpactTransvValid->Fill(ip[0]);
hPrimaryReco_pt->Fill(info[MCinfo.getEventID()][MCinfo.getTrackID()].pt);
}
}
}
}

hFakeTrack_pt->Divide(hNumRecoFake_pt, hDenTrue_pt, 1, 1);
hEfficiency_pt->Divide(hNumRecoValid_pt, hDenTrue_pt, 1, 1);
hFakeTrack_pt->Divide(hNumRecoFake_pt, hDenTrue_pt, 1, 1, "B");
hEfficiency_pt->Divide(hNumRecoValid_pt, hDenTrue_pt, 1, 1, "B");

hFakeTrack_phi->Divide(hNumRecoFake_phi, hDenTrue_phi, 1, 1);
hEfficiency_phi->Divide(hNumRecoValid_phi, hDenTrue_phi, 1, 1);
hFakeTrack_phi->Divide(hNumRecoFake_phi, hDenTrue_phi, 1, 1, "B");
hEfficiency_phi->Divide(hNumRecoValid_phi, hDenTrue_phi, 1, 1, "B");

hFakeTrack_eta->Divide(hNumRecoFake_eta, hDenTrue_eta, 1, 1);
hEfficiency_eta->Divide(hNumRecoValid_eta, hDenTrue_eta, 1, 1);
hFakeTrack_eta->Divide(hNumRecoFake_eta, hDenTrue_eta, 1, 1, "B");
hEfficiency_eta->Divide(hNumRecoValid_eta, hDenTrue_eta, 1, 1, "B");

hFakeTrack_r->Divide(hNumRecoFake_r, hDenTrue_r, 1, 1);
hEfficiency_r->Divide(hNumRecoValid_r, hDenTrue_r, 1, 1);
hFakeTrack_r->Divide(hNumRecoFake_r, hDenTrue_r, 1, 1, "B");
hEfficiency_r->Divide(hNumRecoValid_r, hDenTrue_r, 1, 1, "B");

hFakeTrack_z->Divide(hNumRecoFake_z, hDenTrue_z, 1, 1);
hEfficiency_z->Divide(hNumRecoValid_z, hDenTrue_z, 1, 1);
hFakeTrack_z->Divide(hNumRecoFake_z, hDenTrue_z, 1, 1, "B");
hEfficiency_z->Divide(hNumRecoValid_z, hDenTrue_z, 1, 1, "B");
}

void ITSTrackSimTask::endOfCycle()
{

std::ifstream runNumberFile(mRunNumberPath.c_str());
if (runNumberFile) {

std::string runNumber;
runNumberFile >> runNumber;
if (runNumber != mRunNumber) {
for (unsigned int iObj = 0; iObj < mPublishedObjects.size(); iObj++)
getObjectsManager()->addMetadata(mPublishedObjects.at(iObj)->GetName(), "Run", runNumber);
mRunNumber = runNumber;
}
ILOG(Info, Support) << "endOfCycle" << ENDM;
}
for (unsigned int iObj = 0; iObj < mPublishedObjects.size(); iObj++)
getObjectsManager()->addMetadata(mPublishedObjects.at(iObj)->GetName(), "Run", std::to_string(mRunNumber));
ILOG(Info, Support) << "endOfCycle" << ENDM;
}

void ITSTrackSimTask::endOfActivity(Activity& /*activity*/)
Expand Down Expand Up @@ -322,6 +313,11 @@ void ITSTrackSimTask::reset()

hTrackImpactTransvValid->Reset();
hTrackImpactTransvFake->Reset();

hPrimaryGen_pt->Reset();
hPrimaryReco_pt->Reset();

hAngularDistribution->Reset();
}

void ITSTrackSimTask::createAllHistos()
Expand All @@ -334,11 +330,13 @@ void ITSTrackSimTask::createAllHistos()

hEfficiency_pt = new TH1D("efficiency_pt", ";#it{p}_{T} (GeV/#it{c});t{p}_{T}Efficiency", nb, xbins);
hEfficiency_pt->SetTitle("#it{p}_{T} efficiency of tracking");
hEfficiency_pt->SetBit(TH1::kIsAverage);
addObject(hEfficiency_pt);
formatAxes(hEfficiency_pt, "#it{p}_{T} (GeV/#it{c})", "Efficiency", 1, 1.10);

hFakeTrack_pt = new TH1D("faketrack_pt", ";#it{p}_{T} (GeV/#it{c});Fake-track rate", nb, xbins);
hFakeTrack_pt->SetTitle("#it{p}_{T} fake-track rate");
hFakeTrack_pt->SetBit(TH1::kIsAverage);
addObject(hFakeTrack_pt);
formatAxes(hFakeTrack_pt, "#it{p}_{T} (GeV/#it{c})", "Fake-track rate", 1, 1.10);

Expand All @@ -348,11 +346,13 @@ void ITSTrackSimTask::createAllHistos()

hEfficiency_phi = new TH1D("efficiency_phi", ";#phi;Efficiency", 60, 0, TMath::TwoPi());
hEfficiency_phi->SetTitle("#phi efficiency of tracking");
hEfficiency_phi->SetBit(TH1::kIsAverage);
addObject(hEfficiency_phi);
formatAxes(hEfficiency_phi, "#phi", "Efficiency", 1, 1.10);

hFakeTrack_phi = new TH1D("faketrack_phi", ";#phi;Fake-track rate", 60, 0, TMath::TwoPi());
hFakeTrack_phi->SetTitle("#phi fake-track rate");
hFakeTrack_phi->SetBit(TH1::kIsAverage);
addObject(hFakeTrack_phi);
formatAxes(hFakeTrack_phi, "#phi", "Fake-track rate", 1, 1.10);

Expand All @@ -362,55 +362,77 @@ void ITSTrackSimTask::createAllHistos()

hEfficiency_eta = new TH1D("efficiency_eta", ";#eta;Efficiency", 30, -1.5, 1.5);
hEfficiency_eta->SetTitle("#eta efficiency of tracking");
hEfficiency_eta->SetBit(TH1::kIsAverage);
addObject(hEfficiency_eta);
formatAxes(hEfficiency_eta, "#eta", "Efficiency", 1, 1.10);

hFakeTrack_eta = new TH1D("faketrack_eta", ";#eta;Fake-track rate", 30, -1.5, 1.5);
hFakeTrack_eta->SetTitle("#eta fake-track rate");
hFakeTrack_eta->SetBit(TH1::kIsAverage);
addObject(hFakeTrack_eta);
formatAxes(hFakeTrack_eta, "#eta", "Fake-track rate", 1, 1.10);

hNumRecoValid_eta = new TH1D("NumRecoValid_eta", "", 30, -1.5, 1.5);
hNumRecoFake_eta = new TH1D("NumRecoFake_eta", "", 30, -1.5, 1.5);
hDenTrue_eta = new TH1D("DenTrueMC_eta", "", 30, -1.5, 1.5);

hEfficiency_r = new TH1D("efficiency_r", ";r;Efficiency", 50, 0, 5);
hEfficiency_r = new TH1D("efficiency_r", ";r (cm);Efficiency", 50, 0, 5);
hEfficiency_r->SetTitle("r efficiency of tracking");
hEfficiency_r->SetBit(TH1::kIsAverage);
addObject(hEfficiency_r);
formatAxes(hEfficiency_r, "r", "Efficiency", 1, 1.10);
formatAxes(hEfficiency_r, "r (cm)", "Efficiency", 1, 1.10);

hFakeTrack_r = new TH1D("faketrack_r", ";r;Fake-track rate", 50, 0, 5);
hFakeTrack_r = new TH1D("faketrack_r", ";r (cm);Fake-track rate", 50, 0, 5);
hFakeTrack_r->SetTitle("r fake-track rate");
hFakeTrack_r->SetBit(TH1::kIsAverage);
addObject(hFakeTrack_r);
formatAxes(hFakeTrack_r, "r", "Fake-track rate", 1, 1.10);
formatAxes(hFakeTrack_r, "r (cm)", "Fake-track rate", 1, 1.10);

hNumRecoValid_r = new TH1D("NumRecoValid_r", "", 50, 0, 5);
hNumRecoFake_r = new TH1D("NumRecoFake_r", "", 50, 0, 5);
hDenTrue_r = new TH1D("DenTrueMC_r", "", 50, 0, 5);

hEfficiency_z = new TH1D("efficiency_z", ";z;Efficiency", 100, -5, 5);
hEfficiency_z = new TH1D("efficiency_z", ";z (cm);Efficiency", 100, -5, 5);
hEfficiency_z->SetTitle("z efficiency of tracking");
hEfficiency_z->SetBit(TH1::kIsAverage);
addObject(hEfficiency_z);
formatAxes(hEfficiency_z, "z", "Efficiency", 1, 1.10);
formatAxes(hEfficiency_z, "z (cm)", "Efficiency", 1, 1.10);

hFakeTrack_z = new TH1D("faketrack_z", ";z;Fake-track rate", 100, -5, 5);
hFakeTrack_z = new TH1D("faketrack_z", ";z (cm);Fake-track rate", 100, -5, 5);
hFakeTrack_z->SetTitle("z fake-track rate");
hFakeTrack_z->SetBit(TH1::kIsAverage);
addObject(hFakeTrack_z);
formatAxes(hFakeTrack_z, "z", "Fake-track rate", 1, 1.10);
formatAxes(hFakeTrack_z, "z (cm)", "Fake-track rate", 1, 1.10);

hNumRecoValid_z = new TH1D("NumRecoValid_z", "", 100, -5, 5);
hNumRecoFake_z = new TH1D("NumRecoFake_z", "", 100, -5, 5);
hDenTrue_z = new TH1D("DenTrueMC_z", "", 100, -5, 5);

hTrackImpactTransvValid = new TH1F("ImpactTransvVaild", "Transverse impact parameter for valid tracks; D (cm)", 60, -0.1, 0.1);
hTrackImpactTransvValid->SetTitle("Transverse impact parameter distribution of valid tracks");
hTrackImpactTransvValid = new TH1F("ImpactTransvVaild", "Transverse impact parameter for valid primary tracks; D (cm)", 60, -0.1, 0.1);
hTrackImpactTransvValid->SetTitle("Transverse impact parameter distribution of valid primary tracks");
addObject(hTrackImpactTransvValid);
formatAxes(hTrackImpactTransvValid, "D (cm)", "counts", 1, 1.10);

hTrackImpactTransvFake = new TH1F("ImpactTransvFake", "Transverse impact parameter for fake tracks; D (cm)", 60, -0.1, 0.1);
hTrackImpactTransvFake->SetTitle("Transverse impact parameter distribution of fake tracks");
hTrackImpactTransvFake = new TH1F("ImpactTransvFake", "Transverse impact parameter for fake primary tracks; D (cm)", 60, -0.1, 0.1);
hTrackImpactTransvFake->SetTitle("Transverse impact parameter distribution of fake primary tracks");
addObject(hTrackImpactTransvFake);
formatAxes(hTrackImpactTransvFake, "D (cm)", "counts", 1, 1.10);

hPrimaryGen_pt = new TH1D("PrimaryGen_pt", ";#it{p}_{T} (GeV/#it{c});counts", nb, xbins);
hPrimaryGen_pt->SetTitle("#it{p}_{T} of primary generated particles");
addObject(hPrimaryGen_pt);
formatAxes(hPrimaryGen_pt, "#it{p}_{T} (GeV/#it{c})", "counts", 1, 1.10);

hPrimaryReco_pt = new TH1D("PrimaryReco_pt", ";#it{p}_{T} (GeV/#it{c});counts", nb, xbins);
hPrimaryReco_pt->SetTitle("#it{p}_{T} of primary reconstructed particles");
addObject(hPrimaryReco_pt);
formatAxes(hPrimaryReco_pt, "#it{p}_{T} (GeV/#it{c})", "counts", 1, 1.10);

hAngularDistribution = new TH2D("AngularDistribution", "AngularDistribution", 30, -1.5, 1.5, 60, 0, TMath::TwoPi());
hAngularDistribution->SetTitle("AngularDistribution");
addObject(hAngularDistribution);
formatAxes(hAngularDistribution, "#eta", "#phi", 1, 1.10);
hAngularDistribution->SetStats(0);
}

void ITSTrackSimTask::addObject(TObject* aObject)
Expand Down