BmnCSCDigitizer.cxx

Go to the documentation of this file.

#include "BmnCSCDigitizer.h"
 
#include "BmnCSCPoint.h"
#include "CbmMCTrack.h"
#include "FairRootManager.h"
#include "TSystem.h"
 
static Float_t workTime = 0.0;
static int entrys = 0;
 
BmnCSCDigitizer::BmnCSCDigitizer()
    : fOnlyPrimary(kFALSE)
    , fStripMatching(kTRUE)
{
 
    fInputBranchName = "CSCPoint";
    fOutputDigitsBranchName = "BmnCSCDigit";
    fOutputDigitMatchesBranchName = "BmnCSCDigitMatch";
 
    fVerbose = 1;
 
    fCurrentConfig = BmnCSCConfiguration::None;
    StationSet = nullptr;
    GlobalModTransfSet = nullptr;
    fBmnCSCDigitsArray = nullptr;
    fBmnCSCDigitMatchesArray = nullptr;
}
 
BmnCSCDigitizer::BmnCSCDigitizer(UInt_t run_period)
    : BmnCSCDigitizer()
{
    switch (run_period) {
        case 7:   // BM@N RUN-7
            fCurrentConfig = BmnCSCConfiguration::RunSpring2018;
            break;
        case 8:   // BM@N RUN-8
            fCurrentConfig = BmnCSCConfiguration::FullCSC_Run8;
            break;
        case 9:   // BM@N RUN-9
            fCurrentConfig = BmnCSCConfiguration::Run9;
            break;
 
        default:
            fCurrentConfig = BmnCSCConfiguration::None;
    }
}
 
BmnCSCDigitizer::~BmnCSCDigitizer()
{
    if (StationSet) {
        delete StationSet;
    }
 
    if (GlobalModTransfSet) {
        delete GlobalModTransfSet;
    }
    if (fBmnCSCDigitsArray) {
        fBmnCSCDigitsArray->Delete();
        delete fBmnCSCDigitsArray;
    }
    if (fBmnCSCDigitMatchesArray) {
        fBmnCSCDigitMatchesArray->Delete();
        delete fBmnCSCDigitMatchesArray;
    }
}
 
InitStatus BmnCSCDigitizer::Init()
{
    if (fVerbose)
        cout << "\nBmnCSCDigitizer::Init()\n ";
 
    if (fVerbose && fOnlyPrimary)
        cout << "  Only primary particles are processed!!! " << endl;
 
    if (fVerbose && fStripMatching)
        cout << "  Strip Matching is activated!!! " << endl;
    else
        cout << "  Strip Matching is deactivated!!! " << endl;
 
    FairRootManager* ioman = FairRootManager::Instance();
 
    fBmnCSCPointsArray = (TClonesArray*)ioman->GetObject(fInputBranchName);
    if (!fBmnCSCPointsArray) {
        cout << "BmnCSCDigitizer::Init(): branch " << fInputBranchName << " not found! Task will be deactivated"
             << endl;
        SetActive(kFALSE);
        return kERROR;
    }
 
    fMCTracksArray = (TClonesArray*)ioman->GetObject("MCTrack");
 
    fBmnCSCDigitsArray = new TClonesArray(fOutputDigitsBranchName);
    ioman->Register(fOutputDigitsBranchName, "CSC_DIGIT", fBmnCSCDigitsArray, kTRUE);
 
    if (fStripMatching) {
        fBmnCSCDigitMatchesArray = new TClonesArray("BmnMatch");
        ioman->Register(fOutputDigitMatchesBranchName, "CSC_DIGIT", fBmnCSCDigitMatchesArray, kTRUE);
    }
 
    TString gPathCSCConfig = gSystem->Getenv("VMCWORKDIR");
    gPathCSCConfig += "/parameters/csc/XMLConfigs/";
 
    // Create CSC detector ------------------------------------------------------
    switch (fCurrentConfig) {
        case BmnCSCConfiguration::RunSpring2018:
            XMLConfigFile = gPathCSCConfig + "CSCRunSpring2018.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : RunSpring2018"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::RunSRCSpring2018:
            XMLConfigFile = gPathCSCConfig + "CSCRunSRCSpring2018.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : RunSRCSpring2018"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::Run8:
            XMLConfigFile = gPathCSCConfig + "CSCRun8.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : Run8"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::RunSRC2021:
            XMLConfigFile = gPathCSCConfig + "CSCRunSRC2021.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : RunSRC2021"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::LargeCSC_Run8:
            XMLConfigFile = gPathCSCConfig + "LargeCSCRun8.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : LargeCSCRun8"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::FullCSC_Run8:
            XMLConfigFile = gPathCSCConfig + "FullCSCRun8.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : FullCSCRun8 (4 small + 1 large chambers)"
                     << "\n";
            break;
 
        case BmnCSCConfiguration::Run9:
            XMLConfigFile = gPathCSCConfig + "CSCRun9.xml";
            if (fVerbose)
                cout << "   Current CSC Configuration : Run9"
                     << "\n";
 
        default:;
    }
    if (!gSystem->AccessPathName(XMLConfigFile)) {   // returns FALSE if the file exists
        if (fVerbose)
            cout << "   XMLConfigFile : " << XMLConfigFile << "\n";
        StationSet = new BmnCSCStationSet(XMLConfigFile);
        GlobalModTransfSet = new BmnCSCTransform();
        GlobalModTransfSet->LoadFromXMLFile(XMLConfigFile);
    } else {
        Fatal("BmnCSCDigitizer::Init()", " !!! Current configuration is not set !!! ");
    }
    //--------------------------------------------------------------------------
 
    if (fVerbose)
        cout << "BmnCSCDigitizer::Init() finished\n\n";
    return kSUCCESS;
}
 
void BmnCSCDigitizer::Exec(Option_t* opt)
{
 
    if (!IsActive())
        return;
 
    clock_t tStart = clock();
    fBmnCSCDigitsArray->Delete();
 
    if (fStripMatching) {
        fBmnCSCDigitMatchesArray->Delete();
    }
 
    // if (!fBmnCSCPointsArray) {
    //     Error("BmnCSCDigitizer::Exec()", " !!! Unknown branch name !!! ");
    //     return;
    // }
 
    if (fVerbose) {
        // cout << " BmnCSCDigitizer::Exec(), Number of BmnCSCPoints = " << fBmnCSCPointsArray->GetEntriesFast() <<
        // "\n";
    }
 
    ProcessMCPoints();
 
    // if (fVerbose) cout << " BmnCSCDigitizer::Exec() finished\n\n";
    entrys++;
    clock_t tFinish = clock();
    workTime += ((Float_t)(tFinish - tStart)) / CLOCKS_PER_SEC;
}
 
void BmnCSCDigitizer::ProcessMCPoints()
{
 
    FairMCPoint* CSCPoint;
    Int_t NNotPrimaries = 0;
 
    for (Int_t ipoint = 0; ipoint < fBmnCSCPointsArray->GetEntriesFast(); ipoint++) {
        CSCPoint = (FairMCPoint*)fBmnCSCPointsArray->At(ipoint);
 
        if (fOnlyPrimary) {
            CbmMCTrack* track = (CbmMCTrack*)fMCTracksArray->UncheckedAt(CSCPoint->GetTrackID());
            if (!track)
                continue;
            if (track->GetMotherId() != -1) {
                NNotPrimaries++;
                continue;
            }
        }
 
        Double_t x = -CSCPoint->GetX();   // invert because in the current geometry +x is left, -x is right
        Double_t y = CSCPoint->GetY();
        Double_t z = CSCPoint->GetZ();
 
        Double_t px = -CSCPoint->GetPx();   // invert because in the current geometry +x is left, -x is right
        Double_t py = CSCPoint->GetPy();
        Double_t pz = CSCPoint->GetPz();
 
        Double_t dEloss = CSCPoint->GetEnergyLoss() * 1e6;   // in keV
        Int_t refId = ipoint;
 
        // Information from MC-points
        Int_t mc_station_num = ((BmnCSCPoint*)CSCPoint)->GetStation();
        Int_t mc_module_num = ((BmnCSCPoint*)CSCPoint)->GetModule();
 
        // test output
        // cout << "mc_station_num = " << mc_station_num << "\n";
        // cout << "mc_module_num = " << mc_module_num << "\n";
 
        // Apply inverse global module transforms to MC-point (in right-handed system)
        if (GlobalModTransfSet && mc_station_num < StationSet->GetNStations()) {
            if (mc_module_num < StationSet->GetStation(mc_station_num)->GetNModules()) {
                Plane3D::Point point_transformed =
                    GlobalModTransfSet->ApplyInverseTransforms(Plane3D::Point(-x, y, z), mc_station_num, mc_module_num);
                Plane3D::Point direction_transformed = GlobalModTransfSet->ApplyInverseTransforms(
                    Plane3D::Point(-(px + x), (py + y), (pz + z)), mc_station_num, mc_module_num);
                x = -point_transformed.X();   // convert to detector system (x-inverted)
                y = point_transformed.Y();
                z = point_transformed.Z();
 
                px = -(direction_transformed.X() - point_transformed.X());   // convert to detector system (x-inverted)
                py = direction_transformed.Y() - point_transformed.Y();
                pz = direction_transformed.Z() - point_transformed.Z();
            }
        }
 
        StationSet->AddPointToDetector(x, y, z, px, py, pz, dEloss, refId);
    }
 
    // Int_t NAddedPoints = StationSet->CountNAddedToDetectorPoints();
    // if (fVerbose && fOnlyPrimary) cout << "   Number of not primaries points : " << NNotPrimaries << "\n";
    // if (fVerbose) cout << "   Processed MC points  : " << NAddedPoints << "\n";
 
    for (Int_t iStation = 0; iStation < StationSet->GetNStations(); ++iStation) {
        BmnCSCStation* station = StationSet->GetStation(iStation);
 
        for (Int_t iModule = 0; iModule < station->GetNModules(); ++iModule) {
            BmnCSCModule* module = station->GetModule(iModule);
 
            for (Int_t iLayer = 0; iLayer < module->GetNStripLayers(); ++iLayer) {
                BmnCSCLayer layer = module->GetStripLayer(iLayer);
 
                for (Int_t iStrip = 0; iStrip < layer.GetNStrips(); ++iStrip) {
                    Double_t signal = layer.GetStripSignal(iStrip);
                    if (signal > 0.0) {
                        new ((*fBmnCSCDigitsArray)[fBmnCSCDigitsArray->GetEntriesFast()])
                            BmnCSCDigit(iStation, iModule, iLayer, iStrip, signal);
 
                        if (fStripMatching) {
                            new ((*fBmnCSCDigitMatchesArray)[fBmnCSCDigitMatchesArray->GetEntriesFast()])
                                BmnMatch(layer.GetStripMatch(iStrip));
                        }
                    }
                }
            }
        }
    }
    StationSet->Reset();
}
 
void BmnCSCDigitizer::Finish()
{
    if (StationSet) {
        delete StationSet;
        StationSet = nullptr;
    }
 
    if (GlobalModTransfSet) {
        delete GlobalModTransfSet;
        GlobalModTransfSet = nullptr;
    }
 
    cout << "Work time of the CSC digitizer: " << workTime << endl;
}