TreeWriterDUT.cpp 5.79 KB
Newer Older
1
#include "TreeWriterDUT.h"
2
3
#include "core/utils/file.h"

4
5
6
7
8
#include <vector>

using namespace corryvreckan;
using namespace std;

9
10
TreeWriterDUT::TreeWriterDUT(Configuration config, std::shared_ptr<Detector> detector)
    : Module(std::move(config), detector), m_detector(detector) {
11

12
13
    m_fileName = m_config.get<std::string>("file_name", "outputTuples.root");
    m_treeName = m_config.get<std::string>("tree_name", "tree");
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
}

/*

 This algorithm writes an output file and fills it with trees containing
 the following information for the DUT:

 1)size in x of each cluster,
 2)size in y of each cluster,
 3)the event ID of each cluster,
 4)a list of the X positions of the pixels of each cluster,
 5)a list of the Y positions of the pixels of each cluster,
 6)a list of the ToT of the pixels of each cluster,
 7)a list of the timestamp/ToA of the pixels of each cluster,
 the number of pixels in each cluster (so that in the analysis the lists
 of pixel X positions etc. can be split into their respective clusters).
 8)the intercept with the DUT for each track

 */

34
35
void TreeWriterDUT::initialise() {
    LOG(DEBUG) << "Initialised TreeWriterDUT";
36
37

    // Create output file and directories
38
39
40
    auto path = createOutputFile(add_file_extension(m_fileName, "root"));
    m_outputFile = new TFile(path.c_str(), "RECREATE");
    LOG(DEBUG) << "Made and moved to output file: " << path;
41
42
43
44
45
    gDirectory->Delete("tree;*");
    m_outputTree = new TTree(m_treeName.c_str(), m_treeName.c_str());
    LOG(DEBUG) << "Created tree: " << m_treeName;

    eventID = 0;
46
47
48
49
50
51
52
53
54
55
56
57
58
    filledEvents = 0;

    // Create the output branches
    m_outputTree->Branch("EventID", &v_clusterEventID);
    m_outputTree->Branch("clusterSizeX", &v_clusterSizeX);
    m_outputTree->Branch("clusterSizeY", &v_clusterSizeY);
    m_outputTree->Branch("pixelX", &v_pixelX);
    m_outputTree->Branch("pixelY", &v_pixelY);
    m_outputTree->Branch("pixelToT", &v_pixelToT);
    m_outputTree->Branch("pixelToA", &v_pixelToA);
    m_outputTree->Branch("clusterNumPixels", &v_clusterNumPixels);
    // Branch for track intercepts
    m_outputTree->Branch("intercepts", &v_intercepts);
59
60
}

61
StatusCode TreeWriterDUT::run(std::shared_ptr<Clipboard> clipboard) {
62
63
    // Counter for cluster event ID
    eventID++;
64

65
    // Clear data vectors before storing the cluster information for this event
66
67
68
69
70
71
72
73
74
75
76
    v_intercepts.clear();
    v_clusterSizeX.clear();
    v_clusterSizeY.clear();
    v_clusterEventID.clear();
    v_pixelX.clear();
    v_pixelY.clear();
    v_pixelToT.clear();
    v_pixelToA.clear();
    v_clusterNumPixels.clear();

    // Getting tracks from the clipboard
77
78
    Tracks* tracks = reinterpret_cast<Tracks*>(clipboard->get("tracks"));
    if(tracks == nullptr) {
79
        LOG(DEBUG) << "No tracks on the clipboard";
80
        return StatusCode::Success;
81
82
83
84
    }

    // Iterate through tracks found
    for(auto& track : (*tracks)) {
85
86
87
88
89
90
91
92
        // CHeck if we have associated clusters:
        Clusters associatedClusters = track->associatedClusters();
        if(associatedClusters.empty()) {
            LOG(TRACE) << "No associated clusters, skipping track.";
            continue;
        }
        // Drop events with more than one cluster associated
        else if(associatedClusters.size() > 1) {
93
            LOG(DEBUG) << "More than one associated cluster, dropping track.";
94
95
96
97
98
            continue;
        }

        LOG(DEBUG) << "Found track with associated cluster";

99
        // Get track intercept with DUT in global coordinates
100
        trackIntercept = m_detector->getIntercept(track);
101
102

        // Calculate the intercept in local coordinates
103
        trackInterceptLocal = m_detector->globalToLocal(trackIntercept);
104
        v_intercepts.push_back(trackInterceptLocal);
105
106
107
108

        Cluster* cluster = associatedClusters.front();

        // x size
109
110
        LOG(DEBUG) << "Gets column width = " << cluster->columnWidth();
        v_clusterSizeX.push_back(cluster->columnWidth());
111
112

        // y size
113
114
        LOG(DEBUG) << "Gets row width = " << cluster->rowWidth();
        v_clusterSizeY.push_back(cluster->rowWidth());
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129

        // eventID
        v_clusterEventID.push_back(eventID);
        LOG(DEBUG) << "Gets cluster eventID = " << eventID;

        // Get the pixels in the current cluster
        Pixels* pixels = cluster->pixels();

        // Iterate through all pixels in the cluster
        numPixels = 0;
        for(auto& pixel : (*pixels)) {
            // Increase counter for number of pixels in the cluster
            numPixels++;

            // x position
130
131
            LOG(DEBUG) << "Gets pixel column = " << pixel->column();
            v_pixelX.push_back(pixel->column());
132
133

            // y position
134
135
            LOG(DEBUG) << "Gets pixel row = " << pixel->row();
            v_pixelY.push_back(pixel->row());
136
137

            // ToT
138
139
            LOG(DEBUG) << "Gets pixel tot = " << pixel->adc();
            v_pixelToT.push_back(pixel->adc());
140
141

            // ToA
142
143
            LOG(DEBUG) << "Gets pixel timestamp = " << pixel->timestamp();
            v_pixelToA.push_back(pixel->timestamp());
144
        }
145
146
147
148
        v_clusterNumPixels.push_back(numPixels);
    }

    if(v_intercepts.empty()) {
149
        return StatusCode::NoData;
150
151
152
153
    }

    filledEvents++;
    if(filledEvents % 100 == 0) {
154
        LOG(DEBUG) << "Events with single associated cluster: " << filledEvents;
155
    }
156
157
    // Fill the tree with the information for this event
    m_outputTree->Fill();
158
159

    // Return value telling analysis to keep running
160
    return StatusCode::Success;
161
162
}

163
void TreeWriterDUT::finalise() {
164
    LOG(DEBUG) << "Finalise";
165
166
    auto directory = m_outputFile->mkdir("Directory");
    directory->cd();
167
    LOG(STATUS) << filledEvents << " events written to file " << m_fileName;
168

169
    auto orientation = m_detector->rotation();
170
    directory->WriteObject(&orientation, "DUTorientation");
171
172

    // Writing out outputfile
173
174
175
    m_outputFile->Write();
    delete(m_outputFile);
}