TRestEventTimeSelectionProcess.cxx

/*************************************************************************
 * This file is part of the REST software framework.                     *
 *                                                                       *
 * Copyright (C) 2016 GIFNA/TREX (University of Zaragoza)                *
 * For more information see http://gifna.unizar.es/trex                  *
 *                                                                       *
 * REST is free software: you can redistribute it and/or modify          *
 * it under the terms of the GNU General Public License as published by  *
 * the Free Software Foundation, either version 3 of the License, or     *
 * (at your option) any later version.                                   *
 *                                                                       *
 * REST is distributed in the hope that it will be useful,               *
 * but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the          *
 * GNU General Public License for more details.                          *
 *                                                                       *
 * You should have a copy of the GNU General Public License along with   *
 * REST in $REST_PATH/LICENSE.                                           *
 * If not, see http://www.gnu.org/licenses/.                             *
 * For the list of contributors see $REST_PATH/CREDITS.                  *
 *************************************************************************/
 
 
#include "TRestEventTimeSelectionProcess.h"
 
using namespace std;
 
ClassImp(TRestEventTimeSelectionProcess);
 
TRestEventTimeSelectionProcess::TRestEventTimeSelectionProcess() { Initialize(); }
 
void TRestEventTimeSelectionProcess::Initialize() {
    SetSectionName(this->ClassName());
    fEvent = nullptr;
    fFileWithTimes = "";
    fIsActiveTime = true;
    fDelimiter = ',';
    fStartEndTimes.clear();
    fTimeOffsetInSeconds = 0;
    fTimeStartMarginInSeconds = 0;
    fTimeEndMarginInSeconds = 0;
    fNEventsRejected = 0;
    fNEventsSelected = 0;
    fTotalTimeInSeconds = 0;
}
 
void TRestEventTimeSelectionProcess::InitProcess() {
    // Read the file with the time ranges
    if (!fFileWithTimes.empty()) {
        fStartEndTimes = ReadFileWithTimes(fFileWithTimes, fDelimiter);
    }
    fTotalTimeInSeconds = CalculateTotalTimeInSeconds();
    fNEventsRejected = 0;
    fNEventsSelected = 0;
}
 
std::vector<std::pair<std::string, std::string>> TRestEventTimeSelectionProcess::ReadFileWithTimes(
    std::string fileWithTimes, Char_t delimiter) {
    std::vector<std::pair<std::string, std::string>> startEndTimes;
    string line;
    ifstream file(fileWithTimes);
    if (file.is_open()) {
        while (getline(file, line)) {
            if (line[0] == '#') {  // understand as comment
                continue;
            }
            std::istringstream lineStream(line);
            std::string startDate, endDate;
            if (std::getline(lineStream, startDate, delimiter) &&
                std::getline(lineStream, endDate, delimiter)) {
                // check if the time format is correct. TODO: use better way to check
                // (StringToTimeStamp usually returns a negative big number if not)
                if (StringToTimeStamp(startDate) < 0 || StringToTimeStamp(endDate) < 0) {
                    continue;
                }
 
                startEndTimes.emplace_back(startDate, endDate);
            }
        }
        file.close();
    }
    return startEndTimes;
}
 
Double_t TRestEventTimeSelectionProcess::CalculateTotalTimeInSeconds() {
    Double_t totalTime = 0;
    for (auto id : fStartEndTimes) {
        TTimeStamp startTime = TTimeStamp(StringToTimeStamp(id.first), 0);
        TTimeStamp endTime = TTimeStamp(StringToTimeStamp(id.second), 0);
        // Reduce the time by the margin in both sides
        startTime.Add(TTimeStamp(fTimeStartMarginInSeconds));
        endTime.Add(TTimeStamp(-fTimeEndMarginInSeconds));
        auto timeDiff = endTime.AsDouble() - startTime.AsDouble();
        if (timeDiff < 0) {
            RESTDebug << "End time is before start time in time range: " << id.first << " to " << id.second
                      << RESTendl;
            continue;
        }
        totalTime += endTime.AsDouble() - startTime.AsDouble();
    }
    return totalTime;
}
 
TRestEvent* TRestEventTimeSelectionProcess::ProcessEvent(TRestEvent* inputEvent) {
    fEvent = inputEvent;
 
    TTimeStamp eventTime = fEvent->GetTimeStamp();
    eventTime.Add(TTimeStamp(fTimeOffsetInSeconds));
 
    Bool_t isInsideAnyTimeRange = false;
    for (auto id : fStartEndTimes) {
        TTimeStamp startTime = TTimeStamp(StringToTimeStamp(id.first), 0);
        TTimeStamp endTime = TTimeStamp(StringToTimeStamp(id.second), 0);
        // Reduce the time by the margin in both sides
        startTime.Add(TTimeStamp(fTimeStartMarginInSeconds));
        endTime.Add(TTimeStamp(-fTimeEndMarginInSeconds));
        if (eventTime >= startTime && eventTime <= endTime) {
            isInsideAnyTimeRange = true;
            break;
        }
    }
 
    // Decide if the event is selected or rejected based on the time ranges
    // and their meaning (active or dead periods of time).
    if (fIsActiveTime) {             // time ranges represent active periods of time
        if (isInsideAnyTimeRange) {  // time is inside an active period of time
            fNEventsSelected++;
            return fEvent;
        }
    } else {                          // time ranges represent dead periods of time
        if (!isInsideAnyTimeRange) {  // time is outside all dead period of time
            fNEventsSelected++;
            return fEvent;
        }
    }
 
    // rejected events are not returned
    fNEventsRejected++;
    return nullptr;
}
 
void TRestEventTimeSelectionProcess::EndProcess() {
    // Write here the jobs to do when all the events are processed
}
 
std::string TRestEventTimeSelectionProcess::GetTimeStampCut(std::string timeStampObsName, Bool_t useOffset,
                                                            Bool_t useMargins, Int_t nTimes) {
    std::string timeCut = "";
    std::string timeStampObsNameWithOffset = timeStampObsName;
    if (useOffset && fTimeOffsetInSeconds != 0) {
        timeStampObsNameWithOffset += "+" + to_string(fTimeOffsetInSeconds);
    }
    if (nTimes < 0) nTimes = fStartEndTimes.size();
    Int_t c = 0;
    for (auto id : fStartEndTimes) {
        if (c++ >= nTimes) break;
        auto startTime = StringToTimeStamp(id.first);
        auto endTime = StringToTimeStamp(id.second);
        // Reduce the time by the margin in both sides
        if (useMargins) {
            startTime += fTimeStartMarginInSeconds;
            endTime -= fTimeEndMarginInSeconds;
        }
 
        if (startTime >= endTime) {
            continue;
        }
 
        // Build the cut string
        if (!timeCut.empty()) {
            if (fIsActiveTime)
                timeCut += " || ";  // inside ANY time range
            else
                timeCut += " && ";  // outside ALL time ranges
        }
        if (!fIsActiveTime) timeCut += "!";  // NOT inside the time range
        // inside the time range
        timeCut += "(";
        timeCut += "(" + timeStampObsNameWithOffset + ">=" + to_string(startTime) + ")";
        timeCut += " && ";
        timeCut += "(" + timeStampObsNameWithOffset + "<=" + to_string(endTime) + ")";
        timeCut += ")";
    }
    return timeCut;
}
void TRestEventTimeSelectionProcess::PrintMetadata() {
    BeginPrintProcess();
    std::string typeOfTime = fIsActiveTime ? "Active" : "Dead";
 
    RESTMetadata << "File with times: " << fFileWithTimes << RESTendl;
    // print periods
    RESTMetadata << "Offset time: " << fTimeOffsetInSeconds << " seconds" << RESTendl;
    RESTMetadata << "Start margin time: " << fTimeStartMarginInSeconds << " seconds" << RESTendl;
    RESTMetadata << "End margin time: " << fTimeEndMarginInSeconds << " seconds" << RESTendl;
    RESTMetadata << typeOfTime << " time periods: " << RESTendl;
    for (auto id : fStartEndTimes) {
        RESTMetadata << id.first << " to " << id.second << RESTendl;
        TTimeStamp startTime = TTimeStamp(StringToTimeStamp(id.first), 0);
        TTimeStamp endTime = TTimeStamp(StringToTimeStamp(id.second), 0);
    }
 
    // Get total time in seconds
    TTimeStamp totalTime = TTimeStamp(fTotalTimeInSeconds, 0);
    if (!fStartEndTimes.empty()) {
        TTimeStamp firstTime = TTimeStamp(StringToTimeStamp(fStartEndTimes.front().first), 0);
        TTimeStamp lastTime = TTimeStamp(StringToTimeStamp(fStartEndTimes.back().second), 0);
        totalTime = lastTime - firstTime;
    }
 
    double fractionOfTime = fTotalTimeInSeconds / totalTime.AsDouble() * 100;
    std::string fractionOfTimeStr = StringWithPrecision(fractionOfTime, 4) + " %";
    if ((Int_t)(fTotalTimeInSeconds / 24 / 3600) != 0)  // order of days
        RESTMetadata << "Total " << typeOfTime << " time: " << fTotalTimeInSeconds / 24 / 3600 << " days"
                     << " (" << fractionOfTimeStr << ")" << RESTendl;
    else if ((Int_t)(fTotalTimeInSeconds / 3600) != 0)  // order of hours
        RESTMetadata << "Total " << typeOfTime << " time: " << fTotalTimeInSeconds / 3600 << " hours"
                     << " (" << fractionOfTimeStr << ")" << RESTendl;
    else if ((Int_t)(fTotalTimeInSeconds / 60) != 0)  // order of minutes
        RESTMetadata << "Total " << typeOfTime << " time: " << fTotalTimeInSeconds / 60 << " minutes"
                     << " (" << fractionOfTimeStr << ")" << RESTendl;
    else
        RESTMetadata << "Total " << typeOfTime << " time: " << fTotalTimeInSeconds << " seconds"
                     << " (" << fractionOfTimeStr << ")" << RESTendl;
 
    RESTMetadata << "Number of events rejected: " << fNEventsRejected << " ("
                 << fNEventsRejected * 1. / (fNEventsRejected + fNEventsSelected) * 100 << " %)" << RESTendl;
    RESTMetadata << "Number of events selected: " << fNEventsSelected << " ("
                 << fNEventsSelected * 1. / (fNEventsRejected + fNEventsSelected) * 100 << " %)" << RESTendl;
 
    EndPrintProcess();
}