EventLoaderCLICpix2.cpp 11.8 KB
Newer Older
1
#include "EventLoaderCLICpix2.h"
2

3
4
#include "CLICpix2/clicpix2_pixels.hpp"
#include "CLICpix2/clicpix2_utilities.hpp"
5

6
using namespace corryvreckan;
7
using namespace std;
8
9
using namespace caribou;
using namespace clicpix2_utils;
10

11
EventLoaderCLICpix2::EventLoaderCLICpix2(Configuration config, std::shared_ptr<Detector> detector)
12
    : Module(std::move(config), detector), m_detector(detector) {
13

14
    discardZeroToT = m_config.get<bool>("discard_zero_tot", false);
15
}
16

17
void EventLoaderCLICpix2::initialise() {
18

Simon Spannagel's avatar
Simon Spannagel committed
19
    // Take input directory from global parameters
20
    string inputDirectory = m_config.getPath("input_directory");
21

Simon Spannagel's avatar
Simon Spannagel committed
22
23
    // Open the root directory
    DIR* directory = opendir(inputDirectory.c_str());
24
    if(directory == nullptr) {
25
        throw ModuleError("Directory " + inputDirectory + " does not exist");
Simon Spannagel's avatar
Simon Spannagel committed
26
    }
27

Simon Spannagel's avatar
Simon Spannagel committed
28
29
30
    dirent* entry;

    // Read the entries in the folder
31
    while((entry = readdir(directory))) {
Simon Spannagel's avatar
Simon Spannagel committed
32
33
        // Check for the data file
        string filename = inputDirectory + "/" + entry->d_name;
Simon Spannagel's avatar
Simon Spannagel committed
34
        if(filename.find(".raw") != string::npos) {
CLICdp user's avatar
CLICdp user committed
35
            LOG(INFO) << "Found file " << filename;
Simon Spannagel's avatar
Simon Spannagel committed
36
            m_filename = filename;
37
38
39
40
41
42
43
44
45
46
47
            LOG(INFO) << "Found data file: " << m_filename;
        }
        if(filename.find(".cfg") != string::npos) {
            if(filename.find("matrix") != string::npos) {
                m_matrix = filename;
                LOG(INFO) << "Found matrix file: " << m_matrix;
            }
        }
    }

    if(m_matrix.empty()) {
48
        throw ModuleError("No matrix configuration file found in " + inputDirectory);
49
50
51
52
    }

    // Read the matrix configuration:
    matrix_config = clicpix2_utils::readMatrix(m_matrix);
Simon Spannagel's avatar
Simon Spannagel committed
53
    // If no data was loaded, give a warning
54
    if(m_filename.empty()) {
Simon Spannagel's avatar
Simon Spannagel committed
55
        LOG(WARNING) << "No data file was found for CLICpix2 in " << inputDirectory;
56
    }
57

Simon Spannagel's avatar
Simon Spannagel committed
58
59
    // Open the data file for later
    m_file.open(m_filename.c_str());
60

61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
    // Compression flags
    comp = true;
    sp_comp = true;

    std::streampos oldpos;
    std::string line;

    // Parse the header:
    while(getline(m_file, line)) {

        if(!line.length()) {
            continue;
        }
        if('#' != line.at(0)) {
            break;
        }

        // Replicate header to new file:
        LOG(DEBUG) << "Detected file header: " << line;
        oldpos = m_file.tellg();

        // Search for compression settings:
        std::string::size_type n = line.find(" sp_comp:");
        if(n != std::string::npos) {
            LOG(DEBUG) << "Value read for sp_comp: " << line.substr(n + 9, 1);
            sp_comp = static_cast<bool>(std::stoi(line.substr(n + 9, 1)));
            LOG(INFO) << "Superpixel Compression: " << (sp_comp ? "ON" : "OFF");
        }
        n = line.find(" comp:");
        if(n != std::string::npos) {
            LOG(DEBUG) << "Value read for comp: " << line.substr(n + 6, 1);
            comp = static_cast<bool>(std::stoi(line.substr(n + 6, 1)));
            LOG(INFO) << "     Pixel Compression: " << (comp ? "ON" : "OFF");
        }
    }

    decoder = new clicpix2_frameDecoder(comp, sp_comp, matrix_config);
    LOG(INFO) << "Prepared CLICpix2 frame decoder.";

100
101
102
103
104
105
    // Check if longcounting mode is on, defined by pixel 0,0:
    int maxcounter = 256;
    if(matrix_config[std::make_pair(0, 0)].GetLongCounter()) {
        maxcounter = 8192;
    }

Simon Spannagel's avatar
Simon Spannagel committed
106
    // Make histograms for debugging
107
    std::string title = m_detector->name() + " Hit map;x [px];y [px];pixels";
108
    hHitMap = new TH2F("hitMap",
109
                       title.c_str(),
110
                       m_detector->nPixels().X(),
111
                       0,
112
113
                       m_detector->nPixels().X(),
                       m_detector->nPixels().Y(),
114
                       0,
115
                       m_detector->nPixels().Y());
116
    title = m_detector->name() + " Map of discarded hits;x [px];y [px];pixels";
117
    hHitMapDiscarded = new TH2F("hitMapDiscarded",
118
                                title.c_str(),
119
                                m_detector->nPixels().X(),
120
                                0,
121
122
                                m_detector->nPixels().X(),
                                m_detector->nPixels().Y(),
123
                                0,
124
                                m_detector->nPixels().Y());
125
    title = m_detector->name() + " TOT spectrum;TOT;pixels";
Simon Spannagel's avatar
Simon Spannagel committed
126
    hPixelToT = new TH1F("pixelToT", title.c_str(), 32, 0, 31);
127
    title = m_detector->name() + " TOT map;x [px];y [px];TOT";
128
    hPixelToTMap = new TProfile2D("pixelToTMap",
129
                                  title.c_str(),
130
                                  m_detector->nPixels().X(),
131
                                  0,
132
133
                                  m_detector->nPixels().X(),
                                  m_detector->nPixels().Y(),
134
                                  0,
135
                                  m_detector->nPixels().Y(),
136
137
                                  0,
                                  maxcounter - 1);
138
    title = m_detector->name() + " TOA spectrum;TOA;pixels";
Simon Spannagel's avatar
Simon Spannagel committed
139
    hPixelToA = new TH1F("pixelToA", title.c_str(), maxcounter, 0, maxcounter - 1);
140
    title = m_detector->name() + " CNT spectrum;CNT;pixels";
Simon Spannagel's avatar
Simon Spannagel committed
141
    hPixelCnt = new TH1F("pixelCnt", title.c_str(), maxcounter, 0, maxcounter - 1);
142
    title = m_detector->name() + " Pixel Multiplicity; # pixels; # events";
143
    hPixelMultiplicity = new TH1F("pixelMultiplicity", title.c_str(), 1000, 0, 1000);
144

145
146
147
    title = m_detector->name() + " Timewalk;TOA;TOT;pixels";
    hTimeWalk = new TH2F("timewalk", title.c_str(), maxcounter, 0, maxcounter - 1, 32, 0, 31);

148
    title = m_detector->name() + " Map of masked pixels;x [px];y [px];mask code";
149
    hMaskMap = new TH2F("maskMap",
150
                        title.c_str(),
151
                        m_detector->nPixels().X(),
152
                        0,
153
154
                        m_detector->nPixels().X(),
                        m_detector->nPixels().Y(),
155
                        0,
156
157
158
                        m_detector->nPixels().Y());
    for(int column = 0; column < m_detector->nPixels().X(); column++) {
        for(int row = 0; row < m_detector->nPixels().Y(); row++) {
159
            if(m_detector->masked(column, row)) {
160
                hMaskMap->Fill(column, row, 2);
Simon Spannagel's avatar
Simon Spannagel committed
161
162
            }
            if(matrix_config[std::make_pair(row, column)].GetMask()) {
163
164
165
166
167
                hMaskMap->Fill(column, row, 1);
            }
        }
    }

Simon Spannagel's avatar
Simon Spannagel committed
168
169
    // Initialise member variables
    m_eventNumber = 0;
170
171
}

172
StatusCode EventLoaderCLICpix2::run(std::shared_ptr<Clipboard> clipboard) {
173

Simon Spannagel's avatar
Simon Spannagel committed
174
175
176
    // If have reached the end of file, close it and exit program running
    if(m_file.eof()) {
        m_file.close();
177
        return StatusCode::Failure;
178
    }
179

Simon Spannagel's avatar
Simon Spannagel committed
180
181
    // Pixel container, shutter information
    Pixels* pixels = new Pixels();
182
    long long int shutterStartTimeInt = 0, shutterStopTimeInt = 0;
183
184
    double shutterStartTime, shutterStopTime;
    string datastring;
Simon Spannagel's avatar
Simon Spannagel committed
185
186
    int npixels = 0;
    bool shutterOpen = false;
187
    std::vector<uint32_t> rawData;
Simon Spannagel's avatar
Simon Spannagel committed
188
189

    // Read file and load data
190
    while(getline(m_file, datastring)) {
Simon Spannagel's avatar
Simon Spannagel committed
191
192

        // Check if this is a header
193
194
        if(datastring.find("=====") != string::npos) {
            std::istringstream header(datastring);
195
            std::string guff;
Simon Spannagel's avatar
Simon Spannagel committed
196
197
            int frameNumber;
            header >> guff >> frameNumber >> guff;
198
            LOG(DEBUG) << "Found next header, frame number = " << frameNumber;
Simon Spannagel's avatar
Simon Spannagel committed
199
            break;
200
        }
Simon Spannagel's avatar
Simon Spannagel committed
201
        // If there is a colon, then this is a timestamp
202
203
        else if(datastring.find(":") != string::npos) {
            istringstream timestamp(datastring);
Simon Spannagel's avatar
Simon Spannagel committed
204
205
206
207
            char colon;
            int value;
            long long int time;
            timestamp >> value >> colon >> time;
208
            LOG(DEBUG) << "Found timestamp: " << time;
Simon Spannagel's avatar
Simon Spannagel committed
209
210
211
212
213
214
215
            if(value == 3) {
                shutterOpen = true;
                shutterStartTimeInt = time;
            } else if(value == 1 && shutterOpen) {
                shutterOpen = false;
                shutterStopTimeInt = time;
            }
CLICdp user's avatar
CLICdp user committed
216
        } else {
217
            // Otherwise pixel data
218
            rawData.push_back(static_cast<unsigned int>(std::atoi(datastring.c_str())));
219
        }
220
221
    }

222
    // Now set the event time so that the Timepix3 data is loaded correctly, unit is nanoseconds
223
    // NOTE FPGA clock is always on 100MHz from CaR oscillator, same as chip
224
225
    shutterStartTime = static_cast<double>(shutterStartTimeInt) / 0.1;
    shutterStopTime = static_cast<double>(shutterStopTimeInt) / 0.1;
226

227
    try {
CLICdp user's avatar
CLICdp user committed
228
        LOG(DEBUG) << "Decoding data frame...";
229
230
231
232
233
234
235
236
237
238
        decoder->decode(rawData);
        pearydata data = decoder->getZerosuppressedFrame();

        for(const auto& px : data) {

            auto cp2_pixel = dynamic_cast<caribou::pixelReadout*>(px.second.get());
            int col = px.first.first;
            int row = px.first.second;

            // If this pixel is masked, do not save it
239
            if(m_detector->masked(col, row)) {
240
241
242
                continue;
            }

243
            // Disentangle data types from pixel:
244
            int tot = -1, toa = -1, cnt = -1;
245
246
247
248

            // ToT will throw if longcounter is enabled:
            try {
                tot = cp2_pixel->GetTOT();
249
250
251
252
                if(!discardZeroToT || tot > 0) {
                    hPixelToT->Fill(tot);
                    hPixelToTMap->Fill(col, row, tot);
                }
253
            } catch(caribou::DataException&) {
254
                // Set ToT to one if not defined.
255
                tot = 1;
256
257
            }

258
259
260
            // Time defaults ot rising shutter edge:
            double timestamp = shutterStartTime;

261
262
263
            // Decide whether information is counter of ToA
            if(matrix_config[std::make_pair(row, col)].GetCountingMode()) {
                cnt = cp2_pixel->GetCounter();
264
265
266
                if(!discardZeroToT || tot > 0) {
                    hPixelCnt->Fill(cnt);
                }
267
268
            } else {
                toa = cp2_pixel->GetTOA();
269
270
271
                // Convert ToA form 100MHz clk into ns and sutract from shutterStopTime. Then add configured detector time
                // offset
                timestamp = shutterStopTime - static_cast<double>(toa) / 0.1 + m_detector->timingOffset();
272
273
                if(!discardZeroToT || tot > 0) {
                    hPixelToA->Fill(toa);
274
                    hTimeWalk->Fill(toa, tot);
275
                }
276
277
            }

278
279
280
281
            if(col >= m_detector->nPixels().X() || row >= m_detector->nPixels().Y()) {
                LOG(WARNING) << "Pixel address " << col << ", " << row << " is outside of pixel matrix.";
            }

282
283
            // when calibration is not available, set charge = tot
            Pixel* pixel = new Pixel(m_detector->name(), col, row, tot, tot, timestamp);
284
285
286
287
288
289
290
291

            if(tot == 0 && discardZeroToT) {
                hHitMapDiscarded->Fill(col, row);
            } else {
                pixels->push_back(pixel);
                npixels++;
                hHitMap->Fill(col, row);
            }
292
293
294
        }
    } catch(caribou::DataException& e) {
        LOG(ERROR) << "Caugth DataException: " << e.what() << ", clearing event data.";
295
    }
CLICdp user's avatar
CLICdp user committed
296
    LOG(DEBUG) << "Finished decoding, storing " << pixels->size() << " pixels";
297

CLICdp user's avatar
CLICdp user committed
298
    // Store current frame time and the length of the event:
299
300
    LOG(DEBUG) << "Event time: " << Units::display(shutterStartTime, {"ns", "us", "s"})
               << ", length: " << Units::display((shutterStopTime - shutterStartTime), {"ns", "us", "s"});
Simon Spannagel's avatar
Simon Spannagel committed
301
    clipboard->put_event(std::make_shared<Event>(shutterStartTime, shutterStopTime));
302

Simon Spannagel's avatar
Simon Spannagel committed
303
    // Put the data on the clipboard
CLICdp user's avatar
CLICdp user committed
304
    if(!pixels->empty()) {
305
        clipboard->put(m_detector->name(), "pixels", reinterpret_cast<Objects*>(pixels));
Simon Spannagel's avatar
Simon Spannagel committed
306
307
    } else {
        delete pixels;
308
        return StatusCode::NoData;
CLICdp user's avatar
CLICdp user committed
309
    }
310

Simon Spannagel's avatar
Simon Spannagel committed
311
    // Fill histograms
312
    hPixelMultiplicity->Fill(npixels);
313

Simon Spannagel's avatar
Simon Spannagel committed
314
    // Return value telling analysis to keep running
315
    return StatusCode::Success;
316
}
Simon Spannagel's avatar
Simon Spannagel committed
317
318
319
320

void EventLoaderCLICpix2::finalise() {
    delete decoder;
}