BmnRootIterator.cxx

Go to the documentation of this file.

/*
    BM@N alignment routine
    BM@N experiment at NICA complex, JINR, 2025
 
    Department: Math & Soft Group of HEP lab
    Author:     Igor Polev, polev@jinr.ru
 
    BmnRootIterator class implementation
*/
 
#include "BmnRootIterator.h"
 
#include "BmnAlignDefines.h"
#include "fairlogger/Logger.h"
 
BmnRootIterator::BmnRootIterator(const char* filesPathList, const BmnDetectorModel& detModel, BmnSimpleProgressBar* pb)
    : BmnDataIterator(pb)
    , fpDetModel(&detModel)
    , fpData(new BmnDataReader(filesPathList))
    , fpReader(&fpData->Reader())
{
    // turn off unneeded branches
    // notice: turning off branches works weirdly and does not give any speedup
    // but let it be for now - maybe it will give something later
    fpData->Chain().SetBranchStatus("*", kFALSE);
    fpData->Chain().SetBranchStatus(BMN_TRACKS_BRACH_NAME "*", kTRUE);
    fpData->Chain().SetBranchStatus(BMN_HITS_BRACH_NAME "*", kTRUE);
 
    // Tryed to use TTreeReaderArray but it does not work
    fpRootTracks = new TTreeReaderValue<TClonesArray>{*fpReader, BMN_TRACKS_BRACH_NAME};
    fpRootHits = new TTreeReaderValue<TClonesArray>{*fpReader, BMN_HITS_BRACH_NAME};
 
    fProgressByTracks = kFALSE;   // progress by events until tracks are counted
    fTotalEvents = fpData->Chain().GetEntries();
    SetTickProgress(1.0f / fTotalEvents);
 
    fInitialized = kTRUE;
}
 
BmnRootIterator::~BmnRootIterator()
{
    if (fpRootTracks)
        delete fpRootTracks;
    if (fpRootHits)
        delete fpRootHits;
    if (fpData)
        delete fpData;
}
 
void BmnRootIterator::CountElements()
{
    // Count all tracks and hits that matches minimum number
    // of hits per track and maximum allowed chi2 per hit,
    // Filtration is built in ResetAll() and NextTrack() directly.
    LOG(info) << "    minimum hits per track = " << BMN_MIN_HITS_PER_TRACK;
    LOG(info) << "    maximum track chi2/NDF = " << BMN_MAX_CHI2_PER_NDF;
    LOG(info) << "    counting suitable tracks and other elements...";
 
    auto events{fpData->Chain().GetEntries()};
    SetTickProgress(1.0f / events);
    fProgressByTracks = kFALSE;   // progress by events
    ProgressReset();
    fTotalTracks = 0;
    fTotalHits = 0;
    ResetAll();
    while (!EndOfTracks()) {
        fTotalTracks++;
        fTotalHits += fHitsInTrack;
        NextTrack();
    }
    ProgressClear();
    fProgressByTracks = kTRUE;
    SetTickProgress(1.0f / fTotalTracks);
 
    LOG(info) << "    " << fTotalTracks << " total suitable tracks found";
    LOG(info) << "    " << fTotalHits << " total hits counted";
    LOG(info) << "    " << events << " events processed";
}
 
void BmnRootIterator::ResetHits()
{
    fHitsInTrack = fpCurrentTrack->GetNStsHits();
    if (fHitsInTrack > 0) {
        fpCurrentHit = (Hit_t*)(*fpRootHits)->At(fpCurrentTrack->GetStsHitIndex(0));
        fHitIdx = 0;
        fEndOfHits = kFALSE;
    } else {
        fpCurrentHit = nullptr;
        fEndOfHits = kTRUE;
    }
}
 
void BmnRootIterator::ResetAll()
{
    if (!fInitialized)
        return;
 
    UnsetProperties();
    fpReader->Restart();
    Bool_t event = fpReader->Next();
    if (!event)
        return;   // no events at all
    if (!fProgressByTracks)
        ProgressTick();
 
    // search for event with tracks
    fTracksInEvent = (*fpRootTracks)->GetEntries();
    while (event && fTracksInEvent <= 0) {
        event = fpReader->Next();
        if (!fProgressByTracks)
            ProgressTick();
        fTracksInEvent = (*fpRootTracks)->GetEntries();
    }
    if (!event)
        return;   // no events with tracks
 
    fpCurrentTrack = (Track_t*)(*fpRootTracks)->At(0);
    fTrackIdx = 0;
    fEndOfTracks = kFALSE;
    // skip track with too few hits or too big chi2 per NDF
    if (BMN_MIN_HITS_PER_TRACK > fpCurrentTrack->GetNStsHits()
        || BMN_MAX_CHI2_PER_NDF < fpCurrentTrack->GetChi2() / fpCurrentTrack->GetNDF())
        NextTrack();
    else
        ResetHits();
}
 
void BmnRootIterator::NextTrack()
{
    if (!fInitialized || fEndOfTracks)
        return;
    Bool_t search{kTRUE};   // skip track with too few hits or too big chi2 per hit
    while (search) {
        if (++fTrackIdx >= fTracksInEvent) {
            Bool_t event{fpReader->Next()};
            if (!fProgressByTracks)
                ProgressTick();
            if (event) {
                fTracksInEvent = (*fpRootTracks)->GetEntries();
                while (event && fTracksInEvent <= 0) {
                    event = fpReader->Next();
                    if (!fProgressByTracks)
                        ProgressTick();
                    fTracksInEvent = (*fpRootTracks)->GetEntries();
                }
            }
            // (!event) must be checked againg, don't use else
            if (!event) {   // no more events - no more tracks
                UnsetProperties();
                return;
            }
            fTrackIdx = 0;
        }
        fpCurrentTrack = (Track_t*)(*fpRootTracks)->At(fTrackIdx);
        search = (BMN_MIN_HITS_PER_TRACK > fpCurrentTrack->GetNStsHits()
                  || BMN_MAX_CHI2_PER_NDF < fpCurrentTrack->GetChi2() / fpCurrentTrack->GetNDF());
    }
    ResetHits();
    if (fProgressByTracks)
        ProgressTick();
}
 
void BmnRootIterator::NextHit()
{
    if (!fInitialized || fEndOfHits)
        return;
    if (++fHitIdx >= fHitsInTrack) {
        fpCurrentHit = nullptr;
        fEndOfHits = kTRUE;
    } else
        fpCurrentHit = (Hit_t*)(*fpRootHits)->At(fpCurrentTrack->GetStsHitIndex(fHitIdx));
}
 
Double_t BmnRootIterator::HitX() const
{
    if (!fInitialized || fEndOfHits)
        return 0.0;
    return fpCurrentHit->GetX();
}
 
Double_t BmnRootIterator::HitY() const
{
    if (!fInitialized || fEndOfHits)
        return 0.0;
    return fpCurrentHit->GetY();
}
 
Double_t BmnRootIterator::HitZ() const
{
    if (!fInitialized || fEndOfHits)
        return 0.0;
    return fpCurrentHit->GetZ();
}
 
Double_t BmnRootIterator::HitWx() const
{
    if (!fInitialized || fEndOfHits)
        return 0.0;
    return (fpCurrentHit->GetDx() < BMN_CLOSE_TO_ZERO) ? MAX_WEIGHT
                                                       : 1.0 / (fpCurrentHit->GetDx() * fpCurrentHit->GetDx());
}
 
Double_t BmnRootIterator::HitWy() const
{
    if (!fInitialized || fEndOfHits)
        return 0.0;
    return (fpCurrentHit->GetDy() < BMN_CLOSE_TO_ZERO) ? MAX_WEIGHT
                                                       : 1.0 / (fpCurrentHit->GetDy() * fpCurrentHit->GetDy());
}
 
Int_t BmnRootIterator::HitDetectorID() const
{
    if (!fInitialized || fEndOfHits)
        return 0;
    return fpDetModel->IDfromHit(fpCurrentHit->GetDetectorID());
    // this is the moment when hit detector encoding transformed into alignment encoding
}
 
Int_t BmnRootIterator::HitsInTrack() const
{
    if (!fInitialized || fEndOfTracks)
        return 0;
    return fHitsInTrack;
}
 
void BmnRootIterator::UnsetProperties()
{
    fTracksInEvent = 0;
    fHitsInTrack = 0;
    fpCurrentTrack = nullptr;
    fpCurrentHit = nullptr;
    fEndOfTracks = kTRUE;
    fEndOfHits = kTRUE;
}